Automated optimization of a reduced layer 5 pyramidal cell model based on experimental data

ContentslistsavailableatSciVerseScienceDirect

Journal of Neuroscience Methods

jo u rn al h om 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

Automated optimization of a reduced layer 5 pyramidal cell model based on experimental data

Armin Bahl

, Martin B. Stemmler

, Andreas V.M. Herz

, Arnd Roth

aDepartmentofSystemsandComputationalNeurobiology,MaxPlanckInstituteofNeurobiology,82152Martinsried,Germany

bGraduateSchoolofSystemicNeurosciences,Ludwig-Maximilians-UniversitätMünchen,82152Martinsried,Germany

cDepartmentofBiologyandBernsteinCenterforComputationalNeuroscienceMünchen,Ludwig-Maximilians-UniversitätMünchen,82152Martinsried,Germany

dWolfsonInstituteforBiomedicalResearchandDepartmentofNeuroscience,PhysiologyandPharmacology,UniversityCollegeLondon,LondonWC1E6BT,UK

a r t i c l e i n f o

Articlehistory:

Received14November2011 Receivedinrevisedform4April2012 Accepted5April2012

Keywords:

Pyramidalneuron Compartmentalmodel Evolutionaryalgorithm Multi-objectiveoptimization Automatedfitting

Dendriticgeometry Firingpattern

Dendriticcalciumdynamics

a b s t r a c t

Theconstructionofcompartmentalmodelsofneuronsinvolvestuningasetofparameterstomakethe modelneuronbehaveasrealisticallyaspossible.Whiletheparameterspaceofsingle-compartment modelsorothersimplemodelscanbeexhaustivelysearched,theintroductionofdendriticgeometry causesthenumberofparameterstoballoon.Asparametertuningisadauntingandtime-consuming taskwhenperformedmanually,reliablemethodsforautomaticallyoptimizingcompartmentalmodels aredesperatelyneeded,asonlyoptimizedmodelscancapturethebehaviorofrealneurons.Herewe presentathree-stepstrategytoautomaticallybuildreducedmodelsoflayer5pyramidalneuronsthat closelyreproduceexperimentaldata.First,wereducethepatternofdendriticbranchesofadetailed modeltoasetofequivalentprimarydendrites.Second,theionchanneldensitiesareestimatedusinga multi-objectiveoptimizationstrategytofitthevoltagetracerecordedundertwoconditions–withand withouttheapicaldendriteoccludedbypinching.Finally,wetunedendriticcalciumchannelparameters tomodeltheinitiationofdendriticcalciumspikesandthecouplingbetweensomaanddendrite.More generally,thisnewmethodcanbeappliedtoconstructfamiliesofmodelsofdifferentneurontypes,with applicationsrangingfromthestudyofinformationprocessinginsingleneuronstorealisticsimulations oflarge-scalenetworkdynamics.

© 2012 Elsevier B.V.

1. Introduction

Toincorporaterealismintolarge-scalesimulationsofcortical andothernetworks(Traubetal.,2005;Markram,2006),oneneeds toconstructbiophysicallyrealisticcompartmentalmodelsofthe individualneuronsinthecircuit.Manyparametersofthesemodels havenotbeendirectlymeasuredexperimentally;therefore,these parametersmustbetunedtomatchtheexperimentallyobserved input–outputrelationoftheneuron.Solvingtheresultingnonlinear optimizationproblemisdifficultandrequiresextensivecomput- ingresources,especiallyformodelscomprisingadetailedneuronal morphologyanda largenumber ofcompartments(Traubetal., 2005;AchardandDeSchutter,2006;Markram,2006;Druckmann etal.,2007;Hayetal.,2011).

Herewedevelopreducedmodelsofneocorticallayer5pyra- midalcells withasmall numberof compartmentstorepresent

∗Correspondingauthorat:DepartmentofSystemsandComputationalNeurobi- ology,MaxPlanckInstituteofNeurobiology,82152Martinsried,Germany.

Tel.:+498985783289.

E-mailaddress:arbahl@gmail.com(A.Bahl).

thedendriticgeometry.Comparedtofullydetailedcompartmental models,reducedmodelsconfersignificantspeedadvantagesboth fortheoptimizationofthesingleneuronmodelandforsimulation ofnetworksofsuchneurons.Thereducedmodel’sgeometry,even thoughsimplified,shouldstill incorporatethefactthatsynaptic inputsarrivingatdifferentlayersinthedendritictreeareintegrated differently(Larkumetal.,1999,2004;Schaeferetal.,2003;Branco etal.,2010).Thisimpliesthatasufficientnumberofcompartments mustbeusedforthereduceddendriticmorphology.Theabilityofa neurontolockontofastfluctuationsintheinputdependscritically onhowsharptheactionpotentialonsetisinitsvoltagetime-course (Naundorfetal.,2005;PalmerandStuart,2006;Koleetal.,2007;

Popovicetal.,2011).Astheinitiationsiteoftheactionpotential, locatedintheaxon(StuartandSakmann,1994;PalmerandStuart, 2006;Koleetal.,2007;Popovicetal.,2011),determineshowsteep theactionpotentialonsetis(Yuetal.,2008),areducedmodelmust alsohaveaminimumnumberofcompartmentsfortheaxon.

Weconstructthereducedpyramidalcellmodelsstepbystep, applyingmethodsadaptedtotheproblem(RothandBahl,2009).

Thefirststep,whichdeterminesthereducedcellmorphologyand passivemembraneproperties,followsthetraditionofdefiningcer- tainpassiveelectricalpropertiesofthefullmodelthatarepreserved

0165-0270© 2012 Elsevier B.V.

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

Open access under CC BY license.

Open access under CC BY license.

Konstanzer Online-Publikations-System (KOPS) URL: http://nbn-resolving.de/urn:nbn:de:bsz:352-2-49mjifjo9us20

inthereducedmodel(Stratfordetal.,1989;BushandSejnowski, 1993;Destexhe,2001).Second,weapplyapowerfuloptimization strategy, EvolutionaryMulti-Objective Optimization(EMOO) (Deb etal.,2002;Druckmannetal.,2007).Thismethod,startingfrom afamilyofmodelscharacterizedbymultiplefeatures,generatesa suiteofnewmodelsateachstep,withoutmakinganapriorideter- minationastowhichoneofthemultipledesiredobjectivesismost important.

