Uncertainty quantification of computational coronary stenosis assessment and model based mitigation of image resolution limitations

ContentslistsavailableatScienceDirect

Journal of Computational Science

j o u r n al ho me p a g e :w w w . e l s e v i e r . c o m / l o c a t e / j o c s

Uncertainty quantification of computational coronary stenosis assessment and model based mitigation of image resolution limitations

Jacob Sturdy

^∗

, Johannes Kløve Kjernlie, Hallvard Moian Nydal, Vinzenz G. Eck, Leif R. Hellevik

DepartmentofStructuralEngineering,NorwegianUniversityofScienceandTechnology(NTNU),RichardBirkelandsvei1a,7491Trondheim,Norway

a r t i c l e i n f o

Articlehistory:

Received12November2018 Accepted19January2019 Availableonline2February2019

Keywords:

Coronarystenosis Uncertaintyquantification Sensitivityanalysis Fractionalflowreserve Hemodynamics Modeling

a b s t r a c t

Coronaryarterydiseaseisoneoftheleadingcausesofdeathglobally.Thehallmarkofthisdiseaseisthe occurrenceofstenosedcoronaryarterieswhichreducebloodflowtothemyocardium.Severelystenosed arteriescanbetreatedifdetected,butthediagnosticproceduretoassessfractionalflowreserve(FFR), aquantitativemeasureofstenosisseverity,isinvasive,burdensometothepatient,andcostly.Recent computationalapproachesestimatetheseverityofstenosesfromsimulationsofcoronarybloodflow basedonCTimagery.Thesemethodsallowfordiagnosistobemadenoninvasivelyandusingfewerhos- pitalresources;however,thepredictionsdependonuncertaininputdataandmodelparametersdue totechnicallimitationsandpatientvariability.Toassesstheconsequencesofboundaryconditionand inputuncertaintyonpredictionsofFFR,wedevelopedamodelofcoronarybloodflow.Weperformed uncertaintyquantificationandsensitivityanalysisofthepredictionsbasedonuncertaintiesinboundary conditions,parameters,andgeometricmeasurements.Ourresultsidentifiedthreeinfluentialsourcesof uncertainty:geometricdata,cardiacoutput,andcoronaryresistanceduringhyperemia.Further,uncer- taintyaboutthegeometryofthestenosedcoronarybranchinfluencesestimatesmuchmorethanother geometricaldata.Limitationsofmedicalimagingcontributeuncertaintytopredictionsasvesselsbelowa certainthresholdremainunobserved.Weassessedtheeffectsofunobservedvesselsbycomparingpredic- tionsbasedonbothhighandlowresolutiondata.Moreover,weintroducedanovelmethodthatestimates flowdistributionwhileaccountingforunobservedvessels.ThismethodimprovedFFRpredictionsinthe casesconsideredby50%onaverage.

©2019ElsevierB.V.Allrightsreserved.

1. Introduction

Cardiovasculardiseases(CVDs)accountformorethan17mil- liondeathseachyearglobally,anumberexpectedtogrowto23.6 millionby2030[1].TheWorldHealthOrganizationestimatesthat coronaryarterydisease(CAD)caused7.4milliondeathsworldwide (thelargestsinglecauseofdeath)in2012[2].Inadditiontothe socialandpersonalcostsofCVDs,theepidemic’seconomiccosts aregrowingrapidly.IntheUS,CVDsareprojectedtotriplefrom 273billionUSDto818billionUSDfrom2010to2030[3].Coro- naryarterydisease(CAD)leadingtomyocardialischemiacomposes thelargestportionofthese.CADdescribesabuildupofplaquesin arterialwallsreducingthelumenareaofbloodvessels,effectively

∗Correspondingauthor.

E-mailaddress:jacob.t.sturdy@ntnu.no(J.Sturdy).

resultinginanarrowingofthebloodvesselreferredtoasastenosis.

Thestenosisincreasesworkrequiredtoperfusethemyocardialtis- sueandincreasestheriskofaclotblockingthearteryandleading toaninfarction.GiventheprevalenceandrisksofCAD,significant effortshavebeenmadetoimprovethediagnosisandtreatment ofthisconditionbyidentifyingwhichstenosesaresignificantand intervening,typicallybyinsertionofvascularstent,toensuresuffi- cientbloodflowthroughtothemyocardium.Diagnosisoftenoccurs onlyafteranindividualhasnoticedthesymptomsofCAD-chest painorinmoreseriouscasesheartattackorarrhythmia.Methods toevaluatethesignificanceofastenosisrangefromnon-invasive imagingbasedapproachestoinvasivemeasurementofbloodpres- sures.Ifastenosisisdeemedsignificant,thenthepatienttypically undergoescoronaryangioplasty(percutaneouscoronaryinterven- tion(PCI))toopenthearteryand,ifdeemedbeneficial,toinserta stent.

https://doi.org/10.1016/j.jocs.2019.01.004 1877-7503/©2019ElsevierB.V.Allrightsreserved.

The most widely used method of identifying and assess- ing stenosis severity considers only the visual appearance of the stenosis imaged through coronary angiography via cardiac catheterization.In2010morethan1millioncardiaccatheteriza- tionswereperformedintheUS[4].Amoreobjectivebutinvasive methodisbasedonfractionalflowreserve(FFR),whichisdefined astheratioofmaximalbloodflowthroughthestenosedarteryto themaximalblood flowthroughthesamearteryifthestenosis werenotpresent.FFRmaybeestimatedastheratiobetweenthe pressuredistalandproximaltothestenosismeasuredinvasively viapressurewires[5].DeBruyneetal.[6]foundthatFFRwasinde- pendentofthehemodynamicstateofthepersonshowingthatthe measuremayperformacrossavarietyofphysiologicalconditions.

Additionally,numerousstudiesshowthediagnosticaccuracyof FFR(greaterthan90%)andassociatedbenefitsofimprovedtreat- mentoutcomesandreduced costs[7–9].Despitethesebenefits, FFRistypicallymeasuredonlywhenoperatorsareuncertainabout angiographybasedassessmentofthestenosis[10].Hannawietal.

alsofoundthatnearlythree-fourthsofphysicianssurveyeduseFFR inlessthanone-thirdofcases.Manyphysicians(47%)reportedthat theprocedurewasnotavailableattheirinstitutions,whileanother significantportion(25%)claimedthetimerequiredtosetupthe testdiscouragedbroaderusage.Morrisetal.[11]notethatFFRis measuredinlessthan10%ofindividualswhoundergoPCIforCAD intheUKandclaimmanyoperatorsbelievetheirownjudgment basedonimagedataandnoninvasivetestingissufficientdespite evidencetothecontrary.FailuretoemployFFRmorebroadlylikely resultsinunnecessaryPCIforsomewhilesignificantstenosesgo undiagnosedinothers.

MethodsthatestimateFFRnoninvasivelymayimprovediag- nosis and treatmentof CAD bymaking quantitative evaluation of CAD severity accessible to a wider number of patients and physicians.Computationalfluiddynamics(CFD)maybothimprove understandingofCADandreducethecostoftreatmentbyprovid- ingmoreaccuratephysiologicalinformationthroughsimulation.

Recenteffortsshowpotentialtoimprove diagnosisofCAD[12]

byusingCFDmethods fortheestimation ofFFR [13,14], which mayavoidcostly,uncomfortableandriskyproceduresofsedation, catheterization,andinductionofhyperemia.CTangiographymay identifystenosesanddeterminethegeometryofthecoronaryvas- culature.CFDsimulationsbasedoninputdataconsistingofvascular geometry,andbloodpressureandcardiacoutputmaypredictflow andpressurethroughthestenosedarteryandthusalsoFFR(This approachistypicallyreferredtoasvirtualFFRorvFFR).

Despitethistremendouspotential,thenecessityofappropri- ate boundary conditions (BC) and input data for CFD analysis remainsanobstacle.Further,uncertaintiesaboutinputdatamust beaddressedtoensurethereliabilityandrobustnessofvFFRfor clinicaldecisionmaking.

