A geometric approach to image matching and synthesis of diffeomorphic paths

MISSING IMAGE

Material Information

Title:
A geometric approach to image matching and synthesis of diffeomorphic paths
Physical Description:
1 online resource (102 p.)
Language:
english
Creator:
Seo, Dohyung
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Thesis/Dissertation Information

Degree:
Doctorate ( Ph.D.)
Degree Grantor:
University of Florida
Degree Disciplines:
Electrical and Computer Engineering
Committee Chair:
Vemuri, Baba C
Committee Members:
Rangarajan, Anand
Entezari, Alireza
Ho, Jeffrey Yih Chian

Subjects

Subjects / Keywords:
diffeomorphisms -- image -- registration
Electrical and Computer Engineering -- Dissertations, Academic -- UF
Genre:
Electrical and 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:
Image processing and analysis traditionally deals with methods that do not take the geometry of the image surface into account. Using the geometry of the image surface can however prove to be of immense value in many tasks such as image denoising, image matching to name a few. Earlier work reported in the field has shown the importance of using the geometry of the image surface and this thesis takes it even further by using this rich source of information in the task of image matching and classification. The image graph representation provides one with the local geometry and allows one to use the machinery of differential geometry to characterize images intrinsically which is otherwise not possible. This intrinsic framework is utilized in developing a novel image mathcing/registration formulation that is symmetric and parameterization invariant. The formulation is tested on both scalar-valued image data sets as well as to Gaussian mixture fields estimated from high angular resolution diffusion magnetic resonance scans. In the context of image registration, there are several classes of the registration transformation that one can consider. In this work, the transormations considered belong to the class of diffeomorphisms. Here again, geometry of the space of diffeomorphisms is exploited in developing novel algorithms to estimate a smooth path between two given diffeomorphisms. Several applications of this path computation are then shown illustrating the value of the algorithms developed. These range from filling in of video frames to medical applications of cardiac magenetic resonance image analysis of patients with and without stem cell treatment.
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 Dohyung Seo.
Thesis:
Thesis (Ph.D.)--University of Florida, 2013.
Local:
Adviser: Vemuri, Baba C.
Electronic Access:
RESTRICTED TO UF STUDENTS, STAFF, FACULTY, AND ON-CAMPUS USE UNTIL 2013-11-30

Record Information

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


This item is only available as the following downloads:


Full Text

PAGE 1

AGEOMETRICAPPROACHTOIMAGEMATCHINGANDSYNTHESISOF DIFFEOMORPHICPATHS By DOHYUNGSEO ADISSERTATIONPRESENTEDTOTHEGRADUATESCHOOL OFTHEUNIVERSITYOFFLORIDAINPARTIALFULFILLMENT OFTHEREQUIREMENTSFORTHEDEGREEOF DOCTOROFPHILOSOPHY UNIVERSITYOFFLORIDA 2013

PAGE 2

c 2013DohyungSeo 2

PAGE 3

ToJesus 3

PAGE 4

ACKNOWLEDGMENTS Ioweagreatthankyoutomanypeoplewhohelpedmenishmydissertation.First andforemost,IwouldliketothankDr.BabaVemuri,myadvisor,whointroducedme intotheeldofmedicalimageprocessingandguidedmethroughmyPh.D.studyusing histirelessinspirationandrichacademicexperience.IamalsothankfulthatDr.Vemuri offeredmefullresearchassistantshipwithwhichIcouldputallmyfocusontheresearch andnotworryaboutnancialproblems.Dr.Vemuribuilttheexampleofapassionate andseriousscholarforme. IwouldalsoliketoappreciateDr.JeffreyHo,whoisalwaysavailablefordiscussion whenIneedhistechnicalandtheoreticaladvise.Heselesslysharedhisexperience andknowledgesothatIavoidedmakingmanymistakes.Fromhim,Ilearnedthe importanceofmathematicalthinkingandcarefulnessinsolvingaresearchproblem.Itis mylucktohavehimasapatientfriendandknowledgeablementor IamalsothankfultoDr.AnandRangarajanandDr.AlirezaEntezari,foragreeing toserveonmycommitteeandspendingtheirvaluabletimeinreadingandcritiquingthis dissertation. IfeelveryprivilegedtospendtimewithmanywonderfulpeopleintheLabof ComputerVision,GraphicsandMedicalImaging(CVGMI).Iamheartilygratefulto SanthoshKodipakaandAngelosBarmpoutiswhohavebeenextremelyhelpfuland patientwhenIrstjoinedinCVGMIlab.IwanttotakethisopportunitytothankYuchen Xie,TingChen,MeizhuLiu,WenxingYe,GuangCheng,,SileHu,QiDeng,YanDeng, YuanxiangWang,HesamodinSalehianandotherCVGMIlabmembers.Wehaveshared somuchjoytogether. Finallyandmostimportantly,Ithankmyparentsandmywifefortheirunconditional loveandsupport.MygratitudetothemisbeyondwhatIcanexpressinwords.This dissertationisdedicatedtothem. 4

PAGE 5

GenerousnancialsupportforthisdissertationwasinpartprovidedbytheNSF grantIIS0916001andtheNIHgrantNS066340toDr.BabaC.Vemuri. 5

PAGE 6

TABLEOFCONTENTS page ACKNOWLEDGMENTS .................................. 4 LISTOFTABLES ...................................... 8 LISTOFFIGURES ..................................... 10 ABSTRACT ......................................... 12 CHAPTER 1GENERALINTRODUCTION ............................ 14 1.1AGeometricRepresentationofImages ................... 16 1.1.1RiemannianManifolds ......................... 16 1.1.2ImageGraphs .............................. 17 1.1.3RiemannianMetricofanImageGraph ................ 18 1.2Diffeomorphisms ................................ 19 2MATCHINGANDCLASSIFICATIONOFIMAGESUSINGTHEIMAGEGRAPH REPRESENTATION ................................. 20 2.1Introduction ................................... 20 2.2MatchingImageGraphs ............................ 25 2.2.1CostFunctional ............................. 27 2.2.2CostFunctionalfor3-by-3SPDTensorImageGraphsin3D .... 30 2.2.3CostFunctionalforGaussianMixtureFieldGraphsin3D ...... 30 2.2.4DeformationFields ........................... 34 2.2.5DerivativesofCostFunctional ..................... 35 2.2.5.1Derivativesofcostfunctional:3-by-3SPDtensorimage graphs ............................ 36 2.2.5.2Derivativesofcostfunctional:gaussianmixtureeld graphs ............................ 37 2.2.6NumericalMethod ........................... 38 2.3Experiments .................................. 38 2.4Application ................................... 44 2.4.1DiffusionMaps ............................. 45 2.4.2Application:ClassicationofOASISDataset ............ 46 2.4.2.1Datapreparation ....................... 46 2.4.2.2Experimentalresults ..................... 46 2.4.3Application:ClassicationofGMFDataset ............. 51 2.4.3.1Datapreparation ....................... 52 2.4.3.2Results ............................ 52 2.5Conclusion ................................... 53 6

PAGE 7

3ANOVELFRAMEWORKFORCOMPUTINGDIFFEOMORPHICPATHS ... 56 3.1ModelFormulation ............................... 61 3.1.1StepI.theSpaceofDensitiesandDiffeomorphisms ........ 61 3.1.2Step2.PathLiftingtotheSpaceofDiffeomorphisms ........ 63 3.2NumericalMethod ............................... 64 3.3Experiment ................................... 66 3.3.1TheEstimationofMissingVideoSequences ............. 68 3.3.1.1Bendingmotion ....................... 71 3.3.1.2Dancingmotion ....................... 75 3.3.2CardiacMotionAnalysisandItsApplicationtoClassication .... 78 3.3.2.1Datapreparation ....................... 78 3.3.2.2Classifyingischemiapatients ................ 80 3.4Discussion ................................... 82 4GENERALCONCLUSION ............................. 85 APPENDIX ASUPPLEMENTALMATERIAL1FORCHAPTER2 ................ 88 BSUPPLEMENTALMATERIAL2FORCHAPTER2 ................ 89 CSUPPLEMENTALMATERIAL3FORCHAPTER2 ................ 92 DSUPPLEMENTALMATERIAL1FORCHAPTER3 ................ 93 ESUPPLEMENTALMATERIAL2FORCHAPTER3 ................ 94 REFERENCES ....................................... 98 BIOGRAPHICALSKETCH ................................ 102 7

PAGE 8

LISTOFTABLES Table page 2-1Theangleerror(AE)andtheEuclideandistance(ED)betweentheground truthofthedeformationvectoreldsandtheestimatedonesbyourimage graphapproachandthediffeomorphicDemonsarecomputed. ......... 43 2-2Scoresofleave-one-outandfour-foldcross-validationtests.Foursubgroups havebeenrandomlyselected50times. ...................... 49 2-3Scoresofleave-one-outandfour-foldcross-validationtests.bytheDemons algorithm. ....................................... 49 2-4Scoresofleave-one-outandfour-foldcross-validationtestswithdifferentmetrics ontheambientspace(withandwithoutreorientation). .............. 50 2-5Comparisonofclassicationscoreswith4-foldvalidationbetweenclassication methods ........................................ 51 2-6Scoresfromleave-one-outvalidationtestwithtwodifferentmatchingmethods. 52 2-7Averageclassicationscoresoffour-foldvalidationtestfrom50differentcombinations oftestandtrainingdata. ............................... 53 3-1TheaveragedpSNRscoresbetweengroundtruthsandestimatedimagesequences withtheproposedmethodandLDDMM ...................... 69 3-2Computingtimetohaveresultingdiffeomorphicpathswiththeproposedmethod andLDDMM. ..................................... 69 3-3Thesizeofdeformationineachgroupoftestdata.Sizeismeasuredusing theL2normofdisplacementvectors. ....................... 71 3-4TheDicecoefcientsofthecomputedheartphasesusingourmethodand LDDMM. ....................................... 79 3-5ClassicationscorewithSupportVectorMachine ................. 81 8

PAGE 9

LISTOFFIGURES Figure page 1-1Diagrammaticexpressionoftheconditionslistedinthedenitionofthedifferential manifold[ 14 ]. ..................................... 18 2-1Animagegraphanditsparameterization, X fora2Dimage ........... 24 2-2Aberandsectionof3-by-3symmetricpositivedenitetensoreldsin2D[ 4 ] 27 2-3Matchingimagesasimagegraphscanberealizedviaregistrationmap dened betweenthetwosurfaces. ............................. 29 2-4Plotofthevalueofthecostfunctionalvs.thenumberofiteration. ........ 40 2-5Anexampleofinputincomparisonbetweentheproposedimagegraphapproach andtheDiffeomorphicDemons. .......................... 42 2-6Resultingdeformationvectoreldswithintheimagegraphapproachandthe DiffeomorphicDemons. ............................... 43 2-7Theplotsoferrormeasures. ............................ 45 2-8Slicesof3-DMRIimages. .............................. 48 2-9Slicesof3-Dtensor-valuedimagesofventricleswhicharetherstfundamental formsofthegraphsofFig. 2-8 ........................... 48 2-102-Dplotsofdiffusionmaps. ............................. 49 2-11Datapointsprojectedon2-Dplanefrom R s =10 ofdiffusionmaps. ........ 54 3-1Thetwostepsoftheproposedframeworkfortheinterpolationofdiffeomorphisms. ............................................. 60 3-2Theareaofthepolygonisgivenby 1 / 2 | X 30 X 21 | ................ 65 3-3Thediagrammaticviewofcomputingmissingframesinavideoscene. I ( t ) representanimageateachtimepointt. ...................... 69 3-4Theinputandresultingdiffeomorphicpathsoftherotatingcylinderexample. 71 3-5pSNRscoreplot.Redlinesandbluelinesrepresentscoresfromourmethod andfromLDDMMrespectively. ........................... 72 3-6Keyframesinthebendingmotionexamplesanddeformationvectoreldswhich deformthepreviouskeyframestothenextkeyframes. .............. 74 3-7Theground-truthandcomputedsequencesofmissingframesinthevideoof bendingmotion. ................................... 75 9

PAGE 10

3-8Keyframesinthedancingmotionexamplesanddeformationvectoreldswhich deformthepreviouskeyframestothenextkeyframes. .............. 76 3-9Theground-truthandcomputedsequencesofmissingframesinthevideoof dancingmotion. ................................... 78 3-102-Dsliceimagesoftheleftventricularmyocardium(LVM)atend-systole(ES) andend-diastole(ED). ................................ 80 3-11Thediastolicllingmotionofa2-DsliceofLVMbetweenESandED. ..... 81 3-12TheyellowboxesindicatetheROIchosenforthecardiacmotionanalysis. ... 82 E-1AHexahedralcell ................................... 95 10

PAGE 11

AbstractofDissertationPresentedtotheGraduateSchool oftheUniversityofFloridainPartialFulllmentofthe RequirementsfortheDegreeofDoctorofPhilosophy AGEOMETRICAPPROACHTOIMAGEMATCHINGANDSYNTHESISOF DIFFEOMORPHICPATHS By DohyungSeo May2013 Chair:BabaC.Vemuri Major:ElectricalandComputerEngineering Imageprocessingandanalysistraditionallydealswithmethodsthatdonottake thegeometryoftheimagesurfaceintoaccount.Usingthegeometryoftheimage surfacecanhoweverprovetobeofimmensevalueinmanytaskssuchasimage denoising,imagematchingtonameafew.Earlierworkreportedintheeldhasshown theimportanceofusingthegeometryoftheimagesurfaceandthisthesistakesit evenfurtherbyusingthisrichsourceofinformationinthetaskofimagematchingand classication.Theimagegraphrepresentationprovidesonewiththelocalgeometry andallowsonetousethemachineryofdifferentialgeometrytocharacterizeimages intrinsicallywhichisotherwisenotpossible.Thisintrinsicframeworkisutilizedin developinganovelimagemathcing/registrationformulationthatissymmetricand parameterizationinvariant.Theformulationistestedonbothscalar-valuedimage datasetsaswellastoGaussianmixtureeldsestimatedfromhighangularresolution diffusionmagneticresonancescans.Inthecontextofimageregistration,thereare severalclassesoftheregistrationtransformationthatonecanconsider.Inthiswork, thetransormationsconsideredbelongtotheclassofdiffeomorphisms.Hereagain, geometryofthespaceofdiffeomorphismsisexploitedindevelopingnovelalgorithms toestimateasmoothpathbetweentwogivendiffeomorphisms.Severalapplicationsof thispathcomputationarethenshownillustratingthevalueofthealgorithmsdeveloped. 11

PAGE 12

Theserangefromllinginofvideoframestomedicalapplicationsofcardiacmagenetic resonanceimageanalysisofpatientswithandwithoutstemcelltreatment. 12

PAGE 13

CHAPTER1 GENERALINTRODUCTION Geometricapproachesincomputervisionhavebeenstudiedintensively.Their applicationcoverstopicssuchasimagedenoising,imageenhancement,image segmentation,imageregistration,statisticalanalysisandbeyond.Thereasonwhy weareinterestedingeometricapproachesincomputervisionisthattheycanprovide intrinsicinformationaboutobjectsandallowustousethisinformationtobuildanimage frameworkwhichisfundamentallymoresolid.Forexample,Sochen etal. introduced ageometricimagerepresentationforthepurposeofimagedenoising[ 1 ],anditwas thenappliedtoamulti-channelimageandimagesegmentation[ 1 2 4 6 ].Whilethe usualimagerepresentationisafunctionofimagedomain,Sochen etal. ,Kimmel etal. Gur etal. andSeo etal. [ 1 2 4 6 ]describedanimageasaRiemannianmanifold,ora surfaceembeddedinahigherdimension.Unlikeotherdenoisingframeworks[ 7 ],their frameworksallowustoavoidtuningarticialparametersinordertoachievesatisfactory anisotropicdiffusion:Theyutilizedtheintrinsicgeometricinformationofimages.In addition,Kimmel etal. [ 2 ]discussedhowtoemploytheintrinsicinformationforthe purposeofimagesegmentation.Theimagerepresentationintroducedintheseworksis called"imagegraphrepresentation". Anothergeometricapproacheincomputervisioncanbefoundinregistration frameworks.Themainpurposeofregistrationframeworksistosearchfordiffeomorphisms whichassignadisplacedpositioninatarget(orsource)imagetoeachgridpointor pixel/voxelinasource(ortarget)imageorviceversa[ 8 ].Thespaceandgroupsof diffeomorphismsareessentialsubjects.Thisisnotjustbecausediffeomorphisms canpreservethetopologiesofimageswhilematchingasourceandtargetimages. Understandingthegeometricpropertiesofdiffeomorphismsmayprovideuswitha fundamentallysolidframeworkforcomputingdiffeomorphicpaths.Thisresearchis crucialforustostudythespaceofimagesorstatisticalanalysisoflongitudinaldatasets. 13

PAGE 14

Whileplentyofliteraturediscussesthespaceandgroupsofdiffeomorphismsin imageregistrationframeworks,wecanhardlyndaregistrationframeworkthatmodels animageasageometricalobjectdespiteitsfundamentalusefulness.Additionally, eventhoughavarietyofacademicresourcesdealwithdiffeomorphisms,most formulatediffeomorphicpathswithregistrationinasingleformulation[ 9 13 ].These jointframeworksmaybeinefcientsincethesemethodscarryouttwooptimizations simultaneously.Itmightbedesirabletoformulateregistration,anddiffeomorphismsor diffeomorphicpathseparately.However,fewstudiesseparatelyconsidercomputing diffeomorphicpathsinthespaceofdiffeomorphismsthemselves.Therearetwomain reasons:rst,theinnitedimensionalityofdiffeormorphismsandsecond,thelack ofknowledgeofthetruemetricinthatspace.Thisfacthashinderedeffortstostudy longitudinaldatasetsorconductstatisticalanalysisinthespaceofdiffeomorphisms. Thisdissertationdiscusseshowtointroducegeometricalapproachestoanimage registrationframeworkandtoutilizetheminsolvingtheissuesaddressedabove. Thisdissertationiscomposedoftwomainparts.Chapter 2 concernshowwecan equiparegistrationframeworkwithimagegraphrepresentation.Itturnsoutthatthis approachprovidesasymmetric,parametrization-invariantregistrationanduniversal frameworkregardlessofthespaceofimagevalues.Afterwardchapter 3 introduces anovelframeworkforcomputingadiffeomorphicpath.Inthisframework,instead workingdirectlyonthespaceofdiffeomorphisms,weprojectpointsinthespaceonto thespaceofdensities,whichcanbecanonicallyembeddedintoasphereintheHillbert space,andweliftbackthegeodesicconnectingtheprojectedpointstothespaceof diffeomorphisms.Resultsshowthatthisseparatedframeworkisfarmoreefcientthan thejointframeworksincomputationtimeandaccuracy.Additionally,thediffeomorphic pathscomputedinthismethodyieldhighlyreliableresultsinaclassicationproblem. However,beforeproceedingtothemainsubjects,thefollowingsectionsprovide preliminaryknowledgeaboutdifferentialgeometry.Abriefintroductiontodenitions 14

PAGE 15

relatedwithRiemannianmanifoldsareprovided,andamethodforrepresenting animageasaRiemannianmanifoldisintroduced.Then,abriefintroductionto diffeomorphisms(appliedtoimageregistration)isprovided. 1.1AGeometricRepresentationofImages Thepurposeoffollowingsubsectionsistoprovidepreliminaryknowledgeof Riemannianmanifoldsanddifferentialgeometry. 1.1.1RiemannianManifolds Thereareplentyofsourceswhichprovidedenitionsindifferentialgeometry; however,inthisdissertationwetakedenitionsdirectlyfrom[ 14 ].Alsowefollow notationsusedinthesedenitionsintherestofthisdissertation. Differentialmanifold :Adifferentialmanifoldofdimension n isaset M andafamily ifinjectivemappings X : U # R n $ M ofopensets U of R n intoMsuchthat (1) X ( U )= M (2) Foranypair # $ with X ( U ) X ( U )= W % = ,thesets X 1 ( W ) and X 1 ( W ) areopensetsin R n andthemappings X 1 & X aredifferentiable. (3) Thefamily { ( U X ) } ismaximalrelativetotheconditions(1)and(2). Riemannianmetric :ARiemannianmetriconadifferentiablemanifold M isa correspondencewhichassociatestoeachpoint p on M aninnerproduct < > p on thetangentspace T p M whichvariesdifferentiablyinthefollowingsense:If X : U # R n $ M isasystemofcoordinatesaround p ,with X ( % 1 % 2 ,..., % n )= q X ( U ) and # q #$ i = d X q (0,...,1,...0) ,then # # q #$ i # q #$ j $ q = g ij ( % 1 ,..., % n ) isadifferentiablefunction on U .Thefunction g ij iscalledthelocalrepresentationoftheRiemannianmetricin thecoordinatesystem X : U # R n $ M Riemannianmanifold :AdifferentiablemanifoldwithagivenRiemannianmetricis calledaRiemannianmanifold. Volumeform : Vol = % volumeform = % & det ( g ) d % 1 ... d % n where isanimagedomainwhichiscompact,and g istheRiemannianmetricof animagegraph. 15

