WO2010046573A1 - Method for estimating the concentration of a tracer in a tissue structure assembly, and corresponding storage medium and device - Google Patents

Abstract

The invention relates to a method for estimating the concentration of a tracer in a tissue structure assembly including at least one tissue structure, from a measurement image of the tracer concentration in said tissue structure assembly, which is obtained by an imaging apparatus, wherein said image includes at least one space domain inside which the tracer concentration is homogenous and at least one region of interest in which the tracer concentration is measured. The method includes the following steps: determining a geometric transfer matrix having coefficients representative of the contribution of the space domains in the measurement of the tracer concentration in the regions of interest; optimising the coefficient of the geometric transfer matrix by defining the best regions of interest in terms of errors in order to measure the tracer concentration, the definition of the regions of interest being carried out according to an iterative loop that includes the following steps upon each iteration: modifying the regions of interest, and calculating the coefficients of the geometric transfer matrix from the modified regions of interest; selecting an optimised geometric transfer matrix among the calculated geometric transfer matrices; and estimating the tracer concentration from the optimised geometric transfer matrix.

Description

 Method for estimating the concentration of a tracer in a set of tissue structures, corresponding information medium and device

The present invention relates to a method for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the tracer concentration in said set of tissue structures obtained by a imaging apparatus, the image comprising at least one spatial domain within which the tracer concentration is homogeneous and at least one region of interest within which the tracer concentration is measured. In vivo molecular imaging allows non-invasive local biochemical and pharmacological parameters to be measured in intact animals and humans. Molecular imaging encompasses nuclear techniques, magnetic resonance imaging (MRI) and optical techniques. Among the nuclear techniques, positron emission tomography (PET) is the most sensitive technique and the only one to offer quantitative measurements.

Quantifying measurements with PET for the study of a drug requires the intravenous injection of a radioactive tracer suitable for the study of this drug. An image of the tracer distribution is then reconstructed. Time Activity Curve (TAC), that is, changes in tracer concentration within a volume element, are also strongly spatially correlated not only because of reconstruction but also the existence of physiological regions that respond to the tracer identically, regions that can be called pharmaco-organs.

Post-treatment of the PET image generally requires three steps to access the pharmacological parameter. First, regions of interest (ROI for Region Of Interest) delimiting the organs, more precisely the pharmacokinetic organs, must be defined either on the PET image or on a high-resolution image point-to-point correspondence with the PET image. Then, the TACs of the pharmaco-organs must be extracted and possibly corrected to compensate for the limited spatial resolution of PET. Finally, a physiological model can be defined, based on the TACs and the concentration of the tracer in the plasma, allowing the calculation of the pharmacological parameters of interest to the biologist. Accurate delineation of ROIs is required to extract relevant TACs.

However, the quantification of TACs is hampered by the limited resolution of the PET system, resulting in the so-called Partial Volume Effect (PVE).

A geometrical transfer matrix (GTM) method is known for efficiently correcting the EVP. This method requires the PET image data on one hand, and on the other hand spatial domains delineating the functional organs, regions of interest in which the average TACs of the functional organs will be calculated and a model of resolution of the PET imager.

However, this GTM process is very sensitive not only to delineation errors of said spatial domains and image reconstruction artifacts, but also to image smoothing effects due to periodic physiological motions, such as heartbeat. or breathing movements.

The aim of the invention is to propose a method for limiting the impact of effects due to segmentation errors, image reconstruction artifacts and physiological movements on the efficiency of the GTM process in order to obtain TACs presenting the least bias with the lowest noise uncertainty.

For this purpose, the object of the invention is a process of the aforementioned type, characterized in that it comprises the following steps: - determination of a geometric transfer matrix whose coefficients are representative of the contribution of the spatial domains to the measurement the tracer concentration in the regions of interest;

optimization of the coefficients of the geometric transfer matrix by delimitation of the regions of best interest in terms of errors to measure the concentration of the tracer, the delimitation of the regions of interest being carried out according to an iterative loop comprising at each iteration the steps following: - modification of regions of interest; and calculating the coefficients of the geometric transfer matrix from the modified regions of interest;

selecting an optimized geometric transfer matrix from the calculated geometric transfer matrices; and - estimating the tracer concentration from the optimized geometric transfer matrix.

The method according to the invention may comprise one or more of the following characteristics:

the step of modifying the regions of interest comprises adding and / or excluding at least one picture element;

the iterative loop for the delimitation of the regions of interest comprises, for each iteration, the following steps:

determining a set of elementary modifications each consisting in adding at least one image element to the regions of interest; - estimating, for each elementary modification, the concentration of the tracer in the regions of interest modified by the elementary modification;

- evaluating, for each elementary modification, a criterion making it possible to compare the estimated tracer concentration with the tracer concentrations estimated at at least one previous iteration; - choice of the elementary modification giving an estimate of the tracer concentration closest to tracer concentrations estimated at previous iterations;

