WO2025099779A1 - Method for pulmonary segmentation from a thoracic volumetric scan - Google Patents

Abstract

Computer-implemented method (10) for pulmonary segmentation from a thoracic volumetric scan and computer program suitable to implement such method.

Description

“METHOD FOR PULMONARY SEGMENTATION FROM A THORACIC VOLUMETRIC SCAN”

FIELD OF THE INVENTION The present invention concerns a method for pulmonary segmentation from a patient’s thoracic volumetric scan acquired through computed tomography (CT).

BACKGROUND OF THE INVENTION

In the field of segmentation of thoracic CT scans, various methods are known, generally implemented as computer programs, aimed at the segmentation of the lung parenchyma (or tissue) and/or the segmentation and classification of altered density volumes within the parenchyma. Most currently usable software products allow to achieve segmentation of lung parenchyma volumes with altered density, such as so-called ground glass areas, honeycomb lung and low attenuation areas, while they are usually inadequate in the segmentation of so-called consolidated tissue, the density of which is such that is indistinguishable from the tissues outside the lung parenchyma in a volumetric CT scan in the event this tissue is contiguous with the pulmonary border. Examples of such products are Annotate by ART- PlanTM (Therapanacea TM), Pinnacle3 Treatment Planning (Philips S.p.A.), MIM Maestro® (Tema Sinergie S.p.A.) and Caliper (© 1998-2023 Mayo Foundation for Medical Education and Research).

There are also products, such as Coreline AVIEW COPD and syngo® LungCARE CT (Siemens Healthcare GmbH ©2023), able to provide the segmentation of volumes of consolidated tissue thanks to the use of supervised machine learning algorithms based on atlases of CT volumetric scans, or deep learning algorithms, while there remain limitations in the contact sites of the consolidated tissue with the chest wall, diaphragm and mediastinum. These products require, in addition to the CT volumetric scan to be analyzed, a database of volumetric scans to train the algorithm, which cannot be controlled by the user. As a result, algorithms based on machine learning or deep learning are more complex and prone to reliability and stability problems related to possible incompatibility of the database with the specific volumetric scan in question.

There is therefore the need to perfect a method for pulmonary segmentation from a thoracic volumetric scan that can overcome at least one of the disadvantages of the state of the art.

In order to do this, it is necessary to solve the technical problem of distinguishing consolidated tissue from tissues of similar density contiguous with, but external to, the lung parenchyma, starting from the volumetric CT scan to be analyzed alone.

In particular, one purpose of the present invention is to provide a method for pulmonary segmentation from a thoracic volumetric scan that allows for the segmentation of volumes of the lung parenchyma with altered density, including, if present, consolidated tissue in whatever position it is found. Another purpose of the present invention is to provide a method for pulmonary segmentation from a thoracic volumetric scan that does not require databases or atlases of training scans.

A further purpose of the present invention is to provide a method for pulmonary segmentation from a thoracic volumetric scan that is rapid and efficient. Another purpose of the present invention is to provide a method for pulmonary segmentation from a thoracic volumetric scan whose result can be explained and deduced on the basis of the operations performed in applying the method, avoiding the “black box” effect typical of methods that use artificial intelligence (for example, methods based on machine learning or deep learning). The Applicant has devised, tested and embodied the present invention to overcome the shortcomings of the state of the art and to obtain these and other purposes and advantages.

SUMMARY OF THE INVENTION

The present invention is set forth and characterized in the independent claims. The dependent claims describe other characteristics of the present invention or variants to the main inventive idea.

In accordance with the above purposes and to resolve the technical problem described above in a new and original way, also achieving considerable advantages compared to the state of the prior art, some embodiments described here concern a computer- implemented method for pulmonary segmentation from a thoracic volumetric scan.

In accordance with one aspect of the present invention, the method comprises a step of first densitometric threshold segmentation of the left and right lung, thus dividing a lung volume of the volumetric scan into at least one low-density region and one high-density region, which are delimited overall by a preliminary pulmonary border.

The method also comprises a radiomic analysis step, following or preceding the first segmentation step as above, which provides to obtain, for each pixel of each axial section of the volumetric scan, the value of at least one radiomic feature, thus obtaining, for each axial section, at least one radiomic feature image. By radiomic feature we mean a quantitative characteristic of a medical image.

The method then comprises a lateral and posterior border identification step, which provides to identify a lateral and posterior pulmonary border on the basis of a combination of information deriving from: a geometric modeling of the shape of the lung parenchyma; the preliminary pulmonary border; a threshold densitometric segmentation of the rib cage; and the radiomic analysis step as above. The method also comprises a step of reconstructing a lung volume, which provides to combine at least the results obtained in the previous steps in order to identify a pulmonary border along the mediastinal profile, simultaneously reconstructing the entire volume of the lung parenchyma, thus obtaining a pulmonary segment. The method then comprises a second threshold segmentation step, which provides to segment the above mentioned pulmonary segment further, through the application of at least two densitometric thresholds, delimiting at least one new low-density region, one new high-density region and one very high density region. The method also comprises a vessel removal step, which provides to remove any pulmonary blood vessels at least from the very high density region, thus obtaining a new very high density region. The method disclosed above offers the advantage of not requiring databases or atlases of training volumetric scans.

In accordance with another aspect of the present invention, the method comprises another inferior border identification step aimed at identifying, for each axial section of the thoracic volumetric scan, an inferior pulmonary border on the basis of a combination of information deriving from: the anatomy (relative positioning of the lungs, diaphragm and underlying organs and tissues); the preliminary pulmonary border; a densitometric threshold segmentation of tissues underlying the diaphragm; and the radiomic analysis step as above. In accordance with another aspect of the present invention, the new low-density region corresponds to healthy lung tissue, the new high-density region corresponds to the ground glass and the new very high-density region corresponds to consolidated tissue. The method therefore advantageously gives a segmentation of volumes of the lung parenchyma with altered density, including, if present, the consolidated tissue.

