Functional parcellation of the hippocampus based on its layer-specific connectivity with default mode and dorsal attention networks

(1)

ContentslistsavailableatScienceDirect

NeuroImage

journalhomepage:www.elsevier.com/locate/neuroimage

Functional parcellation of the hippocampus based on its layer-speciﬁc connectivity with default mode and dorsal attention networks

Gopikrishna Deshpande

^a^,^b^,^c^,^d^,^e^,^f^,^g^,^∗

, Xinyu Zhao

^a^,^h

, Jennifer Robinson

^a^,^b^,^c^,^d

aAU MRI Research Center, Department of Electrical and Computer Engineering, Auburn University, 560 Devall Dr, Suite 266D, Auburn, AL 36849, USA

bDepartment of Psychological Sciences, Auburn University, Auburn, AL, USA

cAlabama Advanced Imaging Consortium, Birmingham, AL, USA

dCenter for Neuroscience, Auburn University, Auburn, AL, USA

eKey Laboratory for Learning and Cognition, School of Psychology, Capital Normal University, Beijing, China

fDepartment of Psychiatry, National Institute of Mental Health and Neurosciences, Bangalore, India

gCentre for Brain Research, Indian Institute of Science, Bangalore, India

hQuora Inc., Mountain View, CA, USA

a r t i c le i n f o

Keywords:

Resting state functional connectivity Hippocampus

Parcellation

Layer-speciﬁc connectivity Default mode network Dorsal attention network HERNET model

a b s t r a ct

Recentneuroimagingevidencesuggeststhattheremightbeananterior-posteriorfunctionaldifferentiationofthe hippocampusalongthelong-axis.TheHERNET(hippocampalencoding/retrievalandnetwork)modelproposed anencoding/retrievaldichotomywiththeanteriorhippocampusmoreconnectedtothedorsalattentionnetwork (DAN)duringmemoryencoding,andtheposteriorportionsmoreconnectedtothedefaultmodenetwork(DMN) duringretrieval.EvidencebothforandagainsttheHERNETmodelhasbeenreported.Inthisstudy,wetestthe validityoftheHERNETmodelnon-invasivelyinhumansbycomputingfunctionalconnectivity(FC)inlayer- specificcortico-hippocampalmicrocircuits.Thiswasachievedbyacquiringsub-millimeterfunctionalmagnetic resonanceimaging(fMRI)dataduringencoding/retrievaltasksat7T.Specifically,FCbetweeninfra-granularout- putlayersofDANwithhippocampusduringencodingandFCbetweensupra-granularinputlayersofDMNwith hippocampusduringretrievalwerecomputedtotestthepredictionsoftheHERNETmodel.Ourresultssupport somepredictionsoftheHERNETmodelincludinganterior-posteriorgradientalongthelongaxisofthehippocam- pus.WhilepreferentialrelationshipsbetweentheentirehippocampusandDAN/DMNduringencoding/retrieval, respectively,wereobservedaspredicted,anterior-posteriorspecificityinthesenetworkrelationshipscouldnotbe confirmed.Thestrengthandclarityofevidencefor/againsttheHERNETmodelweresuperiorwithlayer-specific datacomparedtoconventionalvolumedata.

1. Introduction

The hippocampus, located in the medial temporal lobe (MTL), plays important roles in many brain functions including episodic memoryandspatial navigation(Das etal., 2011; Yushkevichet al., 2010; Heckemann et al., 2011). Investigation of functional special- ization within thehippocampus has receivedincreased attention in neuroimaging.Theabnormalitiesofthehippocampushavebeeniden- tiﬁedin many neuropsychiatric disorders,includingAlzheimer’s dis- ease(Lukiw,2007; Zhou et al., 2008; Scheﬀ and Price,2006), ma- jor depression (Campbell and MacQueen, 2004; Stockmeier et al., 2004;Rossoetal.,2005;MacQueenandFrodl,2011),post-traumatic stressdisorder(Asturetal.,2006;JavidiandYadollahie,2012),and schizophrenia(Harrison,2004;Heckers,2001;Grace,2012).Therefore

∗Correspondingauthorat:AUMRIResearchCenter,DepartmentofElectricalandComputerEngineering,AuburnUniversity,560DevallDr,Suite266D,Auburn, AL36849,USA.

E-mailaddress:gopi@auburn.edu(G.Deshpande).

abetterunderstandingofthefunctionalspecializationwithinthehip- pocampushasthepotentialtoleadtoabetterunderstandingofthese disorders.

Recently,manystudieshavepositedthattheremaybeafunctional diﬀerentiationalongthelong-axisofthehippocampus(Robinsonetal., 2015; deWaeletal.,2018).Lepageetal.(1998)performed ameta- analysisofpositronemissiontomography(PET)ofepisodicmemoryand discoveredanorderlyfunctionalanatomicpatterninthehippocampus.

Morespeciﬁcally, theanteriorportionswereprimarilyactivatedwith episodicmemoryencoding,whereastheposteriorportionswereprimar- ilyassociatedwithepisodicmemoryretrieval.Thismodelisreferredto asHIPER(hippocampusencoding/retrieval)model.

TheHIPERmodelhasreceivedsupportfrommanyrecentstudies.

Spaniol et al.(2009) conductedmeta-analysesof event-relatedfMRI studiesofepisodicmemoryandrevealedananterior-posteriorgradient

https://doi.org/10.1016/j.neuroimage.2022.119078.

Received4June2019;Receivedinrevisedform29January2022;Accepted7March2022 Availableonline9March2022.

(2)

G. Deshpande, X. Zhao and J. Robinson NeuroImage 254 (2022) 119078

inthehippocampalactivationsassociatedwithencodingandretrieval.

Nadeletal.(2013)andhiscolleaguescomparedtheanteriorandposte- riorhippocampalactivationsduringretrievalofdifferenttypesofspatial information.Theyfoundthatthereisafunctionaldifferentiationalong thelongitudinalaxisofthehippocampuswiththeposteriorhippocam- pus being crucialfor precise spatialbehavior, andthe anterior hip- pocampusbeinginvolvedincontextcoding.BaumannandMattingley (2013)examinedretrieval-relatedactivityinthetaxidrivers’hippocam- pususingfMRIdataandfoundthattaxidriverswithasmallanterior hippocampushaddifficultyencoding newspatialassociations.More- over,Kim(2015)conductedameta-analysisandrevealedthattheen- codingofsensoryinputinvolvedmainlytheanteriorhippocampusand theexternalattentionnetwork,whereasretrievalengagedmainlythe posteriorhippocampusandtheinternalattentionnetwork.Thismodel wasreferredtoastheHERNET(hippocampalencoding/retrievaland network)model.

Memory encoding is inherently linked with external attention, whereas retrieval is intrinsically related to internal attention (Kim,2010;Chunetal.,2011).Manystudieshaveidentifiedtwobrain networks,i.e.,thedorsalattentionnetwork(DAN)anddefaultmode network(DMN),that areclosely associated withexternal andinter- nal attention,respectively (Buckner etal., 2008; Corbetta andShul- man, 2002). Therefore, the HERNET model predicts that the ante- riorhippocampusandregionsoftheDANco-activateduringencoding whiletheposterior hippocampusandregionsof theDMNco-activate duringretrieval. InfactKim(2015) meta-analysisconfirmsthis pre- diction. However, evidenceconflicting the HERNET model also ex- ists.This includesmeta-analysesofimagingstudieswhichcontradict Kim’sfindings(SchacterandWagner,1999),suggestionsthatthean- teriorhippocampusisactivatedbynovelty(whichispurportedlymis- taken forencoding sinceencoding tasks typicallyuse novel stimuli) (Kumaran andMaguire,2006;Poppenketal.,2008;Zweynertetal., 2011;Poppenketal.,2010a2010b)andalternativemodelsoffunctional specializationwhichattribute“hot” processing(emotion/motivation)to anteriorhippocampusand“cold” processing(cognition)totheposterior part(Murtyetal.,2011;Robinsonetal.,2015;Robinsonetal.,2016).