comparing the criterion corresponding to the elementary modification chosen with a predetermined threshold; and - application or not of the chosen elementary modification or stop of the iterative loop as a function of the result of the comparison of the criterion with the threshold;

the method comprises a step of scheduling the image elements in each spatial domain according to order criteria making it possible to evaluate their risk of introducing non-noise errors in the estimation of the tracer concentration;

the order criteria include at least one of a first criterion for defining whether the concentration contained in the picture element originates from one or more different spatial domains, a second criterion defining whether the image element participates in computing the non-diagonal elements of the geometric transfer matrix, and a third criterion for defining whether the pixel is reliable for estimating the tracer concentration; the first criterion makes it possible to measure the homogeneity of the concentration of the tracer in the neighborhood of the image element, the second criterion makes it possible to measure the contribution of the spatial domains to the measurement of the concentration of the tracer in the element of image, and the third criterion comes from a probabilistic atlas; the iterative loop for the delimitation of the regions of interest is carried out according to the order of the pixels obtained in the scheduling step;

the method comprises a step of estimating the spreading function of the point of the imaging apparatus at any point in the field of view of the imaging apparatus; the step of estimating the spreading function of the point comprises the following steps:

measuring the spreading function of the point at several points of the field of view; and

- interpolation and / or extrapolation of the measurements obtained; the method comprises a step of estimating a region spreading function for each spatial domain by convolution of the point spreading function with each spatial domain;

the method comprises a step of calculating the size of the regions of interest;

the method comprises a step of delimiting the spatial domains so as to obtain regions of interest of maximum size;

the method comprises a step of optimizing the order criteria so as to obtain regions of interest of maximum size;

the measurement image of the tracer concentration in the tissue structure is a sequence of three-dimensional images mapped and comprising at least one three-dimensional image, and the image elements are voxels; and

the three-dimensional image is obtained by positron emission tomography. The invention also relates to an information carrier comprising a code for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the concentration of the tracer in said set of tissue structures obtained by an imaging apparatus, the image comprising at least one spatial domain within which the tracer concentration is homogeneous and at least one region of interest within which the tracer concentration is measured, characterized in that the code comprises instructions for:

determining a geometric transfer matrix whose coefficients are representative of the contribution of the spatial domains to the measurement of the tracer concentration in the regions of interest;

optimizing the coefficients of the geometric transfer matrix by delimiting the regions of best interest in terms of errors for measuring the tracer concentration, the delimitation of the regions of interest being carried out according to an iterative loop comprising at each iteration the steps following:

- modification of regions of interest; and

calculating the coefficients of the geometric transfer matrix from the modified regions of interest;

selecting an optimized geometric transfer matrix from the calculated geometric transfer matrices; and

- estimate the tracer concentration from the optimized geometric transfer matrix, when executed by a data processing system. The invention further relates to a device for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the tracer concentration in said set of tissue structures. obtained by an imaging apparatus, the image comprising at least one spatial domain within which the tracer concentration is homogeneous and at least one region of interest within which the tracer concentration is measured, characterized in that it comprises:

- an imaging device; and

a data processing system comprising: means for determining a geometric transfer matrix whose coefficients are representative of the contribution of the spatial domains to the measurement of the concentration of the tracer in the regions of interest;

means for optimizing the coefficients of the geometric transfer matrix by delimiting the regions of best interest in terms of errors for measuring the tracer concentration, the delimitation of the regions of interest being carried out according to an iterative loop comprising at each iteration the following steps:

- modification of regions of interest; and calculating the coefficients of the geometric transfer matrix from the modified regions of interest;

means for selecting an optimized geometric transfer matrix from the calculated geometric transfer matrices; and

means for estimating the tracer concentration from the optimized geometric transfer matrix.

The invention will be better understood on reading the description which follows, given solely by way of example and with reference to the appended drawings, in which:

- Figure 1 is a flow chart showing the two main phases of the method according to the invention;

FIG. 2 is a diagram illustrating the PET image of a mouse;

FIG. 3 is a flowchart representing the steps of the method according to the invention;

FIGS. 4 and 5 are respectively transaxial and longitudinal sections of the field of view of the PET scanner;

Figure 6 is a flowchart showing in more detail the step of delimiting the regions of interest; and

FIGS. 7A to 11B are diagrams showing the results obtained with an exemplary embodiment of the method according to the invention. The present invention is based on the geometric transfer matrix (GTM) method of correction of the partial volume effect (PVE) proposed by Rousset et al. in the article "Correction for Partial Volume Effects in PET: Principle and Validation," The Journal of Nuclear Medicine, 1998, Vol. 39, No. 5, pages 904 to 911, and implemented in the image space by Frouin et al. in the article "Correction of Partial-Volume Effect for Striatal Imaging PET: Fast Implementation and Study of Robustness", The Journal of Nuclear Medicine, 2002, Vol. 43, No. 12, pages 1715 to 1726. The purpose of the GTM process is to find the true tracer concentration T ₁ within spatial domains D ₁ . Details of the GTM process can be found in the article by Rousset et al. above, but the principles are presented here in order to highlight the features of the present invention. The method according to the invention aims to optimize the calculation parameters of the true concentration T, the tracer by limiting the effects due to segmentation errors, image reconstruction artifacts and physiological movements.