PAGE 16

Figure1-1. Diagrammaticexpressionoftheconditionslistedinthedenitionofthe differentialmanifold[ 14 ]. Fig. 1-1 showstheconditionslistedinthedenitionofdifferentialmanifold.Inthe rstdenition, ( U or X or ) orthemapping X or iscalledaparameterization(ora systemofcoordinates)of M at p X or ( U or ) 1.1.2ImageGraphs Withthepreliminaryknowledgeprovidedintheprevioussection,wecanliftan imagefunctiontoahigherdimension.Inthissection,just2-Dimagesareconsidered. Moregeneralexamplesarepresentedinchapter 2 .Byndimensionalimagefunction, wemeanthat I : $ imagevalues R n (11) or I = I ( x 1 x 2 ,... x n ) .Now,let'ssupposethatwehavea2-Dgrayscaleimage.Inour geometricframeworkofimages,agrayscaleimageisparameterizedbymapping X :( % 1 % 2 ) $ ( x y I ( x y )) where I ( x y ) isassumeddifferentiableandistheheight functionofintensityvalueat ( x = % 1 y = % 2 ) .Thisliftedimageiscalledanimagegraph. Applyingthismappingtomultichannelimagesin2-D/3-Disstraightforwardandmore rigorousdiscussionaboutimagegraphrepresentationisprovidedinSection2. 16

PAGE 17

1.1.3RiemannianMetricofanImageGraph Onceweconsideranimageasadifferentiablemanifoldwithaparameterization, wecanequipthiswithaRiemannianmetric.Inthecaseof2-Dgrayscaleimages,by denition,theRiemannianmetricisgivenasfollows: ( ) 1+ I 2 x I x I y I x I y 1+ I 2 y + (12) InEq.( 12 ),themetrictodeneinnerproductin R 2 R + aregivenasa3-by-3identity matrix,and I x denotes # I # x .Inthecaseof3-DRGBcolorimages,theRiemannianmetric isgivenasfollows: ( ( ( ( ) 1+ i = { R G B } ( I i x ) 2 i = { R G B } I i x I i y i = { R G B } I i x I i z i = { R G B } I i x I i y 1+ i = { R G B } ( I i y ) 2 i = { R G B } I i y I i z i = { R G B } I i x I i z i = { R G B } I i y I i z 1+ i = { R G B } ( I i z ) 2 + + + + (13) Here,themetricin R 3 R 3 + ischosenasa6-by-6identitymatrix.Inthisdissertation, twokindsofimagegraphswillbeconsidered:thegraphsof3-by-3symmetricpositive denite(SPD)tensoreldsandthegraphsoftheGaussianmixtureeldsin3-D.The featurespaceorthespaceofimagevalueof3-by-3SPDtensorsandtheGaussian mixturewithncomponentsaredenotedas P (3) and GM ( n ) respectively.Then,mapping forthesetwoimagegraphsare X : R 3 $ R 3 P (3) and X : R 3 $ R 3 GM (3) respectively. OncetheRiemannianmeticisevaluated,wecandenethevolumeformofthe manifoldasdenedintheprevioussection.Theintegrationinthedenitionofthe volumeformgivestheareaofa2-Dimagegraphorthevolumeofa3-Dimagegraph. Theimportantfeatureofthisdenitionisthattheintegrationisindependentofthe re-parameterizationofsurfacesaslongasthere-parameterizationpreservesthe orientationoftheimagemanifold. 17

PAGE 18

Inchapter 2 ,wewilldiscusshowtoformulateasymmetricandparametrization invariantregistrationframeworkwiththisgeometricrepresentationofimages. 1.2Diffeomorphisms Themainpurposeofaregistrationframeworkistosearchforatransformation, whichistheoptimizerofthefollowingcostfunctional.Ifwehavetwoimages, I 1 and I 2 to beregisteredtoeachother,generally,thecostfunctionalisgivenasfollows: E = dist ( I 1 I 2 & ()) d + E reg ( & ) (14) isthedomainofthetwoimages,and & ( ) isthetransformationdenedas & : $ E reg isaregularizationtermfor & inEq.( 14 ).Thistransformationfunctiongives correspondencesbetween I 1 and I 2 .Therearetwoconditionsthat & shouldfollow.The twoconditionsareprovidedin[ 8 ]: & shouldnotcreateholds:everypoint & ( ) shouldbetheimageofsome pointin ,i.e., & shouldbeonto. Foldsarealsoprohibited:twodistinctpoints x and x in shouldnotcorrespondto thesamepoint & ,i.e., & mustbeone-to-one. Withthesmoothnessrequirementof & ,wecancall & asdiffeomorphismsdenedas following[ 8 ]. Ahomeomorphismof isacontinuousbijection & : $ suchthatits inverse, & 1 iscontinuos.Adiffeomorphismof isacontinuouslydifferentiable homeomorphism & : $ suchthat & 1 iscontinuouslydifferentiable. Diffeomorphismshavebeenstudiedintensivelyinregistrationliterature[ 8 13 15 ]. However,mostly,researchhasbeencarriedoutwithinregistrationframeworks.Thatis, theevaluationofdiffeomorphicpathshasbeenapartofregistrationframework,andit isformulatedwiththevelocityeldswhicharecharacterizedbytheSobolevnorm.As pointedoutpreviously,fewliteraturediscusseshowtointerpolatediffeomorphismsina frameworkseparatedfromregistration.Thisdissertationproposesanovelframeworkto thisaimandchapter 3 concernsthisnovelframework. 18

PAGE 19