Inpreviousapproachestoevolutionarymodeloptimization,up totwentydifferentfeaturesaredefined(Hayetal.,2011),fromthe actionpotentialwidthtothedepthoftheafter-hyperpolarization potential(AHP).Eachfeatureisassociatedwithasinglenumber.

In contrast, we return tothe classical approach of minimizing theleastsquaredifferencebetweentherecorded traceandthe model’sresponse.Wedefinefourleastsquaredifferenceobjective functions,oneforeachofthedifferenttimescalespresentinthe dynamics– fromtheAPonsetdynamicsonthemicrosecondtime scaletotheslowsub-thresholdchargingphaselastingseveraltens ofmilliseconds.

Theexperimentaldatawefitincludetheresponsesofalayer 5pyramidalneurontosomaticcurrentinjectionwhentheapical dendriteisoccludedorpinched(BekkersandHäusser,2007).Such aprocedureyieldsessentialinformationabouttheelectricalprop- ertiesofthecellanditsapicaldendrite.Thisinformationallows parameteroptimizationtonarrowdownthepossiblecombinations ofchanneldensitiesandpropertiesalongthedendriticgeometry thatcouldexplaintheneuron’svoltageresponsetosomaticcurrent injection.

Finally, we show how the evolutionary approach can be extendedtoensurethataneuronalmodelcapturesfeaturesthat donottakeoncontinuousvalues,butarediscrete.Forthispur- pose,we tookanotherdataset fromadifferent experimenton olderlayer5pyramidalneurons.Weoptimizedthedendriticcal- ciumchannelparametersandadjustedthesevaluessuchthatthe modelreproducestheshapeofthedendriticcalciumAPaswellas somato-dendriticcouplingfactorsfoundinexperiments(Schaefer etal.,2003).

Afterpursuingthesethreestepsinoptimizingneuronalmod- els,wepresentafamilyof10reducedmodelsoflayer5pyramidal neuronswhoseinput–outputrelationmatchesarangeofexperi- mentaldata.Thesemodelscouldbeusedinlarge-scalenetwork simulationsoftheneocortex.

2. Methods 2.1. Thecellmodel

Ouraimistocreateareducedpyramidalcellmodelthatissim- pleandfastbutdetailedenoughtoshowcomplexsomato-dendritic interactions.Themodelisbasedonstandardtechniquesfromcom- partmentalandionchannelmodeling(HodgkinandHuxley,1952;

Rall,1962), implementedin NEURON7.1(CarnevaleandHines, 2005)andcontrolledviatheNEURON-Pythoninterface(Hinesetal., 2009).

Toobtaina simplified geometryof thedendrites, we model the functional neuronal sections (soma, basal dendrites, apical dendriteand theapicaldendritictuft) eachbyasinglecylinder whoselengthanddiameterwillbelaterdeterminedbytheopti- mizationalgorithm we describe. Theaxonal geometry is based ona detailedreconstruction(Zhu,2000)and consistsofa coni- calaxonhillock(l=20␮m)whichhasadiameterof3.5␮matthe somaconnectionand tapersto2.0␮m. Theconical axoninitial segment(iseg;l=25␮m)isconnectedtothehillockanditsdiam- etertapersfrom2.0␮mto1.5␮m.Theactualaxon(l=500␮m)is connectedtotheinitialsegmentandhasa uniformdiameterof

1.5␮m.WedidnotmodelnodesofRanvierormyelination.Asthe reducedmodelshouldbefast,thenumberofcompartmentsought tobeassmallaspossible.Wechosethefollowingcompartment numbersforthefunctionalsections:soma=1;basaldendrite=1;

apicaldendrite=5;apicaldendritictuft=2;axonhillock=5;initial segment=5;axon=1.Hence,themodelhasatotalof20compart- ments.

Ion channelswere selected and distributed based onrecent experimentalfindingsandmodelingstudiesandweredownloaded fromModelDB(Hinesetal.,2004):Ahyperpolarization-activated cationchannel(HCN)(Koleetal.,2006)wasinsertedintothebasal dendrite,the apicaldendrite and thedendritic tuft.A transient sodiumchannel(Nat)(Koleetal.,2006)wasplacedintothesoma, theaxonhillock,theinitialsegment,theapicaldendriteandthe dendritictuft.Thevoltagedependencyofthechannelkineticswas shiftedtohighervalues(vshift_Nat=+10mV)inallcompartmentsto yieldhigherthresholds(MainenandSejnowski,1996).Weintro- ducedasecondvoltageshift(vshift2Nat)fortheNatchannelinthe initialsegmenttoaccountfordifferentchannelpropertiesinthis area(ColbertandPan,2002).Natchanneldensitydecayedlinearly intheapicaldendritewithdistancefromthesoma(Mainenetal., 1995;Kerenet al.,2009).Afastpotassium channel(Kfast)(Kole etal.,2006)wasinsertedintothesoma,theapicaldendriteand thetuft.Itsdensitydecayedexponentiallyfromthesomatowards thetuft(Kerenetal.,2009).Aslowpotassiumchannel(Kslow)was inserted intothe soma, the apicaldendrite and tuftand chan- neldensitiesdecayedexponentiallywithdistancefromthesoma (KorngreenandSakmann,2000).Apersistentsodiumchannel(Nap) wasinsertedintothesomatoadjusttheneuron’sexcitability(Traub etal.,2003).Amuscarinicpotassiumchannel(Km)wasinsertedinto thesoma.TheKmchannelisanon-inactivatingvoltage-dependent slowpotassiumchannelwhichisthoughttoplayaroleinspike frequencyadaption(Winogradetal.,2008).Aslowcalciumchan- nel(Cas)wasinsertedintothetuft.The voltagedependencyof itskineticscouldbeshifted(vshift_Cas)toadjustactivationthresh- olds.Finally,acalciumdependentpotassiumchannel(KCa)(Mainen andSejnowski,1996)aswellasacalciumpump(CP)(Koleetal., 2006)wasinsertedintothetuft.The10optimizedreducedmodel neuronsandtheionchannelmodelsareavailablefordownloadat http://senselab.med.yale.edu/ModelDB.

2.2. Freemodelparameters