As shown in FIG. 1, the method according to the invention comprises two main phases: a first phase intended to choose the response function of the PET scanner in an adequate manner; and

a second phase 12 intended to select the optimal regions of interest for estimating the true concentration T ₁ of the tracer.

With reference to FIG. 2, the PET image is a three-dimensional image formed of a set of voxels X (x, y, z) each characterized by a concentration of tracer A (x, y, z), this concentration being susceptible to vary over time.

The PET image is assumed to be composed of I homogeneous spatial domains λPih _≤ i _≤ i in terms of tracer concentration. Each space domain D ₁ has a true tracer concentration denoted T ₁ .

A concentration of tracer if / j _{1 <; KJ} is measured in the PET image at

within regions of interest (ROI) IR ₇ | _1≤7≤J .

For reasons of clarity of the explanation, we will consider in the following that each ROI is included in the corresponding spatial domain and that their number are equal: J = I. Alternatively, J is different from I, preferably with J> I, so that the GTM process makes it possible to recover the true concentration T, of the tracer.

The tracer concentration measured t _j in a region of interest R _j, volume V _j, can be expressed as a linear combination of the concentrations of the tracer V _z Ji _<? < _/ in the different spatial domains:

where RSF, (x) is the contribution of the spatial domain D ₁ to the PET measurement of tracer concentration in voxel x.

This equation can be rewritten in matrix form: t = WT (2) where t is the vector of tracer concentration PET measurements averaged in ROIs R _j , T is the vector of true tracer concentrations in D space domains, and W is the geometric transfer matrix

(GTM) I x J between the spatial domains IA Ji _{≤z≤ /} and the ROIs Wj i _ι≤J≤J . The coefficients w _u of the matrix W are equal to: w; _J = ^ _r \ RSFXx) CbC (3)

V V J R,

The true tracer concentration can then be calculated as: T = W ^"1 .t (4)

As illustrated in FIG. 3, an initial step 14 of phase 12 of the method according to the invention consists in delimiting the spatial domains D ₁ .

The user is free to delimit the D _1s as he wishes, by drawing them manually, semi-automatically or automatically on the PET image, but under the constraint that each spatial domain D is homogeneous in terms of the true concentration of the tracer. this concentration may vary over time in the case where the PET image contains dynamic information.

Alternatively, the Ds are plotted on an image from another imaging modality, such as nuclear NRM or the X-ray scanner.

In parallel, the phase 10 of choice of the response function of the PET scanner is implemented. The scanner response function, or Point Spread Function (PSF) function, is preferably estimated for the GTM process to work effectively.

With reference to FIGS. 4 and 5, the PSF is decomposed in a cylindrical coordinate system (p, θ, z) into:

an axial component along the z axis: _axial PSF;

a radial component in the direction CX: PSF _raC iιaie; and

- a tangential component in the direction XT: PSF _your ngentιeiie-

PSF components _axιa ie, PSF _radιa ι _e and _tan gentιeiie PSF PSF may each comprise several parameters, for example two variances of two Gaussian if the PSF is estimated by two Gaussian, or more values.

Assuming that PSF is invariant according to z, it is necessary to measure the PSF only on the transaxial plane passing through the middle of the PET scanner.

According to this hypothesis, PSF _axιa ι _e does not vary with position in the field of view, and more _radιa ie PSF and PSFT _year gentιeiie are invariant under θ. So just measure PSF _axιa ι _e at a point of view of field and PSF and PSF _radιa ι _{e ta} ngentιeiie in several points of the x axis.

Assuming that the PSF varies as one moves in the direction of the axis of the PET scanner, it is necessary to measure PSF _axιa ι _e at several points along the z axis, and PSF _radιa ι _e and PSF _your ngentιeiie in several points of the x axis.

Referring to Figure 3, a first step 16 of the method 10 of the phase according to the invention to measure PSF _axιa ι _e, PSF _radιa ι _e and TFTP gentιeiie _year at several points the field of view of the PET scanner.

According to the invention, the PSF may have any parametric or non-parametric shape, or more precisely a shape having a good match with the PSF form generated by the PET image reconstruction method, and possibly also with the effect smoothing performed by physiological movements. If the PET image is reconstructed with modeling or PSF compensation, the PSF in the context of this invention will model the residual effects not taken into account during the reconstruction (for example if the reconstruction assumes the uniform PSF in the field of view, then the gap for each point or region of the field of view between the uniform PSF and the actual PSF will be modeled). In a second step 18 of the stage 10, is deduced from these measures the value of PSF _axιa ie, PSF _rac iιaie and PSF _your ngentιeiie anywhere in the field of view by interpolation and / or extrapolation, by a method chosen by the user.