CHAPTER2 MATCHINGANDCLASSIFICATIONOFIMAGESUSINGTHEIMAGEGRAPH REPRESENTATION 2.1Introduction Duringthepastfewdecades,largenumbersofimageregistrationalgorithmshave beendeveloped.Therearemethodsusinglandmarks,contours,surfaces,orvolumes (referencesin[ 16 ]).Morestraightforwardwaysare L 2 basedapproaches.Inthese methods,wesearchforamap whichminimizesanobjectivefunctionalgenerallygiven by % ( I 1 ( x ) ( I 2 ( x & )) 2 d .Themap canbeobtainedbyminimizingtheobjective functionaldirectlywithrespectto withapenaltyon ,[ 16 18 ]orbyevaluatingowsof diffeomorphism[ 9 19 21 ]. However,nomatterwhatkindofregistrationframeworkisconcerned,itisdesirable thatanidealimageregistrationframeworkhasthethreepropertieslistedbelow. Symmetry:Ifweset asaregistrationmapfromasourceimage I 2 toatarget image I 1 ,thecostfunctionalmustsatisfytheequalitygivenby Cost ( I 1 I 2 & )= Cost ( I 2 I 1 & 1 ) ,where Cost ( ) denotesthedesiredcostfunctional.Thismeans thattheestimatedregistrationmapanditsinverseshowsamecorrespondences betweenthetwoimageswhenweswapthesourceandtarget. Parametricinvariance:Thecostfunctionalshouldbeinvariantunderdifferent parameterizationofimagefunctions.Thatis, Cost ( I 1 I 2 )= Cost ( I 1 & I 2 & ) Universality:Itisdesirablethattheframeworkshouldbeapplicabletoany multi-channelimagessuchascolorimagesandtensorialimagesinasingle framework. Symmetryisoneofthemostdesiredpropertiesinimageregistrationalgorithms becauseifthealgorithmisnotsymmetrical,thematchingresultsandresultsfromits application,suchasclassication,maybebiased.Itiseasytoshowthattheusual L 2 basedapproachesarenotsymmetricalinformulationunlesswegivesomeconstraints orredenevolumeform.Symmetricalregistrationframeworksforvarioustypesof imageshavebeeninvestigatedinliteraturesuchas[ 16 18 22 23 ].Christensen etal. havesuggestedasymmetriccostfunctionalgivenby Cost ( I 1 I 2 )= Cost ( I 1 I 2 & 1 )+ 19

PAGE 20

Cost ( I 2 I 1 & 2 ) ,andtheinverseconsistencyconstraintwasimposedtoforce 1 as closeto 1 2 aspossiblein[ 16 ].Whilearegistrationmapestimatedbytheframeworks in[ 16 18 ]isapproximatelysymmetric,theexactsymmetricformulationofthecost functionalhasbeenintroducedin[ 22 ].Intheframework,theLebesguemeasure,which isthestandardvolumeformusedintheusual L 2 basedapproaches,isreplacedwith anovelvolumeformgivenby (1+ J % ) ( .Inthisdenitionofthevolumeform, J % is thedeterminantoftheJacobianmatrixoftheregistrationmap,and ( isthestandard volumeform.Alternatively,Cheng etal. hasintroducedacostfunctionalgivenby Cost ( I 1 & 1 I 2 & 2 ) andregisteredtwoimagesinacommondomain[ 23 ].Themethod issimilartotheatlasconstruction,andthedissimilarityismeasuredinthecommon domain.Inthisframework,theregistrationmapanditsinversearegivenby = 1 & 1 2 and 1 = 2 & 1 Thepropertyofparametricinvarianceisusuallyrequiredinmatching2Dopen/closed surfacesembeddedin3DEuclideanspace.Whileasurfaceisrepresentedaspoints cloud,polygonmeshes,oraparameterizationinsurfacematchingframework,images arerepresentedinimageprocessingasintensityvaluesdenedonpixelsoruniform girdintheimagedomain.Mostimageprocessingframeworksconsideranimageas acontinuousfunctionor I ( x y z ) .However,itisobviousthatresultsfromregistration frameworkswiththisimagefunctionrepresentationarenotconsistentwheninput imagesaredistortedunderre-parametrizationor I ( f x ( r ), f y ( r ), f z ( r )) ,where r =( x y z ) Ingeneral,thisconceptisveryevident,asshownthroughtheseequations. Cost ( I 1 & I 2 & )= F ( I 1 ( ( u )), I 2 ( ( u ))) du = F ( I 1 ( u ), I 2 ( u )) J & 1 d u % = F ( I 1 ( u ), I 2 ( u )) d u = Cost ( I 1 I 2 ) (21) 20

PAGE 21

There-parameterizationisoneofthemust-havepropertiesofshapematching frameworks.Kurtek etal. attacked3Dshapematchingandintroduced q ( map as anovelrepresentationforshapestosolvethere-parametrizationissueresidingin common L 2 based3Dshapematchingframeworks[ 24 25 ].However,iftheambient spaceisnotEuclidean,itisnotclearhowtocompletethecenteringstepinthe constructionofthe q ( map [ 24 ].Moreover,the L 2 distanceof q ( map canonlybe denedinEuclideanspace.Ingeneral,imagevaluesinimageregistrationcanbefound inmorecomplicatedspacethanEuclideanspace.Therefore,the q ( map frameworkis notgoodenoughtobeadoptedforimagematching. Anotherissueinformulatinganimageregistrationframeworkisuniversality.Most imageregistrationframeworksinliteraturearedesignedforspecictypesofimages. Oneofthereasonsisthateachdifferentkindofimagevaluescarriesitsownspecic featuresandasingleframeworkcannottakecareofthedifferences.Forexample, aregistrationframeworkforgrayscaleimagescannotbeappliedtotheregistration problemsoftensorialimages. Consequently,formulatinganimageregistrationframeworkcarryingthesethree propertiesmayrequireadifferentrepresentationofimagesratherthanimagefunction. Inthisdissertation,theimagegraphrepresentationissuggestedastherepresentation ofanimage.Theimagegraphshavebeenusedforanisotropicimagesmoothing,image enhancementandsegmentation[ 1 6 ];however,thereisnoliteraturewhichemploys theimagegraphrepresentationforthepurposeofregistration.Intheimagegraph representation,forexample,a2Dgrayscaleimage, I ,canbeembeddedin R 3 whose maporparameterization, X isgivenby, X :( % 1 % 2 ) $ ( x y I ( x y )) ,where % { 1,2 } are localcoordinatesand ( x y )=( % 1 % 2 ) .Fig. 2-1 showstheparameterization.Wecan extendthistoimageswithmorecomplicatingspaceofimagevaluessuchascolorand tensors. 21

PAGE 22

Figure2-1. Animagegraphanditsparameterization, X fora2Dimage. Inthedifferentialgeometricpointofview,animagegraphisaRiemanniansurface embeddedinaspacewithhigherdimensions.Therefore,theimagegraphcanbe equippedwithaRiemannianmetrictensor.If X denotesaparameterizationofa2D image,itsRiemannianmetrictensorortherstfundamentalformisgivenasfollows. g = ( ) < X $ 1 X $ 1 >< X $ 1 X $ 2 > < X $ 2 X $ 1 >< X $ 2 X $ 2 > + (22) ormoregenerally, g = / ij h ij ) X i ) X j (23) where h ij isthemetricintheembeddingspacewithcomponentindex { i j } ,and { ) ) } arepartialderivativesw.r.t.thelocalcoordinates.OncewehavederivedtheRiemannian tensor,wecandeneavolumeformfromit[ 14 ]. vol = & g ( % 1 % 2 ,..., % n ) d % 1 d % 2 ... d % n (24) where % 'sarelocalcoordinatesand g isthedeterminantoftheRiemannianmetric tensor.Oneoftheimportantpropertiesisthattheintegrationwiththevolumeformis 22

PAGE 23

invariantunderre-parameterization.Thatis, & g ( % 1 ,..., % n ) d % 1 ... d % n = & g ( % 1 ,..., % n ) d % 1 ... d % n (25) Inaddition,thevolumeformisformulatedinaconsistentwayregardlessofthespace ofimagevalues.Therefore,theimagegraphrepresentationischosentosolveallthe issueswhichhavebeendiscussedabove. Chapter ?? ch:part1isorganizedasfollows:Section 2.2 providesabriefdiscussion oftheimagegraphrepresentationandshowshowtoformulatethedesirednovel registrationframeworkwiththeimagerepresentation.Thegeneralformulationofthe costfunctionalisdiscussedinsection 2.2.1 .Asexamples,thissectiondescribes howtoapplythistotheregistrationproblemof3D3-by-3symmetricpositivedenite (SPD)tensorimagesandGaussianmixtureelds(GMF)insection 2.2.2 and 2.2.3 respectively.Inthisdissertation,thedeformationeldsofaregistrationmapare expressedasthecompletesetofcosineandsinefunctionsandnon-linearconjugate gradientmethodisadoptedforthenumericalsolveroftheoptimization.Theseand thederivationofitsgradientdirectionsaredescribedinsection 2.2.4 2.2.5 and 2.2.6 Then,theproposedmethodiscomparedwithDemonsalgorithminsection 2.3 .Inthis comparison,variousdegreesofnoiselevelareaddedtooriginalimages,andknown deformationeldsareappliedtothenoisedimagestodeformtheimages.Evaluated deformationeldsbythetwomethodsarereportedbycomputingtheangleerrorand theEuclideandistance.Insection 2.4 ,thenovelregistrationframeworkisappliedto twoclassicationproblems.Insection 2.4.2 ,3D3-by-3SPDtensoreldsarederived fromOASISbrainMRIdataset[ 34 ],andthesizeofdeformationeldsfrompair-wise registrationbetweensubjectsinthedatasetareusedasfeaturestoclassifyagegroups, andtheclassicationscoresarereported.Insection 2.4.3 ,3DGMFsdataset,which isderivedfromthehighangularresolutiondiffusionweightedMRIs(HARDI)ofrat 23

PAGE 24

brains,isprepared.Inthistest,theratbrainshavebeensimulatedelectricallytoinduce epilepsy.Then,HARDIscanshavebeencarriedoutpriortothestimulation(whichisthe baselinescan)andinregulartimeperiodsafterthestimulation.Thislongitudinaldataset consistsoftwogroups:post-seizureandpre-seizure.Theproposedmethodisapplied toclassifytwogroupsandtheclassicationscoresarereported.Inthetwoapplication problems,thediffusionmapsareadoptedasaclassicationmethod.Abriefintroduction tothediffusionmapsisprovidedinsection 2.4.1 .Finally,chapter 2 isconcludedin section 2.5 2.2MatchingImageGraphs Inthissection,ageneralintroductiontoimagegraphsisprovidedandthe formulationofanovelregistrationframeworkforimagegraphsisalsodiscussed.In thisdissertationtwodifferenttypesofimagesareconsidered:3-by-3symmetricpositive denite(SPD)tensoreldsandtheGaussianmixtureelds(GMF).Imagegraphswere reportedinmanyearlierworks,specicallyforanisotropicimagede-noising/smoothing [ 2 4 ].Gur etal. haveprovidedthedifferentialgeometricpointofviewofimagegraphs indetail[ 4 ].Theyhaveconsideredimagesassectionsinberbundlesandeachsection isamapping f : B $ E where B isabasemanifoldand E isatotalspacedescribed asproductspaceofthebasemanifoldandtheber, B F .Thatis,for2-or3-Dgray valuedorcolorimages,themappingis f : R n $ R n R m + where n is2or3fortwoor threedimensionalimagesrespectivelyand m is1or3forgrayvaluedorcolorimages respectively[ 2 3 ].Iftheimageisa3-by-3SPDtensorialimage,thetotalspaceis R n P (3) [ 4 ].Thebasemanifoldandthetotalspacearenamedasimagedomainand spatial-featuremanifoldrespectively.Fig. 2-2 showsasectionofaberbundleof3-by-3 SPDtensorimagein2D.Wecallthissectionanimagegraph[ 4 ]. Let'ssupposethatthereisanimagegraphwithparameterizationgivenas X :( u v w ) $ ( r I ( r )) (26) 24

PAGE 25

Figure2-2. Aberandsectionof3-by-3symmetricpositivedenitetensoreldsin2D [ 4 ]. where r =( x = u y = v z = w ) isabasemanifoldand I ( r ) isanimagevalue.Thenthe totalspaceisequippedwithitsRiemannianmetric, G whichistheproductmetricofthe Euclideanmetricin R 3 andthemetricinthespaceofberwhichisdenotedby H : G = ( ( ( ( ( ( ( ) 000 0 00 00 0 000 H + + + + + + + (27) ThentheRiemannianmetrictensoroftheimagegraphisgivenasfollows. K = ( ( ( ( ) < X u X u > G < X u X v > G < X u X w > G < X v X u > G < X v X v > G < X v X w > G < X w X u > G < X w X v > G < X w X w > G + + + + (28) whichistherstfundamentalformof X .Then,thevolumeformofaimagegraphis givenby vol = ) + dudvdw (29) where + isthedeterminantof K 25

PAGE 26

2.2.1CostFunctional Let'ssupposethattherearetwoimagegraphstobematched.Theimagegraphs aremappingsdenedas X 1 :( u 1 v 1 w 1 ) $ ( r 1 I 1 ( r 1 )) X 2 :( u 2 v 2 w 2 ) $ ( r 2 I 2 ( r 2 )), (210) whosedomainsare 1 and 2 respectively.If ) + 1 d 1 and ) + 2 d 2 denotethevolume formsofthetwoimagegraphs,wecanwriteasymmetricandre-parametricinvariant costfunctionalfortheimageregistrationinthefollowingway. Cost ( X 1 X 2 ; )= Dist ( X 1 X 2 & ) ( ) + 1 + & + 2 ( 2 & ) J % ) d 1 (211) Registeringtwoimagesmeanssearchingfor suchthat =argmin % ( Cost ( f 1 f 2 ; )). (212) InEq.( 211 ), X 1 and X 2 sharethecommondomain whichis 1 inthiscaseand J % isthedeterminantoftheJacobianof .AsimilarcostfunctionaltoEq.( 211 )was proposedearlierin[ 26 ]whichisdesignedformatching3-Dfacesurfacemeshes. However,theproposedmatchingalgorithmisnotsymmetricnorisitdesignedfor matchingimagegraphs.ThedetailsofderivationofEq.( 211 )arepresentedin AppendixA. Throughtheregistrationmap,thecoordinatessystemsinthetwoimagegraphsare relatedasfollows: ( ( ( ( ) u 2 v 2 w 2 + + + + = ( ( ( ( ) u 1 + U ( u 1 v 1 w 1 ) v 1 + V ( u 1 v 1 w 1 ) w 1 + W ( u 1 v 1 w 1 ) + + + + (213) 26

PAGE 27

Figure2-3. Matchingimagesasimagegraphscanberealizedviaregistrationmap denedbetweenthetwosurfaces.Thetwosurfacesrepresent2Dimage graphsor2Dslicesof3Dimagegraphs. Here, U V ,and W aredeformationeldswhoseJacobinisgivenby J % = ( ( ( ( ) 1+ U u 1 U v 1 U w 1 V u 1+ V v 1 V w 1 W u 1 W v 1 1+ W w 1 + + + + (214) FinallyusingEq.( 210 )throughEq.( 213 ), Dist ( X 1 X 2 & ) inEq.( 211 )canbechosen arbitrarily;however,inthisdissertationthetermisgivenasfollows. Dist ( X 1 X 2 & ) = (( x 1 ( x 2 ) 2 +( y 1 ( y 2 ) 2 +( z 1 ( z 2 ) 2 )+ dist ( I 1 I 2 & ) = ( U 2 + V 2 + W 2 )+ dist ( I 1 I 2 & ), (215) where dist ( I 1 I 2 & ) isvoxelwiseRiemanniandistancebetween I 1 and I 2 & equipped with H 27

PAGE 28

Fig. 2-3 illustratesthismatchingandthecorrespondencesbetweenimagesare givenbythedeformationvectorelds, U V and W whichminimizeEq.( 211 ). Anotherbenetfromusingtheimagegraphrepresentationisthattherepresentation canbeappliedtoanytypeofimageaslongasweknowthespaceofitsimagevalues. Todemonstratethisuniversality,thefollowingtwosectionsdescribehowtoformulate thisnovelframeworkforimagegraphsof3-by-3SPDtensorimagesandGaussian mixtureeldsin3Drespectively. Inaddition,iftheimagevaluescarrydirectionalinformation,reorientationaccording tothedeformationeldsmustbeformulatedin dist ( I 1 I 2 & ) sothatitisguaranteedthat theorientationalinformationispreservedduringtheoptimizationprocess.Specically, forSPDtensorialimages,theorientationalinformationiscontainedintheeigenvectors. Inthiscase,therearemainlytwotypicalreorientationmethods:thePreservationof PrincipleDirection(PPD)strategy[ 21 23 27 ]andtheFiniteStrain(FS)strategy[ 27 ].In thePPDmethod,theprincipledirectionofthe3-by-3tensorislinearlytransformedby J % andnormalized;then,othereigenvectorsareestimatedbythetransformedeigenvector. Ontheotherhand,intheFSmethod, R istheproductoftwoorthogonalmatricesfrom singularvaluedecomposition(SVD)oftheJacobianof .Thechoicebetweenthe twomethodsdependsonhowanobjectistobereorientedanditspropertiestobe preservedthroughreorientation.Forexample,ifwerepresent(orvisualize)SPDtensor asacylinder,thePPDmethodisdesirable.ThisisthecaseofGMFssincetheprinciple axisofellipsoidofGaussianistobereorientedinthiscase.Ontheotherhand,ifwe represent(orvisualize)SPDtensorasanellipsoidandtheanglesbetweenallthree axesoftheellipsoidarerequiredtobepreserved,theFSmethodismoredesirable.In thisdissertation,theFSandPPDmethodswillbeappliedtotheregistrationof3-by-3 SPDtensorimagegraphsandtheGaussianmixtureeldgraphsmatchingrespectively. 28

PAGE 29

2.2.2CostFunctionalfor3-by-3SPDTensorImageGraphsin3D 3-by-3SPDtensorimagegraphsin3Daremappinggivenasfollows: X :( u v w ) $ ( r I ( r )) (216) where I ( r ) P (3) ,andinthiscase,Eq( 28 )iswrittenasfollows[ 4 ]: K SPD = ( ( ( ( ) + tr (( I 1 I u ) 2 ) tr ( I 1 I u ) tr ( I 1 I v ) tr ( I 1 I u ) tr ( I 1 I w ) tr ( I 1 I v ) tr ( I 1 I u ) + tr (( I 1 I v ) 2 ) tr ( I 1 I v ) tr ( I 1 I w ) tr ( I 1 I w ) tr ( I 1 I u ) tr ( I 1 I w ) tr ( I 1 I v ) + tr (( I 1 I w ) 2 ) + + + + (217) dist ( I 1 I 2 & ) inEq.( 215 )canbeFrobeniusdistance,symmetricKLdistance,or Riemanniandistance[ 28 ].Inthisdissertation, dist ( I 1 I 2 & ) isformulatedasthe Riemanniandistancein P (3) anditisequippedwithreorientationtopreservethe orientationinformationoftensors.Then, dist ( I 1 I 2 & ) inEq.( 215 )isre-writtenas follows: dist ( I 1 ,( I 2 & ) R )= Trace (log( I 1 / 2 R   ( I 2 & ) R I 1 / 2 ) 2 ) (218) where R isthereorientationmatrixappliedto I 2 andthesubscript R represents reorientation.InEq.( 218 ),theFSmethodisadoptedforreorientationschemeto preservetheanglesbetweenthethreeprincipleaxesof 3 3 SPDtensorsunder reorientation.Thedetailsofformulationisdiscussedinsection 2.2.5.1 2.2.3CostFunctionalforGaussianMixtureFieldGraphsin3D Inthissection,thenovelframeworkfortheimagegraphsofGaussianmixtureelds (GMF)isformulated. GMFswerestudiedindetailin[ 29 30 ].GMFisonewaytorepresentthehigh angularresolutiondiffusionweightedMRIs(HARDI).Cheng etal. haveintroduceda novelframeworkforsymmetricregistrationofGMFandmoredetailscanbefoundin referencestherein[ 23 ].However,theirframeworkisnotformulatedwithimagegraphs. ThefollowingsectionprovidesashortintroductiontothespaceoftheGaussianmixture 29

PAGE 30

anddiscusseshowtoformulatetheproposalframeworkintermsoftheimagegraphsof GMFs. TheGaussianmixtureisaweightedlinearsumofthediscretesetoftheGaussian componentswithweights.Ifthereare"n"Gaussianswithadifferentcovariance,the Gaussianmixtureat r isgivenasfollows: I ( r )= n / i =1 i ( r ) G ( r ;,0," i ( r )) (219) where G ( r ;,0," i ( r )) istheGaussiandistributionwithzerosmeanandvariance i r is distancevectorfrom r toneighbors,and i ( r ) 'saretheweights.IntheGaussianmixture elds,theweightsandcovariancesareelementsoftheeldsdenedatvoxels.Then, let GM ( n ) denotethespaceofGaussianmixturewith n Gaussiancomponentsandon thisspace,ifwehavetwoGaussianmixturedenedat r ,wedenetheinnerproduct betweentwoGaussianmixturesby < I 1 ( r ), I 2 ( r ) > GM = R 3 I 1 ( r ) I 2 ( r ) d r (220) Eq.( 220 )isintegratedoverneighborof r .Andtheassociatednorm, || I || iswrittenas || I || = 0 R 3 I 2 d r (221) Thenon GM ( n ) wecandenevoxelwise L 2 distancebetweenGaussianmixturesby dist ( I 1 I 2 )= R 3 ( I 1 ( I 2 ) 2 d r (222) Theadvantageof L 2 distanceofGMFisthatweknowthecloseformexpressionfor theinnerproduct.Ifthecovariancesareassumedtohavecylindricalsymmetryalong theprincipleeigenvectorintheGaussianmixturemodel[ 23 ],thecovariancescanbe decomposedas =( 1 ( 2 ) uu # + 2 I (223) 30

PAGE 31

whoseeigenvaluesare 1 2 = 3 andtheprincipleeigenvectoris u whichisaunit vectorand I representsidentitymatrix.ThengiventwoGaussianmixtures, I 1 and I 2 composedofweightelds i and i ,andcovarianceelds i and # i ,respectively,whose eigenvaluesare { i } 3 i =1 and { i } 3 i =1 ,andprincipleeigenvectorsare u and v ,respectively, where i { 1,2,..., n } ,Eq.( 222 )canbeeasilyevaluatedas[ 23 ] dist ( I 1 I 2 )= # A + # B ( 2 # C (224) where =( 1 ,..., n ) and =( 1 ,..., n ) arevectorsofthemixtureweights,and A B and C arematriceswhoseelementsare 1 2 2 2 2 3 2 2 2 2 4 A i 1 i 2 =((2 / ) 3 det ( i 1 + i 2 )) 1 / 2 B j 1 j 2 =((2 / ) 3 det ( # j 1 + # j 2 )) 1 / 2 C i j =((2 / ) 3 det ( i + # j )) 1 / 2 whichareachievedbytheintegrationoftheproductoftwoGaussians.And det ( "+#) thedeterminantofadditionoftwocovariancesisgivenas det ( "+#)= # ( $ ( u # v ) where # =( 1 + 2 )( 2 + 1 )( 2 + 2 ) and $ =( 1 ( 2 )( 1 ( 2 )( 2 + 2 ) .Weuse Eq.( 224 )forthevoxelwisedistancebetweenGaussianmixtureinEq.( 215 ). TocompletecostfunctionalformatchinggraphsofGMFwithmappingdenedin Eq.( 216 )replacing I ( r ) with I ( r ) GM ( n ) ,weneedtoevaluatetherstfundamental form, K GM ofthegraphswhosedeterminantcontributestothevolumeformas & det ( K GM ) d andthereorientationoftheGaussianmixtures,ormoreprecisely thecovariances,accordingtothedeformationelds.First, K GM mustbewrittenas K GM = ( ( ( ( ) < X u X u > G < X u X v > G < X u X w > G < X v X u > G < X v X v > G < X v X w > G < X w X u > G < X w X v > G < X w X w > G + + + + (225) 31

PAGE 32

ortheelementsof K GM aregivenby < X X > = *0 + R 3 ) I ( r ) ) I ( r ) d r (226) where ) ispartialderivativew.r.t.eachvoxel.Thepartialderivativeof I ( r ) isgivenas follows: ) I r ( r )= n / i =1 [( ) i ( r ) G (;0," i ( r ))+ i ( r ) ) G (;0," i ( r ))] (227) However,whenGMFareconstructedfromHARDI,eachcovarianceisconstanteld overvolume.Therefore,Eq.( 227 )canbesimpliedtothefollowingform: ) I r ( r )= n / i =1 ( ) i ( r )) G (;0," i ( r )) (228) whichgivesclosedformofEq.( 226 )asgivenby, < X X > = *0 + n / i j =1 ) i ) j 1 & (2 / ) 3 det ( i + j ) (229) Second,inthecontextofreorientation,eachcovarianceintheGaussianmixture modelisassumedtohavecylindricalsymmetryandonlytheprincipleeigenvectorsare neededtoperformthereorientation.Therefore,thePPDmethodforreorientation[ 23 ]is aptforourgoal.Inthereorientationstep,wecanrotatetheprincipleeigenvectorsof I 1 or I 2 .Inouralgorithm,theprincipleeigenvectorsof I 1 arereorientedasfollows: u iR = J % u i | J % u i | (230) where J % istheJacobianeldofdeformationand i { 1,2,..., n } istheindexof Gaussiancomponents.Thereorientedversionof I 1 ( r ) iswrittenas I 1 ( r ) R = n / i =1 i ( r ) G ( r ;0," i R ( r )) (231) where i R ( r )=( 1 ( 2 ) u iR u # iR + 2 I .Ineveryinterationofoptimizationprocess,the reoriented I 1 ( r ) R mustbecomparedwith I 2 ( r ) inthecostfunctional. 32

PAGE 33

Finally,aregularizationovertheJacobianeld J % isaddedtoEq.( 211 ).Theterm wasproposedin[ 22 31 ]andisgivenasfollows: Cost reg = L r ( J % ( 1)log( J % ) d 1 (232) whichissymmetric.Thentheoverallcostfunctionalcannowbewrittenas Cost total = Cost ( X 1 X 2 ; )+ Cost reg = L r R 3 Dist ( X 1 X 2 & ) R ( ) + 1 + & + 2 ( 2 & ) J % ) d 1 +(1 ( L r ) R 3 ( J % ( 1)log( J % ) d 1 (233) where J % isthedeterminantoftheJacobian, J % .InEq.( 233 ),thesubscript R representsreorientation.Thereorientationcanbeappliedtoeither X 1 or X 2 L r is positivenumberlessthan1. 2.2.4DeformationFields Inthisdissertation,thedeformationeldsareexpressedasthecompletesetof theorthogonaldiscretesineandcosine(DSC)basis.Thereareadvantageswiththis approachsuchthat: IfweconsideronlythelowfrequencycomponentsoftheDSCbasis,wecaneasily achievesmoothnessofdeformationeldoverthevolume. Toavoidtheoptimizationfrombeingtrappedinalocalminimum,coarsegrain methodsarewidelyadoptedinoptimizationprocess.However,insuchascheme, onehastohandledownsamplingoforiginalimages.TheDSCbasisapproach allowsusavoidthisadditionaldownsamplingstep.Instead,wecanstartfromthe lowestfrequencycomponentsofDSCbasisandkeepaddinghigherfrequency componentsafterreachingacertainlevelinoptimization. 33

PAGE 34

TheDSCbasisarewrittenas U ( u v w )= / n =1, m l =0 A nml $ U nml ( u v w ) (234) V ( u v w )= / m =1, n l =0 B nml $ V nml ( u v w ) W ( u v w )= / l =1, n m =0 C nml $ W nml ( u v w ) and 1 2 2 2 2 3 2 2 2 2 4 $ U nml = N u sin( P nh ( u ))cos( Q mw ( v ))cos( Q ld ( w )) $ V nml = N v cos( Q nh ( u ))sin( P mw ( v ))cos( Q ld ( w )) $ W nml = N w cos( Q nh ( u ))cos( Q mw ( v ))sin( P ld ( w )) (235) whereallfunction $ representbasisfunctionsand N { u v w } arenormalizationconstants. Fromhereon,anysubscriptionforlocalcoordinateswillbedroppedunlessitneedsto bespecied. Tounderstandthefunctions P ab ( x ) and Q ab ( x ) ,weneedtoconsiderthediscretization ofthebasespace(domain).Inpractice,thebasespaceisarectangulargriddenedon Z d .Inthisdissertation, d =3 .InEq.( 235 ), u [0, h ( 1] v [0, w ( 1] and w [0, d ( 1] where { h w d } Z d ,and P ab ( x ) and Q ab ( x ) ,aregivenasfollows: P ab ( x )= / ax b ( 1 Q ab ( x )= / a (2 x +1) 2 b (236) DuetothesinebasisinEq.( 235 ),deformationacrossboundariesisprevented. 2.2.5DerivativesofCostFunctional Nowtheregistrationproblemturnsintooptimizingthetotalcost-functionalw.r.t. coefcients, A nml B nml and C nml inEq.( 235 ).Tosearchforanoptimizationsolution,the nonlinearconjugategradientmethod(NCG)isadoptedforanumericsolver,andthe 34

PAGE 35

rst-orderderivativesofEq.( 233 )w.r.t.thecoefcientsaregivenasfollows: ) ) A nml Cost total = ) ) A nml Cost ( X 1 X 2 ; )+ ) ) A nml Cost reg (237) Ateachvoxel,threetermsneedtobeevaluatedinsideofthevolumeintegralin Eq.( 237 ): ) ) A nml Dist ( X 1 X 2 & ) R (238) ) ) A nml & + 2 ( 2 & ) J % (239) ) ) A nml ( J % ( 1)log( J % ) (240) Here,onlyderivativesw.r.t. A nml arepresentedinthefollowingtwosubsections.Other derivativescouldbedoneinasimilarfashion.Eq.( 239 )andEq.( 240 )arecommon tobothexamplesinthisdissertation,anddetailsofthesearepresentedinAppendixC. Eq.( 238 )couldberewrittenasfollows: ) ) A nml Dist ( X 1 X 2 & ) R =2 U $ U nml + ) ) A nml dist ( I 1 I 2 & ) R (241) Intheaboveequation,thersttermonRHSisverystraightforwardanddependssolely onthebasemanifold.However,thesecondtermmustbeconsideredinmoredetail dependingonthespaceofbers.Thefollowingsubsectionsdiscusshowtoevaluate Eq.( 241 )dependingonthespaceofbers(consideredinthisdissertation)whichare P (3) and GM ( n ) 2.2.5.1Derivativesofcostfunctional:3-by-3SPDtensorimagegraphs Forconvenience,let'srewriteEq.( 218 )below: dist ( I 1 ,( I 2 & ) R )= Trace (log( I 1 / 2 1 R # ( I 2 & ) R I 1 / 2 1 ) 2 ) (242) Thewell-knownnitestrain(FS)method[ 27 ]isadoptedhereforthereorientation schemeinEq.( 242 ).Thatis, R = UV # (243) 35

PAGE 36

where U and V # aretheorthogonalmetricobtainedfromthesingularvaluedecomposition(SPD) of J % = UDV # .Setting = I 1 / 2 R # ( I 2 & ) R I 1 / 2 ,therstorderderivativeofEq.( 242 ) canbewrittenas ) ) A nml dist ( I 1 ,( I 2 & ) R = Trace 5 2log( ) 1 ) ) A nml 6 (244) andhere, ) ) A nml = I 1 / 2 1 5 ) R # ) A nml ( I 2 & ) R + R # ) ( I 2 & ) ) A nml R + R # ( I 2 & ) ) R ) A nml 6 I 1 / 2 1 (245) InEq.( 245 ), ) ( I 2 & ) / ) A nml issimply ) ( I 2 & ) ) A nml = ) I 2 ) u 7 7 7 7 u + U $ U nml (246) However,therst-orderderivativeof R w.r.t. A nml requiresmorealgebra.Thetermcould bewrittenas ) R ) A nml = / i j ) R ) J % ij ) J % ij ) A nml (247) where i and j arematrixindexesof J % ,andfromthedenitionoftheJacobianitisclear that ) R ) A nml = ) R ) J % 11 ) $ U nml ) u + ) R ) J % 12 ) $ U nml ) v + ) R ) J % 13 ) $ U nml ) w (248) Theevaluationof ) R / ) J % ij canbeobtainedbyfollowingtheworkin[ 27 ].Wecanthen rewriteEq.( 244 )andEq.( 245 )factoringoutthebasisfunctionsandtheirderivatives. Finally,tosearchforanoptimizationsolution,weneedtoevaluatetherst-order derivativeof Cost total w.r.t.DSCcoefcients,andthiswouldrequiretheintegrationof Eq.( 238 )orEq( 244 ),Eq.( 239 )andEq.( 240 )overdomain. 2.2.5.2Derivativesofcostfunctional:gaussianmixtureeldgraphs Tohavethisreorientationinvolvedinsearchingforaregistrationmap,Eq.( 224 ) wouldberewrittenas dist ( I 1 R I 2 & )= # A R +( & ) # B ( & ) ( 2 # C R ( & ) (249) 36

PAGE 37

Then,itsrst-orderderivativew.r.t A nml couldbewrittenas ) ) A nml dist ( I 1 R I 2 & )= # ) A R ) A nml +2( & ) # B ) ( & ) ) A nml ( 2 # ) C R ) A nml ( & ) ( 2 # C R ) ( & ) ) A nml (250) where ) A R ) A nml 7 7 7 7 i j = $ u # i R u j R [ ) ( u # i R u j R ) / ) A nml ] (2 / det ( i R + j R )) 3 / 2 (251) ) C R ) A nml 7 7 7 7 i j = $ u # i R v j [ ) ( u # i R v j ) / ) A nml ] (2 / det ( i R + # j )) 3 / 2 (252) ) ( & ) ) A nml = ) ) u 7 7 7 7 u + U $ U nml (253) and ) u i R ) A nml = J 1 5 ) J ) A nml 6 u i R ( 8 u # i R J 1 5 ) J ) A nml 6 u i R 9 u i R (254) AppendixBpresentshowtofactoroutthebasisfunctionsfromEq.( 250 ). 2.2.6NumericalMethod ThenumericalalgorithmforoptimizationisdescribedinAlg. 1 .Aswehave discussedinsection 2.2.4 ,optimizationisstartedwithaDSCbasiswithlowfrequencies; higherfrequencymodesareaddedaftertheoptimizationreachesacertaindegree. Fig. 2-4 showshowthevalueofthecost-functionalconvergeswhenhigherfrequency modesareaddedduringoptimization. 2.3Experiments Inthissection,theproposedmethodiscomparedwithamethodinvolving diffeomorphicDemonsregistration.Five2-Dslicesof3-DbrainMRIscans,along withtheirdiffeomorphicallydeformedimages,werepreparedasmovingimagesand targetimagesrespectively.Inaddition,Gaussianwhitenoisewasaddedtotheimages withvedifferentsetsofvariance: 0.0001 0.001 0.005 0.01 and 0.02 .Thenoiselevels 37

PAGE 38

Figure2-4. Plotofthevalueofthecostfunctionalvs.thenumberofiteration.The numberofbasisfunctionsisincreasedafteracertainnumberofiteration fromthelowestfrequencymode.Theredspotsintheplotrepresentwhere higherfrequencymodesadded. addedtotherespectivelevelsofvariancearereportedinTable 2-1 intheunitofpSNR. ThesesamplesareshowninFig. 2-5A andFig. 2-5B Thegoalwiththisexperimentistocomparetheseestimateddeformationelds withinthisproposedframeworkandtheDemonsmethodwiththeground-truth deformationelds.Thesearegray-scalelevelimages,andthereforeformulatinga correspondingenergyfunctionalisstraightforwardwhencomparedwithworkdoneinthe previoussections. Foraquantitativecomparison,theangularerror(AE)andtheend-pointerror (EPE)arecomputedbetweentheground-truthandtheestimateddeformationvector elds.However,deformationeldsdeformingthebackgroundsofbrainimagesare meaninglessintermsofthiscomparison.Hence,onlydeformationeldsinROIare counted.AsampleispresentedinFig. 2-5C 38

PAGE 39

Algorithm1 Registrationofimagegraphs:NCGmethod Initializedeformationleds: = Id Initializefamilyofcoefcients: A =0 Evaluate # # ( A Cost total usingEq.( 237 )throughEq.( 254 ) Set r 0 + # # ( A Cost total k + 0 Set p 0 +( r 0 while Cost total k +1 / Cost total k < (1 ( 1 ) do Compute # k andset A k +1 = A k + # k p k ; EvaluatedeformationeldsusingEq.( 235 ) anditsJacobian; Interpolate I 2 anditsrstorderderivatives; Reorient I 1 ; Evaluate # # ( A Cost total usingEq.( 237 )throughEq.( 254 ); r k +1 + # # ( A Cost total ; $ k +1 + r k +1 r k +1 r k r k ; p k +1 +( r k +1 + $ k +1 p k ; # k +1 = # k r k p k r k +1 p k +1 ; k + k +1 ; endwhile ResultsarethenreportedinFig. 2-6 ,Fig. 2-7 andTable 2-1 .Someoftheseresults arecomparedvisuallyinFig. 2-6 (indicatedbyarangeofcolor-codeddeformation vectorelds).TheresultstabulatedinTable 2-1 andtheerrorplotsinFig. 2-7 showthat theproposedmethodoutperformstheDemonsalgorithmusingtheAEandEPEasthe errormeasures. 39

PAGE 40

A B C D E F Figure2-5. Anexampleofinputincomparisonbetweentheproposedimagegraph approachandtheDiffeomorphicDemons.A)Originalimage.B)Deformed image.C)Thecolor-codeddeformationvectoreld.D-E)TheimageofA andB,respectively,withtheGaussiannoisewithvariance 0.02 .F)The colormapofC. 40

