Stability metrics for optic radiation tractography

Delft University of Technology

Stability metrics for optic radiation tractography

Towards damage prediction after resective surgery

Meesters, Stephan; Ossenblok, Pauly; Wagner, Louis; Schijns, Olaf; Boon, Paul; Florack, Luc; Vilanova

Bartroli, Anna; Duits, Remco

DOI

10.1016/j.jneumeth.2017.05.029

Publication date

2017

Document Version

Final published version

Published in

Journal of Neuroscience Methods

Citation (APA)

Meesters, S., Ossenblok, P., Wagner, L., Schijns, O., Boon, P., Florack, L., Vilanova Bartroli, A., & Duits, R.

(2017). Stability metrics for optic radiation tractography: Towards damage prediction after resective surgery.

Journal of Neuroscience Methods, 288, 34-44. https://doi.org/10.1016/j.jneumeth.2017.05.029

Important note

To cite this publication, please use the final published version (if applicable).

Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons. Takedown policy

Please contact us and provide details if you believe this document breaches copyrights. We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

ContentslistsavailableatScienceDirect

Journal

of

Neuroscience

Methods

j ou rn a l h o m epa g e :w w w . e l s e v i e r . c o m / l o c a t e / j n e u m e t h

Stability

metrics

for

optic

radiation

tractography:

Towards

damage

prediction

after

resective

surgery

Stephan

Meesters

a,b,∗

_,

_Pauly

_Ossenblok

a,c

_,

_Louis

_Wagner

a

_,

_Olaf

_Schijns

a,d

_,

_Paul

_Boon

a

_,

Luc

Florack

b

_,

_Anna

_Vilanova

e,b

_,

_Remco

_Duits

b

a_{AcademicCenterforEpileptologyKempenhaeghe&MaastrichtUniversityMedicalCenter,Netherlands} b_{DepartmentofMathematics&ComputerScience,EindhovenUniversityofTechnology,Netherlands} c_{DepartmentofBiomedicalEngineering,EindhovenUniversityofTechnology,Netherlands} d_{DepartmentofNeurosurgery,MaastrichtUniversityMedicalCenter,Netherlands}

e_{DepartmentofMathematicsandComputerScience,DelftUniversityofTechnology,Netherlands}

h

i

g

h

l

i

g

h

t

s

•Thealignmentofstreamlinesisquantifiedbyfiber-to-bundlecoherencemeasures.

•ReliableML-TPdistancemeasurementbyremovalofspurious(deviating)streamlines.

•ParameterestimationtoremovespuriousstreamlinesandtoretaintheMeyer’sloop.

•ThevalidityofML-TPdistanceisestimatedbypreandpostoperativeORcomparisons.

•ThestabilitymetricsarepromisingtorelateORdamagetoavisualfielddeficit.

a

r

t

i

c

l

e

i

n

f

o

Articlehistory:

Received15February2017

Receivedinrevisedform25April2017 Accepted31May2017

Availableonline23June2017 Keywords:

Opticradiation Meyer’sloop

Diffusionmagneticresonanceimaging Fibertractography

Epilepsy Neurosurgery

a

b

s

t

r

a

c

t

Background:Anaccuratedelineationoftheopticradiation(OR)usingdiffusionMRtractographymay reducetheriskofavisualfielddeficitaftertemporalloberesection.However,tractographyisproneto generatespuriousstreamlines,whichdeviatestronglyfromneighboringstreamlinesandhinderareliable distancemeasurementbetweenthetemporalpoleandtheMeyer’sloop(ML-TPdistance).

Newmethod:Stabilitymetricsareintroducedfortheautomatedremovalofspuriousstreamlinesnear theMeyer’sloop.Firstly,fiber-to-bundlecoherence(FBC)measurescanidentifyspuriousstreamlinesby estimatingtheiralignmentwiththesurroundingstreamlinebundle.Secondly,robustthresholdselection removesspuriousstreamlineswhilepreventinganunderestimationoftheextentoftheMeyer’sloop. Standardizedparameterselectionisrealizedthroughtest–retestevaluationofthevariabilityinML-TP distance.

Results:ThevariabilityinML-TPdistanceafterparameterselectionwasbelow2mmforeachofthehealthy volunteersstudied(N=8).Theimportanceofthestabilitymetricsisillustratedforepilepsysurgery can-didates(N=3)forwhomthedamagetotheMeyer’sloopwasevaluatedbycomparingthepre-and post-operativeORreconstruction.Thedifferencebetweenpredictedandobserveddamageisintheorder ofafewmillimeters,whichistheerrorinmeasuredML-TPdistance.

Comparisonwithexistingmethod(s):Thestabilitymetricsareanovelmethodfortherobustestimateof theML-TPdistance.

Conclusions:Thestabilitymetricsareapromisingtoolforclinicaltrialstudies,inwhichthedamageto theORcanberelatedtothevisualfielddeficitthatmayoccurafterepilepsysurgery.