In accordance with another aspect of the present invention, the radiomic analysis step provides to calculate the radiomic features of skewness, entropy and standard deviation. The Applicants have found that the three radiomic features listed above have greater discriminatory power and are less burdensome in terms of computational time than other known radiomic features.

In accordance with another aspect of the present invention, the radiomic analysis step provides to normalize each radiomic feature image by dividing it by its maximum value, filtering each normalized radiomic feature image with a threshold value and finally combining the radiomic feature images by intersection, thus obtaining a synthetic radiomic image.

In accordance with another aspect of the present invention, the lateral and posterior border identification step provides, separately for each right and left lung, to identify first points of an internal border of the ribs, then interpolate these first points with a parametric three-dimensional cylindroid spline function, obtaining a first border surface. Second points of the preliminary pulmonary border which are distant from the first border surface by less than a certain threshold value are then added to the set of points to be interpolated. These first and second points are then interpolated with a parametric three-dimensional cylindroid spline function, obtaining a second border surface. The identification step then provides to add to the set of points to be interpolated third points of the at least one radiomic feature image, or possibly of the synthetic radiomic image, which are distant from the second border surface by less than a certain threshold value, so as to then interpolate the first, second and third points with a parametric three-dimensional cylindroid spline function, obtaining a third border surface.

In accordance with another aspect of the present invention, the other inferior border identification step described above operates on axial sections of the thoracic volumetric scan, instead of the entire scanned volume. This other identification step can provide to identify a right part of the inferior pulmonary border as the union of a plurality of curves of the right border of the diaphragm corresponding to one border of the liver. These curves can be identified by interpolating points of the at least one radiomic feature image, or possibly of the synthetic radiomic image. The other identification step can also provide to identify a left part of the inferior pulmonary border as the union of a plurality of curves of the left border of the diaphragm corresponding to a border of the spleen and of the perigastric adipose tissue. These curves can be identified by segmenting the spleen from each axial section that contains it, through a region growing technique and taking into account the discontinuities in density which are due to the perisplenic and perigastric adipose tissue, and of the preliminary pulmonary border and the lateral and posterior pulmonary border.

In accordance with another aspect of the present invention, the lung volume reconstruction step described above provides to: join the low and high density regions into a single region; subject this single region to a morphological closing operation with a spherical structuring element, using the lateral and posterior pulmonary borders as above and a possible inferior pulmonary border, which have been previously identified to establish the boundary of the extension of the closing operation toward the exterior of the lung, thus obtaining the above mentioned pulmonary segment.

Advantageously, the morphological closing operation with a spherical structuring element with a suitably chosen radius allows to exclude the mediastinum, which has regular borders, from the pulmonary segment, and to include instead the consolidated tissue, which has irregular borders. In accordance with another aspect of the present invention, the vessel removal step as above provides to apply a vesselness filter to the above mentioned pulmonary segment, thus identifying any portions of the pulmonary segment having high tubular morphology, which correspond to presumed blood vessels. The vessel removal step then provides to associate with each voxel of each presumed blood vessel a distance of the voxel from the mediastinum and a diameter of the presumed blood vessel. The presumed large diameter vessels close to the mediastinum are then identified as blood vessel. The blood vessels thus identified are then removed from the pulmonary segment, in particular from the very high density region mentioned above, obtaining the new very high density region.

In accordance with another aspect of the present invention, the method comprises a step of delimiting a region of interest, prior to the other steps described above, which provides to divide the thoracic volumetric scan to be segmented into a region of interest, the voxels of which are analyzed in the subsequent steps of the method, and an external region, which is not analyzed in the subsequent steps of the method. This advantageously allows to reduce the method’s calculation time.

Some embodiments described here also concern a computer program containing instructions which, if executed by a computer, cause the execution of the method disclosed above.

DESCRIPTION OF THE DRAWINGS

These and other aspects, characteristics and advantages of the present invention will become apparent from the following description of some embodiments, given as a non-restrictive example with reference to the attached drawings wherein:

- fig. 1 is a flowchart of an embodiment of a method for pulmonary segmentation from a thoracic volumetric scan according to the present invention;

- fig. 2 shows an example of a step of delimiting a region of interest of the method of fig. 1; - figs. 3a-b show an example of a first segmentation step of the method of fig. 1 ;

- figs. 4a-d show examples of results of a radiomic analysis step of the method of fig. 1;

- figs. 5a-g show an example of a step of identifying a lateral and posterior pulmonary border of the method of fig. 1 ; - figs. 6a-e show an example of the right diaphragm identification sub-step of an inferior border identification step of the method of fig. 1 ;

- figs. 7a-e show an example of the left diaphragm identification sub-step of an inferior border identification step;

- figs. 8a-d show an example of lung volume reconstruction (figs. 8a-c) and second segmentation (fig. 8d) steps of the method of fig. 1;

- figs. 9a-b show an example of a vesselness filter application in a vessel removal step of the method of fig. 1;

- figs. lOa-b show an example of a vessel removal step of the method of fig. 1. We must clarify that the phraseology and terminology used in the present description, as well as the figures in the attached drawings also in relation as to how described, have the sole function of better illustrating and explaining the present invention, their purpose being to provide a non-limiting example of the invention itself, since the scope of protection is defined by the claims.

To facilitate comprehension, the same reference numbers have been used, where possible, to identify identical common elements in the drawings. It is understood that elements and characteristics of one embodiment can be conveniently combined or incorporated into other embodiments without further clarifications. DESCRIPTION OF SOME EMBODIMENTS OF THE PRESENT INVENTION