PAGE 41

A B C D E F Figure2-6. Resultingdeformationvectoreldswithintheimagegraphapproachandthe DiffeomorphicDemons.A-C)Thecolor-codeddeformationvectoreldswith theimagegraphapproach.D-F)theDiffeomorphicDemons.Theground truthofdeformationvectoreldsispresentedinFig. 2-5C .AandDrepresent estimateddeformationeldsinmatchingimageswithoutnoise.BandE showresultingdeformationeldsfrommatchingnoisedimageswith var =0.0001 .NotethatinimagesCandF,var =0.02 41

PAGE 42

Table2-1. Theangleerror(AE)andtheEuclideandistance(ED)betweenthegroundtruthofthedeformationvectorelds andtheestimatedonesbyourimagegraphapproachandthediffeomorphicDemonsarecomputed. Comparisonsareperformedwithfourdifferentkindsofimages:with/withoutnoisesanddenoisedimagesbythe Laplace-BeltramidiffusionlterandtheGaussianlter. OriginalimagesGaussiannoise GaussiannoiseGaussiannoiseGaussiannoiseGaussiannoise Var =0.0001 Var =0.001 Var =0.005 Var =0.01 Var =0.02 pSNR( dB ) 41.63 41.84 31.70 31.83 24.74 25.02 21.80 21.93 18.98 19.11 Errormethod AEEDAE(avg)ED(avg)AE(avg)ED(avg)AE(avg)ED(avg)AE(avg)ED(avg)AE(avg)ED(avg) ImageGraph0.00640.02760.05880.24450.06500.27260.10310.42170.15990.66790.28401.1880 Diffeomorphic Demons 0.04440.19570.11290.45160.20400.84270.40501.42150.39881.63630.47651.9237 42

PAGE 43

A B Figure2-7. Theplotsoferrormeasures.A)theangleerror(AE)andB)theend-point error(EPE).Thex-axisisthevariancesoftherandomGaussiannoises.The redandbluedotsrepresenterrorsfromthediffeomorphicDemonsandfrom ourimagegraphapproach,respectively.Itisclearthatallerrorsfromour methodaredistributedundererrorsfromthediffeomorphicDemonsmethod. 2.4Application Inthissection,theframeworkinappliedtoclassicationproblemsinordertoverify theproposedframework.First,datasetswithtwoorthreegroupstobeclassied areprepared,andasapre-processingstep,alldataaresubjecttothesimilarity transformationwithintheregionofinterest(ROI).Second,withthesepre-processed data,pairwiseregistrationiscarriedoutbetween X i and X j ,whicharetwodifferent imagegraphsinthedataset; L 2 normsofthedeformationeldsarethenevaluated.The normisgivenasfollows: d ( X i X j )= : | U | 2 + | V | 2 + | W | 2 ; d (255) AdatagraphwithGaussianweights e d (X i ,X j ) / ) isthenbuilt,andthecorresponding Markovmatrixisusedtocomputediffusionmaps[ 32 33 ].Finally,theEuclidean distanceisadoptedinthespaceofdiffusionmapswithreduceddimensionsasfeatures, anditsclassierisnearestneighborclassier. 43

PAGE 44

ThefollowingsectionsserveasabriefintroductiontoDiffusionmapswhile mountingdiscussionofclassicationofa3-by-3SPDtensorimagedatasetandthe Gaussianmixtureeld(GMF)datasetin3-D.TheSPDtensoreldsarederivedfrom OASISdatasetofMRbrainimages[ 34 ],andGMFisconstructedfromadataset consistingofhighangularresolutiondiffusionweightedMRIs(HARDI)ofratbrains. 2.4.1DiffusionMaps Diffusionmapsandtheirapplicationsareintroducedin[ 32 33 ].Itisbenecialto recallthedenitionofDiffusionmaps.Oncewehavebuiltagraphwithaweighting functionandaMarkovmatrix, M ,thefamilyofdiffusionmaps { t } t $ N andthediffusion distance D t ( X i X j ) aredenedasfollows: t ( X i )= ( ( ( ( ( ( ( ) % t 1 & 1 ( X i ) % t 2 & 2 ( X i ) % t s & s ( X i ) + + + + + + + (256) and D t ( X i X j )= < s / l =1 % 2 t l ( & l ( X i ) ( & l ( X j )) 2 = 1 / 2 (257) where { % l } l % 0 and { & l } l % 0 areeigenvaluesandeigenvectorsof M respectivelysuchthat 1=% 0 > | % 1 | *| % 2 | ... and M & l = % l & l [ 32 33 ],and s = max { l N suchthat | % l | t > 0 | % 1 | t } (258) Throughthediffusionmaps,alldataintheset { X i } areembeddedintotheEuclidean spacewithreduceddimension s ; D t reectstheconnectivityamongelementsinthe datasetin R s ,suchthatpointsin R s whichrepresenteach X i ,arecloseriftheyare highlyconnected.Therefore,wecouldusethenearestneighborasaclassier,and D t asthedistancemeasuredin R s 44

PAGE 45

2.4.2Application:ClassicationofOASISDataset 2.4.2.1Datapreparation TheOASISdataset[ 34 ]consistsof3-DMRbrainscansfrom416humansubjects withagesrangingfrom18to92.includingsubjectsdiagnosedwithAlzheimer'sdisease (AD).Inthisexperiment,subjectsaredividedintothreegroups:thoseaged40or below(Young);thoseagedbetween41and60(Middle);andthoseaged61andabove (Old).ThegoalofthisexperimentistoclassifytheMRimagedatasetwithrespect tothedifferentagegroupswithintheproposedframeworkformulatedforthe3-by-3 SPDtensorimagegraphsin3-D.However,sincetheOASISdatacontains16-bit intensityvalues,thequestionthenishowtoconstructthetensor-valuedimagesandthe correspondingimagegraphs.Inthisapplication,thetensor-valuedimagegraphsare constructedasfollows: 1.ThegraphsoftheOASISMRimagesareconstructedby,mapping f : R 3 $ R 3 R + 2.Therstfundamentalformofeachgraphisderived,whichisa3-by-3SPD tensor-valuedimageasafunctionoftheimagedomain. 3.Theimagegraphsofthetensor-valuedimagesareconstructedbymapping: X : R 3 $ R 3 P (3) Inthisapplication,theventricleischosenasaregionofinterest(ROI)sinceitisoneof thestructuresthatcapturesthemostsignicantdifferenceamongagegroups.Fig. 2-8 andFig. 2-9 show2-Dslicesfromthesuperiorviewpoint;mid-sagittalsectionsofOASIS MRimagesandthetensor-valuedimagesarederivedfromthese.Thetensor-valued imagesarevisualizedusingasphericalharmonicbasis. 2.4.2.2Experimentalresults TheclassicationresultsoftheOASISdatasetarereportedinthissection.The resultsconsistofclassicationscoresbetweensubjects:Youngvs.Middle;Youngvs. Old;andMiddlevs.Old,respectively.SubjectswithAlzheimer'sdisease(AD)vs.the controlgrouphavealsobeenclassied.Inordertovalidatethemethodofclassication, 45

PAGE 46

A B C D E F Figure2-8. Slicesof3-DMRIimages.A-C)2-Dslicesfromsuperiorviewpointof ventriclesfrom18,43,and81yearsoldrespectively.D-F)Mid-sagittal sectionsofventriclesinthesameorder. A B C D E F Figure2-9. Slicesof3-Dtensor-valuedimagesofventricleswhicharetherst fundamentalformsofthegraphsofFig. 2-8 .A-C)2-Dslicesfromsuperior viewpointof18,43,and81yearsoldrespectively.D-F)Mid-sagittalsections ofventriclesinthesameorder. leave-one-outvalidationandfour-foldcross-validationmethodsareused.Fig. 2-10 showsclassicationresultswithin2-Dplotof & 1 vs. & 2 (whicharethersttwoaxis ofdiffusionmapsdenedinEq.( 256 )).Inthefour-foldcross-validationtest,each subgrouphasaboutquarteroftotalmembersofitsparentgroup,andthesubgroups arerandomlyselected50times.Theresultsfromthetwovalidationtestsarereportedin Table 2-2 .ToclassifyAD,70subjectsfromOldareselected.Thesesubjectsconsistof twogroups:35aresubjectsdiagnosedasAD,and35serveascontrol.Inthesetests, 46

PAGE 47

theparametersofdiffusionmapsdenedinEq.( 258 )arechosenas 0 =0.07 and t =1 A B C Figure2-10. 2-Dplotsofdiffusionmaps.A)Youngvs.Old,B)youngvs.middle,andC) middlevs.old.Ineachplot,x-axisandy-axisare & 1 and & 2 respectively. Theresultingclassicationscoresarecomparedwithscoresfromothermethods. First,Table 2-3 showstheclassicationscoresbytheDemonsalgorithm,whichis formulatedbasedonmatchingintensities.Thatis,imagesarerepresentedbyimage 47

PAGE 48

Table2-2. Scoresofleave-one-outandfour-foldcross-validationtests.Foursubgroups havebeenrandomlyselected50times.Themetricusedfortheambient space R 3 P (3) istheproductmetricoftheEuclideanmetricon R 3 andthe afne-invariantmetricon P (3) Input:3-by-3SPDimagegraphsOldvs.YoungOldvs.MiddleMiddlevs.YoungADvs.Control Maximum 100% 100% 100% 100% Minimum 97.72% 95.38% 94.54% 75.0% Average 99.54% 98.98% 98.54% 94.87% Standarddeviation 0.79% 1.09% 1.62% 5.55% Leave-one-out 99.43% 98.46% 99.09% 95.32% functionsinthisalgorithm.Theclassicationmethodandvalidationmethodsusedin Table 2-3 areidenticaltothoseinTable 2-2 .Itisobviousthatthenovelframeworkshows superbresultsovertheconventionalregistrationmethod. Table2-3. Scoresofleave-one-outandfour-foldcross-validationtests.bytheDemons algorithm.ThecostfunctionalinthiscaseisSSDofintensityvaluesofthe OASISMRIdataset.Therandomlyselectedfoursubgroupsaresimilartothe subgroupsinTable 2-2 Input:imagefunctionOldvs.YoungOldvs.MiddleMiddlevs.Young Maximum 100% 100% 100% Minimum 96.59% 93.85% 85.45% Average 98.81% 97.35% 94.94% Standarddeviation 0.88% 2.065% 2.96% Leave-one-out 98.3% 97.69% 97.27% DuringthetestrepresentedinTable 2-2 ,themeticusedistheproductmetricof theEuclideanmetricon R 3 andtheafneinvariantmetricon P (3) .Thisisnottheonly choiceonecanmake.Alternatively,wecanconsiderthespaceofberas R 6 ,sothat themappingofthegraphis X : R 3 $ R 3 R 6 andthiswouldreplaceEq.( 218 )with theFrobeniusnorm.However,withthischoice,wewilllosethepositivedenitenessof tensors,whichiscrucialfeatureoftheRiemannianmetrictensor.Anotherchoicewe couldmakeinvolvesthevolumeform.Tagare etal .introducedavolumeformwhich wouldmakeacost-functionalforregistrationsymmetric[ 22 ].However,thisvolumeform isforimagefunctions,notimagegraphs.Classicationscoreswiththischoiceofmetric (withandwithoutreorientation)andthevolumeformproposedinthisdissertationare reportedintherstandsecondcolumnsofTable 2-4 48

PAGE 49

Ontheotherhand,thevolumeformproposedin[ 22 ]andtheRiemanniandistance in P (3) havebeenusedforthetestinthelastcolumnofTable 2-4 .Thegroups comparedinTable 2-4 areYoungvs.Old.IfwecompareresultsinTable 2-2 with resultsinthersttwocolumnsof 2-4 ,wecannoticethatthedistancefunctionofthe tensorialimagegraphsequippedwiththeafneinvariantmetricon P (3) yieldssuperior classicationresultstoonewithFrobeniusmetricon R 6 .Thiscomparisonclearly demonstratesthatpreservingthespaceofberintheformulationofoptimizationis crucialinimageregistrationandclassicationproblems.Also,comparisonbetween Table 2-2 andthethirdcolumnof 2-4 showsthatthesymmetrizedvolumeformwith imagegraphrepresentationworksbetterthanonewithimagefunctionrepresentationin classicationproblems. ItbearsmentioningthatthescoresinthethirdcolumnofTable 2-4 areslightly betterthanthescoresachievedwiththeconventionalintensity-value-matching algorithmreportedinTable 2-3 .Perhapsthetensor-valuedimagesderivedinthis testcarryneighborhoodinformation,notjustinformationateachvoxel,andperhapsthis geometricalinformationhasbolsteredaccuracyinmatchingimagesandinclassication problems. Table2-4. Scoresofleave-one-outandfour-foldcross-validationtestswithdifferent metricsontheambientspace(withandwithoutreorientation).Four subgroupsarerandomlyselected50timesfromYoungandOldgroups.Inthe rstandsecondcolumns,theFrobeniusnormisadoptedasthedistance metricofSPDtensors.Thevolumeformistheoneoftensor-valuedimage graphsinthiscase.Inthelastcolumn,theRiemanniandistanceofSPD tensorsisused,howeverthevolumeformistheoneintroducedin[ 22 ]. Methods Frobenius Frobenius Riemanniandistance withreorientationwithoutreorientationwithreorientation andstandardvolumeform Maximum 100% 100% 100% Minimum 88% 88.63% 96.59% Average 94.96% 96.52% 99.03% Standarddeviation 2.57% 2.15% 0.99% Leave-one-outvalidation 97.03% 98.01% 99.15% 49

PAGE 50

Table2-5. Comparisonofclassicationscoreswith4-foldvalidationbetween classicationmethods Oldvs.Oldvs.Middlevs.ADvs. YoungMiddleYoungControl TensorialImageGraphs 99.54%98.98%98.54%94.87% CAVIAR[ 36 ] 99.14%98.36%97.76%88.0% Adaboost[ 36 ] 98.75%96.80%96.0%84.38% Submanifoldprojection[ 35 ] 96.43%90.23%84.32%88.57% NearestNeighborinPCA[ 35 ] 92.43%87.74%78.42%84.29% InTable 2-5 ,theresultsfromourproposedmethodarecomparedwithpreviously publishedclassicationresultsontheOASISdataset.In[ 35 ],Xie etal .haveregistered eachimagetoacommonatlas;theycomputedthedeformationtenoreldsofeach asmainfeaturesforeachimage.Thesubmanifoldofeachagegroupwasthen constructedfromtrainingsamples;thegeodesicdistancesbetweensubjectsand thesubmanifoldswereusedasthemaindiscriminativefeatureforclassication. Alternatively,histogramsofdeformationvectoreldshavebeenusedasfeatures [ 36 ],andtheCAVIARmethodproposedin[ 36 ]takesaadaboost-likeapproachto integratetheresultsfromacollectionofweakclassiersintoastrongclassication result.Twothingsbearmentioningthatrst,theimagematchingmethodsusedinthese referencesaretheusual L 2 -norm-basedmethods,whichonlycompareimagevalues atcorrespondingvoxels;second,theproposedmethodcomparesfavorablywiththese methodsintermsofclassicationrates.Inparticular,forthemorechallengingproblem ofclassifyingbrainimagesofAlzheimer'sdiseasepatients,theproposedmethod demonstratesasmallbutrealimprovementoverthesetwomethods. 2.4.3Application:ClassicationofGMFDataset Thissectionconcernsthematchingandclassicationofimagegraphsembedded inhigherdimensionsthanthepreviousexample.Here,imagegraphsarethegraphs ofGaussianmixtureeld(GMF).TheformulationofmatchingGMFimagegraphswas discussedpreviouslyinsection3. 50

PAGE 51

Table2-6. Scoresfromleave-one-outvalidationtestwithtwodifferentmatching methods.Classicationscoresofpre-andpost-seizurearereported separatelyandtotalscoreisreported. ImagegraphofGMFImagefunctionofGMF Pre-seizurePost-seizurePre-seizurePost-seizure Scroe 100% 75% 91.18% 50% Totalscore 93.47% 80.43% 2.4.3.1Datapreparation Vivohighangularresolutiondiffusion(HARDI)scansofratbrainswereprepared; GMFsweregeneratedfromtheHARDI(with46Gaussiancomponentsateachvoxel bythereconstructionmethodin[ 30 ]).Themappingofthegraphsisthengivenbythe following: X : R 3 $ R 3 GM (46) (259) Thedatasetconsistsofthebrainsofninedifferentrats;thesebrainshavebeen electricallystimulatedtoinduceepilepsy.HARDIdatahavebeencollectedtemporallyto observewhenseizuresoccurredduringthetimesequenceofscans.Foreachrat,there arevetoeighttemporaldata.Werefertodatabeforeseizurebutpoststimulationas pre-seizure;datacollectedafterseizureisreferredtoaspost-seizure.Theaiminthis sectionistoclassifypre-seizureandpost-seizuredata.Inthistest,thereare34setsof pre-seizuredataand12post-seizureandpair-wisematchinghasbeencarriedoutto constructaMarkovmatrixasinstructedinsection3. 2.4.3.2Results Inthisexperiment, t and 0 aresetequalto1.0and 0.31 inEq.( 258 )respectively. Twovalidationmethodsareadoptedtoverifyclassicationscores:theseinclude four-foldcrossandleave-one-outvalidation.Thescoresfromtheleave-one-out validationtestisreportedinTable 2-6 ;theyarethencomparedwiththescorefromam intensity-basedimageregistrationmethod.Again,itshouldbenotedthattheproposed frameworkdoesoutperformtheusualintensity-matching-basedalgorithm.Fig. 2-11 51

PAGE 52