Giventhisstateofaffairs,wesetouttodirectlytesttheHERNETmodel usingfunctionalconnectivitybetweenthehippocampusandDAN/DMN regionsduringmemoryencodingandretrievaltasks.Unlikeprevious studieswhichemployedvoxel-levelanalysisfromdataobtainedatcon- ventionalfieldstrengths(≤3T),weemployedfMRIwithultra-highspa- tialresolution(sub-millimeter)obtainedat7T.Thisallowedustoin- vestigatelayer-specificmicrocircuitsbetweenDMN/DAN regionsand thehippocampuswiththehypothesisthattheywouldprovideafiner grainedcharacterizationof theconnectivitybetweenthem, whichin turnmayprovidemoredefinitiveevidencefororagainsttheHERNET model.Byfinegrainedcharacterization,wemeanthattheconnections betweenthehippocampusandcorticalregionsin theDMNandDAN areanatomicallylayer-specific.Yet,thefunctionalconnectivityanaly- sesofthesecircuitsunderlyingtheHIPER/HERNETmodelsarecarried attheconventionalvolumelevel,whichignoresthisbiologicalreality andacceptsthelossofspatialspecificity.Ourpointisthatcurrenttech- nologyallowsustoovercomethis limitationby acquiringdatawith higherspatialspecificity,whichcouldprovideevidencefororagainst theHIPER/HERNETmodelwhichisclosertothebiologicalrealityof theselayer-specificcircuits.Wewouldalsoliketonotethatwhilethe ventralattentionnetwork(VAN)alsocontributestotheallocationofat- tentionalresourcestomemoryprocessing,werestrictourselvestothe DANsinceweareinterestedinspecificallytestingtheHERNETmodel whichhypothesizesonlytheDAN.

Previousinvasivestudiesinanimalshaveshownthattheconnec- tionsbetweenthehippocampusandDAN/DMNprimarilyinvolvecor- ticallayersIIandV,withlayerVofthehigherordercortex(frontal andparietalcortices)projectingtothehippocampus,whereaslayerII ofhigherordercortexreceivingthesignalbackfromthehippocampus (SwansonandCowan,1977;ThomsonandBannister,2003).Thislayer-

speciﬁcpathway betweenthehippocampusandDAN/DMNisnotex- clusivesincepathwayswhichoriginate/terminateinotherlayersofthe cortexmayalsocontributetothehippocampalinputoroutput.Thisis notsurprisinggiventhehighlycomplexunderlyingmicrocircuitryand giventhatsignalsbetweenanytwobrainregionscanrelayviamulti- plestructuresincludingthethalamus.However,thepathwaysbetween thehippocampusandlayersIIandVoftheDMN/DANseemtobethe dominantonesbasedonpriorinvasiveanimalliterature(Shepardand Grillner,2010;ThomsonandLamy,2007).

ThepresentstudysoughttoinvestigatetheHIPER/HERNETmodel using the functionalconnectivity between the hippocampus and (1) deeperlayersofDAN/DMNduringanencodingtask,and(2)superﬁ- ciallayers ofDAN/DMN duringaretrievaltask.Speciﬁcally,wehy- pothesizedthatduringamemoryencodingtask,clusteringhippocam- pal voxelsbasedontheirfunctionalconnectivitywithlayerVof the DANmustparcellatethehippocampus inananterior-posteriorgradi- entalongthelongaxis.Similarly,duringamemoryretrievaltask,clus- teringhippocampalvoxelsbasedontheirfunctionalconnectivitywith layerIIoftheDMNmustalsoshowananterior-posteriorsegmentation.

Second, duringanencodingtask,thehippocampalvoxelsmusthave strongerconnectivity withlayerV ofthe DANthanwith layerVof theDMN.Contrarily,duringaretrievaltask,thehippocampalvoxels musthave strongerconnectivity withlayerIIof theDMN thanwith layerIIoftheDAN.Third,consideringthedirectionalityofsignalpro- jection,duringanencodingtask,layerVoftheDANmustshowstronger correlationwithanteriorhippocampalregionsthanwithposteriorhip- pocampalregions,whereasduringretrievaltask,layerIIoftheDMN mustexhibitstrongercorrelationwithposteriorhippocampalregions thanwithanteriorhippocampalregions.Finally,wepredictedthatus- inglayer-specificdatawouldleadtomoredefinitiveresultsthanusing conventionalvolume-leveldatawhileinvestigatingtheconnectionbe- tweenthehippocampusandtheDAN/DMN.Allfourhypotheses,based onpredictionsfromtheHERNETmodel,areillustratedinFig.1.Here itisnoteworthythatalthoughweintendtotestthesehypothesesus- ingfMRI dataextracted from layersII andVof DAN/DMN regions, somepartialvolumeeffectsareexpectedforthespatialresolutionof ourdata.Therefore,signalfromlayerIIbroadlyrepresentsthosefrom supra-granularlayerswithdominantcontributionfromlayerII.Like- wise,signalfromlayerVbroadlyrepresentsthosefrominfra-granular layerswithdominantcontributionfromlayerV.

In order totest these hypotheses, fMRI data was acquired from healthy subjects performing memory encoding and retrieval tasks with faces, pictures of scenes and words in the 7T scanner. Ultra- high resolution anatomical data was also acquired to resolve lay- ersintheDAN/DMN.Unsupervisedclusteringmethodswereapplied, which groupsobjectsin adatadrivenwaywithout using anylabels to guide theresults, to parcellatethe hippocampalvoxels based on theirconnectivitywithdifferentlayersintheDAN/DMN.Threeclus- tering methodswere specifically chosen, i.e.,hierarchical clustering (DasguptaandLong,2005),orderingpointstoidentifytheclustering structure(OPTICS)(Ankerstetal.,1999),anddensitypeakclustering (DPC)(RodriguezandLaio,2014),sincetheydidnotrequireapriori specificationofthenumberofclusters.Sinceclusteringaccuracyisoften lowerinhighdimensionalfeaturespace,ageneticalgorithm(GA)based featureselectionmethodwasalsoemployed,whichislesspronetolo- caloptimum(DyandBrodley,2004),comparingwithotherexistingfea- tureselectionmethods,e.g.,sequentialforwardsearching(DyandBrod- ley,2004),non-linearoptimization(BradleyandMangasarian,1998), etc.

2. Materialsandmethods

In this work theconnectivity between the hippocampus and the DAN/DMNwereinvestigatedbyclusteringhippocampalvoxels,inan unsupervised way, based on their connectivity with (1) layer V of DAN/DMNduringanencodingtask,and(2)layerIIofDAN/DMNdur-

(3)

Fig.1. Illustrationofourhypothesesbasedon predictionsfromtheHERNETmodelaswellas knownanatomicalpathwaybetweendiﬀerent layersofDAN/DMNandhippocampus.

ingaretrievaltask.Theidentiﬁedfunctionalclustersofthehippocam- puswerethencomparedwithananatomicalanterior-posteriorsegmentation(Poppenketal.,2013).Theentireanalysispipelineisillustrated inFig.2andwillbeelaboratedbelow.

2.1. Dataacquisition

Thirty-onehealthyindividuals(26right-handed,12 males,19fe- males,age =21.1± 1.4)were recruitedforthestudy. TheInternal ReviewBoard(IRB)atAuburnUniversityapprovedthestudy,subjects providedinformedconsentandtheexperimentalprocedureswereper- formedinaccordancewithinternationallyacceptedethicalstandards.

Echo-planarimaging(EPI)datawereacquiredontheAuburnUniver- sityMRIResearchCenter(AUMRIRC)Siemens7TMAGNETOMscanner outﬁttedwitha32-channelheadcoilbyNovaMedical(Wilmington, MA).Thesequencewasoptimizedforthehippocampus(37slicesac- quiredparalleltotheAC-PCline,0.85mm×0.85mm×1.4mmvoxels, TR/TE:3000/28 ms,70° ﬂip angle,base/phaseresolution: 234/100, A→Pphaseencodedirection,iPATGRAPPAaccelerationfactor=3,in- terleavedacquisition,123 timepoints,total acquisitiontime6min).

Duringencodingtask,theparticipantswereaskedtoviewaseriesof faces,picturesofscenesandwords.Eachtriallastedfor30sinwhich thesubjectswere presented10 imagesof thesame category for3s each.Thesetrialswereinterspersedwitha6sinter-trialinterval.The paradigmfortheretrievaltaskwasidenticaltotheencodingtaskwith theexceptionthatthesubjectswereprovidedwithanMR-compatible buttonboxtoindicate,viabuttonpresses,whethertheyrecognizedthe imageashavingseenduringtheencodingtaskornot.Awhole-brain high-resolution3DMPRAGEsequencewasusedtoacquireanatomical data(256slices,0.63mmx0.63mmx0.60mm,TR/TE:2200/2.8,7°

ﬂipangle,base/phaseresolution384/100%,collectedinanascending fashion,acquisitiontime=14:06)forextractingdiﬀerentcorticallayers intheDAN/DMNROIsandforregistrationpurposes.

2.2. Preprocessing

(a) Hippocampaldata:Standardpre-processingstepswerecarriedoutus- ingSPM8(Ashburner,2012)includingbrainextraction,slicetiming correction,motioncorrection,regressionofmotionandphysiologi- calartifacts(usingCompCor(Behzadietal.,2007)),registrationto anatomicalspace,andnormalizationtoMNIstandardspace.