As seen previously, PSF _axi aie depends only on the position on the z axis and will be noted PSF _{axi a} (z), while PSF _raC iιaie and PSF _ta ngénieie only depend on p and will therefore be noted PSF _radιa ι _e (p) and PSF _tan gentïeiie (p) -

In step 20 of the method according to the invention, the RSF is estimated (x) at any point in the field of view.

This estimate is obtained by convolution of the convolution mask corresponding to the PSF with the voxels of the spatial domain D ₁ .

The convolution mask corresponding to the PSF of the PET scanner is obtained as follows:

Let C = (x _c , y _c ) and X = (x, y, z) respectively the coordinates of the scanner axis and those of a voxel located in the field of view of the scanner. Let δX = (δx, δy, δz) a small displacement around X.

The value in X + δX of the convolution mask corresponding to the PSF in X is given by:

PSF _xyz (&, δy, δz) = f (δr, PSF _n ((j)), a, PSF ^ _sen ((j)), z, PSF _aa ( ₍ z)) (5)

Where δr = accosθ + δysinθ, = = - & sinθ + δycosθ,

and f is the function representing the pace of the PSF. RSF, (x) is obtained in the same way as the co nvo lution

PSF _xyz {& c, δy, δz) to the mask of D ₁ . This smoothing is non-homogeneous and reproduces the effects produced by the acquisition and reconstruction of the image.

In step 22, a set Λ of voxels whose activity measured in their neighborhood is representative of the concentration of tracer in the spatial domain D 1 is defined in parallel. These voxels can be determined automatically or manually. We also define S as a set of voxels to exclude regions of interest R _j .

A scheduling of voxels in each spatial domain D is then performed in step 24 according to the following principle. At least one voxel x e D, is affected by effects that smooth images

(partial volume effect, tissue fraction in the voxel, physiological movement, etc.).

It is possible to predict, with or without introduction of information a priori, in what proportion this voxel x is affected by these effects from one or more predetermined order criteria.

Among these criteria:

- criterion α (x): the voxel x is affected by the effect of partial volume, the tissue fraction or the periodic physiological movements of the subject, this criterion making it possible to define in other words if the activity contained in the voxel x comes from one or more different spatial domains;

- criterion β (x): the voxel x participates in the computation of the non-diagonal elements of the matrix W, this criterion being able to be easily calculated from the computation of the RSF _k (x) for k ≠ i;

- criterion γ (x): any other criterion or set of criteria that can help to discriminate the most reliable voxels for the estimation of the tracer concentration (reliability of the PET measurement, homogeneity of response within the same organ, etc.), for example from a probabilistic atlas.

The criteria α (x), β (x) and γ (x) are defined in such a way that the smaller they are, the less the aggregation of x to the region of interest R ₁ will introduce correction error due to Segmentation errors or effects not taken into account in the PSF estimate.

Let g, (α (x), β (x), γ (x), Λ _h S ₁ ) be an order of voxels xe D ₁ giving x a rank all the smaller as α (x), β (x ) and γ (x) are small, under the constraint that:

(C1) A voxel x e D, of rank n is connected by voxels of rank lower than one of the points of Λ,; and

(C2) voxels ye S ₁ are excluded from this sort order. In the following, R, ^{(w) will} denote an ROI containing the n voxels of ranks g, (α (x), β (x), γ (x)) the smallest among all the voxels of D ₁ or a part of these voxels.

R, ^(w) may also contain voxels located outside D, but at its periphery.

The automatic optimization of the regions of interest R _j is then performed in step 26 of the method according to the invention in an iterative loop.

Firstly, the parameters useful for the automatic optimization of the regions of interest R _{j are} defined. Let Rj <= Ri, or T 'and T be the concentrations of the tracer estimated in the members respectively with R, and with R ₁ . Let T be the true concentration of the tracer. We can say that:

- Errors due to segmentation errors on the estimation of T using R are less important than those made using R ₁ . In addition, the geometric transfer matrix W estimated using R 1 is closer to the unit matrix than the geometric transfer matrix W estimated using R 1. The conditioning of the matrix W is better than that of W.

- Noise errors made on estimating T using R, are larger than those made using R,. - It is possible to estimate T _{lk τ} , which is the concentration of the tracer in the region of interest R ₁ for the iteration k and the measurement time τ:

or by optimizing R ₁ globally by integrating the information of all the measurement times (referred to as "frames" in the following);

or by optimizing R ₁ separately for each frame τ; or by optimizing R ₁ separately for each frame τ _op t while also integrating information from the other frames τ _mes - The last case corresponds to a general case, the first two cases being special cases. We will therefore consider in the following only the general case. n ( ⁿ ι k)

We define beforehand ^κ ι, s _h ^' _k as the ROI containing at iteration k the n, _ιk voxels having the voxels x of the lowest rank among all the voxels of D ₁ excluding the voxels belonging to the together S ,, k-

Let δ be an elementary modification of the number of voxels of the regions of interest R ₁ and let Δ _{k be} a set of possible elementary modifications to the iteration k.