A B C Figure2-11. Datapointsprojectedon2-Dplanefrom R s =10 ofdiffusionmaps.A) & 1 vs. & 2 B) & 1 vs. & 3 C) & 2 vs. & 3 .Redstarsandbluecrossesrepresent pre-seizureandpost-seizurerespectively. Table2-7. Averageclassicationscoresoffour-foldvalidationtestfrom50different combinationsoftestandtrainingdata. Pre-seizurePost-seizurePre-andPost-seizure Ave.score 96.25%75.6% 90.61% showsdatapointsofthediffusionmapsprojectedon2-Dplanesoftherstthreeaxes, { & l } 3 l =1 inthecaseoftheleave-one-outvalidationtest. Astothefour-foldcrossvalidationtest,8setsofpre-seizuredataand3setsof post-seizuredataarerandomlyselectedoutof34and12setsrespectivelyanddiffusion mapshavebeenbuiltthroughpairwisematching.Thisrandomselectionhasbeen repeated50times;averagescoresarereportedinTable 2-7 2.5Conclusion Chapter ?? hasproposedanovelimagematchingframeworkwiththreeproperties: symmetry,re-parameterizationinvariance,anduniversality.Thesearerealizedby consideringimagesasimagegraphswhichare2-D/3-Dsurfacesembeddedinan ambientspace,whichisthecartesianproductofimagedomainandthespaceofimage values.Theoptimizationprocessintheframeworkistosearchfordeformationelds whichminimizethenovelcostfunctional;ithasbeenassumedthatthedeformation eldsareexpressibleinanorthonormalbasis.Theoptimizationproblemthenfocuses 52

PAGE 53

onasearchforoptimalcoefcientsofthebasisset.Thesmoothnessofthedeformation eldscouldbeachievedwithoutbreakingsymmetryofcostfunctionalonlybydiscarding highfrequencycomponentsinthecompleteset.Inaddition,intheoptimization algorithm,coarse-to-negrainstrategyhasbeenadoptedthroughtheuseofbasis functions.Inthisway,theresultingcostfunctionalcanbesatisfactorilyoptimizedusinga gradientdescent-basedmethod.Toverifytheproposedframework,ithasbeenapplied totheclassicationoftwodifferentkindsofimagedatabases:3-by-3symmetricpositive denitetensorimagesandGaussianmixtureeldsin3-D.Thetensorialimageshave beenconstructedfromOASISMRimagedataset;classicationbetweenagegroupshas beenconducted.Theseclassicationresultshavebeencomparedwithfourcompeting methods;thiscomparisonshowsthattheproposedmethodyieldsresultsthatexceed currentstate-of-the-artresults.Theproposedmethodspecicallyoutperformsother methodsinclassifyingAlzheimer'sdisease(AD).However,amoreconvincingwayto diagnoseADismaybetomeasuretherateofchangeofventriclesizeinindividuals overconsecutiveperiodsoftime.Longitudinalanalysesforthispurposehavenotbeen conductedinthisdissertation. InTable 2-4 ,tensorialimageshavebeenembeddedindifferentambientspaces, andthesameclassicationtestshavebeenconducted.ComparisonsbetweenTable 2-2 andTable 2-4 showthatbesidestheimagegraphrepresentation,thefollowing conditionsmustbesatisedtomakeanimageregistrationframeworkyielding satisfactoryresultsinregistrationandclassicationproblems:1.thespaceofimage valuesmustbepreservedintheprocessofoptimization;2.thereorientationstepis crucialifimagevaluescarryorientationalinformation. Classicationscoresfromtheproposedmethodwerealsocomparedwithscores fromtheDemonsalgorithm,whichisanintensity-basedmethod.Itshouldbenoted thatscoresinTable 2-2 andinthethirdcolumnofTable 2-4 arehigherthaninTable 2-3 Onepossibleexplanationforthisiswasfollows:withinthisclassicationtest,themain 53

PAGE 54

differencesamongagegroupscorrespondwithdifferingsizesofventriclesamong subjects.Ifwecomparetwoagegroupsbymatchingintensityvalues,acomparisonof theareassurroundingtheventriclescouldbecost-prohibitive;atanyrateregistering changeinsurroundingareascouldbemoreimportantthanmeasuringchangesinthe sizeoftheventricle. Ontheotherhand,theSPDtensordatasetcarriesedgeinformationwhichismostly suppliedbyexaminingthechangeofventricleshape;hence,classicationscoresby registeringtheSPDtensordatasetderivedfromimagefunctionsaresuperiortothe intensity-basedregistrationmethod. Forthesecondclassicationproblem,Gaussianmixtureeldswerebuiltfrom HARDIdataacquiredfromratbrains;thedatasetwasdividedintotwogroups: pre-seizureandpost-seizure.Classicationresultsfromwithintheproposedframework werecomparedwiththe L 2 -norm-basedmethod;itisveryclearthattheproposed matchingmethodoutperformstheusual L 2 -norm-basedmethod. Althoughthenovelframeworkisformulatedforimageregistrationor,more precisely,forimagegraphsthisframeworkdoesnotneedtobelimitedtoimage registration.Ifweremindourselvesthatanimagegraphisa2-Dor3-Dsurface embeddedinahigherdimension,andifwenotethattheframeworkmatchesthis surfaceintrinsically,theuniversalitypropertycouldthenbeextendedtoageneral shape-matchingproblem.Withthisextension,thisimagematchingproblembecomesa subsetofashape-matchingproblem;thisisanothernoveltywhichattendstheproposed frameworkinchapter 2 Therefore,thoseconductionresearchinfuturemaybedirectedtoapplythis proposedframeworktomatchingandclassifyingclosedsurfaces(suchas2-Dclosed surfacesembeddedin3-D). 54

PAGE 55