TheHarvard-OxfordStructuralProbabilityAtlasdistributedwiththe FSLneuroimaginganalysissoftwarepackage(http://www.fmrib.ox.

ac.uk/fsl/fslview/atlas-descriptions.html#ho) was used to deﬁne hippocampalROIs.Ourtaskdidnotexplicitlyinvolvespatialmem- orytypicallyusedinnavigation,whichisprimarilyassociatedwith therighthippocampus (Burgess etal.,2002).Therefore,only the lefthippocampus wasconsideredin thisstudy.Inorder todeter- mine a conservative anatomicalrepresentation, the hippocampal ROI was thresholded at 75%.The meanprobability forthe vox- elsinthehippocampalROIsbelongingtothehippocampus(M±SD:

86.41%± 7.10%)andthevoxelcentroid(MNIcoordinates[−26,

−18.8, −17.2]) belonging to the hippocampus was > 86% and 97.1%,respectively.TheidentiﬁedhippocampalROIsareillustrated inFig.3withatotalvolumeof1880mm³.

(b) Corticallayer-specificdataoftheDAN/DMN:Dataextractedfromthe wholebrainwasfirstpreprocessedusingSPM8.Topreservehigh spatialresolution,nospatialfilteringwasapplied.Thecorticallay- erswerethenreconstructedusingFreeSurfer(Fischl,2012).Specifi- cally,twointerfaces,i.e.,thecorticalgraymatterandtheunderlying whitematter(white-grayinterface)andtheinterfacebetweenthe corticalgraymatterandthepialsurface(gray-pialinterface),were automaticallyreconstructedfromtheanatomicalimage.Then,the corticalthicknesswascalculatedastheaverageofthedistancefrom thewhite-grayinterfacetotheclosestpossiblepointonthegray-pial interface,thenfromthatpointbacktotheclosestpointonthewhite- grayinterfaceagain.Toimproveaccuracy,surfacesmoothingand

(4)

Fig.2. Illustrationofproposedanalysispipelineforinvestigatinghippocampalparcellationbasedonitslayer-speciﬁcconnectivitywithDAN/DMNROIs.Thesame processwasrepeatedforencodingandretrievaltasks,separately,aswellaswithconventionalvolumedata(asopposedtolayer-speciﬁcdata).

(5)

Fig.3.Hippocampal ROIusedinthisstudy.

Onlythelefthippocampuswasconsideredin thisstudy.

automatictopologycorrectionwerealsoapplied(Daleetal.,1999; Fischletal.,1999;FischlandDale,2000;Hanetal.,2006).Fromthe corticalthicknessmap,sixcorticallayerswerethenreconstructed withinthecorticalgraymatteratfixedrelativedistancebetweenthe whiteandpialsurfacesdeterminedfromthecorticalthickness,i.e., thefirstlayerwaslocatedat96%ofthecorticalthicknessawayfrom thewhitematter,thesecondlayerat80%,thethirdlayerat64%,the fourthlayerat48%,thefifthlayerat32%,andthesixlayerat16%

(Fig.4).Thelaminarlayerswerederivedfromtheanatomicalimage, soitwasnecessarytoaligntheEPIdatatotheselayers.Aboundary- basedregistrationmethod(GreveandFischl,2009)wasemployed, whichalignstheEPIimagetotheanatomicalimagebymaximizing theintensitygradientacrosswhite-grayinterfaceandgray-pialinterface.Theentirecorticalsurfacewasautomaticallydividedinto 34corticalROIsineachoftheindividualhemispheresbasedonthe Desikan-KillanyatlasinFreesurfer(Desikanetal.,2006).Themajor DANregions,i.e.,frontaleyefield[FEF],inferiortemporalcortex [ITC],inferiorfrontalgyrus[IFG],andsuperiorparietallobe[SPL], andmajorDMNregions,i.e.,anteriorcingulatecortex[ACC],me- dialprefrontalcortex[mPFC],inferiorparietallobe[IPL],posterior cingulatecortex[PCC],andprecuneus,wereidentifiedfromthose 34ROIs,andvertices(voxelsinthevolumebecomeverticesonsur- faces)inlayerIIandlayerVoftheseregionswerethenidentified.

2.3. Connectivitymeasures

FCmeasuresthefunctionalinterrelationshipbetweenpairsofbrain regions byestimating Pearson’s correlationbetween time seriesrep- resentingthose brainregions.Functional connectivity (FC)was esti- matedbetweenthehippocampalvoxelsand(1)verticesinlayerVof theDAN/DMNduringtheencodingtask,and(2)verticesinlayerIIof theDAN/DMNduringtheretrievaltask.

2.4. Layerspeciﬁcclustering

The hippocampal voxels were parcellated using three clustering methodsbasedontheirFCwith(1)verticesinlayerVoftheDAN/DMN duringtheencodingtask,and(2)verticesinlayerIIoftheDAN/DMN during the retrieval task. The same process was repeated on DAN and DMN volume ROIs, separately. Let 𝑌 ={𝑌1,⋯,𝑌_𝑖, ⋯,𝑌_𝑁} represent aset of 𝑁 objects, i.e.,numberof hippocampal voxels. 𝒀_𝑖= (𝑌_𝑖1,𝑌_𝑖2,⋯,𝑌_𝑖𝑑)∈ℝ^𝑑,wheredequalstothenumberofFCfeatures.As- sumetheNobjectsareseparatedintokclusters.Eachclusterisaset ofindexesfrom{1,⋯,𝑁},andeachobject𝑌_𝑖belongstoexactlyone cluster.

(a)HierarchicalClustering(Agglomerative):Asoneofthemostcom- monlyusedconnectivity-basedclusteringmethod,thehierarchicalclus- tering(Liaoetal.,2008;Chengetal.,2006;DasguptaandLong,2005) groupsobjects into diﬀerent clustersby building a hierarchicaltree structure.Theprocedureofthismethodisillustratedbelow:

(1) Initially,eachobject𝑌_𝑖isassignedtoaclusterwithonlyitselfinit.

Fig.4. CorticallayerreconstructionwithFreeSurfer.Sixcorticallayerswere reconstructedwithinthecorticalgraymatteratﬁxedrelativedistancesbetween thewhiteandpialsurfaces.

(2) Distancebetween anytwoclustersismeasured.Then,theclosest pairofclustersaremerged.

(3) Step2–3arerepeateduntilall𝑌_𝑖areinonebigcluster.

The resulting tree structure is usuallyreferred toas thedendrogram (Fig.5).The rootof thedendrogram represents the entire data, eachleaf noderepresentsone singleobject, andtheheight of the dendrogram represents thedistance between each pair of clusters. Differentdata partitionscan be obtained by cutting the dendrogramatdifferentheights.Notethattherearethreelinkage criterions, single-linkage, complete linkage, and average-linkage, whichhavebeenwidelyusedinmeasuringdistancebetweentwo clusters.Thesinglelinkage(Andreopoulosetal.,2008)calculates theshortest distance between twoclusters, the complete linkage (Andreopoulosetal.,2008)calculatesthelongestdistance,andthe averagelinkage(Andreopoulosetal.,2008)calculatesthemeandis- tance.Thesinglelinkagemethodcanhandlenon-ellipticalshapeof clusters,butcanbeaffectedbynoiseandoutliers.Thecompletelink- agemethodislesssensitivetonoiseandoutliersbuttendstobreak largeclusters.Theaveragelinkageisacompromisebetweensingle- linkageandcompletelinkagemethods.Thus,theaveragelinkage methodwasemployedinthiswork.

(6)

Fig.5.Illustrationofhierarchicalclustering.(a)Originalsimulateddataset,(b)Dendrogramderivedfromhierarchicalclustering,and(c)Clusteringresultsobtained withaspecificcuttingheight.Twoclustersthatwereidentifiedaremarkedwithdifferentcolors.

(b) Ordering Points to Identify the Clustering Structure: OPTICS (Ankerstetal., 1999)is oneof themostpopulardensity-basedclusteringmethods(Kriegeletal.,2011)Givenadistancethreshold(𝜀)and theminimumnumberofobjectsrequiredtoformacluster(𝑀𝑖𝑛𝑃𝑡𝑠), objectsinhigh-densityareasaregroupedtogether,whereasobjectsin sparseareas,whicharerequiredtoseparateclusters,areusuallyconsid- eredtobenoiseoroutliers.OPTICScandiscoverclusterswitharbitrary shapesandhastheabilitytoidentifyoutlierobjectsthatdonotbelong toanyoftheclusters.

InOPTICS,twovariablesarecomputedforeachobjectinthedataset:

core-distance and reachability-distance. Let 𝑁_𝜀(𝑌_𝑖) represent the num- berofnearbyobjectswithin𝜀(called𝜀−𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑),and𝑀𝑖𝑛𝑃𝑡𝑠− 𝑑𝑖𝑠𝑡𝑎𝑛𝑐𝑒(𝑌_𝑖)representthedistancefrom𝑌_𝑖toits𝑀𝑖𝑛𝑃𝑡𝑠’neighbor.An object𝑌_𝑖isacoreobjectifatleast𝑀𝑖𝑛𝑃𝑡𝑠objectsarefoundwithits

𝜀−𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑.Thecore-distanceof𝑌_𝑖isdeﬁnedas:

𝑐𝑜𝑟𝑒−𝑑𝑖𝑠𝑡𝜀,𝑀𝑖𝑛𝑃 𝑡𝑠( 𝑌𝑖)

=

{𝑈𝑛𝑑𝑒𝑓𝑖𝑛𝑒𝑑𝑖𝑓 𝑁_𝜀( 𝑌_𝑖)

<𝑀𝑖𝑛𝑃𝑡𝑠 𝑀𝑖𝑛𝑃𝑡𝑠−𝑑𝑖𝑠𝑡𝑎𝑛𝑐𝑒(

𝑌𝑖)

𝑜𝑡ℎ𝑒𝑟𝑤𝑖𝑠𝑒 (1)

which is the smallest distance for 𝑌_𝑖to have 𝑀𝑖𝑛𝑃𝑡𝑠 in its 𝜀− 𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑.

Thereachability-distanceofobject𝑌𝑗withrespecttoobject𝑌𝑖isde- ﬁnedas:

𝑟𝑒𝑎𝑐ℎ𝑎𝑏𝑖𝑙𝑖𝑡𝑦−𝑑𝑖𝑠𝑡_𝜀,_{𝑀𝑖𝑛𝑃 𝑡𝑠}( 𝑌_𝑗, 𝑌_𝑖)

={ 𝑈𝑛𝑑𝑒𝑓𝑖𝑛𝑒𝑑𝑖𝑓𝑁𝜀( 𝑌𝑖)

<𝑀𝑖𝑛𝑃𝑡𝑠 𝑚𝑎𝑥(

𝑐𝑜𝑟𝑒−𝑑𝑖𝑠𝑡𝑎𝑛𝑐𝑒( 𝑌_𝑖)

,𝑑𝑖𝑠𝑡( 𝑌_𝑗,𝑌_𝑖))

𝑜𝑡ℎ𝑒𝑟𝑤𝑖𝑠𝑒 (2)

Where𝑑𝑖𝑠𝑡(𝑌_𝑗,𝑌_𝑖)isthedistancemeasure(e.g.,Euclideandistance)be- tween𝑌_𝑗and𝑌_𝑖.ThecompleteprocedureofOPTICSisdescribedbelow:

(7)

Fig.6.IllustrationofOPTICSclustering.(a)Originalsimulateddataset,(b)reachabilityplotobtainedfromOPTICS,and(c)clusteringresults.Twoclusterswere identiﬁedcorrespondingtovalleysinthereachabilityplot.

(1) Chooseoneobject𝑌_𝑖arbitrarily.

(2) Retrieve the 𝜀−𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑 of 𝑌_𝑖, determine the core-distance of 𝑌_𝑖, and set the reachability-distance of each object 𝑌_𝑗in the 𝜀−𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑of𝑌𝑖toundeﬁned.

(3) If𝑌_𝑖isnotacoreobject,gotostep5.Otherwise,gotostep4.

(4) For each object 𝑌_𝑗 in the 𝜀−𝑛𝑒𝑖𝑔ℎ𝑏𝑜𝑟ℎ𝑜𝑜𝑑 of 𝑌_𝑖, update its reachability-distancefrom𝑌𝑖andinsert𝑌𝑗 intoanOrderSeedslist ifithasnotbeenprocessedyet.

(5) Iftheinput datasetis fully consumed andtheOrderSeeds list is empty,gotostep6.Otherwise,moveontothenextobjectinthe OrderSeedslist(ortheinputlist,iftheOrderSeedslistisempty)and gotostep2.

(6) Outputcore-distance,reachability-distanceofeachobject,andpro- cessedorder.

Thedataobjectsareplotted intheprocessedorder togetherwith theirrespectivereachability-distance(calledreachabilityplot)depicting thehierarchicalstructureoftheclusters.Sinceobjectsbelongingtoa clusterhavealowreachability-distancetotheirnearestneighbor,the clustersshowupasvalleysinthereachabilityplot(seeFig.6).Theﬁnal

datapartitioncanbeobtainedbyusingathresholdonthereachability plot.

(C)DensityPeakClustering:recentlyRodriguezandLaio(2014)proposed anovel density-based clustering method (referredtoas DPC) basedontheideathattheclustercentersarecharacterizedbyahigher densitythantheirneighborsandbyarelativelylargedistancefromob- jectswithhigherdensities.Likeotherdensity-basedclusteringmethods, e.g.,OPTICS,ithasabilitytodetectarbitrarilyshapedclustersandspot outlierobjects.Moreover,DPCoutperformscommonlyusedclustering methods, e.g.,k-meansandhierarchicalclustering, whenthedataset containscomplicatedfeaturessuchasnarrowbridgesbetweenclusters, uneven-sizedclusters,clusterswithhighoverlap,etc.

Foreachobject𝒀_𝒊,twoquantitiesarecomputed:localdensity𝝆(𝒀_𝒊) andminimumdistancewithhigherdensity𝜹(𝒀_𝒊).𝝆(𝒀_𝒊)isdeﬁnedas:

𝝆( 𝒀_𝒊)

= ∑

𝒊 𝝌( 𝒅𝒊𝒔𝒕(

𝒀_𝒊, 𝒀_𝒋)

− 𝒅_𝒄)

(3)

where𝒅_𝒄isacutoﬀ distance,and𝝌(𝒚)canbecomputedby, 𝝌(𝒚)=

{1𝒊𝒇𝒚<0

0𝒐𝒕𝒉𝒆𝑟𝒘𝒊𝒔𝒆 (4)

(8)

Fig.7. IllustrationofDPCclustering.(a)Originalsimulateddataset,(b)Plotof𝛿asafunctionof𝜌foreachobject.Objectswithlarger𝜌and𝛿areclustercenters andobjectswithsmaller𝜌,andlarger𝛿areoutliers.(c)Clusteringresults.Twoclusterswereidentiﬁedcorrespondingtotwoclustercentersinthedecisiongraph.

FromEqs.(3)to(4),itcanbeseenthat𝝆(𝒀_𝒊)equalstothenumber ofobjectswithin𝒅_𝒄withrespecttoobject𝒀_𝒊.𝜹(𝒀_𝒊)ismeasuredby, 𝜹(

𝒀_𝒊)

= min

𝒊∶𝝆( 𝒀_𝒋)

>𝝆(𝒀_𝒊)𝒅𝒊𝒔𝒕( 𝒀_𝒊,𝒀_𝒋)

(5)

Fortheobjectwithhighestdensity,𝜹(𝒀_𝒊)isconventionallysetto, 𝜹(

𝒀𝒊)

=max

𝒊 𝒅𝒊𝒔𝒕( 𝒀𝒊,𝒀𝒋)

(6) Notethatif𝒀_𝒊islocalorglobalmaximainthedensity,𝜹(𝒀_𝒊)willbe muchlargerthanitstypicalnearestneighbor.Thus,objectswithlarger 𝝆and𝜹areconsideredasclustercenters,whereasobjectswithsmaller 𝝆andlarger𝜹areconsideredasoutliers.Otherobjectsareassignedto thesameclusterastheirnearestneighborofhigherdensity(seeFig.7).

(d)InputParameterOptimization: Ineach clusteringmethod,there areseveraluser-speciﬁedinputparameters,whichcansigniﬁcantlyaf- fect theshapeof thecluster andthenumberof thecluster (Fig.8).

Forhierarchicalmethod,thecuttingheightofthedendrogramneeds tobespecifiedandthenumberofclustersvarieswithdifferentcutting heights.ForOPTICS,𝜺cansimplybesettothemaximumpossiblevalue, andAnkerstetal.(1999)showedthatfor𝑴𝒊𝒏𝑷𝒕𝒔usingvaluesbetween 10and20wouldalwaysleadtogoodresults.However,thethresholdfor thereachabilityplot,whichisusedtoextractclusters,stillneedstobe properlydetermined.ForDPC,𝒅_𝒄canbechosenbasedontherulethat theaveragenumberofneighborsisaround1–2%ofthetotalnumber ofobjectsinthedataset(RodriguezandLaio,2014).Athresholdfor 𝝆and𝜹needstobedefinedtodistinguishclustercenters,borders,and outliers.Inthisstudy,theaveragesilhouetteindex(Rousseeuw,1987) wasemployedtodeterminetheoptimalvaluesoftheseparameters.As- sumethedatahavebeenclusteredviaanyclusteringalgorithm,suchas OPTICS,intoKclusters.Foreachobject𝒀_𝒊,let𝒂(𝒀_𝒊)representtheav- eragedistanceof𝒀𝒊withallotherobjectinthesamecluster,and𝒃(𝒀𝒊) representthesmallestaveragedistanceof𝒀_𝒊 toanyothercluster,of

(9)

Fig.8.Illustrationofdependencyofclusteringresultsoninputparameters.Takehierarchicalclusteringasanexample.Witharelativelyhighcuttingheight(a), twoclusters(redandblue)wereidentiﬁed(b).Ascuttingheightwasreduced(c),onebigcluster(red)wasseparatedintotwosmallerclusters(d).

which𝒀_𝒊doesnotbelongto.Thenthesilhouetteindexof𝒀_𝒊isdeﬁned as:

𝒔( 𝑿𝒊)

= 𝒃( 𝒀_𝒊)

−𝒂( 𝒀_𝒊) 𝐦𝐚𝐱{

𝒂( 𝒀_𝒊)

, 𝒃(

𝒀_𝒊)} (7)

FromEq.(7),itcanbeseenthat𝒔(𝒀_𝒊)isboundedbetween−1and 1.If𝒀_𝒊 hasbeenassignedtoa“correct” cluster,𝒔(𝒀_𝒊)willbecloseto 1.Contrarily,if𝒀_𝒊hasbeenassignedtoa“wrong” cluster,𝒔(𝒀_𝒊)willbe closeto−1.𝒔(𝒀_𝒊)willbecloseto0,if𝒀_𝒊islocatedontheborderoftwo naturalclusters.Bycomputingtheaverage𝒔(𝒀_𝒊)overallobjectsinthe entiredataset,theaccuracyoftheclusteringresultscanbequantiﬁed.

Withtheaveragesilhouetteindexastheoptimizationcriterion,a

“gridsearch” (BergstraandBengio,2012)methodwasappliedonde- terminingtheoptimalvalueofeachparameter.Theoptimalnumberof clusterscanbedeterminedsimultaneously.Forexample,inOPTICS,we startedwitharelativelyhighthresholdforthereachability-plot.Ineach iteration,thethresholdwasreducedbyasmallamountandtheaver- age𝒔(𝒀_𝒊)wascomputedandrecordedbasedonthecurrentpartition.

Theiterationcontinuesuntilthethresholdwassmallerthanaspeciﬁed baseline,e.g.,theaveragereachability-distanceofthereachability-plot.

Theoptimalthresholdwasthendeterminedastheonewiththelargest average𝒔(𝒀_𝒊).Thesameiterativeprocedurewasappliedtohierarchical clusteringtodetermineoptimalcuttingheightofdendrogram,andto DPCtodeterminetheoptimalthresholdof𝝆and𝜹.

2.5. Featureselectionandclusteridentiﬁcation

Theclusteringaccuracyisoftenlowerinhighdimensionalfeature space,which maybebecausemostoffeaturesinthedatasetmaybe irrelevant,redundant,orsometimesmayevenmisguideresults.More- over,alargenumberoffeaturesmaketheclusteringresultsdiﬃcultto interpret.Therefore,afeatureselectionmethodisrequiredtoimprove theclusteringaccuracy.Forsupervisedlearning,featureselectioncan betrivial,i.e.,onlythefeaturesthatarerelatedtothegivenclusterla- belsaremaintained.Nevertheless,forunsupervisedlearning,thecluster labelsareunknown.Thus,ﬁndingtherelevantsubsetoffeaturesand clusteringthesubsetofthedatamustbeaccomplishedsimultaneously.

Assumingdtobetheinitialnumberoffeatures,anexhaustivesearch of2^𝒅possiblesubsetsneedtobeexamined,whichiscomputationallyex- pensive.Therefore,inthisstudyanalternativeGAbasedmethodwas appliedtodeterminetheoptimalsubsetoffeatures,aswellastheop- timalclusteringresults.Theaveragesilhouetteindexwasusedasthe optimizationcriteria.

Asoneofthemostpopularsearchheuristicmethods,GAhasbeen widelyusedingeneratingsolutionstooptimizationandsearchingprob- lems(YangandHonavar,1997;ShahamatandPouyan,2015).Diﬀerent fromsingle-statemethods(onlyonesolutionisevaluatedatatime),e.g., simulatedannealing,hillclimbing,etc.,GAisa“population” method thatmaintainsasetofsolutionsevolvingtowardanoptimalsolution.

Theevolutionusuallystartsfromapopulationofrandomlygenerated

(10)

Fig.9.FlowchartofGAforfeatureselection.In the𝑴-by-𝒅 matrix,eachrow representsacan- didatesolution, describinga subsetofselected features.Eachofthe𝒅bitsinarowrepresents whetherafeatureisselected(1)ordiscarded(0).

solutions.Ineachiteration(orsocalled“generation”),“survivorsolu- tions” withlargervaluesofoptimizationcriteria,i.e.,theaveragesilhouetteindex,areselectedtoformanewgenerationofsolutions.These survivalsolutionscanbegeneratedfromthecrossover,whichproduced newsolutionsbyrandomlycombiningtwocurrentsolutions,mutations, whichrandomlychangesnewsolutionswithasmallprobability,orfrom theinitialpopulation.Thenewgenerationofsolutionsisthenusedin thenextiterationofthealgorithm.Conventionally,thealgorithmter- minateswhenthebestsolutioncannotbeimprovedanyfurther.

Inthisstudy,anarrayof𝒅bitswasusedtorepresenttheselected subsetoffeaturesandthepopulationsizeisrepresentedusing𝑴.Each bitinthearrayindicatestheactivationstatusofonespeciﬁcfeature:1 indicatesselectedand0indicateddiscarded.Thecompleteprocedureof GAisdescribedbelow(Fig.9):

(1) Initialization:400candidatesolutionsweregeneratedbyrandomly setting1or0foreachbitinvectors.

(2) Crossover:twocandidatesolutionsAandBwererandomlyselected fromthecurrentpopulation.Avaluevbetween1and𝒅wasran- domlyselected.Thenanewsolutionwasformedbycombiningthe featurebits1tovfromAandfeaturebitsv+1to𝒅fromB.

(3) Mutation:foreachnewgeneratedsolution,amutationwasapplied byreversingbitsinthevectorwithaprobabilityof0.1.

(4) Evaluation:theclusteringmethodwasappliedoneachcandidateso- lution(i.e.,asubsetofselectedfeatures),andtheaveragesilhouette indexwascomputedforeachobtainedpartition.

(5) Selection:280solutionsresultinginhighaveragesilhouetteindex wereselectedalongwith120solutionsrandomlyselectedfromthe restofthesolution(toincreasethediversityofthepopulation).

(6) Iftheresultdidconverge,i.e.,theaveragesilhouetteindexofthebest solutioninthepopulationkeepincreases,weiteratedbacktostep2.

Otherwise,theclusteringresultwiththelargestaveragesilhouette indexandthecorrespondingselectedsubsetoffeaturesweresaved astheoutput.

(11)

Fig.10. Anatomicalanterior-posteriorsegmentation usedinthisstudy.(CoordinatesareinMNIspace).

2.6. Volumelevelclustering

In order to determine whether characterizing layer-specific mi- crocircuits using ultra-highfieldfMRIprovides any advantagesover conventionally computed voxel-level connectivity, the same cluster- ingprocedureenumeratedabovewasalsoperformedontheFCcom- puted between hippocampal and DAN/DMN voxels during encoding/retrieval tasks. Specifically, let 𝒀 ={𝒀1,⋯,𝒀_𝒊, ⋯,𝒀_𝑵} represent a set of 𝑵 objects, i.e., number of hippocampal voxels. 𝒀_𝑖= (𝑌𝑖1,𝑌𝑖2,⋯,𝑌𝑖𝑑)∈ℝ^𝑑,where𝒅equalstothenumberofFCfeaturescom- putedbetweenhippocampalvoxelsandvoxelsintheDAN/DMN,respectively.Subsequently,thesameclusteringandfeatureselectionpro- cesswasrepeatedonencodingandretrievaltasks,DANandDMNROIs, separately.

2.7. Comparisonwithanatomicalanterior-posteriorsegmentation

In this study, clusters were identiﬁed based on the DAN/DMN- hippocampalFCduringencoding/retrievaltasks,respectively.Thispro- videdafunctionalparcellationof thehippocampus,which wascom- paredwithanatomicallydelineatedanteriorandposteriorhippocampal segments.

Thehippocampuscanbe anatomicallyseparatedintohead,body, andtails,with theheadandbody beingconsideredasanterior,and thetailbeingconsideredastheposteriorpart.Variousmethods,such as landmark-based segmentation (Poppenk et al., 2013), percentile- basedaxissegmentation(Hackertetal.,2002;Greiciusetal.,2003),a Talairach/MNIcoordinate-basedsegmentation(Poppenketal.,2013), havebeenusedtodeﬁnetheanteriorandposteriorregionsofthehip- pocampus.WeemployedtheTalairach/MNIcoordinate-basedsegmentation(Poppenketal.,2013)whichchosey=−21inMNIspace(y=−20 inTalairachspace)astheborderbetweenanteriorandposteriorseg- mentationsasshowninFig.10.Thiscoordinatecorrespondstotheun- calapex,whichisconsideredastheendoftheposteriorportionofthe hippocampus(Destrieuxetal.,2013).

Let𝑨={𝑨1, 𝑨2, ⋯,𝑨_𝑴}representmanatomicalparcels,and𝑭= {𝑭1, 𝑭2,⋯, 𝑭_𝑴}denotenfunctionalparcelsidentifiedbyapplying clusteringmethodsonFCfeatures. Thesimilaritybetweenthese two parcellationswasthenquantifiedusingTorres’method(Torresetal., 2008).ThesimilaritymatrixforAandFisan𝒎 × 𝒏matrixdefinedas:

𝑺_𝑨_,_𝑭=

⎡⎢

⎢⎢

⎣ 𝑺11

𝑺_𝒊1

𝑺𝒎1

…

… 𝑺1𝒋

⋮ 𝑺_𝒊_𝒋

⋮ 𝑺𝒎𝒋

…

… 𝑺1𝒏

𝑺_𝒊_𝒏 𝑺𝒎𝒏

⎤⎥

⎥⎥

⎦

(8)

where𝑺_𝒊_𝒋=𝒊∕𝒖,whichisJaccard’sSimilarityCoeﬃcientwith𝒊being thesizeofintersectionand𝒖beingthesizeoftheunionofclustersets 𝑨_𝒊and𝑭_𝒋.Thesimilarityofparcellations𝑨and𝑭isthendeﬁnedas:

𝑺𝒊𝒎(𝑨,𝑭)=

∑𝒊≤𝒎,𝒋≤𝒏𝑺_𝒊_𝒋

𝒎𝒂𝒙(𝒎, 𝒏) (9)

FromEqs.(8)to(9), itcanbeseenthat0≤𝑺𝒊𝒎(𝑨,𝑭) ≤1,and 𝑺𝒊𝒎(𝑨, 𝑭)=1whentwoparcellationsareidentical.

The entireanalysispipeline proposedfor investigating functional diﬀerentiation of the hippocampus and layer-speciﬁc microcircuitry between the hippocampus and the DAN and DMN during encoding/retrievaltasksareillustratedinFig.2.

3. Results

Theoptimalvaluesofeachinputparameterdeterminedforthethree clustering methodsarepresented Tables 1and2for layer-level and volume-levelclustering,respectively.Usingeachclusteringmethod,the hippocampalvoxelswereclusteredintotwodiﬀerentfunctionalparcels basedontheirFCwith(1)layerVoftheDAN/DMNduringtheencod- ingtask,and(2)layerIIoftheDAN/DMNduringtheretrievaltask.This wastrueacrossmethods.Theobtainedclusterswerethenmappedback totheimagespaceandtheresultinghippocampalparcelswereoverlaid ontheanatomicalimageforthevisualization.

SimilarhippocampalparcellationswerediscoveredusingtheirFC withlayerVoftheDAN/DMNduringtheencodingtaskandlayerIIof theDAN/DMNduringtheretrievaltask,andtheaverageclustersimilar- itybetweenDPCandhierarchicalclustering,betweenDPCandOPTICS, andbetweenhierarchicalclusteringandOPTICS,were0.93,0.95,and 0.92,respectively.Forillustration,theclusteringresultsobtainedusing theDPCmethodareshowninFig.11.Theclusteringresultsobtained usinghierarchicalclusteringandOPTICSareshowninSupplementary InformationFigs.S.1andS.2.FromFig.11,itcanbeseenthatthehip- pocampus showedananterior-posterior gradientalongthelong-axis,

(12)

Table1

Estimatedoptimalvaluesofeachinputparameterinlayer-speciﬁcclustering.h:cutting height,ands:

reachabilitythreshold.

Method Name Parameter Left Hippocampus

DAN DMN

Encoding Layer V Retrieval Layer II Encoding Layer V Retrieval Layer II

DPC 𝝆 9.03 7.13 8.29 8.50

𝜹 3.15 2.82 3.14 3.04

Hierarchical h 1.16 1.15 1.15 1.15

OPTICS s 0.55 0.50 0.53 0.25

Fig.11. Clustersofhippocampalvoxelsdetermined(usingtheDPCmethod)basedontheirfunctionalconnectivitywith(1)layerVoftheDMN/DANduringthe encodingtask,and(2)layerIIoftheDMN/DANduringtheretrievaltask.(CoordinatesareinMNIspace).

Table2

Estimatedoptimalvaluesofeachinputparameterinvolume-levelclustering.h:

cuttingheight,ands:reachabilitythreshold.

Method

Name Parameter

Left Hippocampus

DAN DMN

Encoding Retrieval Encoding Retrieval

DPC 𝝆 7.59 6.38 5.61 7.93

𝜹 0.09 0.13 0.18 0.12

Hierarchical h 1.15 1.15 1.16 1.15

OPTICS s 0.39 0.06 0.002 0.03

whichisconsistentwiththeanatomicalanterior-posteriorsegmentation (Fig.10).

To quantitatively characterize the identiﬁed clusters, the cluster similaritybetweenourfunctionallyobtainedhippocampalparcelsand anatomicallydeﬁnedanterior-posteriorparcels(Fig.10)wascomputed.

Themeancorrelationbetweenhippocampalvoxelsandselectedvertices (usingGA-basedfeatureselectionmethod)inlayerVoftheDAN/DMN (duringtheencodingtask)andthelayerIIoftheDAN/DMN(during theretrievaltask)wasalsocomputedwithinanteriorandposteriorhip- pocampalregions(Table3).

Forthelefthippocampus,duringtheencodingtask,theaverageclus- tersimilarity,overdifferentclusteringmethods,betweenfunctionaland anatomicalparcellationsforlayerVoftheDANandDMNwas0.70and 0.70,respectively.TheabsolutecorrelationobservedbetweenlayerV oftheDANandthehippocampuswassignificantlylargerthanthatbe- tweenlayerVoftheDMNandthehippocampusasthep-valuesshown in Table4, whichisinlinewithoursecondhypothesis.Aone-tailed 2-samplet-testwasconductedtotestwhetherthecorrelationoflayer VoftheDANwiththeanteriorhippocampalregionswassignificantly largerthanthatwiththeposteriorhippocampalregions.Thep-values obtainedfordifferentclusteringmethodswereallcloseto1,whichdid notprovideevidencetosupportourthirdhypothesisforDANpart.

Duringtheretrievaltask,theaverageclustersimilaritybetweenfunc- tionalandanatomicalparcellationsforlayerIIoftheDANandDMNwas 0.70and0.77,respectively.Theabsolutecorrelationobtainedbetween layerIIoftheDMNandthehippocampuswassignificantlylargerthan thatbetweenlayerIIoftheDANandthehippocampusasthep-values showninTable5,whichwasconsistentwithoursecondhypothesis.In addition,thecorrelationbetweenlayerIIoftheDMNandtheposterior hippocampalregionswassignificantlylargerthanthatwiththeanterior hippocampalregions(one-tailed2-samplet-test;p<0.001fordifferent clusteringmethods),consideringthesignofthecorrelation.Thisresult wasinlinewithourthirdhypothesisforDMNpart.

(13)

Table3

Clustersimilaritybetweenfunctionalandanatomicalanterior-posteriorparcellationsusingdiﬀerentclusteringmethodsandmeancorrelationobtainedwithineach clusterontheleftsideofthehippocampus.

Clustering Method

Feature Type

Task Name

DAN DMN

Correlation

Sim.

Correlation

Sim.

Anterior Posterior Anterior Posterior

DPC Layer-

speciﬁc

Encoding Layer V) − 0.36 0.36 0.70 − 0.13 0.16 0.70

Retrieval (Layer II) 0.09 − 0.06 0.73 − 0.39 0.31 0.77

Volume- level

Encoding 0.20 0.04 0.69 0.02 0.12 0.69

Retrieval 0.30 0.33 0.74 − 0.01 − 0.02 0.74

Hier. Layer- speciﬁc

Encoding Layer V) − 0.36 0.35 0.70 − 0.13 0.16 0.70

Retrieval (Layer II) − 0.14 0.11 0.63 − 0.41 0.34 0.77

Volume-

level Encoding 0.48 − 0.17 0.68 0.51 − 0.18 0.69

Retrieval 0.08 − 0.05 0.72 0.01 0.02 0.74

OPTICS Layer- speciﬁc

Encoding Layer V) − 0.35 0.34 0.70 − 0.09 0.10 0.70