CFDmethods for vFFR typically solvestandard equations of fluidmechanicsoveradomaingeneratedfromimagingdatawith BCbased ondirect measurements or physiologicallymotivated modelswhichcoupletoCFDmodel.Thisapproachrequiresinput dataconsisting ofgeometric measurements,model parameters, and patient specific characteristics for imposition of BC. We investigatedtheeffectofuncertainties abouttheinputdatafor BC and geometric data for a simplified method of vFFR. This modelintegratesadetailednetworkofporcinecoronaryarteries (downtodiameterofapproximately0.1mm)[15]withasimple stenosismodel [16].To evaluatetheinfluence of uncertainties, weemployedpolynomialchaos,metamodeling,andMonteCarlo methodsinordertoestimatethevariabilityofvFFRestimatesfor assumed uncertaintiesaboutthe modelparameters, inputs and boundaryconditions.We firstdevelopeda simplifiedmodel for coronarybloodflowduringhyperemiabasedonMurray’slawand a stenosismodel proposedbyHuo etal. [16].With this model

weemployedpolynomialchaosquantifytheeffectsofuncertain- ties incardiacoutput,arterialpressure, andtheinternal model parameters. Subsequently, we used hierarchical modeling and MarkovchainMonte Carlotoevaluatethediscrepancy between observed geometrical diameters in our data set and the phys- iological diameters based on Murray’s law. We evaluated the sensitivity of vFFR to such uncertain geometric measurements withmetamodelbasedMonte Carlosensitivityanalysis.Finally, weevaluatedtheconsequencesoflimitedimageresolutiononthe estimationofvFFRanddevelopeda novelapproachtoinferthe flowthroughvesselsthatdonotappearinimagingdatausedfor thepredictionofvFFR.

2. Methods

Toinvestigatetheeffects inputdatauncertaintyonvFFRwe developed a simplified model of coronary blood flow that is amenabletoperforminguncertaintyquantificationandsensitiv- ityanalysis,whilepreservingmanyofthecharacteristicsofmore complexmodels,i.e.dependenceongeometricdataandbasedon similarphysicalprinciples.Weconstructedamodelofhyperemic flow in coronary arteries based onMurray’s law (2), measure- mentsofcoronaryarterygeometry,cardiacoutput,arterialblood pressure,rheologicalparametersandmyocardialflowfraction.The modelconsistsofasystemofnonlinearequations(9)thatdeter- minespressureinthecoronaryarteriesandthusalso,vFFRby(10).

Weinvestigatedtheeffectsofuncertaintyinclinicalmeasure- mentsandmodelparameterswithpolynomialchaosmethod.The impactofgeometricuncertaintiesisassessedbyfirstestimating the measurement error with a BayesianHierarchical modeling approachforthegivendatasetandsubsequentlybyemployinga metamodelingapproachcombiningrandomforestregressionand polynomialchaostoestimatesensitivityindicesforeacharterial segment.Additionally,weassessedtheimpactoflimitedimage resolutionbycomparingsimulationsusingthefulldatasettosim- ulationswhichincludeonlyvesselsdowntoacertainthreshold.

Additionally,weproposedanewmethodtoaccountforbloodflow tovesselsthatdonotappearinimagingdata.

Inthefollowingwewillfirstoutlinethemodelforestimating vFFRgivengeometricdataandpatientdataandthendetaileach component ofthe modelin more detail. Subsequently, wewill presentthemethodsofconceptsusedforUQandSA.Finally,we willpresentthedetailsforeachanalysisconducted.

2.1. Coronarybloodflowmodel

Wepresentahemodynamicalmodelbasedontheassumption oflaminarandsteadyflow,andboundaryconditionsbasedonMur- ray’slaw[17].Bulantetal.[18]recentlycomparedpredictedFFR forsteadyand unsteady3Dsimulationsandfoundlessthan1%

difference.Steadyflowalsoallowstheintramyocardialpressure, whichvariesoverthecardiaccycle,tobeneglectedasitsaverage effectsareaccountedforintheterminalresistancesandassumed geometries.

Basedontheseassumptionswemodeledcoronarybloodflowas alumpedparameternetworkofresistorswithinflow,q,inletarte- rialpressure,pa,andvenouspressure,pv,attheoutlets.Inprinciple thismodelcouldbeappliedtoanysetofappropriategeometric measurementsofcoronarynetworkanatomy,thoughweconsider onlyasinglenetworktopologybasedondetailedgeometricaldata froma porcine heart [15]. Theresistances are calculated using vesseldimensionsandPoiseuille’slaw.Tocompletethemodel,a stenosissubmodelbasedontheworkofHuoetal.[16]isinserted intothecoronarynetwork,and thevFFRiscalculatedbased on thecompletenetwork,boundaryconditions,andflowdistribution

models.Theterminalresistancesdeterminedforthebaselinenet- workarescaledinordertoreflecttheusageofadenosinetoinduce hyperemiainpatientswhenevaluatingFFR.Thefollowingsections describethegeometricdataandsubsequentlyexpandoneachof componentsofthesemathematicalmodelsofcoronarybloodflow.

2.1.1. Geometricdata

Kassabetal.[15] madesiliconeelastomercastsofa porcine coronarynetwork.Fromthesecaststheyrecordeddetaileddescrip- tionsofthethreemainbranchesofthecoronarynetwork,i.e.the rightcoronaryartery(RCA),theleftcoronaryartery(LCA),andthe circumflexartery(CX).The geometryof each branchis charac- terizedbynotingthelengthlanddiameterdofthemaintrunk betweeneachbifurcation.Ateachbifurcationthediametersofthe outletsarealsoreported(seeFig.1B),thusdescribingthemain coronarybranchesasasequenceofvesselsegmentswithoutlets fromthemaintrunklocatedateachjunction.Toexplicitlyrepre- sentthisdata,wedefinethesetN=

(di,li,b_i)

whereiindexes thevesselsfromtheinlet,i=0,totheterminalsegment,i=ns,and b_i=(bi1,...,bin_i)isalistofthen_ioutletdiametersatthejunction ofsegmentsiandi+1.Superscriptnotation,N,denotestheset correspondingtoaparticularbranchofthecoronaryarteries.

Porcinecoronaries are widelyused as a biological model of humancoronariesandtoalargeextenttheanatomiesaresimilar, aswellastherelativesizeincomparisontotalbodysize.Levolas etal.[19]reviewthesimilaritiesanddifferencesofporcineand humanheartsandconcludeporcineheartsarethemostappro- priatebiologicalmodelfor cardiovascularresearchandthatthe hemodynamicsaresufficientlysimilarfortestingandevaluating methodsandequipmentintendedforbiomedicalapplications.The grosssizeoftheheartandcoronariesaswellasthecardiacoutput arequitesimilartothoseofhumans,thuswhilenotanexactmatch ofhemodynamicpatternsthegeneralcharacteristicsanddriving factorsarequitesimilar.

2.1.2. Hyperemicbloodflowmodel

TopredictFFRfromasetofgeometricmeasurementssomesort ofmodel,typicallybasedontheunderlyinghemodynamics,must beemployedtopredictthepressuresandflowsinthecoronaries.

Themodelconsistsofadescriptionofthehemodynamicswithinthe domainwheredataisavailableandthedeterminationofboundary conditionstodescribetheflowsandorpressuresattheinletsand outletsofthenetwork.

First,theboundaryconditionsof therestingstatearedeter- minedbyassumingtheflowtothecoronariesisproportionalto thecardiacoutputqco

qmyo=myoqco. (1)

Todeterminetheboundaryconditionsattheoutlets,weapplied Murray’slawtodeterminetheflowateachoutlet.Murray’slaw [17]predictsthatunderrestingconditionstheflowinagivenblood vesselis

q=ad^c, (2)

whereaisacoefficientdeterminedbytherelationshipbetween viscouslossesandmetabolicdemandsanddisthevesseldiameter.

Applyingthisrelationshiptotheoutletdiametersandenforcing conservationofmassimpliesthat

qmyo=

i

j

ab^c_ij, (3)

andthusthefractionoftheqmyothatflowsthroughthejthoutlet atthejunctionofsegmentsiandi+1is

