Segmentierung der Prostata aus MRT-Bilddaten mittels eines statistischen Modells Stefan Dänzer1,5 , Stefan Freitag2 , Dirk Beyersdorff3 , Markus Scholz4 , Oliver Burgert1 , Jens-Uwe Stolzenburg5 1 ICCAS, Universität Leipzig 2 HTWK, Leipzig 3 Klinik für Radiologie, Charite Berlin 4 IMISE, Leipzig 5 Klinik für Urologie, Universitätsklinikum Leipzig stefan.daenzer@gmail.com Kurzfassung. In dieser Arbeit wird ein semi-automatisches Verfahren zur Segmentierung der Prostata in MRT Bildern vorgestellt. Dieses ba- siert auf der Kombination eines statistischen Formmodells mit einem lokalen statistischen Erscheinungsmodell, welche beide aus handsegmen- tierten Trainingsdaten generiert werden. Eine neue, robuste Kostenfunk- tion auf Basis eines Shrinkage-Schätzers wird an 10 Datensätzen eva- luiert. Dabei zeigt sich für das Segmentierungsergebnis mit Shrinkage- Schätzer eine mittlere Oberflächenabweichung von 1.77 mm, während die verbreitet eingesetzte Mahalanobis Distanz eine Abweichung von 2.67 mm erreicht. 1 Einleitung Die Visualisierung und Modellierung von anatomischen Strukturen findet zuneh- mend Anwendung in der Therapieplanung und der computerassistierten Chirur- gie. Um Strukturen patientenindividuell visualisieren zu können, müssen diese in radiologischen Aufnahmen zunächst segmentiert werden. Eine detailgetreue manuelle Segmentierung ist jedoch sehr zeitaufwändig. Ein akkurates Segmen- tierungsverfahren, welches möglichst wenig Interaktion benötigt stellt hier einen entscheidenden Vorteil dar. In dieser Arbeit wird ein modellbasiertes Verfahren zur Segmentierung vorgestellt, mit welchem anatomische Strukturen mit mini- maler Benutzerinteraktion erkannt werden können. 2 Material und Methoden 2.1 Stand der Forschung In der Literatur finden sich mehrere Publikationen zum Thema der Segmen- tierung der Prostata aus MRT Aufnahmen. Martin verwendet ein Verfahren, welches auf einem statistischen Formmodell (SSM) mit einer Beschränkung der Segmentierung der Prostata aus MRT-Bildern 115 Deformierung basiert [1]. Klein greift auf ein Verfahren basierend auf Atlas Mat- ching mit einer Mutual Information Kostenfunktion zurück [2]. Toth entwicklete ein Verfahren basierend auf SSM kombiniert mit einer combined Mutual Infor- mation Kostenfunktion [3]. 2.2 Beschreibung der verwendeten Methoden Das vorgestellte Verfahren basiert auf einem statistischen Form- und Erschei- nungsmodell, welches aus 87 MRT Aufnahmen von verschiedenen Patienten er- stellt wurde. Ein deformierbares Oberflächenmodell und ein graphbasiertes Such- verfahren in Kombination mit einer Kostenfunktion werden verwendet um das Modell iterativ an die gesuchte Struktur anzupassen. Für die Generierung des Formmodells wurde auf einen bereits implementierten Algorithmus von Hei- mann [4] zurückgegriffen. Die iterative Anpassung erfolgt gemäß der Metho- de in [5], durch ein deformierbares Modell, welches den Lagrange’schen Be- wegungsgleichungen zugrunde liegt. Dabei wird eine Oberflächentriangulation M = (V, E) mit Punkten p, q ∈ V und Kanten [p, q] ∈ E an das Bildvolumen angepasst. Jeder der n Punkte pi ∈ V, (i = 1, . . . , n) korrespondiert hierbei zu einem Punkt p̃i , (i = 1, . . . , n) im statistischen Formmodell. Innere Kräfte Fint und äußere Kräfte Fext verformen die Triangulation iterativ und wirken dabei explizit auf jeden Punkt pi ∈ V pt+1 i = pti + Fint (pti ) + Fext (pti ) (1) Die inneren Kräfte setzen sich aus zwei voneinander unabhängigen Teilkräften zusammen. Die Spannkräfte FT gleichen Unterschiede in den Kantenlängen aus. Die Steifekräfte FR beheben Differenzen bei Winkeln benachbarter Flächen. Bei- de Kräfte werden für jeden Punkt pi ∈ V im geometrischen Modell unabhängig voneinander berechnet und addiert. Äußere Kräfte verformen das geometrische Modell so, dass es möglichst gut an das zugrunde liegende Bildvolumen ange- passt wird. Dazu wird eine Kostenfunktion eingeführt, welche für jeden Punkt pi ∈ V die Abweichung von der Oberfläche bzw. die Anpassungskosten von pi berechnet. Zu diesem Zweck, werden sowohl für jeden Punkt pi ∈ V , als auch an K nach innen und außen verschobenen Positionen die Kosten berechnet. Da- durch erhält man für jeden Punkt 2K + 1 Stellen pik , welche äquidistant mit dem Abstand d entlang des Normalenvektor N (pi ) verteilt sind pik = pi + k ∗ d ∗ N (pi ) k ∈ [−K; K] (2) Unter den 2K + 1 Punkten kann nun der optimale Punkt s(pi ) mit minimalen Kosten bestimmt werden. Die äußeren Kräfte sind als Zugkräfte einer Linearen Feder definiert, welche die pi ∈ V in Richtung ihrer optimalen Position s(pi ) verschieben Fext (pi ) = γ(s(pi ) − pi ) (3) γ stellt einen Parameter zur Regulierung der Stärke der äußeren Kräfte dar. Zur Bestimmung der optimalen Verschiebungen s(pi ) aller Punkte pi des geometri- schen Modells wurde die Methode von Li et al. [6] verwendet. Dabei wird das 116 Dänzer et al. beschriebene Problem in das Problem des Findens eines Minimum s-t-Schnitt in einem gerichteten Graphen umgewandelt. Das Ergebnis des Verfahrens liefert eine Konstellation optimaler Verschiebungen s(p1 ), . . . , s(pn ), welche zusätzlich noch die Bedingungen erfüllt, dass sich benachbarte Punkte pi , pj ∈ V in ihrer optimalen Verschiebung s(pi ), s(pj ) um nicht mehr als ein vorher festgelegtes ∆k unterscheiden dürfen. Als Grundlage zur Modellierung einer Kostenfunktion im Zusammenhang mit statistischen Formmodellen zählen die Grauwertprofile, welche bereits 1993 von Cootes [7] eingeführt wurden. Dazu wird entlang der Normalen des betrachteten Punktes ein Grauwertprofil g abgetastet. Der Grundgedanke bestand darin, dass die Grauwertprofile näherungsweise einer multivariaten Gauß-Verteilung unter- liegen sollten. Um diese zu modellieren wurde zu jeder Landmarke pi in jedem der J Trainingsdatensätze das Profil gij abgetastet, sodass durch den Mittelwert g¯i und die empirische Kovarianzmatrix Σi die Verteilung der Profile für die i-te Landmarke geschätzt werden konnte. Auf dieser Grundlage lässt sich schließlich für jedes weitere Profil g die Wahrscheinlichkeit Pi (g), dass g zur Verteilung der Profile der i-ten Landmarke gehört, bestimmen Pi (g) = e− 2 DM,i (g) DM,i (g) = (g − g¯i )Σi−1 (g − g¯i )T 1 (4) Wobei DM,i (g) die Mahalanobisdistanz des Profils g zur Stichprobe der i-ten Landmarke bestimmt. Um das Verfahren der lokalen Suche stabiler zu machen, erfolgt die Berech- nung anhand eines Multi Resolution Approach. Bei jedem Downsampling des Originalbildes, werden zusätzlich bei jeder Halbierung der Auflösung die Abstän- de d der Proben entlang jedes Profils verdoppelt. Ein bekanntes Problem der Ma- halanobisdistanz besteht darin, dass die Kovarianzmatrix und der Erwartungs- wertvektor für kleine Stichproben stark fehleranfällig für Ausreißer sind. Um diesen Effekt zu mindern und somit auch für nicht so mächtige Trainingsmengen noch eine verlässliche Kostenfunktion zu erhalten wurden Shrinkage-Schätzer aus [8] für die neue Kostenfunktion verwendet. Durch den Einsatz dieser Schät- zer für die Kovarianzmatrix und den Erwartungswertvektor wird die Berech- nung wesentlich weniger anfällig für Ausreißer. Dabei werden die Schätzer für Erwartungswertvektor µ̂shr und Kovarianzmatrix Σ̂shr mit Hilfe der traditionel- len Schätzer µ̂ sowie Σ̂ bestimmt. Der Schätzwert für den Erwartungswertvektor setzt sich als gewichtetes Mittel aus µ̂ und b zusammen   1 1 ∑ . N µ̂shr = (1 − α) ∗ µ̂ + α ∗ b b = µ̂i ∗  . . (5) N i=1 1 In [8] wurde gezeigt, dass sich die optimale Gewichtung α durch 1 N λ̄ − 2λmax α = ∗ (6) T (µ̂ − b)T (µ̂ − b) berechnen lässt, wobei λmax den größten Eigenwert des Kovarianzschätzwertes Σ̂ und λ̄ den Durchschnitt über alle Eigenwerte darstellt. Die Berechnung des Segmentierung der Prostata aus MRT-Bildern 117 Tabelle 1. Mittlerer Fehler ASDmean und Root Mean Square Fehler RMS über alle Segmentierungen mit Mahalanobisdistanz (Mahal) und Shrinkage Schätzer (Shr). Mahal ASDmean Mahal RMS Shr ASDmean Shr RMS Minumum 1.664 1.9864 1.1675 1.4766 Unteres Quartil 2.2702 2.6115 1.549 1.9282 Median 2.6657 3.2821 1.7708 2.1847 Oberes Quartil 2.9124 3.495 2.223 2.7645 Maximum 3.35 3.8769 2.4822 2.9494 Shrinkage-Schätzers für die Kovarianzmatrix erfolgt ebenfalls durch ein gewich- tetes Mittel aus traditionellem Schätzer Σ̂ und der Diagonalmatrix C 1 ∑ N Σ̂shr = (1 − β) ∗ Σ̂ + β ∗ C C = λ̂n ∗ EN (7) N n=1 Mit EN als Einheitsmatrix mit Rang N, lässt sich die optimale Wichtung β, wie in [8] gezeigt, folgender Maßen berechnen ∑J { } ′ 1 1 J j=1 Spur (g j gj − Σ̂)2 β = ∗ { } (8) J Spur (Σ̂ − C)2 Sowohl α als auch β beschränkt man auf das Intervall [0, 1]. Die Shrinkage- Schätzwerte µshr und Σshr werden in Gleichung 4 für ḡi bzw. Σi eingesetzt, zur Berechnung der Anpassungskosten angewandt. 3 Ergebnisse Im Zuge der Validierung wurde die Genauigkeit der Segmentierungsmethode mit den jeweiligen Kostenfunktionen an 10 MRT Bilddatensätzen gemessen. Ein Segmentierungsergebnis ist in Abb. 1 dargestellt. Die Abweichungen der Segmen- tierungen zu den handsegmentierten Referenzen sind in Tabelle 1 dargestellt. 4 Diskussion In dieser Arbeit wurde ein Verfahren zur semi-automatischen Segmentierung von Organkonturen vorgestellt und anhand MRT Bildern der Prostata validiert. Das Verfahren beruht auf einer Kombination von statistischem Formmodell, statisti- schen Erscheinungsmodell und eines deformierbaren Modells. Der entscheidende Forschungsbeitrag stellt hierbei neben dem Segmentierungsergebnis, die Eva- luierung einer neuen, robusten Kostenfunktion dar, welche der herkömmlichen Mahalanobisdistanz überlegen ist. 118 Dänzer et al. Abb. 1. Exempl. Segmentierung unter Verwendung des Shrinkage Schätzers. Literaturverzeichnis 1. Martin S, Troccaz J, Daanenc V. Automated segmentation of the prostate in 3D MR images using a probabilistic atlas and a spatially constrained deformable model. Med Phys. 2010;37(4):1579–90. 2. Klein S, van der Heide UA, Lips IM, et al. Automatic segmentation of the prostate in 3D MR images by atlas matching using localized mutual information. Med Phys. 2008;35(4):1407–17. 3. Toth R, Chappelow J, Rosen M, et al. Multi-attribute non-initializing texture re- construction based active shape model. In: Proc MICCAI. vol. 11; 2008. p. 653–661. 4. Heimann T, Wolf I, Williams T, et al.; Springer. 3D active shape models using gradient descent optimization of description length. Proc Inf Process Med Imaging. 2005; p. 566–77. 5. Heimann T, Münzing S, Meinzer HP, et al.; Springer. A shape-guided deformable model with evolutionary algorithm initialization for 3D soft tissue segmentation. Proc Inf Process Med Imaging. 2007; p. 1–12. 6. Li K, Millington S, Wu X, et al.; Springer. Simultaneous segmentation of multiple closed surfaces using optimal graph searching. Proc Inf Process Med Imaging. 2005; p. 406–417. 7. Cootes T, Taylor C. Active shape model search using local grey-level models: A quantitative evaluation. In: Proc BMVC; 1993. p. 639–48. 8. Meucci A. Risk and Asset Allocation. Springer Verlag; 2005.