Let ' ^{tk τ} ° ^{pt be} the ROI defined at the iteration k for the optimization of the region n of interest R ₁ for the frame τ _op t, containing ^α ^' voxels and excluding the voxels of

the set ^ι ' ^k ' ^τ ° ^pι while remaining consistent with (C1) and (C2).

T ^! ' ^k ' ^{T <} "" the concentration of the tracer that is to be estimated for the

frame τ _opt and let ^α -v. ^ the tracer concentration estimated in ^{ïΛ t} ^ 'at the time τ _m es

Let ^σ x, _{Tmβs be} the variance of the noise at voxel x at time τ _me s- This variance can be determined as a function of the image reconstruction algorithm and possibly of the Voxel x TEP signal at time τ _me s-

Let ^ i, _v , _τ _ be the average tracer concentration measured in the region

of interest R ₁ at the iteration k at time τ _mes for the estimation of

-

Under certain assumptions, the variance ζ ^ _k , _τop ,, _τm of the error due to the noise made on the estimate of

P ^was also determined from

.2

<r _XlTmB and of "_{α> v} .

Knowing this error, the GTM method provides an upper bound of the variance ζl _{k, τo l, τ ^} of the error due to the noise made on the estimate of ^ ι, k, τ,:

where ^ k, τ _opι is the GTM matrix calculated at the iteration k for the estimation of

} \ _ι≤ι≤I -

Is

a quantitative criterion evaluating the difference in noise between the estimate of te. *. w _" L _{rm ≤r,} i _{≤! </} ^at the iteration k and the estimates of fc /, v, r _m J 1 _{un / </ tl≤r} _ _{≤ri≤; ≤ /} at previous iterations. h will be smaller if ψ ,, k, τ _opt , τ _my [ _{≤τ ≤T ι≤i≤I} is

different from fc, _v , r_

-

Let _{hmn be} a threshold below which Vα, _v , ^ j _{1 <7 <r i <} , _<7 is judged

significantly different from Ψi, ι, τ _opt , τ _mes \ _≤ ι _{<kA≤Tmes ≤TA≤i≤I} -

The algorithm of the iterative loop for the optimization of the regions of interest R _j is as follows (FIG. 6). Step 28 corresponds to the initialization of the loop.

The regions of interest R _j are initialized by ^R ι, s ° ^{* °} ", where"_{Î;O; V} > 2 and where S ,, o is an empty set.

At each iteration k, then proceed to a step 30 of determining a set Δ _k of possible elementary modifications δ possible of ψιj ^* _k ', ^ | _<;</ which can result in the addition of voxels to a region of interest R ₁ and / or in the addition of voxels to S, and / or in the modification of A ₁ . At step 32, for ut δ e Δ, fc "ll≤i≤IW is calculated,

> ^τ opt, ^τ mes

<fl ≤ i ≤ I and

In step 36, the elementary modification δ _op t is chosen which gives an estimate of the tracer concentration closest to the tracer concentrations estimated at the previous iterations:

In step 38, the criterion h corresponding to δ _opt is compared with h _r

If, then we apply the

elementary modification δ _opt at R ₁ and S ₁ (step 40) and the next iteration k + 1 is set back to step 30.

Otherwise the algorithm ends (step 42), and we have:)

Opt