CHAPTER3 ANOVELFRAMEWORKFORCOMPUTINGDIFFEOMORPHICPATHS Inchapter 2 ,animagegraphrepresentationhasbeenutilizedinordertobuild anovelregistrationframeworkwhichissymmetric,parametrization-invariant,and universal.Chapter 3 concernsanovelframeworkforcomputingdiffeomorphicpaths betweentwogivendiffeomorphisms.Morespecically,ifweassume denotethe 2-D/3-D(image)domain,andifwelet Di! ( ) denotethespaceofdiffeomorphisms,we thenassume : $ .Giventwodiffeomorphisms 0 1 asmoothpathin Di! ( ) connecting 0 1 wouldofferasmoothmap :[0,1] $ suchthatforeach t [0,1] ( t ,.) Di! ( ) isadiffeomorphismand (0,.)= 0 (1,.)= 1 .Computing ( t ) haslongbeenanimportanttopicinMedical ImageAnalysisforlongitudinalstudies,orforananalysisofshapechangesover time.However,evaluatingdiffeomorphicpathsbetweengiventwodiffeomorphismsisa challengingproblem:itrequiresknowledgeofthemeticinthespaceofdiffeomorphisms (whichisunknown). Oneofthemostintensivelystudiedframeworksforcomputingdiffeomorphicpaths istheLargeDeformationDifferomorphicMetricMappingframework(LDDMM)[ 8 9 ]. Withinthisframework,thediffeomorphicpath isobtainedbysolvingthefollowing equation: d dt ( t ,.)= v ( t ( t ,.)) (31) inwhich v ( t ( t ,.) isavelocityelddenedonthemanifoldof andistheminimizerof thefollowingcostfunctional: E = 1 0 || v ( t ( t ,.) || 2 V dt + || I 1 ( I 2 & ( t =1,.) || 2 L 2 (32) InEq.( 32 ),thenormof v isdenedwithanappropriateSobolevnormon v [ 9 ]. 55

PAGE 56

Beg etal. solvedEq.( 31 )byintegratingEq.( 31 )[ 9 ].Thatis, ( t =1,.)= ( t = 0,.)+ % 1 0 v ( t ( t ,.)) dt Also,Beg etal. haveshownthattheestimatedpathsofdiffeomorphismsby LDDMMstayonthespaceofdiffeomorphisms[ 9 ].ThisLDDMMframeworkhasbeen furtherdevelopedinrecentyearswithinthecontextofevolutionequationsongroupsof diffeomorphismsandmomentumconservationin[ 10 11 15 ].Younes etal. introduced geodesicequationsonthegroupofdiffeomorphismsbasedontheLDDMM[ 11 ]; Sommer etal. introducedamulti-scalekernelbundleforLDDMMinordertoincrease theaccuracyofLDDMM[ 12 ].TheycalledittheLargeDeformationDiffeomorphicKernel BundleMapping(LDDKBM),andtheyintroducedevolutionequationsforthisLDDKBM (shownin[ 13 ]). Alloftheseworksprimarilyexistinordertounderstandthespaceofdiffeomorphisms (oftheimagedomain)denotedby Di! ( ) sothatonemightestimateadiffeomorphic path.Onecanthenutilizetheresultingpathforthepurposeofamedicalimage analysis(suchasanalysesofMRIbraindataorcardiacdata).Speakinggenerally, theseframeworkscouldbeappliedvariously:tomotionestimation;tostatisticsonthe spacesofdiffeomorphismsandimages;andtolongitudinaldataanalysis. Theoretically,frameworksforcomputingadiffeomorphicpathbetweenimages requiretwodifferenttypesofoptimization:theregistrationofimagesandtheinterpolation ofdiffeomorphisms.LDDMM-likeframeworksaredesignedtoachievethesetwo differentoptimizatoinsinasingleformulation.However,therearesomeissuesweneed toconsider. TheseLDDMM-likeapproachesdependonhowweformulatetheSobolev metriconthevelocityeldsshownas V inEq.( 32 )anditbearsmentionthat theseapproachesmaycauseundesirablegeometricdistortioninthespaceof diffeomorphisms. Thereisnoguaranteethattheseformulationscanyieldbestmatchingresults betweeninputimagesregardlessofthespaceofanimagevaluesinceeach 56

PAGE 57

differentspaceofimagevaluesmayrequireauniqueformulationofregistration framework. Thesmoothnessconstraintoveradiffeomorphicpathmayhamperayieldof satisfactorycorrespondencesbetweeninputimages. Tuningthesetofparametersforaconvergentresultisheavilydata-dependent,and iscumbersomeandtedious. Therefore,analternativewaytoovercometheseissuesmightbetoseparate registrationandinterpolationinordertomaximizetheefcienciesofthetwooptimizations. Thatis,wemayrsttrytoachievebestdiffeomorphiccorrespondencesusingaspecic registrationmethod;adiffeomorphicpathwouldbecomputedinaseparateframework, takingthediffeomorphismsfromtheregistrationstageasinput.Moreover,ifweneedto dostatisticalanalysesinthespaceofdiffeomophisms,itmightbedesirabletoworkwith oneofthediffeomorphicpathscomputedintheseparateframework.However,while thereareplentyofregistrationframeworksforvariousimagemodalities,itbearsnoting thatframeworksfortheinterpolationofdiffeomorphismsarefewandfarbetween.Thisis mainlyduetotheinnitedimensionalityofthespaceofdiffeomorphismsandtothelack ofknowledgeofthemetricinthisspace. Inthisstudy,allthesedifcultiesareovercomebyprojectingdiffeomorphsmsontoa knownspaceinwhichweknowtheclosedformofthegeodesicsandwhosemetric. Supposethatwehavetwodiffeomorphismstobeinterpolatedand (0,.)= id and (1,.)= denotethetwodiffeomorphismsrespectively.Thenovelframeworkwould thenconsistoftwosteps.Intherststep,thesquarerootsofthedeterminantsofthe Jacobiansassociatedwiththegivendiffeomorphismsareprojectedontothespaceof densities,andthegeodesicconnectingtheprojectedpointsisevaluatedinthesame space.Inthesecondstep,thegeodesicisliftedbacktothespaceofdiffeomorphisms. Morespecically,wewouldlet M denoteacompact n -dimensionalRiemannian manifold,and Di! ( M ) and Di! ( M ) denote,respectively,theinnite-dimensionalgroup ofdiffeomorphismsof M anditsinnite-dimensionalsubgroupofvolume-preserving 57

PAGE 58

Figure3-1. Thetwostepsoftheproposedframeworkfortheinterpolationof diffeomorphisms. f (0) f ( t f ) and f ( t ) representprojectedpointsoftwo givendiffeomorphismsandthegreatcircleconnectingthetwoprojected pointsonthesphererespectively. diffeomorphisms.Wewouldthenconsiderthequotientspace Di! ( M ) / Di! ( M ) which hasbeenstudiedin[ 37 ]andwhichmaybeidentiedwiththespaceofdensities functions Dens ( M ) on M .Thelatterspacecanthenbecanonicallyembeddedintoa sphereintheHilbertspace(See[ 37 ]fordetails).Thiswouldimplythatthegeodesicof thetwoprojectedpointsisagreatcirclewhoseformiswell-known. Inthesecondstep,thegeodesicin Dens ( M ) isliftedbackto Di! ( M ) bysolvinga quadraticprogrammingproblemwithbilinear/trilinearconstraintsfor2-D/3-Ddomains. Thisoptimizationproblemcouldbesolvedusingthenumericalalgorithmproposed in[ 40 ],whichisbasedonaugmentedLagrangianmethod[ 38 43 ].Itbearsnotingthatin LDDMM,thediffeomorphicpathiscomputeddirectlyin Di! ( M ) withoutreferencingthe quotientspace Dens ( M )= Di! ( M ) / Di! ( M ) ;theproposedframeworkforinterpolating diffeomorphismsiscompletelyseparatedfromanyregistrationframeworks.Fig. 3-1 illustratesthetwostepsinthisnovelapproach. 58

PAGE 59

Theapplicationofthisnovelframeworkmayvary.Itcouldprovideanovel frameworktosynthesizeadiffeomorphicpath.Itmightalsobeappliedtomotion estimation,tovideoanalysis,andtoastatisticalanalysisofshapesandimagesfora longitudinaldataanalysisinmedicalimaging.Inthisdissertation,thismethodisapplied totwotopics:toanestimationofframesinavideo,andtoancardiacmotionanalysisfor thepurposeofclassication. Intheestimationofvideoframes,atestvideofootageistemporallydown-sampled. Thatis,videoframesaredroppedoutremainingkeyframes,andmissingframesare computedwiththeproposedframework.Thisexperimentshowshowtheproposed frameworkcouldsimulatethegroundtruthofadiffeomorphicpath.Inaddition,inthis experiment,theproposedframeworkisalsocomparedwithLDDMM. Intheapplicationoftheframeworktocardiacmotionanalysis,thediffeomorphic pathsofcardiacmotionsarecomputedandcharacterizedinasubject-specic manner.Thissubject-speciccharacterizationcanthenbeusedtoclassifycontrols vs.pathologies.Themainideaistousesmoothpathsinthespaceofdiffeomorphisms asthecharacterizationofthecontinuousvariationinshapesandvolumeofthepatient's leftventricularmyocardium(LVM)acrosscardiaccycles.Thepurposeofthisexperiment istoshowthatvolume-preservingdiffeomorphismsmaybenuisanceparametersin cardiacmotionanalysis.ThemotionofLVMundergoesconstantsmallandrandom perturbationduetopressureextractedfromitssurroundings;thismaycausenoisein computingadiffeomorphicpath.However,themotionofLVMalwaysproducesthesame volumeofejectaineverycardiaccycle;andthismeansthattheserandommovements couldbeviewedasnuisanceparameters.Thisparticularinsightleadsnaturallytothe ideaoftreatingthevolume-preservingdiffeomorphismsasnuisanceparameters;the desireddiffeomorphicpathsshouldbesmoothpathsinanappropriatespacethathas "takenout"theeffectsofthesenuisanceparameters.Thisveriestherststep:the proposedframeworkinwhichdiffeomorphismsareprojectedonthequotientspace. 59

PAGE 60

Therestofthischapterisorganizedasfollows:inthefollowingsection,the theoreticaldetailsofthenovelframeworkareprovided;thisisfollowedbythediscussion ofanumericalschemeandsolverfortheoptimizationproblemoftheproposed frameworkinsection 3.2 .Section 3.3 presentstheexperimentsofthetwotopics. Thedetailsofdatadescriptionandresultsofthevideosequenceexperimentsandthe cardiacmotionanalysisarepresentedinsection 3.3.1 and 3.3.2 ,respectively.Lastly, Thensection 3.4 providesadiscussion. 3.1ModelFormulation Thissectionintroducesanovelframeworktoestimatethepathofdiffeomorphisms. Theestimationofgeodesicsonthespaceofdiffeomorphismsrequiresitsmetric; however,wehavenoknowledgeofitsthegroundtruth.Thereforeratherthanworking onthemanifoldofdiffeomorphismswhosemetricisunknowndiffeomorphismsare mappedtoaknownspace,whichisthespaceofdensitiesinthenovelframework.This frameworkconsistsoftwosteps: Mappingagivenpairofdiffeomorphismstothespaceofdensities, Dens ( M ) ,and estimatingitsgeodesiconthespace. Liftingthegeodesicevaluatedin Dens ( M ) backto Diff ( M ) 3.1.1StepI.theSpaceofDensitiesandDiffeomorphisms Let'ssupposewehaveapairofdiffeomorphisms.Withoutalossofgenerality,we cansetoneofthepairtotheidentitymap id andtheothertoadiffeomorphism .The aimistoestimatediffeomorphisms, ( t ) ,between id and ,suchthat (0)= id and ( t = t f )= where 0 t t f Ifweformulatethediffeomorphismswithvectoreldsdenedinaninnite dimensionalRiemannianmanifoldM,andifweconsiderthemasmembersinthe fulldiffeomorphismgroup Di! ( M ) ,wecanthinkofamap $ whichmaps Di! ( M ) to thespaceofdensities Dens ( M )= Di! ( M ) / Di! ( M ) ,where Di! ( M ) isthegroupof 60

PAGE 61

volume-preservingdiffeomorphisms.Themapisgivenasfollows: $ : $ f = & Jac (33) where Jac isthedeterminantoftheJacobianof .Thebenetofthismappingis thatthespaceofdensitiesisinnitedimensionalspherewithradius ( M ) ,whichisthe volumeof M .Thatis, M $ 2 ( ) d = ( M ) (34) Therefore,wecanadoptknowngeometricalformulationsofinnite(ornite)dimensional spheresinordertogetgeometricdistancesandgeodesicsbetweenpointsonthe sphere.Wecanalsoreducethespheretoaunitspherebyscaling $ with & ( M ) .More detailsaregivenin[ 37 ].Withthesetwopointsmappedonthesphere,wecanthen formulateagreatcirclebetweenthetwopointsonthespaceofdensitiesasfollows[ 44 ]: f ( t )= 1 sin( 2 ) [ sin( 2 ( t ) f 1 +sin( t ) f 2 ] (35) where f 1 = & Jac id f 2 = & Jac Itcouldbenotedthatiftherearevolume-preservingdifferomorphismsdeforming id ,then,bydenition,theseareprojectedtothesamepointin Dens ( M ) .Hence, theirgeodesicin Dens ( M ) isdegenerateditconsistsofjustonepointandthelifted diffeomorphicpathin Di! ( M ) isthenapathconsistingofonlyvolume-preserving diffeomorphisms. Inthisdissertation,diffeomorphismsareexpressedasvectorelds.Inthecaseof 2-D, ( t )=( x + U ( r t ), y + V ( r t )) where r =( x y ) ,andconventionally,wename ( U ( r t ), V ( r t )) thedeformationvectorelds.Then Jac ( t )=(1+ U x )(1+ V y ) ( U y V x =1+ U x + V y + U x V y ( U y V x (36) 61

PAGE 62

andwehave f ( t ) 2 = Jac ( t ) / ( M ) (37) Eq.( 36 )isthedeterminantoftheJacobianof ,and U i and V i arerst-orderderivatives withrespectto i { x y } inEq.( 36 ).The3-DversionofEq.( 36 )isderivedinAppendix D .WehavetoremarkthatEq.( 37 )isabilinearequation,anditwillbeatrilinearin3-D. Thenextstepistorecoverthedeformationelds, ( t ) ,fromEq.( 37 ). 3.1.2Step2.PathLiftingtotheSpaceofDiffeomorphisms Afterageodesiconthespaceofdensitiesisobtained,wehavetoliftthispathto Diff ( M ) .However,solvingEq.( 37 )cannotyieldauniquesolutionalgebraically.This isbecausetheequationisbilinear,andgeometrically,thiscorrespondstothefactthat Eq.( 37 )onlyrequires ( t ) tolieona Di! ( M ) -orbitin Di! ( M ) ,parametrizedbythe point f ( t ) inthesphere. Therefore,inthisstudy,theliftedpathin Di! ( M ) isevaluatedbyintroducing L 2 smoothnessofthedeformationvectoreldsovertimeasthemainregularization criterion.Eventhoughthereareotherchoicesfornormsforregularizationsuchas theSobolevnormthe L 2 normregularizationwaschosenforadoptionbecauseof itscomputationalsimplicity.Theliftingproblemnowturnsouttobethequadratic programming,withbilinearconstraintsfor2Dproblems. min 7 7 7 7 dU ( t ) dt 7 7 7 7 2 + 7 7 7 7 dV ( t ) dt 7 7 7 7 2 d dt s t Jac ( t )= f ( t ) 2 ( M ) (38) Theliftingproblemwillbecomethequadraticprogrammingwithtrilinearconstraintsfor 3Dproblems.ThisisdiscussedinAppendix D ThefollowingchapterconcernshowtooptimizeEq( 38 )numerically. 62

PAGE 63

Figure3-2. Theareaofthepolygonisgivenby 1 / 2 | X 30 X 21 | 3.2NumericalMethod WhileEq.( 38 )isformulatedforcontinuousvariablesof x y and t ,inpractice,we havetoworkwithdiscretepixelpointsandtime;itsdiscreteversioniswrittenasfollows: min U t ij V t ij / i j t > ? 7 7 7 7 7 U t +1 ij ( U t 1 ij 2 x 7 7 7 7 7 2 + 7 7 7 7 7 V t +1 ij ( V t 1 ij 2 y 7 7 7 7 7 2 @ A t x y s t .(1+ U t xij + V t yij + U t xij V t yij ( U t yij V t xij ) ( ( M )( f t ij ) 2 =0 (39) wherethesuperscript t andsubscripts i and j denotethediscretetimeindexandthe pixelpoints,respectively,andthesubscripts x and y intheconstraintsdenotethe rst-orderderivativesw.r.t. x and y respectively. 0 x and 0 y aresettobeone;therefore, thevolume ( M ) istheimagesize. Geometrically,thedeterminantoftheJacobian Jac ( t ) istheratioofthechangeof volumeelementsby ( t ) atdomainpointsandbytime t .Inthisdissertation,thischange isconsideredastheareachangeofpolygonmesheswithfourneighborpixelpoints replacedbythedeformationeldsasverticesforthenumericalscheme.Fig. 3-2 shows adeformedpolygonmesh,andtheareaofthepolygonisgivenby Area = 1 2 7 7 7 X 30 X 21 | 7 7 7 (310) 63

PAGE 64

Ifwesetvertices 0 1 2 and 3 as t ij t ij +1 t i +1 j and t i +1 j +1 with t ij =( i + U t ij j + V t ij ) in Fig. 3-2 respectivelythen X 30 and X 21 arewrittenasfollows: X 30 =(1+ U t i j ( U t ij ,1+ V t i j ( V t ij ) X 21 =( ( 1+ U t i j ( U t ij ,1+ V t i j ( V t ij ) (311) where i and j denote i +1 and j +1 respectively,and ( i j ) [1, H ( 1] [1, W ( 1] with thedomainsizeof H W .Aftervectorizingalldeformationvectorelds, X 30 and X 21 are rewrittenasfollows: X 30 =( 1 + C 30 U ,1+ C 30 V ) X 21 =( ( 1 + C 21 U ,1+ C 21 V ) (312) where U and V denotethevectorversionsof U t ij and V t ij at t ,respectively,and C 30 and C 21 arematriceswhichyieldthenitedifferencesdenedinEq.( 311 ). Jac mu ( t ) can thenbegivenby, Jac ( t )=1+0.5( C 30 ( C 21 ) U +0.5( C 30 + C 21 ) V +0.5( C 30 U 3 C 21 V ( C 30 V 3 C 21 U ) (313) where 3 denotesthecomponent-wiseproduct.Withthisnumericalscheme,thesquare rootoftheJacobiandeterminantisapointonthe ( H ( 1) ( W ( 1) dimensionalsphere. Afterbeingvectorizedandlinearized,Eq.( 39 )isrewrittenasaquadraticprogramming problemwithbilinearconstraints: min 1 2 U RU + 1 2 V RV + B u U + B v V s t h = c + C v V + C u U + D u ( V ) U =0 or h = c + C u U + C v V + D v ( U ) V =0 (314) Here, U =( U ( t =1),..., U ( t = T ( 1)) ,and B u =( U ( t =0), 0 U ( t = T )) where 0 isazerosvectorwiththelengthof H W ( T ( 1) ifwehavetimesequence t [0, T ] includingtwoboundaryconditions. V and B v aredenedinthesameway.InEq.( 314 ), 64

PAGE 65

R isthematrixforthequadraticcomponents,anditssizeis H W ( T ( 1) ; C u U and C v V correspondstothelineartermsof U and V inEq.( 313 )respectively,and D u ( V ) U (or D v ( U ) V )correspondstothebilineartermsinEq.( 313 )afterfactoringout U (or V ). Finally, c isthevectorizationof ( M )( f t ij ) 2 inEq.( 39 )plustheidentityvector.Thedetails ofentriesofthesematrixesareprovidedinAppendix D The3-Dformulationismoreinvolved.Inthiscase,thesquarerootofthedeterminant oftheJacobianisformulatedasthevolumeofhexahedralcells,andEq.( 314 )becomes aquadraticprogrammingwithtrilinearconstraints.DetailsareprovidedinAppendix E Tosolvethisproblem,WeadopttheaugmentedLagrangianmethodwithapenalty [ 38 40 ].Wesolvetheproblemrstwithaxed V (or U ),andwethen,solveiterativelyfor theotherone.Withaxed V ,theproblemisgivenasfollows: min 1 2 U RU + B u U + h + 1 2 c || h || 2 s t U R m (315) where R m istheLagrangianmultiplier, c $. ,and m = H W ( T ( 1) .Thisisan unconstraintedoptimizationproblem,anditcouldberewrittenasfollows: min 1 2 U ( R + c H u H u ) U +( B u + 1 2 c (( G u H u ) + H u G u )+( H u ) ) U s t U R m (316) where H u = D u + C x and G u = c .Thatis, h = H u U + G u =0 Alg. 2 descriptionofthisfulloptimizationmethodwiththebilinearconstraints. 3.3Experiment Thissectionconcernshowtoutilizetheproposedframeworkforthesolutionof practicalproblems:forexample,computingmissingvideosequencesandsolving classicationproblemswithinacardiacdataset.Inbothexamples,imagesare lled-inbetweenasourceandtargetimagesbydiffeomorphicpathscomputedwith theproposedframework.Anoveralldescriptionofthislling-inprocessisasfollows: rst,asourceandtargetimagesareprepared. I s and I t denotethesourceandtarget 65

PAGE 66

Algorithm2 OptimizationoftheaugmentedLagrangianoftheoriginalLagrangianwith bilinearconstraint S et k + 0 S et c k > 0 S et k = 0 R m S et $ > 1 I nitialize V 0 while || h || 2 < 1 do H k u + D U ( V k )+ C x A ( c k ) + R + c k H k U H k U B ( c k k ) + B U +0.5 c k (( G U H k U ) + H k U G U )+( k H k U ) S olve A ( c k ) U k = ( B ( c k k ) for U k k + k + c k ( H k U U k + G U ) c k + $ c k H k v + D V ( U k )+ C y A ( c k ) + R + c k H k V H k V B ( c k k ) + B V +0.5 c k (( G V H k V ) + H k V G V )+( k H k V ) S olve A ( c k ) V k = ( B ( c k k ) for V k k + k +1 k +1 = k + c k ( H k V V k + G V ) c k +1 = $ c k endwhile U + U k V + V k + k images,respectively.Second,aregistrationmapordiffeomorphismdenotedby betweenthetwoimagesiscomputedwitharegistrationframework.Wethenset id as thediffeomorphismwhichmapsthesourceimagetoitself.Inthisway,wecanprepare twodiffeomorphismstobeinterpolated.Next,wecomputeadiffeomorphicpath ( t ) joining id and withtheproposedalgorithm.Finally,withthiscomputeddiffeomorphic path ( t ) ,wecansimulateanimagesequence I ( t ) between I s and I t simplybywarping theimage I ( t )= I s & ( t ) 66

PAGE 67

A B Figure3-3. Thediagrammaticviewofcomputingmissingframesinavideoscene. I ( t ) representanimageateachtimepointt.Imageswithintheblueboxesare missingframes;imagesoutoftheboxesarekeyframes.A)Thegroundtruth. B)Acomputedpath. ( t ) denotesadiffeomorphicpathbetweenkeyframes. Themissingframesarecomputedby I & ( t ) 3.3.1TheEstimationofMissingVideoSequences Thepurposeofthissectionistoshowhowtheproposedframeworkworksforthe problemsofimageinterpolationandvideo-sequenceestimation.Towardthisaim,rst, animagesequenceisextractedfromavideoscene,andframesaredroppedfromthe sequenceremainingkeyframes.Next,missingframesbetweenkeyframesarecomputed asdescribedearlier.Fig. 3-3 illustratesthislling-inprocess. Theselled-inimagesequencesarecomparedwiththegroundtruthofsequences bymeasuringthepSNRofthecomputedonesw.r.t.thegroundtruths.Inaddition, theresultingsequencesarecomparedwithresultsfromtheLargeDeformation DiffeomorphicMetricMapping(LDDMM).Beforemovingtoproblemswithrealvideo datasets,asimpleexampleispresentedinFig. 3-4 .Therearetwocylinders:one ispositionedupright(Fig. 3-4A )andtheotherisa50-degree-rotatedversionofthe uprightcylinder(Fig. 3-4B ).FirstweevaluateadeformationeldsdeformingFig. 3-4A to 67

PAGE 68

Table3-1. TheaveragedpSNRscoresbetweengroundtruthsandestimatedimage sequenceswiththeproposedmethodandLDDMM Method TheproposedmethodLDDMM Rotationtest 46.24 dB 26.59 dB Bendigwoman 48.25 dB 46.2 dB Dancingwoman 45.94 dB Table3-2. Computingtimetohaveresultingdiffeomorphicpathswiththeproposed methodandLDDMM.Thereportedtimesfortheproposedmethodinclude registrationandinterpolationstages.Inthecaseofdancingwoman,itfailed tohaveasetofparametersofLDDMMyieldingaconvergentresult. Method TheproposedmethodLDDMM Rotationtest 24min40sec 42min20sec Bendigwoman 16min7sec 41min57sec Dancingwoman 12min30sec Noconvergence Fig. 3-4B .Inthisstep,boundariesofdomainhavebeenxedwhiledeformingthesource image;therefore,theresultingdeformationmustbenon-rigid.Twentyintermediate stagesoftherotationarethengenerated.Thegroundtruthissimplygeneratedby rotatingtheuprightcylinderby 50 / 19 degreesconsecutively,sothatwethenhave 20sequentialimages.Fig. 3-4 presentstheresultingsequences,providingvisual comparisons.Wecantelltheproposedframeworkyieldsmoreperceptionaland acceptableresultsthandoestheLDDMM.Inaddition,pSNRscoresbetweentheground truthandcomputedimagesaremeasured.ThepSNRscoresvsframesareplottedin Fig. 3-5 ,andtheaveragescoreisreportedinTable 3-1 .Inthistest, c 0 and $ inAlg.1 aresettobe 2 and 1.1 ,respectively.Inaddition,thecomputingtimeandthesizeof deformationarereportedinTable 3-2 andTable 3-3 respectively. Thefollowingsubsectionsconcernmorepracticalexamples:videofootageofa womanbendingherbody,andasecondvideosequenceofdancingwoman.Aswas describedearlier,inthesetestswedropvideoframesexceptkeyframesandwethen computethemissingframeswiththeproposedmethod. 68

PAGE 69

A B C D E F G Figure3-4. Theinputandresultingdiffeomorphicpathsoftherotatingcylinderexample. A)Originalcylinder.B)50degreerotatedoneofAwithdiffeomorphism plottedinC.C)ThecolorplotofdiffeomorphicdeformationeldsrotatingA toB.D)ThecolormapofC.E)Thegroundtruthoftherotatingsequences. F)Thepathofimageswrappedbythepathofdiffeomorphismsestimated withourmethod.G)ThepathofimagesbyLDDMM. 69

PAGE 70

A B C Figure3-5. pSNRscoreplot.Redlinesandbluelinesrepresentscoresfromourmethod andfromLDDMMrespectively.A)Therotatingcylindertest.B)Thebending womantest.C)Thedancingwomantest. Table3-3. Thesizeofdeformationineachgroupoftestdata.Sizeismeasuredusing theL2normofdisplacementvectors.Itsmaximum(Max),minimum(Min), average(Ave)andstandarddeviation(Std)arereportedinpixels. Testdata MaxMinAveStd Rotatingcylinder The35.16017.857.43 Bendingmotionsequence1 13.6903.332.75 Bendingmotionsequence2 5.9 01.4591.3 Dancingwomansequence1 11.1501.381.59 Dancingwomansequence2 20.1802.882.78 Dancingwomansequence3 9.9503.41.77 Dancingwomansequence4 14.5603.42.48 3.3.1.1Bendingmotion Inthissubsection,wepreparevideofootageofawomanwhoisbendingherbody. Thefootagehas29.75fps;sequentialimageshavebeenextractedwitha30-fpsframe rateso,thatthereare34imageframes.Thegroundtruthofthesequenceispresented inthetopandsecondrowsofFig. 3-7 .Therst,17th,andlastframesintheground truthsequencearechosenaskeyframes(sothatthereare15missingframesbetween keyframes).Therefore,inordertorecovertheoriginalmissingframes,weneedtoset T =15 inAlg.1.First,weevaluatediffeomorphicdeformationelds,whichdeformthe rstreferencetothemiddleandthemiddletothelast.Thecolorplotofthedeformation eldsarepresentedinFig. 3-6 .Withthisinputofdiffeomorphisms,themissingframes 70

PAGE 71

arecomputedwithintheproposedframeworkasdescribedintheprevioussectionandin Fig. 3-3 .Theresultingframesarethencomparedwiththegroundtruthandwithresults fromLDDMM,andthesecomparisonsaremadebothvisuallyandquantitatively. Foritsquantitativecomparisons,pSNRscoreswereevaluatedbetweentheground truthandthecomputedframes.Fig. 3-7 presentsvisualcomparisons.Whilewecannot observenoticeabledifferencesbetweenresultsfromtheproposedmethodandLDDMM withthesevisualcomparisons,thequantitativecomparisonsshowthattheproposed frameworkworksslightlybetterthanLDDMMinthislling-in-missing-frameproblem. TheplotsofthepSNRscoresinFig. 3-5B andtheaveragedscoresinTable 3-1 show thattheproposedmethodisslightlymoreaccuratethanLDDMMoverthewhole frames.ItshouldbealsobenotedthatwhileimplementingLDDMMinorderto haveapathbetweentherstandthe17thframes,moresubkeyframeshadtobe introducedinbetween,sincesearchingforasetofparametershasfailedtoproduce convergentLDDMMresults.Therearenomarkedlylargevisualdifferencesbetween theproposedmethodandLDDMM.Thismightbebecausethemotioninthetestmay notbecomplicatedandeasytobeinterpolated.Thecomputingtimeandthesizeof deformationarereportedinTable 3-2 andTable 3-3 ,respectively. 71

PAGE 72

A B C D E Figure3-6. Keyframesinthebendingmotionexamplesanddeformationvectorelds whichdeformthepreviouskeyframestothenextkeyframes.A-C) Keyframes.D-E)ColorplotofdeformationeldswhichdeformAtoBandB toCrespectively. 72

PAGE 73

A B C Figure3-7. Theground-truthandcomputedsequencesofmissingframesinthevideoofbendingmotion.Theimagesof outsideboxesarekeyframes,andtheimagesofinsideofboxescorrespondmissingframes.A)Theground truthofmissingframes.B)Framescomputedbytheproposedmethod.C)FramesestimatedbyLDDMM. 73

PAGE 74

A B C D E Figure3-8. Keyframesinthedancingmotionexamplesanddeformationvectorelds whichdeformthepreviouskeyframestothenextkeyframes.A)Keyframes. B-E)Colorplotsofdeformationeldsamongthereferenceimages 3.3.1.2Dancingmotion Thissubsectioncontainsmorecomplicatedmotionthandoestheprevious example.Asintheprevioussections,29imageframeshavebeensampledfromthe originalvideo.24frameshavedroppedout,leavingvekeyframes.Thesereference imagesandthedeformationeldsamongthemarepresentedinFig. 3-8 ;Fig. 3-9 presentsthewhole-imageframesofgroundtruth,andthecomputedimages.More keyframeshadtobeintroducedinthistestsincemorecomplicatedmultiplemotions areinvolved.Asdoneintheprevioussection,visualandquantitativecomparisonsare providedinFig. 3-9 andTable 3-1 ,respectively.Again,inFig. 3-9 ,theimagesoutsideof theboxesarekeyframes,andtheimagesinsidetheblue-andred-boxedframesare thegroundtruthandthecomputedimagesofmissingframes,respectively.ThepSNR scoresbetweenthegroundtruthandthecomputedframesfromtheproposedmethod areplottedinFig. 3-5C ,andtheaveragedpSNRscoreisreportedinTable 3-1 Inthistest,thecomparisonresultsw.r.tLDDMMarenotreported,sincetheresults failedtosearchforasetofparametersinordertoyieldconvergentresultswithLDDMM. 74

PAGE 75

ThecomputingtimeandthesizeofdeformationarereportedinTable 3-2 andTable 3-3 respectively. 75

PAGE 76

A B Figure3-9. Theground-truthandcomputedsequencesofmissingframesinthevideoofdancingmotion.Theimagesof outsideboxesarekeyframes,andtheimagesofinsideofboxescorrespondmissingframes.A)Theground truthofdancingframes.B)Framescomputedbytheproposedmethod. 76

PAGE 77

3.3.2CardiacMotionAnalysisandItsApplicationtoClassication Thepurposeofthissectionistoshowhowtheproposedframeworkforinterpolating diffeomorphismscouldbeappliedtothecardiac-motionanalysisforclassication purposes.Morespecically,itistoshowthatfeaturesextractedfromdiffeomorphic pathscanbeusefulinclassifyingischemiapatientsundergoingstemcelltreatment. 3.3.2.1Datapreparation CardiacMRIscansofpatientswithacutemyocardialinfraction(AMI)were prepared,andtheywererandomizedintooneoftwogroups:acontrolgroupanda treatmentgroup.Thetreatmentgroupreceivedendothelialprogenitorstemcell(EPCs) therapy;thecontrolgroupreceivedacell-freeinfusion.PriortoEPCtherapy,cardiac MRIexamswerecarriedouttoprovidebaseline;theseexamswererepeatedsix monthspost-treatment[ 45 ].Generally,cardiacanalysisinvolvesameasuringofthethe averagewallmotionorwallthickeningreportedforeachoftheindividualsegments ofhearts.Thissectionproposesapatient-classicationmethod(betweencontroland treatmentgroups)bymeasuringchangesindiffeomorphicpathsoversixmonths.These diffeomorphicpathsarethediastolicllingmotionoftheleftventricularmyocardium (LVM)computedusingtheproposedinterpolationframework.Fortheexperiments, expertshavesegmentedLVM'sfromthecardiacMRIimages.Thediastolicllingmotion ofLVMstartsatastagecalledend-systole(ES)andendsatend-diastole(ED).The examplesof2-Dslicesinthesestages,alongwiththeirsegmentations,areshownin Fig. 3-10 ;awholerangeofmotionbetweentheseESandEDstagesispresentedin Fig. 3-11 Usingacinemulti-slicesequenceusingaSiemensAvanto 1.5 T wholebodyMRI scanner,thecardiacdatawasacquiredalongtheshortaxisoftheheart.Twenty-ve phaseswerecollectedthroughoutthecardiaccycleusingaTrueFISP(abalanced coherentgradientsequenceusingFastImagingwithSteadystate-freePrecession),with anechotimeof1.13msandanapparentrepetitiontimeof 71.82 ms .Datawerecollected 77

PAGE 78

A B C D Figure3-10. 2-Dsliceimagesoftheleftventricularmyocardium(LVM)atend-systole (ES)andend-diastole(ED).A-B)ESandEDrespectively.C-D)LVM segmentationsofAandBrespectively. usingatotalof161phaseencodesintoamatrixof 256 151 ,yieldinganin-plane resolutionof 1.71875 1.71875 mm 2 .Slicethicknesswas 7 mm ,withasliceseparationof 10 mm .throughoutthecardiaccycle,coverageofthemyocardiumextendedfromabove thevalveplaneatthebase,tobelowtheapex. TheproposedmethodisusedtocomputethediastolicllingmotionofLVM betweenESandEDasithasbeendoneinthevideoexperiments.Thatis,The2-D segmentedLVMimagesatESandEDaresettobe I s and I t ,respectively.Theproposed methodthencomputes I ( t )= I s & ( t ) .Throughthiscomputation,wecanachievetwo benets.First,byassignthesamenumberofintermediatephasesbetweenESand EDforallpatients,thecardiaccyclecouldbenormalized.Secondly,wedonotneedto segmentLVMinallphasesbetweenESandED.ThesegmentationofLVMinbetweenis atime-consumingandtediousjob. However,beforewemovetotheclassicationproblem,theinterpolatedpaths betweenESandEDusingtheproposedmethodandLDDMMshouldbecomparedto thegroundtruthofthecardiacmotion.Forthiscomparison,eightintermediatetime pointsarechosenalongtheinterpolatedpath,andatthesepoints,comparisonsare performedagainstthegroundtruth.Thiscomparisonisrepeatedforvedifferent patients.AnexampleofresultinginterpolatedphasesispresentedinFig. 3-11 .For thequantitativecomparison,theDicecoefcientsbetweenthegroundtruthsand 78

PAGE 79

A B Figure3-11. Thediastolicllingmotionofa2-DsliceofLVMbetweenESandED.A) Thegroundtruth.B)Thecomputedmotionusingtheproposedmethod. Table3-4. TheDicecoefcientsofthecomputedheartphasesusingourmethodand LDDMM. AverageMaximumMinimum Ourmethod 0.9180.9960.868 LDDMM 0.880.9150.835 theinterpolatedLVMswhichhavebeenobtainedusingtheproposedmethodand LDDMMarereportedinTable 3-4 .Whencomparedwiththegroundtruth,theresult showsthatthecomputedpathofthesegmentedLVMusingtheproposedmethod arevisuallyacceptable;theproposedmethodoutperformsLDDMMbyusingtheDice coefcientastheevaluationmetric. 3.3.2.2Classifyingischemiapatients Thissectionconcernstheapplicationoftheproposedmethodinordertoclassify thetwogroupsofischemiapatientsdescribedintheprevioussections.Thestem-cell therapyusedinthetreated-patientgroupissupposedtoimprovethediastoliclling motionoftheleftventricularmyocardium(LVM);therefore,anapproachforclassication betweenthetwogroupsistomeasurehowthediastolicllingmotionofLVMchanges oversixmonths(bothwithandwithoutendothelialprogenitorstemcell(EPCs)therapy). TheperiodofLVMdiastolicllingmotionisfromtheend-systole(ES)totheend-diastole (ED),andtheimagesofsegmentedLVMatESandEDareassignedto I s 'sand I t 's 79

PAGE 80

A B Figure3-12. TheyellowboxesindicatetheROIchosenforthecardiacmotionanalysis. A)End-systole.B)End-diastole. respectively.Theproposedmethodthencomputesthediffeomorphicpath ( t ) of eightpatientsinthecontrolgroupand23inthetreatedgroup.Allcardiaccyclesare normalizedsuchthattherearetenintermediatephasesbetweenESandED; c 0 and $ inAlg.1aresettobe 2 and 1.1 ,respectively.Thepathsarethencomputedina slice-by-slicefashionwithinthechosenROI,showninFig. 3-12 Next,thecomputed ( t ) isexpressedindiscretesineandcosine(DSC)bases. Thebasesarechoseninsuchawaythatthe ( t ) x and ( t ) y of ( t )=( ( t ) x ( t ) y ) remainzeroalongboundariesparalleltoy-andx-axes,respectively.Allcoefcientsof ( t ) arethenputintoasinglevectorforeachpatient;themotionchangeofLVMover sixmonthsforeachpatientisthenmeasuredbythe L 2 distanceofDSCcoefcients betweenthebaselinescanandthe6-monthscan.Thisdistanceisusedasafeature forclassication.Inthisexperiment,theSupportVectorMachine(SVM)isadoptedfor theclassierwithapolynomialkernelfunctionofdegreethree.ThenumberofDSC basiselementsforeachcomponentof ( t ) ateachphaseiseleven;theleave-one-out crossvalidationmethodisused.AsshowninTable 3-5 ,theclassicationresultisquite remarkable,indicatingthattheproposedmethodoutperformsLDDMM. 80

PAGE 81

Table3-5. ClassicationscorewithSupportVectorMachine Validationmethod Classicationscore Theproposedmethod LDDMM Leave-one-out ControlTreatmentControlTreatment 88.89%82.62%66.67%69.56% 3.4Discussion Chapter 3 introducesanovelframeworkfortheinterpolationoftwogivendiffeomorphisms whicharetheidentitymapandadiffeomorphismdeformingoneimagetotheother. LDDMM-basedframeworkshavestudiedthisproblemintensively.However,thereare issuesofconcern: TheLDDMM-basedframeworkscontaintwodifferentoptimizationswithinasingle formulation:registrationandinterpolation.Thisjointframeworkmaybeinefcient dependingonimageswewanttoregister. TuningparametersintheLDDMM-basedframeworkisatime-consumingand tediousjob.Attimesitisimpossibletondasetofparametersyieldingconvergent results. Therefore,anidealapproachtosolvetheinterpolationproblemistoseparatethe interpolationtaskfromtheregistrationtask.However,duetotheinnitedimensionality ofdiffeomorphismsandthelackofknowledgeofmetricinthespace,formulatinga separatedframeworkisnoteasy.Inthisdissertation,wehavecircumventedthese difcultiesinthefollowingway.First,ourworkprojectsthetwogivendiffeomorphisms computedfromtheregistrationontothequotientspacethathasbeenshowntobea subsetoftheinnite-dimensionalsphere,anditcomputesthegeodesicconnectingthe twoprojectedpoints.Second,thisworkliftsthegeodesicpathonthespherebackto thespaceofdiffeomorphismsbyusingan L 2 -regularization,resultinginaquadratic programmingproblemwithbilinear/trilinearconstraintsthatcanbeefcientlyoptimized usingaugmentedLagrangianmethod. Theproposedmethodhasbeenappliedtotwodifferentproblems:rstforlling-in missingframesinavideo;second,forananalysisofcardiacmotionforthepurpose 81

PAGE 82

ofclassication.Intherstexperiment,resultsinTable 3-1 andTable 3-2 showthatthe proposedmethodismoreaccurateandfasterthanisLDDMMinestimatingtheground truthofadiffeomorphicpath.Inthevideoexperiments,wehavefoundthatsearching forasetofparametersinordertoyieldconvergentresultsinLDDMMisdifcultand casesensitive;attimes,thoseexperimentsfailed.Wehaveobservedthisthroughthe experimentwiththedancingwomanvideo.Ontheotherhand,theproposedmethodis lesssensitivetoparametersthanisLDDMM.Therearetwoparameterstobeadjusted: c 0 and $ inAlg.1.Ithasbeenlearnedthat 1 c 0 3 and $ / 1.1 workne.Thismeans thatadiffeomorphicpathjoiningtwodiffeomorphismscanbeachievedeasilyaslong aswehaveasatisfactorydiffeomorphismthroughaproperregistrationframework.We maybeabletoincreasetheaccuracyofresultsinLDDMMbyadoptinglandmark-based LDDMMframeworks.However,itrequiresthatwechooselandmarksfromthegiven pairofdiffeomorphismsorassumeallofthegivengridpointstobelandmarkswhich increasesthecomputationalburdenimmenselyforLDDMM,makingitaninfeasible choice.Finally,alltheseresultssupporttheargumentthattheproposedframeworkis moregeneralthanaretheLDDMM-basedframeworks.Inaddition,sincetheproposed frameworkiscompletelyseparatedfromregistrationframeworks,thismustbeanovel frameworkforsynthesizingdiffeomorphisms. Inthesecondexperiment,theproposedalgorithmisappliedtoanischemiapatient classicationproblem.Therearetwogroupsinthedataset:groupswithandwithout endothelialprogenitorstemcell(EPCs)treatment.Thereareeightand23patients inthecontrolgroupandthetreatedgroup,respectively.Thesegroupsareclassied throughanalyzingthediastolicllingmotionoftheleftventricularmyocardium(LVM). Thefeatureforthisclassicationhasbeenobtainedbymeasuringdifferencesbetween baselineand6-monthMRIscansofthemotion.Theproposedframeworkprovides practicalbenetsinthecardiacmotionanalysis,andalsoinitsapplicationtothe classicationproblem.First,thecardiaccyclecouldbenormalized,whichmaydiffer 82

PAGE 83

patient-by-patient.Second,theinterpolationbetweenend-systole(ES)andend-diastole (ED)usingtheproposedframeworkcansavelotsofeffortinsegmentingoutLVMs ineverycardiacphasebetweenESandED.TheresultsinTable 3-4 supportsthis argument.Finally,byquotientingoutthevolume-preservingdiffeomorphisms,wecan takeouttheeffectofthenuisanceparametersinthecardiac-motionanalysis.The classicationscoresinTable 3-5 reectthefactthat1.takingoutthevolume-preserving diffeomorphismsiscrucial,and2.thefeaturesextractedfromthediffoemorphicpaths arerelevantandusefulforthecardiac-motionanalysis(andforitsapplicationtothe classicationproblem). Intheseexperiments, L 2 regularizationhasbeenadoptedintheliftingstep.Even thoughitisnottheonlychoicewecanmake,thishasstillyieldedveryconvincing results. Theapplicationoftheproposedframeworkisextensive.First,wecanimprove thequalityofvideoswithlow-frameratesbyaddingmoreframes,aswehavedonein section 3.3.1 .Wecandothisalsoinanoppositeway.Inotherwords,wecanapply thistovideocompressionbysavingreferenceimagesanddeformationelds.Inthis dissertation,theexperimentsarelimitedin2-Ddata,butwecanextendthisto3-D; allformulationsin3-DareprovidedinAppendix E .Cardiac-motionanalysisin3-Dis especiallycrucialinordertohavemoredependableresults:forexample,withrespect totheclassicationproblem.Therefore,futureworkmustberelatedto3-Dmotion analyses.Anotherapplicationthatwemightconsiderisastudyoflongitudinaldatasets. Yetanotherdirectionfutureresearchcouldtakemightbeanadoptionofthenovel frameworktothestudyofgrowthpatternsofmultiplegroups,andtothestatistical analysiswhichcouldbederivedthroughresultingdiffeomorphicpaths. 83

PAGE 84

CHAPTER4 GENERALCONCLUSION Inthisresearch,wehavediscussedtwogeometricalapproachesforimage registrationandforsynthesizingdiffeomorphisms.Inchapter 2 ,theimagegraphwas introduced;itisutilizedtobuildaregistrationframeworkwiththreeproperties:symmetry, parametricinvarianceanduniversality.Thesepropertiesarecrucialinaregistration framework,sincewithoutthese,resultingimagecorrespondencesandresultsfrom applicationsuchasclassicationmaybebiased.Theimagegraphmakesthethree propertiespossibleinthenovelregistrationframework,becauseimagesarecompared intrinsicallythroughtheimagegraphrepresentation. Theframeworkisveriedthroughtwoexperiments.First,theproposedframework hasbeencomparedwiththeDemonsregistrationalgorithminestimatingtheground truthsofdeformationelds.Second,ithasbeenappliedtoaclassicationproblem. Inthisexperiment,theproposedframeworkwasutilizedtoclassifyYoung,Middle andOldgroupsoutofadatasetoftheMRIbrainscansof416subjects.Fromtheresults oftheseexperiments,wehavelearnedthat1.Theproposedmethodoutperformsthe commonregistrationalgorithminyieldingcorrectcorrespondences;and2.Theimage graphrepresentationiscrucialtoensurebetterclassicationresults.3.Ifmatching imagescontainorientationalinformation,thereorientationstepshouldbeinvolvedinour optimizationprocess.. Inchapter 3 ,wehavediscussedanovelframeworkforinterpolatingtwodiffeomorphisms. UnlikeLDDMM-basedframeworks,theproposedframeworkisformulatedseparately fromregistration.Therefore,theframeworkismoreefcientandgeneralizedthanare theLDDMM-basedframeworks. Ourproposedframeworkconsistsoftwosteps. Intherststep,twodiffeomorphismsareprojectedonthespaceofdensities,and thegeodesicconnectingtheprojectedpointiscomputed. 84

PAGE 85

Inthesecondstep,thegeodesiclineisliftedbacktothespaceofdiffeomorphisms. Thisstepisformulatedasthequadraticprogrammingwithbilinear/trilinear constraints,whichisoptimizedusingtheaugmentedLagrangianmethodwith penalty. Throughtheststep,severalthingsbearnoting:volume-preservingdiffeomorphisms areconsiderednuisanceparameters;thespaceofdensitiesisthecorresponding quotientspace;andtheproposedframeworkcaninduceasmoothpathinthespaceof diffeomorphismstakingouttheeffectofthenuisanceparameters. Tovalidatetheproposedalgorithm,twoexperimentsarecarriedout:lling-in missingvideoframesandcorrespondentananalysisofcardiacmotion.Intherst experiment,theresultingdiffeomorphicpathshaveprovedthattheproposedframework isabletosimulatethegroundtruthswithahighlevelofaccuracy.Inthesecond experiment,itwasarguedthatvolume-preservingdiffeomorphsmswithincardiac motionneedtobeconsideredasnuisanceparameters.Theclassicationresults withinthesecondexperimentsreectthefactthatthefeatureswhichareextracted fromthecardiacmotionsimulatedbytheproposedframeworkmustberelevantinthe classicationproblem,sincetheproposedframeworktakesoutthenuisanceparameters whilecomputingsmoothpathsforcardiacmotion. Theregistrationframeworkusedinchapter 3 isausualintensity-basedmethod. Thatis, + 1 + 2 ,and aresetto 1 0 ,and 0 ,respectively,inEq.( 211 )Eq.( 215 ); deformationeldsareexpressedwiththeDSCbasissothatwecanpreservethetotal volumeofdomainunderdeformation.However,ifwerecallthatthenovelregistration frameworkinPartIcanbeusedforshape-matchingaswell,combiningthetwonovel frameworkswhichareintroducedinthisdissertationmustyieldanovelframeworkfor computingadiffeomorphicpathofshapechange.Thiscanbeappliedtostudyevolution scenariaorsetsoflongitudinaldata.Forexample,tostudyAlzheimer'sdisease,we needtoinvestigatealongitudinaldatasetofbrainMRIscansofpatients.Inthiscase, thecombinedframeworkcouldprovideagrowthscenarioofapatient'sbrain,and 85

PAGE 86

wecanutilizetheresultingscenarioforthepurposeofdiagnosis.Also,wecando statisticalstudiesbyinvestigatingpatternsofgrowthamongdifferentpatients.Therefore, applyingthecombinedframeworktoshapechangewillinformthefuturedirectionofour research. 86

PAGE 87

APPENDIXA SUPPLEMENTALMATERIAL1FORCHAPTER2 Inthissection,weshowhowtoderiveEq.( 211 ).Weadoptinvariantconsistent constraintsinordertomakeourcostfunctionalsymmetricunderinversematching.In theworkthatfollows,weusethesamenotationdenedinEq.( 211 ).First,wewritethe costfunctionalasasummationofmatchingcostfunctionalsintwooppositedirections asfollows: E (( f 1 f 2 ), & )= 1 Dist ( f 1 f 2 & ) ) + 1 d 1 + 2 Dist ( f 2 f 1 & & ) ) + 2 d 2 (A1) where : 1 $ 2 and & : 2 $ 1 ,and ) + d isvolumeform. Ifwechangethecoordinatessystemfrom 2 to 2 inthesecondtermontheleft handsideofEq( A1 ),thetermcanberewrittenasfollows: 1 Dist ( f 2 & & 1 f 1 ) & + 2 ( 2 & & 1 ) J 1 d 1 (A2) andifwerequirethatEq( A1 )issymmetricmatching, & 1 shouldbe and J 1 = J & whichisdeterminantoftheJacobian. Finally,thecostfunctionisgivenasfollows: Cost ( f 1 f 2 ; )= Dist ( f 1 f 2 & )( ) + 1 + & + 2 ( 2 & ) J & ) d 1 (A3) 87