qij

qmyo= b^c_ij

i

jb^c_ij (4)

asthecoefficientacancels(bijisthediameteroftheoutlet).Simi- larly,theflowtoaparticularcoronarybranchis

q=

i∈N n_i

j=1

ijqmyo, (5)

whereisoneofCX,LAD,orRCA,jdenotesthesegmentindex withinN,andij=_qmyo^qij istheflowfractiontoaspecificoutlet.

Giventheflowintoaparticularbranchandtheflowleavingat eachjunction,theflowinaparticularsegmentis

qi+1=qi−

j

ijqmyo. (6)

whereq_i+1denotetheflowfromsegmentitosegmenti+1.

Theseassumptionsfullydeterminetheflowthroughthecoro- naries.Todeterminethepressuredistributionweassumedthatthe pressureatrootofeachcoronarybranchwassimplythesystemic arterialpressurepa.Assumingeachsegmentofcoronaryvascula- tureiscylindricalandbloodflowislaminarandsteady,thepressure dropalongaparticularsegmentofcoronaryvasculatureisgivenby Poiseuille’slaw:

pi+1−pi= 128li

d_i⁴ qi, (7)

whereistheviscosityoftheblood,liisthelengthofthesegment andd_iisthediameterofthesegment.Eqs.(5)–(7)determinethe flowsandpressuresinthecoronarynetworkinrestingconditions.

ThephysiologicalprinciplebehindFFRistoquantifytheimpact ofstenosisastheamountthatitimpedesmaximalphysiological flowrelativetoahealthyvessel.Thus,hyperemicconditionsare inducedbypharmacologicallydilatingthecoronarymicrovascu- latureand decreasingitsresistancetoflow tomeasure FFR.To simulatethis,weassumetheresistancetoflowunderhyperemia R^hyp_ij =˛RijwhereRijistheoutletresistanceinrestingconditions and˛isthetotalcoronaryresistanceindex(TCRI)corresponding totheamounttheresistanceisdecreasedduringhyperemia.The resistanceinrestingconditionsiscalculatedbyapplyingOhm’slaw (p_i+1−pv)=R_ijq_ij,whichimplies

Rij=p_i+1−pv

ijqmyo , (8)

wherepvisthevenouspressure.

Asadenosineprimarilyaffectsthemyocardium,weassumeas othershavethatRsysdoesnotchangeduringinducedhyperemia andfurtherthatqcoisconstant[20,14].Thearterialpressuredur- inghyperemiaispa=qsysRsys,whereqsys=qco−qmyo,andRsysmay bedeterminedfromthecoronaryflow,cardiacoutputandarterial pressureinrestingconditions.

Insteadof(7),thepressurelossacrossastenosedsegmentis treatedasanonlinearfunctionofflow,stenosisgeometry,andrhe- ologicalparametersoftheblood,p_i+1−p_i=f(q_i,d_i,ds,ls,,),as thesimplificationsduetoassumedlaminar andsymmetricflow areviolatedinastenosis,Weemployedastenosismodeldeveloped andvalidatedagainstocclusionsofbothhumancarotidandporcine coronaryarteries(R²=0.9andR²=0.7respectively)byHuoetal.

[16].Thefunctionfispresentedindetailintheappendix(Section A.1).

Thereducedoutletresistances,therelationshipbetweenarte- rialpressureandsystemicflow,andthestenosismodelcomplete

Fig.1. PanelAshowsthecoronaryarteriesontheheart,panelBshowstheextractedleftcoronaryarteryandassociateddimensions,panelCshowstheresultinglumped parametermodelofresistances,astenosismodel,andterminalimpedances.(PanelsAandBmodifiedfromtheworkofPatrickJ.Lynch,medicalillustratorderivativework:

FredtheOysteradaptionandfurtherlabeling:MikaelHäggström–withpermissionCCBY-SA3.0.)Notethefullcoronarycirculationimagedconsistsofthreebrancheslike theoneshowninpanelC.

themodelofcoronaryflowduringhyperemia,whichimposes a pressureconditionattheinletsofthecoronariesandresistancesat theoutlets.Theflowsandpressuresmaybedeterminedbysatis- fyingtheseequationsaswellastheconservationofmassforeach subnetworkN^CX,N^LAD,andN^RCA.Thesystemofequationsisthus

qmyo =qCX+qLAD+qRCA

qsys =qco−qmyo

pa =Rsysqsys

qmyo =

∈CX,LAD,RCA

k∈N n_k

l=1

qkl

qij = p_i+1−pv

˛R_ij qi+1 =qi−

n_i

j

qij

p_i+1 =

⎧ ⎪

⎨

⎪ ⎩

pi−f(qi,di,ds,ls,,) if segmentiis stenosed

pi−128li

d⁴_i qi otherwise.

(9)

Thisnonlinearsystemofequationsintermsoftheflowstoeach branchofthecoronarynetworkandthenodalpressuresmaybe solvedusinganonlinearequationsolversuchasfsolvefromthe Pythonlibraryscipy.optimize.

Insummarygivenanominalcardiacoutputqco,arterialpressure pa,andmyocardialflowfractionmyo,theoutletresistancesmay bedeterminedfromEqs.(4)–(8).Additionally,thesystemicresis- tanceisdeterminedbyRsys= ₍₁₋^p_myo^a _)qco.Oncethesehavebeen calculatedthestenosedmodel(9)issolvedfortheflowsandpres- suresineachsegmentfromwhichvFFRiscalculatedas

vFFR=p_i+1

pa , (10)

wheresegmentiisthestenosedsegment.

2.2. GenerationofsimulatedCT-imagingbaseddata

ThegeometryandtopologyusedforvFFRdependsontheves- selsidentifiedbysegmentationofclinicalimagingdata.Thus,the limitsofclinicalCT-imagingresolutionandvesselsegmentation willpreventobservationofvesselssmallerthansomethreshold, d_th.Consequently,CT-imagingbasedvFFRwillnotaccountforany flowthatgoesthroughvesselssmallerthanthisthreshold.Wesim- ulatedCT-angiographydatabygeneratinganewnetworkNd_thfrom thecastdataconsistingofonlytheoutletsb_ijlargerthanathreshold d_th.TosimulatethepredictionvFFRfromclinicallyderivedimag- ingdata,thesamemathematicalmodel,Eqs.(4)–(8)and(9)was appliedusingonlythedatainNd_th.

2.3. Methodtoestimateflowtounobservedvessels

Asmanyvesselswillnotappearinimagingdata,amethodthat compensatesforthislackofinformationishighlydesirable.We proposeanextensiontothemodelthatestimatestheamountof bloodthat“leaks”outofvisiblevesselsthroughunobservedves- selsinimagingdata.Theleakingbloodflowmaybeaccountedfor byaddinganoutletresistancetoeachobservedvesselsuchthat Murray’slaw,(3),issatisfied,i.e.qi=ad^c_i =

jab^c_ij+q_i+1whereb_ij isthediameterofthejthoutletatthejunction.Theadditionalresis- tanceisaddedbyappendinganadditionaloutletdiameterb_i,leakto b_isuchthat

ad^c_i =

j

ab^c_ij+ab^c_i,leak+ad^c_i+1, (11)

whichsimplifiesto b^c_i,leak=d^c_i−

j

b^c_ij−d^c_i+1. (12)

Onceallb_i,leakaredetermined,thesameprocedureasoutlinedin Section2.1.2maybeappliedtodeterminethenetworkresistances, pressures and subsequently vFFR (vFFR_leak is usedto denote a vFFRcalculatedspecificallyusingtheleakagemodel).Notethatthe valueR_i,leak= ^pi+1

i,leakqmyo representstheequivalentleakage resis- tancepresentatthejunctionbetweensegmentsiandi+1.

2.4. Summaryofmodels

Wehavepresentedamodelofcoronarybloodflowbasedon physicalprinciples,geometricaldataandvaluesforcardiacoutput, myocardialflow,andarterialpressure.Theflowandpressurein thecoronarynetworkmaybepredictedfornominalandhyper- emicconditions.ThesepredictionsalsodirectlyestimatevFFR.We haveproposedmethodtosimulatelimitedclinicallyderivedgeo- metricdata,aswellasamethodthatcouldbeemployedtoinferthe presenceofvesselswhichdonotappearintheexplicitgeometric data.