With reference to fig. 1, a computer-implemented method 10 is described for pulmonary segmentation from a thoracic volumetric scan according to the present invention. The volumetric scan to be segmented can for example be a volumetric scan acquired through computed tomography (CT) and thus represent information on the density of the tissues contained in the scanned volume. Typically, a thoracic CT volumetric scan has voxels of (1 x 1 x 3) mm³, but different formats are also possible.

In some embodiments, the method 10 can be semi-automatic, that is, requiring at input, in addition to the volumetric scan to be segmented, one or more additional pieces of information from a user.

The method 10 can comprise a preliminary step 11 of delimiting a region of interest Rl. The step 11 of delimiting the region of interest R1 provides to divide the thoracic volumetric scan into a region of interest Rl, the voxels of which will be analyzed in the subsequent steps of the method 10, and an external region R2, which will not be analyzed in the subsequent steps of the method 10 (fig. 2). The step 11 of delimiting the region of interest Rl is not essential to the realization of the segmentation, but allows to reduce the calculation time necessary to perform the method 10.

According to possible embodiments, the step 11 of delimiting the region of interest Rl provides to segment the rib cage, calculate its convex hull and consider as a region of interest Rl the region of the thoracic volumetric scan internal to the convex hull (fig. 2). The segmentation of the rib cage can take place, for example, by means of densitometric threshold. In the present document, densitometric threshold values are given in Hounsfield Units (HU), a unit of measurement of radiometric density commonly used for CT volumetric scans. Typical threshold values for different tissues can be found in Lederer D. J. et al., “Cigarette smoking is associated with subclinical parenchymal lung disease: the Multi-Ethnic Study of Atherosclerosis (MESA)-long study”, Am. J. Respir. Crit. Care Med. (2009), 180 (5): 407-14. DOI: 10.1164/rccm.200812- 1966OC; Colombi D. et al., “Well-aerated Lung on Admitting Chest CT to Predict Adverse Outcome in COVID-19 Pneumonia”, Radiology (2020), 296 (2): E86- E96. DOI: 10.1148/radiol.2020201433; Herbert L., “Prostatic Diseases”, W.B. Saunders Company (2000), ISBN 9780721674162 -Page 83; Kuntz, E. and Kuntz, H.-D., “Hepatology, Principles and Practice: History, Morphology, Biochemistry, Diagnostics, Clinic, Therapy”, Springer Science & Business Media (2006), ISBN 9783540289777 - Page 210; Birur, N. P. et al., “Comparison of gray values of cone-beam computed tomography with hounsfield units of multislice computed tomography: An in vitro study”, Indian Journal of Dental Research (2017), 28 (1):

66-70, DOI: 10.4103/ijdr.IJDR_415_16.

For example, to segment bones such as the bones that make up the rib cage, a lower threshold of about 90 HU can be used, that is, classify the voxels having a density greater than about 90 HU as bones. In some embodiments, the step 11 of delimiting a region of interest R1 takes at input, in addition to the thoracic volumetric scan, also an upper fiducial point Pl and a lower fiducial point P2, for example entered by the user, which are able to identify respective upper and lower axial sections of the volumetric scan within which the segmentation will be calculated. A first threshold segmentation step 12, possibly following the step 11 of delimiting a region of interest Rl, if present, provides to perform a densitometric threshold segmentation of the left and right lung. This segmentation can occur in a known manner to distinguish healthy lung tissue from altered density tissues based on appropriate densitometric thresholds. For example, the first segmentation step 12 can divide the lung volume of the volumetric scan into at least two regions, a low density region R3 (figs. 3a-b, dark gray) and a high density region R4 (figs. 3a-b, light gray).

For example, the first segmentation 12 can segment the healthy, low-density lung into the low-density region R3 and the so-called ground glass tissue, characterized by higher densities, into the high-density region R4. For example, the first segmentation step 12 can classify as healthy lung the regions R3 of the volumetric scan having a density comprised between about - 950 HU and - 750 HU, and as ground glass the regions R4 of the volumetric scan having a density comprised between about - 750 HU and - 250 HU.

Some embodiments can allow the user to manually exclude certain regions from the lung volume, such as low-density regions like the trachea, gastrointestinal cavities adjacent to the lungs, cavities caused by pathological conditions, low- density lung regions, or others. In some embodiments, it can be possible to exclude the volume belonging to the trachea starting from a tracheal fiducial point P3, that is, a point inside the trachea. The tracheal fiducial point P3 can be entered manually by the user.

According to possible embodiments, the first segmentation step 12 can provide to apply an erosion technique to the border of the low density region R3, or of the high density region R4, or both. An example of an erosion technique suitable for this purpose can be found in Serra J., “Image Analysis and Mathematical Morphology” (1982), ISBN 0-12-637240-3.

This erosion technique can solve the problem of the so-called partial volume effect, which occurs in correspondence with a separating border between two tissues when the space occupied by the density gradient at the interface between the two tissues is greater than the sizes of a voxel. In this case, the density values of the voxels contiguous with that border do not correspond to the density of either of the two tissues. An example of the resolution of the partial volume effect by applying an erosion technique to the high-density region R4 (fig. 3a) is shown in fig. 3b.

Please note that in the first segmentation step 12, it is not possible to identify any regions of consolidated tissue contiguous with the lung wall (fig. 3b, black arrows), since they have densities comparable to those of tissues outside the lung parenchyma.

The low- and high-density regions R3, R4 are delimited overall by a border that can be considered as a preliminary pulmonary border B 1 , which does not delimit the entire volume of the lung parenchyma, but only the low- and high-density regions R3, R4.