Reliabledataonthebiophysicalpropertiesofionchannelsin pyramidalneuronsarerareandmeasurementsmostlystemfrom differentneuronsfromdifferentanimalsoreven species.More- over,duetoexperimentallimitationsmanyparametersjustcannot bemeasured.Inparticular,informationaboutionchanneldensi- tiesandkineticsindistaldendriticbranchesarenotavailable.It is,therefore,notsufficienttoputtogetherallexistinginformation onpyramidalneuronsandbuildaworkingmodel;alonglistof uncertaintiesremains.Wemadeaselectionofthemostuncertain parametersthatwethoughtcouldbeestimatedbymeansofan optimizationstrategy:Thelengths(l),thediameters(d)andthe axialresistances(Ra)ofthefunctionalsectionsofthereducedmor- phologywerefreeparameters,aswerethevaluesofthemembrane resistance(Rm)andcapacitance(Cm).Weintroducedadendritic scalingfactor(dendscaling)thatallowedtheadjustmentofden- driticmembraneresistanceandcapacitancerelativetothesoma.

Thisshouldaccountfordendriticspines,systematicerrorsinden- driticreconstructionsordifferentratiosbetweenthesomasizeand thedendritictreesizeindifferentcells. Furtherfreeparameters weretheionchanneldensitiesindifferentcompartments,which weredescribedbythechannel’smaximalionicconductanceper membranearea(e.g.soma ¯g_Nat).Itwasshownthatcertainionchan- neldensitiesalongtheapicaldendritecanbedescribedbysimple

Table1

Fulllistofalloptimizedparameterswiththeirloweranduppersearchbounds(LSB,USB)togetherwithasummaryofmodelpropertiesforall10reducedmodels.(a)The morphologicalparametersfoundbyreducingthedetailedmodel(firststep)areshown.Theremainingparametersaregivenbecausethemethodmaintainsthetotalsurface areaofthefunctionalsectionsduringtheoptimization(Asoma=1682␮m²,Abasal=7060␮m²,Aapical=9312␮m²andAtuft=9434␮m²andtherefored^soma=23␮m,L^soma=23␮m, d^basal=8.7␮m,d^apical=5.9␮m,d^tuft=6.0␮m).(b)Themodelparametersareshownthatwerefoundaftertheoptimizationofionicconductances(secondstep).Models1–6 wereoptimizedusingtargetdatabeforeandduringpinching,whilemodels7–10wereoptimizedonlywithtargetdatafromtheintactneuron.(c)Modelparametersfound afteroptimizingthedendriticcalciumdynamics(thirdstep)areshown.Inmodels4,5,9and10noparameterscouldbyfoundsatisfyingtherequirementsforpropercalcium dynamics.(d)Wesummarizedafewpropertiesofthe10resultingmodels:theratioofthesomaticsodiumandpotassiumchanneldensity,theratioofthesodiumchannel densityintheinitialsegmentandsomaaswellastheratioofthetuftHCNandbasalHCNchanneldensity.Itisalsosummarizedthatallmodelsshowback-propagatingAPs andwhichmodelsshowadecayofBAPamplitudewithsomaticdistance(illustratedinFig.7).Finallywesummarizethemodulationoftherestingpotentialalongtheapical dendrite(distalvoltage–somaticvoltageatrest,seeFig.6).Ifpropercalciumdynamicscouldbefoundinthethirdstepweshowtheresultingcouplingfactor.

Parameter Result LSB USB Unit

(a)Parametersfoundinfirststep

somaRa 82 80 200 cm

basalL 257 170 280 ␮m

basalRa 734 700 2000 cm

apicalL 500 500 800 ␮m

apicalRa 261 150 300 cm

tuftL 499 400 600 ␮m

tuftRa 527 500 1200 cm

Parameter Model1 Model2 Model3 Model4 Model5 Model6 Model7 Model8 Model9 Model10 LSB USB Unit (b)Parametersfoundinsecondstep

epas −83.06 −80.40 −80.50 −78.97 −82.55 −85.00 −83.68 −80.74 −84.37 −80.75 −85.00 −60.00 mV

Rm 23,823 20,588 20,514 10,784 17,387 15,159 21,298 11,594 11,081 14,712 10,000 30,000 cm²

Cm 2.30 2.23 2.41 1.79 2.02 2.71 2.37 2.34 1.40 2.51 0.60 3.00 ␮F/cm²

dendscaling 0.86 0.78 0.69 0.50 0.52 0.50 0.55 0.50 0.50 0.50 0.50 2.00 1

soma ¯gNat 284.55 236.62 238.88 182.17 248.75 295.47 447.25 402.17 277.35 371.43 0.00 500.00 pS/␮m²

soma ¯gKfast 50.80 67.20 59.26 45.43 44.77 43.23 43.78 41.35 32.41 48.25 0.00 300.00 pS/␮m²

soma ¯gKslow 361.58 475.82 433.80 467.50 523.27 630.25 190.33 264.79 187.87 621.74 0.00 1000.00 pS/␮m²

soma ¯gNap 0.87 1.44 1.48 3.33 2.29 3.52 0.85 4.18 2.21 4.16 0.00 5.00 pS/␮m²

soma ¯gKm 7.12 10.46 11.12 12.99 14.20 11.91 11.05 14.92 12.22 7.00 0.00 15.00 pS/␮m²

basal ¯gHCN 15.71 11.04 10.72 7.94 13.62 12.87 3.12 11.92 13.90 22.09 0.00 50.00 pS/␮m²

tuft ¯gHCN 17.69 16.19 17.80 18.89 15.73 23.93 40.67 15.27 51.84 3.34 0.00 150.00 pS/␮m²

tuft ¯gNat 6.56 47.82 29.01 76.65 40.38 45.92 0.41 11.56 46.94 87.60 0.00 100.00 pS/␮m²

Kfast 58.52 20.08 55.58 2.15 8.91 65.61 82.07 91.80 73.47 67.52 1.00 100.00 ␮m

_Kslow 42.21 37.71 88.72 55.49 49.67 34.33 65.18 75.61 69.61 83.03 1.00 100.00 ␮m

hillock ¯gNat 8811 9512 8303 5997 4988 9451 8171 8030 4904 8407 0 20,000 pS/␮m²

iseg ¯gNat 13,490 13,327 17,624 12,625 10,730 17,194 19,583 17,591 10,777 15,509 0 20,000 pS/␮m² isegvshift2Nat −9.80 −10.61 −9.57 −10.16 −10.75 −8.92 −5.36 −5.98 −10.26 −8.47 −15.00 0.00 mV

Parameter Model1 Model2 Model3 Model4 Model5 Model6 Model7 Model8 Model9 Model10 LSB USB Unit (c)Parametersfoundinthirdstep