Wewillnowpresentanumberofmethodsforanalysisofthe effectsofuncertaintyintheinputdata.Firstwewillintroducethe basicconceptsofuncertaintyquantificationandsensitivityanaly- sis.Subsequently,wewilldescribethespecificmethodsweemploy inthisstudy.

2.5. Uncertaintyquantificationandsensitivityanalysis

Weperformeduncertaintyquantificationandsensitivityanal- ysis(UQSA)ofthecoronarybloodflowmodelusingpolynomial chaos,metamodelingandMonteCarloalgorithms.Anoverviewof theprimaryconceptsofUQSAmaybefoundintherecentreview byEcketal.[21].Forcaseswithfewerthan10inputparameters weusedpolynomialchaos.Whenconsideringmoreinputparame- ters,weemployedMonteCarlomethodsandfurthermetamodeling inordertoefficientlyevaluatethesensitivityofthemodel.The characterizationofinputuncertaintyisakeypartofperforming UQSA.Wherepossiblewerefertoknowncharacterizationsofmea- surementuncertainty,populationvariability,andreportedranges ofmodelparameterstocharacterizetherangeofinputparame- ters.Insomecasesinferenceofmodelparametersandassociated uncertaintiesisperformedbasedonthenetworkdatausinghier- archicalmodeling.Detaileddescriptionsoftheexactmethodsused andthestatisticalmodelforparameterestimationarereportedin thesubsequentsections.

Asthespecificmethodsarecoveredindetailinthereferenced work,wepresentonlythebasicconceptsandnotation.Thesources ofuncertaintyincludemeasurementerror,intrinsicbiologicalvari- ability,andassumingpopulationoraveragevaluesforunmeasured parameters.Themodelistreated asfunctionoftheseinputs,Z, whichyieldsanoutputofinterestY.Inthispaper,themodelout- putisY=vFFR,while themodelparametersand geometrydata arethe modelinputs. A schematicdiagram of theflow ofdata inuncertainty quantificationand sensitivityanalysisisgiven in Fig.2.ToconductUQSA,firstdeterminetheprobabilitydistribu- tionofZtoreflecttheuncertainties,thencalculatethedistribution ofYbypropagatingZthroughthemodel.TheuncertaintyofYis characterizedbyitsvariance andassociated toleranceintervals,

Y˛/2,Y1−˛/2

,containingto(1−˛)×100percentoftheprobabil- itydensity.FurtherthevarianceofYispartitionedtocalculateSobol sensitivityindicesS_iandS_Tiwhichquantifytheamountofvariance duetotheithcomponentofZ[22,23].S_iisthefirstordersensitivity indexandquantifiesthedirectinfluenceofZiwhile,STiisthetotal sensitivityindexwhichincludestheaverageimpactofZ_iincluding itsinteractionswithothersourcesofuncertainty.

2.5.1. Polynomialchaos

Polynomial chaos expansions approximate a function of stochasticvariables asa sumof multivariate basis polynomials orthogonal withrespect tothe jointprobability distribution of theinputvariables[24].Orthogonalityenablesmoreefficientcal- culationof uncertaintyandsensitivity measuresofthefunction incomparisontoMonteCarlomethodsforsufficientlylownum- bersofinputs.ForanoverviewofUQSAmethodsandpolynomial chaosexpansions,wereferthereadertotherecentpaperonthe

Fig.2.Schematicoftheflowofdataandcomponentsinthemodelinganduncer- taintyquantificationframework.

UQSAforcardiovascularmodels[21].WeusedthePythontool- boxchaospydevelopedbyFeinbergandLangtangentoperform polynomialchaos[25].

2.5.2. MetamodelassistedMonteCarlo

Forcaseswherealargenumberofmodelinputsvaryindepen- dentlyitisnolongercomputationallyefficienttousepolynomial chaosas thenumber of terms –and thus model evaluations– growsrapidlywiththenumberofinputsaccordingtothebinomial coefficient

D+p p

=(D+p)!

p!D! whereDisthenumberofinputs andpisthepolynomialorder.Further,MonteCarlomethodsmay require a prohibitive amount of simulations before converging toaccurateestimatesofsensitivityindices.UQSAofmodelswith largenumbersofinputsthusrequiresalternativemethods.

Oneapproachistoapproximatethemodelwithamoreefficient metamodelandtoperformsensitivityanalysisonthismetamodel asaproxyforthefullmodel.ToestimatetheuncertaintyofvFFRdue touncertaintiesinindividualdiameters,wefirstcreatedameta- modelofthefullmodelusingrandomforestregression[26]and subsequentlyestimatedsensitivityindicesoftheoutletdiameters byMonteCarlomethodsappliedtothemetamodel[27,28].Aran- domforestmetamodelprovidesaveryefficientwaytoevaluatethe largenumberofsamplesrequiredforhighdimensionalMonteCarlo sensitivityanalysis.Withminimaltuningofthefittingprocedure,

randomforestregressionalsotendstoachievea goodaccuracy whileavoidingoverfitting[29].Similarapproachesforsensitivity analysisusingmetamodelingandMonteCarolestimationproved toagreequalitativelyintermsofsensitivityindices[30,31].

Aregressiontreeconsistsofasequenceofbinaryruleswherethe firstbranchisselectedifz>z_candthesecondotherwise.Foreach stageselectthevariablez_iandcutoffzcsuchthatthevalues>=

1 N>

z_i>z_cf(z) and <= _N¹_<

z_i<z_cf(z) minimize

z_i>z_c||>− f(z)||+

z_i<zc||<−f(z)||where ||·|| can bean arbitraryerror metric.<and>neednotbesimplytheaverageeither.Therules foreachofthedaughterbranchesaregeneratedinthesameway, butusingonlythedatacorrespondingtoeachbranch.Thisiscontin- ueddowntoacertainnumberofgenerationsoruntilacertainerror thresholdisachievedforeachleaf,i.e.abranchwithnochildren, ofthetree.

Randomforestregressionaverages,N_f,regressiontreescreated usingdifferentrandomlyselectedsubsetsoftheinputparameters, e.g.onlythevaluesofz₁,z₃,andz₆foronetreeandz₁,z₄,z₅,andz₆ foranother.Further,eachtreeisfitonlytoarandomsubsampleof thedata.Thepredictionisthentheaverageofalloftheindividual predictionstheN_fregressiontrees.

2.5.3. Parameterestimationandinputuncertaintyquantification Sofarwehavepresenteduncertaintyquantificationfromthe perspectiveoferrorpropagation,yetthecharacterizationofthe uncertaintyaboutthedataandparametersisalsoimportant.In thisstudyweemployBayesianhierarchicalmodelingtoestimate theuncertaintyaboutthemeasuredradiiandMurray’sexponent, c.Inthisapproach,thedataareassumedtofollowacertaindistri- butionconditionedonthevaluesofthemodelinputparameters, e.g.independentnormallydistributedvalueswithmeanvaluecor- respondingtothemodel

yi|␪i,e∼N(m(␪i),e) (13) wherey_i denotestheithdatapoint,mthemodel,␪ithemodel parameters’valuesfortheithdata,andetheerrorstandarddevi- ation.Therelationshipbetweenparametersisdescribedthrough ahierarchyofconditionaldistributions,e.g. | ∼D,whereD representsahypothesizeddistributionoftheparameter␾condi- tionalon .Thenumberoflevelscanbechosentoaccommodatethe hypothesizedrelationshipsbetweenthevariablesandparameters.

ForamoredetaileddescriptionofBayesianHierarchicalmodeling, wereferthereadertoGelmanetal.[32].

Giventhehierarchythedistributionofthemodelparameters giventhedataisestimatedbasedonBayes’theorem

p(␪,e|y)= p(␪,e,y)

p(y) = p(y|␪,e)p(␪,e)

p(y) . (14)