According to possible embodiments, the first segmentation step 12 can be followed by a radiomic analysis step 13. In possible alternative embodiments, the radiomic analysis step 13 can precede, rather than follow, the first segmentation step 12.

By radiomic analysis we mean an analysis of medical images aimed at extracting quantitative characteristics, also called radiomic features, from the analyzed images. Generally, for a given pixel of a medical image, a corresponding radiomic feature is computed within a predefined extension region centered on the given pixel. A detailed description of radiomic features in use is available, for example, in Lambin P. et al., “Radiomics: extracting more information from medical images using advanced feature analysis”, European Journal of Cancer (2012), 48 (4): 441 - 6 or in Zwanenburg, A. et al., “Image biomarker standardisation initiative - feature definitions”, eprint arXiv: 1612.07003 (2016). In some embodiments, the radiomic analysis step 13 provides to perform a radiomic analysis on axial sections of the volumetric scan to be segmented, in order to obtain, for each pixel of each axial section, the value of at least one radiomic feature.

For example, the radiomic analysis step 13 can provide to calculate the skewness index, or the entropy, or the standard deviation, or another radiomic feature, or two or more of these or other radiomic features.

The Applicants have selected radiomic features appropriate for pulmonary segmentation from a thoracic volumetric scan, in particular to identify consolidated tissue, performing a t-test study to highlight the radiomic features with the highest discrimination power and choosing among these the radiomic features least burdensome in terms of computational time. Skewness, entropy and standard deviation have therefore been chosen. It is understood that one or more alternative radiomic features, as well as different criteria for choosing radiomic features, are possible and methods that provide for them fall in any case within the scope of the present invention.

For each axial section of the thoracic volumetric scan, as many images can be obtained as there are radiomic features used, each being a so-called radiomic feature image containing, in each pixel, the value of a particular radiomic feature calculated within a predefined extension region centered on that pixel. For example, for each axial section of the thoracic volumetric scan, it is possible to calculate a skewness image, an entropy image and a standard deviation image. Some examples of radiomic feature images for a same axial section of a thoracic volumetric scan are shown in figs. 4a-c, in particular the skewness image (fig. 4a), the standard deviation image (fig. 4b) and the entropy image (fig. 4c).

In some embodiments, the radiomic analysis step 13 can provide to combine, if present, different radiomic feature images relating to a same axial section into a single synthetic radiomic image. As a non-limiting example, each radiomic feature image (figs. 4a-c) can be normalized by dividing it by its maximum value, filtered with an appropriate threshold value so as to better highlight the borders of the various anatomical regions, and finally combined with the other radiomic feature images by intersection (that is, using an AND logical operation). An example of such a synthetic radiomic image is shown in fig. 4d and its non-zero points are called “output points of the features”. The Applicants have found that the synthetic radiomic image highlights anatomical borders particularly effectively, including those of any consolidated tissue.

In those embodiments that provide the step 11 of delimiting the region of interest Rl, it is possible to exclude the external region R2 from the radiomic analysis, thus reducing processing times. In those embodiments in which the radiomic analysis step 13 follows the first segmentation step 12, it is possible to exclude all those areas already segmented in the first segmentation step 12 from the radiomic analysis, for example the low and high density regions R3, R4, thus further reducing processing times. Subsequently, a lateral and posterior border identification step 14 provides to identify a lateral and posterior border of the lung parenchyma, or lateral and posterior pulmonary border B2. The identification step 14 can exploit the shape and position of the rib cage bones (or ribs) that are located very close to the pleura, so that the border of the lung parenchyma can be correctly identified even in the presence of consolidated tissue contiguous with it.

The identification step 14 can be based on a combination of information deriving from: a geometric modeling of the shape of the lung parenchyma; the threshold densitometric segmentation of the lung parenchyma (preliminary pulmonary border Bl identified in the first segmentation step 12) and of the rib cage; and the radiomic analysis (one or more radiomic feature images calculated in the radiomic analysis step 13) in order to identify the lateral and posterior pulmonary border B2. In some embodiments, the identification step 14 can provide to separately identify a right lateral and posterior pulmonary border and a left lateral and posterior pulmonary border. The union of these borders corresponds to the lateral and posterior pulmonary border B2.

The identifications of the right and left lateral and posterior pulmonary borders are independent of each other, therefore it is possible to first perform the identification of the right lateral and posterior pulmonary border and then the identification of the left lateral and posterior pulmonary border, or vice versa, or the two identifications can be performed in parallel, if the specific implementation of method 10 provides for it. The example identification step 14 shown in figs. 5a-g refers to the lateral and posterior border of the right lung parenchyma. Obviously, the same method can also be applied to the left lung parenchyma, thus obtaining the lateral and posterior pulmonary border B2.

In some embodiments, the identification step 14 can comprise at least an interpolation of a set of first points P4 of an internal border of the ribs (fig. 5 a, gray points) with a suitable parametric function.

In some embodiments, the identification step 14 can provide one or more iterations of the interpolation, in each of which the set of points to be interpolated and/or one or more parameters of the interpolating function may be changed in order to refine the interpolation, making the lateral and posterior pulmonary border B2 thus identified more adherent to the actual border of the lung.

For example, the identification step 14 can provide to model the border of the lung parenchyma of each lung with a spline function representing a parametric three-dimensional cylindroid surface. In some embodiments, the identification step 14 can first provide to identify first points P4 of an internal border of the ribs (fig. 5a), then to interpolate these first points P4 with a parametric three-dimensional cylindroid spline function (fig. 5d) in a first interpolation sub-step, obtaining a first border surface SI (fig. 5d, gray line).