apicalRa 454.06 382.22 444.13 – – 445.01 332.92 358.47 – – 250.00 500.00 cm

tuft ¯gCas 3.68 0.45 2.12 – – 0.49 2.81 3.86 – – 0.00 4.00 pS/␮m²

tuftvshiftCas 7.48 7.19 8.35 – – 0.80 2.35 3.79 – – −10.00 10.00 mV

tuft ¯gKCa 9.76 6.15 8.23 – – 9.69 9.55 9.60 – – 0.00 4.00 pS/␮m²

Property Model1 Model2 Model3 Model4 Model5 Model6 Model7 Model8 Model9 Model10 (d)Modelproperties

soma ¯gNat/( ¯gKfast+g¯Kslow) 0.69 0.44 0.48 0.36 0.44 0.44 1.91 1.31 1.26 0.55 iseg ¯gNat/(soma ¯gNat) 47.41 56.32 73.78 69.30 43.14 58.19 43.78 43.74 38.86 41.76 tuft ¯gHCN/(basal ¯gHCN) 1.13 1.47 1.66 2.38 1.16 1.86 13.03 1.28 3.73 0.15

BAPs? Yes Yes Yes Yes Yes Yes Yes Yes Yes Yes

BAPsdecay? Yes Yes Yes No No Yes Yes Yes No No

Mod.ofrest +1.2mV +1.5mV +2mV +1.5mV +1.5mV +2mV +5mV +0.5mV +4mV −2mV

Coupling 0.50 0.50 0.50 – – 0.62 0.50 0.54 – –

functions(Bergeretal.,2001;Koleetal.,2006).Weassumedan exponentialdecayofthepotassiumchanneldensityalongtheapi- caldendritebutallowedthespaceconstant(_Kfast,_Kslow)tobe afree parameter.Natand HCNchanneldensitieswerechanged linearlyalongtheapicaldendritewithslopescalculatedbasedon thechannels’densitiesinthesomaandtuft.Finally,somekinetic parametersofthesameionchanneltypearedistinctindifferent partsoftheneuron(ColbertandPan,2002).Thereforewechosethe voltagedependencyoftheNatchannelintheaxoninitialsegment (isegvshift2Nat)andthevoltagedependencyoftheCas-channelin thetuft(tuftvshiftCas)asfurtherfreeparameters.Thefreemodel parametersusedforoptimizationaswellasloweranduppersearch boundsarelistedinTable1.

2.3. Theoptimizationalgorithm

InallthreeoptimizationstepsweusedtheEvolutionaryMulti- ObjectiveOptimization(EMOO)algorithm(Deb,2001;Debetal., 2002).TheEMOO-algorithmallowsone tosimultaneouslymin- imize multiple and possibly conflicting error functions and is thereforeespeciallywellsuitedfortheoptimizationofspikingneu- ronmodelstoexperimentaldata(Druckmannetal.,2007,2008).

EMOO isa geneticalgorithm and usesmechanismsinspiredby biological evolution, suchas selection, crossoverand mutation, to grow a population of size N to a certain capacity C and to transfergoodindividualsintothenextgeneration.WeusedaSim- ulated BinaryCrossoveroperator (DebandAgrawal,1995)anda

PolynomialMutationoperator(Deb,2001).Theefficacyofanoper- atoriscontrolledbycandmrespectively.Largevaluescandm

meanthattheoperator’seffectisweak.TheEMOO-algorithmwas implementedinPythonandrunparallelizedonanAMDx8664 clusterwith80processorsrunningLinux.ThePython-Framework fortheevolutionarymulti-objectiveoptimizationisavailablefor downloadunderwww.g-node.org/emoo.

3. Results

3.1. Stepwisefittingstrategy

Animportantaspectofagoodmodelisthatitcloselyreproduces datausedtoconstructthemodeland,evenmoreimportantly,pre- dictskeyfeaturesofnewdatathatwerenotusedtoconstrainthe model(Druckmannetal.,2011).Tuningasetofparametersbyhand suchthatthemodelfulfillstheserequirementsisadaunting,time- consumingorevenimpossibletaskandanewsetofdataorminor modelmodifications mightrequire arepetition ofthat process.

Thereforetheprocessofcreatingamodelandtuningitsparam- etersshouldbeasautomatedaspossible.Todividetheparameter spaceintosmallerunitsandtooptimizeeachparametersubset independentlywehavedevelopedathree-stepfittingstrategy:In thefirststepweestimatedthegeometryofthemodel,includingthe axialresistancesofthesections.Thegeometricalparametersofthe modelwerefixedafterthisstep.Inthesecondstepallionchan- nelandmembraneparametersaffectingsomaticspikingbutnot calciumspikingdynamicswereoptimized.Oncethisstepwascar- riedout,theseparameterswerealsofixed.Finally,inthethirdstep weestimatedtheparametersneededfordendriticcalciumspike dynamics.Inanoptimalscenarioallthesedatawouldcomefrom thesameneuronfromthesameanimal.Thiswouldbeanexper- imentalchallengeandcurrentlysuchasetofdatadoesnotexist.

Thereforewehavetocombinedatafromdifferentexperimentsin eachofthesesteps.

3.2. Firststep:optimizingthereducedgeometry

Toobtainarepresentativegeometryforthereducedmodelwe startedwithadetailedmodelreconstructionofalayer5pyramidal neuronfromayoung(≈P21)rat(StuartandSpruston,1998).Todo soweadoptedamodelsimplificationstrategy(Destexhe,2001):

Asweare onlyoptimizing theneuronalgeometry weremoved allactiveconductancesandgloballysetthemembraneresistance (Rm)to15,000cm²,themembranecapacitance(Cm)to1pF/␮m², thereversalpotential(epas)to−70mVinboththedetailedandin thereducedmodel.Thespecificmembrane parameters(R_m,C_m, epas)donotchangewithgeometryandthereforedonotneedto beoptimized in the firststep. However theywill needfurther tuningwhen themodel is matched tospikingdata inthe sec- ondstep. Next, weassigned each compartment in the detailed model to one of the functional sections (basal dendrite,soma, apicaldendrite,tuft),determinedtheirsurfaceareasandsetthe sizeofthefunctionalsectionsinthereducedmodeltothesame values(Asoma=1682␮m²,A_basal=7060␮m²,A_apical=9312␮m²and A_tuft=9434␮m²). Then, in thedetailed model, we setthe axial resistance(Ra)globallytothecommonlyusedvalueof100cm.