- I ⁷ ^ V -L≤ / rU- ⁱ .v> v {≤, ≤ / ^{; and}

- I Pf 'Op, \ f _{≤ ι≤I} = k P' ² * 1 ^ ^ A f≤i≤ / ^'

With reference to FIG. 3, the estimation of the true tracer concentration T ₁ in the spatial domains D ₁ is thus obtained in step 44.

From the optimization step 26 of the regions of interest R _j , the final size of these regions of interest (step 46) is also obtained, which makes it possible to optimize the parameters of the method.

Indeed, the method comprises various parameters determining the order in which the voxels will be aggregated to the ROIs. Among these parameters, we find the A ₁ , the functions α (x), β (x), γ (x) and g, (α (x), β (x), γ (x), A ₁ , S ₁ ). The optimization process makes a compromise between noise estimation errors and estimation errors arising from effects such as segmentation errors and physiological movements.

Noise errors tend to decrease as the size of the ROIs increases.

The final size of the R _j is therefore an indicator of:

- the quality of D ₁ ;

the homogeneity of the PET signal within the D ₁ ; and

the quality of the parameters determining the order of the voxels. AD, and homogeneity of the fixed signal, the final size of the R _j is therefore an indicator of the quality of the choice of A ₁ , the functions α (x), β (x), γ (x) and g, (α (x ), β (x), γ (x), A ₁ , S ₁ ).

We can either:

to launch several optimizations of R _j with D ₁ , A ₁ , α (x), β (x), γ (x) and gι (oc (x), β (x), γ (x), A ₁ , S ₁ ) and choose the D ₁ , A ₁ , α (x), β (x), γ (x) and g, (α (x), β (x), γ (x), A ₁ , S ₁ ) for which the size of the R _j , or any increasing monotonic function of this size or any other function quantitatively evaluating the quality of the estimation of T ₁ , is maximal (step 48);

- give the size of the R _j , or the value of any monotonic function of this size or any other function that quantitatively evaluates the quality of the T ₁ estimate, as an overall indicator of the quality of the D ₁ , A ₁ , α (x), β (x), γ (x) and g, (α (x), β (x), γ (x), A ₁ , S ₁ ) (step 50).

An embodiment of the method according to the invention is as follows:

- RSF | (x) is estimated according to the method of Frouin et al. ; - J = 1 and j = i;

the Ds are defined as the spatial domains delimited by the segmentation using a method based on the local mean analysis (LMA) proposed by Maroy et al. in the article "Segmentation of Rodent Whole-Body Dynamic PET Images: An Unsupervised Method Based on Voxel Dynamics," IEEE Trans. Med. Imaging, 2008; - Λ, is the set of voxels xe D ₁ extracted during the extraction step of local minima of the image _p ² measuring each voxel p the homogeneity of the concentration of the tracer in the neighborhood of p;

α (x) is the order of voxels sorted by increasing _x ² ;

- β (x) is the order of voxels sorted by RSF ₁ (X) decreasing;

- g, (α (x), β (x), γ (x), A ₁ , S ₁ ) is the order of x among voxels xe D ₁ whose voxels of S ₁ are excluded and which are sorted by (axα ² (x) + (1-a) xβ ² (x)) increasing. In this order, the voxels not checking (C1) see their rank shifted until they can check (C1);

- A _k

, where S ₁ is the elementary modification which consists in adding d, voxels to the region of interest R,;

- τ _op t is the frame for which we optimize the {^} _{1≤; ≤7} and τ _mes is a frame where we estimate the {^} i _{≤ * ≤ /} ;

- Let L be a scale level and L _max the maximum scale level, fixed a priori (for example L _max = T / 4). We associate with the frame x _op t at the scale level 0 a set F _L of L + 1 consecutive frames at the level of scale L, these L frames being centered on x _op t; let p _{l | L} the p-value associated with the t statistic

The algorithm is the following.

At the initialization stage: δ, = 0.001 x # D ,, where #D, is the cardinal of the segmented spatial domain D ₁ , ⁿ i, o, τ _opt = δ _h R ₁ is initialized by R ^ '^ and S ,, _o is an empty set. At each iteration k and for all i: we calculate fc "L _,, W, te., _ L_ _ffJM , ^ ZΛ _{^ m} , and

If h _mm , then we add d.

, (», J voxels at ^w _; £ r, redefine ^ s ^' / ^° "' and proceed to the next iteration k + 1.

Otherwise, if d,> 1, we decrease d, and we retry an iteration. Otherwise, if #S, <s _ma χ, add the above d, voxels to S ₁ , increase d, and retry an iteration.

Otherwise, we stop the algorithm and the solution is:

- R, = OIl ':

- K ^, j; _{≤ ≤ ≤ /} Kt-i, _Wop , 1 1 _{î ≤ /} : and

^'

The results were obtained with the following material. For the simulation, we simulated through an analytical simulator fifty PET acquisitions of the MOBY mouse phantom on a PET scanner dedicated to the small animal (the simulated scanner was the Siemens FOCUS 220 with a spatial resolution of 1.3mm), with time frames of one minute each. The images were reconstructed using the Ordered Subsets Expectation Maximization (OSEM) method, which is a conventional PET image reconstruction method, with voxels of 0.5x0.5x0.8 mm ³ .

Experimental data include sixteen mice injected with xenografts of peritoneal human tumors. The mice received an injection of Fluoro-thymidine (18F-FLT), immediately followed by a PET scan of 5400 min consisting of 18 time frames (five one-minute frames, five two-minute frames, three five-minute frames, three ten-minute frames, and two fifteen-minute frames). The images were reconstructed using the OSEM method with voxels of 0.5x0.5x0.8 mm ³ . Mice were sacrificed immediately after acquisition, and tracer concentration was counted in the organs removed. This concentration counted in the organs served as a "gold standard" for the concentration estimated at the end of the acquisition.

The PET images were segmented using the LMA method (Figure 7A for simulation and Figure 7B for experimental data). The estimation of the tracer concentration in the organs is carried out according to three methods:

calculating the average TAC within the segmented spatial domain corresponding to the organ;

partial volume correction performed by the GTM method implemented in the image space as described by Frouin et al. (named simply "GTM" in the following); and

partial volume correction performed by the method proposed in the embodiment of the present invention (hereinafter referred to as "LMA-GTM"). The GTM and LMA-GTM processes were compared on the basis of:

- ETR: which is the difference in absolute value between the corresponding measurement of the gold standard and the measured value, expressed as a percentage of the gold standard measurement; and

- ARC: which is the apparent recovery ratio, in other words the measured value divided by the corresponding value of the gold standard.

Figures 8A and 8B show examples of temporal variation in tracer concentration (TACs) estimated for simulations (Figure 8A) and for experimental data (Figure 8B).

Figures 9A and 9B show that the LMA-GTM method significantly improves the measurement accuracy (the p-value is less than 10 ^{~ 5} for both the simulations and the experimental data) and that the process

GTM improves only TACs extracted for experimental data only. The precision of the measurements estimated using the LMA-GTM method is good (ETR _Slm ulatons = 5.3% ± 8% and ETR _D omex experimental = 4.8% ± 7%) and the apparent contrast is correctly recovered (ARC _Slm s uiatιo _n = 94% ± 10%

Experimental ARC ⁼ 99.98% ± 9%). The results obtained using the LMA-GTM method are significantly better in terms of ARC and ETR than those obtained by the GTM method (the p-value is less than 10 ^{~ 5} ).

The complete TAC extraction process corrected for partial volume takes 10 minutes per image: - segmentation (-15 seconds);

- naming of the organs by the user (~ 5 minutes); and

- partial volume correction (~ 4 minutes).

Figure 10A shows the ARC obtained by the three methods on the simulation images for each organ. For the LMA-GTM method, the ARC varies between 94% ± 5% and 99.3% ± 1% for all organs except for small organs, such as the thalamus and thyroid, and the pancreas because of its form.

ETR (Figure 10B) is less than 10% of the true value for all organs except the thalamus (ETR = 32% ± 5%) and the thyroid (ETR = 21% ± 6%).

The accuracy obtained using the LMA-GTM method is better than that obtained by calculating the average TAC and by the GTM method, and with a p-value less than 10 ^{~ 5} for all the organs.

The results obtained for the experimental data are similar to those obtained for the simulations.

The contrast recovery (Figure 11A) was good (between 98.4% ± 7% and 110% ± 6%).

In addition, the ETR (Figure 11B) was less than 10% for the LMA-GTM process for all organs. In both cases, the LMA-GTM process achieves better results than the GTM process (p-value less than 0.002). For brain studies, GTM is considered as one or even the reference method for correcting the partial volume effect in PET images.

For rodent studies, there is no known generic procedure (applicable to all tracers and pathologies) for partial volume correction.

The LMA-GTM process shows great potential for such studies, but also probably for studies in humans, in whole body or brain imaging. Indeed, the LMA-GTM process:

- uses the automatic LMA delineation of organs in PET images;

- does not require anatomical information or a priori;

- is robust to segmentation errors and errors due to non-modelable effects by PVE correction; and - is fast, easy to use and accurate.

The invention therefore proposes a method for reliably and accurately estimating the concentration of a tracer in a set of tissue structures by optimally defining the regions of interest in which the tracer concentration is measured.

Claims

1. A method for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the tracer concentration in said set of tissue structures obtained by a device imaging, the image comprising at least one spatial domain (D ₁ ) within which the tracer concentration is homogeneous and at least one region of interest (R _j ) within which the tracer concentration is measured, characterized in that it comprises the following steps: - determination of a geometric transfer matrix (W) whose coefficients (w _u ) are representative of the contribution of the spatial domains (D,) to the measurement of the concentration tracer in the regions of interest (R _j );

optimization of the coefficients (w _u ) of the geometric transfer matrix (W) by delimitation of the regions of interest (R _j ) the best in terms of errors to measure the tracer concentration, the delimitation of the regions of interest ( R _j ) being carried out according to an iterative loop comprising at each iteration (k) the following steps:

- modification of the regions of interest (R _j ); and

calculating the coefficients (w _u ) of the geometric transfer matrix (W) from the modified regions of interest (R _j );

selecting an optimized geometric transfer matrix (W) from the calculated geometric transfer matrices (W); and

- Estimation of the tracer concentration from the optimized geometric transfer matrix (W).

2. Method according to claim 1, characterized in that the step of modifying the regions of interest (R _j ) comprises the addition and / or exclusion of at least one image element (x).

3.- Method according to claim 1 or 2, characterized in that the iterative loop for delimiting the regions of interest (R _j ) comprises, for each iteration (k), the following steps:

determining a set of elementary modifications (δ) each consisting in adding at least one image element (x) to the regions of interest (R _j ); estimating, for each elementary modification (δ), the concentration of the tracer in the regions of interest (R _j ) modified by the elementary modification (δ);

- evaluating, for each elementary modification (δ), a criterion (h) making it possible to compare the estimated tracer concentration with the tracer concentrations estimated at at least one previous iteration;

- choice of the elementary modification (δ _op t) giving an estimate of the tracer concentration closest to the tracer concentrations estimated at the previous iterations; - comparison of the criterion (h) corresponding to the selected elementary modification (δ _opt ) with a predetermined threshold (h _mn ); and

- application or not of the selected elementary modification (δ _opt ) or stop of the iterative loop as a function of the result of the comparison of the criterion (h) with the threshold (h _mm ).

4. Method according to any one of claims 1 to 3, characterized in that it comprises a step of ordering the image elements (x) in each spatial domain (D ₁ ) according to order criteria to assess their risk of introducing non-noise errors into tracer concentration estimation.

5. A method according to claim 4, characterized in that the order criteria comprise at least one of a first criterion (α (x)) for defining whether the concentration contained in the picture element (x ) comes from one or more different spatial domains (D ₁ ), a second criterion (β (x)) to define whether the image element (x) participates in the calculation of the non-diagonal elements of the matrix of geometric transfer (W), and a third criterion (γ (x)) for defining whether the pixel (x) is reliable for estimating the tracer concentration.

6. A method according to claim 5, characterized in that the first criterion (α (x)) makes it possible to measure the homogeneity of the concentration of the tracer in the vicinity of the image element (x), the second criterion (β (x)) makes it possible to measure the contribution of the spatial domains (D,) to the measurement of the concentration of the tracer in the image element (x), and the third criterion (γ (x)) comes from a probabilistic atlas.

7.- Method according to any one of claims 4 to 6, characterized in that the iterative loop for the delimitation of the regions of interest (R _j ) is performed according to the order of the image elements (x) obtained u at the scheduling stage.

8. A method according to any one of claims 1 to 7, characterized in that it comprises a step of estimating the spreading point function (PSF) of the imaging apparatus at any point in the field from the view of the imaging apparatus.

9. A method according to claim 8, characterized in that the step of estimating the spreading point function (PSF) comprises the following steps:

- measurement of the point spread function (PSF) at several points of the field of view; and

- interpolation and / or extrapolation of the measurements obtained.

10. A method according to claim 8 or 9, characterized in that it comprises a step of estimating a region spreading function (RSF) for each spatial domain (D ₁ ) by convolution of the function of spreading point (PSF) with each spatial domain (D ₁ ).

11. A method according to any one of claims 1 to 10, characterized in that it comprises a step of calculating the size of the regions of interest (R _j ).

12. A method according to claim 11, characterized in that it comprises a step of delimiting the spatial domains (D,) so as to obtain regions of interest (R _j ) maximum size.

13.- Method according to claim 11 or 12 and any one of claims 4 to 7 taken together, characterized in that it comprises a step of optimizing the order criteria so as to obtain regions of interest ( R _j ) maximum size.

14. A method according to any one of claims 1 to 13, characterized in that the measurement image of the tracer concentration in the tissue structure is a sequence of three-dimensional images mapped and comprising at least one three-dimensional image, and in that and the image elements (x) are voxels.

15.- Method according to claim 14, characterized in that the three-dimensional image is obtained by tomography by emission of positrons.

16. An information medium comprising a code for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the tracer concentration in said set of tissue structures. obtained by an imaging apparatus, the image comprising at least one spatial domain (D,) within which the concentration of the tracer is homogeneous and at least one region of interest (R _j ) within which the concentration of the tracer is measured, characterized in that the code comprises instructions for:

- determining a geometric transfer matrix (W) whose coefficients (W | _j ) are representative of the contribution of the spatial domains (D,) to the measurement of the tracer concentration in the regions of interest (R _j );

optimizing the coefficients (w _u ) of the geometric transfer matrix (W) by delimitation of the regions of interest (R _j ) the best in terms of errors for measuring the tracer concentration, the delimitation of the regions of interest ( R _j ) being carried out according to an iterative loop comprising at each iteration (k) the following steps:

- modification of the regions of interest (R _j ); and

calculating the coefficients (w _u ) of the geometric transfer matrix (W) from the modified regions of interest (R _j );

selecting an optimized geometric transfer matrix (W) from the calculated geometric transfer matrices (W); and

estimating the tracer concentration from the optimized geometric transfer matrix (W) when it is executed by a data processing system.

17.- A device for estimating the concentration of a tracer in a set of tissue structures comprising at least one tissue structure from an image for measuring the tracer concentration in said set of tissue structures obtained by an apparatus of imaging, the image comprising at least one spatial domain (D,) within which the tracer concentration is homogeneous and at least one region of interest (R _j ) within which the concentration of the tracer is measured, characterized in that it comprises:

- an imaging device; and

a data processing system comprising: means for determining a geometric transfer matrix

(W) whose coefficients (w _u ) are representative of the contribution of spatial domains (D ₁ ) to the measurement of tracer concentration in regions of interest

means for optimizing the coefficients (w _u ) of the geometric transfer matrix (W) by delimitation of the regions of interest (R _j ) the best in terms of errors for measuring the tracer concentration, the delimitation of the regions of interest (R _j ) being carried out according to an iterative loop comprising at each iteration (k) the following steps:

- modification of the regions of interest (R _j ); and - calculating the coefficients (w ,,) of the geometric transfer matrix (W) from the modified regions of interest (R _j );

means for selecting an optimized geometric transfer matrix (W) from the calculated geometric transfer matrices (W); and

means for estimating the tracer concentration from the optimized geometric transfer matrix (W).