<%BANNER%>

Tensor Based Representation and Analysis of Diffusion Weighted Magnetic Resonance Images

Permanent Link: http://ufdc.ufl.edu/UFE0024350/00001

Material Information

Title: Tensor Based Representation and Analysis of Diffusion Weighted Magnetic Resonance Images
Physical Description: 1 online resource (125 p.)
Language: english
Creator: Barmpoutis, Angelos
Publisher: University of Florida
Place of Publication: Gainesville, Fla.
Publication Date: 2009

Subjects

Subjects / Keywords: diffusion, dti, estimation, hardi, mri, registration, segmentation, tensor
Computer and Information Science and Engineering -- Dissertations, Academic -- UF
Genre: Computer Engineering thesis, Ph.D.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract: Cartesian tensor bases have been widely used to model spherical functions. In medical imaging, tensors of various orders can approximate the diffusivity function at each voxel of a diffusion-weighted MRI data set. This approximation produces tensor-valued datasets that contain information about the underlying local structure of the scanned tissue. The goal in this dissertation is to use the information provided in the tensor field in an automated system that detects changes in the tensor coefficients and correlates them with different types and levels of injury in spinal cord. In order to achieve this, one has to follow several intermediate steps of tensor field processing. These steps include robust estimation of the various orders' diffusion tensor coefficients, tensor field segmentation, registration, and atlas construction as well as extraction of various features, such as anisotropy measures and fiber orientations. Several methods are presented for solving these problems in this dissertation. The proposed algorithms are either the first ever presented in literature for solving a specific problem, or improved alternative solutions to existing techniques. All methods are validated using synthetic simulated diffusion-weighted MR data and their effectiveness is demonstrated in several real datasets.
General Note: In the series University of Florida Digital Collections.
General Note: Includes vita.
Bibliography: Includes bibliographical references.
Source of Description: Description based on online resource; title from PDF title page.
Source of Description: This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility: by Angelos Barmpoutis.
Thesis: Thesis (Ph.D.)--University of Florida, 2009.
Local: Adviser: Vemuri, Baba C.

Record Information

Source Institution: UFRGP
Rights Management: Applicable rights reserved.
Classification: lcc - LD1780 2009
System ID: UFE0024350:00001

Permanent Link: http://ufdc.ufl.edu/UFE0024350/00001

Material Information

Title: Tensor Based Representation and Analysis of Diffusion Weighted Magnetic Resonance Images
Physical Description: 1 online resource (125 p.)
Language: english
Creator: Barmpoutis, Angelos
Publisher: University of Florida
Place of Publication: Gainesville, Fla.
Publication Date: 2009

Subjects

Subjects / Keywords: diffusion, dti, estimation, hardi, mri, registration, segmentation, tensor
Computer and Information Science and Engineering -- Dissertations, Academic -- UF
Genre: Computer Engineering thesis, Ph.D.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract: Cartesian tensor bases have been widely used to model spherical functions. In medical imaging, tensors of various orders can approximate the diffusivity function at each voxel of a diffusion-weighted MRI data set. This approximation produces tensor-valued datasets that contain information about the underlying local structure of the scanned tissue. The goal in this dissertation is to use the information provided in the tensor field in an automated system that detects changes in the tensor coefficients and correlates them with different types and levels of injury in spinal cord. In order to achieve this, one has to follow several intermediate steps of tensor field processing. These steps include robust estimation of the various orders' diffusion tensor coefficients, tensor field segmentation, registration, and atlas construction as well as extraction of various features, such as anisotropy measures and fiber orientations. Several methods are presented for solving these problems in this dissertation. The proposed algorithms are either the first ever presented in literature for solving a specific problem, or improved alternative solutions to existing techniques. All methods are validated using synthetic simulated diffusion-weighted MR data and their effectiveness is demonstrated in several real datasets.
General Note: In the series University of Florida Digital Collections.
General Note: Includes vita.
Bibliography: Includes bibliographical references.
Source of Description: Description based on online resource; title from PDF title page.
Source of Description: This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility: by Angelos Barmpoutis.
Thesis: Thesis (Ph.D.)--University of Florida, 2009.
Local: Adviser: Vemuri, Baba C.

Record Information

Source Institution: UFRGP
Rights Management: Applicable rights reserved.
Classification: lcc - LD1780 2009
System ID: UFE0024350:00001


This item has the following downloads:


Full Text

PAGE 1

1

PAGE 2

2

PAGE 4

IwouldliketothankmysupervisorDr.BabaC.Vemuri,forhisinexhaustibleinspirationandintellectualsupportduringmydoctoratestudies.Hisinsightfulguidanceaectedpositivelynotonlymyacademicpersonalitybutalsomycharacteringeneral.Itisalsoveryimportantformetoacknowledgethenancialsupportthatmysupervisorgenerouslyoeredmeintheformofresearchassistantshipandtravelgrants.Furthermore,IwouldliketoexpressmysincereappreciationtoDr.JohnR.Forderwhosupportedmyeducationrelatedtothebiomedicalaspectofmyresearch.ThecontinuouscollaborationwithhislabintheMcKnightBrainInstituteattheUniversityofFloridaplayedimportantroleintheclinicalapplicationofthetechniquespresentedinthisdissertation.Moreover,IwouldliketothankDr.ArunavaBanerjee,Dr.JereyHo,andDr.AnandRangarajan,forservingasmembersofmyPhDcommitteeandfordevotingseveralhourstodiscuswithmedetailsofmyresearch.IshouldalsoacknowledgeDr.TimothyShepherd,Dr.DenaHowland,Dr.EvrenOzarslan,Dr.BingJian,andRitwikKumarforourproductivecollaborationanddiscussions.IwouldalsoliketothankSauravB.Chandra,andMinSigHwang,bothPhDcandidatesattheBiomedicalEngineeringDepartment,fortheircontinuousfeedbackandusefultechnicalcomments,whichhelpedmeimprovetheclinicalapplicabilityofthescientictoolsthatIdevelopedduringmystudies.Inaddition,Iwouldliketomentionthefriendlyandproductiveenvironmentthatwasformedbymylab-matesTingChen,SanthoshKodipaka,Dr.NicholasAndrewLord,Dr.AdrianPeter,AjitRajwade,OneilSmith,OzlemSubakan,Dr.FeiWang,YuchenXie,SenihaEsenYukselandmoreotherswhosharedwithmeunforgettablemomentsbothinandoutoflab.Moreover,IwouldliketothankMrs.RoryJeanDeSimoneandDr.GerhardX.Ritterwho,althougharenotdirectlyrelatedtomydissertation,bothhelpedmeto 4

PAGE 5

5

PAGE 6

page ACKNOWLEDGMENTS ................................. 4 LISTOFTABLES ..................................... 9 LISTOFFIGURES .................................... 10 LISTOFSYMBOLS .................................... 12 LISTOFABBREVIATIONS ............................... 13 ABSTRACT ........................................ 14 CHAPTER 1DIFFUSION-WEIGHTEDMRIPROCESSING:INTRODUCTION,PROBLEMS,GOALSANDAPPLICATIONS ........................... 15 1.1Introduction ................................... 15 1.2DiusionTensorFieldProcessing ....................... 19 1.2.1TensorFieldInterpolation ....................... 19 1.2.2TensorFieldSegmentation ....................... 21 1.2.3TensorFieldRegistration ........................ 22 1.2.4HigherOrderDiusionTensors .................... 23 1.3ProblemsinTensorProcessing ......................... 26 1.4GoalsandApplications ............................. 27 1.5RealandSyntheticDW-MRIData ...................... 29 2ESTIMATIONOFREGULARIZED2ND-ORDERDIFFUSIONTENSORFIELDSUSINGTENSORSPLINES ............................. 31 2.1Introduction ................................... 31 2.2MathematicalPreliminaries .......................... 32 2.3TensorSplines .................................. 33 2.3.1B-Splines ................................. 33 2.3.2TensorSplineInterpolation ....................... 34 2.3.3TensorSplineFitting .......................... 36 2.3.4Log-EuclideanSplines .......................... 39 2.3.5RobustTensorSplines ......................... 40 2.4ExperimentalResults .............................. 41 3APPLICATIONOFTENSORSPLINESINSEGMENTING2ND-ORDERTENSORFIELDS ........................................ 47 3.1Introduction ................................... 47 3.2DTISegmentationUsingTensorSplines ................... 48 3.3ExperimentalResults .............................. 50 6

PAGE 7

..... 58 4.1Introduction ................................... 58 4.2DiusionTensorsof4thOrder ......................... 60 4.2.1EstimationFromDWI ......................... 65 4.2.2DistanceMeasure ............................ 67 4.3ExperimentalResults .............................. 69 4.3.1SyntheticDataExperiments ...................... 69 4.3.2RealDataExperiments ......................... 72 5EXPONENTIALTENSORS:AFRAMEWORKFOREFFICIENTHIGHER-ORDERDT-MRICOMPUTATIONS .................. 76 5.1Introduction ................................... 76 5.2ExponentialDiusionTensors ......................... 76 5.2.1DistanceMeasure ............................ 77 5.2.2DistancefromtheClosestIsotropicCase ............... 78 5.2.3EstimationofEDTFieldfromDWI .................. 80 5.2.4DisplacementProbabilityProle .................... 80 5.3ExperimentalResultsandDiscussion ..................... 81 6ESTIMATIONOFDISPLACEMENTPROBABILITYFIELDSUSING4TH-ORDERTENSORS .............................. 85 6.1Introduction ................................... 85 6.2SphericalFunctionTensorialApproximation ................. 86 6.3Higher-OrderBasesforHARDIApproximation ............... 86 6.4FastEstimationfromHARDIDatasets .................... 89 6.5ExperimentalResults .............................. 89 6.5.1SyntheticDataExperiments ...................... 89 6.5.2RealDataExperiments ......................... 91 6.6ApplicationoftheDisplacementProbabilityFieldforTractography .... 92 6.6.1ExtractingTractosemasbyDiusingProbabilityIso-Surfaces .... 94 6.6.2ExperimentalResults .......................... 96 7REGISTRATIONOF4TH-ORDERTENSORFIELDS .............. 99 7.1Introduction ................................... 99 7.2LocallyAneRegistrationof4th-OrderTensorFields ............ 100 7.2.1Distancemeasure ............................ 100 7.2.2Registration ............................... 101 7.2.33DAneTransformationof4th-OrderTensors ............ 103 7.3Group-WiseRegistrationandAtlasConstruction .............. 103 7.3.1RiemannianMetricforPositive-ValuedRealFunctions ........ 103 7.3.2Group-WiseRegistrationof4th-OrderTensorFields ......... 104 7.3.3ImplementationDetails ......................... 106 7.3.4RobustAtlasConstruction ....................... 107 7

PAGE 8

.............................. 108 7.4.1ExperimentsUsingSyntheticDatasets ................ 109 7.4.2RealDataExperiments ......................... 110 8CLINICALAPPLICATIONS ............................ 113 8.1Introduction ................................... 113 8.2ComparisonofControlandInjuredSpinalCordDatasets .......... 113 9CONCLUSIONS ................................... 117 REFERENCES ....................................... 119 BIOGRAPHICALSKETCH ................................ 125 8

PAGE 9

Table page 2-1Tensoreldapproximationerrorsusingvariousalgorithms ............ 43 2-2Tensoreldapproximationerrorsusingrobustfunction .............. 43 2-3Approximationerrorsofrobustsplinesusingdierentmetrics .......... 43 4-1FormulastocomputethetensorcoecientsDi;j;k. ................. 65 5-1ComparisonofDTIframeworks ........................... 81 5-2PropertiesofDTIframeworks ............................ 81 6-1Proposedbasesfunctions ............................... 87 9

PAGE 10

Figure page 1-1IllustrationoftheDW-MRIacquisitionusingvariousdiusionsensitizinggradientdirections. ....................................... 16 1-2ADW-MRIdatasetfromaratbrain. ........................ 17 1-3ThediusiontensorcoecientscomputedfromaDW-MRIdataset. ....... 18 2-1Illustrationofatensorspline ............................ 35 2-2Comparisonoftensoreldapproximationmethods ................ 42 2-3Tensorsplineexampleinrealhumanhippocampusdata. ............. 44 2-4Tensorsplineexampleinisolatedrathippocampus. ................ 46 3-1Tensoreldsegmentationexperiment ........................ 51 3-2Tensoreldsegmentationexampleofsyntheticsmoothlyvaryingtensoreld. .. 52 3-3Illustrationofsegmentationunderpartialvolumingeects ............ 52 3-42DsegmentationoftheDT-MRIeldfromanisolatedrathippocampus .... 55 3-53DsegmentationoftheDT-MRIeldfromanisolatedrathippocampus .... 57 4-1Quantitativecomparisonofberorientationserrorsfrom4th-ordertensors. ... 70 4-2Quantitativecomparisonofthegeneralizedtraceof4th-ordertensors. ...... 71 4-3Fiberorientationerrorsforvariousnoiselevels. .................. 72 4-5Visualcomparisonof4th-ordertensorelds ..................... 74 4-6Visualizationofthe4th-ordertensoreldfromaratspinalcord. ......... 75 5-1ComparisonbetweenFAand4th-orderfisomap. .................. 79 5-2ComparisonofExponentialtensorsandotherhigher-ordertensorframeworks. 82 5-3Exponentialtensoreldfromaratopticchiasm. .................. 83 5-4Illustrationofexponentialtensorprocessing. .................... 84 6-13Dplotsofthreeofthebasesfunctions ....................... 88 6-2SyntheticsimulatedDW-MRIdata ......................... 90 6-3Comparisonofberorientationerrorsforvariousmethods ............ 91 6-4Estimatedprobabilityprolesfromrealdataofarat'sspinalcord ........ 92 10

PAGE 11

.................... 95 6-6Estimatedtractosemasfromreadhippocampaldata. ............... 96 6-7Exampleofbertrackingusingtractosemas. .................... 98 7-1Quantitativecomparisonofregistrationmethods. ................. 108 7-2Comparisonoftensoratlasescomputedbyvariousmetrics. ............ 110 7-3Locallyaneregistrationofhippocampaldata. .................. 111 7-44th-ordertensoreldatlasconstruction. ...................... 112 8-1TheS0imagesfromcontrolandinjuredratspinalcords. ............. 113 8-2Comparisonoftheberorientationsestimatedinthecontrolandinjuredcord. 114 8-3Quantitativecomparisonsofspinalcorddata. ................... 115 8-4ComparisonofspinalcorddatausingRiemanniandistances. ........... 116 8-5Comparisonoftractosemasfromcontrolandinjuredspinalcorddatasets .... 116 11

PAGE 12

Symbols:bThediusionweightingfactor(b-value)d(g)ThediusivityfunctionDThediusiontensorgThedirectionofthemagneticgradient(unitvector)GThemagneticgradientP(r)Thewatermoleculedisplacementprobability,givendisplacementrqThereciprocalspacevectorSTheacquiredDW-MRIsignalThegyromagneticratioThedurationoftheapplieddiusionsensitizinggradients 12

PAGE 13

Abbreviations:DTIDiusionTensorImagingorImage(s)DW-MRIDiusion-WeightedMagneticResonanceImagingorImage(s)EDTExponentialDiusionTensor(s)HARDIHighAngularResolutionDiusionImage(s)MRIMagneticResonanceImagingNMRINuclearMagneticResonanceImagingPSDPositiveSemi-DeniteROIRegionOfInterestSPDSymmetricPositive-Denite 13

PAGE 14

14

PAGE 15

StejskalandTanner ( 1965 )havemodeledtheattenuationoftheacquiredsignalbyamono-exponentialmodelgivenby 15

PAGE 16

ThisgureillustratestheacquisitionofaDW-MRdatasetthatconsistsofNimages.TheimageswereacquiredbyapplyingasetofNmagneticgradientsthatsensitizediusionindierentorientations.Theorientationsaregivenbythereciprocalspacevectorsq1:::qN. ThedisplacementprobabilityP(r)ofthewatermoleculesisrelatedtotheDW-MRIsignalattenuationthroughthefollowingFourierintegralrelationship( Callaghan 1991 ; Ozarslanetal. 2006 ). 1{2 ,risthedisplacementvector,qisthereciprocalspacevector,andS(q)istheDW-MRIsignalvalueassociatedwithvectorq.ThevectorqisdenedasafunctionofthemagneticgradientGgivenbyq=G=2,whereisthegyromagneticratioandisthedurationoftheapplieddiusionsensitizinggradients.Hence,bothvectorsqandGarepointingtothesamedirectiong=q=jqj=G=jGj.Furthermore,therelationbetweenthevectorqandthediusionweightingfactorbisgivenbytheexpressionb=4jqj2t,wheretistheeectivediusiontime.Ifweassumethatthediusionpropagatorfollowsamulti-variate(3-dimensional)Gaussiandiusionmodel,thenitsinverseFouriertransformcanbecomputedanalyticallyandisalsoamulti-variateGaussian( Basseretal. 1994 ).Inthiscase,thesignal 16

PAGE 17

ADW-MRIdatasetfromanexcisedratbrain.TheS0imageisshownonthetopleft.Therestofthe45imageswereacquiredusingasetof45diusiongradientdirections. attenuationmodelinEquation 1{1 canbegeneralizedasfollows 1-1 illustratestheacquisitionofDW-MRimagesusingasetofmagneticeldgradientdirections.Figure 1-2 showsarealratbraindatasetthatconsistsof45diusionweightedimagesandtheS0image(shownonthetopleft).Observationoftheimagesshowsthatthecontrastpatternschangewhenweapplygradientsthatsensitizediusionindierentdirections. 17

PAGE 18

A)ThediusiontensorcoecientscomputedfromtheDW-MRIdatasetinFigure 1-2 .Thelowerthreecoecientsaresymmetrictotheuppercoecientsofthematrixandhenceomitted.B)Eachtensorcanbevisualizedasanellipsoidcenteredatthecorrespondingvoxel. ThecorrespondingdiusiontensorcoecientscanbecomputedbyttingthemodelinEquation 1{3 totheacquiredimagesusingnon-linearminimizationtechniques( Wangetal. 2004 )andareshowninFigure 1-3 .Thisgurealsoshowsanothermoredescriptivewayofpresentingthetensorelddata,byplottingthecorrespondingellipsoidalsurfacesateachvoxel.TheellipsoidscanbecomputedbyevaluatingggTDgforg2S2,whereS2denotesthespaceofunitvectors,i.e.unitsphere.Thisvisualizationofthematrix-valuedeldgivesaverydetailedinsightoftheberstructuresandconnectivityinthescannedtissue.Furtherprocessingandstudyofsuchdatasetshelpsusrelatethevariousbrainfunctionsanddisorderswithspecicpatternsinchangesofthediusivityandnallydevelopmechanismsforautomaticdetectionandclassicationofthoseneurologicalconditions.Suchaframeworkfordiusiontensoreldprocessingwillincludealgorithmsforinterpolation,segmentation,registrationoftensoreldsaswellasmethodsforatlasconstruction.Section 1.2 presentsabriefintroductiononthespaceofdiusiontensorsandareviewofseveralexistingtensoreldprocessingtechniques. 18

PAGE 19

Terras 1988 ),wheretheRiemannianmetricisdenedastheinnerproductassignedtoeachpointofthisspace.Byusingthismetric,onecancomputegeodesic(shortest)distancesbetweenthepoints(diusiontensors)ofthisspaceandcanbeemployedintensoreldinterpolation( FletcherandJoshi 2004 ; WangandVemuri 2004 2005 ; Lengletetal. 2004 ; Pennecetal. 2005 ; San-JoseEsteparetal. 2005 ).Interpolationisrequiredinseveralalgorithmsfortensoreldregistration,segmentation,atlasconstructionandbertacking.Belowwebrieyreviewsomeoftheexistingtensorinterpolationtechniques. Pajevicetal. 2002 ),doesnotpreservemostoftheproperties,e.g.,thevalueofthedeterminantofthediusiontensors.Moreover,anyprocessingoftheseSPDtensorsbytreatingthemcomponentwiseandignoringthepositivedenitenesspropertywillleadtoestimationofunnaturaldiusivityvalues.Smoothinterpolationoforientationeldshasbeenproposedby Barretal. 1992 ).Although2ndordertensorscontainthenotionoforientation(e.gtheorientationofitseigenvectors),theirstructureismuchmorecomplicated.Byinterpolatingsmoothlytheorientationsofatensoreld,weignoreotherinformationprovidedbythetensors,(i.e.,eigenvalues,determinantetc.). WangandVemuri ( 2004 2005 )usedthesymmetrizedKullback-Liebler(KL)divergenceasa\distance"measurebetweentwoSPDtensors.ThesymmetrizedKL,denotedherebyKLswasshowntobeasecondorderapproximationtothetrue 19

PAGE 20