Asitisnotclearwhattheaxialresistanceandthelengthofthe functionalsectionsinthereducedmodelshouldbe,weneededto estimatethesevalues.Onceweknowthelengthofafunctional section,wecanalsocomputeitsdiameter,asthisshouldleadto thesamesurfaceareaasmeasuredinthereconstruction.Destexhe (2001)onlyusedthesteadystatevoltageinresponsetoconstant currentinjectionstofittheseparametersacrossallcompartments.

Weusedthistargetfunction,aswell,butextendedthemethodto

reproduce the detailed neuron’s somatic input impedance and phase shift functions for oscillatory somatic input currents (0–1000Hz).Theseimpedancefunctionswereshowntocontain information about thesub-threshold neuronal dynamicsof the wholemorphology(Fox,1985;BorstandHaag,1996).Wedeter- mined these three functions in the detailed model and, for a givenparametercombination, in thereducedmodel and calcu- latedthesumofsquareddifferences.Thisgaveusthreefeatures tominimizebyEMOO.WeusedapopulationsizeofN=350and acapacity ofC=700individualsandevaluated100generations.

Thecrossoverparameterstartedatc=5andincreasedlinearlyto c=50,whilethemutationparameterincreasedfromm=10to m=500,therebyreducingthestrengthoftheseoperatorsduring evolution.Themutationprobabilityperparameterwasconstantat 10%.Wedroppedtheaxoninthisstepasthedetailedreconstruc- tionlacksanaxon.Theaxonwasappendedagainafterwardsand itsmembranepropertieswerechosentoequalthosefoundforthe soma.Appendingtheaxonreducedthesomaticinputresistance from69Mto63M.

After optimization the passive response properties of the reduced model matched those of the detailed model (Fig. 1).

Onesetofoptimalparametersisgivenin Table1;furtheropti- mization trials have led to different parameter combinations (notshown)that reproducedthepassiveresponsepropertiesof thedetailedmodelsimilarlywell.Tochallengetheoptimization methodweinjectedthesameGaussiannoiseintothesomaofthe reducedandcomplexmodelandmeasuredtheresultingvoltage responsesinthesomaandin thedendrite.Asdemonstrated in Fig.2,thepassivevoltageresponsesinboth modelsarealmost indistinguishable.

3.3. Secondstep:optimizingtheionchannelparameters

In the next step we estimated the ion channel parameters affectingsomaticspiking.Thiswasdoneusingexperimentaldata onthesomaticspikingdynamicsundertwodifferentconditions, firstin theintactneuron,and secondwhiletheapicaldendrite hasbeenoccludedusinga methodcalledPinching(Bekkersand Häusser,2007).Dendritic calciumspikesdevelop morefully in olderpyramidalneurons(Schilleretal.,1997)andareactivated bystrongdendriticlocaldepolarizationorbysomaticinputcom- binedwithdendriticinput(Larkumetal.,1999,2001).Thismeans that theexperimentaldatadoesnot containinformation about calciumdynamics.Thereforethisfeaturecouldnotbeoptimized hereandwe onlytunedtheparametersaffectingsomaticspik- inginthesecondstep.Astherecordingsweremadefromyoung rat(P17-25) pyramidalcells wecan usetheoptimized reduced geometrythatwe obtainedafterthefirststepas itisbased on adetailedgeometryofapyramidalneuronfromaratofsimilar age.

Pinchinghasanumberofeffectsonneuronaldynamics(com- pare black and blue traces in Fig. 3): (1) the input resistance increases (from ≈82M to 131M) and hence the spike fre- quency.By occluding theapical dendrite, one path for current lossis blocked,allowingthe somaticmembrane tobecharged moreeffectively(BekkersandHäusser,2007).(2)Thesomaticrest- ingpotentialbecomesmorehyperpolarized(≈4mV).Ithasbeen shownthatthedensityofI_hincreaseswithdistancefromthesoma intheapicaldendriteoflayer5pyramidalneurons(Bergeretal., 2001;Koleetal.,2006).Whentheapicaldendriteisblockedduring pinching,thedepolarizinginfluenceofdendriticHCNchannelsis reduced,whichmightexplainthehyperpolarizationofthesomatic restingpotential.(3)ThethresholdforAPinitiationbecomeslower (≈3mV)andtheAPpeakvoltageincreases(≈0.8mV),which,at leastinamodel(BekkersandHäusser,2007)couldbeexplainedby thedecreaseinthedendriticcapacitance,reducingtheelectrical

Fig.1. Morphologyandpassiveresponsepropertiesforthecomplexandthereducedcellmodel.(a)Thedetailedreconstructionofalayer5pyramidalneuron(Stuart andSpruston,1998)thatweusedasthestartingpointtocreatethereducedmodel.(b)Anillustrationofthereducedmodel(samescaleasdetailedmodel,geometrical parameterscanbefoundinTable1).Wedividedthecomplexmorphologyintofourfunctionalsections:Thesoma,thebasaldendrites,theapicaldendritesandthetuft.

Theobliquedendritesareconsideredtobepartoftheapicaldendrites.(c)Weinjectedaconstantcurrent(−1nA)intothesomataofbothmodelneuronsandmeasured thesteady-statevoltageatdifferentlocations.(d)Weusedalow-amplitudeoscillatorysomaticinputcurrentandmeasuredtheresultingmembranepotentialoscillationto determinethesomaticfrequency–impedancecurveforbothmodels.(e)Wecalculatedthesomaticphase-shiftbetweentheoscillatoryinputcurrentandresultingmembrane potentialoscillationforbothmodels.Theblackdotsandcurvesin(c–e)describethepassiveresponsepropertiesofthedetailedmodelandservedastargetfunctionsforthe optimizationprocedureinthefirststep.Thereddotsandcurvesshowthecorrespondingresponseofthereducedmodelafteroptimization.

Fig.2. Comparisonofthevoltagetracesinthecomplexandinthereducedmodelinresponsetonoisyinputcurrent.Totestwhetherthereducedmodelisagoodapproximation ofthecomplexmodel,weanalyzedresponsestoGaussiannoisecurrentinjections.Thesamerandomcurrentwasinjectedintothesomataofthemodels(greenelectrodes inaandbandgreentraceine).Forbothmodels,thesomaticvoltageaswellasthevoltagedistally(≈280␮mand≈425␮mfromthesoma)wasrecorded(blackandred electrodesinaandb)andtherecordingsoverlaid(tracesinc–e).