PAGE 88

APPENDIXB SUPPLEMENTALMATERIAL2FORCHAPTER2 Inthefollowingthreesections,weshowhowtofactoroutbasisfunctionsandtheir derivativesinEq.( 238 )throughEq.( 240 ).Hereweshowonlydetailsofderivatives w.r.t. A nml .Derivativesw.r.t. B nml and C nml canbedoneinthesameway,describedas follows.Inthissection,weconsiderEq.( 238 ),andinthefollowing,sectionsEq.( 239 ) andEq.( 240 )willbeconsidered. Eq.( 238 )couldberewrittenas ) ) A nml Dist ( X 1 R X 2 & )= ) ) A nml ( ( U 2 + V 2 + X 2 ) B CD E Part 1 + dist ( I 1 R I 2 & ) B CD E Part 2 ) (B-1) In Part 1 ofEq.( 250 ),factoringoutbasisfunctionsisverystraightforward;therefore,we onlyshowhowtodothisin Part 2 below. Part 2 canbegroupedintotwoparts,asfollows: ) ) A nml dist ( I 1 R I 2 & )= # ) A R ) A nml ( 2 # ) C R ) A nml ( & ) B CD E +2( & ) # B ) ( & ) ) A nml ( 2 # C R ) ( & ) ) A nml B CD E !! (B-2) FactoringoutbasisfunctionsinEq.( B-2 )requiresanevaluationofEq.( 251 ),Eq.( 252 ) andEq.( 253 )specically.First,weneedtoevaluatedderivativesoftheJacobian matrix, J .Given J = ( ( ( ( ) 1+ U u U v U w V u 1+ V v V w W u W v 1+ W w + + + + (B-3) therstorderderivativew.r.t. A nml isgivenas ) J ) A nml = ( ( ( ( ) # # u $ U nml # # v $ U nml # # w $ U nml 0 0 0 0 0 0 + + + + (B-4) 88

PAGE 89

Oncewenoticethat ) J ) A nml u i R = ( ( ( ( ( ( ) u i R 1 ) ) u $ U nml + u i R 2 ) ) v $ U nml + u i R 3 ) ) w $ U nml B CD E i 0 0 + + + + + + (B-5) wherethesubscript i of issamewith i of u i R ,andoncewerewrite J 1 as J 1 = ( ( ( ( ) J 1 11 J 1 12 J 1 13 J 1 21 J 1 22 J 1 23 J 1 31 J 1 32 J 1 33 + + + + = ( ( ( ( ) J 1 1 J 1 2 J 1 3 + + + + (B-6) thenwecanfactoroutbasisfunctionsinEq.( 254 )as ) u i R ) A nml =( J 1 1 ( ( u # i R J 1 1 ) u i R ) i (B-7) andinthesameway, ) u i R ) B nml =( J 1 2 ( ( u # i R J 1 2 ) u i R ) i (B-8) ) u i R ) C nml =( J 1 3 ( ( u # i R J 1 3 ) u i R ) i (B-9) WithEq.( B-7 ),wecanrewriteEq.( 251 )andEq.( 252 )asfollows: ) A R ) A nml 7 7 7 7 i j = a ij u # i R u j R F u # i R ( J 1 1 ( ( u # j R J 1 1 ) u j R ) j +( J 1 1 ( ( u # i R J 1 1 ) u i R ) # u # j R i G = J 1 ij j + J 2 ij i (B-10) and ) C R ) A nml 7 7 7 7 i j = b ij u # i R v j ( J 1 1 ( ( u # i R J 1 1 ) u i R ) # v j i = Z ij i (B-11) 89

PAGE 90

where a ij = $ (2 / det ( i R + j R )) 3 / 2 b ij = $ (2 / det ( i R + # j R )) 3 / 2 Finally,wefactoroutbasisfunctionsandtheirderivativesinEq.( B-2 )asfollows: = i 5 ) A R ) A nml 6 ij j ( 2 i 5 ) C R ) A nml 6 ij ( & ) j = i ( J 1 ij j + J 2 ij i ) j ( 2 i Z ij ( & ) j i = i [( J 1 ij u j R 1 + J 2 ij u i R 1 ) j ( 2 Z ij ( & ) j u i R 1 ] ) $ U nml ) u + i [( J 1 ij u j R 2 + J 2 ij u i R 2 ) j ( 2 Z ij ( & ) j u i R 2 ] ) $ U nml ) v + i [( J 1 ij u j R 3 + J 2 ij u i R 3 ) j ( 2 Z ij ( & ) j u i R 3 ] ) $ U nml ) w (B-12) !! =2(( & ) # B ( # C R ) ) ) u 7 7 7 7 u + U $ U nml (B-13) 90

PAGE 91

APPENDIXC SUPPLEMENTALMATERIAL3FORCHAPTER2 Inthissection,weshowhowtoevaluateEq( 239 )andEq.( 240 ).Ifwenotate thedeterminantoftheJacobian J % denedinEq.( 214 )as J % ,itsderivativew.r.t. A nml denedinEq.( 235 ),iswrittenas ) ) A nml J % =(1+ V v + W w + V v W w ( V w W v ) ) $ U nml ) u +( V + wW + u ( V u W w ( V u ) ) $ U nml ) v +( V u W v ( V v W u ( W u ) ) $ U nml ) w (C-1) Then,evaluationofEq( 239 )andEq.( 240 )becomesstraightforwardasfollows: ) ) A nml & + 2 ( 2 & ) J % = 1 2 ) + 2 )+ 2 ) u 7 7 7 7 u + U $ U nml J % + ) + 2 ) J % ) A nml (C-2) and ) ) A nml ( J % ( 1)log( J % )= 8 log( J % )+1 ( 1 J % 9 ) ) A nml J % (C-3) 91