Thehierarchicalapproachprovidesaconceptuallysimplewayto buildupp(y|␪,e)p(␪,e)asp(y|␪,e)isdeterminedby(13)and p(␪,e)bytheD.Theposteriordensityof␪ande(thelefthand sideof(14))quantifiesthemostcharacteristicvaluesoftheparame- tersbasedonthedataaswellastheuncertaintyaboutthesevalues, i.e.thewidth,varianceorquantilesofthisdensity.Thedenominator oftheequationisdefinedby

p(y) =

_␪×e

p(y,␪,e)d␪de

=

_␪×e

p(y|␪,e)p(␪,e)dde

(15)

whereandedenotetheregionswhere␪andehavenonzero probabilitydensity.Directevaluationof(14)requirescalculationof this(potentially)highdimensionalintegral.Further,calculatingthe

Table1

Tableofassumedinputuncertaintiesforevaluationoftheeffectsofparametric inputuncertaintyonvFFR.UniformdistributionsaredenotedbyU(min, max)and normaldistributionsN(mean, std. dev.).

Input Symbol Nominalvalue Distribution

Cardiacoutput qco 6L/min N(6,1.0)

Arterialpressure pa 93.0mmHg N(93,4.77)

Myocardialflowfraction myo 0.045 U(0.04,0.05)

TCRI ˛ 0.23 U(0.15,0.30)

Murray’slawexponent c 3.0 U(2.4,3.0)

Hematocrit H 0.45 N(0.45,0.08)

mean,median,standarddeviationorpercentilesofandrequires evaluatingsimilarintegrals,e.g.E(␪)=

_␪×e␪p(␪,e|y)d␪de. Thedirectevaluationof(15)maybeavoidedbysamplingthe posteriordistributionand calculatingthecorrespondingsample statistics.MarkovchainMonteCarlo(MCMC)methodsallowsam- plingtheunknownposteriordensitybyproduceachainofsamples basedonthelikelihood function,p(y|␪,e)p(␪,e).Thedetails ofthis process and a specificimplementationarepresented by Salvatieretal.[33]whodevelopedthePythontoolboxPYMC3for performingMCMC.

2.6. Analyses

Weappliedtheaforementionedmethodsofuncertaintyquan- tification and sensitivity analysis to the vFFR model for three specificanalyses.First, we consideredtheimpactof inputdata uncertaintyintheabsenceofgeometricuncertainty.Second,we employed MCMC to quantify the uncertainty about the geo- metric data and the parameter c in Murray’s law, and then performeduncertainty quantificationand sensitivity analysisof vFFRgivengeometricuncertaintiesbasedontheestimatedparam- eters.Finally,weassessedtheperformanceofthemodelforlimited imageresolutiondataandevaluatedtheperformanceofournewly proposedmissingvesselcorrectionmethod.Eachanalysisispre- sentedindetailinthefollowingsections.

2.6.1. Impactofnon-geometricinputuncertainties

Toaccountforuncertaintiesaboutthemodelparametersand input values, we performed UQSA via the polynomial chaos approach.Usingthewholecoronarynetworkmodel(seeSection 2.1.2),wefirstconsideredtheeffectsofparametricuncertaintyon vFFRfora75%areastenosis1mmlongtoreplacesegment10ofthe LAD,theentrancelengthwastakentobethelengthofthesegment.

Thevaluesofthearterialpressure,cardiacoutput,myocardial flow split,TCRI,and hematocrit, which determinesdensity and viscosity,weretreatedasuncertaininputs.WeusedHammersley samplingtogeneratetwiceasmanysamples(924)ascoefficients fora5thorderpolynomialapproximationin6inputvariablesand estimatedthecoefficientsusingregression.Hammersleysampling wasusedduetoitsnestednaturewhichisconvenientforassess- ingconvergence.Thepolynomialorderwaschosensuchthatwe foundtheestimatesofsensitivityindiceswerelessthan1%differ- entbetweensuccessiveorders.Theprimaryquantitiesofinterest werea95%predictionintervalforvFFR,andthesensitivityindices ofvFFRforeachinputparameter.Thesewerecalculatedfromthe polynomialchaosexpansion.

WeappliedUQSAtounderstandtheconsequencesofvariability ininputparametersandboundaryconditionsonthepredictionof vFFRinclinicalsettings.Tothisend,weattemptedtocharacter- izethevariationoftheseinputsaccordingtoclinicalmeasurement modalities,and assumed a normalpatientasidefromthepres- enceofthestenosisforcharacterizingthevariabilityofunmeasured parameters.Inthefollowingsection,wepresentthebasisforthe assumedinputuncertaintiessummarized inTable1.Inputdis-

tributionswerecharacterizedaseithernormalrandomvariables whenmean andstandard deviationweregiven andas uniform randomvariablesovertheplausiblerangeotherwise.

Arterialpressure.Thesphygmomanometeristhecurrentclin- icalstandardfornoninvasivebloodpressuremeasurement.This methodwasfoundtohavestandarddeviationsof3.3mmHgand 5.5mmHginsystolicanddiastolicpressuremeasurementrespec- tively[34].Themeanpressureisfurtherapproximatedusingthe widelyusedformulapa= ²3p_d+¹3ps[35],which,assumingperfect positivecorrelationbetweensystolicanddiastolicmeasurements hasstandarddeviation4.77mmHg.Wemodelthemeasurement ofmeanpressureasnormalrandomvariablewithmean93mmHg andstandarddeviation4.77mmHg.

Cardiac output. The gold standard measurement of cardiac outputrequires invasivecatheterization[36,37];however, non- invasivemethodsderivedfromultrasoundshowgoodagreement withinvasivemeasurements[38].Studiescomparingultrasound measurementstothepulmonary arterycatheterizationthermal dilutionmethodhavefoundthedifferenceinmeasurementshad astandarddeviationofbetween0.82L/minand1L/minwithmean biaseswerebetween0.03L/minand0.18L/min[39,40].Basedon thesestudiestheclinicalmeasurementofcardiacoutputwasrep- resentedasnormallydistributedwithameanof6L/minastandard deviationof1L/min.

Myocardialflowfraction.Ideallytheflow tothemyocardium couldbemeasureddirectlyasrecentstudieshaveshownpromise towardsthisend[41];however,mostofthesemethodsrequire bothdyeinjectionandcatheterizationtoaccuratelymeasurethe volumetricflow.Thus,weassumetheflowtothecoronaryarteries isafractionofcardiacoutput,myo,whichisbetween4%to5%on average[42].

Bloodcharacteristics. Thedensity and viscosity of blood areimportant parametersin describingthemechanicsof blood flow and both dependon hematocrit levels, H. Sankaran et al.

[43]describethisusingtherelationship=_(1−H)^p2.5 andnotethat apopulationaveragehematocritof0.45 witha standarddevia- tionof0.08. Blooddensitydepends onhematocritaccordingto =eH+(1−H)p,wheretheplasmadensity,p,is1018kg/m⁻³ anddensityoferythrocytes,e,is1085kg/m⁻³[44].Wemodel andaccordingtotheserelationshipsandhematocritasanor- malrandomvariablewiththeaforementionedmeanandstandard deviation.

2.6.2. Estimationofgeometricuncertaintyanditsassociated effectonvFFR

Giventheavailabledataitisofinteresttoestimatewhatval- uesoftheexponentcareplausible.Additionally,thesmallsizeof manyarteriesaswellasthehighpressureinjectionprocessusedto producingacastofthecoronaryarteries[15,45]couldleadtodis- crepanciesbetweenthemeasureddiametersandthephysiological diameters.WeemployedBayesianhierarchicalmodelingtoformu- lateastatisticalmodelthatcanaccountforuncertaintyaboutbothc andthetruephysiologicaldiametersandestimateposteriorprob- abilitydistributionsofcandthegeometricuncertaintygiventhe measureddatausingaMarkovchainMonteCarlomethod.

Themeasureddiametersd_iwereassumedtorelatetothephys- iologicaldiameters ˜dias