Fig.3.Fittingresultsforthereducedmodel2afterevolvingtheionchannelconductancesusingEMOOinthesecondstepofoptimization.Foursub-thresholdstepcurrents andonesupra-thresholdstepcurrentweredelivered,andtheneuron’sandmodel’svoltageresponsewasrecordedinthesomabeforeandduringpinchingoftheapical dendrite.Thesevoltagetraceswereusedfortheoptimization.Theleftpartofthefigurecomparestheexperimentalrecordings(black)withthemodelresponses(red)forthe intactneuronbeforepinching.Thetracesintherightpart(blueandorange)showthecorrespondingresponsesduringpinching.Acomparisonoftheblackandbluetraces (a–c,leftandrightparts)revealstheeffectspinchinghasontheneuronalresponseproperties.Weusedthefollowingfourobjectivesfortheoptimization:(1)theshapeof theAPonset(a,lefttraces);(2)theshapeoftheAPoffset(a,righttraces);(3)Theinterspikeintervaltimes(ISIs)ofthespiketrain(b,onlythefirst600msareshownfor bettervisualization).(4)Thefoursub-thresholdtraces(c).ThefivestepcurrentinjectionswereIamp=−0.1nA,−0.05nA,0nA,0.05nAand0.4nA(d,greentraces).Thefour objectives(a–c)weredeterminedbeforeandduringpinchingandcomparedwiththeexperimentaldata.Theresultingdistancesweresummedupyieldingfourdistance functionsthatweminimizedbyEMOO.Nocalciumchannelswerepresentinthisstepoftheoptimization.

loadon the AP initiation site in the axon. (4) The spike after- hyperpolarizingpotential(AHP)becomesstrongerduringpinching (≈3mV).ItisunknownwhatcausestheincreaseintheAHP.Under normalconditionsintheintactcell,back-propagatingactionpoten- tials(BAPs) (StuartandSakmann,1994)arecarriedbydendritic voltage-dependentcurrents,andthisinturnleadstodendriticcur- rents flowing backtothesomathat counteract theAHP. Upon pinching,thissourceofdepolarization isremoved,whichcould resultinanapparentincreaseoftheAHP.

Forourtargetdata,wechosea pyramidallayer5cell’s volt- agetracesin responsetofourweakcurrentinjections(−0.1nA,

−0.05nA,0nA,0.05nA)andonestrongercurrentinjection(0.4nA) that led tospiking responses under both conditions, before as wellasduringpinching.Thegoalofparameteroptimizationusing EMOOwastoreproducethesedataandtheeffects ofpinching.

EMOO automaticallydistributedthe activeconductances inthe axon,somaanddendrites.

Aspike was definedas a voltage excursionabove a thresh- old(=−20mV).Thespiketimewasgivenbythetimeatwhich the voltage reaches its maximum, whereas the spike width is measuredbyvoltagecrossingthethresholdbeforeandafterthe spike.

Inresponsetosupra-thresholdstimulationof0.4nA,themodel hadtorespondwithatleastsixspikes.Furtherprerequisitesforthe model’sresponseatthiscurrentwere:nospikewidthcouldexceed 3ms;theabsolutespikeheights(voltagepeaks)fromthethirdto

thepenultimatespikecouldnotchangebymorethan20%;thevolt- ageminimumbetweenthethirdandthefourthspikecomparedto thevoltageminimumbetweenthepenultimateandthelastspike shouldnotchangebymorethan10%;and,finally,thereshouldnot beanyinterspikeinterval(ISI)below15ms.Forallothercurrents, themodelwasrequirednottospike.

Iftheseprerequisiteswerenotmet,forboth theintactneu- ronaswellasforthemodelneuroninwhichtheapicaldendrite waspinched,theEMOOalgorithmseverelypunishedthemodelby assigningitanextremelyhigherrorvalue.Withtheprerequisites fulfilled,foursquareddistancevaluesweremeasuredbetweenthe modelresponseandtheexperimentaldata:

(1) We determined the distances between model and data foreachofthefoursub-thresholdtracesbetweent=−50msand 300ms (relativetothe currentinjectiononset, t0=100ms)and summedthesedistancesup:

E₁=

k

−50

(vkModel(t)−vkExp(t))²dt

(2)WedeterminedthesquareddistancebetweentheaverageAP onsetinthemodelandtheAPonsetinthedata.Forthispurpose, thetimesegmentbetween−0.5and−0.1msbeforetheAPpeak wasconsidered,andboththevoltageanditstimederivativewere

Fig.4.Modelgeneralization.Toillustratehowwelltheoptimizedmodelgeneralizes,wecomparethemodelresponses(redandorangetraces)toanewinputcurrent(c, 0.6nA)withthecorrespondingexperimentaldata(blackandbluetraces).Theseexperimentaltraceshavenotbeenusedduringtheoptimization.Shownarethemodel predictionsoftheAPshape(a)andofthespiketrain(b)beforeandduringpinching.

used:

E₂=

(v^Model(t)−v^Exp(t))²

+0.01·

dtv^Model(t)− d dtv^Exp(t)

dt

(3)WealsodeterminedthedistancebetweentheaverageAP offsets,includingtheirfirstderivatives(timewindowrangedfrom 0.1msto14msafterthespikepeaks):

E₃=

(v^Model(t)−v^Exp(t))²

+0.01·

dtv^Model(t)− d dtv^Exp(t)

dt

As we included the time derivative for feature extraction our approach resembles the idea of taking thephase-plane of a spike train for theoptimization (LeMasson and Maex, 2001).

Thatapproach,however,putsmoreweightonthesub-threshold responsesthanonthespikeshape.Thiscanleadtoimperfectfit- tingresults,especiallywhenexperimentaldataisusedasatarget (Druckmannetal.,2008).Weinsteadfocusontheaveragespike ratherthanonthewholespiketrainandhenceovercomethatlim- itation.Theexperimentaldataweusedfortheoptimizationwas recordedat50kHz(t=20␮s)anddoesnotprovidesufficienttime resolutionforproperspikealignmentandhenceforcalculatingthe averagespikeanditstimederivatives.Thereforeweappliedacubic splineinterpolation(newt=5␮s)toeachspikeandthenaligned andaveragedtheseinterpolatedspikes(WheelerandSmith,1988).

(4) Finally we also determined a distance for each ISI and summedthesedifferencesup.Thisdistancefunctionhasalready