FletcherandJoshi ( 2004 ); Batcheloretal. 2005 ); Lengletetal. 2006 ).Havingatensoreld,e.g.,avolumetricDT-MRI,wecanusethegeodesiccurvesbetweenspatiallyconsecutivetensorsinordertointerpolateandapproximatethedataset.However,noneoftheearlierreportedmethodsongeodesiccurvecomputationbetweentensors( FletcherandJoshi 2004 ; Batcheloretal. 2005 ; Lengletetal. 2006 )usehigherordersmoothnessconstraintsinachievingtheinterpolation/approximation.Thus,althoughthereiscontinuityoftheinterpolateddataset,higherordercontinuityandhencesmoothnessislacking.Recently,aLog-Euclideanmetricwasproposedby Arsignyetal. 2005 )forcomputingwithtensors.Inthiswork,theelementsfromthespaceofpositivedenitediusiontensors,P(3),aremappedtotheirtangentspace,denotedbySym(3),usingthematrixlogarithmmap.ThetangentspaceofP(3)formsavectorspaceofdimension6.Therefore,onecanusetheEuclideannormforcomputationsinthistangentspaceandnallybyusingtheinversemapping,theinterpolateddataaremappedbacktothespaceofpositivedenitematricesP(3).Approximationofmatrix-valuedimagescanbeachievedviavariousregularizationmethods.Forinstance,aPDE-basedapproachaswasproposedby WeickertandWelk ( 2005 ).Althoughthisisafeasibleapproachfordenoising/smoothingandinterpolationofmatrix-valuedimages,thelackofuseofanappropriatemetricdenedonthespaceofSPDmatricescouldleadtoundesirablelimitationsinthesolution,e.g.,lackofaneinvariancepropertyetc.Anothertensoreldregularizationmethodwasproposedby Westinetal. 2006 )usingnormalizedconvolutionandMarkovRandomFields(MRF)inaBayesianframework.TheSPDtensorsaretreatedasvectorsin6Dandtheircomponents 20

PAGE 21

1.2.2 Zhukovetal. 2003 ),proposedalevelsetsegmentationmethodthatsegmentsthescalaranisotropicpropertycomputedfromthediusiontensor.Byusingsuchascalareld,thedirectioninformationcontainedinthetensoreldhasbeenignored.Thus,suchamethodwillfailtocorrectlysegmenttwohomogeneousregionsofatensoreldthathavethesamescalaranisotropypropertybutareorientedindierentdirections. Feddernetal. 2003 )extendedthemeancurvatureowandself-snakesmodelstomatrix-valueddata.However,theirmethodemploystheEuclideanmetrictomeasuredistancebetweentensorsandnottheRiemannianmetricdiscussedearlier.Thus,itdoesnotpossesssomeoftheinterestingpropertiesthataccruefromtheuseofaRiemannianframework(e.g.,theaneinvariancepropertyetc.). WangandVemuri ( 2004 ),developedaregion-basedactivecontourmodelfortensoreldsegmentation.Theygeneralizedthewellknownregion-basedactivecontourmodelforscalareldsegmentationtothatoftensorelds,anddevelopedavariationalprincipleusingtheForbeniusnormofthedierenceoftensorsasadiscriminantinthedataterm.Inalaterarticlebythem,( WangandVemuri 2005 ),ananeinvarianttensordissimilaritymeasurewasproposedbyusingtheJ-divergencebetweenorientedGaussiansrepresentingtheSPDtensorsinthetensoreld.Thiswasincorporatedinaregion-basedactivecontourmodelforthesegmentationofDT-MRIdataintoeitherpiece-wiseconstantorpiece-wisesmoothregions. 21

PAGE 22

Lengletetal. 2006 )developedastatisticalsurfaceevolutionframeworkusingtheFisher-Raometric.Theyemployedtheprinciplethatwithinaregionthediusiontensorscanbemodeledbyusingstatisticsanddistributionsofdiusiontensors.Thesurfaceevolutionframeworkandtheactivecontourmodelscanbeextendedtocopewithmultipletypesofregions.Segmentationcanalsobeachievedbymeansofregistration.Inthisapproachweestimatetheunknowndeformationthattransformsadatasetintoanothergivendatasetwhosesegmentationisknown.ThereareseveralmethodsforregisteringdiusiontensoreldsandarebrieyreviewedinSection 1.2.3 Alexanderetal. 2001 ).Inthisworkatensorre-orientationoperationwasdenedasasignicantpartofthediusiontensoreldtransformationprocedure.Thisre-orientationoperationrotatesappropriatelythediusiontensorsduringtheregistrationprocedure,inordertomaintainthelocaltopologyofthediusivitypatterns.Aframeworkfornon-rigidregistrationofmulti-modalmedicalimageswasproposedby Ruiz-Azolaetal. 2002 ).Thistechniqueperformsregistrationbasedonextractionofhighlystructuredfeaturesfromthediusiontensorelddatasets.RegistrationofDTIusingquantitieswhichareinvarianttorigidtransformationsandcomputedfromthediusiontensorswasproposedby Guimondetal. 2002 ).Byregisteringtherigid-tranformationinvariantmaps,oneavoidsthere-orientationstepandthuscanreducethetimecomplexity.Thelocallyanemulti-resolutionscalarimageregistrationproposedby Juetal. 1996 )wasextendedtoDTIimagesby Zhangetal. 2004 ).Inthismethodtheimagedomainoftheimagebeingregisteredissubdivided(usingamulti-resolutionframework)intosmallerregions,andeachregionisregisteredusinganetransformation.Theanetransformationisparametrizedusingatransformationvector,arotation,andanSPD 22

PAGE 23

Zhangetal. 2007 ).Recently,variousmethodshavebeenproposedfornon-rigidregistrationoftwotensorimagesemployingcontinuousmodelsofdeformation(incontrasttothediscontinuouslocallyanemodels)andtensormetricsthatmakeuseofthegeometryofsymmetricpositivedenitetensors.Themethodproposedby Caoetal. 2006 )registersthetensoreldsusinglargedeformationdieomorphicmetricmappingoftensorelds,resultinginoptimizingforgeodesicsonthespaceofdieomorphismsconnectingtwodiusiontensorimages.Adeformablediusiontensorregistrationalgorithmwasproposedby Irfanogluetal. 2008 )whichcomputesaB-splinedeormationbetweentwotensoreldbyminimizingtheGeodesic-Loxodromesdistancesbetweenthetensors.However,inalltheseapproachesoneoftheregisteredimagesisxed,whichisbiasingtheobtainedresult.Alltheabovemethods(inSections 1.2.1 1.2.2 ,and 1.2.3 )performprocessing(interpolation,segmentation,andregistrationrespectively)ofDW-MRIdatasetsbasedonscalarimagesordiusiontensorimages.However,thediusiontensorsare2nd-ordertensorialapproximationsofthelocaldiusivity,whichfailstorepresentcomplexlocaltissuestructures,suchasbercrossings,andthereforeDTIprocessingofdatasetcontainingsuchcrossingsleadstoinaccurateresults.Wecanovercomethislimitationbyusingalgorithmsthatemployhigherorderapproximationsofthediusivityfunction,andtheyarebrieyreviewedintheSection 1.2.4 23

PAGE 24

BasserandPierpaoli ( 1996 )).Thesignicanceofthesemeasuresisthattheyprovideaneasytool,sincetheyarescalarquantities,formonitoringchangesintheanisotropicpropertiesofwhitematterberbundlescausedeithernaturallybyneurologicaldisordersorarticiallypost-surgeryorinresponsetoatreatment.Fractionalanisotropyhasbeensuccessfullyusedclinically,andareductionofitsvaluehasbeenreportedinpatientswithtrauma Ptaketal. 2003 ),multiplesclerosis Cercignaniaetal. 2001 ); Hasanetal. 2005 ),Hypoxic-IschemicEncephalopathy Wardetal. 2006 ),multiplesystematrophy Itoetal. 2007 ),meningitis Nathetal. 2007 )andotherpathologies.Similarly,reductionofthemeandiusivity,anotherscalarquantity,hasbeenalsoreportedinthecaseofischemictissue vanGelderenetal. 1994 ).However,theabovescalarmeasuresarecomputedusingthe2nd-ordertensorialapproximationofthediusivity,whichalthoughworkswellforsimpletissuestructures,failstoapproximatemorecomplextissuegeometrywithmulti-lobeddiusivityproles.Forinstance,thevalueoffractionalanisotropydropssignicantlyinareasofbercrossings,althoughtheselocationsareanisotropic.Inordertoovercometheabovelimitation,higherordertensorswereintroducedby OzarslanandMareci ( 2003 )torepresentmorecomplexdiusivityproleswhichbetterapproximatethelocaldiusivityfunction.Generalizedscalarquantitiessuchasthevarianceofdiusivityandthegeneralizedanisotropywerederivedasfunctionsofthehigherordertensorcoecients( EvrenOzarslan 2005 ).However,inalloftheseworksthehigher-ordertensorsareestimatedwithoutimposingthepositivityofthediusivityfunctionapproximation,whichissignicantlyimportantsincenegativediusivityvaluesarenonphysical.Theuseofa4th-ordercovariancetensorwasproposedby BasserandPajevic ( 2003 ).ThiscovariancetensorisemployedindeningaNormaldistributionof2ndorderdiusiontensors.Thisdistributionfunctionhasbeenemployedby BasserandPajevic ( 2007 ) 24

PAGE 25

MoakherandNorris ( 2006 )presentedawaytominimizethedistancebetweenagivenhigher-ordertensorandtheclosestelasticitytensorofhighersymmetry,usingdierentmetrics.Thisworkwasrecentlyextendedby Moakher ( 2008 )whereinhedevelopedaninterestingframeworktodescribethegeometryof4th-ordertensors.A4th-ordersymmetricpositivedenitetensorinthreedimensionsisrepresentedbya2nd-ordersymmetricpositivedenitetensorinsixdimensionsandtherefore,onecanusetheRiemannianmetricofthespaceof66SPDmatricesfortheSPD4th-ordertensorcomputations.However,thisframeworkfailstoparameterizethefullspaceofSPD4th-ordertensors.ForexamplethegroupsofSPDtensorsT(g)=ag41+bg42+cg43andT(g)=(ag21+bg22)2+cg43,wherea;b;carepositivevaluedcoecients,cannotbewrittenintheformof66SPDmatrices.Furthermore,theRiemannianmetricofsymmetricpositive-denite2nd-ordertensorsinsixdimensionsassignsinnitedistancebetween66SPDand66semi-denitematrices.TheSPDtensorsintheaforementionedexamplescorrespondto66semi-denitematricesandtherefore,theRiemannianmetricassignsinnitedistancebetweenthemalthoughtheircorrespondingcoecientscanbearbitrarilyclose.Thisleadstoarbitrarilylargeerrorsintensorprocessingandoneshouldthereforestrivetoseekamoreappropriateframeworktocharacterizethefullspaceofsymmetricpositivedenite4th-ordertensors.Morethanonedistinctbertractstructurewithinavoxelcanbealsoestimatedbyemployingothernon-tensor-basedmodelsforreconstructionofthediusion-weightedMRsignal.AlthoughthemaintopicofthisdissertationisthetensorbasedrepresentationofDW-MRIdatasets,webrieyreviewherethemajornon-tensor-basedtechniques.Someofthemodelsthathavebeenproposedinliteratureincludediscreteandcontinuousmixture 25

PAGE 26

Tuchetal. 2003 )and JianandVemuri ( 2007b )respectively,andthesphericalharmonictransformationby Frank ( 2002 ).Afterreconstructionofthesignal,onehastocomputeitsFouriertransforminordertoobtainthedisplacementprobabilitywhosepeakscorrespondtodistinctberorientations.Thedisplacementprobabilityprolescanalsobecomputedbytransformingthediusivityprolesusingthediusionorientationtransform(DOT)( Ozarslanetal. 2006 ).Multipleberorientationscanalsobeestimatedbyreconstructingtheorientationdistributionfunction(ODF)( Descoteauxetal. 2007 )usingthesocalledQ-ballimaging( Tuch 2004 ).Mostoftheabovetechniques( Basseretal. 1994 ; Tuchetal. 2003 ; Descoteauxetal. 2007 ; JianandVemuri 2007b )canbeexpressedasaspecialcaseofamoregeneralizedmethodinwhichtheDW-MRsignalcanbeexpressedastheconvolutionoverthesphereofaberbundleresponsefunctionwiththeODF( Tournieretal. 2004 ; JianandVemuri 2007c ).Inthissphericaldeconvolutionapproachthereisnolimitationregardingthenumberoftheestimateddistinctberpopulations. 26

PAGE 27

27

PAGE 28

2 ,TensorSplinesareproposedasamethodforestimatingacontinuousandsmoothlyvaryingdiusiontensoreldfromagivendiusion-weightedMRIdataset.Thismethodperformshigh-ordercontinuousinterpolationofatensoreldincontrasttotheexistingmethodsthatemploylinearinterpolationusingtheRiemanniangeometryofthespaceofSPDmatrices.TensorSplinesareemployedinChapter 3 asamoduleinaprobabilisticdiusiontensoreldsegmentationalgorithm.InChapter 4 Ipresenttherstevermethodforestimatingpositive4th-orderdiusiontensorsfromdiusion-weightedMRimages.Thismethodguaranteesthepositivityoftheestimateddiusivityvalues,whichisaphysicalpropertythatisnotbeingimposedbytheexistingmethods.Asaresultourmethodoutperformstheexistingtechniquesintermsofaccuracy.Anotherhigher-orderparametrizationofthediusivityispresentedinChapter 5 usingExponentialTensors.Inthismethodwemodelthediusivityfunctionusingtheexponentialofhigherorderpolynomials.Thisframeworkguaranteestheestimationofpositivevaluesandalsoallowsfasttensor-eldcomputationsduetothesimplegeometryoftheparametricspace.Toourknowledgethisisthefastestframeworkforhigher-ordertensoreldprocessingthatimposesthepositivityconstraintontheestimates.Tensorsof4thorderarealsousedinChapter 6 forcomputingthewatermoleculedisplacementprobability,whichisthenemployedforestimatingsymmetricaswellasasymmetriclocalberstructures.Weshouldemphasizethattheexistingmethodsestimateberorientationdistributionsthatarenaturallysymmetricandconsequentlycannotmodelasymmetriessuchassplayingbers.InChapter 7 wepresentamethodforgroupwiseregistrationandatlasconstructionof4thordertensorelds.Thismethodextendsexistingtechniquesthatemployscalarimagesor2nd-ordermatrixvaluedimageswhichcannothandlebercrossingdiusivityproles.Finally,inChapter 8 thetensoreldframeworkisdemonstratedinarealapplicationforclassifyingdierentlevelsofspinalcordinjuries.Thismethodinvolvesregistrationand 28

PAGE 29

SodermanandJonsson 1995 ):E(q)=1Xn=01Xk=11Xm=02Knm2(2q)4sin2(2)km2 L2#D0!(1{4)Inthisexpression,JmisthemthorderBesselfunction,kmisthekthsolutionoftheequation~Jm()=0withtheconvention10=0,andKnm=n0m0+2[(1n0)+(1m0)]. 29

PAGE 30

1{4 tosimulatemultiplebercrossings.Althoughthissystemissimplerthanthediusionalprocesseswithinrealneuraltissue,itisasuitabletestbedfortheproceduresdevelopedforberorientationestimationandithasbeenusedby Ozarslanetal. 2006 ).Similartothatin vondemHagenandHenkelman ( 2002 ),inoursimulationsweusedthefollowingparameters:L=5mm,=5m,D0=2:02103mm2/s,=20:8ms,=2:4ms,b=1500s/mm2andtheinniteseriesinEquation 1{4 wereterminatedatn=1000andk,m=10. 30

PAGE 31

2.2 ,wepresentabriefnoteonthegeometryofthespaceofdiusiontensorsandtherelatedmathematicalbackground.ThisisfollowedbythepresentationofanovelandrobustalgorithmforcomputingtensorsplinesinSection 2.3 .Finally,inSection 2.4 ,comparisonsofourtensorsplinealgorithmwithexistingmethodsofinterpolationandapproximationofDT-MRIdataarepresented.Validationresultsusingsyntheticallygeneratednoisytensorelddatawithoutliersarealsopresented. 31

PAGE 32

Pennecetal. 2005 ); Lengletetal. 2004 ); FletcherandJoshi ( 2004 ).ThespaceofranktwodiusiontensorscanbeviewedasaRiemannianSymmetricspace( Helgason 2001 ),whereaRiemannianmetricassignsaninnerproducttoeachpointofthisspace.Byusingthismetricwecancomputegeodesicdistancesbetweendiusionstensorsandcalculatestatisticsonthisspace( FletcherandJoshi 2004 ; Lengletetal. 2004 ; Pennecetal. 2005 ).Forexample,themeantensorofasetofdiusiontensors,cannowbecomputedasthattensorwhichminimizesthesumofsquaredRiemanniandistancesbetweenitselfandthegivensetoftensors.Themeantensormaybeemployedas,aninterpolant,forperformingprincipalgeodesicanalysisetc.IntheaforementionedRiemannianframework,thedistancebetweentwotensorsD1andD2isgivenbydist2(D1;D2)=trace(log(D11=2D2D11=2)2),wherelogisthematrixlogarithmoperation.Byusingthisdistancemeasure,thegeodesiccurve(shortestpath)betweenD1andD2isdeneduniquely.Thetangent,speciedatthersttensorwithrespecttotheotheronealongtheuniquegeodesicbetweenthem,.isa33symmetricmatrixandisgivenbytheRiemannian-logmap,LogD1(D2)=D11=2log(D11=2D2D11=2)D11=2.TheinverseoperationisgivenbytheRiemannian-expmap,ExpD1(D2)=D11=2exp(D11=2DD11=2)D11=2,whereexpisthematrixexponentialoperationandDisa33symmetricmatrix.WewillusethisRiemanniandistancebetweenSPDtensorsincomputingthedistancebetweenthegivendataandthetensorsplineapproximationofthedataaswellasincomputingtheweightedaveragefordeningthetensorspline.InthefollowingsectiontheRiemannianexponentialandlogarithmicmaps,andtheexpressionforageodesicbetweentwodiusiontensorswillbeusedinordertodeneandcomputethetensorsplines. 32

PAGE 33

1{1 by StejskalandTanner ( 1965 ),andusingalinearleastsquaresmodeltoestimatethediusiontensorsoverthe3Dlattice.Thissectionhasthreesubsections.FirstweprovideabriefreviewofB-splines.Next,wepresentanovelalgorithmforcomputingsplinesonagivenSPDtensoreld.AfterthatwepresenttensorsplinesusingtheLog-Euclideanmetricasanimprovementoverrecentworkby Arsignyetal. 2005 )(describedearlierinSection 1.2.1 ).Then,wepresentourrobusttensorsplineapproximationalgorithm.Finallyarobusttensorsplineapproximation(tting)techniqueispresented. 33

PAGE 34

ti+kti+1(2{3)Ni;k(t)functionsarepolynomialsofdegreek-1.CubicbasisfunctionsNi;4canbeusedfora3rddegreeB-spline.Knotsmustbeseriesofmonotonicallyincreasingnumbers.AmoredetaileddiscussiononB-splinescanbefoundinthearticleby deBoor ( 1972 ).OneusefulpropertyofthefunctionsNi;k(t)isthatNi;k(t)0,foralliandXi=0Ni;k(t)=1.Consideringtheaboveproperties,functionsNi;k(t)behaveasblendingfunctionsandEquation 2{1 isaweightedaverageofthecontrolpointsci. Pajevicetal. 2002 ),oncontinuoustensoreldapproximationachievessmoothness,however,aRiemannianframeworkisnotemployedfortensorcalculations.Inthissectionwedenetensorsplineswhicharecurvesinterpolatingorapproximatingmatrixvaluedfunctions,constructedusingthegeometryofthespaceofSPDtensors.Notethatwearedeningtensor-splinesbydoingweightedintrinsicaveragesonP(n)andchoosingtheweightfunctionstobeB-splines.Asanillustrationofinterpolationona1Dgridoftensors,gure( 2-1B )depictstheideaofusingweightfunctions(B-splineshere)toperformweightedaverageoftensorsusingtheRiemannianmetric.Thisweightedaveragingleadstothedesireddegreesplineinterpolant(approximantwhenusedinattingproblem)ofthediusiontensordata. 34

PAGE 35

BFigure2-1. A)TangentspaceofthemanifoldMofdiusiontensorsatpointp1.ThetangentvectorXpointstothedirectionofgeodesic(t)betweenthepointsp1andp2.B)AcubictensorsplineS(t),thatapproximatespisofa1-Dtensoreld.ThegivenpointspiandthepointsofthetensorsplineS(t)areSPDmatrices,elementsoftheRiemannianmanifoldM.However,thetensorsplineisinP(n)
PAGE 36