di=d˜i(1+Ei) (16) where E_i are independently and identically distributed normal randomvariableswithmean0andstandarddeviationerepre- sentingaproportionalerrorofd_irelativeto ˜di.Thelikelihoodof d_igiventhephysiologicaldiameter ˜d_iandeissimplyanormal distributionwithmeani=d˜iandstandarddeviationi=e·i. Murray’s law (3) implies that for physiological diameters ˜di=

n_s j=i

n_j k=1(_1+E^bjk

jk)^c

¹_c

where iindicatestheinletsegmenttoa subnetworkN_i consistingofthesegmentsdistaltoandinclud- ingsegmentiinbranch.Thus,theconditionaldistributionofd_i giventheerrorsoftheoutlets,E_jk,andtheMurrayexponent,c,was specifiedforeachsubnetwork,N_i .

Tocompletethehierarchicalmodel,weassumeduniformprior distributionsforeandcaslittlepriorknowledgeaboutthephys- iologicaldiameterswasavailable.Theerrorswereassumedtobe lessthan50%fore.Forctheintervalwas(2,4)corresponding totherangeofvaluesproposedintheliterature.(Weconsidered widerpriorsforcande,butitdidnotchangethefinalinferences significantly.)Thus,thecompletehierarchicalmodelis

i =

⎧ ⎨

⎩

k

j=i n_j

l=1

( bjl

1+Ejl

)

c

+d˜^c_k+1

⎫ ⎬

⎭

1 c

di|Eij,c,e ∼N(i,e·i) Eij ∼N(0,e)

c ∼U(2,4) e ∼U(0,0.5)

(17)

Toestimatetheposteriordistributions ofe and cgiven the measureddiametersweperformedaMarkovchainMonteCarol simulation as described in Section 2.5.3 withthe Python tool- boxPYMC3[33].Thechainwasrunforonemillionstepsandthe posteriordistributionsfortheaveragevalueofeandcwerechar- acterizednumericallyfromtheMarkovchain.

2.6.3. Impactofgeometricinputuncertainties

Duetotheevidenceofgeometricuncertaintyaboutthemea- sured diameters inferred from the analysis of Murray’s law presentedintheprevioussection(Section2.5.3), weconsidered the effect this uncertainty would have on vFFR. To assess the effectsofthesediametricuncertainties,weanalyzedthesame75%

areastenosiscaseasthepreviousanalysis.A1mmlongstenosis replacedsegment10oftheLAD.Thepre-stenoticlengthwastaken tobethelengthofthesegment.Thevaluesoftheoutletdiameters wereassumedtofollow

b˜ij=bij(1+Eij), (18) whereE_ijareindependentlyandidenticallydistributednormalran- domvariableswithmeanzeroandstandarddeviationbasedonthe meanofe(0.073)fromtheMarkovchainMonteCarloanalysis.

Toconstructthemetamodel,weevaluatedvFFRwiththecoro- narymodelpresentedinSection2.1.2for5000differentsamples oftheoutletdiametersdeterminedbyLatinHypercubesampling overthejointdistributionofoutletdiameters. (Latinhypercube samplingwaschosenasitremainsaLatinHypercubesampleas parameterdimensionsareeliminated,i.e.wellsuitedforeachtree intherandomforest.)Wethenfittedarandomforestregression modeltothevaluesgeneratedatthesesamplepointswithRandom- ForestRegressorfromthePythonpackagescikit-learnusing thedefaultparameters.Thenumberofsubtrees,100,wasselected toensuresufficientaccuracyofthemeta-modelsuchthattheaver- ageout-of-bagerrorwaslessthan0.01andthemaximumabsolute errorwaslessthan0.02.

Toestimatethesensitivitytoparticularoutletdiameters,we computedthevFFRviatherandomforesttodeterminesensitivity indicesofvFFRtoeachoutletbasedontheMonteCarloalgorithms proposedbySaltelli(forS_i)andSobol(forS_T)withsamplematrices AandBof50,000samplessuchthatadifferenceoflessthan10⁻³ wasachievedbetweenSTforsuccessivesamplesof10,000[46].As

Fig.3.Sobolsensitivityindices(firstordergraywithhatchesandtotalblack)forthefullcoronaryarterymodelsubjectedtouncertaintyacrossallinputparameters.Cardiac output,Q,andTCRI,˛,arethemostsignificantinputsintheresultinguncertaintyofvFFR.

atestofthisapproach, wecomparedpolynomialchaosandthe metamodelapproachonasmallproblemofonly6inputs,such thatPCwasfeasible,andfoundgoodagreementbetweenthetwo approaches.

2.6.4. Analysisoftheimpactsoflimitedresolutionimaging

Clinicalassessmentofcoronaryarteriesislimitedbythereso- lutionofimagingmodalitiesaswellastheprocessofsegmenting andidentifyingvesselsfromimages.CTangiographycanresolve objectsaround0.5mminoptimalconditions,buttheresolution maybereducedtocompensatefornoiseandinclinicalsettingsis closerto1mm[47,48].Theprocessofsegmentationrequiresgood qualityimagesoverthelengthofthevesselandthusinmoststudies thevesselsanalyzedaretypicallylargerthat1.5mm[49–58].

TocharacterizetheimpactvesselresolutionmayhaveonvFFR, weemployedtheprocedureoutlinedinSection2.2tosimulateclin- icalimagingofthecoronarynetworkfromwhichthecastdatawas obtained.Astheimpactofimageresolutionlikelydependsonthe locationofthestenosedsegment,weevaluatedvFFRfromtheimag- ingbasednetwork Nd_th for eachpossiblelocationofa75%area stenosis1mmlong.Theerrorduetothelimitedresolutionwas quantifiedby

RMSE(dth)=

1 Ns

i

vFFR(Nd_th,,i)−vFFR(N,,i) (19)

wherevFFR(N,,i)indicatesthepredictedvalueofvFFRgiventhe networkdataNwithastenosislocatedinsegmentiofbranch andN_s=94isthetotalnumberofsegmentsinallbranches.Asthe thresholddiameterforobservedvesselsisnotanexactvalue,this errorwascalculatedatvariousvaluesofd_thfrom0mmto2mm toprovideamorecompletepictureofthepotentialimpactimage resolutioncouldhaveonpredictedvFFR.

Justastheimpactsofunobservedvesselswerequantifiedby (19),theperformanceoftheproposedleakagemodel,seeSection 2.3,atagivenvisibilitythresholdwasassessedinthesamemanner bycomparingvFFR_leak(Nd_th)andvFFR(N)forallpossiblestenosis locationsofa75%areastenosis1mmlong.

2.7. Summaryofmethods

Insummarywehavepresentedamodelofhyperemicflowin thecoronaryarteriesbasedonMurray’slaw(2),measurementsof coronaryarterygeometry,cardiacoutput,arterialbloodpressure, rheologicalparametersandmyocardialflow fraction.Ourmodel consistsofasystemofnonlinearequations(9)whichwhensolved forthenodalpressuresmaybeusedtocalculatevFFRby(10).We

thenpresentedtwomethodsforpropagatinguncertaintythrough thismodel.Thepolynomialchaosapproachwasappliedtoquan- tifyuncertaintyandsensitivitywhenassuminggeometricaldata washighlyaccurate.Anapproachusingrandomforestregression togenerateametamodelwaspresentedforquantificationofuncer- taintyandsensitivityindiceswhenassuminguncertaintyaboutall outletdiameters.Thisapproachismoresuitablethanpolynomial chaosduetothelargenumberofinputparameterswhenconsid- ering114outlets.We alsopresent ahierarchicalmodel (17)to estimatethevariabilityoftheMurraycoefficientcandthephys- iologicaldiameterswhenassumingthatMurray’slawaccurately representsthetruephysiology.We furtherpresentedamethod forsimulatingtheconsequencesofreducedimageresolutionby removinggeometricaldatabelowacertainthreshold,aswellasa methodthatcorrectsformissinggeometricaldatabyestimatinga leakageflowateachjunctiontoaccountforflowthroughinvisible vessels.

3. Results