proventobewellsuitedfortheoptimizationofconductance-based models(Kerenetal.,2005):

E4=

i

(ISI^Model_i −ISI^Exp_i )²

Foreachofthefourdistances,thevaluewasdeterminedbefore and during pinching and summed, yielding four final distance values. If all four distances are minimized, then the resulting modelreproducesexperimentalsub-thresholdresponses,spiking responsesaswellasdetailedAPshapebeforeandduringpinching.

Theeffectofpinchingwasintroducedintothemodelbyincreasing theaxialresistanceoftheapicaldendritetoahighvalue(10⁶cm).

Tominimizethefourdistancevalues,weusedEMOOwithapop- ulationsizeofN=1000andacapacityofC=2000individuals.1000 generationswereevaluated.Thecrossoverandmutationparam- etersremainedconstantatc=10andm=20respectively.The probabilityofamutationperparameterwas20%.Bydesign,EMOO endswithapopulationcontainingmodelsthatperformwellon oneofthedistancemeasuresalone,butcouldbefaroffthemarkin theothers.Butthefinalpopulationalsocontainsmodelsthathave intermediatefitnessineachofthedistancefunctions.Itisupto themodelertochooseanindividualfromthepopulation.Inorder todothisstepautomaticallyaswell,wenormalizedeachdistance functionbythelowestvaluefoundintherespectivedistancefunc- tionofthelastgeneration.Thisnaturalnormalizationthenallowed ustogothroughallgenerationsandpicktheindividualforwhich thesum(=totalerror,ET)ofthesenormalizeddistancevalueswas minimal.

Weperformedsixindependentoptimizationtrials,witheach trialrequiring approximatelytwo daysof runtime. Therelative improvementofthetotalerrorduringevolutionwas86±5%(mea- suredas(E_T(0)−E_T(i))/E_T(0);E_T(0)istheminimaltotalerrorinthe initialpopulationandET(i)istheminimaltotalerrorinthefinal

Fig.5.Modelpredictionoffiringfrequency.Wecomparetheexperimentallymea- suredwiththepredictedfiringfrequenciesbeforeandduringpinching(blackand reddotsandblueandorangedotsrespectively)forcurrentinjectionsrangingfrom 0nAto1.1nA.Thecurrentinjection’amplitudetakenformodeloptimizationwas 0.4nAwhichisillustratedwiththeblackdashedline.Errorbarsfortheexperimental datashowthestandarddeviationofthetwomeasurementsforeachdatapoint.

population)indicating a significantimprovement of thefit.We obtainedsixmodelsthatcloselyreproducetheexperimentalsub- thresholdresponses,thespiketrainaswellastheAPshapebefore andduringpinching,buthavedifferentsetsofparameters(models 1–6,Table1).Fig.3showsthefittingresultsformodel2.Theother fivemodelsreproducethedatasimilarlywell(notshown).

3.3.1. Optimizationwithoutthepinchingdata

Thedatasetwehaveusedforoptimizationisunusual,giventhat theneuron’sresponsewithandwithouttheapicaldendritewas measured.Wewishedtoquantifyhowmuchisgainedbyusingsuch data.Therefore,werepeatedtheoptimizationstrategydescribed abovebutexcludedthedistancesobtainedduringpinching and onlyusedmodelresponsesanddatafromtheintactneuron.This wasdoneforfourindependentruns(models7–10,Table1).The qualityofthefitimprovedsignificantlyduringevolutionasshown byarelativeimprovementofthetotalerrorof82±10%.Thesemod- elsreproducetheexperimentalrecordingsfortheintactneuron well,includingtheAPshape.Howeverthepredictionoftherecord- ingsduringpinchingispoorinallthesemodelsandalsovariable betweenmodels(seeFigs.S1andS2).Thisshowsthatthediffer- enceinneuronalresponsesbeforeandduringpinchingcontains usefulinformationaboutdendriticparameters.

3.3.2. Modelevaluation

Wecheckedhowwelltheoptimizedreducedmodelsgeneralize toinputcurrentsthatwerenotusedintheoptimizationprocess.

Forthiswecomparedtheresponsesofmodel2withtheexperi- mentaldataforanothercurrentinjection(0.6nA),whichthemodel predictswell(Fig.4,theothermodels’predictionsweresimilarly good).Moreover,themodel’sspikefrequencybeforeandduring pinchingisqualitativelycorrectforabroadsetofcurrentinjections (between0and1.1nA),withaslightlyhigherfiringfrequencyfor strongerinputcurrents(Fig.5).

We werealso interestedin how theoptimized models pre- dictotherexperimentalfindingsinpyramidalneurons.Models1–6 showedarestingpotentialmodulationwithdistancefromthesoma ofapproximately+2mV(Fig.6AandTable1)whichisqualitatively alsoseen in experiments (Stuart etal., 1997).Followingprevi- ousstudies(Kerenetal.,2009)weusedmodel2toevaluatethe conductancesactive atrestalong theapical dendrite.It canbe

seenthattheenhanceddendriticdepolarizationatrestisdueto anincreaseddendriticdepolarizingHCNcurrent(Fig.6B)anddue toalackofhyperpolarizingdendriticpotassiumcurrents(Fig.6D andE).Themodulationoftherestingpotentialin models7–10 wasvariable,andinmodel10theneuron’svoltagebecameeven morehyperpolarizedthefartherawayfromthesoma(Table1).

Thisshowsthatthepinchingdatais beneficialtoproperly con- straindendriticparameters.Next,wetestedifthemodelpredicts realisticBAPs(StuartandSakmann,1994;Stuartetal.,1997).In all10modelssomaticAPsactivelypropagatedintotheapicalden- dritewhiletheAPhalf-width(widthathalfwayfrom−60mVto theAPpeak)increased.Inmodels1–3and6–8theAPamplitude decreasedwithsomaticdistance,whileinmodels4,5,9and10the APamplituderemainednearlyconstantalongtheapicaldendrite (notshown).Fig.7illustratesBAPsformodel2.

Finallywetestedhowwellanaverageoftheoptimizedmodels 1–6wouldfittheexperimentaldata.Theresultingmodelpresented asurprisinglygoodfittothedata(Fig.S4)andthepredictionofthe experimentalIF–curveappearedevenbetterthaninanyofthe6 optimizedmodels(Fig.S5).Itmightbepossiblethattheparameter rangeleadingtogoodfittingresultsisbroadand thattheaver- agemodelstilllieswithinthisrange.Totestforthiswereplaced onlyasingleparameter(likeRm)permodelwithitsaverage,which producedamodelthatfailedreproducingtheexperimentaldata.

