TAVITIAN, Bertrand (186 rue de Vaugirard, Paris, Paris, F-75015, FR)
MAROY, Renaud (22 rue de Chambéry, Paris, F-75015, FR)
TAVITIAN, Bertrand (186 rue de Vaugirard, Paris, Paris, F-75015, FR)
REVENDICATIONS
1.- Procédé d'estimation de la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial (D 1 ) à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt (R j ) à l'intérieur de laquelle la concentration du traceur est mesurée, caractérisé en ce qu'il comprend les étapes suivantes : - détermination d'une matrice de transfert géométrique (W) dont les coefficients (w u ) sont représentatifs de la contribution des domaines spatiaux (D,) à la mesure de la concentration du traceur dans les régions d'intérêt (R j ) ;
- optimisation des coefficients (w u ) de la matrice de transfert géométrique (W) par délimitation des régions d'intérêt (R j ) les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt (R j ) étant réalisée suivant une boucle itérative comprenant à chaque itération (k) les étapes suivantes :
- modification des régions d'intérêt (R j ) ; et
- calcul des coefficients (w u ) de la matrice de transfert géométrique (W) à partir des régions d'intérêt (R j ) modifiées ;
- choix d'une matrice de transfert géométrique (W) optimisée parmi les matrices de transfert géométrique (W) calculées ; et
- estimation de la concentration du traceur à partir de la matrice de transfert géométrique (W) optimisée. 2.- Procédé selon la revendication 1 , caractérisé en ce que l'étape de modification des régions d'intérêt (R j ) comprend l'ajout et/ou l'exclusion d'au moins un élément d'image (x).
3.- Procédé selon la revendication 1 ou 2, caractérisé en ce que la boucle itérative pour la délimitation des régions d'intérêt (R j ) comprend, pour chaque itération (k), les étapes suivantes :
- détermination d'un ensemble de modifications élémentaires (δ) consistant chacune à ajouter au moins un élément d'image (x) aux régions d'intérêt (R j ) ; - estimation, pour chaque modification élémentaire (δ), de la concentration du traceur dans les régions d'intérêt (R j ) modifiées par la modification élémentaire (δ) ;
- évaluation, pour chaque modification élémentaire (δ), d'un critère (h) permettant de comparer la concentration du traceur estimée avec les concentrations du traceur estimées à au moins une itération précédente ;
- choix de la modification élémentaire (δ op t) donnant une estimation de la concentration du traceur la plus proche des concentrations du traceur estimées aux itérations précédentes ; - comparaison du critère (h) correspondant à la modification élémentaire choisie (δ opt ) avec un seuil prédéterminé (h mιn ) ; et
- application ou non de la modification élémentaire choisie (δ opt ) ou arrêt de la boucle itérative en fonction du résultat de la comparaison du critère (h) avec le seuil (h mm ). 4.- Procédé selon l'une quelconque des revendications 1 à 3, caractérisé en ce qu'il comprend une étape d'ordonnancement des éléments d'image (x) dans chaque domaine spatial (D 1 ) en fonction de critères d'ordre permettant d'évaluer leur risque d'introduire des erreurs non dues au bruit dans l'estimation de la concentration du traceur. 5.- Procédé selon la revendication 4, caractérisé en ce que les critères d'ordre comprennent l'un au moins parmi un premier critère (α(x)) permettant de définir si la concentration contenue dans l'élément d'image (x) provient d'un seul ou de plusieurs domaines spatiaux (D 1 ) différents, un deuxième critère (β(x)) permettant de définir si l'élément d'image (x) participe au calcul des éléments non diagonaux de la matrice de transfert géométrique (W), et un troisième critère (γ(x)) permettant de définir si l'élément d'image (x) est fiable pour l'estimation de la concentration du traceur.
6.- Procédé selon la revendication 5, caractérisé en ce que le premier critère (α(x)) permet de mesurer l'homogénéité de la concentration du traceur dans le voisinage de l'élément d'image (x), le deuxième critère (β(x)) permet de mesurer la contribution des domaines spatiaux (D,) à la mesure de la concentration du traceur dans l'élément d'image (x), et le troisième critère (γ(x)) est issu d'un atlas probabiliste.
7.- Procédé selon l'une quelconque des revendications 4 à 6, caractérisé en ce que la boucle itérative pour la délimitation des régions d'intérêt (R j ) est réalisée suivant l'ordre des éléments d'image (x) obten u à l 'éta pe d'ordonnancement.
8.- Procédé selon l'une quelconque des revendications 1 à 7, caractérisé en ce qu'il comprend une étape d'estimation de la fonction d'étalement du point (PSF) de l'appareil d'imagerie en tout point du champ de vue de l'appareil d'imagerie.
9.- Procédé selon la revendication 8, caractérisé en ce que l'étape d'estimation de la fonction d'étalement du point (PSF) comprend les étapes suivantes :
- mesure de la fonction d'étalement du point (PSF) en plusieurs points du champ de vue ; et
- interpolation et/ou extrapolation des mesures obtenues.
10.- Procédé selon la revendication 8 ou 9, caractérisé en ce qu'il comprend une étape d'estimation d'une fonction d'étalement de région (RSF) pour chaque domaine spatial (D 1 ) par convolution de la fonction d'étalement du point (PSF) avec chaque domaine spatial (D 1 ).
11.- Procédé selon l'une quelconque des revendications 1 à 10, caractérisé en ce qu'il comprend une étape de calcul de la taille des régions d'intérêt (R j ).
12.- Procédé selon la revendication 11 , caractérisé en ce qu'il comprend une étape de délimitation des domaines spatiaux (D,) de manière à obtenir des régions d'intérêt (R j ) de taille maximale.
13.- Procédé selon la revendication 11 ou 12 et l'une quelconque des revendications 4 à 7 prises ensemble, caractérisé en ce qu'il comprend une étape d'optimisation des critères d'ordre de manière à obtenir des régions d'intérêt (R j ) dé taille maximale. 14.- Procédé selon l'une quelconque des revendications 1 à 13, caractérisé en ce que l'image de mesure de la concentration du traceur dans la structure tissulaire est une séquence d'images tridimensionnelles mises en correspondance et comprenant au moins une image tridimensionnelle, et en ce que et les éléments d'image (x) sont des voxels.
15.- Procédé selon la revendication 14, caractérisé en ce que l'image tridimensionnelle est obtenue par tomographie par émission de positrons. 16.- Support d'information comprenant un code pou r estimer la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial (D,) à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt (R j ) à l'intérieur de laquelle la concentration du traceur est mesurée, caractérisé en ce que le code comprend des instructions pour :
- déterminer une matrice de transfert géométrique (W) dont les coefficients (W| j ) sont représentatifs de la contribution des domaines spatiaux (D,) à la mesure de la concentration du traceur dans les régions d'intérêt (R j ) ;
- optimiser les coefficients (w u ) de la matrice de transfert géométrique (W) par délimitation des régions d'intérêt (R j ) les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt (R j ) étant réalisée suivant une boucle itérative comprenant à chaque itération (k) les étapes suivantes :
- modification des régions d'intérêt (R j ); et
- calcul des coefficients (w u ) de la matrice de transfert géométrique (W) à partir des régions d'intérêt (R j ) modifiées ;
- choisir une matrice de transfert géométrique (W) optimisée parmi les matrices de transfert géométrique (W) calculées ; et
- estimer la concentration du traceur à partir de la matrice de transfert géométrique (W) optimisée, lorsqu'il est exécuté par un système de traitement de données.
17.- Dispositif destiné à estimer la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial (D,) à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt (R j ) à l'intérieur de laquelle la concentration du traceur est mesurée, caractérisé en ce qu'il comprend :
- un appareil d'imagerie ; et
- un système de traitement de données comprenant : - des moyens pour déterminer une matrice de transfert géométrique
(W) dont les coefficients (w u ) sont représentatifs de la contribution des domaines spatiaux (D 1 ) à la mesure de la concentration du traceur dans les régions d'intérêt
- des moyens pour optimiser les coefficients (w u ) de la matrice de transfert géométrique (W) par délimitation des régions d'intérêt (R j ) les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt (R j ) étant réalisée suivant une boucle itérative comprenant à chaque itération (k) les étapes suivantes :
- modification des régions d'intérêt (R j ) ; et - calcul des coefficients (w,,) de la matrice de transfert géométrique (W) à partir des régions d'intérêt (R j ) modifiées ;
- des moyens pour choisir une matrice de transfert géométrique (W) optimisée parmi les matrices de transfert géométrique (W) calculées ; et
- des moyens pour estimer la concentration du traceur à partir de la matrice de transfert géométrique (W) optimisée. |
Procédé d'estimation de la concentration d'un traceur dans un ensemble de structures tissulaires, support d'information et dispositif correspondants
La présente invention concerne un procédé d'estimation de la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt à l'intérieur de laquelle la concentration du traceur est mesurée. L'imagerie moléculaire in vivo permet de mesurer des paramètres biochimiques et pharmacologiques locaux de manière non invasive dans des animaux et êtres humains intacts. L'imagerie moléculaire englobe des techniques nucléaires, l'imagerie par résonance magnétique (IRM) et des techniques optiques. Parmi les techniques nucléaires, la tomographie par émission de positrons (TEP) est la technique la plus sensible et la seule à offrir des mesures quantitatives.
La quantification des mesures avec la TEP pour l'étude d'un médicament nécessite l'injection intraveineuse d'un traceur radioactif adapté à l'étude de ce médicament. Une image de la distribution du traceur est alors reconstruite. Les courbes d'activité temporelle (TAC pour Time Activity Curve), c'est-à-dire les variations de la concentration du traceur à l'intérieur d'un élément de volume, sont également fortement corrélées spatialement du fait non seulement du procédé de reconstruction mais également de l'existence de régions physiologiques qui répondent au traceur de façon identique, régions qui peuvent être appelées des pharmaco-organes.
Le post-traitement de l'image TEP requiert en général trois étapes pour accéder au paramètre pharmacologique. Tout d'abord, des régions d'intérêt (ROI pour Région Of Interest) délimitant les organes, plus précisément les organes pharmacocinétiques, doivent être définis soit sur l'image TEP soit sur une image de modalité haute résolution mise en correspondance point à point avec l'image TEP. Ensuite, les TACs des pharmaco-organes doivent être extraites et éventuellement corrigées pour compenser la résolution spatiale limitée de la TEP. Enfin, un modèle physiologique peut être défini, sur la base des TACs et de la
concentration du traceur dans le plasma, permettant le calcul des paramètres pharmacologiques intéressants pour le biologiste. Une délinéation précise des ROIs est nécessaire afin d'extraire des TACs pertinentes.
Cependant, la quantification des TACs est entravée par la résolution limitée du système TEP, résultant en ce que l'on appelle l'effet de volume partiel (PVE pour Partial Volume Effect).
On connaît un procédé par matrice de transfert géométrique (GTM pour Géométrie Matrix Transfer) permettant de corriger le PVE de manière efficace. Ce procédé nécessite la donnée de l'image TEP d'une part, et d'autre part de domaines spatiaux délinéant les organes fonctionnels, de régions d'intérêt au sein desquelles les TACs moyennes des organes fonctionnels seront calculées et d'un modèle de résolution de l'imageur TEP.
Cependant, ce procédé GTM est très sensible non seulement aux erreurs de délinéation desdits domaines spatiaux et aux artefacts de reconstruction de l'image, mais également aux effets de lissage de l'image dus à des mouvements physiologiques périodiques, tels que les battements du cœur ou encore les mouvements respiratoires.
L'invention a pour but de proposer un procédé permettant de limiter l'impact des effets dus aux erreurs de segmentation, aux artefacts de reconstruction de l'image et aux mouvements physiologiques sur l'efficacité du procédé GTM afin d'obtenir des TACs présentant le moins de biais avec la plus faible incertitude due au bruit.
A cet effet, l'invention a pour objet un procédé du type précité, caractérisé en ce qu'il comprend les étapes suivantes : - détermination d'une matrice de transfert géométrique dont les coefficients sont représentatifs de la contribution des domaines spatiaux à la mesure de la concentration du traceur dans les régions d'intérêt ;
- optimisation des coefficients de la matrice de transfert géométrique par délimitation des régions d'intérêt les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt étant réalisée suivant une boucle itérative comprenant à chaque itération les étapes suivantes : - modification des régions d'intérêt ; et
- calcul des coefficients de la matrice de transfert géométrique à partir des régions d'intérêt modifiées ;
- choix d'une matrice de transfert géométrique optimisée parmi les matrices de transfert géométrique calculées ; et - estimation de la concentration du traceur à partir de la matrice de transfert géométrique optimisée.
Le procédé selon l'invention peut comporter une ou plusieurs des caractéristiques suivantes :
- l'étape de modification des régions d'intérêt comprend l'ajout et/ou l'exclusion d'au moins un élément d'image ;
- la boucle itérative pour la délimitation des régions d'intérêt comprend, pour chaque itération, les étapes suivantes :
- détermination d'un ensemble de modifications élémentaires consistant chacune à ajouter au moins un élément d'image aux régions d'intérêt ; - estimation, pour chaque modification élémentaire, de la concentration du traceur dans les régions d'intérêt modifiées par la modification élémentaire ;
- évaluation, pour chaque modification élémentaire, d'un critère permettant de comparer la concentration du traceur estimée avec les concentrations du traceur estimées à au moins une itération précédente ; - choix de la modification élémentaire donnant une estimation de la concentration du traceur la plus proche des concentrations du traceur estimées aux itérations précédentes ;
- comparaison du critère correspondant à la modification élémentaire choisie avec un seuil prédéterminé ; et - application ou non de la modification élémentaire choisie ou arrêt de la boucle itérative en fonction du résultat de la comparaison du critère avec le seuil ;
- le procédé comprend une étape d'ordonnancement des éléments d'image dans chaque domaine spatial en fonction de critères d'ordre permettant d'évaluer leur risque d'introduire des erreurs non dues au bruit dans l'estimation de la concentration du traceur ;
- les critères d'ordre comprennent l'un au moins parmi un premier critère permettant de définir si la concentration contenue dans l'élément d'image provient d'un seul ou de plusieurs domaines spatiaux différents, un deuxième critère
permettant de définir si l'élément d'image participe au calcul des éléments non diagonaux de la matrice de transfert géométrique, et un troisième critère permettant de définir si l'élément d'image est fiable pour l'estimation de la concentration du traceur ; - le premier critère permet de mesurer l'homogénéité de la concentration du traceur dans le voisinage de l'élément d'image, le deuxième critère permet de mesurer la contribution des domaines spatiaux à la mesure de la concentration du traceur dans l'élément d'image, et le troisième critère est issu d'un atlas probabiliste ; - la boucle itérative pour la délimitation des régions d'intérêt est réalisée suivant l'ordre des éléments d'image obtenu à l'étape d'ordonnancement ;
- le procédé comprend une étape d'estimation de la fonction d'étalement du point de l'appareil d'imagerie en tout point du champ de vue de l'appareil d'imagerie ; - l'étape d'estimation de la fonction d'étalement du point comprend les étapes suivantes :
- mesure de la fonction d'étalement du point en plusieurs points du champ de vue ; et
- interpolation et/ou extrapolation des mesures obtenues ; - le procédé comprend une étape d'estimation d'une fonction d'étalement de région pour chaque domaine spatial par convolution de la fonction d'étalement du point avec chaque domaine spatial ;
- le procédé comprend une étape de calcul de la taille des régions d'intérêt ;
- le procédé comprend une étape de délimitation des domaines spatiaux de manière à obtenir des régions d'intérêt de taille maximale ;
- le procédé comprend une étape d'optimisation des critères d'ordre de manière à obtenir des régions d'intérêt de taille maximale ;
- l'image de mesure de la concentration du traceur dans la structure tissulaire est une séquence d'images tridimensionnelles mises en correspondance et comprenant au moins une image tridimensionnelle, et les éléments d'image sont des voxels ; et
- l'image tridimensionnelle est obtenue par tomographie par émission de positrons.
L'invention a également pour objet un support d'information comprenant un code pour estimer la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt à l'intérieur de laquelle la concentration du traceur est mesurée, caractérisé en ce que le code comprend des instructions pour :
- déterminer une matrice de transfert géométrique dont les coefficients sont représentatifs de la contribution des domaines spatiaux à la mesure de la concentration du traceur dans les régions d'intérêt ;
- optimiser les coefficients de la matrice de transfert géométrique par délimitation des régions d'intérêt les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt étant réalisée suivant une boucle itérative comprenant à chaque itération les étapes suivantes :
- modification des régions d'intérêt ; et
- calcul des coefficients de la matrice de transfert géométrique à partir des régions d'intérêt modifiées ;
- choisir une matrice de transfert géométrique optimisée parmi les matrices de transfert géométrique calculées ; et
- estimer la concentration du traceur à partir de la matrice de transfert géométrique optimisée, lorsqu'il est exécuté par un système de traitement de données. L'invention a en outre pour objet un dispositif destiné à estimer la concentration d'un traceur dans un ensemble de structures tissulaires comprenant au moins une structure tissulaire à partir d'une image de mesure de la concentration du traceur dans ledit ensemble de structures tissulaires obtenue par un appareil d'imagerie, l'image comprenant au moins un domaine spatial à l'intérieur duquel la concentration du traceur est homogène et au moins une région d'intérêt à l'intérieur de laquelle la concentration du traceur est mesurée, caractérisé en ce qu'il comprend :
- un appareil d'imagerie ; et
- un système de traitement de données comprenant :
- des moyens pour déterminer une matrice de transfert géométrique dont les coefficients sont représentatifs de la contribution des domaines spatiaux à la mesure de la concentration du traceur dans les régions d'intérêt ;
- des moyens pour optimiser les coefficients de la matrice de transfert géométrique par délimitation des régions d'intérêt les meilleures en termes d'erreurs pour mesurer la concentration du traceur, la délimitation des régions d'intérêt étant réalisée suivant une boucle itérative comprenant à chaque itération les étapes suivantes :
- modification des régions d'intérêt ; et - calcul des coefficients de la matrice de transfert géométrique à partir des régions d'intérêt modifiées ;
- des moyens pour choisir une matrice de transfert géométrique optimisée parmi les matrices de transfert géométrique calculées ; et
- des moyens pour estimer la concentration du traceur à partir de la matrice de transfert géométrique optimisée.
L'invention sera mieux comprise à la lecture de la description qui va suivre, donnée uniquement à titre d'exemple et faite en se référant aux dessins annexés, sur lesquels :
- la Figure 1 est un organigramme représentant les deux phases principales du procédé selon l'invention ;
- la Figure 2 est un schéma illustrant l'image TEP d'une souris ;
- la Figure 3 est un organigramme représentant les étapes du procédé selon l'invention ;
- les Figures 4 et 5 sont des coupes respectivement transaxiale et longitudinale du champ de vue du scanner TEP ;
- la Figure 6 est un organigramme représentant plus en détail l'étape de délimitation des régions d'intérêt ; et
- les Figures 7A à 11 B sont des schémas présentant les résultats obtenus avec un exemple de réalisation du procédé selon l'invention. La présente invention se base sur le procédé par matrice de transfert géométrique (GTM) de correction de l'effet de volume partiel (PVE) proposé par Rousset et al. dans l'article « Correction for Partial Volume Effects in PET : Principle and Validation », The Journal of Nuclear Medicine, 1998, Vol. 39, No. 5,
pages 904 à 911 , et implémenté dans l'espace image par Frouin et al. dans l'article « Correction of Partial-Volume Effect for PET Striatal Imaging : Fast Implementation and Study of Robustness », The Journal of Nuclear Medicine, 2002, Vol. 43, No. 12, pages 1715 à 1726. Le but du procédé GTM est de retrouver la vraie concentration T 1 du traceur au sein de domaines spatiaux D 1 . Les détails du procédé GTM peuvent être trouvés dans l'article de Rousset et al. précité, mais les principes en sont présentés ici afin de mettre en avant les caractéristiques de la présente invention. Le procédé selon l'invention a pour but d'optimiser les paramètres de calcul de la vraie concentration T, du traceur en limitant les effets dus aux erreurs de segmentation, aux artefacts de reconstruction de l'image et aux mouvements physiologiques.
Comme représenté sur la Figure 1 , le procédé selon l'invention comprend deux phases principales : - une première phase 10 destinée à choisir la fonction de réponse du scanner TEP de manière adéquate ; et
- une deuxième phase 12 destinée à choisir les régions d'intérêt optimales pour l'estimation de la vraie concentration T 1 du traceur.
En référence à la Figure 2, l'image TEP est une image tridimensionnelle formée d'un ensemble de voxels X(x,y,z) chacun caractérisé par une concentration de traceur A(x,y,z), cette concentration étant susceptible de varier au cours du temps.
L'image TEP est supposée être composée de I domaines spatiaux homogènes λPih ≤ i ≤ i en termes de concentration de traceur. Chaque domaine spatial D 1 a une vraie concentration de traceur notée T 1 .
Une concentration de traceur if/j 1<; KJ est mesurée dans l'image TEP au
sein de régions d'intérêt (ROI) IR 7 | 1≤7≤J .
Pour des raisons de clarté de l'explication, nous considérerons dans la suite que chaque ROI est incluse dans le domaine spatial correspondant et que leur nombre sont égaux : J = I.
En variante, J est différent de I, avec de préférence J > I pour que le procédé GTM permette de retrouver la vraie concentration T, du traceur.
La concentration de traceur mesurée t j
dans une région d'intérêt R j
, de volume V j
, peut être exprimée comme une combinaison linéaire des concentrations du traceur V z
Ji< ?
< /
dans les différents domaines spatiaux :
Cette équation peut être réécrite sous forme matricielle : t = W.T (2) où t est le vecteur des mesures TEP de concentration du traceur moyennées dans les ROIs R j , T est le vecteur des vraies concentrations du traceur dans les domaines spatiaux D, et W est la matrice de transfert géométrique
(GTM) I x J entre les domaines spatiaux IA Ji ≤z≤/ et les ROIs Wj i ι≤J≤J . Les coefficients w u de la matrice W sont égaux à : w ;. J = ^ r \ RSFXx)CbC (3)
V V J R,
La vraie concentration du traceur peut alors être calculée comme : T = W "1 .t (4)
Comme illustré sur la Figure 3, une étape initiale 14 de la phase 12 du procédé selon l'invention consiste à délimiter les domaines spatiaux D 1 .
L'utilisateur est libre de délimiter les D 1 comme il le souhaite, en les traçant manuellement, semi-automatiquement ou automatiquement sur l'image TEP, mais sous la contrainte que chaque domaine spatial D, soit homogène en termes de vraie concentration du traceur, cette concentration pouvant varier au cours du temps dans le cas où l'image TEP contient une information dynamique.
En variante, les D, sont tracés sur une image provenant d'une autre modalité d'imagerie, telle que NRM nucléaire ou encore le scanner X.
Parallèlement, la phase 10 de choix de la fonction de réponse du scanner TEP est mise en œuvre.
La fonction de réponse du scanner, ou fonction d'étalement du point (PSF pour Point Spread Function), est de préférence estimée pour que le procédé GTM fonctionne efficacement.
En référence aux Figures 4 et 5, la PSF est décomposée dans un système de coordonnées cylindriques (p,θ,z) en :
- une composante axiale le long de l'axe z : PSF axιa ie ;
- une composante radiale dans la direction CX : PSF raC iιaie ; et
- une composante tangentielle dans la direction XT : PSF ta ngentιeiie-
Les composantes PSF axιa ie, PSF radιa ι e et PSF tan gentιeiie de la PSF peuvent comprendre chacune plusieurs paramètres, par exemple deux variances de deux gaussiennes si la PSF est estimée par deux gaussiennes, ou plusieurs valeurs.
Si l'on suppose que la PSF est invariante selon z, il n'est nécessaire de mesurer la PSF que sur le plan transaxial passant par le milieu du scanner TEP.
Selon cette hypothèse, PSF axιa ι e ne varie pas avec la position dans le champ de vue, et de plus PSF radιa ie et PSFt an gentιeiie sont invariantes selon θ . Il suffit donc de mesurer PSF axιa ι e en un point du champ de vue, et PSF radιa ι e et PSF ta ngentιeiie en plusieurs points de l'axe x.
Si l'on suppose que la PSF varie lorsqu'on se déplace dans la direction de l'axe du scanner TEP, il est nécessaire de mesurer PSF axιa ι e en plusieurs points de l'axe z, et PSF radιa ι e et PSF ta ngentιeiie en plusieurs points de l'axe x.
En référence à la Figure 3, une première étape 16 de la phase 10 du procédé selon l'invention consiste à mesurer PSF axιa ι e , PSF radιa ι e et PSFt an gentιeiie en plusieurs points du champ de vue du scanner TEP.
Selon l'invention, la PSF peut avoir une forme paramétrique ou non paramétrique quelconque, ou plus précisément une forme présentant une bonne adéquation avec la forme de PSF générée par le procédé de reconstruction de l'image TEP, et éventuellement également avec l'effet de lissage opéré par les mouvements physiologiques. Si l'image TEP est reconstruite avec modélisation ou compensation de la PSF, la PSF dans le cadre de cette présente invention modélisera les effets résiduels non pris en compte lors de la reconstruction (par exemple si la reconstruction suppose la PSF uniforme dans le champ de vue, alors l'écart pour chaque point ou région du champ de vue entre la PSF uniforme et la PSF réelle sera modélisé).
Dans une deuxième étape 18 de la phase 10, on déduit de ces mesures la valeur de PSF axιa ie, PSF raC iιaie et PSF ta ngentιeiie en tout point du champ de vue par interpolation et/ou par extrapolation, par une méthode au choix de l'utilisateur.
Comme vu précédemment, PSF axi aie ne dépend que de la position sur l'axe z et sera noté PSF axi aie(z), alors que PSF raC iιaie et PSF ta ngentιeiie ne dépendent que de p et seront donc notées PSF radιa ι e (p) et PSF tan gentιeiie(p)-
A l'étape 20 du procédé selon l'invention, on estime la RSF,(x) en tout point du champ de vue.
Cette estimation est obtenue par convolution du masque de convolution correspondant à la PSF avec les voxels du domaine spatial D 1 .
Le masque de convolution correspondant à la PSF du scanner TEP est obtenu de la façon suivante :
Soient C = (x c ,y c ) et X = (x,y,z) respectivement les coordonnées de l'axe du scanner et celles d'un voxel situé dans le champ de vue du scanner. Soit δX = (δx,δy,δz) un petit déplacement autour de X.
La valeur en X + δX du masque de convolution correspondant à la PSF en X est donnée par :
PSF x y z (&,δy,δz) = f(δr,PSF n ^(j)),a,PSF^ sen ^(j)),z,PSF aaale ( < z)) (5)
Où δr = âccosθ + δysinθ , â = -&sïnθ + δycosθ ,
PSF x y z {&c,δy,δz) au masque de D 1 . Ce lissage est non homogène et reproduit les effets produits par l'acquisition et la reconstruction de l'image.
On définit parallèlement, à l'étape 22, un ensemble λ, de voxels dont l'activité mesurée dans leur voisinage est représentative de la concentration de traceur dans le domaine spatial D,. Ces voxels peuvent être déterminés de manière automatique ou de manière manuelle.
On définit également S, comme un ensemble de voxels à exclure des régions d'intérêt R j .
Un ordonnancement des voxels dans chaque domaine spatial D, est alors réalisé à l'étape 24 selon le principe suivant. Au moins un voxel x e D, est affecté par les effets qui lissent les images
(effet de volume partiel, fraction des tissus dans le voxel, mouvement physiologiques, etc.).
Il est possible de prédire, avec ou sans introduction d'information a priori, dans quelle proportion ce voxel x est affecté par ces effets à partir d'un ou plusieurs critères d'ordre prédéterminés.
Parmi ces critères :
- critère α(x) : le voxel x est affecté par l'effet de volume partiel, la fraction tissulaire ou les mouvements physiologiques périodiques du sujet, ce critère permettant de définir en d'autres termes si l'activité contenue dans le voxel x provient d'un seul ou de plusieurs domaines spatiaux différents ;
- critère β(x) : le voxel x participe au calcul des éléments non diagonaux de la matrice W, ce critère pouvant aisément être calculé à partir du calcul des RSF k (x) pour k ≠ i ;
- critère γ(x) : tout autre critère ou ensemble de critères pouvant aider à discriminer les voxels les plus fiables pour l'estimation de la concentration du traceur (fiabilité de la mesure TEP, homogénéité de réponse au sein d'un même organe, etc.), par exemple issu d'un atlas probabiliste.
Les critères α(x), β(x) et γ(x) sont définis de telle manière que plus ils sont petits et moins l'agrégation de x à la région d'intérêt R 1 va introduire d'erreur de correction due aux erreurs de segmentation ou aux effets non pris en compte dans l'estimation de la PSF.
Soit g,(α(x), β(x), γ(x), λ h S 1 ) un ordre des voxels x e D 1 donnant à x un rang d'autant plus faible que α(x), β(x) et γ(x) sont petits, sous la contrainte que :
(C1 ) Un voxel x e D, de rang n est connecté par des voxels de rang inférieur à l'un des points de λ, ; et
(C2) Les voxels y e S 1 sont exclus de cet ordre de tri.
Dans la suite, R, (w) désignera une ROI contenant les n voxels de rangs g,(α(x), β(x), γ(x)) les plus petits parmi tous les voxels de D 1 ou une partie de ces voxels.
R, (w) peut également contenir des voxels situés en dehors de D, mais à sa périphérie.
L'optimisation automatique des régions d'intérêt R j est ensuite effectuée à l'étape 26 du procédé selon l'invention suivant une boucle itérative.
On définit tout d'abord les paramètres utiles pour l'optimisation automatique des régions d'intérêt R j . Soit Rj <= Ri , soit T' et T les concentrations du traceur estimées dans les organes respectivement avec R, et avec R 1 . Soit T la vraie concentration du traceur. On peut avancer que :
- Les erreurs dues aux erreurs de segmentation sur l'estimation de T en utilisant R, sont moins importantes que celles commises en utilisant R 1 . De plus, la matrice de transfert géométrique W estimée en utilisant R, est plus proche de la matrice unité que la matrice de transfert géométrique W estimée en utilisant R, . Le conditionnement de la matrice W est meilleur que celui de W.
- Les erreurs dues au bruit commises sur l'estimation de T en utilisant R, sont plus importantes que celles commises en utilisant R, . - II est possible d'estimer T l k τ , qui est la concentration du traceur dans la région d'intérêt R 1 pour l'itération k et le temps de mesure τ :
- soit en optimisant R 1 globalement en intégrant l'information de tous les temps de mesure (dénommés « frames » dans la suite) ;
- soit en optimisant R 1 séparément pour chaque frame τ ; - soit en optimisant R 1 séparément pour chaque frame τ op t tout en intégrant également une information provenant des autres frames τ mes -
Le dernier cas correspond à un cas général dont les deux premiers cas sont des cas particuliers. Nous ne considérerons donc dans la suite que le cas général. n ( n ι k )
On définit au préalable κ ι,s h ' k comme la ROI contenant à l'itération k les n, ιk voxels présentant les voxels x de plus petit rang parmi tous les voxels de D 1 à l'exclusion des voxels appartenant à l'ensemble S,,k-
Soit δ une modification élémentaire du nombre de voxels des régions d'intérêt R 1 et soit δ k un ensemble de modifications élémentaires possibles à l'itération k.
Soit ' t k τ ° pt la ROI définie à l'itération k pour l'optimisation de la région n d'intérêt R 1 pour la frame τ op t, contenant α ^' voxels et excluant les voxels de
l'ensemble ι ' k ' τ ° pι tout en restant conforme à (C1 ) et (C2).
T Soit ! ' k ' T< "" la concentration du traceur que l'on cherche à estimer pour la
frame τ opt et soit α -v.^ la concentration du traceur estimée dans ïλ t ^' au temps τ m es-
Soit σ x, Tmβs la variance du bruit au voxel x au temps τ me s- Cette variance peut être déterminée en fonction de l'algorithme de reconstruction de l'image et éventuellement du signal TEP mesuré au voxel x au temps τ me s-
Soit ^i, v , τ _ la concentration de traceur moyenne mesurée dans la région
d'intérêt R 1
à l'itération k au temps τ mes
pour l'estimation de
Sous certaines hypothèses, la variance ζ^ k
, τop
,, τm
de l'erreur due au bruit commise sur l'estimation de
.2
<r XlTmB et de » α>v .
Connaissant cette erreur, le procédé GTM fournit une borne supérieure de la variance ζl k,τo l,τ^ de l'erreur due au bruit commise sur l'estimation de ^ι,k,τ , :
Soit
différent de fc, v
,r_
Soit h mιn un seuil en dessous duquel Vα, v ,^j 1<7 <r i< , <7 est jugé
significativement différent de ψi,ι,τ opt ,τ mes \ ≤ ι <kA≤Tmes ≤TA≤i≤I -
L'algorithme de la boucle itérative pour l'optimisation des régions d'intérêt R j est le suivant (Figure 6). L'étape 28 correspond à l'initialisation de la boucle.
Les régions d'intérêt R j sont initialisées par R ι,s ° *° " , où « î;O;V > 2 et où S,,o est un ensemble vide.
A chaque itération k, on passe ensuite à une étape 30 de détermination d'un ensemble δ k de modifications élémentaires δ possibles des ψιj * k ',^ | <;</ qui peuvent résulter en l'addition de voxels à une région d'intérêt R 1 et/ou en l'addition de voxels à S, et/ou en la modification de A 1 .
A l ' ét a p e 32 , p o u r to u t δ e δ, on calcule fc"l l≤i≤I W,
A l'étape 36, on choisit la modification élémentaire δ op t qui donne une estimation de la concentration du traceur la plus proche des concentrations du traceur estimées aux itérations précédentes :
A l'étape 38, on compare le critère h correspondant à δ opt avec h r
Si , alors on applique la
Sinon l'algorithme se termine (étape 42), et on a : )
- I 7 ^V -L≤/ rU- i .v>v {≤,≤/ ; et
- I Pf ' Op, \ f ≤ ι≤I =k P ' 2 * 1 ^ ^ A f≤i≤/ '
En référence à la Figure 3, on obtient ainsi à l'étape 44 l'estimation de la vraie concentration T 1 du traceur dans les domaines spatiaux D 1 .
On obtient également, à partir de l'étape 26 d'optimisation des régions d'intérêt R j , la taille finale de ces régions d'intérêt (étape 46) qui permet d'optimiser les paramètres du procédé.
En effet, le procédé comprend divers paramètres déterminant l'ordre dans lequel les voxels vont être agrégés aux ROIs. Parmi ces paramètres, on trouve les A 1 , les fonctions α(x), β(x), γ(x) et g,(α(x), β(x), γ(x), A 1 , S 1 ).
Le procédé d'optimisation réalise un compromis entre les erreurs d'estimation dues au bruit et les erreurs d'estimation provenant des effets tels que les erreurs de segmentation et les mouvements physiologiques.
Les erreurs dues au bruit ont tendance à décroître quand la taille des ROIs augmente.
La taille finale des R j est donc un indicateur de :
- la qualité des D 1 ;
- l'homogénéité du signal TEP au sein des D 1 ; et
- la qualité des paramètres déterminant l'ordre des voxels. A D, et homogénéité du signal fixes, la taille finale des R j est donc un indicateur de la qualité du choix des A 1 , les fonctions α(x), β(x), γ(x) et g,(α(x), β(x), γ(x), A 1 , S 1 ).
On peut donc soit :
- lancer plusieurs optimisations des R j avec des D 1 , A 1 , α(x), β(x), γ(x) et gι(oc(x), β(x), γ(x), A 1 , S 1 ) différents et choisir les D 1 , A 1 , α(x), β(x), γ(x) et g,(α(x), β(x), γ(x), A 1 , S 1 ) pour lesquels la taille des R j , ou toute fonction monotone croissante de cette taille ou toute autre fonction évaluant de manière quantitative la qualité de l'estimation des T 1 , est maximale (étape 48) ;
- donner la taille des R j , ou la valeur de toute fonction monotone de cette taille ou toute autre fonction évaluant de manière quantitative la qualité de l'estimation des T 1 , comme indicateur global de la qualité des D 1 , A 1 , α(x), β(x), γ(x) et g,(α(x), β(x), γ(x), A 1 , S 1 ) (étape 50).
Un exemple de réalisation du procédé selon l'invention est le suivant :
- RSF|(x) est estimée selon la méthode de Frouin et al. ; - J=l et j=i ;
- les D, sont définis comme les domaines spatiaux délimités par la segmentation en utilisant un procédé basé sur l'analyse de moyennes locales (LMA pour Local Mean Analysis) proposé par Maroy et al. dans l'article « Segmentation of Rodent Whole-Body Dynamic PET Images : An Unsupervised Method Based on Voxel Dynamics », IEEE Trans. Med. Imaging, 2008 ;
- λ, est l'ensemble des voxels x e D 1 extraits lors de l'étape d'extraction de minima locaux de l'image â p 2 mesurant à chaque voxel p l'homogénéité de la concentration du traceur dans le voisinage de p ;
- α(x) est l'ordre des voxels triés par â x 2 croissant ;
- β(x) est l'odre des voxels triés par RSF 1 (X) décroissant ;
- g,(α(x), β(x), γ(x), A 1 , S 1 ) est l'ordre de x parmi les voxels x e D 1 dont les voxels de S 1 sont exclus et qui sont triés par (axα 2 (x)+ (1-a)xβ 2 (x)) croissant. Dans cet ordre, les voxels ne vérifiant pas (C1) voient leur rang décalé jusqu'à ce qu'ils puissent vérifier (C1 ) ;
- A k
- τ op t est la frame pour laquelle on optimise les {^} 1≤;≤7 et τ mes est une frame où l'on estime les {^}i ≤*≤/ ;
- soit L un niveau d'échelle et L max le niveau d'échelle maximum, fixé a priori (par exemple L max = T/4). On associe à la frame x op t au niveau d'échelle 0 un ensemble F L de L+1 frames consécutives au niveau d'échelle L, ces L frames étant centrées sur x op t ; soit p l|L la p-value associée à la statistique t
L'algorithme est le suivant.
A l'étape d'initialisation : δ, = 0,001 x#D,, où #D, est le cardinal du domaine spatial segmenté D 1 , n i,o,τ opt = δ h R 1 est initialisée par R^'^ et S,, o est un ensemble vide. A chaque itération k et pour tout i :
- on calcule fc"L, , W, te.,_L_ ffJM , ^Zλ ^m , et
Si h mm
, alors on ajoute d.
,(», J voxels à w ; £ r , on redétermine ^ s ' / ° "' et on passe à l'itération suivante k+1.
Sinon, si d, > 1 , on diminue d, et on retente une itération. Sinon, si #S, < s ma χ, on ajoute les d, voxels ci-dessus à S 1 , on augmente d, et on retente une itération.
Sinon, on stoppe l'algorithme et la solution est :
- R, = Oïl' :
- K^, j; ≤î≤/
Kt-i, Wop
, j 1≤î≤/
: et
Les résultats ont été obtenus avec le matériel suivant. Pour la simulation, nous avons simulé au moyen d'un simulateur analytique cinquante acquisitions TEP du fantôme de souris MOBY sur un scanner TEP dédié au petit animal (le scanner simulé était le FOCUS 220 de Siemens avec une résolution spatiale de 1 ,3mm), avec des frames temporelles d'une minute chacune. Les images ont été reconstruites en utilisant le procédé OSEM (pour Ordered Subsets Expectation Maximization), qui est un procédé classique de reconstruction d'image TEP, avec des voxels de 0,5x0,5x0,8 mm 3 .
Les données expérimentales incluent seize souris injectées avec des xénogreffes de tumeurs humaines péritonéales. Les souris ont reçu une injection de Fluoro-thymidine (18F-FLT), immédiatement suivie par une acquisition TEP de
5400 min composée de 18 frames temporelles (cinq frames d'une minute, cinq frames de deux minutes, trois frames de cinq minutes, trois frames de dix minutes, et deux frames de quinze minutes). Les images ont été reconstruites en utilisant le procédé OSEM avec des voxels de 0,5x0,5x0,8 mm 3 . Les souris ont été sacrifiées immédiatement après l'acquisition, et la concentration du traceur a été comptée dans les organes prélevés. Cette concentration comptée dans les organes a servi de « gold standard » pour la concentration estimée à la fin de l'acquisition.
Les images TEP ont été segmentées en utilisant le procédé LMA (Figure 7A pour la simulation et Figure 7B pour les données expérimentales). L'estimation de la concentration du traceur dans les organes est effectuée selon trois procédés :
- calcul de la TAC moyenne au sein du domaine spatial segmenté correspondant à l'organe ;
- correction de volume partiel opérée par le procédé GTM implémenté dans l'espace image tel que décrit par Frouin et al. (nommé simplement « GTM » dans la suite) ; et
- correction de volume partiel opérée par le procédé proposé dans l'exemple de réalisation de la présente invention (nommé « LMA-GTM » dans la suite). Les procédés GTM et LMA-GTM ont été comparés sur la base de :
- ETR : qui est l'écart en valeur absolu entre la mesure correspondante du gold standard et la valeur mesurée, exprimé en pourcentage de la mesure gold standard ; et
- ARC : qui est le ratio de recouvrement apparent, en d'autres termes la valeur mesurée divisée par la valeur correspondante du gold standard.
Les Figures 8A et 8B montrent des exemples de variation temporelle de la concentration de traceur (TACs) estimée pour les simulations (Figure 8A) et pour les données expérimentales (Figure 8B).
Les Figures 9A et 9B montrent que le procédé LMA-GTM améliore de manière significative la précision des mesures (la p-value est inférieure à 10 ~5 tant pour les simulations que pour les données expérimentales) et que le procédé
GTM n'améliore que les TACs extraites pour les données expérimentales seulement.
La précision des mesures estimées au moyen du procédé LMA-GTM est bonne (ETR Slm ulatιons = 5,3% ± 8% et ETR Dθ nnées expérimentales = 4,8% ± 7%) et le contraste apparent est correctement récupéré (ARC Slm uiatιo n s = 94% ± 10% et
ARCoonnées expérimentales = 99,98% ± 9%). Les résultats obtenus en utilisant le procédé LMA-GTM sont significativement meilleurs en termes d'ARC et d'ETR que ceux obtenus par le procédé GTM (la p-value est inférieure à 10 ~5 ).
Le processus complet d'extraction de TAC corrigée du volume partiel prend 10 minutes par image : - segmentation (-15 secondes) ;
- nommage des organes par l'utilisateur (~5 minutes) ; et
- correction de volume partiel (~4 minutes).
La Figure 10A montre l'ARC obtenu par les trois procédés sur les images de simulations pour chaque organe. Pour le procédé LMA-GTM, l'ARC varie entre 94% ± 5% et 99,3% ± 1 % pour tous les organes à l'exception des petits organes, tels que le thalamus et la thyroïde, et le pancréas du fait de sa forme.
L'ETR (Figure 10B) est inférieur à 10% de la valeur vraie pour tous les organes à l'exception du thalamus (ETR = 32% ± 5%) et de la thyroïde (ETR = 21 % ± 6%).
La précision obtenue au moyen du procédé LMA-GTM est meilleure que celle obtenue par le calcul de la TAC moyenne et par la procédé GTM, et ce avec une p-value inférieure à 10 ~5 pour tous les organes.
Les résultats obtenus pour les données expérimentales sont similaires à ceux obtenus pour les simulations.
La récupération de contraste (Figure 11A) était bonne (entre 98,4% ± 7% et 110% ± 6%).
De plus, l'ETR (Figure 11 B) était inférieur à 10% pour le procédé LMA-GTM pour tous les organes. Dans les deux cas, le procédé LMA-GTM obtient des meilleurs résultats que le procédé GTM (p-value inférieure à 0,002).
Pour les études du cerveau, le procédé GTM est considéré comme un voire le procédé de référence pour la correction de l'effet de volume partiel dans les images TEP.
Pour les études chez le rongeur, il n'existe pas à notre connaissance de procédé générique (applicable à tous les traceurs et à toutes les pathologies) de correction de volume partiel.
Le procédé LMA-GTM montre un fort potentiel pour de telles études, mais aussi probablement pour des études chez l'homme, dans des image corps entier ou cerveau. En effet, le procédé LMA-GTM :
- utilise la délinéation automatique LMA d'organes dans les images TEP ;
- ne nécessite pas d'information anatomique ni d'à priori ;
- est robuste aux erreurs de segmentation et aux erreurs dues aux effets non modélisables par la correction de PVE ; et - est rapide, facile à utiliser et précise.
L'invention propose donc un procédé permettant d'estimer de manière fiable et précise la concentration d'un traceur dans un ensemble de structures tissulaires en délimitant de manière optimale les régions d'intérêt dans lesquelles la concentration du traceur est mesurée.