Weappliedthemethodspresentedintheprevioussectionto performthreeanalysesoftheconsequencesofuncertaintiesindata requiredforcomputingvFFR.First,weinvestigatedthecasewhere onlytheinternalmodelparametersandinletboundaryconditions weretreatedasuncertain.Second,wesimulatedtheconsequences oflowresolutionimagingbycomparingthevFFRmodel(Section 2.1.2)forfullgeometricdatatothecasewhereonlygeometricdata correspondingtodiameterslargerthanathresholdd_thwhichvar- iedfrom0mmto2mm.Finally,weconsidereduncertaintyabout thegeometryitselfbyevaluatingthesensitivityofvFFRtoerrorsin measurementofoutletdiametersperformedhierarchicalmodel- ingtoestimatethediscrepancybetweenmeasureddiametersand physiologicaldiametersbasedonMurray’slaw.Theresultsofeach oftheseinvestigationsarepresentedinthefollowingsection.

3.1. Simplestenosismodel

TheuncertaintyquantificationofvFFRforasinglestenosisin theLADconsideredthecardiacoutput(qco),aorticpressure(pa), myocardialflowfraction(myo),Murray’slawexponent(c),hyper- emicresistancecoefficient(˛),and hematocrit(H)assourcesof uncertainty.Thesewerebasedontherangeofvaluesreportedin theliteraturereviewedinthemethodssectionandsummarizedin Table1.TheresultingestimateofvFFRhadmean0.78and95%

predictionintervalof(0.59,0.90),notethethresholdforclinical significanceof0.75.Fig.3showsfirstorderandtotalSobolsensi-

Fig.4.PosteriordistributionsofMurray’slawcoefficientcandstandarddeviationoftheerrorbetweenmeasuredandphysiologicaldiametersestimatedusingMarkov ChainMonteCarlo(see2.5.3).

tivityindicescalculatedfromthepolynomialexpansionofwhich themeasurementsofcardiacoutput(qco)andtheestimateofTCRI (˛)werethelargestsourcesofuncertainty.Interestingisthefact thatarterialpressureandtheMurray’slawcoefficientccontribute almostnouncertainty,suggestingthattheexactvalueoftheMur- ray’slawcoefficientisnotsoimportantindeterminingvFFR.

3.2. EstimationofvariabilitybasedonMurray’slaw

SinceMurray’slawprovidesthebasisofdeterminingoutflows for thesimulationof blood flow throughthecoronary vascula- ture,itisofinteresttoevaluatetheaccuracyofthisprinciplewith respecttothegeometricdata,and inparticulartoestimatethe bestfitoftheexponentc,aswellastoquantifytheuncertainty aboutthisestimate.WeemployedMarkovChainMonteCarloto performBayesianparameterestimationforthehierarchicalmodel ofMurray’slawanddiameteruncertainty(see(17)andSection 2.5.3).Fig.4showstheresultingposteriordistributionswherechas mean2.47andstandarddeviation0.040,andehasmean0.073 and standarddeviation 0.005.These resultssuggestthat differ- encesbetweenphysiologicaldiametersandmeasureddiameters aresomewhatuncertain,whereastheexponentvalueofcisfairly certain(standarddeviationoflessthan2%ofthemeanvalue)given thisdata.

3.3. UncertaintyquantificationandSensitivityAnalysisgiven uncertainradii

Toevaluatetheimpactofgeometricmeasurementerrors,we performedmetamodelbasedsensitivityanalysisforanLADsteno- sissubject touncertainmeasurementsof outletdiameters. The sensitivityofvFFRtoeachoutletisshowninFig.5,whereallsen- sitivitieshavebeengroupedintofivebinsaccordingtosensitivity fromverylowsensitivity(thelowest20%)toveryhigh(thehigh- est20%).TheabsolutesensitivitiesarealsoshowninFig.6with respecttothediameteroftheoutletandpositionrelativetothe stenosisintheLAD.Thelargerthearterythemoreinfluentialit isonthepredictionofvFFR,whileitalsoseemstherelationship betweendiameterandsensitivityissteeperforvesselsinthesame coronarybranchasthestenosis.Finally,thevesselsdistaltothe stenosiswereonaveragealsomarkedlymoreinfluentialonvFFR thanthoseproximal.

3.4. Effectofunobservedvesselsduetolimitedresolution

The impact of imaging limitations of clinical imaging was assessedbysimulatingreducedresolutiondataasdescribedinSec- tions2.2and 2.6.4.Foravisibilitythresholdofd_th=1.5mmand

Fig.5. CoronarytreewithastenosisintheupperLADandoutletsgroupedafterthe relativeimportanceaccordingtothetotalSobolindicescalculatedusingameta- modelandMonteCarloestimation.Eachcircleindicatesanoutletofthenetwork, andthelargeranddarkerthecircleindicatevFFRismoresensitivetouncertaintyin thediameterofthatparticularoutlet.Thesmallestcirclesindicatethesensitivityto thoseoutletswasamongthebottom20%ofalloutlets,thelargestcirclesindicates thoseoutletsrankedinthetop20%ofoutlets.Thesizesinbetweenthuscorrespond tovesselsfallingbetweenthe20thand40th,the40thand60th,andthe60thand 80thpercentilesrespectively.

Fig.6. TotalSobolsensitivityindicesofvFFRtouncertaintyintheradiiofoutletsinthecoronarytreewithastenosisintheupperLAD.Sensitivityindicesareplottedaccording toreportedvesseldiameterfrom[15].

Fig.7. TheRMSEofvFFRrelativetofullvisibilityisshownfortheconventional approach(dashedorange)andforthecorrectivemodel(solidblue)atdifferentvalues ofthevisibilitythreshold,dth.

stenosesof75% area, theRMSE(19)was0.131.TheRMSEwas furthercharacterizedoverarangeofd_thasshowninFig.7.

Theleakagemodeltoaccountforunobservedoutlets(seeSec- tion2.3)wasassessedonthesamecasessimulatingclinicalimaging data. The RMSE of vFFR_leak calculated with this correction for d_th=1.5mmwas0.067.TheRMSEwasalsocalculatedoverarange ofvaluesofd_thasshowninFig.7.

Inthis analysisthevFFRvalues predictedbyvFFR(Nd_th)and vFFR_leak(Nd_th)werecomparedtothereferencecasevFFR(N).The 1.5mmthresholdismotivatedbyclinicallimitations,butitisnota precisecutoffbelowwhichnovesselsareobserved.Consequently, thedependenceoftheerrorinvFFRontheexactvalueofthethresh- oldisofinterest.Toaddressthis,thesameanalysiswasappliedat numerousvaluesofd_thfrom0mmto2mm,correspondingtofull visibilityandnooutletsvisible,respectively.TheresultingRMSE ofvFFR(Nd_th)andvFFR_leak(Nd_th)areshowninFig.7.Furtherthe differencebetweentheconventionalapproachandthecorrective, leakyvesselapproachwascharacterizedonapervesselbasisby comparingtheerrorofthetwoapproachesforstenoseslocatedin eachparticularvessel(seeFig.8).

4. Discussion

Thisstudyconsideredthreecategoriesofuncertaintiesrelevant tosimulationsofcoronarybloodflow:parametricuncertaintyin boundaryconditions,geometricuncertaintyofbloodvesselsinthe domain,andtopologicaluncertaintyduetolimitedimagingreso-

lution.Thefocusofthisanalysiswastoquantifyhowuncertainties abouttheinflowsandoutflowsofthewholenetworkaswellasthe uncertaintiesaboutunderlyinghemodynamicalparametersmay affectthepredictionofvFFR ingeneral.Assuch, weperformed sensitivityanalysisonasimplifiedmodelofcoronarybloodflow thatallowedalargenumberofsimulationstobeperformed,while alsoretainingdependenceonsimilarinputdataandassumptions asmoresophisticatedmodels.Someotherstudiesofuncertainty quantificationofvFFRhavebeenperformedbySankaranetal.[43]

andEcketal.[21].Thesestudieshaveconsistentlyfoundthecross- sectionofthestenosisasthemostimportantparameterasopposed torheologicalparametersorstenosislength.Sankaranetal[43]