Some embodiments provide to place a limit on the chi-square of the interpolation operation normalized to the number of points. The limit value, also called the smoothness parameter of the spline function, is associated with a spline function that is all the “smoother” and more regular in shape, the greater the parameter. Preferred embodiments provide to perform the first interpolation substep with a first relatively high smoothness parameter, in order to avoid possible artifacts related to the identification of an erroneous pulmonary border too far away from the true one. For example, the first smoothness parameter can be comprised between about 4 and 5, for example it can be 4.5.

In some embodiments, the identification step 14 can then provide to add to the set of points to be interpolated, for example the first points P4 of an internal border of the ribs, one or more second points P5 to be interpolated (fig. 5b), for example points of the preliminary pulmonary border B 1 obtained in the first segmentation step 12.

In some embodiments, these points of the preliminary pulmonary border Bl can be added to the set of points to be interpolated only if their distance from the first border surface SI is less than a certain threshold value. This threshold value can for example be of the order of the larger linear dimension of a voxel. For example, for a typical thoracic CT volumetric scan, this threshold value can be approximately 3 mm. A second interpolation sub-step then provides to interpolate the first points P4 of an internal border of the ribs and the second points P5 (fig. 5b, gray points) selected from the preliminary pulmonary border Bl with a parametric three-dimensional cylindroid spline function, obtaining a second border surface S2 (fig. 5e, gray line).

Preferred embodiments provide to perform the second interpolation sub-step with a second smoothness parameter, lower than the first smoothness parameter, in order to bring the second border surface S2 closer to the real pulmonary border. For example, the second smoothness parameter can be equal to about 0.8 - 0.9 times the first smoothness parameter, for example 0.82 times the first smoothness parameter.

In some embodiments, the identification step 14 can then provide to add to the set of points to be interpolated, for example the first points P4 and the second points P5 (fig. 5b), one or more third points P6 to be interpolated (fig. 5c), for example points of the at least one radiomic feature image obtained in the radiomic analysis step 13 (figs. 4a-c), or points of the synthetic radiomic image (fig. 4d), if this has been calculated. In some embodiments, a point of the radiomic feature image or of the synthetic radiomic image can be added to the set of points to be interpolated only if its distance from the second border surface S2 is less than a certain threshold value. For example, this threshold value can be about 3 mm. A third interpolation substep then provides to interpolate the first points P4 of an internal border of the ribs, the second points P5 selected from the preliminary pulmonary border Bl and the third points P6 selected from the radiomic feature image or from the synthetic radiomic image (fig. 5c, gray points) with a parametric three-dimensional cylindroid spline function, obtaining a third border surface S3 (fig. 5f, gray line). Preferred embodiments provide to perform the third interpolation sub-step with a third smoothness parameter lower than the second smoothness parameter, in order to bring the third border surface S3 even closer to the real pulmonary border. For example, the third smoothness parameter can be equal to about 0.8 - 0.9 times the second smoothness parameter, for example 0.84 times the second smoothness parameter. Fig. 5g shows an example of a third border surface S3 as a light gray grid and the corresponding points P4, P5 and P6 interpolated to obtain it as dark gray points. The third border surface S3 of figs. 5f-g exemplifies a right part of the lateral and posterior pulmonary border B2. The black arrows of figs. 5e-f illustrate how the subsequent interpolation sub-steps described above allow to achieve border surfaces SI, S2, S3 gradually closer to the real pulmonary border.

According to possible embodiments, the identification step 14 can optionally be followed by another inferior border identification step 16.

The other inferior border identification step 16 is aimed at identifying a border of the diaphragm below the lung parenchyma, that is, an inferior pulmonary border B3. The other identification step 16 is necessary when there is consolidated tissue adjacent to the diaphragm; otherwise, the first segmentation step 12 is sufficient to identify the inferior pulmonary border B3.

The other identification step 16 can be based on a combination of information deriving from: the anatomy (relative positioning of lungs, diaphragm, and underlying organs or tissues); threshold densitometric segmentation of the lung parenchyma (preliminary pulmonary border Bl identified in the first segmentation step 12) and of the tissues underlying the diaphragm; and the radiomic analysis (one or more radiomic feature images calculated in the radiomic analysis step 13) to identify the inferior pulmonary border B3.

In some embodiments, the other identification step 16 does not operate on the entire scanned volume, but on axial sections of the thoracic volumetric scan.

The other identification step 16 can comprise a sub-step of identifying the right diaphragm (figs. 6a-e) and a sub-step of identifying the left diaphragm (figs. 7a-e), which are diversified according to the different anatomy of the right and left part of the body below the lungs. The sub-steps of identifying the right and left diaphragm are independent of each other, therefore it is possible to perform the right diaphragm identification sub-step first and then the left diaphragm identification sub-step, or vice versa, or the two sub-step can be performed in parallel, if the specific implementation of the method 10 allows parallel executions.

According to possible embodiments, the right diaphragm identification sub-step takes as input a central hepatic fiducial point P7, corresponding to the center of the liver’s apex, and a border hepatic fiducial point P8, corresponding to a border point of the liver’s apex (fig. 6a). The hepatic fiducial points P7, P8 can be indicated by the user in the first axial section of the volumetric scan, in the cranio-caudal direction, in which the liver appears. In other words, the hepatic fiducial points P7, P8 are chosen in correspondence with the apex of the hepatic dome.

In some embodiments, the right diaphragm identification sub-step provides to calculate a circumference centered on the central hepatic fiducial point P7 and with a radius equal to the distance between the hepatic fiducial points P7, P8. This circumference is used to select a set of points to be interpolated in order to obtain a first preliminary curve Cl of the right border of the diaphragm. For example, this circumference can be used to define along it a circular crown with an appropriate width, for example a width equal to the larger linear dimension of a voxel of a typical thoracic volumetric scan, for example 3 mm.