Thisshowsthatthequalityoftheaveragedmodelisratheranindi- cationthatcertainparameterratios(thataremaintainedduring averaging)determinethequalityofthefit.

3.4. Thirdstep:optimizingthecalciumspikedynamics

Olderpyramidalneuronsshowelaboratecalciumspikedynam- ics,i.e.astrongdendriticcurrentinputinducesalocaldendritic calcium spike, which can in turn depolarize the soma suffi- ciently to evoke axosomatic APs (Schiller et al., 1997; Larkum etal.,1999).Moreover,axosomaticallyinitiatedspikescanactively back-propagatealongtheapicaldendriteandreducethecurrent thresholdfordendriticcalciumspikeinitiation(Larkumetal.,1999, 2001).Thisthresholdreduction translatesintoacouplingfactor betweensomaanddendrite,estimatedtobearound0.5forpyra- midalneurons(Schaeferetal.,2003).

We wanted ourreduced model to reproduce thesecalcium spikedynamicsand showhow themethodcanbeextendedto matchrather differentdata sets.In general,fourparameters in themodeldominatethecalciumdynamics.Thestrengthandini- tiationthresholdofadendriticcalciumspikearedeterminedby the calcium channel density in the tuft ( ¯gCas) and the voltage shift(vshift_Cas)ofthischannel.Thecalcium-dependentpotassium channel curtailsthelength of thecalciumplateau and thereby thenumberofsomaticaction potentialsina burst.Thedensity ofthis channel ( ¯gKCa)wasalsooptimized. Theapicalresistivity (apicalRa),eventhoughithadbeenpreviouslyoptimizedinthe firststep,hadtobeleftasafreeparameteraswell.Allremain- ing parameters that had been previously optimized were not modified.

To illustrate the power of the evolutionary optimization approach,wedonotusequantitativeexperimentaldatatocon- structthedistancefunction.Instead,wetookgeneralexperimental observationsofhowdendriticcalciumdynamicsdependonthe couplingof thesomaandthedendrite(seeexperimentaltraces inFig.8)and constructeda discrete step-likedistancefunction describing which of these interactions the model qualitatively reproduced.Geneticalgorithmsfortheoptimization,asopposedto gradientdescent,forinstance,canhandlesuchastep-likedistance function.Foragivenparametercombinationwesetthedistance function(orerror)Easfollows:

-70 -69 -68

0 200 400 600 800 1000

Resting potential (mV)

Distance to soma (µm)

a

0 100 200 300 400

0 200 400 600 800 1000 gHCN (fS/µm2 )

Distance to soma (µm)

b

0 0.05 0.1

0 200 400 600 800 1000 gNat (fS/µm2 )

Distance to soma (µm)

c

0 20 40 60

0 200 400 600 800 1000 gKfast (fS/µm2 )

Distance to soma (µm)

d

0 100 200 300 400 500

0 200 400 600 800 1000 gKslow (fS/µm2 )

Distance to soma (µm)

e

Fig.6. Modelpredictionofdendriticpropertiesandchanneldensities.Therestingpotentialismoredepolarizedinthedistalregionsthanintheproximityofthesoma(a).

TheHCNchanneldensityincreasesalongtheapicaldendrite,addingtothedepolarizationatrest(b).Thesodiumconductanceatrestisnegligibleanddecayslinearlywith somaticdistance(c).Bothpotassiumchannelconductancesdecayrapidlywithsomaticdistance(d,e).Atrest,theseconductancesthusonlyhyperpolarizethesoma.

(1) We determined the threshold somatic current pulse required for a somatic spike. Observing more than a single somaticspikeatthresholdimpliesthatasingleback-propagating AP already elicits a dendritic calcium spike and thereby fur- thersomatic spikes. Thisis not seen in experiments, however.

Hence, a model exhibiting multiple spikes was penalized and associatedwiththehighesterrorvalue(E=5000)duringoptimiza- tion.

(2)Now,ifatthresholdonly asinglesomaticspikewasini- tiated, we nexttested whethera strongEPSP shapeddendritic currentinjectionalonecouldinducealocaldendriticcalciumspike.

Wealsocheckedwhetherthisalsoresultedinaquicklyforward spreading Caspike andeventually multiplesomatic spikes.We searchedforthedendriticcurrentamplitudethreshold(=ICA)that ledto this behavior. Based onexperimental observations, such a currentshouldelicit aburstof 2–4spikesinthesoma. Yeta somaticspikecanoccursimplyduetodepolarizationofthesoma, withoutadendriticspike.Sowedouble-checkedwhetherthefirst somaticspikewas,infact,duetoasomaticdepolarizationresult- ingfromadendriticCaspike.Forthispurpose,weintegratedthe voltageinthetuftfrom−10msto0msbeforethefirstsomatic spikepeak.Alargevalueforthisintegralindicatesthatthecalcium spikewastriggeredlocally;aslongasthetuftvoltageintegralwas largerthan500mVmsweconsideredthecalciumspiketobelocal.

Ifnocurrentamplitudewasfoundtoproducealocallyinitiated andforwardpropagatingcalciumspikewesettheerrorvalueto E=4000.

(3)Thenexttestwastodeterminehowaback-propagatingNa- spikefromthesomainfluencesthatthreshold.Tomeasurethis, weinjectedabriefsomaticcurrentpulseintothesomatoiniti- ateaback-propagatingAPandsearchedfortheminimaldendritic currentthreshold(=IBAC)neededtoinitiateacalciumspike.The resultingcalciumspikehadtofulfillthefollowingrequirements:It shouldproduceasomaticburstof2–4spikesanditshouldconsist ofprolongeddendriticdepolarizationwithafastshutoff.Toquan- tifytheserequirementswedetermined2voltageintegralsinthe tuft(I₁from0msto50msandI₂from100msto150msafterthe firstsomaticspikerespectively):I1hadtobelargerthan500mVms whileI2hadtobesmallerthan1000mVms.Ifnodendriticcurrent couldbefoundthatinitiatedacalciumspikewesettheerrorvalue toE=3000.

(4)Ifthemodelpassedalltheprecedingtests,theerrorwas determinedcompletelybythesomato-dendriticdegreeofcoupling (C)(Schaeferetal.,2003):

C= I_CA−I_BAC I_CA