2{1 tothespaceoftensors.WecancomputethevalueS(t)ofthek-1degreeB-splineoftensorsforaparticulartvalue,bycalculatingaweightedintrinsicaverage,fP,ofthecontroltensorsci,wheretheweightsequaltothebasisfunctionswi=Ni;k(t),discussedearlier. 2{4 )oftensorsisdenedusingtheRiemanniandistanceinsteadoftheEuclideandistance,anditistheminimizerofthefunction min2P(n)()=min2P(n)1 2nXi=0widist(;ci)2(2{5)wheredist(:;:)istheRiemanniangeodesicdistance.Theweightedaveragecanbecomputedusingagradientdescentalgorithmwhichisanextensionofthealgorithmdescribedby FletcherandJoshi ( 2004 )forcomputingthemeanoftensors.Thegradientof()isgivenby Algorithm1: IntrinsicWeightedMeanofTensors 36

PAGE 37

2NN1Xi=0dist(S(ti);Di)2(2{7)InEquation( 2{7 ),theRiemannianmetricshouldbeusedforthedistancecalculation,sincethetensorspace,wherethedataandcontrolpointslive,isacurvedmanifold(convexcone).Weneedtondasetofcontrolpoints(c0;c1;:::;cN1+k2)thatformthesplineS(t)whichminimizestheenergyE.ThegradientofEwithrespecttocjisthengivenby, 2NXN1i=0rS(ti)dist(S(ti);Di)2rcjS(ti)(2{8)ThegradientofthesquaredistancebetweenS(ti)andDiwithrespecttoS(ti)equals, 2{8 )iswithrespecttocj,weneedtoexpressthegradientinEquation 2{9 byusingtangentvectorsatpointcj.Takingthisintoconsideration,Equation 2{9 canbeapproximatedbytheformulacj(Di;S(ti))=Logcj(Di)Logcj(S(ti)),soweobtain 2{8 is 2{11 and 2{10 inEquation 2{8 weobtain, 37

PAGE 38

~cj=Expcj1 Algorithm2: ControlTensorsestimation Thetimecomplexityforasingleiterationofalgorithm 2 isoforderO((kdc)N)wherekisthedegreeofthespline,disthedimensionalityofthedataset(for3Ddatad=3)andcisthenumberofiterationsofalgorithm 1 andNisthegiveninputdatasize(numberoftensorstobeapproximated).Intheexperimentsthatweperformedwefoundthatalgorithm 1 convergesinatmostin5iterations(c5),aCPUtimeof9.37secperiterationonaPentium2.4GHzprocessorforttingacubic(k=3)tensorsplineinadatasetofsize128128.Asexpected,thetimecomplexityofalgorithm 2 increasesasweincreasethedegreeofthesplinekorthedimensionalityofthedatasetd.Notethatwe 38

PAGE 39

2{10 canbelarge,ifthetensorsplineapproximationS(ti)isfarfromthetargetDi.WhenS(ti)tends(movescloser)toDiduringthesplinettingprocedure,theerrorintroducedbytheapproximationofEquation 2{10 tendstowardzero.Bysettingasmallthresholdeonthedierencebetweenconsecutiveiterates,theouterloopofalgorithm 2 willbeiteratedsucienttimesinorderfortheerrorofEquation 2{10 tobeassmallasneeded.Thus,thecontroltensorscj,whichareobtainedastheoutputofalgorithm 2 ,areestimatedbytakingthetruegeometryofP(n)intoaccountandnotdoingapproximationsthatundermineit.TensorSplinescanbeeasilyextendedtohigherdimensionaltensorelds.Forexampleconsiderthecaseofa2DNMtensoreld.InthiscasewearedoinganapproximationonP(n)<2.A(k1)degreetensorsplinethattstoourdatarequires(N+k2)x(M+k2)controltensorsand(N+2(k1))x(M+2(k1))monotonicallyincreasing(inboththedimensions)knots(tk+1;k+1;:::;tN+k2;M+k2).Notethatinthiscasetheknotsarevectorswith2elements,oneforeachparametricdimension.Finallythenewbasisfunctionsareformedbythetensorproductof1-dimensionalbasisfunctionsNi;j;k([t1t2])=Ni;k(t1)Nj;k(t2). Arsignyetal. 2005 )proposedanewLog-Euclideanmetricfortensorcalculations.Inthisframework,thediusiontensorsarerstmappedusingthematrixlogarithmicmaptothespaceofthesymmetricmatricesSym(n).Thereafter,theEuclideannormisusedinallcalculationsinthisspace.Finallybyusingthematrixexponentialmapping,thecomputedvaluesaremappedbacktothemanifoldP(n).UsingtheLog-Euclideanmetricwecanalsotasplinetothelogarithmicallymappeddatainthespaceofsymmetrictensors,andafterthatwecanmaptheinterpolated/approximated 39

PAGE 40

2.4 ,weprovideaquantitativecomparisonbetweentensorsplinesandsplinesusingtheLog-Euclideanmetricthatwewillcall\Log-EuclideanSplines". Wangetal. 2003 ).Arobustalgorithmshouldrejecttheseoutliersfromfurtherconsiderationinanyprocessingalgorithmsappliedtothedataset.Arobustfunctioncanbeappliedtotheenergyfunction,inordertoweightthegivendataDiappropriately.Wecanusearobustfunctionthatassignsweightsintheinterval[0;1],whereweightswhicharealmostzeroimplyrejectionofthecorrespondingdatapoint.Furthermore,highweightsshouldbeassignedtothedatapointswhosedistancefromtheunknownsplinecurveissmallandontheotherhandlowerweightsshouldbeassignedtothedatapointswhosedistancefromtheunknownsplineislarger.Letusconsiderthefollowingfunction(x)=ex2 2NXN1i=0(dist(S(ti);Di))(2{14)ThegradientofthisenergywithrespecttothecontroltensorsnowbecomesrcjE=1 40

PAGE 41

2-2 ).Forthegenerationofthiseld,wesimulatedthediusion-weightedMRsignalusingtheequationpresentedby SodermanandJonsson ( 1995 )(seeSection 1.5 fordetailsonformulasforgeneratingtheMRsignal).UsingthisprocessateachvoxeltheMRsignalwassimulatedasafunctionoftheanglebetweentheapplieddiusiongradientandtheorientationoftheber.Ateachvoxelofthe2Dsyntheticeldtheorientationoftheberwasassumedtobetangenttocirclescenteredatthelowerleftcorneroftheeld(Figure 2-2 ).Usingthisberstructurethediusion-weightedMRsignalattenuationwassimulatedfor21orientationsthatcorrespondtothesecond-ordertessellationoftheicosahedrononaunithemisphere,usingb-value=1500s/mm2.Formoredetailsonothersimulationparameters,seeAppendixB.Thediusiontensoreldwasestimatedfromthe21diusion-weightedimagesusingalinearleastsquarestechniqueappliedtotheloglinearizedStejskal-TannerEquation 1{1 .TheprincipaleigenvectorsoftheestimateddiusiontensoreldareshowninFigure 2-2 .Thisdiusiontensoreldwillbeconsideredasthegroundtrutheldfortheexperimentdescribedbelow.GaussiannoisewasaddedtotherealandimaginarypartsofthesimulateddiusionMRsignalandthenthemagnitudesignalcomputedfromthisnoisycomplex-valueddata.Fromthissignal,weestimatedthediusiontensorsasbeforeandthensubsampleditbyafactorof4.Thisprocesswasrepeatedfordierentamountsofsignaltonoise 41

PAGE 42

2-2 depictstheprimaryeigenvectoreldofthe99subsamplednoisytensoreldscorrespondingtosignaltonoiseratiosof5.0(top)and3.0(bottom).Ourgoalnowistocompareourtensorsplineapproximationandinterpolationmethodagainstexistingmethodsinliterature,aswellasourownmodicationsofthesetechniques. Figure2-2. ComparisonofapproximationmethodsusingaSNR=5.0(top)andaSNR=3.0(bottom).A)Primaryeigenvectorsofthenoisytensorelds.TherestofthecolumnsshowstheerrorinrobustapproximationusingB)Euclideanspline,C)PDEinterpolation,D)Log-Euclideanspline,E)Tensorspline.TheRiemannianmetricwasemployedforcomputingtheseerrors. Werstapproximated(tted)thenoisytensoreldsbyusingfourdierenttechniquesincludingoursandtheninterpolatedtheapproximation(tted)resultsbyafactorof4.Thefourmethodsthatweemployedwere,A)LinearapproximationoftheelementsoftheSPDtensors,B)Log-Euclideangeodesicapproximation( Arsignyetal. 2005 ),C)Riemanniangeodesicapproximation( FletcherandJoshi 2004 ),D)PDE-basedanisotropicnon-lineardiusion( WeickertandWelk 2005 )(table 2-1 ).Followingthis,wepresentatableofresultscomparingourmethodwithstatisticallyrobustimplementationofallthemethods(exceptthePDE-baseddiusionlterasitsnotasimplemattertoimplementthislterinarobustframework).Afterthis,wepresentanothertablewhichdepictsresultsofcomparisonofasplineversionofothers'workwithourown,allina 42

PAGE 43

Tensoreldapproximationerrorsusingvariousalgorithms Riem.Riem.Fro.Fro. Table2-2. Tensoreldapproximationerrorsusingrobustfunction Riem.Riem.Fro.Fro. robustframework.Inallthesecomparisons,asexpected,ourtensorsplinesalgorithmsout-performedthecompetingmethods.Weusetwomethodstomeasurethedistanceoftheestimatedtensoreldsfromthegroundtruthtensoreldofgure 2-2 ;a)theRiemannianmetricandb)theFrobeniusnormdenedasq 2-1 2-2 and 2-3 fortwonoisydatasetscorrespondingtoasignaltonoiseratioof3.0.Asevident,theerrorismuchsmallerforouralgorithmincomparisontotheothers.Theseresultsdemonstratethesuperiorperformanceofouralgorithmoverotherexistingmethods. Table2-3. Approximationerrorsofrobustsplinesusingdierentmetrics Riem.Riem.Fro.Fro. 43

PAGE 44

BFigure2-3. Zero-gradientimageofaslicefrom3Dvolumeofahumanhippocampusautopsyspecimen.Primaryeigenvectoreldbefore(A)andafter(B)cubictensorsplineapproximationfromthemarkedregionofimage. FortheexperimentpresentedinFigure 2-3 ,an8-mmcoronalsegmentofthehumanhippocampuswasdissectedfromanautopsyspecimenwithshortpostmortemintervaltoxation(21.2hrs)andnohistologicalevidenceofneuropathology.Diusiontensormicroscopydatawerecollectedfromthehippocampalspecimenusinga14.1-Tmagnetwithaprotocolthatincluded21uniquediusiongradientorientations,diusiontime=17ms,andb=1250s/mm2( Shepherdetal. 2007 ).Adiusiontensoreldwasestimatedfromthediusion-weightedimagesusinglineartting.Figure 2-3A showsathediusion-weightedMRimageobtainedusingazerogradienteldofaslicefromthe3Dvolume.Aredboxhasbeenusedtoidentifytheregionofinterest(ROI)containingareasonablynoisysectionintheDT-MRIvolume.Figures 2-3A and 2-3B showtheprimaryeigenvectorsofthediusiontensoreldinthe3DROIbeforeandaftercubictensorsplineapproximationrespectively.Noticethatalloftheoutliershavebeenrejectedandtheeldhasbeenapproximatedsmoothly. 44

PAGE 45

2-4 showsanotherrealdataexamplefromanisolatedrathippocampus.ThediusionweightedMRimagesforthisexamplewereacquiredusingthefollowingprotocol.Thisprotocolincludedacquisitionof22imagesusingapulsedgradientspinechopulsesequencewithrepetitiontime(TR)=1.5s,echotime(TE)=28.3ms,bandwidth=35kHz,eld-of-view(FOV)=4:54:5mm,matrix=9090with20-30continuous200-m-thickaxialslices(orientedtransversetothesepto-temporalaxisoftheisolatedhippocampus).Aftertherstimagesetwascollectedwithoutdiusionweighting(b0s/mm2),21diusion-weightedimagesetswithgradientstrength(G)=415mT/m,gradientduration()=2.4ms,gradientseparation()=17.8msanddiusiontime(T)=17mswerecollected.Eachoftheseimagesetsuseddierentdiusiongradients(withapproximatebvaluesof1250s/mm2)whoseorientationsweredeterminedfromthe2ndordertessellationofanicosahedronprojectedontothesurfaceofaunithemisphere.Theimagewithoutdiusionweightinghad36signalaverages(time=81min),andeachdiusion-weightedimagehad12averages(time=27minperdiusiongradientorientation)togiveatotalimagingtimeof10.8hperhippocampus.Temperaturewasmaintainedat200:2oCthroughouttheexperimentsusingthetemperaturecontrolunitofthemagnetpreviouslycalibratedbymethanolspectroscopy.Ingure 2-4 ,theproposedtensorsplineapproximationalgorithmiscomparedwiththerecentlyproposedLog-Euclideanmetricbasedapproximationalgorithmby Arsignyetal. 2005 ),discussedearlier.Figure 2-4B showsa3DviewoftheresultsusingLog-Euclideangeodesicapproximationandnon-robustTensorSplinealgorithms.NoticethatintheTensorSplineapproximationresultsthenoisehasbeenconsiderablyreduced.Thismaybeattributedtothehigher-ordersmoothnessimposedbytheTensorSplinedevelopedinthiswork.Figure 2-4A depictstheFAmapoftheoriginal(noisy)andtheapproximateddata.NotethattheFAmapaftertheapproximationismuchsmootherthattheFAmappriortheapproximation(see 2-4A ). 45

PAGE 46

BFigure2-4. RealDTIfromanisolatedrathippocampus.A)FAmapsbefore(top)andafter(bottom)tensorsplineapproximation.B)TheprincipaleigenvectoreldafterLog-Euclideangeodesicapproximation,andnon-robusttensorsplineapproximation. 46

PAGE 47

2 )asanapproximationmodule.Thesegmentationisachievedbyjointlyestimatingthelabel(assignedtoeachtensorresidingatavoxel)eldandtheparameters(controlpoints)ofeachsmoothtensorsplinemodelrepresentingthelabeledregions.ThelabeleldismodeledbyaGaussMarkovMeasureField(GMMF)andthesegmentationalgorithmveryecientlycomputestheposteriormarginalprobabilitydistributionofthelabeleld(giventheparametersdescribingtheregionmodel)astheglobalminimizerofalinearlyconstrainedquadraticenergy.Furthermorewepresentseveral3DDT-MRIsegmentationresultsusingsyntheticandrealdatasetfromanisolatedrathippocampus.Themotivationforsegmentingandanalyzingthehippocampusisduetoitsimportanceinsemanticandepisodicformationthatisparticularlyvulnerabletoacuteorchronicinjury( AmaralandWitter 1995 ; Squireetal. 2004 ).Incurrentclinicalpractice,weonlylookatthewholehippocampusanddescribeatrophyforepilepsy(hippocampalsclerosis),schizophrenia,depression,hypoxia-ischemia,traumaandAlzheimer'sdiseaseandotherdementias.Obviouslythesecannotbedistinguishedfromeachotherandhippocampalatrophyisalateimagingsignofpathology.However,weknowfrommorethan100-yearsofliteraturethatthehippocampusismadeofmanydierentcytoarchitecturalregionsandthattheseregionsareselectivelyvulnerabletotheaforementioneddiseases(e.g.theCA1andsubicularregionsareaectedbyAlzheimer'sdiseasewhilethedentategyrusisaectedbymedialtemporallobeepilepsy).TheseregionscanbedistinguishedbydiusiontensorMRI.Thus,thesegmentingtechniquesbeingdevelopedherecouldproveusefultoimprovingthesensitivityandspecicity 47

PAGE 48

Shepherdetal. 2006 ).Inourexperimentsweusethesestructuralresults( Shepherdetal. 2006 )forvalidationoftheobtainedsegmentation. Riveraetal. 2005 ).WeassumethatthegivenDT-MRIdatasetDiconsistsofKregions,where"i"isthelatticeindex.FurthermoreweassumethateachregionisrepresentedbyamodelSk(ti)wherek=1::K.Therearedierentchoicesforthemodel,whichcanbeeitherpiece-wiseconstantorasmoothlyvaryingtensoreldmodel.Inthespaceofdiusiontensorsthepiece-wiseconstantmodelhasaparameterkwhichisa33SPDmatrix.ThereforeforthismodeltheequationSk(ti)=kholdsforallti.InthecaseofsmoothlyvaryingtensoreldswechoseatensorproductoftensorsplinesasourmodelSk(ti),whoseparametersarethecontroltensorscidenedinSection 2.3.1 .TherelationbetweenthecontroltensorsandthetensorsplineisgivenbytheEquation 2{4 .Letbkibethelabelmap,wherebki=1indicatesthatthediusiontensoratthei-thlatticepointbelongstothek-thregionclassandbki=0otherwise.Consideringtheabovenotation,adiusiontensordatasetcanbemodeledasbeinggeneratedusingagenerativemodelgiveninthefollowingequation: 48

PAGE 49

Riveraetal. 2005 ).Ourformulationheresignicantlyextendstheirformulationtocopewithtensorelds.LetvkibetheprobabilitythatthediusiontensorDiwasgeneratedbythek-thmodel,thelikelihoodofthelabeleldisthengivenby 22(3{3)wheredist(:;:)istheRiemanniangeodesicdistancebetweentwoSPDtensors.Ifwecanestimatethemarginalprobabilitiespki,thelabeleldforahardsegmentation(assignslabelsinayes/nofashionandnotwithaprobability))canbeestimatedbythemaximumposteriormarginals(MPM)estimator( Marroqunetal. 2001 )whichisdenedasbki=8>><>>:1ifpkipli,fork6=l0otherwiseByusingtheabovelabeleldestimatortheunknownvariablesofourproblemaretheparametersofthemodelsandthemarginalprobabilitiespki.Therearetwowaystoestimatethemarginals:a)themeaneldapproximation( Zhang 1992 )andb)GaussMarkovMeasureField(GMMF)model( Marroqunetal. 2001 ).Themeaneldapproximationleadstoalgorithmsthatareratherslowandsensitivetonoiseandthe 49

PAGE 50

Marroqunetal. 2001 )leadstosignicantlydierent(inthesenseofentropy)distributionsfromthetrueones. Riveraetal. 2005 )developedaclevertechniquethatcontrolstheentropyofthesolutiondistributionandconstraintsittobeclosertothetruedistribution.Usingtheformulationpresentedby Riveraetal. 2005 ),wecanecientlyestimatetheunknownparametersoftheabovementionedGMMFmodelbyminimizingthefollowingenergyfunctionE(pjc)=XkXi(pki2(logvki)+Xs2Ni(pkipks)2)Xii(1Xkpki)wheretocontrolstheentropyofthemarginals,controlsthesmoothnessofthelabeleldandiareLagrangemultipliersthatwereintroducedinordertoenforcethecondition Riveraetal. 2005 ).Theenergycanbeminimizedusingtheexpectationmaximization(EM)procedureorageneralizedEM( NealandBarry 1998 ).InthissegmentationalgorithmthenumberoflabelsKisnotahiddenvariableandispredened.Thisnumbercanbesetequaltoorgreaterthanthenumberofregionsthataneuroanatomistexpectstondinaparticulardataset. 50

PAGE 51

A)FAmapofasynthetictensoreldsimilartothesynthetictensoreldin 3-2D butwithsmoothlyvaryingtensorswithineachregionanddierentamountsofnoise.FromlefttorighttheSNRis:1)nonoise,2)8.8,3)3.8,and4)2.8.B)Thecorrespondingsegmentationusingpiece-wiseconstantrepresentationoftheregions,C)usingsmoothlyvaryingrepresentations. resultsarealsopresentedtodemonstratetheperformanceofourproposedalgorithmfordiusiontensoreldsegmentationunderdierentamountsofnoiseinthedata.Wesynthesizedseveral2Dsynthetictensoreldsofsize32x32,withdierentshapesofregionsanddierentanatomyofthediusiontensoreldineachregion.Allsynthetictensoreldsweregeneratedbysimulatingthediusion-weightedMRsignal(Equation 1{4 by SodermanandJonsson ( 1995 ))asdescribedearlierinSection 1.5 .Figure 3-2B presentsatensoreldconsistingoftwopiece-wiseconstantregions:1)asmallsquareregioninthecenterofthetensoreldand2)theregionformingtherestofthetensoreld.TheFAis0.6inbothregions.Theprincipaleigenvectorofthediusiontensorsintherstregionhashorizontaldirectionandinthesecondregionhasverticaldirection.Theregionsofthetensoreldingure 3-2C arethesameasingure 3-2C .However,nowthetensoreldwithintheregionsissmoothlyvaryinginsteadofpiece-wiseconstant.TheFAisagain0.6inbothregions.Theprincipaleigenvectorsofthediusiontensorsintherstregion(smallsquare)aretangenttocircleswithcenterintheupperleftcorneroftheimage.Similarly,theprincipaleigenvectorsoftheotherregionaretangenttocircleswith 51