Retrieval (Layer II) − 0.14 0.07 0.73 − 0.36 0.27 0.77

Volume- level

Encoding 0.51 − 0.18 0.69 0.02 0.03 0.71

Retrieval − 0.33 0.41 0.57 0.21 − 0.09 0.74

Table4

ComparisonoftheabsolutecorrelationsobtainedbetweentheDANwiththehippocampusandtheDMNwiththehippocampusduringtheencodingtaskby conductingone-tailedtwo-samplet-test(𝑯0:correlationofDANandthehippocampus≤correlationoftheDMNandthehippocampus).

Feature Type

Clustering Method

Left Hippocampus Mean Absolute Correlation

P-value

DAN DMN

Anterior Posterior Anterior Posterior Anterior Posterior

Layer- speciﬁc (Layer V)

DPC 0.39 0.39 0.34 0.35 < 0.001 0.005

Hierarchical 0.38 0.39 0.34 0.35 0.004 0.009

OPTICS 0.37 0.38 0.34 0.35 0.002 0.007

Volume-

level DPC 0.51 0.18 0.51 0.18 0.500 0.500

Hierarchical 0.48 0.17 0.51 0.18 0.750 0.564

OPTICS 0.51 0.18 0.45 0.19 0.170 0.507

Table5

ComparisonoftheabsolutecorrelationsobtainedbetweentheDMNwiththehippocampusandtheDANwiththehippocampusduringtheretrievaltaskby conductingone-tailedtwo-samplet-test(𝑯0:correlationoftheDMNwiththehippocampus≤correlationoftheDANandthehippocampus).