In some embodiments, the set of points to be interpolated in order to obtain the first preliminary curve Cl can be selected from the at least one radiomic feature image calculated in the radiomic analysis step 13 (figs. 4a-c), or from the synthetic radiomic image, if this was calculated in the radiomic analysis step 13 (fig. 4d), by selecting the points inside the circular crown described above. For example, the points to be interpolated in order to obtain the first preliminary curve Cl can be the points of the skewness image (fig. 4a) located inside the circular crown. This set of points can be interpolated with a spline function, thus obtaining the first preliminary curve Cl (fig. 6b). In some embodiments, the right diaphragm identification sub-step can provide to make the first preliminary curve Cl adhere along the preliminary pulmonary border B 1 at its points that are distant from the first preliminary curve Cl by no more than a certain threshold value, for example

3 mm. In this way, a first definitive curve Cl ’ of the right border of the diaphragm is obtained (fig. 6c).

The right diaphragm identification sub-step then provides to calculate additional curves of the right border of the diaphragm, one for each axial section of the next volumetric scan in the cranio-caudal direction. The calculation is done starting from the previous curve of the right border of the diaphragm. For example, for the axial section (figs. 6d-e) immediately below or immediately following, in the cranio-caudal direction, the axial section (figs. 6a-c) containing the apex of the hepatic dome, the first definitive curve Cl ’ is replicated on this inferior axial section (fig. 6d). The replicated first definitive curve Cl ’ is expanded radially, in a non-uniform manner, based on the density gradient present in the radial expansion direction, thus obtaining a second preliminary curve C2 of the right border of the diaphragm (fig. 6d).

In some embodiments, a set of points to be interpolated in order to obtain a second definitive curve C2’ of the right border of the diaphragm can be selected from the at least one radiomic feature image calculated in the radiomic analysis step 13 (figs. 4a-c), or from the synthetic radiomic image, if this was calculated in the radiomic analysis step 13 (fig. 4d), by selecting its points in a vicinity of the second preliminary curve C2, for example the points which are distant no more than 1.5 mm from the second preliminary curve C2. For example, the points to be interpolated in order to obtain the second definitive curve C2’ can be the points of the skewness image (fig. 4a) which are distant no more than 1.5 mm from the second preliminary curve C2. This set of points can be interpolated with a spline function, thus obtaining another curve of the right border of the diaphragm, which can be made to adhere to the preliminary pulmonary border Bl at its points that are distant from the curve itself by no more than a certain threshold value, for example 3 mm. In this case, the second definitive curve C2’ is obtained (fig. 6e).

The right diaphragm identification sub-step then provides to repeat the method described above for each following axial section of the thoracic volumetric scan in the cranio-caudal direction. A surface formed by the union of the definitive curves of the right border of the diaphragm thus calculated represents a superior border of the liver, or equivalently the right inferior border of the lung parenchyma, therefore each definitive curve of the right border of the diaphragm is part of a right part of the inferior pulmonary border B3. The gray arrows in figs. 6c-d indicate regions of the consolidated tissue correctly included in the lung volume.

According to possible embodiments, the left diaphragm identification sub-step first provides to apply a densitometric threshold segmentation in order to identify the border portion of the lung parenchyma contiguous with the perigastric adipose tissue. The left diaphragm identification sub-step then provides to identify the border portion of the lung parenchyma contiguous with the spleen (perisplenic diaphragmatic profile). In some embodiments, the left diaphragm identification sub-step can provide that the user manually selects the axial section of thoracic volumetric scan immediately prior cranially, in the cranio-caudal direction, to the axial section in which the adjacency of the consolidated tissue to the perisplenic diaphragmatic profile is detected (fig. 7a, gray arrow), thus indicating in this axial section a splenic fiducial point P9 inside the spleen. The left diaphragm identification sub-step can provide to segment the spleen starting from the splenic fiducial point P9, from the preliminary pulmonary border Bl, from the left lateral and posterior pulmonary border and from the pulmonary border contiguous with the perigastric adipose tissue. In some embodiments, the segmentation of the spleen is performed through a region growing technique (white region in fig. 7b). This segmentation can be used to calculate a spline function corresponding to a curve of the left border of the diaphragm, which can be made to adhere to the preliminary pulmonary border B 1 at the points of the latter whose distance from the curve of the left border itself is no larger than the larger linear dimension of a voxel of a typical pulmonary volumetric scan, for example no larger than 3 mm (fig. 7c), thus obtaining a third definitive curve C3 of the left border of the diaphragm.

The left diaphragm identification sub-step can then provide to segment the spleen in the subsequent, or underlying, axial sections of the thoracic volumetric scan in the caudal cranial direction. For each subsequent axial section (figs. 7d-e), the left diaphragm identification sub-step can provide to replicate the segmentation of the previous axial section (figs. 7a-c), then expand it (fig. 7d) in order to then identify inside it a discontinuity in density caused by the perisplenic and perigastric adipose tissue, also taking into account the preliminary pulmonary border B 1 , the left lateral and posterior pulmonary border and the border contiguous with the perigastric adipose tissue. Similarly to what described for the third definitive curve C3, the new segmentation can be used to calculate a spline function corresponding to a curve of the left border of the diaphragm, which can be made to adhere to the preliminary pulmonary border Bl at the points of the latter whose distance from the curve of the left border is no larger than the larger linear dimension of a voxel of a typical pulmonary volumetric scan, for example no larger than 3 mm (fig. 7e), thus obtaining a fourth definitive curve C4 of the left border of the diaphragm.

