Patent Application
20230260109
The present invention relates to the field of medical imaging, and more particularly to tumour imaging. Tumour imaging encompasses methods for non-invasive analysis of tumours. Analyses may be aimed at determining the extent of cancer cells within the tumour. Tumour imaging is carried out by magnetic resonance imaging (MRI), and in particular by diffusion-weighted magnetic resonance imaging. However, these two imaging techniques do not currently provide precise information about a tumour in a patient's body, such as the absolute distribution of cancer cells or the tumour load. By tumour load, it is meant here the total number of tumour cells of all types or each type individually, in particular the cancer cell type.
Biopsies are known to be used as an invasive method to supplement or replace tumour imaging. A biopsy consists in locally sampling tissue by a biopsy needle. Biopsy provides local information about the tumour that is not accessible by non-invasive techniques. However, biopsy is a local sampling, and the tissue sampled is usually not representative of the tumour being studied. Furthermore, biopsy is an invasive surgical procedure and poses medical risks such as tumour enlargement, spread of cancer cells into healthy tissue, bleeding, infection or nerve damage.
There is no method to determine a representative location to perform a biopsy. As a result, it is currently necessary to carry out many biopsies, sometimes as many as ten, to increase the chances of obtaining enough information to determine information about the absolute distribution sought or the tumour load.
The invention improves the situation. To this end, the invention provides an imaging device comprising:
This imaging device is very advantageous, as it makes it possible to automatically determine a reduced number of biopsy locations ensuring relevant information. Indeed, the punctures determined by the imaging device are much more representative of the entire tumour than randomly taken punctures, while at the same time being fewer in number.
In various alternatives, the device may have one or more of the following characteristics:
Dj=Dini+Delta(j−1),
where Dini satisfies the equation:
Dini=M−min(|M−Dmin|,|M−Dmax|)
and Delta satisfies the equation:
Delta=2(M−Dini)/(N−1),
Dmin being equal to the value of the first diffusion parameter determined by the selector and Dmax being equal to the value of the second diffusion parameter of the subset of diffusion parameters determined by the selector,
Further characteristics and advantages of the invention will be set out in detail in the following description, made with reference to the attached drawings, in which:
The appended drawings contain, for the most part, elements of certainty. They may therefore not only serve to better understand the present invention, but also contribute to its definition, where appropriate.
The description is followed by an Appendix A comprising mathematical formulae. This appendix forms an integral part of the description.
Reference is made to
An imaging device 1 according to the invention comprises a memory 10, an upsampler 20, a slicer 30, a selector 40 and a guide 50.
The memory 10 may be any type of data storage suitable for receiving digital data: hard disk, flash memory hard disk (SSD), flash memory in any form, random access memory, magnetic disk, locally or cloud distributed storage, etc. The data calculated by the device 1 can be stored on any type of memory similar to or on the memory 10. This data can be erased after the device has performed its tasks or retained.
The upsampler 20, the slicer 30, the selector 40 and the guide 50 are here programs run by the computer processor. Alternatively, one or more of these elements could be implemented in a different way by means of a dedicated processor. By processor, it has to be understood any processor adapted to data processing as described below. Such a processor may be realised in any known way, in the form of a microprocessor for a personal computer, a dedicated chip of the FPGA or SoC type, a computing resource on a grid, a microcontroller, or any other form suitable to provide the computing power necessary for the embodiment described below. One or more of these elements may also be realised in the form of specialised electronic circuits such as an ASIC. A combination of processors and electronic circuits may also be contemplated.
The memory 10 receives imaging data 12, comprising at least one raster image 14 from a cross-sectional diffusion-weighted magnetic resonance imaging (diffusion MRI or DWI) image of a tumour 200. The raster image 14 is also known as a “D-map”. Diffusion-weighted magnetic resonance is a tissue imaging technique consisting in measuring the diffusion of water molecules within the tissue. Each pixel of the raster image 14 is thus associated with a diffusion parameter whose value represents the mobility of water molecules within the tissue of the tumour 200. The higher the value of the diffusion parameter, the greater the diffusion of water molecules, or, in other words, the freer the water molecules are to move within the tissue.
Diffusion parameter values are highly dependent on the diffusion-weighted magnetic resonance imaging machine that generates them, which generally makes it impossible to analyse them without additional information. In the example described here, the raster image 14 represented in
The applicant's work has shown that the diffusion value is correlated with the cell density. Thus, starting from a raster image 14 derived from a cross-sectional diffusion-weighted magnetic resonance image, it would in theory be possible to observe cellular heterogeneities within the tumour, at least qualitatively. However, the resolution of the raster images provided by current diffusion-weighted magnetic resonance imaging devices is too low to represent such heterogeneities in a usable way, so that they cannot be used in practice, except to carry out a full tumourectomy to see the correspondence. A systematic tumourectomy is of course out of the question (and in most cases impossible), especially since if this were possible, imaging and biopsies would be useless.
In one embodiment, the imaging data 12 comprises a plurality of raster images 14 representing successive cross-sectional images of the tumour 200. This provides a three-dimensional view of the tumour 200. In the embodiment described here, the imaging data 12 comprises 9 raster images 14. In another embodiment, the imaging data 12 comprises 40 raster images 14 from which 9 raster images 14 representing sections closest to a median plane of the tumour 200 are derived. One of these raster images 14 is represented in
The memory 10 also receives biopsy data 16, comprising needle data 160 and procedure data 162.
The needle data 160 defines the dimensions of a puncture carried out with a biopsy needle. In the embodiment described here, the dimensions of the puncture comprise a puncture diameter and a puncture length, and the sampling with the biopsy needle is a substantially cylindrical puncture zone defined by that diameter and length. In the example described here, the puncture diameter is substantially 0.84 mm, and the puncture length is 4.55±1.45 mm.
The procedure data 162 comprises a number of punctures to be performed in a biopsy procedure. The number of punctures to be performed is an integer greater than or equal to 2. Indeed, it is necessary to perform at least two punctures to obtain sufficient information about the tumour 200. In the embodiment described here, the number of punctures is comprised between 2 and 4, inclusive.
The procedure data 162 also comprises a set of constraints 300. The set of constraints defines interventional constraints to be met in a biopsy procedure to carry out one or more punctures. The constraints in the sets of constraints may arise from the type of organ in which the tumour 200 is located, the size of the tumour 200, the presence of necrotic tissue 500 within the tumour 200, or any other parameter constraining the performance of the biopsies.
In the example described here, the tumour 200 is a lung tumour, and the biopsy constraints associated thereto, represented in
In the example described here, the maximum puncture angle 310 is 20° and the maximum puncture depth 312 is 2.2 cm.
Constraints may also include
As discussed above, the pixels each correspond to a region of the cross-sectional image that is substantially square in shape and 2.1 mm on a side, and the puncture zone for the puncture needle has a diameter of substantially 0.84 mm and a length of 4.55±1.45 mm. Thus, there are two obstacles to an analysis of the raster image 14:
The upsampler 20 refines the raster image to overcome the first obstacle.
To do so, the upsampler 20 receives the raster image 14 and upsamples it into a processing image 22. The processing image 22 is an image of much higher spatial resolution than the raster image 14, as can be seen by comparing
In the embodiment described here, the upsampler 20 applies an interpolation technique to upsample the raster image 14. The interpolation techniques used may be bicubic interpolation, bilinear interpolation, Bell filter, Hermite filter, Mitchell filter, Lanczos filter, or any other conventional interpolation technique adapted to images. In the example described here, the upsampler 20 applies the bicubic interpolation technique, described for example in Keys, R. (1981). “Cubic convolution interpolation for digital image processing.” IEEE transactions on Signal Processing. In Acoustics, Speech, and Signal Processing (Vol. 29, p. 1153).
Alternatively, the upsampler 20 may apply a Discrete Cosine Transform (DCT) technique to upsample the raster image 14, such as those described in Dugad, R., & Ahuja, N. (2001). “A fast scheme for image size change in the compressed domain.” IEEE Transactions on Circuits and Systems for Video Technology, 11(4), 461-474. and in Park, H., Park, Y., & Oh, S. K. (2003). “L/M-fold image resizing in block-DCT domain using symmetric convolution.” IEEE Transactions on Image Processing, 12(9), 1016-1034. Alternatively, a machine learning-based technique, in particular deep learning, such as that described in Wang, Z., Chen, J., & Hoi, S. C. (2019). “Deep learning for image super-resolution: A survey.” arXiv preprint arXiv:1902.06068, could be used.
The processing image 22 has pixels whose side corresponds to a substantially square zone with a side at least 5 times smaller than the side of the zone of a pixel of the raster image 14, preferably 20 times smaller, as in the example described here.
In another embodiment, the imaging device 1 works in three dimensions, and the imaging data 12 comprises 9 successive cross-sectional raster images 14, separated in twos by a distance of 6 mm. The raster images 14 thus define voxels of dimension 2.1 mm×2.1 mm×6 mm. The processing images 22 define voxels at least 5 times smaller in each dimension, preferably 20 times smaller. In the example described here, the voxels of the processing images 22 each define a zone of 0.105×0.105×0.3 mm.
The processing image 22 has a much higher resolution than the raster image 14. The regions corresponding to the pixels in the processing image are much smaller than the heterogeneities in the tumour. The pixels of the processing image 22 now have to be grouped together to slice the cross-sectional image into a plurality of regions sufficiently homogeneous to overcome the first obstacle and sufficiently large to overcome the second obstacle.
To do so, the slicer 30 receives and partitions the processing image 22 into a set of regions 32. Each of the regions 34 of the set of regions 32 comprises pixels of the processing image 22 that are directly adjoining, at least in twos.
The region 34 is associated with a diffusion parameter value substantially equal to the mean of the diffusion parameter values of its constituent pixels. This diffusion parameter value of the region 34 may be chosen, for example, as the mean value of the diffusion parameter values associated with the pixels of the region 34, or as the median value.
The slicer 30 determines the set of regions 32 such that, for each region 34, the variance of the diffusion parameter values associated with the pixels of that region 34 is less than a threshold value. The regions 34 of the set of regions 32 are thus sets of contiguous pixels of the processing image 22 whose associated diffusion parameter values are substantially homogeneous, or, in other words, whose variance is low.
The slicer 30 further determines the set of regions 32 such that the region 34 corresponds to a portion of the cross-sectional diffusion-weighted magnetic resonance image whose dimensions are greater than or substantially equal to the dimensions of a puncture carried out by means of a biopsy needle of the needle data 160. Here, the portions to which the regions 34 correspond have substantially similar dimensions, with for example substantially equal respective areas. The diffusion parameters associated with each region 34 of the set of regions 32 form a set of diffusion parameters 36.
In one embodiment, the slicer 30 partitions the processing image 22 into the set of regions 32 by a super-pixel method with constraints on the size and low variance of the diffusion parameter values associated with the pixels of each region 34. In the example described here, the super-pixel method used is that described in Ren, X., & Malik, J. (2003, October). “Learning a classification model for segmentation.” In Proceedings Ninth IEEE International Conference on Computer Vision. (p. 10). IEEE. In the example described here, represented in
In
Alternatively, the slicer 30 may partition the processing image 22 into the set of regions 32 by a convex polygon partitioning method as described in Duan, L., & Lafarge, F. (2015). “Image partitioning into convex polygons.” In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (pp. 3119-3127). This method is based on the construction of a Voronoi diagram whose cells are convex polygons by conforming to straight lines already drawn within the processing image 22. Here, regions 34 can be obtained whose dimensions can be controlled, and which preserve the edges of shapes or objects within the processing image 22 at a subpixel scale. The set of regions 32 here provides a representation of the heterogeneities within the processing image 22.
In the three-dimensional case, where several processing images 22 are produced by the upsampler 20, the slicer 30 partitions the processing images 22 using a super-voxel method, such as one described in the paper Xu, C., & Corso, J. J. (2012, June). “Evaluation of super-voxel methods for early video processing.” In 2012 IEEE conference on computer vision and pattern recognition (pp. 1202-1209). IEEE.
The regions 34 determined by the slicer 30 are thus
The selector 40 determines, from the regions 34, which are relevant to a sample, so that the diffusion parameter values in those regions are representative of all diffusion parameters.
To do so, the selector 40 determines within the set of diffusion parameters 36 a subset of diffusion parameters 42 which comprises as many diffusion parameters as the number of punctures to be performed of the procedure data 162. Since the number of punctures to be performed is at least two, the subset of diffusion parameters 42 comprises at least two diffusion parameters. A first diffusion parameter 420 of the diffusion parameters is selected from the first decile of the diffusion parameter values of the set of diffusion parameters 36. A second diffusion parameter 426 of the diffusion parameters is selected from the last decile of the diffusion parameter values of the set of diffusion parameters 36. The first diffusion parameter 420 and the second diffusion parameter 426 may for example each be the median value of their respective decile. This ensures that the subset of diffusion parameters 42 determined by the selector 40 is representative of the diffusion parameter values of the regions 34 of the set of regions 32, in particular the extreme values of the diffusion parameters of the set of diffusion parameters 36.
In one embodiment, the first diffusion parameter 420 and the second diffusion parameter 426 are selected from the 5-th and 95-th percentile of the set of diffusion parameters, respectively, preferably from the 2-th and 98-th percentile. Here, the first diffusion parameter 420 and the second diffusion parameter 426 are each the median value of their respective quantile.
When the number of punctures to be performed of the procedure data 162 is greater than 2, the selector 40 determines additional diffusion parameters 424 in addition to the first diffusion parameter 420 and the second diffusion parameter 426. The respective values of these additional diffusion parameters 424 are between the value of the first diffusion parameter 420 and the value of the second diffusion parameter 426.
In one embodiment, the selector 40 chooses the diffusion parameters 424 comprised between the first diffusion parameter 420 and the second diffusion parameter 426 such that the diffusion parameters of the subset of diffusion parameters 42 have values that are equidistant from each other in twos. To do so, the selector 40 applies the equation (1) in Appendix A, where Dmin is the value of the first diffusion parameter 420, Dmax is the value of the second diffusion parameter 426, N is the number of punctures to be performed in a biopsy procedure of the procedure data 162, j is an integer between 1 and the number of punctures to be performed inclusive, and Dj is the value of the j-th diffusion parameter of the subset of diffusion parameters 42. The subset of diffusion parameters 42 thus obtained is represented in
Alternatively, when the set of diffusion parameters 36 has a distribution of its diffusion parameters around the median value M of the diffusion parameters of the set of diffusion parameters 36 that is highly skewed, the selector 40 determines the diffusion parameters 424 according to formula (2) of Appendix A, where Dini is determined according to formula (3) of Appendix A and Delta is determined according to formula (4) of Appendix A. This alternative allows less weight to be given to the longest tail of the diffusion parameter distribution of the set of diffusion parameters 36, such that the subset of diffusion parameters 42 is representative of the diffusion parameter distribution of the set of diffusion parameters 36. This skewness may, for example, occur if there is too much noise in the raster image 12, which disturbs small diffusion parameter values more than high diffusion parameter values. The determination of this skewness, and thus the choice of this skewed method or the equidistance method described above, can be carried out manually, or automatically by the selector 40.
In one embodiment, the selector 40 calculates the skewness coefficient of the set of diffusion parameters 36. If this skewness coefficient exceeds a certain threshold, the selector 40 determines the diffusion parameters 424 according to formula (2). Otherwise, the selector 40 applies formula (1).
The selector 40 returns a subset of regions 44 comprising regions 46 derived from the set of regions 32. Each region 46 of the subset of regions 44 is associated with a diffusion parameter substantially equal to one of the diffusion parameters of the subset of diffusion parameters 42. Furthermore, in the example described here, all of the diffusion parameters are each associated with at least one of the regions 46. Alternatively, for each diffusion parameter of the subset of diffusion parameters 42, the subset of regions 44 may comprise all of the regions 34 of the set of regions 32 associated with said diffusion parameter of the subset of diffusion parameters 42. In the embodiment where the set of diffusion parameters is sliced into quantiles, the diffusion parameter of each region 46 of the subset of regions 44 is included in the same quantile as the diffusion parameter of the associated subset of diffusion parameters 42.
The subset of regions 44 comprises regions 46 having diffusion parameter values well distributed in the set of diffusion parameters 36. In fact, a sampling from a plurality of these regions 46 with diffusion parameter values being different in twos is representative of all regions 34.
However, the tumour 200 has a number of features, such as the shape of the set of constraints 300 or the presence of necrotic regions 500. Based on these features, the guide 50 determines which regions of the subset of regions 44 should be sampled, and how to perform the sampling.
To do so, the guide 50 determines from the subset of regions 44, puncture parameter set data 52, equal in number to the number of punctures to be performed of the procedure data 162. Each set of puncture parameters of the puncture parameter set data 52 comprises a puncture depth 520, a puncture orientation 522 and a puncture entry point 524, represented in the example described here in
In the example represented in
In some cases, the guide 50 may determine that the subset of regions 44 is incompatible with the set of constraints 300. This means that it is not possible to perform as many punctures as originally intended meeting the sets of constraints 300. In this case, the selector 40 is run to determine a new subset of diffusion parameters 42 comprising one less diffusion parameter than previously, and then a subset of regions 44. This is repeated until the guide 50 determines a subset of regions 44 compatible with the set of constraints 300, or it is determined that the biopsy procedure is not possible.
In another example, represented in
The puncture parameter set data 52 determined by the guide 50 is provided to a surgeon or radiologist who is then able to carry out much more relevant biopsies than if he or she were to select them by hand or randomly. Furthermore, as the number of biopsies performed is small, less than 4 in the example described here, the medical risks are greatly reduced.
Reference is now made to
In the embodiment described here, the imaging device 1 retrieves the result of a biopsy procedure carried out according to the puncture parameter set data 52. This result is stored in the memory as puncture data 58. From the puncture data 58, the imaging device 1 reconstructs a density map of the tumour 200.
To do so, the imaging device 1 comprises a correlator 60 and a constructor 70. The correlator 60 derives from the puncture data 58 a correlation between the diffusion parameter values of the set of diffusion parameters 36 and cell density values in the tumour 200. The cell density may be a total cell density or a cancer cell density. The constructor 70 uses this correlation to reconstruct a cell density map within the tumour 200.
The puncture data 58 comprises a set of cell densities. Each cell density is associated with one of the puncture zones 526 of the puncture parameter set data 52. Each density has a value representing a cell density of the puncture zone 526 associated thereto.
This puncture data 58 may be obtained after a biopsy procedure by a surgeon or radiologist and after analysis of samples taken from the puncture zones 526 of the puncture parameter set data 52. The analysis of this puncture data 58 can be carried out by a histology method as described in Yin, Y., Sedlaczek, O., Müller, B., Warth, A., González-Vallinas, M., Lahrmann, B., Grabe, N., Kauczor, HU., Breuhahn, K., Vignon-Clementel, IE., Drasdo, D. (2018). “Tumour cell load and heterogeneity estimation from diffusion-weighted MRI calibrated with histological data: an example from lung cancer.” IEEE transactions on medical imaging, 37(1), 35-46.
The correlator 60 receives the puncture data 58 and associates with each cell density value associated with one of the puncture zones 526 the diffusion parameter value of the region 46 in which said puncture zone 526 is substantially included. The correlator 60 determines correlation data 62 that associates cell density values with diffusion parameter values of the set of diffusion parameters 36.
In the embodiment described here, the correlator 60 determines pairs each associating a diffusion parameter value and a cell density value associated with one of the puncture zones of the puncture parameter set data 52. Based on these pairs, the correlator 60 performs a linear interpolation, the result of which relates the diffusion parameter values of the set of diffusion parameters 36 to cell density values. The correlator 60 determines correlation data 62 based on the result of the linear interpolation.
Alternatively, the interpolation carried out by the correlator 60 may be a strictly monotonic polynomial interpolation, or other strictly monotonic function.
The constructor 70 determines density data 72 from the processing image 22 and the correlation data 62. The density data 72 comprises a density image. Each pixel of the density image corresponds to a region of the cross-sectional image to which a pixel of the processing image 22 also corresponds. The constructor 70 associates with each pixel of the density image the cell density value associated in the correlation data 62 with the diffusion parameter value of said pixel of the processing image 22.
Thus, from the raster image 14, the imaging device 1 determines here a reduced number of samplings to be carried out, and then, based on the result of these samplings, is able to reconstruct an accurate map of the cell density of the tumour 200.
Reference is made to
In one embodiment, not represented in the figures, the imaging device 1 further comprises an estimator arranged to determine a tumour load of the tumour 200. Here, the estimator carries out a numerical area integration over the entire area of the density image from which it derives the tumour load.
Here, the imaging device 1 is thus able to determine the tumour load of the tumour 200 from the result of samplings according to the puncture parameter set data 52, and to visualise and/or quantify the tumour heterogeneity.
Reference is made to
In order to show the interest of the invention compared to the existing methods, a histogram derived from biopsies selected by means of the device according to the invention is represented in
These figures show that the puncture parameter set data 52 determined by the imaging device 1 allows punctures to be carried out for determining a tumour load that is much closer to the actual tumour load than random punctures. These figures also show that the tumour load measurements made on the basis of the puncture parameter set data 52 have much less scatter than the measurements made with random punctures.
The use of puncture parameter set data 52 determined by the imaging device 1 therefore allows the correlator 60, the constructor 70 and the estimator to determine a much more reliable tumour load, that is, closer to the reference value 170.
The imaging device 1 allows a surgeon or radiologist to accurately and robustly measure a tumour load from a very small number of punctures, typically three to five times fewer than conventional procedures.
The correlator 60, the constructor 70 and the estimator are programs run by the computer processor. The alternatives described for the upsampler 20, the slicer 30, the selector 40 and the guide 50 can also be applied to the correlator 60, the constructor 70 and the estimator.
In the foregoing, the slicer 30 determines regions 34 having a size greater than or substantially equal to the puncture zone of a biopsy needle. Alternatively, the slicer 30 could determine regions 34 having dimensions less than the puncture length, but greater than the puncture diameter. For biopsies carried out on this basis, histological measurements would then be carried out only on the part of the samplings corresponding to the regions 540, 542, 544.
Alternatively, the device may also receive in memory 10 data from other types of non-invasive imaging whose pixels of the raster images 14 are associated, instead of the diffusion parameter in the case of DWI, with a parameter that reflects cell density.
Dj=Dmin+(j−1)(Dmax−Dmin)/(N−1) (1)
Dj=Dini+Delta(j−1) (2)
Dini=M−min(|M−Dmin|,|M−Dmax|) (3)
Delta=2(M−Dini)/(N−1) (4)
Although the present disclosure has been described with reference to one or more examples, workers skilled in the art will recognize that changes may be made in form and detail without departing from the scope of the disclosure and/or the appended claims.
| Number | Date | Country | Kind |
|---|---|---|---|
| FR2006763 | Jun 2020 | FR | national |
This Application is a Section 371 National Stage Application of International Application No. PCT/FR2021/051168, filed Jun. 24, 2021, the content of which is incorporated herein by reference in its entirety, and published as WO 2021/260334 on Dec. 30, 2021, not in English.
| Filing Document | Filing Date | Country | Kind |
|---|---|---|---|
| PCT/FR2021/051168 | 6/24/2021 | WO |