Feature Type

Clustering Method

Left Hippocampus Mean Absolute Correlation

P-value

DAN DMN

Anterior Posterior Anterior Posterior Anterior Posterior

Layer- speciﬁc (Layer II)

DPC 0.37 0.37 0.51 0.41 < 0.001 0.002

Hierarchical 0.34 0.38 0.51 0.41 < 0.001 0.023

OPTICS 0.45 0.28 0.50 0.40 0.015 < 0.001

Volume- level

DPC 0.32 0.30 0.35 0.24 0.298 0.932

Hierarchical 0.23 0.27 0.31 0.28 0.053 0.444

OPTICS 0.33 0.41 0.32 0.26 0.540 0.973

Forvolumelevelanalysis,i.e.clusteringofhippocampalvoxelsbased ontheirFCwithvoxelsinDAN/DMNvolume(asopposedtolayerIIand layerVoftheDAN/DMNasbefore),similarparcellationswerediscov- eredusingdiﬀerentclusteringmethods,andtheaverageclustersimilar- itybetweenDPCandhierarchicalclustering,betweenDPCandOPTICS, andbetweenhierarchicalclusteringandOPTICS,were0.96,0.88,and 0.87,respectively.Forillustration,theclusteringresultsusingDPCare showninFig.12,andtheclusteringresultsusinghierarchicalcluster- ingandOPTICSareshowninsupplementaryinformation(Figs.S.3and S.4).