PAGE 92

APPENDIXD SUPPLEMENTALMATERIAL1FORCHAPTER3 ThissectionconcernsthedetailsofEq.( 314 ).Werecallthedeterminantof Jacobianandthequadraticprogrammingproblemwithbilinearconstraintsinbelow: Jac ( t )=1+0.5( C 30 ( C 21 ) U +0.5( C 30 + C 21 ) V +0.5( C 30 U 3 C 21 V ( C 30 V 3 C 21 U ) (C-1) and min 1 2 U RU + 1 2 V RV + B u U + B v V s t h = c + C v V + C u U + D u ( V ) U =0 or h = c + C u U + C v V + D v ( U ) V =0 (C-2) Ifthesizeofanimagedomainis M N ,and i [1, M ] and j [1, N ] arepixelpoints,the entriesof C 21 and C 30 aregivenasfollows: C 21 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p + N 1 q = p +1 0 otherwise C 30 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p 1 q = p +1+ N 0 otherwise (C-3) where p = j +( i ( 1) N D U ( V ) and D V ( U ) arewrittenasfollows: D U ( V )= C 21 V C 30 ( C 30 V C 21 D V ( U )= C 30 U C 21 ( C 21 U C 30 (C-4) Theproductsymbol inEq.( C-4 )isdenedasfollows: ( Q C !" ) p q =( Q ) p ( C !" ) p q (C-5) Here, Q isavectorandtherepeatedindex, p ,meansjustcomponentproductnot summation. 92

PAGE 93

APPENDIXE SUPPLEMENTALMATERIAL2FORCHAPTER3 FigureE-1. AHexahedralcell Thissectionconcernshowtoformulatethe3DversionofEq.( 314 ).Forthe3D formulation,weneedtoconsider Jac asthevolumeofahexahedralcellwhichis presentedinFig. E-1 andexpressthisasatrilinearequationsothatthisservesas trilinearconstraints.Inthissection,weuseaformulationforthevolumeofatetrakis hexahedron(TH)zone.Moredetailscanbefoundin[ 46 ].First,wediscusshowto formulatethetrilinearconstraints.Thevolume, Vol isgivenby 12 Vol = H X 71 + X 60 7 7 7 X 72 7 7 7 X 30 I + H X 60 7 7 7 X 72 + X 50 7 7 7 X 74 I + H X 71 7 7 7 X 50 7 7 7 X 74 + X 30 I (C-1) InEq. C-1 X !" denotes X ( X ,and [ A | B | C ] isthetripleproduct H A 7 7 7 B 7 7 7 C I = 7 7 7 7 7 7 7 7 7 7 AxBxCx AyByCy AzBzCz 7 7 7 7 7 7 7 7 7 7 (C-2) 93

PAGE 94

InordertouseEq.( C-1 )fortheformulationof Jac ateachgridpoint, ( i j k ) ,let's assignthevertexesinFig. E-1 asfollows: 0 =( i j k ), 1 =( i +1, j k ), 2 =( i j +1, k ), 3 =( i +1, j +1, k +1), 4 =( i j k +1), 5 =( i +1, j k +1), 6 =( i j +1, k +1), 7 =( i +1, j +1, k +1). (C-3) Then, X 0 canbewrittenasafunctionofdeformationvectoreldsdenotedby ( U i j k V i j k W i j k ) asfollows. X 0 = X i j k =( i + U i j k j + V i j k k + W i j k ) =( i + U 0 j + V 0 k + W 0 ). (C-4) Inthissection, i j and k areintegerswhosedomainsare [1, M ] [1, N ] and [1, L ] respectively,foranimagewiththesizeof M N L Therestof X 'scanbewritteninthesamefashionandall X !" inEq. C-1 isrewritten asfollows. X 71 =( U 71 ,1+ V 71 ,1+ W 71 ), X 72 =(1+ U 72 V 72 ,1+ W 72 ), X 74 =(1+ U 74 ,1+ V 74 W 74 ), X 60 =( U 60 ,1+ V 60 ,1+ W 60 ), X 50 =(1+ U 50 V 50 ,1+ W 50 ), X 30 =(1+ U 30 ,1+ V 30 W 30 ), (C-5) Inordertouse Jac astrilinearconstraint,weneedtolinearizeallentitiesandexpress Vol asatrilinearequation.Therefore,werewriteEq.( C-5 )asfollows. X 71 =( C 71 U 1 + C 71 V 1 + C 71 W ), X 72 =( 1 + C 72 U C 72 V 1 + C 72 W ), X 74 =( 1 + C 74 U 1 + C 74 V C 74 W ), X 60 =( C 60 U 1 + C 60 V 1 + C 60 W ), X 50 =( 1 + C 50 U C 50 V 1 + C 50 W ), X 30 =( 1 + C 30 U 1 + C 30 V C 30 W ), (C-6) InEq.( C-6 ),theboldtextof U V and W denotethevectorexpressiionsof U i j k V i j k and W i j k respectivelywhoseindexingisgivenby U p = k +( j +( i 1) N 1) L = U i j k ,and C !" isa 94

PAGE 95

matrixwiththesizeof ( M ( 1)( N ( 1)( L ( 1) MNL whichequates U !" to C !" U C !" is asparseuppertriangularmatrixanditsentriesaregivenasfollows. C 71 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p + NL 1 q = p +( N +1) L +1 0 otherwise C 72 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p + L 1 q = p +( N +1) L +1 0 otherwise C 74 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p +1 1 q = p +( N +1) L +1 0 otherwise C 60 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p 1 q = p + L +1 0 otherwise C 50 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p 1 q = p + NL +1 0 otherwise C 30 pq = 1 2 2 2 2 3 2 2 2 2 4 ( 1 q = p 1 q = p +( N +1) L +1 0 otherwise (C-7) Ouraimistoexpress Vol asatrilinearequationof U V and W asfollows. Vol = 1 + C V V + C W W +( C U + D U ( V )+ D U ( W )+ E U ( V W )) U = 1 + C U U + C W W +( C V + D V ( U )+ D V ( W )+ E V ( U W )) V = 1 + C U U + C V V +( C W + D W ( U )+ D W ( V )+ E W ( U V )) W (C-8) Afteratediousalgebraicwork,allmetricsinEq.( C-8 )aregivenby C U =0.25( ( C 71 + C 72 + C 74 ( C 60 + C 50 + C 30 ) (C-9) C V =0.25( C 71 ( C 72 + C 74 + C 60 ( C 50 + C 30 ) 95

PAGE 96

(C-11) C W =0.25( C 71 + C 72 ( C 74 + C 60 + C 50 ( C 30 ) (C-12) 12 D U ( V )= ( (2 C 30 + C 74 ) V C 71 +(2 C 30 + C 74 ) V C 72 (C-13) +(2 C 60 ( C 72 ( 2 C 50 + C 71 ) V C 74 ( ( C 30 +2 C 74 ) V C 60 +(2 C 74 + C 30 ) V C 50 +(2 C 71 + C 60 ( 2 C 72 ( C 50 ) V C 30 12 D u ( W )= ( (2 C 50 + C 72 ) W C 71 +(2 C 50 + C 72 ) W C 74 (C-14) +(2 C 60 + C 71 ( 2 C 30 ( C 74 ) W C 72 ( (2 C 72 + C 50 ) W C 60 +(2 C 72 + C 50 ) WC 30 +(2 C 71 + C 60 ( 2 C 74 ( C 30 ) W C 50 12 D V ( W )= ( (2 C 60 + C 71 ) W C 72 ++(2 C 60 + C 71 ) W C 74 (C-15) +(2 C 50 + C 72 ( 2 C 30 + C 72 ) W C 71 ( (2 C 71 + C 60 ) W C 50 +(2 C 71 + C 60 ) W C 30 +(2 C 72 + C 50 ( 2 C 74 ( C 30 ) W C 60 Theproductsymbol intheaboveisdenedasfollows. ( Q C !" ) p q =( Q ) p ( C !" ) p q (C-16) Here, Q isavectorandtherepeatedindex, p ,meansjustcomponentproductnot summation. D V ( U ) D W ( U ) and D W ( V ) canbeobtainedbyequatingtheseto ( D U ( U ) ( D U ( U ) 96

PAGE 97

and ( D V ( V ) respectively.ThethirdorderternsinEq.( C-8 )aremoreinvolved 12 E U ( V W )= ( V 72 3 W 30 ( V 30 3 W 72 + V 50 3 W 74 ( V 74 3 W 50 + V 50 3 W 30 ( V 30 3 W 50 ) C 71 +( V 74 3 W 60 ( V 60 3 W 74 + V 30 3 W 71 ( V 71 3 W 30 + V 30 3 W 60 ( V 60 3 W 30 ) C 72 +( V 60 3 W 72 ( V 72 3 W 60 + V 60 3 W 50 ( V 50 3 W 60 + V 71 3 W 50 ( V 50 3 W 71 ) C 74 +( V 72 3 W 74 ( V 74 3 W 72 + V 50 3 W 74 ( V 74 3 W 50 + V 72 3 W 30 ( V 30 3 W 72 ) C 60 +( V 74 3 W 60 ( V 60 3 W 74 + V 74 3 W 71 ( V 71 3 W 74 + V 30 3 W 71 ( V 71 3 W 30 ) C 50 +( V 71 3 W 72 ( V 72 3 W 71 + V 60 3 W 72 ( V 72 3 W 60 + V 71 3 W 50 ( V 50 3 W 71 ) C 30 (C-17) InEq.( C-17 ), 3 denotesthecomponent-wiseproduct,andothertwothirdorderterms inEq.( C-8 )areobtainedbyrotating U V and W incyclicorder.Thatis, E V ( W U )= E U ( W U ) and E W ( U V )= E U ( U V ) Thenwecanformulatethequadraticprogrammingwithtrilinearconstraints. min 1 2 U RU + 1 2 V RV + 1 2 W RW + B u U + B v V + B w W s t h = c + C V V + C W W +( C U + D U ( V )+ D U ( W )+ E U ( V W )) U =0 or h = c + C U U + C W W +( C V + D V ( U )+ D V ( W )+ E V ( U W )) V =0 or h = c + C U U + C V V +( C W + D W ( U )+ D W ( V )+ E W ( U V )) W =0, (C-18) HowtoconstructvectorsandmatrixesinEq.( C-18 )areprovidedinsection 3.2 InEq.( C-18 ), c issettobe 1 ( f ( t ) 2 ( ) ,where f ( t ) isthegeodesicinthespaceof densitiesand ( ) isthevolumeofdomain.Inthisnumericalscheme, ( ) isequalto ( M ( 1)( N ( 1)( L ( 1) 97

PAGE 98

REFERENCES [1] N.A.Sochen,R.Kimmel,andR.Malladi,Ageneralframeworkforlowlevelvision, IEEETrans.Imag.Proc. 7 (3),310-318(1998) [2] R.Kimmel,andN.A.Sochen,Geometric-variationalapproachforcolorimage enhancementandsegmetation,Scale-Space,294-305(1999) [3] R.Kimmel,R.Malladi,andN.A.Sochen,Imagesasembeddedmapsandminimal surfaces:movies,color,texture,andvolumetricmedicalimages,IJCV 39 (2),111-129 (2000) [4] Y.GurandN.Sochen,"Coordinate-baseddiffusionoverthespaceofsymmetric positive-denitematrices",in"VisualizationandProcessingofTensorFIeldsAdvancesandPerspectives",pp325-340,Springer-Verlag,BerlinHeidelberg,(2009). [5] D.Seo,andB.C.Vemuri,Complexdiffusiononimagegraphs,ICIP2969-2972 (2009) [6] D.SeoandB.C.Vemuri,D.SeoandB.Vemuri,Complexdiffusiononscalarand vectorvaluedimagegraphs,EMMCVPR98-111(2009) [7] P.PeronaandJ.Malik,Scale-spaceandedgedetectionusingansiotropicdiffusion, IEEETrans.PatternAnal.MachineIntell.12:629639(1990) [8] L.Younes,Shapesanddiffeomorphisms,AppliedMathematicalSciencesVol.171, Springer,Heisenberg,(2010) [9] M.F.Beg,M.I.Miller,A.Trouv eandL.Younes,Computinglargedeformationmetric mappingsvisgeodesicowsofdiffeomorphisms,IJCV 61 (2),139-157,(2005) [10] M.I.Miller,A.Trouv eandL.Younes,Geodesicshootingforcomputationalanatomy, J.Math.ImagingVis., 24 (2),209-228,(2006) [11] L.Younes,F.Arrate,andM.I.Miller,Evolutionsequationsincomputational anatomy,Neuroimages 45 ,S40-S50,(Mar2009) [12] S.Sommer,M,Nielsen,F.Lauze,andX,Pennec,Amulti-scalekernelbundlefor LDDMM:towardssparsedeformationdescriptionacrossspaceandscales, IPMI 624-635,(2011) [13] S.Sommer,F.Lauze,M,Nielsen,andX,Pennec,KernelBundleEPDiff:Evolution EquationsforMulti-scaleDiffeomorphicImageRegistration, SSVM ,677-688,(2011) [14] M.P.DoCarmo RiemannianGeometry 1992. [15] L.Younes,A.Qiu,R.L.WinslowandM.I.Miller,Transportofrelationalstructures ingroupsofdiffeomorphisms,J.Math.ImagingVis., 32 (1),41-56,(2008) 99

PAGE 99

[16] G.E.ChristensenandH.J.Johnson,Consistentimageregistration,IEEETans. Med.Imag. 20 (7),568-582(2001) [17] H.J.JohnsonandG.E.Christensen,Consistentlandmarkandintensity-based imageregistration,IEEETrans.Med.Imag 21 (5),450-461(2002) [18] G.E.ChristensenandH.J.Johnson,Invertibilityandtransitivityanalysisfor non-rigidimageregistration.J.Electron.Imag. 12 (1),106-117(2003) [19] G.E.Christensen,R.D.Rabbitt,andM.I.Miller,Deformabletemplatesusinglarge deformationkinematics,IEEETrans.Imag.Proc. 5 (10),1435-1447,(1996) [20] S.Joshi,B.Davis,M.Jomier,andG.Gerig,Unbiaseddiffeomorphicatlas constructionforcomputationalanatomy,NeuroImage 23 ,151-160,(2004) [21] Y.Cao,M.I.Miller,S.Mori,R.L.WinslowandL.Younes,Diffeomorphicmatching ofdiffusiontensorimages Proc .CVPR,17-22(2006) [22] H.D.Tagare,D.Groisser,andO.Skrinjar,Symmetricnon-rigidregistration:a geometrictheoryandsomenumericaltechniques,J.MathImagingVis 34 ,61-88 (2009) [23] G.Cheng,B.C.Vemuri,P.R.Carney,andT.H.Marcei,Non-rigidregistrationof highangulardiffusionimagesrepresentedbygaussianmixtureelds.MICCAI(1), 190-197(2009) [24] S.Kurtek,E.Klassen,Z.DingandA.Srivastava,Anovelriemannianframeworkfor shapeanalysisof3Dobjects, Proc. CVPR1625-1632,(2010). [25] S.Kurtek,E.Klassen,Z.Ding,S.Jacobson,J.Jacobson,M.AvisonandA. Srivastava,Parameterization-InvariantShapeComparisonsofAnatomicalSurfaces, IEEETrans.Med.Imaging, 30 ,3,849-858,(2011) [26] N.Litke,M.Droske,M.RumpfandP.Schr 'oder,Animageprocessingapproachto surfacematching,SymposiumonGeometryProcessing,207-216,(2005) [27] D.C.Alexander,C.Pierpaoli,P.J.Basser,andJ.C.Gee,SpatialTransformations ofdiffusiontensormagneticresonanceimages, IEEETrans. Med.Imag.2011, 1131-1139,2001 [28] M.MoakherandP.G.Batchelor,SymmetricPositive-DeniteMatrices:From GeometrytoApplicationsandVisualization,in"VisualizationandProcessingof TensorFIelds,Springer-Verlag,NewYork,285-298,2006. [29] D.S.Tuch,DiffusionMRIofcomplextissuestructure.PhDthesis,MIT,2002 [30] B.JianandB.C.Vemuri,Auniedcomputationalframeworkfordeconvolutionto reconstructmultiplebersformdiffusionweightedMRI, IEEETrans.Med.Imaging26 11 ,1464-1471(2007) 100

PAGE 100

[31] I.Yanovsky,P.M.Thompson,S.Osher,A.D.Leow,Topologypreserving Log-unbiasednonlinearimageregistration:Theoryandimplementation, IEEE CVPR,1-8(2007) [32] R.R.Coifman,S.Kafon,A.B.Lee,M.Maggioni,B.Nadler,F.Warner,andS.W. Zucker,Geometricdiffusionsasatoolforharmonicanalysisandstructuredenitionof data:Diffusionmaps,PNAS 102 (21),7426-7432(2005) [33] R.R.CoifmanandS.Lafon,DiffusionMaps,Appl.Comput.Harmon.Anal.textbf21 ,5-30(2006) [34] D.S.Marcus,T.H.Wang,J.parker,J.G.Csernansky,J.C.Moris,R.L.Buckner, OpenAccessSeriesofImagingStudies(OASIS):Cross-SectionalMRIDatain Young,MiddleAged,Nondemented,andDementedOlderAdults.JournalofCognitive Neuroscience(2007) [35] Y.Xie,B.C.Vemuri,andJ.Ho,"StatisticalAnalysisofTensorFields",InInt.Conf. onMedicalImageComputingandComputerAssistedIntervention(MICCAI),2010. [36] T.Chen,A.Rangarajan,andB.C.Vemuri,"CAVIAR:ClassicationviaAggregated RegressionandItsApplicationinClassifyingtheOASISBrainDatabase",InIEEE InternationalSymposiumonBiomedicalImaging,pp.1337-1340,2010. [37] Khesin,B.,Lenells,J.,Misiolek,G.,Preston,S.C.:GeometryofDiffeomorphism Groups,CompleteIntegrabilityandGeometricStatistics.arXiv:1105.0643(2011) [38] D.P.Bertsekas,ComputationalOptimizationandLagrangeMultiplierMethods, AcademicPress,NewYork,1982 [39] D.P.Bertsekas,NonlinearProgramming,seconded.,AthenaScientic,Belmont, Massachusetts,1999 [40] L.CarotenutoandG.Raiconi,Ontheminimizationofquadraticfunctionswith bilinearconstraintsviaaugmentedLagrangians,Journalofoptimizationtheoryand applications, 55 ,1,23-36,(Oct.1987) [41] G.D.Pillo,G.Liuzzi,S.LucidiandL.Palagi,Anexactaugmentedlagrangian functionfornonlinearprogrammingwithtwo-sidedconstraints,Computational OptimizationandApplication, 25 ,57-83,(2003) [42] X.Du,L.ZhangandY.Gao,AclassofaugmentedLagrangiansforequality constraintsinnonlinearprogrammingproblems,Appliedmathematicsand computation, 172 ,644-663,(2006) [43] P.E.GillandD.P.Robinson,Aprimal-dualaugmentedLagrangian,Computational OptimizationandApplication, 47 ,1-15,(2010) [44] A.Srivastava,I.JermynandS.H.JoshiRiemanniananalysisofprobabilitydensity functionswithapplicationinvision, Proc .CVPR,18-23(2007) 101

PAGE 101

[45] Petersen,J.W.,Forder,J.R.,Thomas,J.D.,Moye,L.A.,Lawson,M.,Loghin,C., Traverse,J.H.,Baraniuk,S.,Silva,G.,Pepine,C.J.:QuanticationofMyocardial SegmentalFunctioninAcuteandChronicIschemicHeartDiseaseandImplications forCardiovascularCellTherapyTrials:aReviewfromtheNHLBI-CardiovascularCell TherapyResearchNetwork,JACCCardiovascImaging, 4 (6),71-9,(2011) [46] J.Grandy,Efcientcomputationofvolumeofhexahedralcells.No. UCRL-ID128886.LawrenceLivermoreNationalLab.,CA(UnitedStates),1997. 102

PAGE 102

BIOGRAPHICALSKETCH DohyungSeoreceivedhisM.S.degreefromtheDepartmentofPhysicsatthe UniversityofFloridain2006.HethentransferredtotheDepartmentofElectricaland ComputerEngineeringin2006HereceivedhisM.S.degreefromtheDepartment ofElectricalandComputerEngineeringattheUniversityofFloridain2009.Hethen chosetopursuehisPh.Dinthesamedepartment.In2008,SeojoinedtheCenterfor ComputerVision,Graphics,andMedicalImaging,andheconductresearchonimage denoising,imageregistration,anddiffeomorphisms.Seohopestomakesignicant contributionstosocietyintheeldsofcomputervisionandmedicalimaging,among others,throughinnovativeandgroundbreakingresearch. 103