FR 2952452 A1 20110513 FR 0905759 A 20091130 Procédé pour estimer des paramètres hémodynamiques par estimation conjointe des paramètres d'un modèle global de perfusion contraint L'invention concerne un système et un procédé pour estimer des paramètres hémodynamiques de perfusion. L'invention s'appuie notamment sur des techniques d'imagerie de Perfusion par Résonance Magnétique (Perfusion Weighted Magnetic Resonance Imaging (PW-MRI) en langue anglaise) ou par Tomodensitométrie (Computed Tomography (CT) en langue anglaise). Ces techniques permettent d'obtenir rapidement des informations précieuses sur l'hémodynamique d'organes tels que le cerveau ou le coeur. Ces informations sont particulièrement cruciales pour un praticien cherchant à établir un diagnostic et à prendre une décision thérapeutique dans le traitement en urgence de pathologies telles que les accidents vasculaires cérébraux. Pour mettre en oeuvre de telles techniques, on utilise un appareil d'imagerie par Résonance Magnétique Nucléaire ou par Tomodensitométrie. Celui-ci délivre une pluralité de séquences d'images numériques d'une partie du corps, notamment du cerveau. Ledit appareil applique pour cela une combinaison d'ondes électromagnétiques à haute fréquence sur la partie du corps considérée et mesure le signal réémis par certains atomes. L'appareil permet ainsi de déterminer la composition chimique et donc la nature des tissus biologiques en chaque point (ou voxel) du volume imagé. Des séquences d'images sont analysées au moyen d'une unité de traitement dédiée. Cette unité de traitement délivre in fine à un praticien, une estimation des paramètres hémodynamiques à partir des images de perfusion, au moyen d'une interface homme-machine adaptée. Le praticien peut ainsi réaliser un diagnostic et décider de l'action thérapeutique qu'il jugera adéquate. Des images de perfusion par Résonance Magnétique Nucléaire ou par Tomodensitométrie sont obtenues en injectant un agent de contraste (par exemple un sel de gadolinium pour l'Imagerie par Résonance Magnétique) par voie intraveineuse et en enregistrant son bol au cours du temps au niveau de chaque voxel de l'image. Par souci de concision, nous omettrons les indices x,y,z pour identifier des voxels. Par exemple, au lieu de noter SxY z (t) le signal pour un voxel de coordonnées x,y,z, nous le noterons simplement S(t). Il est entendu que les opérations et les calculs décrits dans la suite sont généralement effectués pour chaque voxel d'intérêt, de sorte à obtenir au final des images ou des cartes représentatives des paramètres hémodynamiques que l'on cherche à estimer. Un modèle standard permet de relier l'intensité des signaux S(t) mesurée au cours du temps t, à la concentration C(t) dudit agent de contraste. Par exemple, en Tomodensitométrie de Perfusion, le signal pour chaque voxel est directement proportionnel à la concentration : S(t)=a.C(t), a étant une constante non nulle. En Imagerie de Perfusion par Résonance Magnétique Nucléaire, il existe une relation exponentielle S(t)=So.e-k.TE.C(t) où So est l'intensité moyenne du signal avant l'arrivée de l'agent de contraste, TE est le temps d'écho (echo time en langue anglaise) et k est une constante dépendant de la relation entre la susceptibilité paramagnétique et la concentration de l'agent de contraste dans le tissu. La valeur de la constante k pour chaque voxel étant inconnue, celle-ci est fixée à une valeur arbitraire pour tous les voxels d'intérêt. On obtient ainsi des estimations relatives et non pas absolues. Lesdites informations relatives restent pertinentes puisque l'on est intéressé principalement par la variation de ces valeurs relatives dans l'espace. La concentration C(t) peut alors être exprimée elle-même par le modèle de convolution standard de perfusion C(t)=BF Ca(t) R(t) où Ca(t) est la concentration de l'agent de contraste dans l'artère alimentant le volume de tissu dans un voxel (Fonction d'Entrée Artérielle ou Arterial Input Function (AIF) en langue anglaise), BF est le Flux Sanguin dans le volume de tissu (Blood Flow en langue anglaise), R(t) est la Fonction de Répartition Complémentaire du temps de transit dans le volume de tissu (residue function en langue anglaise) et 0 désigne le produit de convolution. Partant d'un signal expérimental S(t), des procédés connus consistent à estimer tout d'abord So en prenant, par exemple, sa moyenne avant que l'agent de contraste n'arrive. On obtient ainsi une estimation de la concentration de l'agent de contraste C(t) grâce à la relation C(t)=- 1 ln k•TE Connaissant C(t) et la fonction d'entrée artérielle Ca(t), on peut alors obtenir par déconvolution numérique la fonction BF-R(t) sous le modèle standard de perfusion. Par suite, on peut obtenir une estimation du paramètre BF puisque R(0)=1 par définition. On obtient en outre une estimation du Temps de Transit Moyen dans le tissu (Mean Transit Time (MTT) en langue anglaise), S (t) So +ao puisque MTT = J R(t)dt. On peut obtenir également une 0 estimation du volume sanguin dans le tissu (Blood Volume (BV) en langue anglaise) par la relation BV =BF.MTT . Selon l'état de la technique, il ne semble pas possible de concevoir une unité de traitement apte à mettre en oeuvre un procédé de déconvolution du modèle standard de perfusion C(t) = BF•Ca (t) ©R(t) sans devoir fournir une fonction d'entrée artérielle Ca(t). En effet ledit modèle est exprimé sous la forme d'une seule 10 équation et il existe deux inconnues. Des techniques ont donc été envisagées pour fixer une ou plusieurs Fonctions d'Entrée Artérielle Ca(t) et rendre ainsi possible la déconvolution de courbes de 15 concentration C(t) afin d'estimer les paramètres hémodynamiques. Pour estimer des paramètres hémodynamiques d'un organe tel que le cerveau humain, selon une première 20 technique, une fonction d'entrée artérielle globale est choisie manuellement par un praticien. Le choix peut porter par exemple sur l'artère sylvienne contralatérale à l'hémisphère pathologique dans le cas de l'imagerie de perfusion du cerveau. On peut obtenir en variante une 25 fonction d'entrée artérielle globale au moyen de mesures supplémentaires, par exemple optiques. Si elle permet d'obtenir des signaux avec hauts rapports signal sur bruit, cette technique induit néanmoins de nombreux inconvénients. Tout d'abord, elle 30 nécessite une intervention humaine et/ou des mesures supplémentaires. Ceci peut être très pénalisant en situation d'urgence clinique. En outre, les procédures et les résultats finaux sont difficilement reproductibles. Ensuite et surtout, la fonction d'entrée artérielle globale ne correspond pas aux fonctions d'entrée artérielle locales pour les tissus d'intérêt. Elle en diffère en terme de délai car les fonctions d'entrée artérielle locales sont le plus souvent en retard par rapport à la fonction d'entrée artérielle globale prise en amont. Elle en diffère en outre en terme de dispersion car la propagation de l'agent de contraste est plus lente en aval qu'en amont. Or il est connu que de tels phénomènes ont une influence considérable sur les estimations des paramètres hémodynamiques. En particulier, selon cette première technique, on n'obtient pas une estimation du véritable Temps de Transit Moyen MTT entre la fonction d'entrée artérielle locale et la fonction de sortie veineuse. On obtient seulement un Temps de Transit Moyen entre la fonction d'entrée artérielle globale et la fonction de sortie veineuse. Pour tenter de quantifier ces discordances par rapport au modèle standard de perfusion, certains ont introduit de nouveaux paramètres descriptifs, tels qu'un paramètre TMAX =argmax R(t), quand bien même lesdits paramètres t n'interviennent pas dans le modèle de perfusion standard. D'autres méthodes tendent à minimiser l'influence de ces discordances sur l'estimation des paramètres hémodynamiques mais introduisent de nouvelles inconnues dans le problème global. Selon une deuxième technique, une fonction d'entrée artérielle globale est obtenue automatiquement à partir d'images de perfusion au moyen de techniques de traitement du signal telles que le Partitionnement de Données (Clustering en langue anglaise) ou l'Analyse en Composantes Indépendantes (Independent Component Analysis (ICA) en langue anglaise). Cette deuxième technique permet de s'affranchir de l'intervention humaine. Toutefois cette approche ne résout pas les problèmes de délai et de dispersion inhérents aux fonctions d'entrée artérielle globales. Selon une troisième approche, des fonctions d'entrée artérielle locales sont produites automatiquement à partir d'images de perfusion à l'aide de techniques de traitement du signal et de critères de sélection. On recherche par exemple la « meilleure » fonction dans un voisinage immédiat du voxel courant pour lequel on souhaite estimer les paramètres hémodynamiques. L'objet de cette troisième technique est d'obtenir au final des estimations moins biaisées et plus précises. On s'affranchit, du moins en partie, des problèmes de délai et de dispersion. Cependant, rien ne garantit a priori et a posteriori que les fonctions d'entrée artérielle locales ainsi obtenues sont des approximations pertinentes de la « véritable fonction locale » pour le voxel d'intérêt. Par exemple, cette « véritable » fonction pourrait ne pas se situer dans le voisinage considéré si celui-ci est trop petit. Ladite « véritable » fonction pourrait également être confondue avec une autre fonction d'entrée artérielle, si le voisinage considéré est trop grand. Par conséquent, même si les résultats peuvent être meilleurs qu'avec une fonction d'entrée artérielle globale, les incertitudes pesant sur les fonctions d'entrée artérielles locales (et a fortiori sur les paramètres hémodynamiques) demeurent en grande partie. Quelle que soit la technique jusqu'alors employée, la déconvolution de la courbe de concentration par la fonction d'entrée artérielle, que celle-ci soit globale ou locale (et par conséquent l'estimation des paramètres hémodynamiques) est effectuée en utilisant, par exemple, des méthodes de déconvolution non-paramétriques telles que la Décomposition en Valeurs Singulières (Singular Value Decomposition en langue anglaise), la Déconvolution de Hunt dans le domaine fréquentiel, la régularisation de Tikhonov, etc. Des méthodes, utilisant des modèles physiologiques paramétriques et semi-paramétriques Ca(t,Oa) pour la fonction d'entrée artérielle et R(t,OR) pour la fonction de répartition complémentaire du temps de transit dans le tissu, ont été toutefois proposées dans la littérature. L'utilisation connue de tels modèles n'a pas pour objet de supplanter l'opération de déconvolution mais de tenter d'améliorer certains résultats obtenus par les méthodes non paramétriques telles que décrites précédemment. Par exemple, un modèle physiologique paramétrique Ca(t,0a) de la fonction d'entrée artérielle locale, modèle « mono-Gamma » à une composante et quatre paramètres ao-1 _(t-t0)110 a t-t0 ) e flOa°F (a0 Oa =(a,a0,fl0,t0) I' ( ) fonction Gamma d'Euler a été introduit pour : - synthétiser des données de simulation pour l'étude 25 et la validation des méthodes de déconvolution non-paramétriques ; - réduire le bruit de mesure et lisser des données en remplaçant la ou les fonctions d'entrée artérielle expérimentales par un ajustement de ces modèles lisses auxdites courbes avant d'effectuer la déconvolution numérique et l'estimation des paramètres ; - pallier les problèmes liés à l'échantillonnage temporel des fonctions d'entrée artérielle en sur-5 échantillonnant artificiellement un ajustement de ces modèles auxdites courbes avant d'effectuer la déconvolution numérique ; - pallier les problèmes dus au phénomène de recirculation de l'agent de contraste dans le voxel 10 d'intérêt, lorsqu'elle recouvre temporellement le premier passage de l'agent de contraste. De même, des modèles paramétriques ou semiparamétriques R(t,eR) de la fonction de répartition 15 complémentaire du temps de transit dans le tissu, ont été proposés, parmi lesquels nous pouvons citer : - le modèle paramétrique « créneau » (box-shaped) : R(t,OR)=x[o 7T](t) t>0 OR=MTT XE ( ) fonction indicatrice - le modèle paramétrique « triangle » R (t,®R)=(1 2•tMTT)x[o,2•MTT](t) t>_0 OR =MTT 20 [ZE ( ) fonction indicatrice - le modèle paramétrique « monoexponentiel » (monoexponential) : _ t R(t,®R)=e MIT t>0 OR=MTT - le modèle paramétrique « Gamma Intégral » à deux 25 paramètres : +ao R(t,OR)= J h(r,a,/3)dr t>_0 t ta-le-tlt3 h(r,a,/3)= aI (a) a>0,/3>0 OR = (a, fl) I' ( ) fonction Gamma d'Euler Ces modèles ont été introduits pour : - synthétiser des données de simulation pour l'étude et la validation des méthodes de déconvolution non-paramétriques ; c'est par exemple le cas pour les modèles « créneau » et « triangle » qui représentent des cas extrêmes qui ne peuvent pas être rencontrés dans la nature ; - obtenir, réciproquement par ajustement de ces modèles aux données expérimentales, des fonctions de répartition complémentaire du temps de transit dans le tissu R(t) plus « physiologiques » que celles obtenues par les méthodes non-paramétriques (en particulier, des fonctions positives, décroissantes et valant un à l'origine, telles que définies par le modèle standard de perfusion) ; - obtenir analytiquement d'autres fonctions et paramètres d'intérêt dépendant de la fonction de répartition complémentaire, telle que la fonction de densité de probabilité h(t) du temps de transit dans le tissu (réponse impulsionnelle) donnée par h(t)=8(t)-dR(t) dt où 8(t) est la fonction généralisée de Dirac. L'utilisation de tels modèles a permis d'améliorer peu ou prou certaines estimations de paramètres. Toutefois, les inconvénients inhérents au fait de fixer la fonction d'entrée artérielle, que celle-ci soit globale ou locale, restent présents. L'invention permet de répondre à l'ensemble des inconvénients soulevés par les solutions connues. Elle offre des procédures rapides, objectives et reproductibles. L'invention ne porte pas sur une méthode de diagnostic en tant que telle car un procédé selon l'invention n'englobe pas toutes les étapes permettant un tel diagnostic. En particulier, l'invention ne concerne pas des phases d'investigation permettant un recueil de données réalisées par exemple au moyen d'un appareil d'imagerie par Résonance Magnétique Nucléaire. L'invention ne concerne pas en outre d'étapes d'interaction avec un corps humain ou animal permettant à un praticien de mettre en oeuvre un diagnostic. Toutefois, les avantages apportés par l'invention sont nombreux et significatifs. L'invention peut se révéler particulièrement appréciable en situation d'urgence clinique, offrant ainsi à un praticien une aide pour parfaire un diagnostic et prendre une décision thérapeutique adéquate. La mise en oeuvre de l'invention permet en outre à un chercheur de progresser dans la compréhension physiologique d'un organe d'un humain ou d'un animal. Parmi les nombreux avantages, nous pouvons mentionner: - la suppression de toute intervention humaine et/ou de mesures expérimentales additionnelles pour tenter de choisir manuellement ou automatiquement une fonction d'entrée artérielle globale ou des fonctions d'entrée artérielle locales ; - obtenir des estimations des paramètres hémodynamiques d'intérêt moins biaisées et plus précises en s'affranchissant de problèmes de déconvolution « mal-posés » et numériquement instables ainsi que par construction, des problèmes de délai et de dispersion inhérents aux fonctions d'entrée artérielle globales ou locales ; - fournir optionnellement des estimations des fonctions d'entrée artérielle locales ainsi que des fonctions de répartition complémentaire pour chaque voxel de tissu afin, par exemple, de comparer et de valider l'apport de l'invention par rapport à l'état de l'art ; - fournir optionnellement des intervalles de confiance sur l'estimation de chacun des paramètres et même des paris sur lesdits intervalles de confiance ; - fournir optionnellement une mesure quantitative et objective de l'adéquation du modèle de perfusion global aux données expérimentales et, par conséquent, permettre la comparaison et la sélection de modèles les plus adéquats. A cette fin, il est prévu un procédé pour estimer des paramètres hémodynamiques de perfusion d'un volume élémentaire - dit voxel - d'un organe, ledit procédé étant mis en oeuvre par une unité de traitement d'un système d'analyse d'imagerie de perfusion, et comportant une étape pour estimer des paramètres hémodynamiques à partir de signaux de perfusion S(t), caractérisé en ce que l'étape pour estimer lesdits paramètres hémodynamiques consiste en une estimation conjointe des paramètres O d'un modèle global de perfusion contraint comportant : - une première relation entre le signal de perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans ledit voxel au cours du temps t ; une deuxième relation entre la concentration C(t) d'un agent de contraste, un paramètre physiologique B, un modèle paramétrique ou semi-paramétrique Ca(t,Oa) d'une fonction d'entrée artérielle et un modèle paramétrique ou semi-paramétrique R(t,OR) d'une fonction de répartition complémentaire du temps de transit dans le voxel, Oa et OR étant respectivement les paramètres desdits modèles Ca(t,Oa) et R(t,OR) ; une troisième relation Y(Oa)=0 entre les paramètres Oa du modèle de la fonction d'entrée artérielle Ca(t,Oa). Selon un premier exemple d'application de l'invention, un tel procédé peut porter sur : un signal de perfusion S(t) préalablement obtenu partir d'images numériques délivrées par un appareil d'imagerie de perfusion par résonance magnétique ; un modèle global de perfusion contraint dont : la première relation peut être exprimée par -k•rE•c(r) S(t) = Soe où s0 est l'intensité moyenne du signal avant l'arrivée de l'agent de contraste dans le voxel, TE est un temps d'écho et k une constante non nulle ; la deuxième relation peut être exprimée par C(t)=77.0•Ca(t,Oa)1 R(t,OR) où désigne le produit de convolution, 77 est une constante 30 non nulle et le paramètre physiologique 0 est le flux sanguin BF dans le voxel d'intérêt. Selon un deuxième exemple d'application de l'invention , un tel procédé peut porter sur :25 un signal de perfusion S(t) préalablement obtenu à partir d'images numériques délivrées par un appareil d'imagerie de perfusion par tomodensitométrie ; un modèle global de perfusion contraint dont : - la première relation peut être une relation de proportionnalité S(t)=a.C(t) où a est une constante non nulle; - la deuxième relation peut s'exprimer par C(t)=17.8-Ca(t,Oa)r3R(t,OR) où r:+ désigne le produit de convolution, 17 est une constante non nulle et le paramètre physiologique 8 est le flux sanguin BF dans le voxel. Selon un premier mode de réalisation, un procédé conforme à l'invention peut être mis en oeuvre par itérations successives pour une pluralité de voxels considérés. A titre d'exemple, l'invention prévoit que l'estimation conjointe des paramètres O du modèle global de perfusion contraint mise en oeuvre par l'unité de traitement peut être réalisée au moyen de la méthode de Bayes ou des méthodes du Maximum de Vraisemblance ou encore des méthodes des Moindres Carrés non linéaires. Selon un mode de réalisation préféré, le procédé peut comporter une étape pour quantifier l'adéquation du modèle global de perfusion contraint à des données de perfusion D=[S(tl),...,S(tN)] expérimentales par un calcul d'une probabilité p(DIM)= jp(OIM)p(DIO,M)dO de ces 0 données sachant ledit modèle. Il peut être en outre prévu que le procédé selon l'invention puisse comporter une étape préalable pour choisir un modèle global de perfusion contraint parmi une pluralité. Ainsi il est possible qu'un tel procédé puisse être mis en oeuvre itérativement par l'unité de traitement pour chaque modèle global de perfusion contraint connu de ladite unité de traitement. Selon le mode de réalisation préféré précédemment décrit, il peut être alors prévu que sont délivrés les paramètres estimés selon le modèle de perfusion dont la probabilité sachant ledit modèle est la plus grande. Outre, une estimation de paramètres hémodynamiques, l'invention prévoit que le procédé puisse comporter une étape pour calculer une information complémentaire sous la forme : d'une fonction d'entrée artérielle estimée Ca(té) a ; d'une fonction de répartition complémentaire estimée R(t,OR) ; d'un intervalle de confiance paramètre du modèle global contraint ; ou encore d'un taux de pari paramètre du modèle global contraint. associé à un de perfusion associé à un de perfusion L'invention prévoit un mode de réalisation du procédé pour qu'il puisse comporter une étape pour délivrer des paramètres hémodynamiques estimés, voire toutes autres informations complémentaires qui leur seraient associées, à une interface homme-machine apte à restituer à un utilisateur lesdits paramètres estimés. L'invention concerne en outre une unité de traitement 5 comportant des moyens de mémorisation, des moyens pour communiquer avec le monde extérieur et des moyens de traitement, de sorte que : - les moyens pour communiquer sont aptes à recevoir du monde extérieur un signal S(t) obtenu par 10 imagerie de perfusion ; - les moyens de mémorisation comportent un modèle global de perfusion contraint, ledit modèle comportant : - une première relation entre le signal de 15 perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans ledit voxel au cours du temps t ; - une deuxième relation entre un paramètre physiologique 8, la concentration C(t), un 20 modèle paramétrique ou semi-paramétrique Ca(t,Oa) d'une fonction d'entrée artérielle et un modèle paramétrique ou semi-paramétrique R(t,OR) d'une fonction de répartition complémentaire du temps de transit dans le 25 voxel, Oa et OR étant respectivement les paramètres desdits modèles ; - une troisième relation YOOa)=0 entre les paramètres Oa du modèle de la fonction d'entrée artérielle Ca(t,0a) ; 30 - les moyens de traitement sont adaptés pour mettre en oeuvre un procédé pour estimer des paramètres hémodynamiques de perfusion conforme à l'invention. L'invention concerne en variante une unité de traitement dont les moyens pour communiquer sont aptes à recevoir du monde extérieur des images numériques de perfusion à partir desquelles les moyens de traitement déterminent un signal de perfusion S(t) et mettent en œuvre un procédé pour estimer des paramètres hémodynamiques de perfusion conforme à l'invention. L'invention prévoit en outre que les moyens pour communiquer d'une unité de traitement conforme à l'invention, peuvent délivrer des paramètres estimés selon un format approprié à une interface homme-machine apte à restituer à un utilisateur lesdits paramètres estimés. L'invention concerne également tout système d'analyse d'imagerie de perfusion comportant une unité de traitement telle que décrite précédemment et une interface homme-machine apte à restituer à un utilisateur des paramètres estimés selon un procédé conforme à l'invention et mis en œuvre par ladite unité de traitement. D'autres caractéristiques et avantages apparaîtront plus clairement à la lecture de la description qui suit et à l'examen des figures qui l'accompagnent parmi lesquelles : les figures 1 et 2 présentent deux variantes de 30 réalisation d'un système d'analyse d'imagerie de perfusion ; Les figures 3 et 4 présentent respectivement une image de perfusion, obtenue par un appareil d'imagerie par Résonance Magnétique Nucléaire, 35 d'une tranche d'un cerveau humain avant l'injection d'un agent de contraste et pendant la circulation de celui-ci dans les tissus dudit cerveau ; les figures 5a et 5b présentent un signal de perfusion S(t) par Résonance Magnétique Nucléaire typique relatif à un voxel d'un cerveau humain ; la figure 6 présente une courbe de concentration CO') typique d'un agent de contraste circulant au sein d'un voxel d'un cerveau humain; la figure 7 présente une fonction d'entrée artérielle Ca(t) typique ; la figure 8 présente un procédé conforme à l'invention ; les figures 9 à 14 présentent respectivement une carte relative à un paramètre hémodynamique estimé conformément à l'invention ; les figures 15 et 16 concernent respectivement des écart-types des flux sanguins et du temps de transit moyen estimés conformément à l'invention ; la figure 17 concerne une carte relative à une probabilité qu'un flux sanguin cérébral soit dans un intervalle de confiance. La figure 1 permet de présenter un système d'analyse d'images de perfusion. Un appareil 1 d'imagerie par Résonance Magnétique Nucléaire ou par Tomodensitométrie est commandé à l'aide d'une console 2. Un utilisateur peut ainsi choisir des paramètres 11 pour piloter l'appareil 1. A partir d'informations 10 produites par l'appareil 1, on obtient une pluralité de séquences d'images numériques 12 d'une partie d'un corps d'un humain ou d'un animal. A titre d'exemple préféré, nous illustrerons les solutions issues de l'art antérieur ainsi que l'invention à l'aide d'images numériques issues de l'observation d'un cerveau humain. D'autres organes pourraient aussi être considérés. Les séquences d'images 12 peuvent optionnellement être stockées au sein d'un serveur 3 et constituer un dossier médical 13 d'un patient. Un tel dossier 13 peut comprendre des images de différents types, telles que des images de perfusion ou de diffusion. Les séquences d'images 12 sont analysées au moyen d'une unité de traitement 4 dédiée. Ladite unité de traitement comporte des moyens pour communiquer avec le monde extérieur pour recueillir les images. Lesdits moyens pour communiquer permettent en outre à l'unité de traitement de délivrer in fine à un praticien 6 ou à un chercheur, une estimation des paramètres hémodynamiques à partir des images de perfusion 12, au moyen d'une interface homme-machine adaptée 5. L'utilisateur 6 du système d'analyse peut ainsi confirmer ou infirmer un diagnostic, décider une action thérapeutique qu'il jugera adéquate, approfondir des travaux de recherche... Optionnellement, cet utilisateur peut paramétrer le fonctionnement de l'unité de traitement 4 au moyen de paramètres 16. Par exemple, il peut ainsi définir des seuils d'affichage ou choisir les paramètres estimés qu'il souhaite visualiser. La figure 2 illustre une variante de réalisation d'un système d'analyse pour lequel une unité de prétraitement 7 analyse des séquences d'images 12 pour en déduire des signaux de perfusion 15 par voxel. L'unité de traitement 4 chargée d'estimer les paramètres hémodynamiques 14 est ainsi déchargée de cette action et met en œuvre un procédé d'estimation à partir de signaux de perfusion 15 réceptionnés par ses moyens pour communiquer avec le monde extérieur. La figure 3 illustre un exemple d'image typique 12 d'une tranche de 5 millimètres d'épaisseur d'un cerveau humain. Cette image est obtenue par Résonance Magnétique Nucléaire. A l'aide de cette technique, on peut obtenir, pour chaque tranche, une matrice de 128 x 128 voxels dont les dimensions sont de 1,5 x 1,5 x 5 millimètres. A l'aide d'une interpolation bilinéaire on peut produire une image à plat de 458 x 458 pixels telle que l'image 20. La figure 4 illustre une image 20 similaire à celle présentée en liaison avec la figure 3. Toutefois cette image est obtenue après une injection d'un agent de contraste. Cette image est un exemple d'image de perfusion typique d'un cerveau. Les artères apparaissent ainsi clairement contrairement à la même image décrite en figure 3. Selon des techniques connues, il est possible de choisir une ou plusieurs fonctions d'entrée artérielle 21 dans l'hémisphère contralatérale à l'hémisphère pathologique pour estimer des paramètres hémodynamiques. La figure 5b permet d'illustrer un exemple de signal de perfusion S(t) par Résonance Magnétique Nucléaire tel que les signaux 15 délivrés par l'unité de prétraitement 7 décrite en liaison avec la figure 2. Le signal de perfusion est ainsi représentatif de l'évolution d'un voxel au cours du temps t à la suite d'une injection d'un agent de contraste. A titre d'exemple, la figure 5b décrit un tel signal sur une durée de 50 secondes. L'axe des ordonnées décrit l'intensité du signal dont l'unité est arbitraire. Pour obtenir un tel signal, l'unité de traitement 4 selon la figure 1 (ou en variante l'unité de prétraitement 7 selon la figure 2), analyse une séquence de n images de perfusion par Résonance Magnétique Nucléaire Il, I2, , Ii, In à des instants tl, t2 , ..., tà , ..., tn comme le décrit, à titre d'exemple, la figure 5a. Pour un voxel donné, par exemple pour le voxel V, on détermine un signal de perfusion S(t) représentatif de l'évolution du voxel au cours du temps t à la suite d'une injection d'un agent de contraste. La figure 6 présente une courbe de concentration déduite d'un signal de perfusion tel que celui-ci décrit en figure 5b. Comme déjà évoqué précédemment, il existe une relation entre un signal de perfusion et une courbe de concentration associée. Ainsi en Imagerie de Perfusion par Résonance Magnétique Nucléaire, il existe une -k.TE.C(t) relation exponentielle S(t) = Sp.e où SD est l'intensité moyenne du signal avant l'arrivée de l'agent de contraste, TE est le temps d'écho (echo time) et k est une constante dépendant de la relation entre la susceptibilité paramagnétique et la concentration de l'agent de contraste dans le tissu. La valeur de la constante k pour chaque voxel étant inconnue, celle-ci est fixée à une valeur arbitraire pour tous les voxels d'intérêt. La figure 6 permet ainsi de visualiser au cours du temps, l'évolution de la concentration d'un agent de contraste au sein d'un voxel. On note un pic de forte amplitude lors du premier passage (first pass, en langue anglaise) de l'agent de contraste dans le voxel suivi de pics d'amplitudes plus faibles lié à un phénomène de recirculation (second pass, en langue anglaise) dudit agent de contraste. La figure 7 illustre quant à elle, une fonction d'entrée artérielle typique Ca(t) représentative de la circulation d'un agent de contraste au sein d'un voxel artériel tel que le voxel 21 présenté en liaison avec la figure 4. La figure 7 permet de constater notamment que le phénomène de recirculation après un premier passage de l'agent de contraste est très faible.35 Pour mettre en oeuvre l'invention, on introduit un modèle physiologique paramétrique ou semi-paramétrique Ca(t,Oa) d'une fonction d'entrée artérielle. Un tel modèle est stocké par des moyens de mémorisation ou programmé dans l'unité de traitement 4 d'un système d'analyse tel que décrit en liaison avec les figures 1 et 2. Selon un mode de réalisation préféré de l'invention, 10 ce modèle est un modèle dit « tri-Gamma » à douze paramètres défini par : )a0-1 -(t-t0)')30 )al-1 -(t-tl)//31 )a2-1 -(t-t2)/,a2 a tùt0 ((e b tùt1 ((e l e tùt2 (e ll ,oaoF(ao) + filalr ) la1 + fl2a2F(a2) Oa =(a,b,c,a0,fip,tp,al,/3i,tl,a2,/32't2) F ( ) : fonction Gamma dBuler 15 En outre, on introduit un modèle physiologique paramétrique ou semi-paramétrique R(t,OR) de la fonction de répartition complémentaire du temps de transit dans le tissu. Au même titre que le modèle de fonction d'entrée 20 artérielle, un tel modèle est stocké ou programmé dans l'unité de traitement 4 d'un système d'analyse tel que décrit en liaison avec les figures 1 ou 2. A titre d'exemple, ce modèle de fonction de 25 répartition complémentaire du temps de transit peut être un modèle dit « Gamma Intégral Corrigé » à deux paramètres : t R(t,®R) = H(t)- Jh(r,MTT,/3)dr 0 MTT_1 t e-tl/3 h(r,MTT,/3)= MTT MTT>0,/3>0 f r MTT OR = (MTT, /3) H( ) : fonction généralisée créneau de Heaviside t ( ) : fonction Gamma d'Euler Avec la seule définition des modèles paramétriques ou semi-paramétriques pour la fonction d'entrée artérielle et pour la fonction de répartition complémentaire, le problème de l'estimation des paramètres de perfusion reste sous-déterminé. On peut en effet écrire un modèle standard de perfusion (en imagerie par Résonance Magnétique Nucléaire) de la manière équivalente suivante : S t = S0e-k•TE•C(t) 2 étant une constante non nulle. Ainsi, si BF et Ca(t,Oa) sont deux estimations permettant d'ajuster le modèle théorique de perfusion aux Ça (t,Oa ) données expérimentales, alors 2BF et 2 sont deux autres estimations ajustant tout aussi bien le modèle aux données expérimentales pour tout 2#0. Autrement dit, BF et Ca(t,Oa) sont déterminés seulement à une constante multiplicative près sous un modèle standard de perfusion. C'est-à-dire qu'ils ne sont pas du tout déterminés à moins que précisément, on ne C (t) _ (2.BF) Ca t,0 R(t,eR) M: fixe préalablement, globalement ou localement, la fonction d'entrée artérielle par des moyens tels que ceux décrits précédemment ou qu'on ajoute une contrainte supplémentaire. L'invention prévoit donc d'introduire une contrainte supplémentaire lP(Oa) sur le modèle de la fonction d'entrée artérielle Ca(t,Oa). Cette troisième relation peut exprimer par exemple la conservation d'une quantité physique sur toutes les fonctions d'entrée artérielle d'intérêt, telle que la masse totale de l'agent de contraste circulant dans un voxel artériel au cours du temps. De manière particulièrement avantageuse, cette contrainte est indépendante du modèle Ca(t,Oa) considéré. Cette contrainte peut s'exprimer, de manière générale, par une relation entre les paramètres Oa du modèle de la fonction d'entrée artérielle telle que lP(Oa)=0 Cette relation supplémentaire sous la forme d'une contrainte sur le modèle de la fonction d'entrée artérielle Ca(t,Oa) permet de fixer, une fois pour toute et pour tous les voxels considérés, la constante 2. Cette troisième relation est stockée ou programmée dans l'unité de traitement 4 d'un système d'analyse tel que décrit en liaison avec les figures 1 ou 2. L'ensemble des trois relations constituent un modèle global de perfusion contraint. Le problème de l'estimation conjointe des paramètres hémodynamiques, en particulier de BF et Oa, devient alors bien déterminé et bien posé. Selon un premier mode de réalisation de l'invention, cette contrainte W(Oa)=O peut s'exprimer comme : +0o f Ca(t,Oa)dt=CO où Co est une constante arbitraire 0 non nulle identique pour tous les voxels d'intérêt. 5 Dans le cas du modèle « tri-Gamma » pour une fonction d'entrée artérielle de paramètres Oa = (a,b,c,a0, f30,t0,a1,l1,t1,a2, f32,t2) , cette contrainte sur les paramètres peut être exprimée comme a+b+c=Co. 10 On a donc la relation : 1I'(Oa)=a+b+c-00. Sous cette contrainte, on peut ainsi éliminer l'un des trois paramètres, par exemple c en posant c=C0-a-b, et se ramener à un modèle Ca(t,Oa) contraint à onze 15 paramètres libres, tel que : ao-1 -(t-to)1io al-1 it-qA a2-1 a(t-to) e +b(t-tl) e +~Co-a-b)(t-t2) Platr(al) ,2a2r(a2 Oa =(a,b,a0,p0,t0,a1,131,t1,a2,fi2't2) F ( ) fonction Gamma d'Euler Selon un deuxième mode de réalisation, la contrainte +0o peut s'exprimer comme f Cca(t,Oca)dt=C0 où Cca(t,@ca) est 0 20 une composante du modèle Ca(t,0a) modélisant la concentration de l'agent de contraste lors de son premier passage (first pass). Dans le cas du modèle « tri-Gamma », il s'agit simplement de la première de ses trois composantes : f0a°r (ao de sorte que la contrainte devient a = C0 soit W(Oa)=a-Co . On peut donc à nouveau se ramener à un modèle 5 Ca(t,0a) contraint à onze paramètres libres tel que : rca(t,®a)= Co (- o )aoe-(t-to~30 + b (t - t, )a11 e-(t-trl+ c (t - t2 )a2 i e-(t-t2)l/32 ~oa-r(ao) ,f31a1T (a1) a2r(a2) L'invention permet ainsi de se ramener à un problème d'estimation conjointe des paramètres O={0a,OR,BF,S0,0B} 10 où OB est un vecteur de paramètres caractérisant les bruits de mesure sur le signal de perfusion S(t), sans avoir à fournir ou à estimer préalablement une quelconque fonction d'entrée d'artérielle. 15 L'invention se fonde ainsi sur un modèle global de perfusion contraint comprenant : une première relation entre un signal de perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans un voxel de tissu au cours 20 du temps t ; - une deuxième relation entre un paramètre physiologique 6, tel que par exemple le flux sanguin BF, ladite concentration C(t), un modèle paramétrique ou semi-paramétrique Ca(t,0a) d'une 25 fonction d'entrée artérielle et un modèle paramétrique ou semi-paramétrique R(t,OR) d'une ao-1 -(t-to)/,60 a(t-t0) e Cca(t,Oca)- floaortaol Oca ={a,ao,fl0,to} Oa F ( ) fonction Gamma d'Euler fonction de répartition complémentaire du temps de transit dans le voxel, Oa et OR étant respectivement les paramètres desdits modèles ; - une troisième relation LP(Oa)=0 entre les paramètres 0a du modèle de la fonction d'entrée artérielle Ça (t, Oa) On pourra dénommer les modèles respectivement d'une fonction d'entrée artérielle Ca(t,0a) et d'une fonction de répartition complémentaire du temps de transit R(t,OR) comme étant des sous-modèles d'un modèle global de perfusion contraint par ladite troisième relation 1Y(Da)=0 entre les paramètres Oa dudit modèle de la fonction d'entrée artérielle Ca(t,0a). A titre d'exemple en imagerie par Résonance Magnétique Nucléaire, un tel modèle global de perfusion contraint M peut être exprimé sous une forme telle que : = -k•TE•C(t) S (t) S0e M: C(t)=t7.0•Ca(t,Oa)O R(t,OR) 1Y(Oa)=0 où -3 désigne le produit de convolution, r7 est une constante non nulle fixée arbitrairement et 0 est le flux sanguin BF circulant dans le voxel considéré. De même, en imagerie par tomodensitométrie un modèle global de perfusion contraint M peut être exprimé sous une forme telle que : {s (t) = a.0 (t) M: C(t)=t7.O•Ca(t,Oa LI'(Oa)= R où ® désigne le produit de convolution, rJ est une constante non nulle fixée arbitrairement, 8 est le flux sanguin BFcirculant dans le voxel considéré et a une constante non nulle. En liaison avec la figure 8, un procédé pour estimer des paramètres hémodynamiques, conforme à l'invention, comporte une étape 51 pour estimer conjointement les paramètres O=}Oa,OR,B,So,OB} d'un modèle M à partir de signaux de perfusion S(t). Un tel procédé est mis en oeuvre par une unité de traitement 4, telle que décrite en figures 1 ou 2, adaptée pour pouvoir notamment mettre en oeuvre des techniques mathématiques et numériques pour résoudre un problème d'estimation conjointe de paramètres. A titre d'exemples, une telle unité de traitement 4 est apte à mettre en oeuvre les Méthodes du Maximum de Vraisemblance, les Méthodes des Moindres Carrés Non- linéaires ou la Méthode de Bayes. Afin d'illustrer un mode de réalisation préféré de l'invention, nous allons considérer que : - l'unité de traitement 4 met en oeuvre la méthode de 25 Bayes ; - le modèle de la Fonction d'Entrée Artérielle, stocké ou programmé dans l'unité de traitement 4, est le modèle « tri-Gamma » ; - le modèle la Fonction de Répartition Complémentaire, 30 stocké ou programmé dans l'unité de traitement 4, est le modèle « Gamma Intégral Corrigé » ; la contrainte, stockée ou programmée dans l'unité de traitement 4, est Y(Oa)=Co-a-b-c. Pour une application à l'imagerie de perfusion par Résonance Magnétique Nucléaire, le modèle Ca(t,Oa) contraint s'écrit donc : Ca t,Oa= a(t-t0 )a°-1 e-(t-t°)/fl° + b(t-t1~a1-1 e-(t-tly fil + (Co -a-b)(t-t2)a2_1 e-(t-t2)/f32 floc°r(ao) fer(al) r2a2r(a2) = (a,b,ap,/30,t0,a1, j31,t1,a2,f32,t2 ) F ( ) fonction Gamma d'Euler Soit ti,i=1,N des instants d'échantillonnage et S(ti),i=1,N l'intensité du signal de perfusion mesurée en ces instants. Les données expérimentales sont donc 10 D=[S(tl),...,S(tN)] . Si l'on suppose un bruit de mesure additif, blanc, Gaussien, stationnaire et d'écart-type 6 = OB, alors la fonction de vraisemblance ou distribution directe des 15 données de perfusion D sachant tous les paramètres s'écrit : N ir N p(DIO,M)=(2 r) 2 ff-N - e 20 C(t)=r1.O Ca(t,Oa)R(t,OR) est calculée numériquement après discrétisation temporelle soit, de préférence, analytiquement en utilisant éventuellement une approximation pour le produit de convolution de deux distributions de probabilité F. 25 Etant donnée une distribution conjointe a priori de tous les paramètres sachant le modèle S(tt)-Soe k•TE.C(tz)'2 26-2 p(OI M)=p(Oa,OR,0,SO,6I M) représentant l'information dont on dispose avant l'expérience sur ces paramètres, on obtient une distribution conjointe a posteriori par la règle de Bayes : r l l p( I jP(elm)p(Dio,m)de p(OIM)p(DIO,M)"I (" OD,M- acpOMpDO,M A partir de cette dist ribution conjointe, le procédé selon l'invention permet d'obtenir 52 la distribution marginale a posteriori de chacun des paramètres 10 hémodynamiques d'intérêt en marginalisant tous les autres. A partir de cette distribution marginale, on peut enfin obtenir un estimateur du paramètre, par exemple l'estimateur de Bayes sous fonction de coût quadratique en prenant l'espérance mathématique de cette distribution 15 ou encore la valeur la plus probable du paramètre (estimateur du maximum a posteriori). Ainsi l'invention permet par exemple de pouvoir estimer pour tout voxel d'intérêt, une estimation du flux sanguin CBF d'un cerveau ainsi qu'une estimation du Temps de Transit Moyen 20 AM'. De manière optionnelle, un procédé conforme à l'invention permet d'obtenir 52a des informations complémentaires sous la forme d'intervalles de confiance 25 sur des paramètres estimés et même des paris 52b sur lesdits intervalles de confiance. Ces informations complémentaires peuvent être particulièrement utiles, pour un utilisateur tel qu'un praticien, et se révéler être une aide supplémentaire dans l'élaboration d'un 30 diagnostic. A l'aide des techniques de l'art antérieur à l'invention, un praticien n'a aucune idée de la précision des estimations dont il dispose pour chaque voxel d'intérêt.5 Par exemple, pour le Flux Cérébral Sanguin CBF, l'invention permet d'obtenir une distribution marginale a posteriori : p(CBFI D,M) = f p(OID,M)d (O\CBF) _ O\CBF f f f f p(OQ,OR,CBF,Sp,oI D,M)do-dSodORdOa ®a %R so u puis un Estimateur de Bayes CBF de CBF sous fonction de _ )=(CBF-CBF _ 12 coût qua dratique L(CBF-CBFQ Q CBFQ = f CBF. p (CBFI D, M) dCBF CBF ou la valeur la plus probable du paramètre : CBFP = arg max p (CBFI D, M) CBF L'unité 4 peut aussi déduire optionnellement : - la variance de l'estimateur de Bayes sous fonction de coût quadratique : (CBF-CBF 2 d'CBF =CQ .p(CBFID,M)dCBF - un intervalle de confiance, par exemple l'intervalle de confiance à :Ela: [CBF - a-CBF , CBF + &CBF ] 20 Selon un mode de réalisation, l'unité de traitement 4 mettant en œuvre un procédé selon l'invention, permet de calculer 52b une probabilité/taux de pari que le paramètre soit dans l'intervalle de confiance, par exemple : CBF+&CBF p (CBF [CBF - 6_cBF, CBF +6_CBF ]) = f p (CBFI D, M) dCBF 25 CBFù&cBF 10 15 En variante, l'invention prévoit que le procédé permette d'obtenir 53 une information complémentaire correspondant à la distribution marginale conjointe a posteriori des paramètres de la fonction d'entrée artérielle locale pour un volume de tissu considéré, telle que : p(OaID,M)= f p(OID,M)d(O\Oa)= o\ea f f f f p(Oa,OR,CBF,S0,6l D,M)dadSodCBFd®R OR CBF So 0 L'invention permet ainsi obtenir par exemple une 10 fonction d'entrée artérielle « moyenne » Ca(t,0a) sous fonction de coût quadratique LIOa- OaI= O -O Q a a de 2 paramètres conjoints : OaQ = f Oa • p(OaI D,M)dda Oa ou encore une fonction d'entrée artérielle « la plus 15 probable » de paramètres : --P Oa =argmaxp(OaID,M) Oa De la même manière, l'invention prévoit une variante de réalisation pour laquelle un procédé conforme à 20 l'invention comporte une étape 54 pour calculer une fonction de répartition complémentaire « moyenne » ou « la plus probable » R(t,OR). On peut ainsi en outre, obtenir également une information complémentaire sous la forme d'estimations de 25 la concentration de l'agent de contraste par C(t,0). Selon un mode de réalisation, l'invention prévoit que ces fonctions d'entrée artérielle ou fonctions de répartition complémentaire estimées, « moyennes » ou « plus probables », puissent être comparées à celle(s) obtenue(s) par les approches antérieures à l'invention ou encore à celles obtenues à partir de modèles de fonction d'entrée artérielle ou de fonction de répartition complémentaire distincts dans le but de permettre la sélection des modèles les plus pertinents. En effet, il peut être intéressant de comparer les fonctions d'entrée artérielle estimées à la véritable fonction d'entrée artérielle lorsque celle-ci est connue et tâcher de minimiser leur écart en introduisant des modèles plus appropriés, si besoin est. En liaison avec la figure 8, l'invention prévoit en outre la possibilité qu'un procédé conforme à l'invention puisse comporter une ou plusieurs étapes complémentaires pour obtenir 53a des intervalles de confiance voire pour obtenir 53b des paris sur ces intervalles pour les paramètres des fonctions d'entrée artérielle. Il peut en être de même, pour obtenir 54a des intervalles de confiance voire pour obtenir 54b des paris sur ces intervalles pour les paramètres des fonctions répartition complémentaire. Les dimensions des modèles paramétriques et semiparamétriques typiques pour les fonctions d'entrée artérielle Ca(t,Oa) et les fonctions de répartition complémentaire R(t,OR) sont suffisamment petites pour que puissent s'appliquer des méthodes d'échantillonnage ou d'approximation de la distribution conjointe a posteriori p(OID,M) et des distributions marginales a posteriori telles que les méthodes de Monte-Carlo à Chaînes de Markov ou les Méthodes de Bayes Variationnelles. L'invention permet d'obtenir dès lors des approximations des estimations souhaitées pour les paramètres d'intérêt. Selon un mode de réalisation particulier, l'invention prévoit que le procédé, décrit en liaison avec la figure 8 et mis en oeuvre par une unité de traitement 4 adaptée en conséquence, puisse quantifier l'adéquation du modèle global de perfusion contraint aux données expérimentales D par le calcul 55 de la probabilité de ces données sachant le modèle p(DI M) = fp(OIM) p(DIO,M)dO, 0={Oa,OR,0,S0,0g} . Le calcul 55 de cette quantité n'était pas concevable avec les techniques de l'art antérieur à l'invention puisque selon ces dernières, on ne dispose pas de modèles globaux de signaux bien déterminés. L'invention prévoit que cette quantité puisse permettre de comparer de manière objective plusieurs modèles de perfusion globaux contraints obtenus à partir de différents sous-modèles pour la fonction d'entrée d'artérielle Ca(t,Oa) et/ou de différents sous-modèles pour la fonction de répartition complémentaire R(t,@R) et/ou de différentes contraintes supplémentaires TODa)= 0 tels que décrits précédemment. Ainsi, l'invention prévoit que le procédé puisse être mis en oeuvre par l'unité de traitement par itérations successives en utilisant des modèles globaux de perfusion contraints distincts. L'unité de traitement peut conserver le résultat du calcul 55 pour chaque modèle utilisé. L'invention permet de sélectionner au final le meilleur modèle global de perfusion contraint. Ce modèle est celui pour lequel la probabilité des signaux de perfusion sachant le modèle est la plus grande en moyenne sur l'ensemble des voxels d'intérêt. En variante, l'invention prévoit que le procédé puisse comporter des étapes 56a, 56b et 56c pour choisir lors d'une itération respectivement un sous-modèle Ca(t,ea), un sous-modèle R(t,OR), une contrainte `P(Oa) = 0 parmi une pluralité et ainsi ajuster 57 un modèle global de perfusion contraint. Cette sélection peut être réalisée manuellement par un utilisateur 6 au moyen de paramètres 16 ou automatiquement par le procédé lui-même. Par un processus itératif et essais successifs, il devient possible grâce à l'invention de progresser dans la modélisation et la compréhension des phénomènes de perfusion et à terme d'identifier et de valider les meilleurs modèles en fonction des patients, des pathologies, des types de tissus, etc. A titre d'exemple d'application de l'invention, nous pouvons citer les principales étapes de mise en oeuvre de l'invention au moyen d'un système d'analyse d'imagerie de perfusion adapté tel que celui décrit en figures 1 ou 2 : - ouverture d'un dossier patient ou prise en compte de séquences d'images par l'unité de traitement 4 (ou de prétraitement 7) pour sélectionner des séquences d'images d'intérêt - en particulier, sélection des images Il à In de perfusion au cours du temps à partir desquelles sont obtenus les signaux de perfusion S(t) pour chaque voxel, tel qu'illustré en figure 5a ; - prévisualisation au moyen d'une interface homme-machine 5 des images pour permettre à un utilisateur 6 d'identifier des tranches ou des zones d'intérêt ; - estimation conjointe par l'unité de traitement 4 des paramètres hémodynamiques 14, tels que CBF ou MTT, pour un organe tel que le cerveau humain et délivrance desdits paramètres 14 à l'interface homme-machine 5 pour que celle-ci les présente in fine sous la forme de cartes où l'intensité ou la couleur de chaque pixel dépend de la valeur calculée, par exemple de manière linéaire ; cette étape peut inclure également la délivrance d'autres quantités estimées tels que des intervalles de confiance, des paris, ou autres informations complémentaires associées aux paramètres du modèle de perfusion global contraint ; - affichage desdites cartes pour en restituer la teneur à l'utilisateur 6 ; - seuillage et filtrage sous l'impulsion 16 de l'utilisateur 6 pour ne conserver qu'une région digne d'intérêt et/ou supprimer certaines aberrations ; - sélection assistée de ladite zone pathologique d'intérêt par l'utilisateur, caractérisée par une anomalie de la distribution d'un ou de plusieurs paramètres hémodynamiques 14 ; - estimation, par l'unité de traitement 4, du volume de la zone de tissus anormalement perfusée et éventuellement de certaines quantités telles que le ratio des volumes des zones lésées et anormalement perfusées, sur lesquelles un praticien pourra peaufiner son diagnostic et sa prise de décision thérapeutique (thrombolyse intraveineuse afin de résorber un caillot sanguin par exemple). Les figures 9 à 17 permettent d'illustrer un mode d'affichage sous la forme de cartes, des paramètres hémodynamiques 14 estimés conformément à l'invention voire des écarts-types ou probabilité qui leur sont associés. Ainsi pour un cerveau humain analysé à l'aide d'imagerie par Résonance Magnétique Nucléaire, la figure 9 permet de visualiser une estimation des volumes sanguins cérébraux CBV. Une telle carte permet de mettre en évidence une zone ischémique probable 80. En effet il est possible de constater à l'aide d'une interface 6 adaptée, une nette augmentation du paramètre CBV dans le territoire de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale. Une vasodilatation consécutive à l'ischémie peut être révélée par une lecture de la carte telle qu'illustrée en figure 9. La figure 10 permet d'illustrer une carte relative à l'estimation des flux sanguins cérébraux en cas d'ischémie cérébrale. On peut constater en analysant la carte, une légère diminution du paramètre CBF dans le territoire de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale consécutive à l'ischémie. La figure 11 permet d'illustrer une carte relative à l'estimation des temps de transit moyen MTT. On peut constater en analysant la carte, une nette augmentation des MTT dans le territoire 80 de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale consécutive à l'ischémie. La figure 12 permet de décrire une carte relative à l'estimation du paramètre a d'un modèle de fonction d'entrée artérielle « tri-Gamma » à 12 paramètres conforme à l'invention. On peut constater une nette diminution de a dans le territoire 80 de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale. Cette information permet de constater une diminution relative de la quantité d'agent de contraste lors de la circulation par rapport à la quantité pendant la recirculation. Cette information n'est pas disponible avec les techniques selon l'art antérieur à l'invention. La figure 13 permet de décrire une carte relative à l'estimation du paramètre )62 du modèle de fonction d'entrée artérielle « tri-Gamma » à 12 paramètres conforme à l'invention. On peut constater en analysant la carte, une nette augmentation de A2 dans le territoire 80 de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale. Cette information n'est pas disponible avec les techniques selon l'art antérieur à l'invention. La figure 14 permet de décrire une carte relative à l'estimation du paramètre tl du modèle de fonction d'entrée artérielle « tri-Gamma » à 12 paramètres conforme à l'invention. On peut constater en analysant la carte, une nette augmentation de Il dans le territoire de l'artère cérébrale postérieure droite par rapport à l'hémisphère contralatérale. Cette information n'est pas disponible avec les techniques selon l'art antérieur à l'invention. La figure 15 permet de décrire une carte relative à l'estimation de l'écart-type des flux sanguins cérébraux estimés o.cBF à partir du modèle de fonction d'entrée artérielle « tri-Gamma » à douze paramètres conforme à l'invention. On peut constater en analysant la carte, que les écarts-types sont les plus élevés dans des zones 81 ne contenant pas de tissu mais des artères. Les écarts- types ne sont pas plus élevés dans la zone ischémique que dans la zone contralatérale, confortant l'estimation des CBF dans cette zone, telle que décrite par la figure 10. Cette information n'est pas disponible avec les techniques selon l'art antérieur à l'invention. La figure 16 permet de décrire une carte relative à 5 l'estimation de l'écart-type des temps de transit moyen estimés a T à partir du modèle de fonction de répartition complémentaire « intégral Gamma corrigé » à deux paramètres conforme à l'invention. On peut constater en analysant la carte, que les écarts-types sont les plus 10 élevés dans des zones 82 ne contenant pas de tissu mais des artères. Il s'agit précisément de zones similaires à celles où les écarts-types des flux cérébraux sanguins sont également les plus élevés. Les écarts-types ne sont pas plus élevés dans la zone ischémique que dans la zone 15 contralatérale, confortant l'estimation des MTT dans cette zone, estimation telle que décrite par la figure 11. Cette information n'est pas disponible avec les techniques selon l'art antérieur à l'invention. 20 La figure 17 permet de décrire une carte relative à l'estimation de la probabilité que le flux sanguin cérébral soit dans l'intervalle de confiance [CBF-cr-cBF,CBF+crcBF] sous le modèle de fonction d'entrée artérielle « tri-Gamma » à douze paramètres conforme à 25 l'invention. On peut constater en analysant la carte que, mis à part des ajustements aberrants, les probabilités sont centrées autour de 0,68. C'est une valeur à laquelle on est en droit de s'attendre si la distribution de probabilité a posteriori de CBF suit une loi normale. 30 Grâce aux cartes présentées précédemment, l'invention permet de mettre à la disposition d'un utilisateur tout un ensemble d'informations utiles, informations qui ne pouvaient être disponibles à l'aide des techniques 35 connues de l'état de la technique. Cette mise à disposition est rendue possible par une adaptation de l'unité de traitement 4 selon les figures 1 ou 2 en ce que ses moyens pour communiquer avec le monde extérieur de l'unité 4 sont aptes à délivrer les paramètres estimés 14 selon un format approprié à une interface homme-machine 5 apte à restituer à un utilisateur 6 lesdits paramètres estimés sous la forme par exemple de cartes telles qu'illustrées par les figures 9 à 17. Nous pouvons remarquer en outre que l'invention permet de supprimer des étapes nécessaires selon les techniques de l'art antérieur à l'invention. Parmi diverses étapes nous pouvons citer la suppression de : - la sélection manuelle ou automatique ou détermination d'une fonction d'entrée artérielle globale ou de fonctions d'entrée artérielle locales ; - la mise en oeuvre par l'unité de traitement d'un algorithme de déconvolution du modèle standard de perfusion C(t)=BF•Ca(t)OR(t) pour chaque voxel, sur la base de ou desdites fonction(s) d'entrée artérielle(s) sélectionnée(s), afin d'obtenir une estimation de BF.R (t) . Nous pouvons remarquer en outre que selon l'invention, l'unité de traitement est adaptée pour mettre en oeuvre un algorithme d'estimation conjointe des paramètres du modèle de perfusion global contraint en lieu et place d'un algorithme de déconvolution conformément aux techniques antérieures. Grâce à l'invention, les informations délivrées sont plus nombreuses et justes. Les informations dont dispose un chercheur ou un praticien sont ainsi de nature à accroître la confiance de ce dernier dans la détermination d'un diagnostic et d'une décision thérapeutique. Pour améliorer les performances du système selon l'invention, celle-ci prévoit que l'unité de traitement 4 puisse être dotée de moyens pour paralléliser des calculs sur les voxels de l'image pour lesquels l'estimation des paramètres est requise. Une telle unité de traitement peut comporter des moyens tels que des microprocesseurs graphiques (Graphical Processor Unit (GPU) en langue anglaise) ou mettre en oeuvre des grappes de calcul (clusters, en langue anglaise). En variante une telle unité de traitement peut mettre en oeuvre des moyens programmés tels que les Méthodes de Monte-Carlo parallèles. En variante, l'unité de traitement conforme à l'invention peut s'appuyer sur des moyens de calcul distants. Les temps de calculs peuvent ainsi être considérablement réduits. En outre, une unité de traitement conforme à l'invention peut comporter des moyens pour mémoriser un ou plusieurs modèle(s) global(aux) de perfusion contraint(s) ou une pluralité de sous-modèles de fonction d'entrée artérielle ou de fonction de répartition complémentaire du temps de transit ou une pluralité de contraintes basées sur les paramètres desdits sous-modèles de fonction d'entrée artérielle. En variante une telle unité de traitement peut autoriser un chargement de tels modèles globaux de perfusion contraints ou de sous-modèles voire de contraintes par l'intermédiaire de ses moyens pour communiquer avec le monde extérieur. Ces derniers peuvent être aptes à mettre en oeuvre une communication filaire ou distante selon des techniques connues de l'état de l'art. L'invention a été décrite et illustrée dans le domaine de l'Imagerie de Perfusion par Susceptibilité de Contraste Dynamique (Dynamic Susceptibility Contrast Magnetic Resonance Imaging (DSC-MRI) en langue anglaise). L'invention ne saurait être limitée à ce seul exemple d'application mais peut s'appliquer également à d'autres technologies médicales telles que l'imagerie de perfusion par Tomodensitométrie notamment.5 REVENDICATIONS1. Procédé pour estimer des paramètres hémodynamiques de perfusion (14) d'un volume élémentaire - dit voxel - d'un organe, ledit procédé étant mis en oeuvre par une unité de traitement (4) d'un système d'analyse d'imagerie de perfusion, et comportant une étape pour estimer des paramètres hémodynamiques à partir de signaux (15) de perfusion S(t), caractérisé en ce que l'étape pour estimer lesdits paramètres hémodynamiques consiste en une estimation conjointe (51) des paramètres e d'un modèle global de perfusion contraint comportant : - une première relation entre le signal de perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans ledit voxel au cours du temps t ; - une deuxième relation entre la concentration C(t) d'un agent de contraste, un paramètre physiologique 8, un modèle paramétrique ou semi-paramétrique Ca(t,0a) d'une fonction d'entrée artérielle et un modèle paramétrique ou semi-paramétrique R(t,OR) d'une fonction de répartition complémentaire du temps de transit dans le voxel, Oa et OR étant respectivement les paramètres desdits modèles Ca(t,Oa) et R(t,OR) ; - une troisième relation 11(0a)=0 entre les paramètres Oa du modèle de la fonction d'entrée artérielle Ca(t,ea).30 2. Procédé selon la revendication 1 caractérisé en ce que . le signal de perfusion S(t) a été préalablement obtenu partir d'images numériques (12, 13) délivrées par un appareil (1) d'imagerie de perfusion par résonance magnétique ; la première relation du modèle global de -k-TE-C(t) perfusion contraint s'exprime par S(t) = Soe où S0 est l'intensité moyenne du signal avant l'arrivée de l'agent de contraste dans le voxel, TE est un temps d'écho et k une constante non nulle ; la deuxième relation du modèle global de perfusion contraint s'exprime par C(t)=77.6•Ca(t,Oa)f R(t,OR) où désigne le produit de convolution, i7 est une constante non nulle et le paramètre physiologique B est le flux sanguin BF dans le voxel. 3. Procédé selon la revendication 1 caractérisé en ce que : le signal de perfusion S(t) a été préalablement obtenu à partir d'images numériques (12, 13) délivrées par un appareil (1) d'imagerie de perfusion par tomodensitométrie ; la première relation du modèle global de perfusion contraint est une relation de proportionnalité S(t)=a.C(t) où a est une constante non nulle; la deuxième relation du modèle global de perfusion contraint s'exprime par C(t)=q.O.Ca(t,Oa)OR(t,OR) où 0 désigne le produit de convolution et r7 est une constante non nulle et le paramètre physiologique 6 est le flux sanguin BF dans le voxel. 4. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il est mis en oeuvre par itérations successives pour une pluralité de voxels considérés. 5. Procédé selon l'une des revendications précédentes caractérisé en ce que l'estimation (51) conjointe des paramètres O du modèle global de perfusion contraint mise en oeuvre par l'unité de traitement (4) est réalisée au moyen de la méthode de Bayes ou des méthodes du Maximum de Vraisemblance ou encore des méthodes des Moindres Carrés non linéaires. 6. Procédé selon l'une des revendications précédentes caractérisé en ce qu'il comporte une étape (55) pour quantifier l'adéquation du modèle global de perfusion contraint à des données de perfusion D=[S(tl),...,S(tN)] expérimentales par un calcul d'une probabilité p(DI M) = f p(OI M) p(DIO,M)d0 de ces o données sachant ledit modèle. 7. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte une étape préalable (56a, 56b, 56c, 57) pour choisir un modèle global de perfusion contraint parmi une pluralité. 8. Procédé selon la revendication précédente caractérisé en ce qu'il est mis en oeuvre itérativement par l'unité de traitement (4) pour chaque modèle global de perfusion contraint connu de ladite unité de traitement. 9. Procédé selon la revendication précédente lorsqu'elle dépend de la revendication 6 caractérisé en ce que sont (55a) délivrés les paramètres estimés selon le modèle de perfusion dont la probabilité sachant ledit modèle est la plus grande. 10. Procédé selon l'une des revendications précédentes caractérisé en ce que le modèle Ca(t,Oa) d'une fonction d'entrée artérielle est un modèle dit « tri-Gamma » à douze paramètres défini par : C t O = a (t - t0 )a0-1 e-(t-t0)/e0 + b (t - t1 )a' 1 e-(t-tl,/ c (t - t2 )a2-1 (, "a) ~3oaaT (ao) fia'F(al) fi2a2r(a2) Oa =(a,b,c,an,Qp,t0, d'Eulepa2,fi2,t2) t ( ) : fonction Gamma d'Euler + 11. Procédé selon l'une des revendications précédentes caractérisé en ce que le modèle R(t,®R) d'une fonction de répartition complémentaire du temps de transit dans le voxel est un modèle dit « Gamma Intégral Corrigé » à deux paramètres défini par : t R(t,OR)=H(t)- fh(r,MTT, f)dr 0 MTT h(r,MTT,f)= t /3 e-tl/3 OR = (MTT,fl) H( ) : fonction généralisée créneau de Heaviside F ( ) : fonction Gamma d'Euler 12. Procédé selon l'une des revendications précédentes caractérisé en ce que la troisième relation Y(Oa)=0 du modèle global de perfusion MTT fi I MTT MTT>0,fi>0+0o contraint s' exprime comme f Ça (t,Oa)dt = CO où Co est 0 une constante arbitraire non nulle, identique pour tout voxel d'intérêt. 13. Procédé selon l'une des revendications précédentes caractérisé en ce que la troisième relation du modèle global de perfusion contraint +00 l(Oa)=O s'exprime comme Cca(t,Oca)dt=CO où 0 Cca (t,Oca) est une composante du modèle Ça (t,Oa) modélisant la concentration de l'agent de contraste lors de son premier passage. 14. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte une étape (53) pour calculer une information complémentaire sous la forme d'une fonction d'entrée artérielle estimée Ca(40a). 15. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte une étape (54) pour calculer une information complémentaire sous la forme d'une fonction de répartition complémentaire estimée R(t,OR). 16. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte une étape (52a, 53a, 54a) pour calculer une information complémentaire sous la forme d'un intervalle de confiance associé à un paramètre du modèle global de perfusion contraint. 17. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte uneétape (52b, 53b, 53c) pour calculer une information complémentaire sous la forme d'un taux de pari associé à un paramètre du modèle global de perfusion contraint. 18. Procédé selon l'une quelconque des revendications précédentes caractérisé en ce qu'il comporte une étape pour délivrer (52) des paramètres hémodynamiques estimés (14) à une interface homme- machine (5) apte à restituer à un utilisateur (6) lesdits paramètres estimés. 19. Procédé selon l'une quelconque des revendications 14 à 18 caractérisé en ce qu'il comporte une étape pour délivrer (52) toute information complémentaire à une interface homme-machine (5) apte à restituer à un utilisateur (6) les paramètres estimés ainsi que ladite information complémentaire. 20. Unité de traitement (4) comportant des moyens de mémorisation, des moyens pour communiquer avec le monde extérieur et des moyens de traitement, caractérisée en ce que : - les moyens pour communiquer sont aptes à recevoir du monde extérieur un signal S(t) (15) obtenu par imagerie de perfusion ; - les moyens de mémorisation comportent un modèle global de perfusion contraint, ledit modèle 30 comportant : - une première relation entre le signal de perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans ledit voxel au cours du temps t ; 35 - une deuxième relation entre un paramètre physiologique B, la concentration C(t), un modèle paramétrique ou semi-paramétrique Ca(t,Oa) d'une fonction d'entrée artérielle et un modèle paramétrique ou semiparamétrique R(t,OR) d'une fonction de répartition complémentaire du temps de transit dans le voxel, Oa et OR étant respectivement les paramètres desdits modèles ; - une troisième relation Y(Da)=0 entre les paramètres Oa du modèle de la fonction d'entrée artérielle Ca(t,Oa) ; - les moyens de traitement sont adaptés pour mettre en œuvre un procédé pour estimer des paramètres hémodynamiques de perfusion (14) selon l'une quelconques des revendications 1 à 19. 21. Unité de traitement comportant des moyens de mémorisation, des moyens pour communiquer avec le monde extérieur, des moyens de traitement, caractérisée en ce que : - les moyens pour communiquer sont aptes à recevoir du monde extérieur des images numériques de perfusion (12, 13) ; - les moyens de mémorisation comportent un modèle global de perfusion contraint, ledit modèle comportant : - une première relation entre un signal de perfusion S(t) et une concentration C(t) d'un agent de contraste circulant dans ledit voxel au cours du temps t ; - une deuxième relation entre un paramètre physiologique 9, la concentration C(t), un modèle paramétrique ou semi-paramétrique Ca(t,Oa) d'une fonction d'entrée artérielle et un modèle paramétrique ou semi-paramétrique R(t,OR) d'une fonction de répartition complémentaire du temps de transit dans le voxel, Oa et OR étant respectivement les paramètres desdits modèles ; une troisième relation 1«Oa)=0 entre les paramètres Oa du modèle de la fonction d'entrée artérielle Ca(t,Oa) ; les moyens de traitement sont adaptés pour déterminer un signal de perfusion S(t) à partir d'images de perfusion (12, 13) reçues et pour mettre en oeuvre un procédé pour estimer des paramètres hémodynamiques de perfusion (14) selon l'une quelconques des revendications 1 à 19. 22. Unité de traitement selon les revendications 20 ou 21 caractérisé en que les moyens pour communiquer délivrent (52) des paramètres estimés (14) selon un format approprié à une interface homme-machine (5) apte à restituer à un utilisateur (6) lesdits paramètres estimés. 23. Système d'analyse d'imagerie de perfusion comportant une unité de traitement (4) selon l'une quelconques des revendications 20 à 22 et une interface homme-machine (5) apte à restituer à un utilisateur (6) des paramètres estimés (14) selon un procédé conforme à l'une quelconque des revendications 1 à 19 et mis en oeuvre par ladite unité de traitement (4).