The union of the definitive curves of the left border of the diaphragm thus calculated corresponds to a border of the spleen. A surface formed by the border contiguous with the perigastric adipose tissue and by the union of the definitive curves of the left border represents the left inferior border of the lung parenchyma, therefore each definitive curve of the left border of the diaphragm is part of a left part of the inferior pulmonary border B3. Subsequently, a lung volume reconstruction step 17 (figs. 8a-d) provides to combine the results obtained in the previous steps in order to identify a pulmonary border along the mediastinal profile, simultaneously reconstructing the entire volume of the lung parenchyma. In some embodiments, the low-density R3 and high-density R4 regions obtained in the first segmentation step 12 (fig. 8a) are joined in a single segment, or single region R5 (fig. 8b).

The single region R5 can be subjected to a morphological closing operation, or closing with a structuring element (or kernel) suitably chosen to obtain a single connected lung volume including therein the cavities present in the single region R5 which are due to the presence of consolidated tissue (examples of such cavities, or regions of consolidated tissue which are not included, are indicated in fig. 8b by black arrows), simultaneously also delimiting the boundary between lung and mediastinum, thus obtaining a pulmonary region or segment R6 (fig. 8c). As known for example from Serra J., “Image Analysis and Mathematical Morphology” (1982), ISBN 0-12-637240-3, the morphological closing operation is capable of attaching to the region subjected to the morphological closing operation cavities present within that region with a diameter up to about twice the size of the associated structuring element. Therefore, a spherical structuring element having a diameter smaller than the radii of curvature typical of the mediastinal surface allows to exclude the mediastinum from the pulmonary segment R6. A structuring element having sizes comparable with the typical sizes of the irregularities of the borders of the regions of consolidated tissue instead allows, through morphological closure, to attach these regions to the pulmonary segment R6.

The Applicants have verified that, for a typical thoracic CT volumetric scan, a spherical structuring element with a radius of 10 mm allows to not lose the morphological details of the border of the mediastinum, such as the profile of the aorta for example, while simultaneously allowing to include the consolidated tissue in the pulmonary segment R6.

In some embodiments, the lateral and posterior pulmonary border B2 and the possible inferior pulmonary border B3, identified in the lateral and inferior border identification step 14 and in the possible other inferior border identification step 16, respectively, can be used as the boundary of the extension of the closing operation toward the outside of the lung, that is, the points that the closing operation could include in the pulmonary segment R6 but that are outside the aforementioned borders are not included in the pulmonary segment R6 (figs. 8b- c).

A subsequent second threshold segmentation step 18 can provide to segment the lung segment R6, through the application of suitable densitometric thresholds, into several regions, such as for example a new low-density region R7, a new high- density region R8 and a very high-density region R9 (fig. 8d).

For example, the new low density region R7 could include the voxels of the pulmonary segment R6 which are associated with density values comprised between about - 950 HU and about - 750 HU; the new high density region R8 could include the voxels of the pulmonary segment R6 associated with density values comprised between about - 750 HU and about - 250 HU; the very high density region R9 could include the voxels of the pulmonary segment R6 associated with density values greater than about - 250 HU.

When the aforementioned densitometric thresholds are chosen appropriately, the new low-density region R7 corresponds to healthy lung tissue; the new high- density region R8 corresponds to ground glass; the very high-density region R9 corresponds to consolidated tissue, but can also include pulmonary blood vessels, since the latter have a density similar to that of consolidated tissue.

Some embodiments provide a vessel removal step 19 (figs. lOa-b), aimed at removing pulmonary blood vessels from the pulmonary segment R6, at least from the very high density region R9. The vessel removal step 19 can first provide to identify any portions of the pulmonary segment R6 having tubular morphology. In order to identify such tubular morphology portions, a vesselness filter can be applied to the pulmonary segment R6 (Piccinelli, M., et al., “A framework for geometric analysis of vascular structures: application to cerebral aneurysms”, IEEE transactions on medical imaging (2009), 28 (8), 1141 - 1155).

The vesselness filter assigns to each voxel of the volume to which it is applied, for example to each voxel of the pulmonary segment R6, a value that quantifies the tubular morphology of that voxel. For example, the portions of the volumetric scan section shown in fig. 9a that exhibit high tubular morphology based on the application of the vesselness filter are shown in fig. 9b in light gray tones or in white. This allows to identify portions of the thoracic volumetric scan with high tubular morphology, corresponding, within the lung volume, to presumed blood vessels. In some embodiments, with each voxel of each presumed blood vessel thus identified there is associated a distance of that voxel from the mediastinum and a diameter of the presumed blood vessel to which the voxel belongs. Each voxel can then be classified as belonging or not belonging to a blood vessel based on these parameters. In fact, the diameter of the blood vessels decreases as the distance from the mediastinum increases (and as the distance from the rib cage decreases). The vessel removal step 19 can then provide to identify as blood vessels the presumed blood vessels with a large diameter and close to the mediastinum. In contrast, the presumed blood vessels with a large diameter but close to the rib cage can be identified as consolidated tissue. Once the blood vessels have been identified, the vessel removal step 19 can provide to remove the blood vessels identified by the pulmonary segment R6, in particular by the very high density region R9 (fig. 10a), thus obtaining a new very high density region RIO (fig. 10b), corresponding to only the consolidated tissue.

It is clear that modifications and/or additions of steps may be made to the method 10 as described heretofore, without thereby departing from the field and scope of the present invention, as defined by the claims.

It is also clear that, although the present invention has been described with reference to some specific examples, a person of skill in the art will be able to achieve other equivalent forms of method for pulmonary segmentation from a thoracic volumetric scan, having the characteristics as set forth in the claims and hence all coming within the field of protection defined thereby.

In the following claims, the sole purpose of the references in brackets is to facilitate their reading and they must not be considered as restrictive factors with regard to the field of protection defined by the claims.

Claims

1. Computer- implemented method (10) for pulmonary segmentation from a thoracic volumetric scan, characterized in that it comprises:

- a first segmentation step (12), which provides to perform a densitometric threshold segmentation of the left and right lung, thus dividing the lung volume of said scan into at least one low-density region (R3) and one high-density region (R4), which are delimited overall by a preliminary pulmonary border (Bl);

- a radiomic analysis step (13), following or preceding said first segmentation step (12), which provides to obtain, for each pixel of each axial section of said thoracic volumetric scan, the value of at least one radiomic feature, thus obtaining, for each of said axial sections, at least one radiomic feature image;

- an identification step (14), which provides to identify a lateral and posterior pulmonary border (B2) on the basis of a combination of information deriving from: a geometric modeling of the shape of the lung parenchyma, said preliminary pulmonary border (Bl), a threshold densitometric segmentation of the rib cage, and said radiomic analysis step (13);

- another identification step (16) of the inferior border aimed at identifying, for each axial section of said thoracic volumetric scan, an inferior pulmonary border (B3) on the basis of a combination of information deriving from: the anatomy; said preliminary pulmonary border (Bl); a densitometric threshold segmentation of tissues underlying the diaphragm; and said radiomic analysis step (13);

- a lung volume reconstruction step (17), which provides to combine at least the results obtained in said previous steps (12, 13, 14, 16) in order to identify a border of the lung along the mediastinal profile, simultaneously reconstructing the entire volume of the lung parenchyma, thus obtaining a pulmonary segment (R6);

- a second segmentation step (18), which provides to segment said pulmonary segment (R6), through the application of at least two densitometric thresholds, into at least one new low-density region (R7), one new high-density region (R8) and one very high density region (R9); - a vessel removal step (19), which provides to remove any pulmonary blood vessels at least from said very high density region (R9), thus obtaining a new very high density region (RIO).

2. Method (10) as in claim 1, characterized in that said new low-density region (R7) corresponds to healthy lung tissue, said new high-density region (R8) corresponds to the ground glass and said new very high-density region (RIO) corresponds to consolidated tissue.

3. Method (10) as in any claim hereinbefore, characterized in that said radiomic analysis step (13) provides to calculate the radiomic features of skewness, entropy and standard deviation.

4. Method (10) as in claim 3, characterized in that said radiomic analysis step (13) comprises normalizing each radiomic feature image by dividing it by its maximum value, filtering each normalized radiomic feature image with a threshold value and finally combining said radiomic feature images by intersection, thus obtaining a synthetic radiomic image.

5. Method (10) as in any claim hereinbefore, characterized in that said lateral and posterior border identification step (14) provides, separately for each right and left lung, to: identify first points (P4) of an internal border of the ribs; interpolate these first points (P4) with a parametric three-dimensional cylindroid spline function, obtaining a first border surface (SI); add to the set of points to be interpolated (P4) second points (P5) of said preliminary pulmonary border (Bl) which are distant from said first border surface (SI) by less than a certain threshold value; interpolate said first and second points (P4, P5) with a parametric three- dimensional cylindroid spline function, obtaining a second border surface (S2); add to the set of points to be interpolated (P4, P5) third points (P6) of said at least one radiomic feature image, or possibly of said synthetic radiomic image, which are distant from said second border surface (S2) by less than a certain threshold value; interpolate said first, second and third points (P4, P5, P6) with a parametric three-dimensional cylindroid spline function, obtaining a third border surface (S3).

6. Method (10) as in any claim hereinbefore, characterized in that said other inferior border identification step (16) operates on axial sections of said thoracic volumetric scan and provides to:

- identify by interpolation of points of said at least one radiomic feature image, or possibly of said synthetic radiomic image, a right part of said inferior pulmonary border (B3) as the union of a plurality of curves (Cl ’, C2’) of the right border of the diaphragm corresponding to one border of the liver;

- identify, through a region growing technique and taking into account the discontinuities in density due to the perisplenic and perigastric adipose tissue and of said preliminary pulmonary border (Bl) and said lateral and posterior pulmonary border (B2), a left part of said inferior pulmonary border (B3) as the union of a plurality of curves (C3, C4) of the left border of the diaphragm corresponding to a border of the spleen and of the perigastric adipose tissue.

7. Method (10) as in any claim hereinbefore, characterized in that said lung volume reconstruction step (17) provides to: join said low and high density regions (R3, R4) into a single region (R5); subject said single region (R5) to a morphological closing operation with a spherical structuring element, using said lateral and posterior (B2) and inferior (B3) pulmonary borders previously identified as the boundary of the extension of said closing operation toward the exterior of the lung, thus obtaining said pulmonary segment (R6).

8. Method (10) as in any claim hereinbefore, characterized in that said vessel removal step (19) provides to: apply a vesselness filter to said pulmonary segment (R6), thus identifying any portions of said pulmonary segment (R6) having high tubular morphology, corresponding to presumed blood vessels; associate with each voxel of each presumed blood vessel a distance of said voxel from the mediastinum and a diameter of said presumed blood vessel; identify as blood vessels the presumed large diameter vessels close to the mediastinum; remove said blood vessels thus identified from said pulmonary segment (R6), in particular from said very high density region (R9), obtaining said new very high density region (RIO).

9. Method (10) as in any claim hereinbefore, characterized in that it comprises a step (11) of delimiting a region of interest (Rl), prior to said other steps (12, 13, 14, 16, 17, 18, 19), which provides to divide said thoracic volumetric scan to be segmented into a region of interest (Rl), the voxels of which are analyzed in the subsequent steps (12, 13, 14, 16, 17, 18, 19) of the method (10), and an external region (R2), which is not analyzed in the subsequent steps (12, 13, 14, 16, 17, 18, 19) of the method (10).

10. Computer program containing instructions which, if executed by a computer, cause the execution of said method (10) as in any claim hereinbefore.