AsshowninFig.12,ananterior-posteriorgradientwasobtainedas well.Theclustersimilaritybetweenfunctionalandanatomicalparcel- lations,andthemeancorrelationwithinanteriorandposteriorregions werealsocomputedforvolumelevelclustering.AsshowninTable3, theaverageclustersimilarityfortheDAN(duringtheencodingtask) andtheDMN(during theretrievaltask)were0.69and0.74,respec-

tively,whichwerequalitativelylessthanthevalueobtainedusingthe layer-specificdata.Duringtheencodingtask,thecorrelationbetween DANandthehippocampuswasnotsignificantlygreaterthanthecor- relationbetween DMNandthehippocampus(p>0.05)forallthree clustering methods(Table 4),whereas during theretrievaltask, the correlationbetweenDMNandthehippocampuswasnotsignificantly greaterthanthecorrelationbetweenDANandthehippocampusforall threemethods(Table5).Thisresultwasdifferentfromtheresultob- tainedusinglayer-specificdata.Itwasalsocontradictorytooursecond hypothesis.Thecorrelationsbetweenanteriorandposteriorhippocam- palregionswithDAN/DMNforencoding/retrievaltaskswerealsocom- paredusingone-tailed2-samplet-test.Duringtheencodingtask,thep- valuesobtainedforDPC,hierarchical,andOTPICSwere0.019,<0.001, and0.02,respectively,whereasduringtheretrievaltask,thep-values obtained forthesethree methodswere1,1, and0.436,respectively.

This result provided evidence for our third hypothesis for only the

(14)

Fig.12. Clustersofhippocampalvoxelsdetermined(usingtheDPCmethod)basedontheirfunctionalconnectivityDMN/DANvolumeduringencodingandretrieval tasks.(CoordinatesareinMNIspace).

DAN partconsidering thesign of the correlation, not for the DMN part.

4. Discussion

In this work, we investigated functional differentiation of the hippocampus and the connectivity between the hippocampus and DAN/DMNregionsduringencoding/retrievaltasks.Giventhepredic- tionsobtainedfromtheHIPER/HERNETmodelandknownanatomical pathwaysbetweenthesenetworksandthehippocampus,wetestedfour hypothesesasdescribedintheIntroductionsection.Inordertodoso, weinvestigatedtheconnectivitybetweenthehippocampusandlayerV oftheDANduringtheencodingtask,andbetweenthehippocampusand layerIIoftheDMNduringtheretrievaltask.Theproposedfirst,second andfourthhypotheseswereconfirmed,whereas thethirdhypothesis waspartiallyconfirmedforDMNpart(consideringsignofthecorrelation),notforDANpart.Thediscussionofresultsisorganizedasfollows.

First,wediscusstheresultsobtainedbythevalidationoftheanterior- posteriordifferentiationofthehippocampususingtheFCbetween(1) hippocampalvoxelsandlayerVoftheDANduringtheencodingtask, and(2)hippocampalvoxelsandlayerIIoftheDMNduringtheretrieval task.Second,wediscusstheconnectivitybetweenthehippocampusand theDAN/DMNduringencoding/retrievaltasks.Third,wediscussthe layer-specificfunctionalmicrocircuitsbetweenthelayerVoftheDAN andthehippocampusduringtheencodingtask,andbetweenthelayerII oftheDMNandthehippocampusduringtheretrievaltask.Fourth,we discusstheimportanceofultra-highfieldfunctionalneuroimaginginde- velopingaccurateandrobustmodelsoffunctionalconnectivity.Finally, wediscussthefunctionaldifferentiationofthehippocampusalongthe long-axis.

4.1. Anterior-Posteriorfunctionaldiﬀerentiationofthehippocampus

Theideathattheanteriorandposteriorpartsofthehippocampus mayservediﬀerentfunctionsemergedhalfacenturyago(Nadel,1968).

More recently, Robinson et al.(2015) conducted meta-analyses and foundananterior-posteriorlong-axissegmentationonboththeleftand righthippocampi. Duarteetal. (2014) studiedfunctionalspecializa- tionofthehippocampususingfMRIandavirtualreality3Dparadigm.

Theyfoundafunctionaldichotomywherebytheanterior/posteriorhippocampus shows antagonistic processing patterns for spatial encod- ingandretrievalof3D spatialinformation. Princeetal.(2005) also provided evidence for an anterior-posterior parcellation that corre- spondedtoencoding/retrievalprocesses.Notjustinthehippocampus, butWangetal.(2016)showedthistobetrueinotherregionsoftheme- dialtemporallobe,suchastheperirhinalcortex(PRC)aswell.Specif- ically,theyshowedpreferentialfunctionalconnectionoftheanterior PRCwiththeDANandtheposteriorPRCwiththeDMN.

Inthisstudy,thefunctionalspecializationofthehippocampushas beeninvestigatedbyusingunsupervisedclusteringoffunctionalcon- nectivity betweenthehippocampus andlayerII/VoftheDAN/DMN duringencoding/retrievalprocesses,respectively.Ourresultsyieldeda consistentanterior-posteriorlong-axisparcellation,whichshowedhigh clustersimilarityacrossclusteringmethodsbasedoncompletelydiﬀer- entprinciples.Thisresultindicatesthatduringtheencoding/retrieval processes,thereisarobustanterior-posteriorfunctionaldiﬀerentiation alongthelong-axisofthehippocampus.

4.2. Layer-speciﬁcconnectivitybetweentheDAN/DMNandhippocampus