alsoconsidereduncertaintyabouttheoutletresistanceandfound ittohavesimilarthoughslightlysmallerinfluenceonvFFRincom- parisonto thecross-section of thestenosis; however,they did notaccountforuncertaintiesintheinflowboundarycondition.As previousstudiesprimarilyfocusedonthecharacterizationofthe stenosisitself,andlessonthedeterminationofboundarycondi- tionsforpredictionofvFFR,ourresultscontributeamissingpiece ofanalysistothefield.

Theuncertaintiesofboundarymodelparameterswerebasedon stateoftheartclinicalmeasurementmodalitieswhereapplicable, andliteraturesurveyofpopulationandstudy-to-studyvariability otherwise.Subjecttotheseuncertaintiesthemodel’spredictionsof vFFRhadsignificantuncertaintiesandwereunabletoconfidently determineifthe75%areareductionisclinicallysignificantaccord- ingtothecutoffvalueof0.75identifiedbythestudyofToninoetal.

[7]asthemodel’s95%predictionintervalspannedfrom0.59to0.90.

Sensitivityanalysisofthiscase indicatedthatthemeasurement ofcardiacoutput,which determinesmyocardialinflow,andthe degreeofhyperemiainducedinfluencedvFFRtothegreatestextent oftheallparametersconsideredinouranalysisofboundarycondi- tionparameteruncertainty(seeFig.3).SincethesensitivityofvFFR tothemyocardialflowwasquitelarge,itwouldbeadvantageous toimprovedeterminationofthetotalinflowtothemyocardium, andifpossibletobettercharacterizethedeterminantsofadeno- sine’sefficacyatinducinghyperemia.Perhapsimprovedimaging modalitiessuchasPET-MRIwillbecomeclinicallyfeasibleinthe future[59],allowingdirectmeasurementofcoronarybloodflow;

however,theyarecurrentlynotwidelyusedduetotheincreased radiationexposureandexpense.Basicinvestigationofmyocardial responsetoadenosinecouldelucidateprocedurestobetterchar- acterizeindividualpatient’scoronaryvasodilation.

ThevalueofvFFRpredictedwasmuchlesssensitivetounder- lying rheological parameters or the exponent c from Murray’s law,even thoughthis exponentisquitevariedintheliterature.

Since the sensitivity of vFFR to these parameters is quite low,

Fig.8. DifferenceintheerroroftheconventionalandcorrectedvFFRestimate,

vFFR(N_dth,,i)−vFFR(N,,i)

−

vFFRleak(N_dth,,i)−vFFR(N,,i)

,shownonpervessel basisvsvisibilitythreshold.RedindicatestheconventionalestimateofvFFRisworsethanvFFRleakbasedoncorrectionforinvisiblevessels.Bluehatchedregionsindicate theconventionalvFFRwasmoreaccuratethanvFFRleak.

effortsforimprovingtheconfidenceofvFFRshouldnotfocuson characterizingthesephysicalmodelparametersmoreaccurately, ratherthefocusshouldbeonprovidingaccurateinflowconditions anddeterminingtheperipheralvasculature’sresponsetoinduced hyperemia.

OuranalysisofthefitofMurray’slawtothegeometricdatafor theporcinearterysuggestedthatthemismatchbetweenMurray’s lawandthereportedgeometrymaybeexplainedbyerrorsinmea- surementofthediametersofvessels,whiletheuncertaintyabout theestimatedexponentcisquitelowrelativetotherangeofvalues foundintheliterature.Themeanvalueof2.47lieswithin,buton thelowerendof,therangereportedbyvariouspreviousanalyses [17,60–63].

However, as our analysis shows c only influences vFFR marginally(seeFig.3),themoreinterestingresultsisthatMur- ray’slaw suggestsdiscrepancy betweentheobserved values of vesseldiametersand thevaluesinphysiologicalconditions.The discrepancymaybeduetodeviationsfromphysiologicalcondi- tionsinducedbythecreationoftheelastomercastasnotedby Suwaetal.[45],orsimplyduetomeasurementerrors.Inanycase thefindingsaboveemphasizetheimportanceofaccuratemeasure- mentofthecross-sectionaldiameterofthevasculatureasthiscan greatlyinfluencetheflowpredictedthroughagivenvessel.

Giventhelikelyuncertaintyinoutletdiametersevaluationof oursensitivityanalysisofvFFRtosuchmeasurementsisparticu- larlyrelevant.ThisanalysisofvFFRwithrespecttooutletdiameters showstwotrends.Firstuncertaintyinlargervesselsismoreinflu- entialrelativetosmallervessels(seeFig.6).Secondtheinfluence ofvesselsdistaltothestenosisaremarkedlymoreinfluentialthan vesselsof similardiameter proximaltothestenosis orinother branchesofthenetwork.Theseresultsemphasizetheimportance ofidentifyingandsegmentinglargervesselsaccurately,andthat measurementandcharacterizationofthevasculaturedistaltothe stenosisseemstobeoffargreaterimportancethancharacterizing theproximalvasculatureorthevasculatureinotherbranchesof thecoronarynetwork.

Analysisofgeometricuncertaintynaturallyleadstoconsidera- tionofdifferentlayersofuncertaintyingeometricdata,i.e.physical limitationsofimaging,segmentationprocessing,andsmoothingof computationaldomains.Identifyinghowmuch uncertaintyeach layercontributeswouldbeavaluablecontribution,albeitdifficult.

Thisisparticularlytrueregardingimagesegmentationandgenera-

tionofcomputationaldomains,astherearenotgenerallyabsolute measurementsavailableforestimationoftheerror.Toourknowl- edgenostudieshavebeenconductedthatsystematicallydetermine thecontributionsofeachofthesestagesofimageprocessingfor computationalsimulations.Somestudieshaveusedinterobserver variabilityasametric forcomparingperformanceofautomated algorithms[64],andBulantetal.[18]comparedthevFFRresult- ingfromCTandIVUSdatainthesamepatients.Theyfoundthe resultswerequitesimilarforbothimagingmodalitiessuggesting thatuncertaintiesduetophysicallimitationsofimagingmodalities arenotdrasticallyaffectingtheresultingsimulations.

Since theresolution of clinicalimaging ofcoronary vascula- tureislimited,clinicallyreliablemethodsofestimatingvFFRmust accountfor theeffectof this thresholdtoensuretheestimates arerobusttovariationsinimagequality.Imagingmodalitiesand subsequentprocessing,suchassegmentation,areencumberedby limitationsofthelevelofdetailthatmayberesolved.Thus,many vesselsbelowacertaindiameterwillnotberepresentedinthesim- ulations.Toevaluatetheeffectofthislimitation,wecomparedthe resultsofourestimationofvFFRusingthehighresolutiondatafrom theelastomercasttothepredictedvFFRbasedonlowerresolution dataconsistingofonlyvesselsofdiameterlargerthanathreshold d_thbetween0mmand2mm(seeFig.7).Theresultsdemonstrate thatlackofdetaileddiameterobservationsresultsinanRMSEthat substantiallyincreasesford_th>1mmandattainsavalueof0.131if d_th=1.5mm,aclinicallyrelevantcutoffvalue.Theseresultssuggest thattheimpactsofunobservedoutletsmaybesubstantialformod- elsthatdistributeflowaccordingtoobservedoutletsandefforts shouldbemadetomitigatetheimpact.

Weproposedonesuchmethod(seeSection2.3)whichapplies theprincipleofMurray’sLaw(3)toaccountforinvisiblevessels.

Ourmethodestimatesanexpectedleakageflowbasedonthedis- agreementbetweentheimagedvesselsandMurray’slaw(12).We appliedourmethodtothesamecasesconsideredforevaluatingthe effectsoflimitedresolution.TheRMSEoftheproposedmethodis showninFig.7atvariousvaluesofd_thbetween0mmand2mm.

Thecorrectivemethodperformsmarkedlybetterifd_thissomewhat largerthan1mm.Ford_th=1.5mmtheRMSEwas0.067or48%lower thanfortheuncorrectedapproach.Comparingtheconventional andcorrectiveapproachonapervesselbasisshowsanincreased benefitoftheleakyvesselapproachforstenoseslocatedinthedistal regionsofthecoronaryvasculature(seeFig.8).Thispatterncould