©2017TheAuthors.PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBYlicense (http://creativecommons.org/licenses/by/4.0/).

∗ Correspondingauthorat:AcademicCenterforEpileptologyKempenhaeghe& MaastrichtUniversityMedicalCenter,Netherlands.

E-mailaddress:s.p.l.meesters@tue.nl(S.Meesters).

1. Introduction

Withdiffusiontensorimaging(DTI)themorphologyofbrain

tissue,andespeciallythewhitematterfiberbundles,canbe

inves-tigated in vivo (Mori, 2007), offering new possibilities for the

evaluation of brain disorders and preoperative counseling. The

opticradiation(OR)isacollectionofwhitematterfiberbundles

http://dx.doi.org/10.1016/j.jneumeth.2017.05.029

whichcarriesvisualinformationfromthethalamustothevisual

cortex(Rubinoetal.,2005).Numerousstudies(Yogarajahetal.,

2009;Taokaetal.,2005;Chenetal.,2009;Winstonetal.,2012; Boriusetal.,2014;Jamesetal.,2015)haveaccomplishedto

recon-structtheORwithDTI,bytrackingpathwaysbetweenthelateral

geniculate nucleus (LGN) and theprimary visual cortex.In the

curvedregionoftheOR,configurationswithmultiplefiber

orienta-tionsappear,suchascrossings,becausewhitemattertractsofthe

temporalsteminterminglewiththefibersoftheMeyer’sloop(Kier

etal.,2004).Therefore,itisespeciallychallengingtoreconstructthe

Meyer’sloop,whichisthemostvulnerablebundleoftheORincase

ofsurgicaltreatmentofepilepsyinwhichpartofthetemporallobe

isremoved(Jamesetal.,2015).However,alimitationofDTIisthat

itcanextractonlyasinglefiberdirectionfromthediffusionMRI

data.

Withtheadventofmulti-fiberdiffusionmodelsithasbecome

possibletodescriberegionsofcrossingfiberssuchasthehighly

curvedMeyer’sloop.Tractographybasedonconstrained

spheri-caldeconvolution(CSD)(Tournieretal.,2007;Descoteauxetal.,

2009)hasbeenshowntohavegoodfiberdetectionrates(Wilkins

etal.,2015)andhasbeenappliedinseveralstudiestoreconstruct

theOR(Limetal.,2015;Martínez-Herasetal.,2015).Furthermore,

probabilistictractographyisconsideredsuperiorincomparisonto

deterministictractographyforresolvingtheproblemofcrossing

fibersintheMeyer’sloop(LiljaandNilsson,2015).The

probabilis-tictrackingresultsbetweentheLGNandthevisualcortexfora

healthyvolunteerareillustratedinFig.1.Thetrackingresultsare

showninacompositeimagealongwithotherbrainstructuressuch

astheventricularsystem.

However,acommonoccurrenceintractogramsobtainedfrom

probabilistic tractography are spurious (deviating) streamlines.

Spuriousstreamlinesarebydefinitionnotwell-alignedwith

neigh-boring streamlines and may hinder the measurement of the

distancebetweenthetemporalpoletothetipoftheMeyer’sloop

(ML-TPdistance).AnaccuratemeasurementoftheML-TPdistance

isrequiredforestimatingthepotentialdamagetotheORafter

tem-poralloberesection(TLR).Methodshavebeenproposedforthe

identificationandremoval ofspurious streamlines,forexample

basedonoutlierdetection(Yeatmanetal.,2012;Martínez-Heraset

al.,2015;Khatamietal.,2016),basedonthepredictionofdiffusion

measurementsbywhole-brainconnectomics(Pestillietal.,2014),

orbasedontheuncertaintyinthemaineigenvectorofthediffusion

tensor(Parkeretal.,2003).Mostofthesemethodsforreducing

spu-riousstreamlinesarebasedondensityestimationinR3_.Incontrast,

inthecurrentstudyfiber-to-bundlecoherence(FBC)tractometry

measuresareemployedthat arebasedondensityestimationin

thespaceofpositionsandorientationsR3_×S2_{.Thestability}

met-ricsintroducedinthisstudyarebasedontheFBCmeasures.These

metricsprovideareliableORreconstructionthatisrobustunder

stochasticrealizationsofprobabilistictractography.Toachievea

reliablereconstructionofthefullextentoftheMeyer’sloop,an

appropriateselectionofstreamlinesisrequiredsuchthat

spuri-ousstreamlinesareremovedwhilepreservingstreamlinesthatare

anatomicallymorelikelytoexist.ForthispurposetheFBC

param-eter isestimatedbased onthemeasuredvariability in ML-TP

distance.Herewerespectana-prioriconstraintonthemaximal

ML-TPdistancevariabilityforatest–retestprocedureon

stream-linetrackinganddeterminethecorrespondingminimalthreshold

selected ontheFBCmeasures.Thisthresholdremovesaminimal

amountofspuriousstreamlineswhileallowingforastable

estima-tionoftheML-TPdistance.

Inthecurrentstudythevalidityofthedistancemeasurements

isevaluatedbasedonpre-andpost-operativecomparisonsofthe

reconstructed ORof patientswhounderwent aTLR.It is

inves-tigated whetheritis feasibletoassess pre-operativelyfor each

individualpatientthepotentialdamagetotheORasanadverse

eventoftheplannedTLR.Thedeviationbetweentheprediction

of thedamage totheORand themeasureddamage in a

post-operativeimageiscompared,givinganindicationoftheoverall

errorindistancemeasurement.

Themaincontributionsofthispaperare:

• Quantificationofspuriousstreamlines.WeprovideFBCmeasures

thatquantifyhowwell-alignedastreamlineiswithrespectto

neighboringstreamlines.

• Stability metrics for the standardized removal of spurious

streamlinesneartheanteriortipoftheMeyer’sloop.

• Robust estimation of the variability in ML-TP distance by a

test–retestevaluation.

• DemonstrationoftheimportanceoftheFBCmeasuresby

ret-rospective predictionof thedamageto theORbasedon

pre-andpost-operativereconstructionsoftheORofepilepsysurgery

candidates.

Fig.1. Left:AnexampleofthereconstructionresultoftheORusingprobabilistictractographyfromanaxialview.Asinserts,close-upsareshownoftheanteriortipsofthe reconstructionsoftheORfromacoronalview.Right:Thetrackingresultsareshownforthesamevolunteerinacompositeimagealongwithotherbrainstructuressuchas theventricularsystem.TheML-TPdistancemeasurementisindicated.

2. Materialsandmethods

2.1. Subjects

Eighthealthy volunteerswithoutanyhistory ofneurological

orpsychiatric disorders were includedin ourstudy. All

volun-teersweremaleandintheagerangeof21–25years.Furthermore,

threepatientswereincludedwhowerecandidatesfortemporal

lobeepilepsysurgery.Foreachpatientastandardpre-and

post-operativeT1-weightedanatomical3D-MRIwasacquired.Patient

1(46/F)wasdiagnosedwitharightmesiotemporalsclerosisand

hadarightTLR,includinganamygdalohippocampectomy.Patient

2(23/F)wasdiagnosed witha leftmesiotemporalsclerosisand

hadanextendedresectionofthelefttemporalpole.Lastly,Patient

3(38/M)wasdiagnosedwithacavernomalocatedinthebasal,

anteriorpartofthelefttemporallobeandhadanextended

lesionec-tomy.Allpatientshadpre-andpost-operativeperimetrycarried

outbyconsultantophthalmologists.Thestudywasapprovedby

theMedicalEthicalCommitteeofKempenhaeghe,andinformed

writtenconsentwasobtainedfromallsubjects.

2.2. Dataacquisition

Datawasacquiredona3.0Tmagneticresonance(MR)scanner,

usinganeight-elementSENSEheadcoil(Achieva,PhilipsHealth

Care,Best,TheNetherlands).AT1-weightedscanwasobtainedfor

anatomicalreferenceusingaTurboFieldEcho(TFE)sequencewith

timingparametersforechotime(TE=3.7ms)andrepetitiontime

(TR=8.1ms).Atotalof160sliceswerescannedwithanacquisition

matrixof224×224withisotropicvoxelsof1×1×1mm,leadingto

afieldofviewof224×224×160mm.Diffusion-weightedimaging

(DWI)wasperformedusingtheSingle-ShotSpin-EchoEcho-Planar

Imaging(SE-EPI)sequence.Diffusion sensitizing gradientswere

applied,accordingtotheDTIprotocol,in32directionswitha

b-valueof 1000s/mm2 _{in additiontoan imagewithoutdiffusion}

weighting.Atotalof60sliceswerescannedwithanacquisition

matrixof112×112withisotropicvoxelsof2×2×2mm,leading

toafieldofviewof224×224×120mm.ASENSEfactorof2and

ahalfscanfactorof0.678wereused.Acquisitiontimewasabout

8minfortheDWIscanand5minfortheT1-weightedscan.The

maximaltotalstudytimeincludingsurveyimageswas20min.

2.3. Datapreprocessing

ThepreprocessingoftheT1-weightedscanandDWIdatais

out-linedinFig.2(top-leftbox).Alldatapreprocessingisperformed

usinga pipelinecreatedwithNiPype(Gorgolewskietal.,2011),

whichallowsforlarge-scalebatchprocessingandprovides

inter-facestoneuroimagingpackages(FSL,MRtrix). TheT1-weighted

scanwasfirstalignedtotheAC-PCaxisbyaffinecoregistration

(12degrees-of-freedom)totheMNI152templateusingtheFMRIB

SoftwareLibraryv5.0(FSL)(Jenkinsonetal.,2012).Secondly,affine

coregistration,consideredsuitableforwithin-subjectimage

reg-istration,was appliedbetweentheDWI volumestocorrect for

motion.Eddycurrentinduceddistortionswerecorrectedwithin

thePhilipsAchievascanningsoftwareanddidnotrequirefurther

post-processing.TheDWIb=0volumewassubsequentlyaffinely

coregisteredtotheaxis-alignedT1-weightedscanusing

normal-izedmutual information, and theresulting transformation was

appliedtotheotherDWIvolumes.TheDWIvolumeswere

resam-pledusinglinearinterpolation.Aftercoregistration,thediffusion

orientationswerereorientedusingthecorresponding

transforma-tionmatrices(LeemansandJones,2009).

2.4. Probabilistictractography

Probabilistictractography of theOR(outlinedin Fig.2,

top-middle box) is based on the Fiber Orientation Density (FOD)

function,firstdescribedbyDescoteauxetal.(2009).With

prob-abilistic tractography, streamlines are generated between two

regionsofinterest(ROIs):theLGN,locatedinthethalamus,andthe

primaryvisualcortex(seeFig.1).TheLGNwasdefinedmanuallyon

theaxialT1-weightedimageusinganatomicalreferences(lateral

andcaudaltothepulvinarofthethalamus)(Fujitaetal.,2001)using

asphereof4mmradius,correspondingtoavolumeof268mm3_.

Theipsilateralprimaryvisualcortexwasmanuallydelineatedon

theaxialandcoronalT1-weightedimage.Theprimaryvisualcortex

ROI’susedinthisstudyhaveanaveragevolumeof1844mm3_.

TheFODfunctiondescribestheprobabilityoffindingafiberata

certainpositionandorientation(Tuch,2004).Inthecurrentstudy

theFODfunctionisestimatedusingCSD,whichisimplementedin

theMRtrixsoftwarepackage(Tournieretal.,2012).During

track-ing,thelocalfiberorientationisestimatedbyrandomsamplingof

theFODfunction.IntheMRtrixsoftwarepackage,rejection

sam-plingisusedtosampletheFODfunctioninarangeofdirections

restrictedbyacurvatureconstraintimposedonthestreamlines.

StreamlinesareiterativelygrownuntilnoFODfunctionpeakcan

beidentifiedwithanamplitudeof10%ofthemaximumamplitude

oftheFODfunction(Jeurissenetal.,2011;Tournieretal.,2012).In

MRtrixtracking,20,000streamlinesaregenerated,whichprovides

agoodbalancebetweencomputationtimeandreconstruction

abil-ity.Astepsizeof0.2mmandaradiusofcurvatureof1mmwere

used.Thesesettingsarereasonableforourapplicationof

recon-structingtheORandarerecommendedbyTournieretal.(2012).

TheFODfunctionwasfittedwithsixsphericalharmonic

coeffi-cients,whichissuitablefortheDTIscanningprotocolusedinthis

study.

AnatomicalconstraintsareappliedwhenreconstructingtheOR

inordertopreventtheneedformanualpruningofstreamlines

andtoreduceasubjectivebias.Firstly,streamlinesarerestricted

withintheipsilateralhemisphere.Secondly,fibersoftheORare

expectedtopassoverthetemporalhornoftheventricularsystem

(Sincoffetal.,2004).Theventricularsystemismanuallydelineated

usingITK-SNAPimagesegmentationsoftware(Yushkevichetal.,

2006).Streamlinesthatcrossthroughtheareasuperior-laterally

tothetemporalhornareretained.Thirdly,anexclusionROIis

cre-atedmanuallyofthefornixtoremovestreamlinesthatcrossthis

region,whichisincloseproximitytotheLGNandMeyer’sloop.

Furthermore,in ordertoremove long anatomicallyimplausible

streamlines,themaximumlengthofthestreamlinesissetto114

mmbasedonafiber-dissectionstudyoftheORbyPeltieretal.

(2006).

2.5. Quantificationofspuriousstreamlines

Thestabilitymetricstoidentifyspuriousstreamlinesare

out-linedinFig.2,top-rightbox.Thesemetricsareusedtoprovidea

reconstructionoftheORthatisrobustagainstthepresenceof

spu-riousstreamlines,whichoccurespeciallyneartheanteriortipof

theMeyer’sloopasshowninFig.1(left).Theapplicationofthese

metricsisimportanttoobtainastablemeasurementoftheML-TP

distanceasindicatedinFig.1(right).

TheFiber-to-BundleCoherenceFBCmeasure,providingthebasis

ofthestability metrics,is a quantitativemeasureof streamline

alignmentandisusedforremovingspuriousstreamlines.Spurious

streamlinesare(partially)poorlyalignedwithsurrounding

stream-linesinthestreamlinebundle,whichisillustratedschematicallyin

Fig.2. Aschematicoverviewoftheanalysisproceduresfollowedtoreconstructtheopticradiation(OR)andtoacquirearobustestimateoftheMeyer’sLooptoTemporal Pole(ML-TP)distance.Thestagesinwhichdataareprocessedareindicatedbythedashedboxes.Thereddashedboxesindicatethenewcontributionsofthestudy.The varioussoftwarepackagesarecolor-coded.Theinputsofthepipelineareadiffusion-weightedimaging(DWI)datasetandananatomicalT1-weightedMRIimage.Outputs ofthepipelineareshownwithdouble-headedarrows.Abbreviations:FOD,fiberorientationdensity;CSD,constrainedsphericaldeconvolution;ROI,regionofinterest;FBC, fibertobundlecoherence;RFBC,relativeFBC.(Forinterpretationofthereferencestocolorinthisfigurelegend,thereaderisreferredtothewebversionofthisarticle.)

Fig.3. Top:Thefiber-to-bundlecoherenceFBCmeasureisdeterminedviakernel den-sityestimation.ABrownianmotionkernelisused(shownleft),whichisdefinedon thespaceofpositionsandorientations.Thestreamlinesarecolor-codedaccording totheirFBCmeasure,scaledfromhigh(blue)tolow(white).Bottom:TheRFBCis computedusingaslidingwindowofsize˛andproducesasinglevalueforeach streamline.(Forinterpretationofthereferencestocolorinthisfigurelegend,the readerisreferredtothewebversionofthisarticle.)

to5Dcurvesbyincludingthelocalorientationofthetangenttothe

streamline.Aliftedstreamlineicanbewrittenas

i=



(yk_i,nk_i)∈R3_× S2 | k=1,...,Ni



, (1)

whereyandnarethepositionandorientationofastreamline

ele-ment,Niisthenumberofpointsinthestreamlineandidenotes

theindexwithinthestreamlinebundle=∪N

i=1{i}.Toincludea

notionofalignmentbetweenneighboringstreamlinetangents,we

embedtheliftedstreamlinesintothedifferentiablemanifoldofthe

rigid-bodymotionLiegroupSE(3).Withinthisdifferential

struc-ture,ameasureisdefinedthatquantifiesthealignmentofanytwo

liftedstreamlinepointswithrespecttoeachotherinthespaceof

positionsandorientationsR3_×S2_{(Mumford,1994;CittiandSarti,}

2006;Duitsetal.,2007).Inordertocomputethismeasure,

ker-neldensityestimationisappliedusinga(hypo-elliptic)Brownian

motionkernel(seeFig.3,topleft).Thekernelsusedinthekernel

densityestimationhaveaprobabilisticinterpretation:theyarethe

limitingdistributionofrandomwalkersinR3_×S2_thatrandomly

moveforwardorbackward,randomlychangetheirorientation,but

cannotmovesideways(DuitsandFranken,2011;Portegiesetal.,

2015).TheFBCmeasureresultsfromevaluatingthekerneldensity

estimatoralongeachelementofallliftedstreamlines,shownin

Fig.3(topright)wheretheFBCiscolor-codedforeachstreamline.

AspuriousstreamlinecanbeidentifiedbyalowFBCthatoccurs

anywherealong itspath.For thispurpose, a scalarmeasurefor

theentirestreamlineisintroduced,calledtherelativeFBC(RFBC),

whichcomputestheminimumaverageFBCinaslidingwindow

alongthestreamlinei∈relativetothebundle.TheRFBCfor

astreamlineiiscalculatedaccordingto

RFBC˛(i,)=

AFBC˛(i,)

AFBC() . (2)

ThenumeratorAFBC˛_(

i,)givestheminimumaverageFBCof

anysegmentoflength˛alongthestreamlinei.Thedenominator

AFBC()isusedfornormalizationandistheaverageFBCofallthe

streamlinesinthebundle,computedovertheentirelengthofeach

streamline.Thesegmentlength˛wasdeterminedempiricallyas

2mm(correspondingto10streamlinepointswhenusinga

step-sizeof0.2mm),whichisconsideredsmallenoughtocharacterize

localdeviationsofthestreamlinebutcontainsenoughstreamline

pointsforstablequantificationoflocalFBC.Fora formal

defini-tionofthenumeratoranddenominatorinEq.(2),seeEqs.(A.5)

and(A.6)inAppendixA,respectively.Furtherdetailsregardingthe

implementationofFBCmeasures,whichincludesseveral

optimiza-tionstepssuchaspre-computedlookuptablesfortheBrownian

2.6. Standardizedparameterselection

Tocontroltheremovalofspuriousstreamlinesthethreshold

parameterisintroduced,whichisdefinedasthelowerbound

criteriononRFBCthatretainsastreamline.Moreprecisely,every

streamlineithatmeetstheconditionRFBC˛(i,)≥isretained.

However,acarefulselectionofthisthresholdisrequiredinorderto

preventanunderestimationofthefullextentoftheMeyer’sloop.A

methodisintroducedforthestandardizedselectionoftheminimal

thresholdselectedthroughtest–retestevaluationofthevariability

inML-TPdistance.Tothisend,probabilistictractographyoftheOR

isperformedmultipletimes,followedbythecomputationofthe

RFBCmeasureineachrepetition.Subsequently,aparametersweep

isperformedinwhichisvariedbetween0≤≤maxwheremax

correspondstothestatewhereallstreamlinesareremovedfrom

.Duringeverystepoftheparametersweep,theML-TPdistance

iscalculatedforalltest–retestrepetitionsbycomputingthe

Haus-dorffdistance(RockafellarandWets,2005)betweenthetemporal

poleandtheOR.Using thesedistancemeasurements,themean

andthestandarddeviation(variability)oftheML-TPdistanceare

determinedforeachvalueof.

TheprocedureisillustratedforahealthysubjectinFig.4,

show-ingthemeanand standarddeviationof theML-TPdistancefor

increasingvaluesof.Initially,ahighvariabilityisseenat=0,

indicatingthepresenceofspuriousstreamlinesneartheanterior

tipoftheMeyer’sloop.At=0.075mostspuriousstreamlinesare

removedandavariabilityintheorderofseveralmillimetersisseen.

Thevariabilityrisesandfallsduring0.1≤≤0.3.Astableregion

isobtainedat≈0.3,howeveratthispointtoomanystreamlines

havebeendiscardedaccordingtotheconditionRFBC˛_(

i,)≥and

therebytheML-TPdistancewillbeoverestimated.Inorderto

esti-matetheminimalthresholdselected,inwhichtheML-TPdistance

isneitherunder-noroverestimated,amaximumissetforthe

vari-abilityof2mm.Thismaximumisbasedonthemaximalaccuracy

of2–5mmthatmaybeachievedduringresectivesurgery.Inthe

selectionprocedure,issetatthefirstoccurrenceoflowvariability,

i.e.

selected=min{>0 | ()≤2mm,()=0,()>0} (3)

where()denotesthestandarddeviationinML-TPforthechosen

.Aftercrossingthe2mmthresholdonvariability,selectedisplaced

onthelocalminimumof().Usingthisprocedure,inthe

exam-pleshowninFig.4theML-TPisestimatedfor=0.075at36mm.

ThisML-TPdistanceiswithintherangeof22–37mmasreported

byEbelingetal.(1988),whoperformedadissectionstudyon25

humancadavers.

Forthepatientsstudied,thedistancemeasurementoutcomes

arecomparedtothepredicteddamageoftheORaftersurgery,as

outlinedinFig.2(bottomrow,reddashedbox).Theresectionarea

ismanuallydelineatedinthepost-operativeT1-weightedimage

usingITK-SNAP(Yushkevichetal.,2006).Theresectionlengthis

measuredfromthetemporalpole,attheanteriortipofthemiddle

sphenoidfossa,uptotheposteriormarginoftheresection.The

predicteddamageisdeterminedbythedistancebetweenthe

pre-operativeML-TPdistanceandtheresectionlength.Thedifference

betweenthepredicteddamageandtheobserveddamage,given

bythedistancebetweenpre-andpost-operativeML-TPdistances,

isnamedthemarginoferror.Themarginoferrorindicatesthe

maximalerrorindistancemeasurements,whichincludesboththe

variabilityinprobabilistictractographyandunaccountedsources

oferrorsuchasbrainshiftordistortions.

2.7. Opensourcesoftware

The methodology for the robust reconstruction of the OR

(outlined in Fig. 2) is available as an open source software

Table1

ListedaretheML-TPdistancesestimatedfortheleftandrighthemispheresofthe healthyvolunteersstudied(N=8)andthecorrespondingselectedvaluesfortheFBC thresholdingparameter.

Volunteer ML-TPdistance 

Left(mm) Right(mm) Left(–) Right(–)

1 36.4±1.5 32.1±1.3 0.075 0.15 2 30.0±0.6 27.8±1.0 0.13 0.14 3 33.4±1.5 23.5±0.9 0.2 0.35 4 34.9±1.7 31.4±0.2 0.45 0.1 5 36.8±1.4 32.2±1.0 0.075 0.33 6 28.3±0.3 25.8±0.6 0.025 0.28 7 32.3±0.4 23.4±1.1 0.15 0.05 8 22.5±0.5 30.7±1.0 0.125 0.18

package. TheNiPype basedpipeline for thebasicprocessing of DW-MRI data, tractography, and FBC measures is available at

https://github.com/stephanmeesters/DWI-preprocessing.Anopen

source implementation of the FBC measures for the reduction

of spurious streamlines described in Appendices A and B is

available inthe DIPY (Diffusion Imagingin Python) framework

(Garyfallidis et al., 2014) or as a C++ stand-alone application at https://github.com/stephanmeesters/spuriousfibers.

Visualiza-tion wasperformedin the opensource vIST/etool (Eindhoven

UniversityofTechnology,ImagingScience&TechnologyGroup,

http://sourceforge.net/projects/viste/).

3. Results

3.1. RobustestimationofML-TPdistance

Theeffectoftheremovalofspuriousstreamlinesonthe

ML-TPdistancemeasurementusingtheFBCmeasuresisdemonstrated

foreighthealthyvolunteers.ForeachvolunteerthemeanML-TP

distanceanditsstandarddeviationarelistedinTable1fortheleft

andrighthemisphere,togetherwithitscorrespondingtest–retest

variability.TheadditionalvalueoftheFBCmeasuresforarobust

ML-TPdistancemeasurementisfurtherevaluatedforthreepatients

whounderwentaTLR.

Theparameter estimationbased ontest–retestevaluation is

illustratedinFig.5forthereconstructedORofthelefthemisphere

fortheeighthealthyvolunteersstudied,showingfor arangeof

parameter(0–0.6)thestandard deviation(left) andthemean

(right)oftheestimatedML-TPdistance.Thetest–retestevaluation

wasperformedwith10repeatedtractogramsoftheOR,whichwas

empiricallydeterminedtobeagoodbalancebetweengroupsize

andcomputationtime.Forallvolunteersevaluated,ahigh

stan-darddeviationoftheML-TPdistance(over2mm)wasobservedat

lowvaluesof(0.0–0.05),whichindicatesthepresenceofspurious

streamlineswithaverylowRFBC.ThecorrespondingmeanML-TP

distancereflectslargejumpsforanincreaseofthevalueoffrom

0to0.05,showinganaverageincreasefortheeighthealthy

volun-teersof8mm.Foreachhealthyvolunteertheselectedisselected

accordingtoEq.(3).TheselectedcorrespondstoameanML-TP

dis-tancethatisdepictedbythearrowsinFig.5(right)fortheeight

healthyvolunteersstudied.Aftertheinitialhighvariabilityofthe

ML-TPdistance,astableregionoccurredforallhealthyvolunteers

inwhichthestandarddeviationwasbelow2mm.Thehealthy

vol-unteers1,5and4indicatedregionsofinstabilityforrelativelyhigh

valuesof.Thiscanbeattributedtogapswithinthereconstructed

ORwitha lowernumber of streamlinescompared tothemain

streamlinebundle.Lastly,itcanbeobservedthatforvolunteer4

theselectedislargecomparedtotheotherhealthyvolunteers.

However,forthisvolunteerthemeanML-TPdistanceisstablefrom

=0.15onwardandthereforedoesnotreflectanoverestimationof

Fig.4.BoxplotshowingthemeanandstandarddeviationoftheestimatedML-TPdistancesfortest–retestevaluationofthereconstructionoftheORforanexamplehealthy volunteer.AsweepfromlowtohighisperformedtoevaluatetheeffectofremovingstreamlinesonthestabilityoftheestimatedML-TPdistance.

Fig.5. ShownistheparameterestimationforthereconstructedleftORoftheeighthealthyvolunteersstudied.Left:ThestandarddeviationoftheML-TPdistanceisshown asafunctionof.Foreachhealthyvolunteerasuitablechoiceofismadeatthepointwherethestandarddeviationfirstdropsbelowthethresholdof2mmandreaches alocalminimum,shownbytheblackdottedline.Right:TheestimatedmeanML-TPdistanceisshownasafunctionof.Theselectedforeachvolunteerisindicatedbyan upwardspointedarrow,indicatedalongwiththevaluesoftheassociatedestimateoftheML-TPdistance.

OnthegroupleveltheML-TPdistanceslistedinTable1areon

average31.7±4.7mmforthelefthemisphereand28.4±3.8mm

fortherighthemisphere.Themeanvariabilityinprobabilistic

trac-tographyontheindividuallevelforthegroupofhealthyvolunteers

is1.0mmand0.9mmfortheleftandrighthemispheres,

respec-tively.LargedeviationsinML-TPdistancewereobservedbetween

theleftandrighthemispheres,especially,forvolunteers3,7and8.

3.2. Pre-andpost-operativecomparisons

TheimportanceoftherobustML-TPdistancemeasurementis

illustratedfor threepatientswho underwent resectiveepilepsy

surgery.Fig.6displaysthepre-operative(firstandlastcolumns)

andpost-operativereconstructions(secondandthirdcolumns)of

theORandindicatesforbothhemispherestheestimatedML-TP

distances (firstandsecond column).Givenis alsotheresection

length(thirdcolumn)andthepre-operativereconstructionofthe

ORalongwiththepredicteddamage,indicatedbytheredcolored

streamlines(fourthcolumn).Thepre-andpost-operativedistance

measurementsandthecorrespondingvaluesofarelistedforboth

theleftandrighthemisphereinTable2.Furthermore,thepredicted

damageislistedinTable2andreflectsthedistancebetweenthe

Fig.6.Tractographyanddistancemeasurementresultsforthethreepatientsincludedinthestudy.ThefirstandsecondcolumnsshowthereconstructionsoftheORbefore andaftersurgery,respectively.ForeachreconstructiontheML-TPdistanceandassociatedvariabilityaredisplayed.Thethirdandfourthcolumnsshowa3Dviewofthe reconstructionoftheORintheaffectedhemisphereafterandbeforesurgery,respectively.Theresectionareaisdisplayedinredandthepredicteddamageisindicatedby color-codedredstreamlines.(Forinterpretationofthereferencestocolorinthisfigurelegend,thereaderisreferredtothewebversionofthisarticle.)

Table2

Theresultslistedforthepre-andpost-operativecomparisonofthereconstructionoftheORforbothhemispheresofthethreeepilepsysurgerycandidatesincludedinour study.DistancemeasurementsoftheanteriorextentoftheORtothetemporalpole(ML-TP)aredisplayedalongwiththevariabilityinprobabilistictractographyforthe correspondingselected.Furthermore,theresectionlengths,predictedandobserveddamages,andthemeasuredmarginsoferrorarelistedfortheaffectedhemispheres.

Patient/hemisphere ML-TPdistance  Resectionlength

(mm) Predicted damage(mm) Observed damage(mm) Marginoferror (mm)

Pre-op(mm) Post-op(mm) Pre-op

(–)/Post-op(–) Patient1 Left 30.2±0.4 28.5±1.0 0.10/0.20 – – – – Right 30.1±0.6 42.1±2.0 0.48/0.18 41.0 10.9±0.6 12.0±2.6 4.3 Patient2 Left 28.7±0.4 48.2±1.6 0.13/0.4 45.0 16.3±0.4 19.5±2.0 5.6 Right 32.0±0.7 30.5±0.5 0.13/0.2 – – – – Patient3 Left 35.3±0.7 36.2±0.9 0.10/0.18 21.0 0.0 0.0 1.6 Right 29.2±0.9 30.1±1.2 0.22/0.18 – – – –

marginoferrorisindicated,definedasthedifferencebetweenthe predicteddamageandtheobserveddamage.

Thetractographyresultsindicatethatforpatients1and2the ORisdamaged,likelyresultinginadisruptedMeyer’sloopforboth patients.Theperimetryresultsofthesepatientsindicatedavisual fielddeficit(VFD)of60degreesforpatient2,whichwassmaller thantheVFDmeasuredforpatient1 at90 degreesdespitethe largerresectionofpatient2(seeTable2).Note,thatforpatient

3,forwhomtherewasnodamagetotheOR,thereconstruction

oftheORiswellreproducibleforbothhemispheres,witha

dif-ferenceof maximally3.0mmincludingthevariability inML-TP

distance.Thedifferencebetweenthepredicteddamageandthe

observeddamagewassmallforthesepatients,indicatingan

max-imumerrorofthepredicteddamageoftheORof5.6mmorless.

Thereproducibilityofthereconstructionresultsobtained

follow-ingtheproceduresasheredescribedisfurtherconfirmedbythe

unaffectedhemispheresofeachindividualpatient,whichshowa

similaranteriorextentforbothpre-andpost-operative

reconstruc-tionsoftheOR.TheML-TPdistanceoftheORreconstructedforthe

ORofthenon-pathologichemisphereshoweddeviationsforthe

twodifferentscansofmaximally3.1mm,2.7mmand3.0mmfor

Patient1,Patient2andPatient3,respectively,includingthe

vari-abilitymeasure.TheoverallmeanML-TPdistancepre-operatively

is31.4±3.5mmforthelefthemisphereand30.4±1.4mmforthe

righthemisphere. Themean variability inprobabilistic

tractog-raphyis0.5mmand 0.7mmfortheleftandrighthemispheres,

respectively.

4. Discussion

Stabilitymetricswereintroducedforarobustestimationofthe

distancebetweenthetipoftheMeyer’sloopandtemporalpole.

Standardizedremovalofspuriousfiberswasachieved,firstlyby

quantificationofspuriousstreamlinesusingtheFBCmeasures,and

secondlybyaprocedurefortheautomaticselectionofthe

mini-malthresholdselectedontheFBCmeasures.Theresultspresented

indicatethatareliablelocalizationofthetipoftheMeyer’sloopis

possibleandthatitisfeasibletopredictthedamagetotheORas

resultofaTLRperformedtorenderpatientsseizurefree.

4.1. ProceduresforthereconstructionoftheOR

FortheestimationoftheFODfunction,CSDwasappliedon

dif-fusiondataobtainedwiththeprevalentDTIacquisitionscheme,

thusallowingforabroadclinicalapplicability.Inthecurrentstudy,

theDTIacquisitionscheme(b=1000,32directions)hasarelatively

lownumberofdirectionsofdiffusion.SincethetipoftheMeyer’s

loophasahighcurvature,itsreconstructioncouldespecially

bene-fitfromtheHARDIacquisitionscheme(Tuchetal.,2002),which

measuresa largernumber ofdirections ofdiffusionsuchas 64

or128directions.However,unlikeDTI,HARDIisnotcommonly

appliedwithinamedicalMRIdiagnosis.Instead,theDTIdatamay

beimprovedbyapplyingcontextualenhancement(Taxetal.,2014;

Portegiesetal.,2015),suchastheoneavailableintheDIPY frame-work(http://dipy.org).Additionally,inordertoimprovetheimage

qualityofthediffusionmeasurementsitmaybebeneficialtoapply

denoising.Thismay,forexample,beachievedbyarecently

pro-poseddenoisingapproachbasedonnon-localprincipalcomponent

analysis(PCA)(Manjonetal.,2015).

TheMRtrixsoftwarepackagewasemployedfortheestimation

of theFODfunction and for performingprobabilistic

tractogra-phy.Asanalternativetotherejectionsamplingmethodthat is

implementedinMRtrixforsamplingtheFODduringtracking,the

importancesamplingmethodasintroducedinFrimanetal.(2006)

couldbeused.Incontrasttothehardconstraintsusedinrejection

sampling,theimportancesamplingmethodprovidesasoft

con-straintonthespaceofpositionsandorientations,whichisinline

withthemathematicalframeworkintroducedinthis paper(see

AppendixA).

TheseedregionsoftheLGNandvisualcortexarehighly

influen-tialforthetractographyresults(Liljaetal.,2016).Itmaybepossible

toimprovethefiberorientationestimationatthewhitematterto

graymatterinterface,suchasneartheLGNandvisualcortexROIs,

byapplyingtherecentlyintroducedinformedconstrained

spher-icaldeconvolution(iCSD)(Roineetal.,2015).iCSDimprovesthe

FODbymodifyingtheresponsefunctiontoaccountfornon-white

matterpartialvolumeeffects,whichmayimprovethe

reconstruc-tionoftheOR.Inthecurrentstudy,theLGNwasidentifiedmanually

andcouldpossiblybeimprovedbyusingasemi-automaticmethod

suchaspresentedbyWinstonetal.(2011).Anotherapproach

pro-posedbyBenjaminetal.(2014)istoplacedifferentROIsaround

theLGNandwithinthesagittalstratum,orbyseedingfromthe

opticchiasm(Kammenetal.,2016).Arecentstudysuggestedusing

seedingaroundtheMeyer’sloopwithana-priorifiberorientation

(Chamberlandetal.,2017).

4.2. Applicationofthestabilitymetrics

The FBC measures are used for the quantification of

spuri-ousstreamlines.TheseFBCmeasuresarebasedontheestimation

ofstreamline densityin thespaceofpositionsand orientations

R3_×S2_{.Anadvantage oftheFBC methodisthat it isgenerally}

applicable,regardlessofthetypeofdiffusionmodelandthe

track-ingalgorithmbeingused,sinceitdependsonlyontheoutcome

oftractography.ApossiblelimitationoftheFBCmeasuresarethe

number ofstreamlinesthat canbeprocessed,since fordensely

populatedregionsofstreamlinesthemethodiscomputationally

expensive.However,throughtheuseofseveraloptimizationsteps

suchaspre-computedlookuptablesfortheBrownianmotion

ker-nel,multi-threadedprocessing,subsamplingofstreamlines,and

theexclusionoffar-awaystreamlinepoints,thecomputationtimes

maintainmanageable.DetailsareavailableinAppendixB.

Inordertoremovespuriousfiberswhilepreventingan

under-estimationofthefullextentoftheMeyer’sloop,aprocedurefor

estimatingselectedwasintroducedbasedonthetest–retest

evalu-ationofthevariabilityinML-TPdistance.Usingthismethodology,a

robustmeasurementoftheML-TPdistancewasachievedintheleft

andrighthemispheresofeighthealthyvolunteers.Thevariabilityin

thereconstructionresultsoftheORstemsmostlyfromdata

acquisi-tion(e.g.SNR,partialvolumeeffects,andpatientmotion)(Wakana

etal.,2007).Therefore,selectedmayvarybetweenpre-and

post-operativescansinthenon-affectedhemisphere(seeTable2).The

meanML-TPdistancesforbothbrainhemispheres,measuredtobe

30.0±4.5mmforthehealthyvolunteergroupand30.9±2.4mm

for thepatientgroup(pre-operatively), arewithintherange of

theML-TPdistancereportedonbyEbelingetal.(1988)and

out-comesfromotherORreconstructionmethodologies.Forexample,

ConTrack(Sherbondyetal.,2008)showing28±3.0mm,

Stream-lines Tracer technique (STT) showing 37±2.5mm (Yamamoto

et al., 2005)and 44±4.9mm (Nilsson et al., 2007), Probability

IndexofConnectivity(PICo)showing36.2±0.7mm(Dayanetal.,

2015), tractographyonHumanConnectiveProject (HCP)

multi-shell data showing 30.7±4.0mm (Kammen et al., 2016), and

MAGNET showing 36.0±3.8mm (Chamberland et al., 2017). It

appeared,furthermore,thatthemeanML-TPforboththehealthy

volunteersandthepatientswaslargerinthelefthemisphere

com-paredtotherighthemisphere,whichisnotconsistentwitharecent

studybyJamesetal.(2015)thatindicatedasignificantlyhigher

ML-TPintherighthemisphere.

Apossiblelimitationoftheparameterestimationprocedureis

FBCmeasures,whichcanbeusedforanytractogram,the

parame-terestimationproceduremaynotbegenerallyapplicableforother

fiberbundlessinceadistancemeasurementbetweenwell-defined

landmarksis required.However,a possibleapproachfor

gener-alizedparameter selection is to fit thestreamline bundle on a

manifoldsuchasusedbyBundleMAP(Khatamietal.,2016)and

optimizeselectedbyminimizingthespreadonthemanifold.

4.3. Towardsdamagepredictionforepilepsysurgery

ThemethodologyfortheestimationoftheML-TPdistanceis

appliedforthesurgicalcandidates,firstlytoassessthevalidityof

thedistancemeasurements,andsecondlytoindicateitsadditional

valueforresectiveepilepsysurgery.Anindicationofthevalidity

ofdistancemeasurementswasgivenbythemarginoferror,which

wasthelargestforpatient2amountingto5.6mm.Themarginof

errorobservedforthethreepatientscanbelowered,e.g.by

correct-ingforbrainshiftsthatoccurduetoresectionandCSFloss(Warfield

etal.,2005)andbycorrectingfordistortionspresentinMR

echo-planarimaging(JezzardandBalaban,1995;Hollandetal.,2010).

ThemeasurementoftheML-TPdistancemaybefurther

compli-cateddue toashiftedlocationofthetemporalpole,orevenits

completeabsence.However,thereproducibilityof thepre- and

post-operativereconstructionsoftheORinthenon-pathological

hemisphereindicatesthattheeffectsofbrainshiftand imaging

distortionsmaybelimited.SmalldeviationsintheML-TPdistance

wereseen(seeTable2), which suggestsa goodreproducibility,

albeitforalimitednumberofpatients.

Inthestandardizedestimationprocedureofselectedthe

max-imalvariabilitywassetat2mm,bothfortheORreconstructions

ofthehealthyvolunteersandthepatients,whichisbasedonthe

maximalsurgicalaccuracythatcanbeachievedduringstandardor

tailoredanteriortemporallobectomybeforetheleakageof

cere-brospinalfluid(CSF). Asurgicalaccuracy below2mmhasbeen

reported(Tibbals,2010)ifastereotacticframeisusedorrobotic

assistanceisinvolved.AftertheleakageofCSFhowever,cortical

displacementupto24mmmaybeseen(Hastreiteretal.,2004),

whileothersourcesofinaccuracyarelikelypresentsuchas

echo-planarimagingdistortion,partialvolumeeffects,andimagenoise.

However,despitetheseinaccuraciesthepre-andpost-operative

comparisonoftheORreconstructionsindicatesthattheprocedures

developedinthisstudyareavalidtooltoassesstherobustnessof

thedistancemeasurements.

ItappearedthattherobustestimationoftheML-TPdistance

enabledtopredictthedamageoftheORaftersurgery,whichwas

concordantwiththeactualdamageforthethreepatientsstudied.

Basedonthedamagepredictionthemarginoferrorwasestimated,

givinganindicationoftheoverallerrorindistancemeasurement.

The perimetry results of two of thepatients studied indicated

damageofeithertheleftorrightvisualfield,correspondingtoa

dis-ruptionoftheMeyer’sloop.ArelativelysmallVFDwasindicatedfor

patient2despitethelargetemporalloberesection.Thisresultmay

beindicativeofthelargeinter-patientvariabilityinORanatomy

andfunction,butmayalsobetheresultofthenon-standardized

proceduresforvisualfieldtestingin-betweenhospitals.Itis

rec-ommendedtoevaluatethedevelopedmethodologyfurtherin a

clinicaltrialincludingasizablegroupofpatientswhoare

candi-dateforaTLRinordertobeabletoassesswhattherelationis

betweenaVFDandthedamagetotheORafteraTLR.

5. Conclusion

Itwasshownforagroupofhealthyvolunteersincludedinthis

studythatstandardizedremovalofspuriousstreamlinesprovides

areliableestimationofthedistancefromthetipoftheMeyer’sloop

tothetemporalpolethatisstableunderthestochasticrealizations

ofprobabilistictractography.Pre-andpost-operativecomparisons

ofthereconstructedORindicated,furthermore,(1)thevalidityofa

robustML-TPdistancemeasurementtopredictthedamagetothe

ORasresultofresectivesurgery,and(2)thehighreproducibility

ofthereconstructionsofthenon-pathologicalhemisphere.In

con-clusion,thedevelopedmethodologybasedondiffusion-weighted

MRItractographyisasteptowardsapplyingopticradiation

trac-tographyforpre-operativeplanningofresectivesurgeryandfor

providinginsightinthepossibleadverseeventsrelatedtothistype

ofsurgery.

Conflictsofinterest

Theauthorsdeclarethattheresearchwasconductedinabsence

ofanycommercialorfinancialrelationshipsthatcouldbeconstrued

asapossibleconflictofinterest.

Acknowledgments

WewouldliketothankBartterHaar Romenyfor

contribut-ingtotheresearchcollaborationbetweentheAcademicCenter

for Epileptology, Kempenhaeghe& MUMC+ and theEindhoven

Universityof Technology.Furthermore,we thankJan Verwoerd

(Philips HealthcareBenelux)for contributions totheMRI

scan-ningprotocols,JorgPortegiesforthehelpindevelopingtheFBC

measures,andRemcoBertingforassistanceduringMRIscanning.

Thepatientswere evaluatedand discussed presurgically inthe

localpresurgicalworkgroupofKempenhaegheandMaastrichtUMC+

(AWEC) and in the national Dutchepilepsy surgery workgroup

(LWEC).Theresearchleadingtotheseresultshasreceivedfunding

fromtheEuropeanResearchCouncilundertheEuropean

Commu-nity’sSeventhFrameworkProgramme(FP7/2007-2014)/ERCgrant

agreementno.335555.

AppendixA. Mathematicalbackground

The fiber-to-bundle coherence(FBC)measures are based on

kerneldensityestimationinthenon-flat5Dposition-orientation

domain.Firstof all,each equidistantlysampledstreamline i=



k i | k=1,...,Ni



with k i =(yki,nki)∈R 3_×S2 _are

repre-sentedbydeltadistributionsı_(yj

i,n j i)

.Here,Nidenotesthenumberof

streamlinepointsofstreamlineiandyk_iandnk_idenotetheposition

andtangentorientationofthestreamlinepointk

i,respectively.The

fullliftedoutputofthetractographyisgivenby

F(y,n)= 1 Ntot 2



=1 Ntot



i=1 N_i



j=1 ı_(yj i,(−1)n j i) (y,n), (A.1) where=∪Ntot

i=1{i}denotesthestreamlinebundleandNtot

indi-catesthenumberofstreamlinesinthebundle.Thesummation

overisusedtoincludeantipodalsymmetry(whereweidentify

nj_i∼−nj_i)ofeachtangentorientationnj_i.

ThekerneldensityestimatorisdefinedbyFokker-Planck

dif-fusion equations, which describe Brownian motion on R3_S2

followingevolutionprocessisusedwhereF=F_servesastheinitial condition,



∂

tWF(y,n,t)=(Dspat(n·

∇

y)2+Dang_S2)W_F(y,n,t),

WF(y,n,0)=F(y,n).

(A.2)

Here,t≥0istheevolutiontime,Dspat>0isthecoefficientfor

spa-tialsmoothingstrictlyinthedirectionofn.Dangisthecoefficient

forangularsmoothing(S2istheLaplace-Beltramioperatoronthe

sphereS2_{).Inthisevolutionprocess,W}

F(y,n,t)representsthe

tran-sitiondensityofamovingparticlewithpositionyandorientation

natthetimet≥0,giventhatitstartedwithinitialdistributionF(y,

n)att=0.Then,theLocalFBC(LFBC)istheresultofevaluatingthe

Brownianmotionkernelpt(seeFig.3,left)alongeachelementof

theliftedstreamline

LFBC(y,n,)=(pt∗R3_×S2F) ( · ) =



R3



S2 pt(RTn(y−y),RTnn)F(y,n)d␴(n)dy. (A.3)

Here,pt(y,n)denotestheGreen’sfunctionoftheevolution

pro-cessinEq.(A.2),whichequalstheprobabilitydensityoffindinga

randomorientedparticleatpositiony,withorientationn,attime

t∈R+_{giventhatitstartedatposition0andorientatione}

z ∈S2at

time0.Likewise,pt(RTn(y−y),RTnn)istheprobabilitydensityof

findingarandomorientedparticleatpositionyandorientationn

giventhatitstartedatpositionyandorientationnatt=0.Here,

distheusualmeasureonthesphereS2_{.Asaresult,by}

super-positionin(A.3),LFBC(y,n,)denotestheprobabilitydensityof

findingarandomorientedparticleatyandpointingatorientation

nattimet>0giventhatitstartedatsomepointofthebundle

att=0.Forexactformulasforthekernelpt(y,n),andtheGaussian

approximationsthatweusedforourcomputations,seePortegies

etal.(2015).

Awholestreamline measure,therelative FBC(RFBC),is

cal-culatedbytheminimumofthemovingaverageLFBCalongthe

streamlinei

RFBC˛(i,)=

AFBC˛(i,)

AFBC() , (A.4)

whereAFBC˛_(

i,)indicatestheminimalaverageLFBCoverasmall

segmentofthestreamlinewithlength˛,givenby

AFBC˛(i,)=mina∈[0,l_i−a]_˛1



a+˛

a

LFBC(i(s),)ds, (A.5)

andAFBC()istheaverageFBCoftheentirestreamlinebundle

AFBC()= 1 N N



i=1 FBC(i,), (A.6) where FBC(i,)= 1 li



l_i 0 LFBC(i(s),)ds. (A.7)

Here,liisthetotallengthofthespatiallyprojectedcurvexi(·)of

fiberi=(xi,ni).FurtherdetailsofthecomputationofFBCcanbe

foundinPortegiesetal.(2015).

AppendixB. Computationaloptimization

TheFBCmeasuresareimplementedinsideDIPY(Garyfallidiset

al.,2014)usingthehigh-speedCython(C++inPython)language.

Thekerneldensityestimationisexecutedwithmultithreadingvia

theOpenMPlibrary,whichespeciallyforclustercomputing

pro-videsasignificantspeedup.Tofurtheracceleratethekerneldensity

estimation, lookup-tablesarecomputed containingrotated

ver-sionsof thekernelpt rotatedovera discretesetoforientations

(Rodriguesetal.,2010).Therotatedversionsareequallydistributed

over aspheretoensurerotationallyinvariantprocessing. Tobe

abletousethelookuptableduringkerneldensityestimation,each

(continuous)streamlinetangentorientationismatchedwiththe

closest(discrete)orientationonthesphere.Forefficient

imple-mentationoforientationmatching,aKD-treeisused,whichisa

multi-dimensional(K=3)binaryspacepartitioning,tominimize

thenumberofangulardistancecomputations.

References

Benjamin,C.F.,Singh,J.M.,Prabhu,S.P.,Warfield,S.K.,2014.Optimizationof tractographyoftheopticradiations.Hum.BrainMapp.35(2),683–697, http://dx.doi.org/10.1002/hbm.22204http://www.ncbi.nlm.nih.gov/pmc/ articles/PMC4083652/.

Borius,P.-Y.,Roux,F.-E.,Valton,L.,Sol,J.-C.,Lotterie,J.-A.,Berry,I.,2014.CanDTI fibertrackingoftheopticradiationspredictvisualdeficitaftersurgery?Clin. Neurol.Neurosurg.122,87–91,

http://dx.doi.org/10.1016/j.clineuro.2014.04.017.

Chamberland,M.,Scherrer,B.,Prabhu,S.P.,Madsen,J.,Fortin,D.,Whittingstall,K., Descoteaux,M.,Warfield,S.K.,2017.ActivedelineationofMeyer’sloopusing orientedpriorsthroughMAGNEtictractography(MAGNET).Hum.BrainMapp., http://dx.doi.org/10.1002/hbm.23399.

Chen,X.,Weigel,D.,Ganslandt,O.,Buchfelder,M.,Nimsky,C.,2009.Predictionof visualfielddeficitsbydiffusiontensorimagingintemporallobeepilepsy surgery.Neuroimage45(2),286–297,

http://dx.doi.org/10.1016/j.neuroimage.2008.11.038http://www.ncbi.nlm.nih. gov/pubmed/19135156.

Citti,G.,Sarti,A.,2006.Acorticalbasedmodelofperceptualcompletioninthe roto-translationspace.J.Math.ImagingVis.24(3),307–326,

http://dx.doi.org/10.1007/s10851-005-3630-2.

Dayan,M.,Munoz,M.,Jentschke,S.,Chadwick,M.J.,Cooper,J.M.,Riney,K., Vargha-Khadem,F.,Clark,C.A.,2015.Opticradiationstructureandanatomyin thenormallydevelopingbraindeterminedusingdiffusionMRIand

tractography.BrainStruct.Funct.220(1),291–306, http://dx.doi.org/10.1007/s00429-013-0655-y.

Descoteaux,M.,Deriche,R.,Knösche,T.R.,Anwander,A.,2009.Deterministicand probabilistictractographybasedoncomplexfibreorientationdistributions. IEEETrans.Med.Imaging28(2),269–286,

http://dx.doi.org/10.1109/TMI.2008.2004424http://www.ncbi.nlm.nih.gov/ pubmed/19188114.

Duits,R.,Franken,E.,2011.Left-invariantdiffusionsonthespaceofpositionsand orientationsandtheirapplicationtocrossing-preservingsmoothingofHARDI images.Int.J.Comput.Vis.92(3),231–264.

Duits,R.,Felsberg,M.,Granlund,G.,Romeny,B.,2007.Imageanalysisand reconstructionusingawavelettransformconstructedfromareducible representationoftheEuclideanmotiongroup.Int.J.Comput.Vis.72(1), 79–102,http://dx.doi.org/10.1007/s11263-006-8894-5.

Duits,R.,DelaHaije,T.,Creusen,E.,2013.Morphologicalandlinearscalespacesfor fiberenhancementinDW-MRI.J.Math.ImagingVis.46(3),326–368. Ebeling,U.,Reulen,H.,Universitfitsklinik,N.,1988.Neurosurgicaltopographyof

theopticradiationinthetemporallobe.ActaNeurochir.92(1–4),29–36. Friman,O.,Farnebäck,G.,Westin,C.-F.,2006.Abayesianapproachforstochastic

whitemattertractography.IEEETrans.Med.Imaging25,965–978.

Fujita,N.,Tanaka,H.,Takanashi,M.,Hirabuki,N.,Abe,K.,Yoshimura,H.,Nakamura, H.,2001.Lateralgeniculatenucleus:anatomicandfunctionalidentificationby useofMRimaging.Am.J.Neuroradiol.22(9),1719–1726

http://www.ncbi.nlm.nih.gov/pubmed/11673167.

Garyfallidis,E.,Brett,M.,Amirbekian,B.,Rokem,A.,vanderWalt,S.,Descoteaux, M.,Nimmo-Smith,I.,2014.DipyContributorsDipyalibraryfortheanalysisof diffusionMRIdata.Front.Neuroinform.8,8,

http://dx.doi.org/10.3389/fninf.2014.00008.

Gorgolewski,K.,Burns,C.D.,Madison,C.,Clark,D.,Halchenko,Y.O.,Waskom,M.L., Ghosh,S.S.,2011.Nipype:aflexible,lightweightandextensibleneuroimaging dataprocessingframeworkinPython.Front.Neuroinform.5,

http://dx.doi.org/10.3389/fninf.2011.00013.

Hastreiter,P.,Rezk-Salama,C.,Soza,G.,Bauer,M.,Greiner,G.,Fahlbusch,R., Ganslandt,O.,Nimsky,C.,2004.Strategiesforbrainshiftevaluation.Med. ImageAnal.8(4),447–464,

http://dx.doi.org/10.1016/j.media.2004.02.001http://www.ncbi.nlm. nih.gov/pubmed/15567708.

Holland,D.,Kuperman,J.M.,Dale,A.M.,2010.Efficientcorrectionof inhomogeneousstaticmagneticfield-induceddistortioninEchoPlanar Imaging.Neuroimage50(1),175–183,

http://dx.doi.org/10.1016/j.neuroimage.2009.11.044, http://www.pubmedcentral.nih.gov/articlerender.fcgi? art%id=2819607&tool=pmcentrez&rendertype=abstract.

James,J.S.,Radhakrishnan,A.,Thomas,B.,Madhusoodanan,M.,Kesavadas,C., Abraham,M.,Menon,R.,Rathore,C.,Vilanilam,G.,2015.Diffusiontensor imagingtractographyofMeyer’sloopinplanningresectivesurgeryfor

drug-resistanttemporallobeepilepsy.EpilepsyRes.110,95–104, http://dx.doi.org/10.1016/j.eplepsyres.2014.11.020.

Jenkinson,M.,Beckmann,C.,Behrens,T.,Woolrich,M.,Smith,S.,2012.FSL. Neuroimage62(2),782–790.

Jeurissen,B.,Leemans,A.,Jones,D.K.,Tournier,J.-D.,Sijbers,J.,2011.Probabilistic fibertrackingusingtheresidualbootstrapwithconstrainedspherical deconvolution.Hum.BrainMapp.32(3),461–479,

http://dx.doi.org/10.1002/hbm.21032http://www.ncbi.nlm.nih.gov/ pubmed/21319270.

Jezzard,P.,Balaban,R.,1995.Correctionforgeometricdistortioninechoplanar imagesfromB0fieldvariations.Magn.Reson.Med.34(1),65–73. Kammen,A.,Law,M.,Tjan,B.S.,Toga,A.W.,Shi,Y.,2016.Automatedretinofugal

visualpathwayreconstructionwithmulti-shellHARDIandFOD-based analysis.Neuroimage125,767–779,

http://dx.doi.org/10.1016/j.neuroimage.2015.11.005.

Khatami,M.,Schmidt-Wilcke,T.,Sundgren,P.C.,Abbasloo,A.,Schölkopf,B., Schultz,T.,2016.BundleMAP:anatomicallylocalizedclassification,regression, andhypothesistestingindiffusionMRI.PatternRecognit.,

http://dx.doi.org/10.1016/j.patcog.2016.09.020http://www.sciencedirect.com/ science/article/pii/S00313%20316302849.

Kier,E.L.,Staib,L.H.,Davis,L.M.,Bronen,R.A.,2004.MRimagingofthetemporal stem:anatomicdissectiontractographyoftheuncinatefasciculusinferior occipitofrontalfasciculusandMeyer’sloopoftheopticradiation.Am.J. Neuroradiol.25(5),677–691.

Leemans,A.,Jones,D.K.,2009.TheB-matrixmustberotatedwhencorrectingfor subjectmotioninDTIdata.Magn.Reson.Med.61(6),1336–1349, http://dx.doi.org/10.1002/mrm.21890.

Lilja,Y.,Nilsson,D.T.,2015.Strengthsandlimitationsoftractographymethodsto identifytheopticradiationforepilepsysurgery.Quant.ImagingMed.Surg.5 (2),288–299,http://dx.doi.org/10.3978/j.issn.2223-4292.2015.01.08. Lilja,Y.,Gustafsson,O.,Ljungberg,M.,Nilsson,D.,Starck,G.,2016.Impactof

region-of-interestmethodonquantitativeanalysisofDTIdataintheoptic tracts.BMCMed.Imaging16(1),42,

http://dx.doi.org/10.1186/s12880-016-0145-9.

Lim,J.C.,Phal,P.M.,Desmond,P.M.,Nichols,A.D.,Kokkinos,C.,Danesh-Meyer,H.V., Kaye,A.H.,Moffat,B.A.,2015.ProbabilisticMRItractographyoftheoptic radiationusingconstrainedsphericaldeconvolution:afeasibilitystudy.PLOS ONE10(3),e0118948,

http://dx.doi.org/10.1371/journal.pone.0118948http://www.ncbi.nlm. nih.gov/pmc/articles/PMC4351098/.

Manjon,J.V.,Coupe,P.,Buades,A.,2015.MRInoiseestimationanddenoisingusing non-localPCA.Med.ImageAnal.22(1),35–47,

http://dx.doi.org/10.1016/j.media.2015.01.004.

Martínez-Heras,E.,Varriano,F.,Prˇckovska,V.,Laredo,C.,Andorrà,M., Martínez-Lapiscina,E.H.,Calvo,A.,Lampert,E.,Villoslada,P.,Saiz,A., Prats-Galino,A.,Llufriu,S.,2015.Improvedframeworkfortractography reconstructionoftheopticradiation.PLOSONE10(9),e0137064, http://dx.doi.org/10.1371/journal.pone.0137064.

Mori,S.,2007.IntroductiontoDiffusionTensorImaging.ElsevierScience. Mumford,D.,1994.AlgebraicGeometryandItsApplications.Springer-Verlag,pp.

491–506(Ch.Elasticaandcomputervision).

Nilsson,D.,Starck,G.,Ljungberg,M.,Ribbelin,S.,Jönsson,L.,Malmgren,K., Rydenhag,B.,2007.Intersubjectvariabilityintheanteriorextentoftheoptic radiationassessedbytractography.EpilepsyRes.77(1),11–16,

http://dx.doi.org/10.1016/j.eplepsyres.2007.07.012http://www.ncbi.nlm. nih.gov/pubmed/17851037.

Øksendal,B.,2003.StochasticDifferentialEquations,6thed.Springer-Verlag. Parker,G.J.M.,Haroon,H.A.,Wheeler-Kingshott,C.A.M.,2003.Aframeworkfora

streamline-basedprobabilisticindexofconnectivity(PICo)usingastructural interpretationofMRIdiffusionmeasurements.J.Magn.Reson.Imaging18(2), 242–254,http://dx.doi.org/10.1002/jmri.10350.

Peltier,J.,Travers,N.,Destrieux,C.,Velut,S.,2006.Opticradiations:amicrosurgical anatomicalstudy.J.Neurosurg.105(2),294–300,

http://dx.doi.org/10.3171/jns.2006.105.2.294http://www.ncbi.nlm.nih.gov/ pubmed/17219837.

Pestilli,F.,Yeatman,J.D.,Rokem,A.,Kay,K.N.,Wandell,B.A.,2014.Evaluationand statisticalinferenceforhumanconnectomes.Nat.Methods11(10), 1058–1063,http://dx.doi.org/10.1038/nmeth.3098.

Portegies,J.M.,Fick,R.H.J.,Sanguinetti,G.R.,Meesters,S.P.L.,Girard,G.,Duits,R., 2015.ImprovingfiberalignmentinHARDIbycombiningcontextualPDEflow withconstrainedsphericaldeconvolution.PLOSONE10(10),1–33, http://dx.doi.org/10.1371/journal.pone.0138122.

Rockafellar,R.T.,Wets,R.J.-B.,2005.VariationalAnalysis.Springer-Verlag,pp.117. Rodrigues,P.,Duits,R.,Romeny,B.M.t.H.,Vilanova,A.,2010.Accelerateddiffusion

operatorsforenhancingdw-mri.In:Bartz,D.,Botha,C.,Hornegger,J., Machiraju,R.,Wiebel,A.,Preim,B.(Eds.),EurographicsWorkshoponVisual ComputingforBiologyandMedicine.TheEurographicsAssociation, http://dx.doi.org/10.2312/VCBM/VCBM10/049-056.

Roine,T.,Jeurissen,B.,Perrone,D.,Aelterman,J.,Philips,W.,Leemans,A.,Sijbers,J., 2015.Informedconstrainedsphericaldeconvolution(iCSD).Med.ImageAnal., http://dx.doi.org/10.1016/j.media.2015.01.001.

Rubino,P.A.,Rhoton,A.L.,Tong,X.,deOliveira,E.,2005.Three-dimensional relationshipsoftheopticradiation.Neurosurgery57(Suppl.4),219–227, http://dx.doi.org/10.1227/01.NEU.0000176415.83417.16http://content. wkhealth.com/linkback/openurl?sid=WKPTLP:landingpage&an=00006123-200510004-00001.

Sherbondy,A.J.,Dougherty,R.F.,Wandell,B.A.,2008.Identifyingthehumanoptic radiationusingdiffusionimagingandfibertractography.J.Vis.8(10),12, http://dx.doi.org/10.1167/8.10.12.Introduction,1–11.

Sincoff,E.H.,Tan,Y.,Abdulrauf,S.I.,2004.Whitematterfiberdissectionoftheoptic radiationsofthetemporallobeandimplicationsforsurgicalapproachestothe temporalhorn.J.Neurosurg.101(5),739–746,

http://dx.doi.org/10.3171/jns.2004.101.5.0739http://www.ncbi.nlm.nih. gov/pubmed/15540910.

Taoka,T.,Sakamoto,M.,Iwasaki,S.,Nakagawa,H.,Fukusumi,A.,Hirohashi,S., Taoka,K.,Kichikawa,K.,Hoshida,T.,Sakaki,T.,2005.Diffusiontensorimaging incaseswithvisualfielddefectafteranteriortemporallobectomy.Am.J. Neuroradiol.26(4),797–803

http://www.ncbi.nlm.nih.gov/pubmed/15814923.

Tax,C.M.,Duits,R.,Vilanova,A.,terHaarRomeny,B.M.,Hofman,P.,Wagner,L., Leemans,A.,Ossenblok,P.,2014.Evaluatingcontextualprocessingindiffusion MRI:Applicationtoopticradiationreconstructionforepilepsysurgery.PLOS ONE9(7).

Tibbals,H.F.,2010.MedicalNanotechnologyandNanomedicine.CRCPress. Tournier,J.-D.,Calamante,F.,Connelly,A.,2007.Robustdeterminationofthefibre

orientationdistributionindiffusionMRI:non-negativityconstrained super-resolvedsphericaldeconvolution.Neuroimage35(4),1459–1472, http://dx.doi.org/10.1016/j.neuroimage.2007.02.016http://www.ncbi.nlm. nih.gov/pubmed/17379540.

Tournier,J.-D.,Calamante,F.,Connelly,A.,2012.MRtrix:diffusiontractographyin crossingfiberregions.Int.J.ImagingSyst.Technol.22(1),53–66,

http://dx.doi.org/10.1002/ima.22005.

Tuch,D.S.,Reese,T.G.,Wiegell,M.R.,Makris,N.,Belliveau,J.W.,Wedeen,V.J.,2002. Highangularresolutiondiffusionimagingrevealsintravoxelwhitematter fiberheterogeneity.Magn.Reson.Med.48(4),577–582,

http://dx.doi.org/10.1002/mrm.10268.

Tuch,D.S.,2004.Q-ballimaging.Magn.Reson.Med.52(6),1358–1372, http://dx.doi.org/10.1002/mrm.20279.

Wakana,S.,Caprihan,A.,Panzenboeck,M.M.,Fallon,J.H.,Perry,M.,Gollub,R.L., Hua,K.,Zhang,J.,Jiang,H.,Dubey,P.,Blitz,A.,vanZijl,P.,Mori,S.,2007. Reproducibilityofquantitativetractographymethodsappliedtocerebral whitematter.Neuroimage36(3),630–644,

http://dx.doi.org/10.1016/j.neuroimage.2007.02.049http://www.ncbi.nlm.nih. gov/pmc/articles/PMC2350213/.

Warfield,S.K.,Haker,S.J.,Talos,I.-F.,Kemper,C.A.,Weisenfeld,N.,Mewes,A.U.J., Goldberg-Zimring,D.,Zou,K.H.,Westin,C.-F.,Wells,W.M.,Tempany,C.M.C., Golby,A.,Black,P.M.,Jolesz,F.A.,Kikinis,R.,2005.Capturingintraoperative deformations:researchexperienceatbrighamandwomen’shospital.Med. ImageAnal.9(2),145–162,http://dx.doi.org/10.1016/j.media.2004.11.005. Wilkins,B.,Lee,N.,Gajawelli,N.,Law,M.,Lepore,N.,2015.Fiberestimationand

tractographyindiffusionMRI:developmentofsimulatedbrainimagesand comparisonofmulti-fiberanalysismethodsatclinicalb-values.Neuroimage 109,341–356,http://dx.doi.org/10.1016/j.neuroimage.2014.12.060. Winston,G.P.,Mancini,L.,Stretton,J.,Ashmore,J.,Symms,M.R.,Duncan,J.S.,

Yousry,T.A.,2011.Diffusiontensorimagingtractographyoftheopticradiation forepilepsysurgicalplanning:acomparisonoftwomethods.EpilepsyRes.97 (1–2),124–132,

http://dx.doi.org/10.1016/j.eplepsyres.2011.07.019http://www.ncbi. nlm.nih.gov/pmc/articles/PMC3223565/.

Winston,G.P.,Daga,P.,Stretton,J.,Modat,M.,Symms,M.R.,McEvoy,A.W., Ourselin,S.,Duncan,J.S.,2012.Opticradiationtractographyandvisionin anteriortemporalloberesection.Ann.Neurol.71(3),334–341, http://dx.doi.org/10.1002/ana.22619.

Yamamoto,T.,Yamada,K.,Nishimura,T.,Kinoshita,S.,2005.Tractographyto depictthreelayersofvisualfieldtrajectoriestothecalcarinegyri.Am.J. Ophthalmol.140(5),781–785,

http://dx.doi.org/10.1016/j.ajo.2005.05.018http://www.ncbi.nlm.nih.gov/ pubmed/16310456.

Yeatman,J.D.,Dougherty,R.F.,Myall,N.J.,Wandell,B.A.,Feldman,H.M.,2012.Tract profilesofwhitematterproperties:automatingfiber-tractquantification. PLOSONE7(11),1–15,http://dx.doi.org/10.1371/journal.pone.0049790. Yogarajah,M.,Focke,N.K.,Bonelli,S.,Cercignani,M.,Acheson,J.,Parker,G.J.M.,

Alexander,D.C.,McEvoy,a.W.,Symms,M.R.,Koepp,M.J.,Duncan,J.S.,2009. DefiningMeyer’sloop-temporalloberesections,visualfielddeficitsand diffusiontensortractography.Brain132(6),1656–1668,

http://dx.doi.org/10.1093/brain/awp114http://www.pubmedcentral.nih.gov/ articlerender.fcgi?artid=2685925&tool=pmcentrez&rendertype=abstract. Yushkevich,P.A.,Piven,J.,CodyHazlett,H.,GimpelSmith,R.,Ho,S.,Gee,J.C.,Gerig,

G.,2006.User-guided3Dactivecontoursegmentationofanatomical structures:significantlyimprovedefficiencyandreliability.Neuroimage31 (3),1116–1128.