TheanatomicalpathwaysbetweenthehippocampusandDAN/DMN havebeenstudiedusinginvasiveanimalmodelsbymanyresearchers (Sugar et al., 2011; Thomson and Bannister, 2003; Thomson and Lamy, 2007). During an encoding task, deeper infra-granular layers (specifically, layer V) of the DAN project to the hippocampus, whereas the hippocampaloutput reaches primarily more superficial supra-granular layers (specifically,layer II) of theDMN. This layer- specificorganizationhasalsobeenreportedbymanyanatomicalstud- ies.However,ithasneverbeendirectlyinvestigatedusingconnectivity basednon-invasivemethodsinhumans.

(15)

This layer-speciﬁc pathway between the hippocampus and DAN/DMNisnotexclusivesincepathwayswhichoriginate/terminate inotherlayersofthecortexmayalsocontributetothehippocampal input or output. This is not surprising given the highly complex underlying microcircuitry and given that signals between any two brainregionscanrelayviamultiplestructuresincludingthethalamus.

However,thepathways between thehippocampus andlayers IIand V of the DMN/DAN seemto be thedominant ones based on prior invasiveanimalliterature(ShepardandGrillner,2010;Thomsonand Lamy,2007)andthereforewehavechosentotesttheminthisstudy.

We investigated the connectivity between the hippocampus and theDAN/DMN usingFC betweenlayerV ofthe DAN/DMNandthe hippocampus duringthe encoding task,andbetween layerII of the DAN/DMNandthehippocampusduringtheretrievaltask.Duringthe encodingtask,thecorrelationobservedbetweenlayerVoftheDANand thehippocampuswassignificantlylargerthanthatbetweenlayerVof theDMNandthehippocampus.Contrarily,duringtheretrievaltask,the correlationobtainedbetweenlayerIIoftheDMNandthehippocampus wassignificantlylargerthanthatbetweenlayerIIoftheDANandthe hippocampus.Together,thesefindingsprovideevidencethattheoutput fromlayerVoftheDANlikelydrivesthehippocampusduringencoding, whereasduringretrieval,theoutputfromthehippocampusislikelyre- layedasinputstolayerIIoftheDMN.GiventhatPearson’scorrelation isadirectionlessquantity,ourinferenceofdirectionalityisindirectat best,basedonknowndirectionalityoftheunderlyinganatomicalpro- jections.

4.3. Layer-speciﬁcfunctionalpathwaybetweentheDAN/DMNand hippocampus

Considering the directionality of signal projectionduring encoding/retrievalprocessesandtheHERNETmodel,theinformationflow betweendifferentlayersoftheDAN/DMNanddifferentregionsofthe hippocampusshouldfollowthepattern:layerVoftheDAN→anterior hippocampus→posteriorhippocampus→layerIIoftheDMN.Thus, wehypothesizedthatduringanencodingtask,layerVoftheDANmush showstrongercorrelationwithanteriorhippocampalregionsthanwith posteriorhippocampalregions,whereasduringretrievaltask,layerIIof theDMNmustexhibitstrongercorrelationwithposteriorhippocampal regionsthanwithanteriorhippocampalregions.Ourresultsprovidepar- tialsupportforthishypothesis,i.e.,onlytheDMNpartwasconfirmed duringretrievaltaskconsideringthesignofthecorrelation(negativein anterior,positiveinposterior),butnotforDANpartduringencoding task.However,ifweonlyconsiderthemagnitudeofthecorrelationig- noringitssign,thethirdhypothesiswasnottrueforboththeDANand DMNpart.

4.4. High-resolutionfunctionalimaging:layervsvolumedata

Recent advances in ultra-high field fMRI have provided a non- invasivewayofinvestigatingcorticalcolumns.Thistechniqueprovides severaladvantagesoverconventionalfieldstrengths,e.g.,improvedspa- tialresolution,increasedsignaltonoiseratio,etc.Moreimportantly, thistechniquemakesitfeasibletoexaminelayer-specificbrainactiva- tionacrossdifferentbrainareas.Severalrecentstudieshaveshowedthat investigatingchangesinfMRIactivationasafunctionoflaminardepth canleadtomorepreciseresults(Olmanetal.,2012;Koketal.,2016).In thisstudy,weinvestigatedfunctionaldifferentiationofthehippocam- pusalongthelong-axisusingunsupervisedclusteringoflayer-specific functionalconnectivitybetweenthehippocampusandtheDAN/DMN regions. The same clustering processwas also applied on the func- tionalconnectivitybetweenthehippocampusandtheDAN/DMNvol- ume.Asweexpected,thelayer-specificdataledtomoredefinitivere- sults,i.e.,theproposed first,second andforth hypotheseswerecon- firmed,whereasthethirdhypothesiswaspartiallyconfirmedforDMN part(consideringsignofthecorrelation),notforDANpart.However,

usingvolume-leveldata,onlythefirsthypothesiswasconfirmedandthe thirdhypothesis(onlyforDANpartwithsignofthecorrelationbeing considered)waspartiallyconfirmed.Therefore,itisimportanttonote therelevanceofhigh-resolutionfunctionalneuroimaging asthefield progressestowarddevelopingmoreaccurateandrobustnetworkmod- elsofbrainfunctioningeneralandhippocampalfunctioninspecific.

4.5. Long-axisdiﬀerentiationofthehippocampus

Recent evidencehassuggested that there isan anterior-posterior functionaldiﬀerentiation ofthehippocampus alongthelong-axis. In Lepage’sHIPERmodel(Lepageetal.,1998)andmorerecentlyKim’s HERNET mode(Kim,2015),theanterior andposterior hippocampus are posited to be more associated with encoding and retrieval processes,respectively.Thisencoding/retrievaldichotomyfacedconﬂict- ingevidencefromsomemeta-analysesandstudies.SchacterandWagner (1999)revieweddatafromdiversefMRIstudiesandobservedthatboth theanteriorandposteriorhippocampalregions wereassociatedwith theencodingactivation.KumaranandMaguire(2006),Poppenketal.

(2008),andZweynertetal.(2011)foundthatmostencodingstudies usednovelstimuli,whichwereassociatedwiththeanteriorhippocam- pus.Ontheother hand,Poppenketal., 2010a2010b) observedthat familiarstimuliwereassociatedwiththeposteriorhippocampus and superiorsourcememory.Together,thesestudiessuggestthattheen- coding/retrievaldiﬀerentiationcannotrelyonﬁndingsthatlinkthean- terior/posteriorhippocampustonovel/familiarstimuli.

Althoughtheencoding/retrievaldichotomyis predominantinhu- mans, many alternative specializations have also been proposed. A motivational processing model has been proposed with the anterior hippocampus being mainly engaged in “hot” processing (emotion/motivation), whereas the posterior hippocampus being mainly associated with “cold” processing (cognition) (Murty et al., 2011).

Robinson etal.(2015)(2016). alsofoundsupportfor thishypothesis basedonhippocampalparcellationsobtainedfrommeta-analyses,rest- ingstatefMRIconnectivityanddiﬀusiontensorimaging().However, Wolosinetal.(2012)foundthattheposteriorhippocampusalsocon- tributedtonegativeemotionalmemory.Someotherstudieshavepro- posedthattheposteriorhippocampusisespeciallyimportantforspa- tialprocessing,whereastheanteriorhippocampusmaybeimportantfor episodicmemoryorotherfunctions(Ryanetal.,2010;Hirshhornetal., 2012).Thismodelwasundercutbyevidencethattheanteriorhippocam- pusalsoplaysaspatialrole(WoollettandMaguire,2012).

In this study, our results provided partial support for HIPER/HERNET model. The strength of evidence in favor of the HERNET model was superior with layer-speciﬁc data compared to conventionalvolumedata.Aconsistentanteriortoposteriorlong-axis segmentationwasfoundduringencoding/retrievaltasks.Wealsofound that duringanencodingtask,thehippocampus wasmore correlated with the layer V of the DAN, whereas during a retrieval task, the hippocampuswasmoreassociatedwithlayerIIoftheDMN.However, wedidnotﬁndsupportforstrongercorrelationoflayerVoftheDAN andthe anteriorportionsduring theencoding task.Our resultsalso did not support theprediction of stronger correlation of layerII of theDMNandtheposteriorhippocampalsegmentsduringtheretrieval task (ignoring the sign of the correlation). Together, these results suggest that the underlying neurophysiology of the hippocampus and its interactionwith theneocortex under variousneurocognitive contextsisfarmorecomplexthanrelativelysimplisticmodelscurrently available.Ourstudydemonstratesthatitwilltakebetterqualitydata intermsofspatial-temporalresolutioninordertobuildbettermodels ofhippocampalfunctionandspecialization.

4.6. Limitations

Thepresentstudyhasafewlimitationsthatmaybeaddressedin futureresearch.First,inthisstudy,thecorticallayersweredelineated