PAGE 52

B C DFigure3-2. A)Theprimaryeigenvectorsofasmoothlyvaryingsynthetictensoreld.B,C,D)Segmentationofdierentsynthetictensorelds. Figure3-3. Illustrationofsegmentationunderpartialvolumingeects.Thetoprowshowsaveragingkernelsatdierentlocationsofthediusion-weightedimages.Thebottomrowdepictsthecorrespondingestimatedtensorelds. centerintheupperrightcorneroftheimage.Finally,gure 3-2D consistsofthefollowingregions:1)aringwithprincipaleigenvectorstangenttocirclescenteredinthelowerleftcorneroftheimage,2)twotriangularregionswithhorizontalprincipaleigenvectorsand3)twotriangularregionswithverticalprincipaleigenvectors.ThelasttworegionshavethesameFAequalto0.6,whilewithintheringtheFAis0.8.Wesegmentedthethreetensoreldsofgure 3-2 usingtheproposedentropycontrolledsegmentationalgorithmdescribedinSection 3.2 .Theparametersofthealgorithmthatweusedare=1and=0:1.Thesegmentationresultsarepresentedinthebackgroundofthevectoreldsofgure 3-2 ,indierentgraylevelsforeachsegmentedregion.Inordertodemonstratethesegmentationperformanceofourproposedalgorithmunderpartialvolumingeectswesynthesizeda2Deldconsistingoftworectangular 52

PAGE 53

3-3 ).Diusion-weightedimagesforthis2DeldwereobtainedbysimulatingtheMRsignal,similarlytopreviousexperiments.Thenthediusion-weightedimageswereaveragedusinga1010kernel,whichproducesasingleaveragemeasurementperimageforevery1010pixels.Thisprocesswasrepeatedfordierentlocationsoftheaveragingkernelandeachtimethecorrespondingtensoreldwasestimatedfromtheaveragedimages(seeillustrationinFigure 3-3 ).Theestimatedtensoreldswerethensegmentedusingourproposedalgorithm.Thesegmentationperformancedenedastheratioofthecorrectlyclassiedareaoverthetotalareawascomputedforthepixelsontheboundarybetweenthetworegions,andfoundtobe0.75.Fortherestofthepixels(non-boundarypixels)thisratiowas1.0.Finallywesynthesizedatensoreldsimilartotheoneingure 3-2D butwithapiece-wisesmoothvariationoftensorswithineachregionandaresolutionof320320,andthenbyfollowingtheaboveaveragingprocesswelteredthetensorelddowntoa)3232andb)6464.Afterthatthetwoobtainedtensoreldsweresegmentedusingouralgorithmandthesegmentationperformancefoundtobe0.92and0.96respectively.Asexpected,inthecaseofthehigherresolution(6464)weobservedbettersegmentationperformance.Inthefollowingexperimentsweusedtherealdiusiontensordatasetfromanisolatedrathippocampus( Shepherdetal. 2006 ),shownearlieringure 2-4 .Figure 3-4 presentssegmentationresultsona2Dsliceselectedarbitrarilyfromthe3Ddataset.Figure 3-4A depictstheFA(fractionalanisotropy)mapsegmentationobtainedusingthescalareldentropy-controlledsegmentationalgorithmpresentedby Riveraetal. 2005 ).SincetheFAmapdoesnotcontainanyinformationabouttheorientationofdiusion,asexpected,itssegmentationyieldserroneousresultswhencomparedtotheexpert'shandsegmentationingure 3-4D ,takenasthegroundtruth.Thisexampledemonstratesthatthetensoreldsegmentationobtainedbysegmentingascalar-valuedfunction 53

PAGE 54

3-4B depictsthesegmentationresultusingapiece-wiseconstantmodelinourtensoreldsegmentationalgorithmdescribedearlierandgure 3-4C depictsthesegmentationresultusingapiece-wisesmoothsegmentationmodel.Finallygure 3-4D isamanuallylabeledsegmentation,basedonexpertknowledgeofhippocampalanatomy.Ingure 3-4E wecompareourresults(showninthebackground)withthemanualsegmentation(shownasanoverlay).Thereweresomedierencesbetweenexpertmanualsegmentation(Figure 3-4D )andtheautomaticsegmentationsoftherathippocampus,althoughautomaticsegmentationbasedonthesmoothlyvaryingtensorsplinemodelbettermatchedtheexpertmanualsegmentationthansegmentationbasedonthepiecewiseconstantmodel.Expertmanualsegmentationoftherathippocampuswasbasedonpreviousknowledgeofhippocampalcytoarchitectureobtainedfrom2-dimensionalvisualizationofthehippocampuswithcontrastgeneratedbyimmunochemistry,varioushistologicalstainingmethodsandtracerstudies( Shepherdetal. 2006 ).Thedierencesappreciatedandboundariesdenotedbytheseoldertechniquesoeronlyindirectorinferentialinformationabouttheorientationsofneuronsandgliawithinthehippocampus.Thus,manualsegmentationmaydiersignicantlyincertainhippocampalregionsfromautomaticsegmentationbasedonthe6-dimensionaltensormodelofbercoherencesintherathippocampus.Itmustbefurtherconsideredthatstudyofthehippocampususingthesenewmethodsmayproducenovelinsightintothecytoarchitectureofdierenthippocampalregionsnotpreviouslyknownbasedonprevioustechniques.SeveralexamplesareshowninFigure 3-4C wheretheautomaticsegmentationusingthesmoothlyvaryingtensorsplinemodelrecognizedsignicanttensordierencesthatsubdividehippocampalregionsintoregionsthatmakesenseinthecontextoftohippocampalconnectivity(pleaseseetheworkby Shepherdetal. 2006 )forindepthdetails),althoughthoseregionsareusuallynotdistinctinhistologicalstains.InFigure 54

PAGE 55

B C D EFigure3-4. 2DsegmentationofanisolatedrathippocampusDT-MRI.A)FAmapsegmentationusingalgorithmby Riveraetal. 2005 ).B)tensoreldusingpiece-wiseconstantmodels.C)tensoreldusingasmoothlyvaryingrepresentationoftheregions.D)Manuallylabeledimagebasedonknowledgeofhippocampalanatomy.Theindexofthelabelscorrespondsto:1)dorsalhippocampalcommissure,2)subiculum,3)alveus,4)stratumoriens,5)stratumradiatum,6)stratumlacunosum-moleculare,7)molecularlayer,8)hilus,X)mixtureofCA3stratumpyramidaleandstratumlucidum,Y)stratumoriensbutambiguous,12)mbria.E)comparisonofourresults(showninthebackground)withthemanualsegmentation(shownasanoverlay). 55

PAGE 56

,themolecularlayerwasfurtherdividedintoinfra-andsuprapyramidalblades('e'end'i'respectively)aswellasthecrest('h')thatcontainsperforantbersenteringthedentategyrusfromtheentorrhinalcortex.Theseregionsdonotdierinhistologicalstainingcontrast,butareoftendividedbasedonspatialinformationbecausetheydohavedierentfunctionalconnectivities.Automaticsegmentationusingthesmoothlyvaryingtensorsplinemodelalsodividedthesubicularcomplexlabeled5and2inFigure 3-4D intosomeofitscomponentpartsincludingpresubiculum('k'),prosubiculum('j')andsubiculumproper('l').ThetensorsplinemodelalsodividedthestratumradiatumintoproximalanddistalCA1regions('c'and'f'respectively)thatreectdierencesinsepto-temporalconnectivitybetweenproximalanddistalCA1pyramidalneurons,aswellasdistinguishingaregionofipsilateralCA3-CA3associationalconnections('a').Incontrast,thesegmentationbasedonthe6-dimensionaltensorhasdicultydistinguishingregionssometimesrecognizedbyhistologystainsthatareconnectedbyacommonberpathway(andthus,waterdiusioncoherence).Forexample,theperforantpathwaycontinuesfromthepresubiculum('k')tothestratumlacunosum-moleculare('g')tosynapseontheterminaldendritesofpyramidalneurons.Theautomaticsegmentationdoesnotdistinguishtheseareasbecausetheysharethesamebercoherenceastheperforantbers.Asdescribedabove,thisperforantpathwayalsobranchesintothecrestofthemolecularlayer('h'),butthereisbettercontrastinthetensorsplinemodelrepresentationoftheorthogonally-orientedbasaldendritesofgranulecells.Anotherexampleofthismergingphenomenaiscreatedbymossybersextendingthroughthehilus('d')ontotheproximalapicaldendritesofCA3pyramidalneuronsinthestratumlucidum('b').Similarregionalmergencebasedonanatomicalconnectivityisnotedbetweentheprosubiculum('j')andthedistalCA1regionofstratumradiatum('f')whicharefunctionallyconnected.Theseinterestingobservationsserveasareasforfutureresearchmorerelatedtodiusiontensorcontrastgenerationandarenotproblemsinherenttothesmoothlyvaryingtensorsplinemodelorautomaticsegmentationemployedherein. 56

PAGE 57

BFigure3-5. A)Aviewofthe3Dsegmentationofanisolatedrathippocampus.B)DierentviewsofthemolecularlayerfromthesegmentationinA. Finally,wepresentresultsofsegmentingthewhole3Dvolumeoftheisolatedrathippocampususingoursegmentationalgorithmwiththesmoothlyvaryingtensorsplinerepresentationoftheregions.Figure 3-5A presentsa3Dviewofmostoftheregionsdetectedbytheproposedalgorithm.Finallygure 3-5B showsdierentviewsofthemolecularlayercontainedinthesegmentationresultshowningure 3-5A .Tothebestofourknowledgethisistherstreportonautomaticsegmentationofanisolatedhippocapmus.Theclinicalsignicanceofanalgorithmforautomaticsegmentationofanisolatedhippocapmushasalreadybeendiscussedearlier. 57

PAGE 58

1 2 and 3 ,a2ndordertensorhascommonlybeenusedtoapproximatethediusivityproleateachimagelatticepointinaDW-MRI( Basseretal. 1994 ).Theapproximateddiusivityfunctionisgivenby Basseretal. 2000 ).Inordertoovercometheabovelimitation,higherordertensorswereintroducedby OzarslanandMareci ( 2003 )torepresentmorecomplexdiusivityproleswhichbetterapproximatethelocaldiusivityfunction.Generalizedscalarquantitiessuchasthevarianceofdiusivityandthegeneralizedanisotropywerederivedasfunctionsofthehigherordertensorcoecients( EvrenOzarslan 2005 ).However,inalloftheseworksthehigher-ordertensorsareestimatedwithoutimposingthepositivityofthediusivityfunctionapproximation,whichissignicantlyimportantsincenegativediusivityvaluesarenonphysical.Inthischapterweproposeanovelparametrizationofthe4th-ordertensorsforimposingthepositivityoftheestimateddiusivity.Inthisparametrizationweexpressthe66matrixG=CCTusingitsIwasawacoordinates( Iwasawa 1949 ).BycombiningtheIwasawaformulationwiththetheoremonternaryquarticsby Hilbert ( 1888 ),weareinterestedinestimatingIwasawadecomposable66symmetricpositivesemi-denitematricesofatmostrank-3.Wederiveformulasforuniquelycomputingthetensorcoecientsasfunctionsoftheparametersofourmodel.Theunknownparametersare 58

PAGE 59

PowersandReznick ( 2000 )andb)allowsthesimultaneousestimationandregularizationofthe4th-ordertensoreld.Themotivationforusinganestimatedpositivehigher-ordertensormodelisthatonecanderivethegeneralizedhigher-orderextensionsofthecommonlyused2nd-orderscalarmeasuresformonitoringvariousbraindiseases.Inourexperimentalresultsweshowthatsuchscalarmeasurescanbemoreaccuratelyestimatedbyusingtheproposedframeworkcomparedtoothermethods.Weshowthatourmodelisrobusttonoiseinestimatingberorientationsaswellasscalarmeasuresderivedfromthe4th-ordertensorcoecientssuchasthegeneralizedtracedenedby EvrenOzarslan ( 2005 ).Thisdemonstratesthattheproposedframeworkisadirectlyapplicabletoolthatextendsecientlytheclinicallyusedmeasuresbyovercomingthelimitationsofthecurrentmodelsinliterature.Therestofthechapterisorganizedasfollows:InSection 4.2 ,wepresentanovelparametrizationofthe4th-ordertensorsthatisusedtoenforcethepositivityoftheestimatedtensorsusingtheIwasawadecompositionoftheGrammatrix.InSection 4.2.1 ,wepresentafunctionalminimizationmethodtoestimate4th-ordertensorsfromdiusion-weightedMRimages.Furthermore,inSection 4.2.2 weproposeadistancemeasureforthespaceof4th-ordertensors,andemployitforestimatingasmoothlyvaryingtensoreldaswellasestimatingthegeneralizedvarianceofthetensors.Section 4.3 containstheexperimentalresultsandcomparisonswithothermethodsusingsimulateddiusionMRIdataandrealMRdatabaseofexcisedratspinalcords. 59

PAGE 60

Ozarslanetal. 2006 ))haveshownthat2nd-ordertensorsfailtomodelcomplexlocalstructuresofthediusivityinrealtissuesandhigher-ordertensorsorothermethodsformulti-berapproximationmustinsteadbeemployed(seediscussioninSection 1.2.4 ).Inthecaseof4th-ordertensors,thediusivityfunctioncanbeexpressedusingthestandardnotationofquarticforms(alsoemployedby PowersandReznick ( 2000 ))asfollows: 4{2 thereisnoguaranteethattheestimatedcoecientsDi;j;kformapositivetensor.Therefore,thereisneedforthedevelopmentofanewparametrizationofthe4th-ordertensor,whichguaranteesthepositivityoftheestimatedtensor.Lettingg1,g2andg3in 4{2 bethevariables,theequivalencebetweensymmetrictensorsandhomogeneouspolynomialsisstraightforwardtoestablish.Moreover,ifasymmetrictensorispositive,thenitscorrespondingpolynomialmustbepositiveforallreal-valuedvariables.Hence,hereweareconcernedwiththepositivedenitenessofhomogenouspolynomialsofdegree4in3variables,orthesocalledternaryquartics.Hilbertin1888provedthefollowingtheoremonternaryquartics(see Rudin ( 2000 ),and Powersetal. 2004 )foramodernexposition)thatwewillemployinthiswork:Theorem1.Thereexistsanidentityd=p21+p22+p33inwhichdisarealternaryformd=d(g1;g2;g3)ofdegreefourwhichispositivesemi-deniteandthepiarequadraticformswithrealcoecients.

PAGE 61

4{2 canbeexpressedasasumof3squaresofquadraticformsas 4{3 ispositivesemi-deniteandhasatmostrank3.HereweshouldemphasizethatGbeingsemi-denitedoesnotimplythatthe4th-ordertensordisalsosemi-denite,sincethevectorspaceoftheformeris6-dimensional(i.e.vectorsv)whilethecorrespondingvectorspaceofthelatteris3-dimensional(i.e.vectorsg).Infact,usingtheparametrizationofEquation 4{3 thewholespaceofthestrictlypositivedeniteternaryquarticsaswellasthosewhicharesemi-deniteisspanned.GivenaGrammatrix,theCCTparametrizationinEquation 4{3 isnotunique,i.e.thereexistdierentmatricesCwhichparameterizethesameGrammatrix.ForexampleCparameterizesthesameGrammatrixasitsantipodalpairC.Furthermore,thereareinnitelymanyorthogonalmatricesthatyieldthesameGwhenappropriatelyappliedintheaforementionedparametrization.Thisisduetotheorthogonalityproperty(OOT=I)oftherotationmatricesO,whereIistheidentitymatrix.Thus,inthecasethatCisofsize63,forany33orthogonalmatrixOwehave(CO)(CO)T=CCT.Thisproducesaninnitesolutionspace,whichtheoreticallycannotbehandledbytheoptimizationtechniques.WeshouldalsonotethatinEquation 4{3 thereareintotal18parametersinC,whichare3morecoecientsthanthenumberofuniquetensorcoecientsinEquation 4{2 .ThenonuniquenessissuesoftheGrammatrixwerealsodiscussedby PowersandReznick ( 2000 ),investigatinghowmanyfundamentallydiernentGrammatricesparametrizethesameternaryquartic. 61

PAGE 62

Iwasawa 1949 ; JianandVemuri 2007a ).EverynnpositivedenitematrixGcanbeuniquelyexpressedusingitsIwasawacomponentsasfollows. 4{4 bysettingrank(W)+rank(V)k.BycomputingthematrixmultiplicationsinEquation 4{4 wederivethefollowingparametrizationofpositivesemi-denitematrices 4{5 andthatinEquation 4{3 bydeningBtobethematrixthatsatisestheequationBAT=XTW.ByusingtheabovewederivethefollowingparametrizationfortheGrammatrix 62

PAGE 63

4{6 isunique.NotethatinEquation 4{6 thereareintotal15parameters,6inAand9inB,whichisequaltothenumberofuniquetensorcoecientsinEquation 4{2 .AccordingtoTheorem1,Equation 4{6 spansthewholespaceofpositivesemi-deniteternaryquarticsbutcanbeusedtodeneaparametrizationforthecaseofstrictlypositivedenite4th-orderternaryquarticsasfollows:Lemma1.Foreveryrealpositivedeniteternaryquarticdthereexistsanarbitrarilysmallpositiverealnumbercsuchthatdcanbewrittenasthesumoftheternaryquarticc(x2+y2+z2)2andthreesquaresofquadraticforms.Proof.Letd(g)beastrictlypositivedeniteternaryquartic.Therefore,d(g)>08g2S2,where,S2istheunitsphere,i.e.thespaceofunitvectorsg.Wedenethefollowingternaryquarticc(g21+g22+g23)2,wherecisanyrealnumberintheinterval00thisternaryquarticisalsopositivedenite.Letusnowdenef(g)=d(g)c(g21+g22+g23)2.Sinceming2S2f(g)=ming2S2d(g)c0,f(g)isapositivesemi-deniteternaryquarticandthereforecanbeexpressedasasumofthreesquaresofquadraticforms(byTheorem1).Thus,everypositivedeniteternaryquarticd(g)canbewrittenasd(g)=f(g)+c(g21+g22+g23)2wherecisanarbitrarilysmallpositiverealnumber,whichprovesthelemma.ThecorrespondingdiusivityfunctioncanbeexpressedusingtheIwasawacoordinatesasfollows: 63

PAGE 64

4{5 ,andAandBaredenedasinEquation 4{6 .Hereweshouldnotethatinpracticeduetoniteprecisioncomputationsccanbesettothesmallestpossiblevalueofaniteprecisionmachineandtherefore,isconsideredasaknownvariable.Inthedouble-precisionoating-pointIEEE754-1985standardthesmallestpositivevalueisapproximatelyc=10308.Furthermore,weshouldemphasizethatthesemi-denitepropertyofGdoesnotnecessarilyimplysemi-denitenessofd(g)sincetheformerissemi-denitein<6whilethelattermaybesemi-denitein<3.Inparticular,usingtheproposedparametrization(Equation 4{7 ),wecomputesemi-denitematricesGthatcorrespondtostrictlypositive-denited(g).GivenaGrammatrix,whichinourcaseisparameterizedusingthematricesAandBinEquation 4{7 ,wecanuniquelycomputethetensorcoecientsDi;j;kbyusingtheformulasintable 4-1 .IntheseformulasthetensorcoecientsareexpressedasfunctionsofthecomponentsofA,Bandthexedparameterc.Therefore,wecanemploytheparametrizationinEquation 4{7 fortheestimationofthecoecientsDi;j;kofthediusiontensorfromMRimagesusingthefollowingtwosteps:1)rstestimatetheAandBmatricesfromthegivenimagesbyusingafunctionalminimizationmethod(minimizationofEquation 4{8 usingtheLavenberg-Marquardtnonlinearleast-squaresmethod),andthen2)computetheuniquecoecientsDi;j;kofthe4th-ordertensorbyusingformulasintable 4-1 .InthefollowingsectionwewillemploythismethodtoenforcethepositivedenitepropertyoftheestimatedfourthorderdiusiontensorsfromthediusionweightedMRimages. 64

PAGE 65

FormulastocomputethetensorcoecientsDi;j;kgiventhecomponentsa11,a22,a33,a21,a31,a32,b11,b12,b13,b21,b22,b23,b31,b32,b33andthexedparameterc=10308. 4{7 ,Misthenumberofthediusionweightedimagesassociatedwithgradientvectorsgiandb-valuesbi;SiisthecorrespondingacquiredsignalandS0isthezerogradientsignal.Usingthemagneticeldgradientdirectionsgiweconstructthe6-dimensionalvectorsvi=[g2i1g2i2g2i3gi1gi2gi1gi3gi2gi3]T.InEquation 4{8 ,the4thorderdiusiontensorisparameterizedusingthe33matricesAandBwhichformtogethertheGrammatrixG.HavingestimatedthematricesAandBthatminimizeEquation 4{8 ,thecoecientsDi;j;kcanbecomputeddirectlyusingtheformulasinTable 4-1 .S0caneitherbeassumedtobeknownorestimatedsimultaneouslywiththecoecientsDi;j;kbyminimizingEquation 4{8 65

PAGE 66

Wangetal. 2004 ).Inourexperimentsweusedtheexponentialmappinga1;1=exp(~a1;1),a2;2=exp(~a2;2)anda3;3=exp(~a3;3)andtherefore,intheminimizationwesolvefor~a1;1,~a2;2and~a3;3insteadofa1;1,a2;2anda3;3.ThetotalnumberofunknownparametersinGinEquation 4{8 are15andtheyarethefollowing:~a1;1,~a2;2,~a3;3,a2;1,a3;1,a3;2andthe9elementsofmatrixB.StartingwithaninitialguessforS0,AandBwecanuseanygradientbasedoptimizationmethodinordertominimizeEquation 4{8 .WeshouldnoteherethattheexponentinEquation 4{8 isintheformofapolynomialandthereforeitsgradientwithrespecttotheunknowncoecientsiseasilyderivedanalytically.GivenAandBateachiterationoftheoptimizationalgorithmwecanupdateS0byagainminimizingEquation 4{8 .ThederivativeofthisequationwithrespecttotheunknownS0is 4{9 equaltozero,wederivethefollowingupdateformulaforS0 66

PAGE 67

4.2.1 wediscussedhowtoestimatethepositive-denite4th-ordertensorsfromDW-MRIdatabyminimizingEquation 4{8 .The4th-ordertensorsarebydenition(Equation 4{2 )smoothlyvaryingfunctions(withinavoxel).Inordertoimposesmoothnessacrosstheimagelatticewecanaddtotheenergyfunctionthefollowingregularizationterm 4{11 weneedtoemployanappropriatedistancemeasurebetweenthetensorsDiandDj.HereweusethenotationDtodenotethesetof15uniquecoecientsDi;j;kofa4th-ordertensor.Wecandeneadistancemeasurebetweenthe4th-orderdiusiontensorsD1andD2bycomputingthenormalizedL2distancebetweenthecorrespondingdiusivityfunctionsd1(g)andd2(g)leadingtotheequation, 4ZS2[d1(g)d2(g)]2dg(4{12) 67

PAGE 68

315[(4;0;0+0;4;0+0;0;4+2;2;0+0;2;2+2;0;2)2+4[(4;0;0+2;2;0)2+(4;0;0+2;0;2)2+(0;4;0+2;2;0)2+(0;4;0+0;2;2)2+(0;0;4+0;2;2)2+(0;0;4+2;0;2)2]+24(24;0;0+20;4;0+20;0;4)6(22;2;0+20;2;2+22;0;2)+2(4;0;0+0;4;0+0;0;4)2+(2;1;1+0;3;1+0;1;3)2+(1;2;1+3;0;1+1;0;3)2+(1;1;2+3;1;0+1;3;0)2+2(3;1;0+1;3;0)2+(3;0;1+1;0;3)2+(0;3;1+0;1;3)2+2(23;1;0+23;0;1+21;3;0+20;3;1+21;0;3+20;1;3)]where,theintegralofEquation 4{12 isoverallunitvectorsg,i.e.,theunitsphereS2andthecoecientsi;j;karecomputedbysubtractingthecoecientsDi;j;kofthetensorD1fromthecorrespondingcoecientsofthetensorD2.Asshownabove,theintegralofEquation 4{12 canbecomputedanalyticallyandtheresultcanbeexpressedasasumofsquaresofthetermsi;j;k.Inthissimpleform,thisdistancemeasurebetween4th-ordertensorscanbeimplementedveryeciently.Notethatthisdistancemeasureisinvarianttorotationsin3-dimensionalspacebecauseitistheL2norm,whichiswellknownforitsinvariancewithrespecttorigidmotions.Furthermore,byusingtheformulasfromTable 4-1 wecanwritetheabovedistanceasafunctionoftheelementsoftheestimatedmatricesAandB.WeshouldnotethattheobtainedregularizationtermisintheformofpolynomialandthereforeitsgradientwithrespecttotheunknownsAandBofEquation 4{8 canbealsocomputedanalytically.Anotherpropertyoftheabovedistancemeasureisthattheaverageelement(meantensor)^DofasetofNtensorsDi,i=1:::NisdenedastheEuclideanaverageofthecorrespondingDi;j;kcoecientsofthetensors.Thispropertycanbeprovedbyverifyingthat^DminimizesthesumofsquaresofdistancesPdist(^D;Di)2.Similarly,itcanbeshownthatgeodesics(shortestpaths)between4th-ordertensorsaredenedasEuclideangeodesics. 68

PAGE 69

EvrenOzarslan ( 2005 ))ofasingledisplacementprobabilityfunctionparameterizedasa4th-ordertensorD.Thevarianceisgivenby 9dist(~D;0)2 EvrenOzarslan ( 2005 )andweuseitinourexperimentalresultsaswell. 1.5 by SodermanandJonsson ( 1995 )(b-value=1250s=mm2,21gradientdirections).Then,weaddeddierentamountsofRicciannoisetothesimulateddatasetandestimatedthe4th-ordertensorsfromthenoisydataby:a)minimizingPMi=1(SiS0exp(bid(gi)))2withoutusingtheproposedparametrizationtoenforcetheSPDconstraint,byemployingthemethodemployedby Ozarslanetal. 2004 )andb)ourmethod,whichguaranteestheSPDpropertyofthetensors.(SiistheMRsignaloftheithimageandS0isthezero-gradientsignal).Studiesonestimatingberorientationsfromthediusivityprolehaveshownthatthepeaksofthediusivityproledonotnecessarilyyieldtheorientationsofthedistinctberbundles( Ozarslanetal. 2006 ).OneshouldinsteademploythedisplacementprobabilityfunctionsgivenbyEquation 1{2 forcomputingtheberorientations.Inour 69

PAGE 70

Comparisonoftheberorientationerrorsfordierentamountsofnoiseinthedata,obtainedbyusing:a)ourparametrizationtoenforcepositivityandb)withoutenforcingpositivityoftheestimatedtensors. experiments,weestimatedthedisplacementprobabilityprolesfromthe4th-ordertensorsusingthemethodthatispresentedinChapter 6 .Then,wecomputedthedisplacementprobabilityfunctionsfromthe4th-ordertensorsestimatedearlierusingthetwodierentmethods,andafterthatwecomputedtheberorientationsfromthemaximaoftheseprobabilityfunctions.Theerrorangles(meanandstandarddeviation)betweengroundtruthandestimatedorientationsfromthetwomethodsfordierentamountofnoiseinthedataareplottedinFigure 4-1 .Here,weshouldemphasizethattheintensityofthesimulatedsignalwassignicantlysmallerwhenthediusiongradientorientationformedasmallanglewiththeberorientationthanthatcorrespondingtogradientorientationsperpendiculartotheber.Thesesmallerintensityvaluesaremostlikelytobecomenegativewhenwecorruptthesignalwithnoise,andasaresultanymethodthatdoesnotimposetheSPDpropertyisaecteddrastically.Asexpected,ourmethodyieldssmallererrorsincomparisonwiththemethodthatdoesnotenforcetheSPDpropertyofthetensors.Whenweincreasetheamountofnoiseinthedata,theerrorsobservedbythelattermethodaresignicantlyincreased,whileourmethoddepictslesssensitivity.ThisdemonstratestheneedforenforcingtheSPDpropertyoftheestimatedtensorsandvalidatestheaccuracyofourproposedmethod. 70

PAGE 71

Comparisonofthestandarddeviationofthegeneralizedtracedenedby EvrenOzarslan ( 2005 ),obtainedbyusing:a)ourparametrizationthatenforcespositivityandb)withoutenforcingpositivityoftheestimatedtensors. Aswasexpected,similarresultswereobservedbycomparingthegeneralizedtrace(denedby EvrenOzarslan ( 2005 ))ofthe4th-orderdiusiontensorsestimatedinthepreviousexperiment.Inthecaseofestimatingthe4th-ordertensorswithoutimposingtheSPDproperty,44%oftheestimatedtensorsyieldednegativegeneralizedtracevalues,whilethevaluescomputedusingtheproposedmethodwereallpositive.Inbothcaseswecomputedthestandarddeviationoftheestimatedgeneralizedtracefordierentamountsofnoiseinthedata(excludingfromourcalculationsthenegativevaluescomputedbythenonSPDmethod).Figure 4-2 showsthestandarddeviationofthegeneralizedtracecomputedusingbothmethods.Byobservingthegureweconcludethattheproposedmethodproducessignicantlymoreaccurateresultseveninthehighernoisecases.Finally,inordertocompareourproposedmethodwithotherexistingtechniquesthatdonotemploy4th-ordertensors,weperformedanotherexperimentusingsyntheticdata.ThedataweregeneratedfordierentamountsofnoisebyfollowingthesamemethodaspreviouslyusingthesimulatedMRsignalofa2-bercrossing.WeestimatedSPD4th-ordertensorsfromthecorruptedsimulatedMRsignalusingourmethodandthencomputedtheberorientationsfromthecorrespondingprobabilityfunctions.For 71

PAGE 72

FiberorientationerrorsfordierentSNRinthedatausingourmethodfortheestimationofpositive4th-ordertensorsandtwootherexistingmethods:1)DOTand2)ODF. comparisonwealsoestimatedtheberorientationsusingtheDOTmethoddescribedby Ozarslanetal. 2006 )andtheODFmethodpresentedby Tuch ( 2004 ).Forallthreemethodswecomputedtheestimatedberorientationerrorsfordierentamountsofnoiseinthedata(showninFigure 4.3.1 ).Theresultsconclusivelydemonstratetheaccuracyofourmethod,showingsmallberorientationerrors('5o)fortypicalamountofnoisewithstd.dev.'0:50:8.Furthermore,byobservingtheplot,wealsoconcludethattheaccuracyofourproposedmethodisveryclosetothatoftheDOTmethodandissignicantlybetterthantheODFmethod.Forlargeramountsofnoiseourmethodyieldedsmallererrorsthanalltheothermethods. 72

PAGE 73

TheelementsofmatricesAandBestimatedbytheproposedmethodaswelltheestimatedS0andgeneralizedanisotropydenedby EvrenOzarslan ( 2005 ). withoutdiusionweightinghad8signalaverages,andeachdiusion-weightedimagehad2averages.First,weestimatedanSPD4th-ordertensoreldbyapplyingtheproposedmethodtotherealDW-MRIdata.Figure 4-4 showstheestimatedelementsofmatricesAandB,whichparameterizetheGrammatrixdiscussedinSection 4.2 .IntherstrowtheestimatedS0andgeneralizedanisotropydenedby EvrenOzarslan ( 2005 )arealsoshown.Theredcolorintheimagesdenotesnegativevalues.Byobservingtheimageswecanseethateachcoecientshowsdierentdetailsoftheunderlyingtissue. 73

PAGE 74

4th-ordertensorsestimatedA)withoutimposingtheSPDpropertyandB)byusingtheproposedmethod.OnthetoprowthecorrespondingestimatedS0imagesareshowncoloredbymappingtheX;Y;ZcoordinatesofthelargestdiusivityorientationtotheR;G;Bcolorcomponents.Inthisregionofinterestweexpectedsingle-lobeddiusivitieswithpeakspredominantlyintheaxialdirection(showninblue). Furthermore,wecomputedthe4th-ordertensoreldrstwithoutimposingthepositive-deniteconstraint.Inordertocomparetheobtainedresultswiththatestimatedbytheproposedalgorithm,inFigure 4-5 weplotthecorrespondingtensorsfromaregionofinterestinthewhitematter.Inthisregionofinterestweexpectedsingle-lobeddiusivitieswithpeakspredominantlyintheaxialdirectionwhichisrepresentedbythebluecolor.TheX;Y;ZcomponentsofthedominantorientationofeachprobabilityproleareassignedtoR;G;B(red,green,blue)componentsofthecolorofeachtensor.Byobservingthisgure,wecansaythatthetensoreldisincoherentifwedonotenforcetheSPDconstraint(Figure 4-5 left)andtheexpectedsingle-lobednatureofthiswhite 74

PAGE 75

Visualizationofthe4th-ordertensoreldestimatedbyapplyingproposedmethodtoarealDW-MRIdatasetfromanexcisedrat'sspinalcord. matterregionislost.Ontheotherhandthetensorsobtainedbyourmethod(Figure 4-5 right)aremorecoherent.NotethatthisisaresultofenforcingtheSPDconstraint,sinceinthisexperimentwedidnotuseanyregularization.Similarilytothesimulateddataexamples(Section 4.3.1 ),theROIinFigure 4-5 correspondstohighlyanisotropicbersinthewhitematterofthespinalcord,whicharemostlikelytoyieldnegativediusivitieswhentheSPDpropertyisnotimposed.ThisdemonstratesthesuperiorperformanceofouralgorithmandmotivatestheuseoftheproposedSPDconstraint.Figure 4-6 depictsavisualizationofthetensoreldestimatedbytheproposedtechnique.Inthisguremulti-lobeddiusivityprolesandothercomplexstructurescanbeobserved.Finally,azoomedinregionofinterestshowssmoothtransitionsfromsingle-lobeddiusivitiestomulti-lobedstructuresintheestimated4th-ordertensoreld. 75

PAGE 76

4 weparametrized4th-orderdiusiontensorsusingtheIwasawadecompositionofpositivedenitematricesandHilbert'stheoremonpositiveternaryquartics.Thismethodguaranteestheestimationofpositivedenite4th-orderdiusivityfunctions.Theestimatedcoecientsofthesehigher-ordertensorsaredirectlyrelatedtousefullyphysicalquantitiessuchaspartialderivativesofthediusivity.Furthermore,thegeneralizedtraceofthetensorscorrespondstothemeandiusivityvalue,whichhasbeenusedinclinicalapplications(seediscussioninSection 1.2.4 ).However,itisnottrivialtoextendthisparametrizationtoordershigherthan4.Todate,noneofthereportedmethodsinliteraturefortheestimationofthecoecientsofhigherordertensors(e.g.6th,8th)preservethepositivedenitenessofthediusivityfunction.Inthischapterwepresentanovelparameterizationofthediusivityfunctionthatguaranteesthepositivedenitepropertyusingthescalarexponentialofahigher-orderhomogeneouspolynomial.Wepresentanecientframeworkforcomputingdistancesandgeodesicsinthespaceofthecoecientsofourproposeddiusivityfunction.Thekeycontributionofthischapteristhatweemploythisframeworkforestimatinghigher(4th,6thand8nd)orderdiusivityapproximationsfromdiusion-weightedMRimages.Notethatweareonlyinterestedforsymmetrictensorsandthereforeweconsideronlyevenorders.WecompareourmethodwithotherexistingDTImethodsshowinghigheciencyoftheproposedmethod.Finally,wevalidateourframeworkusingrealdiusion-weightedMRdatafromexcised,perfusion-xedratopticchiasm. 76

PAGE 77

5{1 ispositiveforanysymmetricmatrix.Forexample,inthecasethatEisthe33zeromatrix,wehaved(g)=e0=18g.IfweusethestandarddiusivityfunctiongTDg(Equation 4{1 ),thepreviousexamplecorrespondstothediusiontensorD=I.InthiscasewehavegTDg=18g.InEquation 5{1 thediusivityfunctionwasdenedbyusinga2ndorderexponentialtensorE.Thisfunctiond(gjE)canbegeneralizedbyusinghigherordertensors.Inthecaseofa4thordersymmetrictensorwehave15uniquecoecientscollectedintoavector~E=(E4;0;0,E0;4;0,E0;4;0,E2;2;0,E0;2;2,E2;0;2,E2;1;1,E1;2;1,E1;1;2,E3;1;0,E3;0;1,E1;3;0,E0;3;1,E1;0;3,E0;1;3).InthecaseofhigherordertensorswewillusethenotationEp1;p2;p3toindicatethatitisthecoecientofthetermg1p1g2p2g3p3.ByusingthisnotationEquation 5{1 canbegeneralizedas 77

PAGE 78

4Z[log(d(gj~E1))log(d(gj~E2))]2dg(5{3)Byanalyticallycomputingtheintegral,Equation 5{3 canbewritteninthefromofsumofsquares,whichisveryfasttocompute.Asanexampleinthecaseof2ndorderEDTsEquation 5{3 canbeevaluatedusing8additionsand10multiplications,inthe4thordercaseusing47add.and33mul..Theanalyticexpressionforthedistanceinthe4th-ordercaseisgiveninEquation 4{13 ,whichhasbeenalsousedintheexperimentspresentedinChapter 4 .Notethatthemetricdenedaboveisrotationinvariantinthecaseofanyorderexponentialtensors.Furthermore,byusingthisdistancemeasureitiseasytoprovethatthemeanelement~EisdenedastheEuclideanaverage(~E1+:::+~EN)=N(orgeometricmeanofd1,...,dN)andthegeodesic(shortestpath)betweentwoelements~E1and~E2isdenedasEuclideangeodesic(t)=(1t)~E1+t~E2,t2[0;1].Inthefollowingsectionweemploythisdistancemeasuretodeneananisotropymapof2;4;6&8th-orderEDTs. 78

PAGE 79

ComparisonbetweenFAand4th-orderfisomapusingratopticchiasmdata.A)FA.B)fiso.C)fisomappedto[0,1]usingthettedfunctionshowninD.D)plotoffisovs.FA. GivenanarbitraryKth-order~E,wecancomputetheclosestisotropytensorcoecients~Eisobyndingthescalarcthatminimizesthedistanceof~Efromtheisotropiccase.Inthe2nd-ordercasec=(E1;1+E2;2+E3;3)=3,inthe4th-ordercasec=(E4;0;0+E0;4;0+E0;0;4+E2;2;0+E0;2;2+E2;0;2)=9,inthe6th-ordercasec=(E6;0;0+E0;6;0+E0;0;6+E4;2;0+E4;0;2+E2;4;0+E0;4;2+E2;0;4+E0;2;4+E2;2;2)=27andinthe8th-ordercasec=(E8;0;0+E0;8;0+E0;0;8+E6;2;0+E6;0;2+E2;6;0+E0;6;2+E2;0;6+E0;2;6+E4;4;0+E0;4;4+E4;0;4+E4;2;2+E2;4;2+E2;2;4)=81.Thefunctionfiso(~E)=dist(~E;~Eiso)mapsthespaceof~Etothespaceofnon-negativerealnumbers.Thesmallerthevalueofthefunction,thecloseris~Etothe~Eiso.Thebehaviorofthefunctionfisoissimilartothatofthewell-knownfractionalanisotropy(FA)mapof2ndorderdiusiontensors.Toillustratethis,weestimatedtheDTeldandtheEDTeldfromarealdataset,andthenwecomputedtheFAandthefisomaprespectively(Figure 5-1 Aand 5-1 B).Furthermore,weplotthefisoasafunctionofFAinFigure 5-1 D.Thesamegurealsocontainstheplot(inred)ofthettedfunction(c)log(1FA)foranestimatedc=0:318.Theinverseoftheabovefunctionis1exp((1=c)fiso(~E))andcanbeusedtomapfisotovaluesintheinterval[0,1]aswasdoneinFigure 5-1 C. 79

PAGE 80

4.3.1 ,thepeaksofthediusivityproledoesnotnecessarilyyieldtheorientationsofthedistinctberpopulations.OneshouldinsteademploythedisplacementprobabilityP(R)whichisgivenbytheFourierintegralinEquation 1{2 .Inthecaseof2ndorderexponentialdiusiontensorsthepeakofthediusivityprolecoincideswiththepeakofthedisplacementprobability.ThereforetheberorientationcanbeestimatedbycomputingtheeigenvectorassociatedwiththelargesteigenvalueofmatrixE.Inhigherordercase,insteadofndingthemaximaofP(R)wecancomputethemaximaoftheexpression:NXi=1Z1E(qqi)exp(2iqqiR)4q2dq,whichapproximatesthereciprocalspaceusingtheicosahedraltessellationoftheunithemisphere,whereqisascalarandqiareunitvectorscorrespondingtothetessellation.Inthecaseofthird-ordertessellationwehaveN=81,whichistheapproximationthatweusedinourexperiments.Bycomputingtheintegralanalyticallyintheaboveapproximationwehave 3=2(5{4) 80

PAGE 81

ComparisonofDTIframeworks Frameworks:Aneinvar.Log-Euc.EDTframework Distancemap0.19sec0.45sec0.04secSmoothing5.65sec0.64sec0.10sec Table5-2. PropertiesofDTIframeworks Properties/FrameworksAneinv.Log-Euc.EDT Ane.InvarianceXRotationInvarianceXXXFastDTIprocessingXXUnconstrainedestimationXUseofhigherordertensorsX where=(2qiR)2,=4td(qi;~E)andtistheeectivediusiontime. 1{4 by SodermanandJonsson ( 1995 ).Inordertocomparethetimeperformanceoftheproposedframeworkwithotherexistingframeworks(aneinvariantby Pennecetal. 2005 ),andLog-Euclideanby Arsignyetal. 2005 ))forprocessingtensoreldswesynthesizedasinglerowofa2nd-orderDTIandEDTeldofsize10000andthenweappliedtwosimplecalculationsoneverypairoftensors:a)computingtheirdistanceandb)smoothingbyndingtheiraverageusingtheabovementionedframeworks.AccordingtothetimesreportedinTable 5-1 ,ourframeworkisthefastest.Inthecaseofsmoothing,itissignicantlyfasterthantheAneinvariantframeworkandasymptoticallyfasterthantheLog-Euclidean.AcomparisonoftheirpropertiesispresentedinTable 5-2 .NotethatonlytheEDTframeworkcanbeusedforhigher-orderapproximations.Forthespecialcaseof4th-orderapproximations,wecompareda)themethodintroducedby Ozarslanetal. 2004 )forcomputinggeneralizeddiusiontensors,b)themethodpresentedinChapter 4 ,andc)theexponentialtensorframeworkdescribed 81

PAGE 82

Comparisonof4th-orderDTIusingpositivityconstraint(seeChapter 4 ),withoutimposingthepositivity( Ozarslanetal. 2004 )andEDTinestimatingberorientationsfordierentSNRinthedata. inthischapter.Figure 5-2 presentscomparisonofthesetechniquesinestimatingberorientationsusingsimulatedMRsignal(Equation 1{4 )fordierentamountsofRicciannoiseinthedata.TheerrorsobservedbyusingtheexponentialtensorframeworkorthemethodinChapter 4 aresignicantlysmallerthanthoseof4th-ordergeneralizedDTI( Ozarslanetal. 2004 ),whichconclusivelyvalidatestheaccuracyofourproposedmethods.Asexpected,thetwoalgorithmspresentedinthisdissertation,havesimilarperformanceregardingtheaccuracyoftheestimates,sincetheybothguaranteethepositive-denitepropertyoftheestimateddiusivityfunctions.TheadvantageofthemethodpresentedinChapter 4 isthatitcomputesthecoecientsofapositivedenite4th-ordertensor,whicharedirectlyrelatedtousefulphysicalquantities,e.g.meandiusivity,whiletheexponentialtensorcoecientsdonothavesuchproperties.Ontheotherhand,theadvantageoftheexponentialtensorframeworkisthatitcanbeusedinordershigherthan4,anditscoecientslieintheEuclideanspace.Furthermore,weestimatedthe4th-orderexponentialtensoreldofadatasetacquiredfromexcised,perfusion-xedratopticchiasm Ozarslanetal. 2006 ).Figure 82

PAGE 83

Displacementprobabilityprolesofa4th-orderEDTeldfromaratopticchiasmdataset( Ozarslanetal. 2006 ).Inthebackgroundthedistancefromtheclosestisotropyfisoisshown. 5-3 showsthedisplacementprobabilityprolescomputedfromtheestimatedeld.Theprobabilityprolesdemonstratethedistinctberorientationsinthecentralregionoftheopticchiasmwheremyelinatedaxonsfromthetwoopticnervescrossoneanothertoreachtheirrespectivecontralateraloptictracts.Theseorientationmapsareconsistentwithotherstudiesonthisanatomicalregionoftheratopticchiasm( Ozarslanetal. 2006 ).Furthermore,thefisomap(Figure 5-1 c)hasslightlybrighterintensitiesinthecentralregion,comparedtotheFAmap(Figure 5-1 a).ThisisbecauseFAuses2nd-orderapproximation,whichfailsinapproximatingthebercrossingsinthisregionandproducesestimationsclosetotheisotropy(lowerintensities).Inourframework,afterhavingestimatedthecoecientvectors~E,wecanusealgorithmsdevelopedforvectoreldprocessinginordertocomputestatistics(average,principalcomponents),interpolateEDTeldsetc.Figure 5-4 showssomeexamplesofprocessesforresolvingbercrossings,interpolatingandcomputingprincipalcomponentsusingtheproposedframework,fortheentireimageofFigure 5-3 83

PAGE 84

Applicationsoftheproposedframeworkfor:resolvingbercrossings(8th-orderexample),geodesicinterpolationandcalculatingPCA(4th-orderexample)fromdatasetofFigure 5-3 84

PAGE 85

Ozarslanetal. 2006 ).Hence,oneshouldinsteademploythedisplacementprobabilityprolesgivenbytheFourierintegral(Equation 1{2 ).However,thecomputationoftheintegralcannotbeperformedanalyticallyandhenceisataskthatinvolvesnumericalintegrationorapproximations.Inordertoavoidtheaforementionedcomputationaleortandthepossibleinaccuraciesintroducedbythisstep,onecandirectlyestimatethedisplacementprobabilityfromthegivenDW-MRIdata.Inthischapterweproposeanovelrepresentationofthedisplacementprobabilityprolebyusingthe4th-orderCartesiantensorbases.4th-ordertensorshavebeenstudiedinChapter 4 showingtheircapabilitytomodelmulti-lobeddiusivityfunctions.Inthemodelthatwewillstudyinthischapter,thelobesoftheprobabilityprole,whichisexpressedintheformofa4th-ordertensor,corresponddirectlytoorientationsofdistinctberdistributionsandthusthereisnoneedtoevaluatetheFourierintegralinEquation 1{2 .Wealsopresentanovelmethodforecientlyestimatingthe15unknowntensorcoecientsofthedisplacementprobabilityfromagivenHARDIdataset.Wecomparetheperformanceofourmethodincomputingaccurateberorientationswithseveralotherexistingtechniques,demonstratingtheeciencyandaccuracyofourmodel.Finally,weemploytheestimateddisplacementprobabilitiesinanasymmetricdiusionframeworkinordertoestimateasymmetriccomplexberstructuresbyinferringinter-voxelinformation. 85

PAGE 86

6{1 areantipodalsymmetric(T(v)=T(v))forevenorders.TheabilityofaCartesiantensorstoapproximatethecomplexgeometryofasphericalfunctionwithmanylobesincreaseswiththeorder.A2nd-ordertensor(usedinDT-MRIprocessingpresentedinChapters 2 and 3 )canonlybeusedforapproximatingantipodalsymmetricsphericalfunctionswithasinglelobe.Forapproximatingfunctionswithmorelobes,higher-ordertensorsarerequired(see4th-orderapproximationofthediusivityfunctioninChapter 4 ).Inthefollowingsections,weemployhigher-ordertensorstomodelthedisplacementprobabilityproleinHARDIdatasets. 6{2 isasphericalfunctionsincetheonlyargumentistheunitvectorr.Inthecaseof4th-ordertensorsthereare15uniqueunknowncoecientsci;j;k(r0)thatneedtobeestimated.However,sinceinourapplication,thegivendataistheDW-MRIsignalandnotthedisplacementprobability,weneedtodeneasetoffunctionsthatapproximatethesignal,andhavethefollowingproperties:a) 86

PAGE 87

Proposedbasesfunctions theirFourierintegralgivesthe4th-ordertensorialbasesfunctions,andb)canbecomputedanalytically.Asetoffunctionsthathavetheabovepropertiesispresentedintherightcolumnoftable1.Thissetoffunctionswasobtainedbytakingallthepossiblecombinationsofpartialderivatives(@=@q1)i(@=@q2)j(@=@q3)kB(q),wherei+j+k=4andB(q)isareal-valuedfunctionin<3denedas 6{3 isexpressedas 6{4 isaCartesiantensorbasisfunction.Hencewecanemploythefunctions 87

PAGE 88

3DplotsofthreeofthebasesfunctionsfromTable1.Thefunctionswereevaluatedforvaryingkqkoveraunitcircleofdirectionsq.Thecirclewasdenedbyxingtheelevationsphericalcoordinateto=3andvaryingazimuth. presentedintable 6-1 tomodelthedisplacementprobabilityproleas 6{5 ,theapproximatedDW-MRIsignalhasasymptoticbehavior(whenkqk!1)similartothatofthecommonlyusedStejskal-TannerEquation 1{1 .Furthermore,bysubstitutingEquation 6{4 into 6{5 thedisplacementprobabilityproleismodeledasa4th-ordertensor(Equation 6{2 ).ByusingknownpropertiesoftheFouriertransform,theintegralinEquation 6{5 canbere-writteninthefollowingmoreconvenientformthatwewilluselaterinSection 6.4 6-1 shows3Dplotsofthreeofthebasesfunctions(@=@q1)i(@=@q2)j(@=@q3)kB(q).Byobservingthegure,wecanseethevariabilityintheshapeofthebasesfunctions(e.g.crossesandpeanut-likeshapes).Thesecharacteristicsofthefunctionsprovidetheabilitytoourmodeltoapproximatecomplexgeometriessuchasbercrossings.InSection 6.4 weemployourproposedmodel(givenbyEquation 6{6 )toapproximatethedisplacementprobabilityprolesfromagivensetofhighangularresolutiondiusion-weightedimages. 88

PAGE 89

6{6 byminimizingthefollowingenergy 6{7 canbeminimizedbysolvinganover-determinedlinearsystem.ThissystemisformedbyconstructinganN-dimensionalvectorSthatconsistsofthesignalvaluesSnandaN15matrixAwhoseentriesarethevaluesofour15dimensionalbasis(@=@q1)i(@=@q2)j(@=@q3)kB(gn)8n=1:::N,where=kqnk=r0.Wenotethatdependsonlyonr0ifthediusion-weightedimageswereacquiredwithaxedb-valuesinceinthiscasekqnkisthesameconstant8n=1:::N.Inourexperimentsweused=0:5.Theestimatedcoecients~ci;j;karethecomponentsofthevectorxinthefollowingover-determinedlinearsystemAx=S,whichcanbesolvedveryeciently.Aftersolvingtheabovesystem,wecancomputethedirectionsofthedistinctberpopulationsbyndingthepeaksofEquation 6{2 .TheprobabilityprolecanbevisualizedbyplottingEquation 6{2 asasphericalfunction(i.e.forallunitvectorsr)asshownintheexperimentalresultsection. 1{4 (b-value=1250s=mm2,81gradient 89

PAGE 90

Probabilityprolesestimatedbyapplyingourmethodtosimulateddataof:A)2-bercrossingbundleandB)corruptedcrossingsfordierentamountsofRicciannoise. directions).TheprobabilityprolesthatwereestimatedfromthisdatasetbyusingourmethodarepresentedinFigure 6-2 Ashowingcorrectberorientations.Furthermore,inordertocompareourproposedmethodwithotherexistingtechniques,weperformedanotherexperimentusingsimulatednoisyMRsignalof2-bercrossingswithdierentamountsofRicciannoise(Figure 6-2 B).Weestimatedthedisplacementprobabilityprolesfromthecorruptedsignalusingourproposedmethodandthefollowingexistingmethods:a)theDOTmethoddescribedby Ozarslanetal. 2006 ),b)theODFmethodpresentedby Tuch ( 2004 )andc)thepositive4th-orderdiusiontensormodelinChapter 4 .Forallmethodswecomputedtheestimatedberorientationerrorsfordierentamountofnoiseinthedata(showninFigure 6-3 ).Theresultsconclusivelydemonstratetheaccuracyofourmethod,showingsmallberorientationerrors(6o)fortypicalamountofnoisewithsignaltonoiseratios(SNR):12.5-16.6.Furthermore,byobservingtheplot,wealsoconcludethattheaccuracyofourproposedmethodisveryclosetothatoftheDOTmethodandisbetterthantheODFmethodforhighernoisecases.Here,weshouldnotethatourmethodestimatedeachprobabilityproleinapproximately2ms(onaPentium2.4GHz)whichdemonstratestheeciencyofourmethod. 90

PAGE 91

FiberorientationerrorsfordierentSNRinthedatausingourmethod(P4)andthreeotherexistingmethods:1)DOT,2)ODFand3)4th-orderDT.IntheexperimentweusedsimulatedMRsignalofa2-bercrossing,whoseprobabilityproleisshowninFigure 6-2 (right). 6-4 showsarealdataexamplefromaxedratspinalcord.Theprotocolthatusedinthisexperimentincludedacquisitionof22imagesusingapulsedgradientspinechopulsesequencewithrepetitiontime(TR)=1.5s,echotime(TE)=27.2ms,bandwidth=30kHz,eld-of-view(FOV)=4:34:3mm.Aftertherstimagesetwascollectedwithoutdiusionweighting(b0s/mm2),21diusion-weightedimagesetswithgradientstrength(G)=664mT/m,gradientduration()=1.5ms,gradientseparation()=17.5msanddiusiontime(T)=17mswerecollected.Theimagewithoutdiusionweightinghad8signalaverages,andeachdiusion-weightedimagehad2averages.Byobservinggure 6-4 itisclearthatthewhitematterprobabilityprolesshowspeaksthatcorrespondtobertractswhich,asexpected,arepredominantlyintheaxialdirectionwhichisrepresentedbythebluecolor.ThisindicatesthatourproposedmethodestimatescorrectberorientationsinrealHARDIdatasets.Thesamegurealsopresents(inthe 91

PAGE 92

Estimatedprobabilityprolesfromrealdataofarat'sxedspinalcord.ThezoomedROIshowssingleberdistributionsinwhitematterandothermorecomplextissuestructures. zoomedplate)regionswheremorecomplexprobabilityproleswereestimatedthatshowtheunderlyingcomplexityofthetissuestructures. Basseretal. 2000 ).Todatetherearenoexistingmethodsinliterature

PAGE 93

Melonakosetal. 2008 ; DericheandDescoteaux 2007 ; Basseretal. 2000 ),inordertoinferthepresenceofasproutingoranti-symmetriccrossingstructures.Inthissectionwepresentanovelmethodforestimatinganintra-voxelasymmetricsphericalfunctionthatcanmodelcomplexlocalberstructuresusinginter-voxelinformation.Thepeaksoftheestimatedsphericalfunctioncorrespondtodirectionsthatpointtodistinctlocalbertractsandareappropriatelydubbedtractosemas.Tractosemaisapointer/signofusedhereforneuraltractsandhasitsrootsintheGreekwordsema(sign).Inourworkhere,weextractaeldoftractosemasfromagiveneldofODFsordisplacementprobabilitiesbyfollowingasymmetricandorientationdependeddiusionofsphericalfunctions.Thekernelthatcontrolsthediusionprocessbetweentwoelements(inourcasesphericalfunctions)isdenedasafunctionoverthespatiallocation(<3)andthedomain(S2unitsphere)ofthetwoelements,whichleadsustothespace(<3S2)(<3S2).WeconstructthediusionkernelasatensorproductofthevonMisesandGaussianprobabilitydistributionsandbyusingitwederiveanupdateformulafortheeldoftractosemaswhichisexpressedintheformofadiscretekernelconvolution.Themaincontributionoftractosemasisthattheycandepictcomplexasymmetricberstructureswithouttheneedforbertracking.Tothebestofourknowledge,itistherstmethodthatestimatesaeldofasymmetricsphericalfunctionsformodelingsplayingbersandotherasymmetricaswellassymmetricstructures.Furthermore,theestimatedeldoftractosemascanbeusedasinputbyanyexistingbertrackingalgorithmforndingberjunctionsandbrancheswithouttheneedformultipleseeds(acommonrequirementinmanyexistingmethods( Melonakosetal. 2008 ; Frimanetal. 2006 ; McGrawetal. 2004 ; Basseretal. 2000 )).Finally,theexperimentalresultsdemonstratetherobustnessandaccuracyofourmodelinestimatingberorientationsinthepresence 93

PAGE 94

1{4 ). 6{8 isexpressedintheformofakernelintegration,wheredist(:)canbeanynormor\edge-stopping"function Blacketal. 1998 ),thekernelK(:)isafunctionofx;y;r;v,andtheintegrationisoverallvectorsyandunitvectorsv.Inourparticularapplication,thekernelisaprobabilityfunctionexpressingtheprobabilityofdiusionbetweentheelementspx(r)andpy(v).Thekernelweseekshouldexhibitthefollowingproperties:a)theprobabilityofdiusionbetweenlocationsxandydecreaseswiththeirdistance,b)theprobabilityofdiusionbetweenorientationsrandvdecreaseswiththeanglebetweenthem,andc)theprobabilityofdiusionislargeratthelocationsalongthemaximaofpx(r).Thesepropertiesaresatisedbysinglepeakeddistributions.Onesuchfunctionusedhereis, (2)3=2ekyxk2 23(6{10) 94

PAGE 95

Syntheticdataexample.A)Simulateddata.B)Theeldofcomputedtractosemas.C)TractosemasinROIundervaryingnoise.D)Plotofberorientationerrors. ThemostnaturalwaytoimposethelasttwopropertiesistoemploythesinglepeakedvonMisesdistributionforbothKorientandKfiber,givenby, 6{10 and 6{11 respectivelycontrolthesharpnessofthekernel.Havingadiscretelatticeofprobabilitiespx(r)theintegralover<3inEquation 6{8 becomessummationoverthelattice.Furthermore,sincetheGaussianpartofthekerneltakesitslargestvaluesintheregionarounditscenter(atlocationx),wecandeneasetN(x)thatcontainsthelatticeindicesintheneighborhoodofx.Furthermore,wediscretizethespaceofunitvectorsbyusinga4thordersubdivisionoftheicosahedraltessellationoftheunitsphere.Byusingtheabovediscretization,Equation 6{8 canbewritteninthefollowingform 6{12 withrespecttopx(r)andsettingitequaltozero,wederivethefollowingupdate 95

PAGE 96

Realhippocampaldata.A)Thedatasetshownin3D(top)andtheregionofinterestshownenlarged(bottom).Therestoftheplatesdepictthedisplacementprobabilityproles(bottom)andtheorientationscorrespondingtotheirmaximashownastubes(top)obtainedbyusing:B)DTI,C)fourthordertensors,andD)tractosemas. formulafortheeldofsphericalfunctions(tractosemas) 6{13 isexpressedintheformofadiscretekernelconvolutionanditisappliediterativelytoallindicesxandvectorsronthediscretizedS2.ThismethodproducesveryecientimplementationssinceonlykernelmultiplicationsareinvolvedintheevaluationofEquation 6{13 ,whichisafullyparallelizableprocess.Furthermore,onlyfewiterations(2to3)arerequiredtoobservevisuallythediusedasymmetrictractosemas.Finally,choosingadierentdist(e.g.L1norm),wouldleadtomoreanisotropicsolutions. 96

PAGE 97

1{4 (b-value=1250s=mm2,81gradientdirections).Afterthat,weestimatedthedisplacementprobabilityeld(Figure 6-5 A)fromthesimulatedsignalbyusingthemethodin Barmpoutisetal. 2008 )(onecanalsouseanyothermethod).Theaboveobtainedeldofprobabilityfunctionswastheninputtoourproposedmethodforextractingtractosemas(=1;=10,3iterations).Figure 6-5 Bshowstheeldoftractosemascomputedbyourtechnique.Byobservingthegure,wecanseethatourmethodestimatedcorrectlysingleberdistributionsinthelowerpartoftheimageandsplayingbersinthecentralregionoftheeld,whichdemonstratestheeectivenessofourtechnique.NotethesmoothtransitionfromsinglebertosplayingstructureintheROI,andtheexpectedanti-aliasingeectobservedinthevoxelsclosetothesplayingbers.Furthermore,toquantitativelytesttheperformanceofourmethodinestimatingberorientationsweaddedvaryingamountsofRicciannoise(SNRbetween20:1and3.3:1)tothedata.Weappliedourmethodtothesenoisecorrupteddatasetsandthencomputedtheestimatedberorientationerrors.Figure 6-5 Ddepictsaplotofthemeanandthestandarddeviationoftheangleerrorbetweencomputedandgroundtruthorientations(indegrees).Theseresultsvalidatetheaccuracyofourmodelanddemonstrateitsrobustnesstonoise.TheproposedmethodwasalsoappliedtoarealDW-MRIfromanisolatedrathippocampus(Figure 6-6 A).Thedatasetconsistsof22imagesacquiredusingapulsedgradientspinechopulsesequencewithTR=1.5s,TE=28.3ms,G=415mT/m,=2.4ms,=17.8ms,T=17msandb'1250s=mm2. 97

PAGE 98

Theeldoftractosemasestimatedfromthehippocampaldataset.A)Threezoomedvoxelsdepictingthevariabilityintheestimatedstructures.B)Thebersproutingwiththeestimatedtractosemassuperimposed. Figure 6-6 showsaregionofinterest(ROI)inthehippocampuscontainingmixtureofCA3stratumpyramidale,stratumlucidumandpartofthehilus.TherestoftheimagesinthisgureshowacomparisonoftheestimatedlocalberstructuresusingDiusiontensors(order-2DTs),fourthordertensors,andtractosemas.IntheDTeldwecanobservetwodominantorientationsonepointingtotheupperleftandtheothertotheupperrightcorneroftheROI,however,thestructureatthejunctionislost.Thejunctionwasrecoveredusingthefourthordertensorshowever,theydepictthetwoaforementionedberorientationsassymmetricstructures.Thecomplicatedjunctionstructureiscorrectlycapturedintheestimatedeldoftractosemaswithasymmetricstructuresthatdepictsplayingbers.Figure 6-7 depictsbertracksestimatedfromthehippocampaldatasetbyfollowingthepeaksoftractosemas.ThecapabilityoftractosemasincapturingvariousstructuresisdemonstratedontheplateBofthisgure. 98

PAGE 99

4 5 and 6 weemployedsymmetricpositive4th-orderCartesiantensorstoapproximatethelocaldiusivityfunction(seeChapters 4 and 5 )orthelocalwatermoleculedisplacementprobability(seeChapter 6 ).Allthesemethodsproduce4th-ordertensoreldsfromgivenDW-MRIdatasets.Hence,giventwodierentDW-MRIdatasetsdepictingthesameordierentsubjects,onecanregisterthembyusingtheinformationprovidedbythecoecientsDi;j;kofthecorresponding4th-ordertensorelds.Inthischapterwepresenttwodierentmethodsfor4th-ordertensoreldregistration.Therstmethodcomputesalocallyanetransformthatregistersonedatasetwithaxedtargetdataset.Theadvantageofalocallyanetransformationisthatitstoresthelocaltranslationvectorandrotationmatrixofthetransformation.Therotationmatricesareusedtoreorientappropriatelythe4th-ordertensorcoecientsinordertopreservethelocalgeometryofthemicrostructuralpropertiesinthetransformedimage.Thesecondmethodperformsgroup-wiseregistrationofseveraldatasetsandsimultaneouslyestimatestheaveragetensoreld(atlas).Thisalgorithmcomputesdieomorphicdeformationsthatmapeachdatasettotheestimatedatlasusingafunctionalminimizationmethod.Theadvantageofthismethodisthatthecomputedresultisunbiased,i.e.thereisnoxedtargetdatasetandthereforetheresultisnotdependedontheorderofthegivendatasets.Thechapterisstructuredasfollows.Section 7.2 presentsanalgorithmforlocallyaneregistrationof4th-ordertensorelds.InSection 7.3 wereviewtheRiemanniangeometryofpositive-valuedfunctionsandweuseitforgroup-wiseregistrationandatlasconstructionofhigherordertensorelds.BothmethodsaredemonstratedinsyntheticandrealdatasetsandtheresultsarediscussedinSection 7.4 99

PAGE 100

Watson ( 1983 ).ThefamilyofangularcentralGaussiandistributionsontheq-dimensionalsphereSqwithradiusoneisgivenbyp(g)=1 2wheregisa(q+1)-dimensionalunitvector,Tisasymmetricpositive-denitematrixandZq(T)isanormalizingfactor.IntheS2case,gisa3dimensionalunitvector,Z2(T)=4p 4{2 andtheintegralisoverS2(i.e.overallunitvectorsg).TheintegralinEquation( 7{1 )canbeanalyticallycomputedanditcanbewritteninasum-of-squaresformasithasbeenshowninSection 4.2.2 100

PAGE 101

RS2(d1(g))2d2(g) RS2(d2(g))2)2(7{2)HereweusethenotationDtodenotethe15-dimensionalvectorconsistingoftheuniquecoecientsDi;j;kofa4th-ordertensor.Equation 7{2 canalsobeanalyticallyexpressedinasum-of-squaresformanditisinvarianttoscaleandrotationsofthe3Dspace,i.e.thedistancebetweenp1(g)andp2(g)isequaltothedistancebetweenp1(sRg)andp2(sRg),whereRisa33rotationmatrixandsisascaleparameter.Inthefollowingsectionweusetheabovedistancemeasureforregisteringapairofmisaligned4th-ordertensorelds. 7.2.1 ,andtheintegralisoverthe3Dimagedomain.A1I2(x)denotessomehigher-ordertensorre-transformationoperation.Thisoperationappliestheinversetransformationto 101

PAGE 102

Alexanderetal. 2001 )thattheunknowntransformationparameterscanbesuccessfullyestimatedbyapplyingonlytherotationcomponentofthetransformationtothedatasetI2.Thishappensbecauseofthefactthat2nd-ordertensorscanapproximateonlysingleberdistributions,whoseprincipaldirectiontransformationcanbeadequatelyperformedbyapplyingrotationsonlytothetensors.Inthecaseof4th-ordertensors,multipleberdistributionscanberesolvedbyasingletensor,whoserelativeorientationscanalsobeaectedbythedeformationpartoftheappliedtransformation.Therefore,tensorre-orientationisnotmeaningfulforhigher-ordertensorsandinthiscaseatensorre-transformationoperationmustbeperformedinstead,usingthefullanematrixA,whichisdenedinsection 7.2.3 .AnetransformationhasbeenalsousedinDTI;formoredetailsandjusticationofthescheme,thereaderisreferredtotheworkby WangandVemuri ( 2005 ).Equation 7{3 canbeextendedfornon-rigidregistrationof4th-ordertensoreldsbydividingthedomainofimageI1intoNsmallerregionsandthenregisteringeachsmallerregionbyusinganetransformations.Similarmethodhasbeenusedforscalarimageregistrationby Juetal. 1996 )andDTIregistrationby Zhangetal. 2004 ).Theunknowntransformationparameterscanbeestimatedbyminimizing 7{4 canbeecientlyminimizedbyaconjugategradientalgorithmusedinamulti-resolutionframework,similartothatusedby Juetal. 1996 )and Zhangetal. 2004 ). 102

PAGE 103

4{2 canbewrittenasP15n=1Dng1ing2jng3kn.Ifweapplyananetransformationdenedbythe33matrixAtothe3Dspace,thepreviousequationbecomesP15n=1Dn(a1g)in(a2g)jn(a3g)kn,where(a1g)in(a2g)jn(a3g)knisapolynomialoforder4in3variablesg1,g2,g3,andaiistheithrawofA.Inthissummationthereare15suchpolynomialsandeachofthemcanbeexpandedasP15m=1Cm;ng1img2jmg3km,bycomputingthecorrespondingcoecientsCm;nasfunctionsofmatrixA.Forexampleifweusethesamevectorizationaswedidinthepreviousexample,wehaveC1;1=(A1;1)4,etc.Therefore,wecanconstructthe1515matrixC,whoseelementsCm;naresimplefunctionsofA,anduseittodenetheoperationoftransforminga4th-ordertensorDbya3DanetransformationAas 7.2 wepresentedanalgorithmforregisteringtwodatasetsbyxingoneofthemanddeformingtheotheronebyminimizinganappropriatelydeneddissimilaritymeasure.Herewepresentamethodforunbiasedgroup-wiseregistrationofseveral4th-ordertensoreldsandsimultaneouslyestimationoftheatlaseld.ThecalculationsinthisalgorithmarebeingperformedintheRiemannianspaceofrealpositive-valuedfunctions(sincesymmetricpositive4th-ordertensorbelongtothisspace)whichisstudiedinSection 7.3.1 103

PAGE 104

b(7{6)whichsatisesscaleinvariance,i.e.dist(sa;sb)=dist(a;b)8a;b;s2R+,additionallytothepropertiesofdistancemeasures.TheRiemannianmetricinR+canalsobeusedtodenedistancesbetweenn-tupletswhoseelementsarepositiverealnumbers.ThedistancebetweenA=fa1;a2;:::;angandB=fb1;b2;:::;bngai;bi2R+canbedenedasdist2(A;B)=PNi=1dist2(ai;bi).Similarly,thedistancebetweenpositive-valuedfunctionsfa(x)andfb(x)x2canbedenedusingtheRiemannianmetricinR+asdist2(fa;fb)=Rdist2(fa(x);fb(x))dx.Intheparticularcaseofparametricsphericalfunctionsd(gjD1)andd(gjD2),whereg2S2andD1andD2arethecorrespondingparametervectors,thedistanceisgivenby 7{7 isoverS2,i.e.thespaceofunitvectorsg.Theabovedistancefunctionisinvariantwithrespectto3Drotationsandscale,i.e.dist(sRD1;sRD2)=dist(D1;D2)8s2R+andR2SO3.InSection 7.3.2 weemploythedistancemeasureinEquation 7{7 inordertoachievesimultaneousgroup-wiseregistrationandatlasconstructionofeldsofsphericalfunctionsmodeledusingCartesiantensorbases. 6 andthediusivityfunctioninChapter 4 andthekurtosiscomponentofthediusionindiusionkurtosisimages(studiedby Jensenetal. 2005 )). 104

PAGE 105

4 .Thisproduceseldsofpositive-valuedsphericalfunctionswhoseprocessingcanbeachievedusingtheRiemannianmetricpresentedinSection 7.3.1 .Theproblemofgroup-wiseregistrationofNtensor-eldsandsimultaneousatlasestimationcanbeformulatedintoafunctionalminimizationproblem.ByusingEquation 7{7 theenergyfunctiontobeminimizedisgivenby 7{8 canbedenedasR10kvn(x;t)k2dt.Wewillminimizetheenergyfunction(Equation 7{8 )byevolvingthedeformationeldsnusingagreedyiterativeschemewhichapproximatesthesolutiontotheaboveminimizationproblem,similarlytotheworkby Joshietal. 2004 ).ForthispurposewewillconstructaeldofforcesbycomputingthevariationoftherstterminEquation 7{8 105

PAGE 106

7{9 )andrrotisrelatedwiththelocalrotation(i.e.variationof(Rxg)i1(Rxg)j2(Rxg)k3inEquation 7{9 ).ThecomputationofthesetermsisdiscussedinSection 7.3.3 .AftertheestimationoftheeldsofforcesFn,n=1:::Nwecomputetheupdatevectoreldsvn=RK(x)Fn(x)dx,whereKisakernelappliedtotheeldofforces.InourexperimentsweemployedthekernelK(x)=(x)G(x),whereGisaGaussiankernelcenteredatxandisasmoothfunctionthattakeszerovalueattheboundariesandthereforeimposeszeroboundaryconditionstothekernelKsimilarlytotheworkby Caoetal. 2006 ).NotethattheintegrationofKwithFnisaconvolutionthatbecomesmultiplicationinthefrequencydomain,henceitcanbeecientlycomputedusingthediscreteFouriertransform(seetheworkby Joshietal. 1997 )).Then,thedeformationeldsareupdatedasnewn=oldn(x+vn)usingasmallstep.Finally,thetensorcoecientsoftheatlascanbeupdatedbyalsominimizingtherstterminEquation 7{8 withrespecttotheparametersofapositivedenite4th-ordertensorusingtheparametrizationinChapter 4 7{10 cannotbecomputedanalyticallywhentheCartesiantensorparametrizationisusedformodelingthediusivityfunction.Thereforeweapproximatetheintegrationoverthespherebyusingasumoverasetofunitvectorsgmm=1:::Muniformlydistributedonthesphere.Thissetofvectorscanbeconstructedbytessellatingtheicosahedronandthenprojectingthevectorsontheunithemisphere(weconsideronlyahemisphereduetoantipodalsymmetryofdiusivityfunctions).Weusethissetofvectorsinordertoevaluatethespherical 106

PAGE 107

7{11 duringtheminimizationprocedure,sincetheRiemannianmetricemployedinEquation 7{10 involvesonlycomputationsinthelog-space.ThecorrespondingdrivingforcesinEquation 7{10 arenowcomputedasfollows 4 7{8 isintheformofsum(andintegral)ofsquaredistances(L2).Howeveritiswellknownthatthesumofabsolutedistances(L1),ismorerobustto 107

PAGE 108

A)SyntheticallygenerateddatasetbysimulatingtheMRsignalusingEquation 1{4 .B)Datasetgeneratedbyapplyinganon-rigidtransformationtoA.C)Quantitativecomparisonoftheregistrationerrors.TheerrorsweremeasuredbyEquation 7{2 forthewholeeld. variationsinthedata.Inthiscasetheeldofforcesiscomputedasfollows: 7.3.1 asfollows:d(gm)=exp(mediann=1:Nflog(gmjDnn)g).InSection 7.4 wecomparetheatlasescomputedusingtheRiemannianmeanandmedianintermsofaccuracyinestimatingberorientationsunderdierentlevelsofoutliersandvariationsinthedata. 108

PAGE 109

1{4 .Adatasetofsize128128wasgeneratedbysimulatingtwoberbundlescrossingeachother(Figure 7-1 A).Then,anon-rigiddeformationwasrandomlygeneratedasab-splinedisplacementeldandthenappliedtotheoriginaldataset.TheobtaineddatasetisshowninFigure 7-1 B.Inordertocomparetheaccuracyofour4th-ordertensoreldregistrationmethodswithothermethodsthatperformDTIregitrationorregistrationofscalarquantitiescomputedfromtensors(e.g.GA),weregisteredthedatasetofFigure 7-1 AwiththatofFigure 7-1 Bbyperforming:a)GeneralAnisotropy(GA)mapregistrationusingthemethodby Juetal. 1996 ),b)DTIregistrationusingthealgorithmby Zhangetal. 2004 )andc)4th-ordertensorregistrationusingourproposedmethods.Figure 7-1 Cshowsaquantitativecomparisonoftheaboveresultsbymeasuringthedistancebetweenthecorrespondingmisalignedtensorsbyusingthemeasuredenedinsection 7.2.1 .Theresultsconclusivelyvalidatetheaccuracyofourmethodsanddemonstratetheirsuperiorperformancecomparedtotheotherexistingmethods.Wecanalsoseethattheunbiasedregistrationmethodperformsslightlybetterthanthelocallyanetransformationmethod.Furthermore,inordertocomparetheRiemannianmetricspresentedinSection 7.3 intermsofberorientationaccuracyoftheatlasestimatedbyeachmetric,weperformedthefollowingexperiment.Wesynthesizeda2-bercrossingDW-MRIdataset(inasinglevoxel)usingthesamesimulationmethodasabove.Then,wecomputeda4th-ordertensorfromthesyntheticdatasetusingthemethodinChapter 6 anditisshowninFigure 7-2 .Thenwegenerated100moredatasetsbyapplyingsmallrotationstothesimulatedcrossing(fewofthemareshowninFigure 7-2 ).Thearbitraryrotationsweregeneratedbycomputingthematrixexponentialof33skewsymmetricmatriceswhoseelementsweregeneratedusingaGaussiandistributionwithzeromeanandvaryingstandarddeviation. 109

PAGE 110

Comparisonofthe4th-ordertensoratlasescomputedbyvariousmetrics:a)Euclideanmean,b)Riemannianmeanandc)Riemannianmedian(thelasttwowerecomputedinthespacepresentedinSection 7.3.1 ). Finally,wealsoaddedvariousquantitiesof1-beroutliersinthedatasetinordertotesttherobustnessoftheproposedmetrics.Afterthatwecomputedtheatlastensorbyminimizinga)thesumofEuclideansquaredistances,b)thesumofRiemanniansquaredistances,andc)thesumofRiemanniandistances.TheresultsareshowninthebarchartofFigure 7-2 .Asexpected,theRiemannianmeanoutperformstheEuclideanmeansincethephysicalspaceofthedataisthespaceofpositivereal-valuedfunctions,wheretheRiemannianmetricwasconstructed.Furthermore,asitiswidelyknown,themedianismorerobusttothepresenceofoutliersinthedatacomparedtothemean. 7-3 showstwoS0imagesfromthe3Dvolumesoftwohippocampi(A,B),andtwo"checkers"imagesshowingthedatasetsbefore(C)andafterlocallyaneregistration(D).A"checkers"imageisawaytodisplaytwoimagesatthesametime,presentingoneimageinthehalfboxes,andtheotherin 110

PAGE 111

AandB)S0imagesfromtwoHARDIvolumesofhumanhippocampi.C,D)datasetsbeforeandafterlocallyaneregistration.TensorsfromtheROIinDshowingcrossings. therestoftheboxes.Basedonknowledgeofhippocampalanatomy,bercrossingsareobservedinseveralhippocampalregionssuchasCA3stratumpyramidaleandstratumlucidum.Therefore,oneshouldemployoneofour4th-ordertensormethodsinsteadofDTIregistration.ByobservingFigure 7-3 dallthehippocampalregionsweresuccessfullyallignedbyourmethod,transformingappropriatelythebercrossingsFigure 7-3 (right).Finallyweappliedourgroup-wise4th-ordertensoreldregistrationalgorithminordertoco-registerallthe4datasetsandcomputetheatlasusingtheRiemannianmedian.Figure 7-4 showstheinitialmisalignmentoftheS0imagesandthedirectionoflargestdiusionbefore(toprow)andaftergroup-wiseregistrationusingourmethod(bottomrow).Byobservingtheimageswecanseethatalldatasetswerecorrectlyalignedandtheprimarydirectionsofdiusionsbecomemorecoherent. 111

PAGE 112

OverlayedS0imagesandthecorrespondingRiemannianmedianatlas(directionoflargestdiusionisshown)before(toprow)andafter(bottomrow)group-wise4th-ordertensoreldregistration. 112

PAGE 113

8-1 .Thegoalinthisclinicalapplicationistoperformunsupervisedclusteringofthedatasetsandcomparequantitativecomparisons. 4 totheDW-MRimages.Then,weregisteredtheestimatedtensoreldsusingthe4th-ordertensor-eldgroup-wiseregistrationmethodpresentedinChapter 7 .Figure 8-2 showstwocorrespondingslicesofthecontrolandaninjuredcorddataset.Inthisgure,theorientationofthepeaksintheestimated4th-ordertensorsareshownastubes.BycomparingthecorrespondingS0images(shownonthetopofFigure 8-2 )onecannoteasilyobservechangesintheunderlyingberstructuresduetotheinjury.However,bycomparingthecorrespondingberorientationswecanobservesignicantchangesinthewhitematterregionbetweenthehorns(markedwitharedcircle).Figure 8-3 (left)depictsa3Dvisualizationoftheregionofinterest(ROI)inthewhitematterregionbetweenthehorns(showninpink).InordertoperfromvariousquantitativecomparisonsintheROI,wecomputedtheaveragenormalizeddiusivityastheminimizeroftheHelinger'sdistancebetweenthe4th-ordertensorsintheROI.Inthisregion,we Figure8-1. TheacquiredS0imageofacontrol(left)andthreeinjuredratspinalcords. 113

PAGE 114

Comparisonoftheberorientationsestimatedinthecontrolandthecorrespondingregisteredinjuredcorddataset.TheS0imagesareshownonthetopofthegure. expectedtondanisotropicdiusivitieswithberorientationspredominantintheaxialdirection.ByobservingtheorientationsinFigure 8-2 ,onecanseethattheestimateddiusivitiesfromtheinjureddatasetarelesscoherentthanthoseinthecontrolspinalcord.InordertocomparetheseorientationplotswecomputedtheanglebetweenthepeaksoftheaveragediusivityintheROIandtheresultsareshowninFigure 8-3 .Inthisplot,thereisaclearreductionoftheberorientationangleinallinjuredspinalcordcasesduetothechangesintheunderlyingstructurescausedbytheinjury.Theintra-voxelvarianceswerealsoestimatedintheROIandthehistogramsoftheobtainedvaluesareshownontherightplotofthesamegure.Thisplotalsoshowsareductionofthevarianceaftertheinjurywhichcorrespondstoadropoftheanisotropycausedbytheinjury.Inalltheaboveexperimentswederivedeitherascalaroranorientationalquantityinordertocomparethespinalcorddatasets.However,onecanusealltheinformationincludedinthe15real-valuedparametersofourmodel,whichfullycharacterizethe 114

PAGE 115

VisualizationoftheROI(showninpink)in3D.TheplotsshowcomparisonsbetweentheberorientationangleoftheaveragediusivityintheROIandthehistogramofvariancesintheROI. corresponding4th-ordertensors.IntheROI,wetreatedthe15coecientsaselementsof<15andweconstructeda1515covariancematrixforeachdataset.Thesematricesaresymmetricandpositive-deniteandtherefore,wecanemploytheRiemannianmetricofP15 2005 )inordertocomputedistancesbetweenthedatasetsandvariousstatistics(e.g.PrincipalGeodesicanalysis FletcherandJoshi ( 2004 )).Figure 8-4 showsthetableofRiemanniandistancesbetweenthedatasets.ThistablewasthenusedintheAglomerativeclusteringtechnique Dudaetal. 2001 ); Gil-Garciaetal. 2006 )toconstructthehierarhicaltreeshownontherightofthesamegure.Inthisplot,thedistancebetweenthebranchesshowstheanitiesbetweenthegivendatasets.InthiscasethethreeinjureddatasetswereclusteredtogetherwhoseRiemanniandistancefromthecontroldatasetissignicantlylargerthanthedistancebetweentheinjureddatasets.Finally,weextractedtractosemasfromthecontrolandinjuredratspinalcorddatasetsusingthetechniquepresentedinChapter 6 .Figure 8-5 showstheCornuPosteriusregioninoneofthecontrolandoneoftheinjuredspinalcords.Avarietyofdierentberstructuresareshown(singlebundles,crossings,branchings).Inordertocomparetheestimatedstructuresweplottedtheaveragepercentageoftractosemaswith1,2,3...peaksfoundinallcontrolandinjuredsets.Byobservingtheplot,onecanseequantitativedierencesincontrolandinjuredcordsintermsofthedistributionofthenumber 115

PAGE 116

QuantitativecomparisonoftheratspinalcorddatasetusingtheRiemannianmetricofthe1515positivedenitematrices.ThetableshowstheRiemanniandistancesbetweenthecovariancematrices.ThecorrespondinghierarchicaldendrogramwascomputedusingtheRiemanniandistances. Figure8-5. Tractosemasextractedfromacontrol(A)andaninjured(B)cat'sspinalcorddataset.C)comparisonofthenumberofpeaksintheestimatedtractosemas. ofpeaks.Thisresultisconsistentwiththendingsinthepreviousexperimentsanddemonstratethatthereisclearreductioninthediusionobservedinthedatasetsforminjuredspinalcords. 116

PAGE 117

4 and 5 )thatthetensoreldestimationmethodswhichimposethepositivityconstraintontheestimatedtensorsproducemoreaccurateresultsintermsofberorientationerrors.InChapter 4 anovelparameterizationthathandlesthepositivitypropertywaspresentedusingHilbert'stheoreminternaryquartics.AdierentapproachthathandlesthisissuewaspresentedinChapter 5 .TheadvantageoftheformermethodisthattheestimatedtensorcoecientsaredirectlyrelatedtophysicalquantitiesinDW-MRI(e.g.partialderivativesofthediusivity),whilethecoecientsestimatedbythelatterdoesnothavesuchaphysicalinterpretation.Ontheotherhand,thelattercanalsobeusedforhigherorders(6th,8thetc.)anditscoecientsbelongtoahigh-dimensionalEuclideanspaceandthereforecanbeprocessedusingstandardvectorprocessingtechniques.Furthermore,inChapter 6 ,4th-ordertensorbaseswereemployed 117

PAGE 118

7 .Inthismethod,eachtensoreldisdeformedinordertomatchwiththerestofthetensorelds.AttheendofthisprocesstheaveragetensoreldisalsoestimatedusingaRiemannianmetriconthespaceofpositive-valuedsphericalfunctions.InChapter 8 aclinicalapplicationofthistensoreldprocessingframeworkwaspresented.TheresultsdemonstratethattherearesignicantchangesinthediusivitypatternsobservedininjuredratspinalcordDW-MRIdatasetscomparedtothoseincontroldatasets.Thesechangesaredepictedintheestimatedtensorcoecientsandaectthecorrespondingberorientationsintheinjuredcorddatasets.Anunsupervisedclassicationalgorithmcanusethesefeaturesinordertogroupthedatasetsintocategorieswithdierenttypesanddegreesofinjury.InmyfutureworkonDW-MRIprocessing,Iwillemploythemethodspresentedinthisdissertationinordertoconstructanagespecichumanbrainatlasthatcanbeusedindetectingstructuralabnormalitiescausedduetoadiseaseorinjury.Furthermore,thisframeworkcanassistinmonitoringandtreatingvariousdiseasessuchasepilepsy,Alzheimer'sdiseaseandotherdementia. 118

PAGE 119

Alexander,D.C.,Pierpaoli,C.,Basser,P.J.,andGee,J.C.(2001).Spatialtransformationsofdiusiontensormagneticresonanceimages.TMI,20(11),1131{1139. Amaral,D.andWitter,M.(1995).Hippocampalformation.TheRatNervousSystem,pp.443{493. Arsigny,V.,Fillard,P.,Pennec,X.,andAyache,N.(2005).Fastandsimplecalculusontensorsinthelog-Euclideanframework.InProceedingsofMICCAI,pp.259{267. Barmpoutis,A.,Vemuri,B.C.,andForder,J.R.(2008).Fastdisplacementprobabilityproleapproximationfromhardiusing4th-ordertensors.InIEEEInternationalSymposiumonBiomedicalImaging,pp.911{914. Barr,A.H.,Currin,B.,Gabriel,S.,andHughes,J.F.(1992).Smoothinterpolationoforientationswithangularvelocityconstraintsusingquaternions.SIGGRAPH,26(2),313{320. Basser,P.,Mattiello,J.,andLebihan,D.(1994).Estimationoftheeectiveself-diusiontensorfromtheNMRspinecho.J.Magn.Reson.B,103,247{254. Basser,P.J.andPajevic,S.(2003).Anormaldistributionfortensor-valuedrandomvariables:ApplicationstodiusiontensorMRI.IEEETrans.onMedicalImaging,22,785{794. Basser,P.J.andPajevic,S.(2007).Spectraldecompositionofa4th-ordercovariancetensor:ApplicationstodiusiontensorMRI.SignalProcessing,87,220{236. Basser,P.J.andPierpaoli,C.(1996).Microstructuralandphysiologicalfeaturesoftissueselucidatedbyquantitative-diusion-tensormri.JMagnResonB.,111(3),209{19. Basser,P.J.,Pajevic,S.,Pierpaoli,C.,Duda,J.,andAldroubi,A.(2000).Invivobertractographyusingdt-mridata.MagneticResonanceinMedicine,44(4),625{632. Batchelor,P.G.,Moakher,M.,Atkinson,D.,Calamante,F.,andConnelly,A.(2005).Arigorousframeworkfordiusiontensorcalculus.MRM,53,221{225. Black,M.J.,Sapiro,G.,Marimont,D.H.,andHeeger,D.(1998).Robustanisotropicdiusion.IEEETransactionsonImageProcessing,7(3),421{432. Callaghan,P.T.(1991).PrinciplesofNuclearMagneticResonanceMicroscopy.ClarendonPress,Oxford. Cao,Y.,Miller,M.,Mori,S.,Winslow,R.,andYounes,L.(2006).Dieomorphicmatchingofdiusiontensorimages.InComputerVisionandPatternRecognitionWorkshop,2006.CVPRW'06.Conferenceon,pp.67{67. 119

PAGE 120

deBoor,C.(1972).Oncalculatingwithb-splines.J.Approx.Theory,6,50{62. Deriche,R.andDescoteaux,M.(2007).Splittingtrackingthroughcrossingbers:Multidirectionalq-balltracking.InISBI,pp.756{759. Descoteaux,M.,Angelino,E.,Fitzgibbons,S.,andDeriche,R.(2007).Regularized,fastandrobustanalyticalq-ballimaging.MRM,58,497{510. Duda,R.O.,Hart,P.E.,andStork,D.G.(2001).PatternClassication.Wiley. EvrenOzarslan,BabaCVemuri,T.H.M.(2005).Generalizedscalarmeasuresfordiusionmriusingtrace,variance,andentropy.MagnResonMed.,53(4),866{76. Feddern,C.,Weickert,J.,andBurgeth,B.(2003).Level-setmethodsfortensor-valuedimages.InIEEEWorkshoponVariational,GeometricandLevelSetMethodsinComputerVision,pp.65{72. Fletcher,P.andJoshi,S.(2004).Principalgeodesicanalysisonsymmetricspaces:Statisticsofdiusiontensors.InComputerVisionandMathematicalMethodsinMedicalandBiomedicalImageAnalysis,pp.87{98. Frank,L.R.(2002).Characterizationofanisotropyinhighangularresolutiondiusion-weightedMRI.MagnResonMed,47(6),1083{1099. Friman,O.,Farneback,G.,andWestin,C.-F.(2006).ABayesianapproachforstochasticwhitemattertractography.TMI,25(8),965{978. Gil-Garcia,R.J.,Badia-Contelles,J.M.,andPons-Porrata,A.(2006).Ageneralframeworkforagglomerativehierarchicalclusteringalgorithms.PatternRecognition,2006.ICPR2006.18thInternationalConferenceon,2,569{572. Guimond,A.,Guttmann,C.R.G.,Wareld,S.K.,andWestin,S.F.(2002).DeformableregistrationofDT-MRIdatabasedontransformationinvarianttensorcharacteristics.InISBI,pp.1{4. Hasan,K.M.,Gupta,R.K.,Santos,R.M.,Wolinsky,J.S.,andNarayana,P.A.(2005).Diusiontensorfractionalanisotropyofthenormal-appearingsevensegmentsofthecorpuscallosuminhealthyadultsandrelapsing-remittingmultiplesclerosispatients.JournalofMagneticResonanceImaging,21(6),735{743. Helgason,S.(2001).Dierentialgeometry,Liegroups,andsymmetricspaces.AmericanMathematicalSociety. Hilbert,D.(1888).Uberdiedarstellungdeniterformenalssummevonformenquadraten.Math.Ann.,32,342{350. 120

PAGE 121

Ito,M.,Watanabe,H.,Kawai,Y.,Atsuta,N.,Tanaka,F.,Naganawa,S.,Fukatsu,H.,andSobue,G.(2007).Usefulnessofcombinedfractionalanisotropyandapparentdiusioncoecientvaluesfordetectionofinvolvementinmultiplesystematrophy.JournalofNeurology,Neurosurgery,andPsychiatry,78,722{728. Iwasawa,K.(1949).Onsometypesoftopologicalgroups.TheAnnalsofMathematics,50(3),507{558. Jensen,J.H.,Helpern,J.A.,Ramani,A.,Lu,H.,andK.,K.(2005).Diusionalkurtosisimaging:Thequanticationofnon-gaussianwaterdiusionbymeansofmagneticresonanceimaging.MRM,53(6),1432{1440. Jian,B.andVemuri,B.C.(2007a).Metriclearningusingiwasawadecomposition.InIEEEInternationalConferenceonComputerVision,pp.1{6. Jian,B.andVemuri,B.C.(2007b).Multi-berreconstructionfromdiusionMRIusingmixtureofWishartsandsparsedeconvolution.InIPMI,pp.384{395. Jian,B.andVemuri,B.C.(2007c).AuniedcomputationalframeworkfordeconvolutiontoreconstructmultiplebersfromdiusionweightedMRI.TMI,26(11),1464{1471. Joshi,S.,Grenander,U.,andMiller,M.(1997).Onthegeometryandshapeofbrainsub-manifolds.InternationalJournalofPatternRecognitionandArticialIntelligence:SpecialIssueonProcessingofMRImagesoftheHuman,11(8),13171343. Joshi,S.,Davis,B.,Jomier,A.,andG.,G.(2004).Unbiaseddieomorphicatlasconstructionforcomputationalanatomy.NeuroImage,23,151{160. Ju,S.X.,Black,M.J.,andJepson,A.D.(1996).Skinandbones:Multi-layer,locallyane,opticalowregularizationwithtransparency.InCVPR,pp.307{314. Lenglet,C.,Rousson,M.,Deriche,R.,andFaugeras,O.(2004).Statisticsonmultivariatenormaldistributions:AgeometricapproachanditsapplicationtodiusiontensorMRI.INRIAResearchReport5242. Lenglet,C.,Rousson,M.,andDeriche,R.(2006).DTIsegmentationbystatisticalsurfaceevolution.IEEETrans.Med.Imaging,25(6),685{700. Marroqun,J.L.,Velasco,F.A.,Rivera,M.,andNakamura,M.(2001).Gauss-Markovmeasureeldmodelsforlow-levelvision.IEEETrans.PatternAnal.Mach.Intell.,23(4),337{348. McGraw,T.,Vemuri,B.C.,Chen,Y.,Rao,M.,andMareci,T.(2004).DT-MRIdenoisingandneuronalbertracking.MedIA,8,95{111. 121

PAGE 122

Moakher,M.(2008).Fourth-ordercartesiantensorsoldandnewfacts,notionsapplications.QUARTERLYJOURNALOFMECHANICSANDAPPLIEDMATHE-MATICS,61,181{203. Moakher,M.andNorris,A.N.(2006).Theclosestelastictensorofarbitrarysymmetrytoanelasticitytensoroflowersymmetry.JournalofElasticity,85(3),215{263. Nath,K.,Husain,M.,Trivedi,R.,Kumar,R.,Prasad,K.N.,Rathore,R.,andGupta,R.(2007).Clinicalimplicationsofincreasedfractionalanisotropyinmeningitisassociatedwithbrainabscess.JournalofComputerAssistedTomography,31(6),888{893. Neal,R.andBarry,R.(1998).AviewoftheEMalgorithmthatjustiesincremental,sparse,andothersvariants.InLearninginGraphicalModels,pp.355{368.Dordrecht:KluwerAcademicPublishers. Ozarslan,E.andMareci,T.H.(2003).GeneralizeddiusiontensorimagingandanalyticalrelationshipsbetweenDTIandHARDI.MRM,50(5),955{965. Ozarslan,E.,Vemuri,B.C.,andMareci,T.(2004).Fiberorientationmappingusinggeneralizeddiusiontensorimaging.InISBI,pp.1036{1038. Ozarslan,E.,Shepherd,T.M.,Vemuri,B.C.,Blackband,S.J.,andMareci,T.H.(2006).Resolutionofcomplextissuemicroarchitectureusingthediusionorientationtransform(DOT).NeuroImage,31,1086{1103. Pajevic,S.,Aldroubi,A.,andBasser,P.J.(2002).AcontinuoustensoreldapproximationofdiscreteDT-MRIdataforextractingmicrostructuralandarchitecturalfeaturesoftissue.JMagnReson.,154(1),85{100. Pennec,X.,Fillard,P.,andAyache,N.(2005).ARiemannianframeworkfortensorcomputing.InternationalJournalofComputerVision,65. Powers,V.andReznick,B.(2000).NotestowardsaconstructiveproofofHilbert'stheoremonternaryquartics.Contemp.Math.,272,209{227. Powers,V.,Reznick,B.,Scheiderer,C.,andSottile,F.(2004).AnewapproachtoHilbert'stheoremonternaryquartics.C.R.Acadsci.Paris,339,617{620. Ptak,T.,Sheridan,R.L.,Rhea,J.T.,Gervasini,A.A.,Yun,J.H.,Curran,M.A.,Borszuk,P.,Petrovick,L.,andNovelline,R.A.(2003).Cerebralfractionalanisotropyscoreintraumapatients:anewindicatorofwhitematterinjuryaftertrauma.AJRAmJRoentgenol.,181(5),1401{7. Rivera,M.,Ocegueda,O.,andMarroqun,J.L.(2005).EntropycontrolledGauss-Markovrandommeasureeldmodelsforearlyvision.VLSM,pp.137{148. 122

PAGE 123

Ruiz-Azola,J.,Westin,C.F.,Wareld,S.K.,Alberola,C.,Maier,S.,andKikinis,R.(2002).Nonrigidregistrationof3dtensormedicaldata.Med.Im.Anal.,6,143{161. San-JoseEstepar,R.,Haker,S.,andWestin,C.-F.(2005).Riemannianmeancurvatureow.ISVC,3804,613{620. Shepherd,T.M.,Ozarslan,E.,King,M.A.,Mareci,T.H.,andBlackband,S.J.(2006).Structuralinsightsfromhigh-resolutiondiusiontensorimagingandtractographyoftheisolatedrathippocampus.NeuroImage,32(4),1499{1509. Shepherd,T.M.,King,M.A.,Ozarslan,E.,andBlackband,S.J.(2007).Diusiontensormicroscopyindicatesthecytoarchitecturalbasisfordiusionanisotropyinthehumanhippocampus.AmericanJournalofNeuroradiology,28,958{964. Soderman,O.andJonsson,B.(1995).Restricteddiusionincylindricalgeometry.J.Magn.Reson.,A(117),94{97. Squire,L.,Stark,C.,andClark,R.(2004).Themedialtemporallobe.Annu.Rev.Neurosci.,27,279{306. Stejskal,E.O.andTanner,J.E.(1965).Spindiusionmeasurements:Spinechoesinthepresenceofatime-dependenteldgradient,.ournalofChemicalPhysics,42,288{292. Terras,A.(1988).HarmonicAnalysisonSymmetricSpacesandApplications,volumeII.Springer-Verlag. Tournier,J.D.,Calamante,F.,Gadian,D.G.,andConnelly,A.(2004).DirectestimationoftheberorientationdensityfunctionfromDW-MRIdatausingsphericaldeconvolution.NeuroImage,23(3),1176{1185. Tuch,D.(2004).Q-ballimaging.MagnResonMed,52,1358{1372. Tuch,D.S.,Reese,T.G.,Wiegell,M.R.,andWedeen,V.G.(2003).DiusionMRIofcomplexneuralarchitecture.Neuron,40,885{895. vanGelderen,P.,deVleeschouwer,M.H.M.,DesPres,D.,Pekar,J.,vanZijl,P.C.M.,andMoonen,C.T.W.(1994).Waterdiusionandacutestroke.MagneticResonanceinMedicine,31(2),154{163. vondemHagen,E.A.H.andHenkelman,R.M.(2002).Orientationaldiusionreectsberstructurewithinavoxel.Magn.Reson.Med.,48(3),454{459. Wang,Z.andVemuri,B.C.(2004).Ananeinvarianttensordissimilaritymeasureanditsapplicationstotensor-valuedimagesegmentation.CVPR(1),pp.228{233. Wang,Z.andVemuri,B.C.(2005).DTIsegmentationusinganinformationtheoretictensordissimilaritymeasure.IEEETransactionsonMedicalImaging,24(10),1267{1277. 123

PAGE 124

Wang,Z.,Vemuri,B.C.,Chen,Y.,andMareci,T.H.(2004).Aconstrainedvariationalprinciplefordirectestimationandsmoothingofthediusiontensoreldfromcomplexdwi.IEEETrans.Med.Imaging,23(8),930{939. Ward,P.,Counsell,S.,Allsop,J.,Cowan,F.,Shen,Y.,Edwards,D.,Scia,M.,andRutherford,M.(2006).Reducedfractionalanisotropyondiusiontensormagneticresonanceimagingafterhypoxic-ischemicencephalopathy.Pediatrics,117(4),619{630. Watson,G.S.(1983).StatisticsonSpheres.NewYork:Wiley. Weickert,J.andWelk,M.(2005).TensoreldinterpolationwithPDEs.Preprint142,SaarlandUniversity. Westin,C.-F.,Martin-Fernandez,M.,Alberola-Lopez,C.,Ruiz-Alzola,J.,andKnutsson,H.(2006).Tensoreldregularizationusingnormalizedconvolutionandmarkovrandomeldsinabayesianframework.InJ.WeickertandH.Hagen,editors,VisualizationandImageProcessingofTensorFields,pp.381{398.Springer. Zhang,H.,Yushkevich,P.A.,andGee,J.C.(2004).Registrationofdiusiontensorimages.CVPR,1,842{847. Zhang,H.,Yushkevich,P.A.,Rueckert,D.,andGee,J.C.(2007).Unbiasedwhitematteratlasconstructionusingdiusiontensorimages.InMedicalImageComputingandComputer-AssistedIntervention,pp.211{218. Zhang,J.(1992).ThemeaneldtheoryinEMprocedureforMarkovrandomelds.IEEETransactionsonSignalProcessing,40(10),2570{2583. Zhukov,L.,Museth,K.,Breen,D.,H.,B.A.,andWhitaker,R.(2003).LevelsetsegmentationandmodelingofDT-MRIhumanbraindata.J.Electron.Imag.,12(1),125{133. 124

PAGE 125

AngelosBarmpoutisreceivedhisB.S.degreeincomputerscience(withthehighesthonor)fromtheAristotleUniversityofThessalonikiin2003,M.S.inelectronicsandelectricalengineeringfromtheUniversityofGlasgowin2004,andPh.D.incomputerengineeringfromtheUniversityofFloridain2009.Hehascoauthoredmorethan30journalandconferencepublications.Hiscurrentresearchinterestslieintheareasofbiomedicalimageprocessing,machinelearningandcomputervision.Dr.Barmpoutishasreceivedvariousscholarshipsandawardsduringhisundergraduate,master'sanddoctoralstudies,includingtheAlumniFellowshipfromtheUniversityofFlorida. 125