Sub-cellular population imaging tools reveal stable apical dendrites in hippocampal area CA3

Anatomically segregated apical and basal dendrites of pyramidal neurons receive functionally distinct inputs, but it is unknown if this results in compartment-level functional diversity during behavior. Here we imaged calcium signals from apical dendrites, soma, and basal dendrites of pyramidal neurons in area CA3 of mouse hippocampus during head-fixed navigation. To examine dendritic population activity, we developed computational tools to identify dendritic regions of interest and extract accurate fluorescence traces. We identified robust spatial tuning in apical and basal dendrites, similar to soma, though basal dendrites had reduced activity rates and place field widths. Across days, apical dendrites were more stable than soma or basal dendrites, resulting in better decoding of the animal’s position. These population-level dendritic differences may reflect functionally distinct input streams leading to different dendritic computations in CA3. These tools will facilitate future studies of signal transformations between cellular compartments and their relation to behavior.


34
Neurons receive distributed inputs onto their dendritic trees. These inputs are often arranged 35 in a stratified manner, in which synapses from specific afferent brain areas arrive onto different 36 dendritic compartments of a pyramidal cell. Broadly, dendrites of a prototypical pyramidal 37 neuron are divided into apical and basal compartments. Beyond differences in the inputs they 38 receive, apical and basal dendrites may differ in their molecular composition, structural and 39 biophysical properties, and function 1 . It is known that dendrites support non-linear processing 40 in the form of dendritic spikes 2-8 and other ionic conductances 9,10 , which can tremendously 41 enhance the computational capacity of individual neurons [11][12][13][14][15] and support ensemble coding 42 and plasticity [16][17][18] . However, it is unclear if the population activity of apical and basal dendrites 43 differ in terms of content or stability across time.

45
Dendritic properties and computations have primarily been studied using single neuron 46 recordings in brain slices or with calcium imaging 19,20 in sparsely labeled in vivo preparations 21-47 25 . While these approaches allow faithful tracking of the activity of individual dendrites, they 48 restrict experimental throughput, are often confined to only one dendritic compartment (apical 49 or basal) and provide limited insight into the population dynamics of large numbers of 50 dendrites. However, in vivo two-photon calcium imaging does allow simultaneous recording 51 from large populations of neurons over long time periods. With stable implants and genetically 52 expressed Ca 2+ indicators, the same field of view can be revisited across several days. This 53 allows experimenters to track the activity of the same neurons over time and quantify the 54 stability or flexibility of response properties. Densely labeled preparations can yield hundreds to 55 thousands of neurons, which enables population vector decoding of behavior 26,27 and can yield 56 insights into network properties and circuit-level function. Though calcium imaging provides 57 subcellular resolution, thus far studies have not fully utilized the ability to track large 58 populations of sub-cellular processes like dendrites over time.

60
A large factor limiting progress in understanding dendritic function in vivo is technical, 61 particularly when recording from dense populations. Regions of interest (ROIs) belonging to 62 individual neurons or sub-compartments of neurons must be identified, which can be tedious 63 and time-consuming when done manually. While several software suites have been developed 64 to automate this process when recording from fields of view composed of cell bodies 28-31 , 65 dendrites' diverse shapes makes morphology-based methods of dendritic ROI detection less 66 reliable. To relate neural activity to behavior, the time-varying fluorescence values from each 67 ROI must be estimated and transient periods of elevated activity identified. These calcium 68 transients, or "events," can then be used to quantify tuning properties. A simple approach of 69 averaging values within ROIs is prone to contamination from overlapping ROIs, more so for 70 dendrites compared to soma due to the smaller number of pixels per ROI. Human screening of 71 detected transients can address these problems, but the amount of necessary manual labor 72 and time will scale with the number of dendrites in a field of view, making manual checking 73 infeasible. Overcoming such difficulties will tremendously expand the field's ability to fully 74 investigate sub-compartment level population activity and behaviorally relevant functional 75 properties.

77
To address these challenges, we developed an automated detection algorithm flexible enough 78 to identify dendritic and somatic ROIs in dense datasets. The algorithm identifies initial ROI 79 estimates using minimal morphological assumptions and refines them using constrained non-80 negative matrix factorization (CNMF) 28,32 . It then screens putative calcium transients with an 81 additional goodness-of-fit measure to eliminate spurious activity from undetected ROIs. To 82 demonstrate the efficacy and utility of this approach, we applied it to several sets of data 83 recorded from area CA3 of the mouse hippocampus 33 , an area containing "place cells" crucial 84 for spatial navigation. Using this preparation, we characterized the spatial tuning properties of 85 soma, apical dendrites, and basal dendrites during a simple navigation task. We find that while 86 basal dendrites show lower event rates than soma, apical and basal dendrites have comparable 87 spatial selectivity. However, by comparing activity across days in a familiar environment, we 88 show that the apical dendrites have higher day-to-day stability than either basal dendrites or 89 soma. This results in a better population decoding accuracy of position across days when 90 recording from CA3 apical dendrites, indicating a functional divergence between sub-cellular 91 compartments of hippocampal pyramidal neurons for long-term spatial representation. 92 93

94
Dense dendritic fields of view are highly overlapping 95 To generate dendritic imaging data sets for developing and testing our algorithm, we expressed 96 GCaMP6f or 7b in mice using a range of virus titers to control the labeling density ( Fig. 1,  97 Extended Data Fig. 1, see Methods) 19,21,[33][34][35] . This allowed us to express GCaMP extremely 98 sparsely to validate the basic mechanism of our algorithm, and then apply the same approaches 99 to densely labeled datasets suitable for population level analysis. At low titers, single neurons 100 were labeled in a single field of view ( Fig. 1b). At higher titers the density of neurons and 101 dendrites was so high that completely isolated soma or dendrites were quite rare. Most ROIs 102 overlapped with at least one other ROI, and most pixels were covered by more than one ROI. To 103 quantify the degree of overlap, human experts manually labeled these datasets, identifying 104 groups of contiguous pixels with coordinated activity (Fig. 1c, see Methods). In the 8 mice with 105 the highest titer virus injected, an average of 1300 ROIs/FOV were identified, covering 35% of 106 the field of view (Fig. 1d). Manual labeling of these datasets took approximately 50 hours per 107 field of view. On average, any given ROI at least partially overlapped with 3.6 other ROIs, and 108 only 40% of pixels belonged to a single ROI, demonstrating a high degree of spatial overlap.

110
Dendritic NMF algorithm enables efficient ROI extraction 111 The effort required to manually curate these datasets is massive, and the large amount of 112 spatial overlap complicates efforts to accurately estimate activity. Hence, we sought a method 113 to identify ROIs in an automated and unbiased way. CNMF 28,32 has previously been used to 114 identify ROIs and activities from partially overlapping neurons. The goal of CNMF is to identify 115 ROIs and corresponding fluorescence traces that best approximate the original image stack. The 116 objective function is the following: 117 argmin , ‖ − ‖ 2 + ‖ ‖ 2 + ‖ ( )‖ 2 ; ≥ 0, ≥ 0 Here, Y represents the raw calcium imaging sequence, A represents the spatial footprints 118 (ROIs), and C represents the fluorescence traces for each ROI. Background fluorescence from 119 out-of-plane neuropil is accounted for by additional columns and rows in A and C, respectively. 120 The η and β terms are penalties to encourage sparsity in both the ROIs and fluorescence traces. 121 This equation is solved by iteratively solving for A and C, fixing one while optimizing the other 122 (see Methods). This process can identify large neural populations but has not been 123 systematically applied to densely labeled dendritic preparations.

125
We developed a pipeline of unique initialization and refinement steps to be amenable to 126 detecting soma and dendrites, which we call "dendritic NMF," or d-NMF. We utilized temporal 127 downsampling and patch-wise processing to reduce the memory and time requirements of the 128 algorithm to enable processing on standard desktop computers (Fig. 2a, see Methods).

129
Processing the image in patches had the added benefit of incorporating multiple background 130 sources for each image patch, to better capture regional differences in neuropil signal.

131
ROI initialization was performed by identifying contiguous regions of pixels active above a 132 threshold (see Methods). Spatial footprints and temporal traces were then estimated by 133 iterating through CNMF until convergence (Fig. 2b). ROIs were then thresholded and split into 134 connected components (see Methods). These final morphological operations ensured that all 135 identified ROIs were relatively compact and easily identifiable as single dendritic branches or 136 soma. 137

138
The raw output of d-NMF contains any set of contiguous pixels that had simultaneously high 139 activity at least once, which can result in identification of ROIs that do not clearly correspond to 140 neural processes (false positives). In many brain regions, including the hippocampus, excitatory 141 neural activity is characterized by long periods of quiescence punctuated by short periods of 142 activity, resulting in a high skewness value. Hence, we quantified the skewness of the extracted 143 calcium traces for each ROI (Fig. 2c) 173 3c). The percentage of the FOV covered by ROIs (coverage) using d-NMF (27%) was also higher 174 than CNMF (13%) and suite2p (21%) (Fig. 3d, Right). In particular, CaImAn missed several ROIs 175 in dense regions of the FOV. While suite2p performed comparably with d-NMF in sparse 176 preparations (Fig. 3a, Left), dense preparations were more completely labeled using d-NMF 177 (Fig. 3a, Right). d-NMF performed well on both sparse and dense datasets with a single set of 178 parameters, whereas other methods may need careful parameter tuning depending on the 179 density of the preparation.

181
A fitness trace to refine activity estimates and eliminate cross contamination 182 Densely labeled datasets present added difficulties in estimating the fluorescence activity of 183 each ROI due to overlapping sources. Given a complete set of ROIs, CNMF accurately de-mixes 184 these signals and properly assigns activity to each ROI. However, despite the best efforts of 185 manual labeling or automated labeling, some ROIs will go undetected. This will lead to 186 "contamination" of signals from identified ROIs by unidentified ones. Such contamination can 187 affect estimates of activity rates, tuning properties, stability, or correlation with related ROIs. 188 While simply setting a high threshold for transient detection may reliably screen out false 189 transients, this runs the risk of excluding low-amplitude calcium transients which may represent 190 different firing modes.

192
To both evaluate the degree of contamination and to provide an additional tool to correct this, 193 we developed the concept of a "Fitness Trace" for each ROI (Fig. 4a), defined as the frame-by-194 frame correlation between the spatial footprint of that ROI with the video (see Methods). The 195 intuition behind this trace is that true activity in a given ROI occurs when all pixels show 196 elevated activity, not just some. Significant calcium transients were detected by thresholding 197 the Fitness trace along with the original calcium trace (Fig. 4b). Overall classifier performance 198 was evaluated by computing the Jaccard index (the size of the intersection divided by the size 199 of the union) 36 between detected events and manually labeled true events (see Methods).

200
Across 9 different FOVs there was a relatively consistent optimal parameter range, at a ΔF/F 201 threshold of 2.3z and a Fitness threshold of 0.28.

203
We compared detection accuracy and the rate of detected transients from our optimized 204 Fitness method to those of other event detection methods (Fig. 4c). A widely used approach is 205 to set a threshold of 2 standard deviations above the mean on the raw ΔF/F traces, with no de-206 mixing (with an additional false positive rejection criterion, see Methods). We call this method 207 "2z." We call our method using optimal values for the ΔF/F threshold and Fitness threshold 208 "D Fit ." Event classification performance was significantly higher using D Fit (0.65) compared to 2z 209 (0.48). Other detection methods incorporating different elements of these two approaches 210 achieved intermediate results (Extended Data Fig. 4). The use of the Fitness Trace also was 211 comparable to or exceeded other related measures 37 (Extended Data Fig. 4e-h). More striking, 212 the estimated calcium transient rate differed by nearly a factor of 2 when comparing the two 213 methods (Fig. 4c). This demonstrates that accurate event detection and screening is a vital step 214 in processing data from densely labeled FOVs.

216
Spatial coding properties of dendrites in CA3 217 To demonstrate the utility of d-NMF, we applied it to dense fields of view obtained in area CA3 218 of the mouse hippocampus, as mice performed a random foraging task on a textured belt ( Fig.  219 5, see Methods). Somata of pyramidal cells in area CA3 demonstrate place coding, where they 220 have reliable elevated activity levels in restricted regions of space. Using d-NMF, across 8 mice 221 we identified 794 total somatic ROIs, 2910 apical dendritic ROIs, and 4577 basal dendritic ROIs. 222 This allowed us to examine and compare the activity and spatial tuning properties of 223 populations of apical and basal dendrites.

225
While spatial tuning in CA3 PN somata has been well described through electrophysiology 226 approaches 38-40 and more recently with Ca 2+ imaging 41,42 , we know little about how reliably CA3 227 dendrites encode place and their spatial tuning properties. Using our automated approaches, 228 we identified spatially tuned CA3 dendrites ("place dendrites") both in the apical and basal 229 compartments ( Fig. 5a-d). Overall activity rates and spatial tuning properties were quite similar 230 between soma and dendrites, with basal dendrites having a lower event rate than soma ( Fig.  231 5e). A similar percentage of ROIs were significantly spatially tuned (Fig. 5f). The place field 232 width of tuned apical dendrites and soma were largely similar (Fig. 5g), while basal dendrite 233 place field widths were significantly smaller than soma. Other spatial tuning properties such as 234 the number of place fields, information content, and sparsity were also similar across 235 compartments (Extended Data Fig. 5). Basic spatial tuning characteristics were not 236 systematically related to the ROI size (Extended Data Figs. 6, 7) or the number of overlapping 237 ROIs (Extended Data Fig. 8).

239
The above tuning properties were computed using events detected using the D Fit procedure. 240 We re-computed spatial tuning properties using the more traditional 2z detection method, 241 which overestimated the event rate, place field width, and number of place fields (Extended 242 Data Fig. 5), and reported a significantly lower firing rate in the soma compared to dendrites 243 (Fig. 5e). This underscores the importance of proper event detection and illustrates that 244 scientific findings can be affected by the choice of analysis method.

246
Apical dendrites are more stable across days and provide better place decoding 247 An advantage of optical recording techniques is the ability to track the same ROIs across 248 multiple days. Thus, we quantified the short-term (within-day) and long-term (across-day) 249 stability of apical dendrites, soma, and basal dendrites in a familiar environment (Fig. 6). 250 Within-day stability was not significantly different across compartments whether quantified 251 using either Tuning Curve (TC) correlation or Population Vector (PV) correlation (see Methods). 252 However, apical dendrites were more stable across days compared to soma (either using TC or 253 PV correlation) and basal dendrites (using PV correlation). The differences in across-day stability 254 were not due to fields of view distorting across days (Extended Data Fig. 9), and could not be 255 attributed to differences in ROI size or overlap (Extended Data Fig. 10).

257
Because apical dendrites showed better across-day stability compared to soma and basal 258 dendrites, we predicted that apical dendrites would also outperform the other compartments 259 in population vector decoding. Thus, we constructed population vector decoders using only 260 apical dendrites, soma, or basal dendrites, and tested their ability to decode position within a 261 session ( Fig. 7a) or across sessions (Fig. 7b). There was no significant difference between 262 decoding accuracy within-day. However, apical dendrites had significantly lower decoding error 263 compared to basal dendrites across days (Fig. 7c). Population vector decoding is known to be 264 highly dependent on the number of sources used to train the model, so we repeated the 265 analysis using random subsets of data. Consistently, we found no differences across groups for 266 within-day decoding (Fig. 7d), but a significant effect of compartment group on across-day 267 decoding, with apical dendrites showing significantly better decoding accuracy compared to 268 basal dendrites (Fig. 7e).

271
In summary, we developed an algorithm and toolkit that facilitates analysis of densely labeled 272 fields of view containing both apical and basal dendrites and cell bodies ( Fig. 1), enabling large-273 throughput in vivo investigation of dendritic activity. The algorithm utilizes the well-validated 274 mathematical framework of CNMF while providing automated initialization flexible enough to 275 identify neural processes of varied shape and size (Figs. 2,3). The addition of the Fitness Trace 276 post-hoc helps screen out false calcium transients to more accurately estimate dendritic tuning 277 properties (Fig. 4). Using these tools, we demonstrate spatial tuning in apical and basal 278 dendrites of pyramidal neurons in CA3 (Fig. 5). While the proportion of tuned ROIs and within-279 day tuning properties were largely comparable across all compartments, basal dendrites 280 showed lower event rates and place field width than soma. By comparing longer term 281 population dynamics across days, we demonstrate that apical dendrites are more stable than 282 soma or basal dendrites across-days (Fig. 6). This leads to better accuracy when constructing a 283 population vector decoder to predict the animal's position (Fig. 7). These population-level 284 analyses would not be possible with sparsely labeled preparations, and we were only able to 285 take advantage of our densely labeled preparations using the analytical tools we have 286 developed.

288
Densely labeled dendritic fields of view present a number of challenges for existing automated 289 detection methods, particularly the initialization step. Initialization is often based on 290 morphology, with parameters to be set for the expected diameter of an ROI. This assumes a 291 uniform size of circular ROIs throughout the FOV, which is not the case for dendrites. Random 292 initialization does not suffer from the drawbacks of morphology-based approaches, but there is 293 no guarantee of complete labeling or consistent labeling across iterations. Additionally, the 294 number of components must be specified for initialization, which may vary over the field of 295 view. One approach is to over-specify the model and then discard unusable components 296 afterwards 32 , but this comes at a higher computational cost.

298
Another drawback other automated methods have when dealing with overlapping dendrites is 299 a reliance on correlation. Seed pixels for initialization can be chosen from the peaks of the 300 correlation map, constructed by computing the average correlation of a pixel with its 301 neighbors 43 . However, thin dendritic segments a few pixels wide may skew towards low 302 correlation values and be systematically missed compared to wider somatic ROIs. Furthermore, 303 when so many dendrites cross a given region ( Fig. 1), the activity of those pixels will be a 304 combination of all dendrites, which may obscure correlations within individual branches. d-305 NMF does not use correlation to initialize ROIs per se, but rather detects instances of 306 synchronous activity across neighboring pixels to define initial spatial estimates (Fig. 2). This 307 approach works in our densely labeled datasets because soma and dendrites tended to be 308 sparsely active in time, such that only a few ROIs are active on any given frame. This mirrors the 309 method by which manual segmentation is done, and allows identification of ROIs that may 310 overlap with many others (Fig. 3).

312
Different aspects of d-NMF could be used in conjunction with pre-existing methods. Other 313 segmentation software could be initialized with the initial ROI cores or the final ROIs obtained 314 after de-mixing. The Fitness Trace (Fig. 4) could be used to screen any activity trace regardless 315 of the method used to obtain it. Ultimately, different tools may be better suited to different 316 data sets, so it is worthwhile for investigators to be able to combine the strengths of all 317 available approaches.

319
Large population recordings allow the use of fewer animals to record from a given number of 320 neurons, saving on time and resources. For characterizing single dendrite properties, this is 321 particularly advantageous if a relatively small percentage of units show tuning to a given 322 stimulus or behavioral feature. This is the case in CA3, where only ~20% of dendrites or soma 323 show spatial selectivity (Fig. 5). Methods of inducing selectivity through artificial means 44,45 can 324 somewhat circumvent this, but the effects of these manipulations may be limited or 325 constrained 46 . Larger population data not only adds more power to the analysis but also 326 provide insights into functional heterogeneity and dynamics [47][48][49] . Simultaneous recording of 327 multiple dendrites in the same animal provides within-subject controls on the effect of 328 behavior on dendritic activity, and large enough population recording enables analysis of 329 population-level stability (Fig. 6) and decoding of behavior ( Fig. 7)  Our results demonstrate that both apical and basal dendrites show spatial selectivity on par 342 with nearby soma (Fig. 5), and the apical dendritic signal is more stable than the somatic or 343 basal dendritic signal (Fig. 6). This suggests an interesting arrangement where the input arriving 344 on to apical dendrites is stable across days yet the basal input changes from day to day. This 345 may be related to the varying functionality of CA3 as performing pattern completion or pattern 346 separation 61-65 . Such a switch may be the result of the network biasing propagation of inputs 347 from apical or basal dendrites. A similar mechanism may be at work in the neocortex, where 348 anatomy 66-68 and theory [69][70][71] suggests that bottom-up input arrives onto the basal dendrites 349 while top-down input arrives on to the apical dendrites, providing a spatial segregation 350 matching the functional segregation.

352
What could drive the difference in stability across dendritic compartments? These properties 353 may be inherited from their upstream sources. Both entorhinal cortex and dentate gyrus, which 354 target the apical dendrites of CA3, show a high degree of stability 42,60 . In contrast, the basal 355 dendrites primarily receive recurrent excitation from CA3 soma, which themselves are largely 356 unstable across days. Apical dendrites may also be more stable due to integrating signals from 357 several different inputs. Conjunctive activity of entorhinal cortex, dentate gyrus, and recurrent 358 collaterals may trigger supra-linearities such as dendritic spikes 7,18 , triggering long-term 359 synaptic plasticity 16,18,44,72,73 which could contribute to increased stability 21 . Differences in task 360 demands also affect representational stability of soma in CA1 50 , so the stability of dendritic 361 compartments under different behavioral tasks should be explored in the future.

363
It is also important to consider the contribution of forward-propagating versus back-364 propagating signals in defining dendritic activity. Basal dendrites are electrotonically close to 365 the soma, so their activity may be largely back-propagated. However, the lower calcium 366 transient rates in the basal dendrites suggest this propagation may not be entirely faithful. 367 Basal dendrites had smaller place fields than soma, which may be a result of decreased event 368 rate or could reflect more selective firing within a single session. Distal apical dendrites may 369 enjoy a larger degree of independence because back-propagating action potentials often decay 370 with distance 74 . Using calcium signals as an indicator of across-compartment correlation should 371 be done with caution though, as calcium signals may reflect large bursts of activity more 372 robustly than single action potentials 19 . Future studies using voltage-sensitive indicators or 373 more sensitive calcium indicators in large populations should provide higher temporal 374 resolution into these questions.

376
This study examined the activity of ensembles of apical dendrites, soma, and basal dendrites 377 separately. A more in-depth investigation of input-output transformations at a sub-cellular level 378 would require simultaneous measurement of connected dendrites and soma. By comparing the 379 tuning properties of the apical dendrite, soma, and basal dendrite of the same neuron within 380 densely labeled fields of view, one could learn about branch-specific 6,22,25,75 versus branch 381 prevalent activity 21,23,76 and soma-dendrite coupling dynamics while sampling large populations.

382
Such information has so far been limited to very sparse preparations, which may overlook 383 cellular heterogeneity within regions 47-49 .

385
Our approach already merges together ROIs of neighboring dendrites and soma which have 386 substantial spatial overlap and highly correlated activity. Such ROIs obtained from merging the 387 soma and dendrites indicate a high degree of coupling between the soma and its proximal 388 dendritic branches, as would be expected from the diffusion of calcium over short distances.

389
But the possibility of dendrites exhibiting semi-independent activity 6,22,25,75 from their cell 390 bodies poses a challenge in linking ROIs based on correlations alone. Furthermore, many distal 391 dendrites are connected to soma that are not in the focal plane, which could only be resolved 392 with multi-planar imaging. These issues of connecting dendrites with parent soma will be 393 important problems to address in future studies.

395
We developed and tested these tools primarily using viral expression of GCaMP6f. In principle, 396 our approaches should work with different genetically encoded calcium indicators 77 or voltage 397 indicators 78,79 . These different indicators have varying signal-to-noise ratios and kinetics, so a 398 reasonable parameter search may be necessary to obtain optimal performance of our tools. 399 The rapid kinetics of voltage indicators with respect to the speed of propagation within 400 dendrites may also result in the identification of finer-grained ROIs, as an assumption of our 401 mathematical model is that the signal in a given ROI is uniform.

403
The methods presented here enable investigation into dendrites and soma, but are broadly 404 applicable to any cell or cell process exhibiting dynamic fluorescence activity, including 405 axons 18,80,81 and astrocytes 82,83 . Any contiguous region of interest is identifiable using our core-406 finding approach. In our datasets we simultaneously detect soma, dendrites, and even putative 407 axons without needing to change parameters. Characterizing the activity of afferent axons to a 408 field of view is another tool in the circuit-mapper's toolbox, providing a measure of the input to 409 a brain region. Clever intersectional genetic approaches should be able to provide access to 410 axonal input, dendritic processing, and somatic output in a single preparation, and our tools will 411 facilitate the rigorous analysis of such datasets.

635
Animals 636 All experiments were conducted in accordance with the National Institute of Health guidelines 637 and with the approval of the New York University School of Medicine Institutional Animal Care 638 and Use Committee (IACUC). Imaging experiments used C57BL/6J mice as well as Ai9, Ai14, SST-639 Cre, and PV-Cre transgenic mice on a C57BL/6J background, from both sexes, 15-25 weeks old. 640 641 Cranial window surgery 642 Surgery procedures were similar as described previously 8 . Mice were anesthetized using 643 isoflurane (1.5-2.5%) and a 0.5 mm hole was drilled in the skull above dorsal area CA3 of 644 hippocampus (1.6mm lateral, and 1.4 mm caudal of Bregma). For densely expressing 645 preparations, injections of AAV1.CamKII. GCaMP6f.WPRE.SV40 (commercially generated at 646 Penn Vector Core, titer: 2.76x10 13 GC/ml, 23 nL per site) were made at two sites, two depth 647 levels each (1.5 mm lateral, 1.3 mm caudal of Bregma; 1.7 mm lateral, 1.5 mm caudal of 648 Bregma; depths of 1.8 and 2 mm below the dural surface). For sparser expression, we used a 649 mixture of AAV1.CamKII-Cre (commercially generated at Penn Vector Core, diluted with ACSF) 650 and AAV1.Syn.Flex.GCaMP6f.WPRE.SV40 or AAV1.Syn.Flex.jGCaMP7b.WPRE (commercially 651 generated at Penn Vector Core) at 23 nL per site, at different ratios 8 (see Fig. S1 for specific 652 titers). After injection, a 3 mm craniotomy was made with the injection site at the center. The 653 skull was removed, and a vacuum system was used to gently remove the overlying cortex and 654 external capsule. Ice-cold ACSF was used to irrigate the area throughout the duration of the 655 procedure. A cranial window (3 mm diameter, 1.7 mm length stainless steel cannula attached 656 to 3 mm diameter glass coverslip) was then implanted over the area. The window was sealed to 657 the skull using Vetbond, and a custom designed 3D-printed plastic headpost was cemented 658 over the skull.

660
Two-photon imaging 661 We used the same imaging system as described previously 8,50,84,85 . Mice were head-fixed on a 662 treadmill belt and trained to run for 5% sucrose water delivered at random locations. In vivo 663 two-photon imaging was performed using a dual galvanometric and resonant laser scanning 664 two-photon microscope (Ultima, Bruker), coupled to a tunable Ti:Sapphire laser (MaiTai eHP 665 DeepSee, Spectraphysics) pulsed at a 80 MHz repetition rates and <70 fs pulse width along with 666 dispersion compensation. GCaMP fluorophore was excited at 920 nm, using a resonant 667 scanning X-galvanometer (8 kHz) paired with a 6 mm standard scanning Y-galvanometer. The 668 scanning system was mounted on movable objective Ultima microscope, equipped with an 669 orbital nosepiece coupled to a 16X, 0.8NA, 3 mm water immersion objective (Nikon) and a 670 piezo drive for angled imaging and ultrafast volumetric scanning. Imaging was performed at a 671 scan speed of 29 fps, using 512x512 frame size (1.085 mm/pixel resolution). Fluorescence signal 672 was detected using high-sensitivity GaAsP photomultiplier tubes (model 7422PA-40 PMTs, 673 Hamamatsu). Recording sessions were ~10 minutes long, yielding datasets of ~10 GB.

675
Head-fixed spatial navigation 676 Mice were trained to run head-fixed on a 2-meter linear treadmill belt similar to that described 677 previously 18,50 . Sugar-water (5% sucrose) rewards were available from a lick spout in front of 678 the mouse at random positions uniformly distributed on the track. The position of the mouse 679 was measured using an optical rotary encoder (S5-720, US Digital). Lap onset was detected by 680 reading RFID tags with an RFID reader mounted below the animal (ID-20LA, SparkFun 681 Electronics). Behavioral programs were controlled with an Arduino Mega 2560 microcontroller. 682 All behavioral data was acquired at a sampling rate of 10 kHz synchronized to the two-photon 683 imaging. Mice ran on the same treadmill belt for 3 days to become acclimated to the textures 684 and cues. The following days (Days 4-7) were designated as "familiar days" for the purposes of 685 stability analysis. We always compared activity across consecutive familiar days for each animal. 686 Taking data recording and template matching artifacts into account, Days 5 & 6 were used for 687 all animals, except for two for which we compared Days 6 & 7 for one mouse and Days 4 & 5 for 688 another.

690
Statistics 691 Unless otherwise stated, all data are reported as the median and 95% confidence interval of the 692 median, in the form M, [L, U], where M is the median, L is the lower bound of the 95% 693 confidence interval, and U is the upper bound of the 95% confidence interval. The confidence 694 interval was found through a resampling procedure previously described 6 . 695 All statistical comparisons were done using a Wilcoxon sign-rank (for paired data) or rank-sum 696 (for unpaired data) unless otherwise specified (See section Subsampled Decoder). 697 698 All analyses were performed using custom-written codes in MATLAB R2018b (MathWorks). 699 700 d-NMF 701 702 Pre-processsing 703 Image stacks were motion corrected for XY-motion using NoRMCorre 86 , and temporally 704 downsampled by a factor of 20, resulting in image stacks with a time resolution of ~1.5 Hz. Each 705 downsampled frame was computed as the mean of the surrounding frames, such that frame N 706 of the downsampled image stack was computed from frames (N-1)*20+1 to N*20 of the 707 original stack. These values were rounded and stored as unsigned 16-bit integers. 708 709 Division into Patches 710 Image stacks were divided up into patches of size 64x64, with an overlap of 8 pixels, yielding 81 711 patches.

713
ROI Core Detection 714 Within each patch, the image sequence was first de-trended and transformed into ΔF/F values 715 by subtracting and dividing by the minimum value in a 45 frame (30 second) moving window. 716 The threshold for a pixel to be "active" was set to be 3 times the median ΔF/F for each pixel. 717 The thresholded binary image stack was then passed through a median filter of size 3 in the 718 time dimension. Connected components of active pixels were then detected in the resulting 3-719 dimensional (X-Y-Time) image stack. Components with a minimum membership size of 30 pixels 720 were projected into 2 (X-Y) dimensions. Those projected components with a minimum size of 721 15 pixels were kept, and the rest discarded. See Algorithm 1.

723
Components were then merged based on the Jaccard index (size of the intersection divided by 724 the size of the union) between components, with a merge threshold of 0.5. Note that the same 725 ROI active in distinct time spans will be detected multiple times in the previous step. Refinement and Activity Extraction 733 The following steps are standard fitting procedure for constrained non-negative matrix 734 factorization 28,32 (See Algorithm 2). The following steps are performed on the original image 735 sequence patch, without de-trending or ΔF/F transformation. An image patch with t time steps 736 and a d x d pixel size is stored as a matrix Y ∈ × × . This is reshaped into the matrix Y * 737 ∈ 2 × . The objective function is the following: 738 739 argmin , ‖ − ‖ 2 + ‖ ‖ 2 + ‖ ( )‖ 2 ; ≥ 0, ≥ 0 where ‖•‖ denotes the Frobenius norm. For k components, A ∈ 2 × represents the spatial 740 footprints (ROIs), and C ∈ × represents the time-varying activity for each ROI. η and β are 741 regularization terms to enforce sparsity of the temporal components (C) and spatial 742 components (A), respectively. For results presented in this manuscript, η = 0.01, and β = 0.5.

744
The algorithm then iterates over solving for C and A while enforcing non-negativity. The optimal 745 solutions for C and A are as follows: 746 747 C * = max((A T A + ηI)(Y T A) † , 0) 748 A * = max(YC T (CC T + β1) † , 0), 749 750 Where I represents the identify matrix, and 1 represents a matrix with each element equal to 1. 751 At each iteration n, C and A are then updated from their current values with a learning rate of 752 0.  Fig. 3c-d).

791
Signal to noise ratio 792 We also tested the efficacy of the signal-to-noise ratio as a classifier threshold (Extended Data 793 Fig. 3). For a detrended calcium trace, the signal-to-noise ratio was defined as SNR = P/M 794 Where P is the 99.9 th percentile of the trace, and M is the median absolute deviation of the 795 trace (median(|xmedian(x)|)). We defined SNR this way to be robust to noise and to be able 796 to define the noise levels without contamination from large transients.

798
CaImAn Parameters 799 The python implementation of CaImAn 28 was run on the motion corrected, downsampled 800 image stacks using default parameters using the "sparse_nmf" initialization procedure. Patch 801 sizes and overlaps were chosen to be 64 x 64 and 8 pixels to match those used by d-NMF. K = 15 802 components per patch were chosen to be estimated. Only those ROIs labeled as "valid" were 803 used for analyses (Fig. 3).
suite2p Parameters 806 Similarly, suite2p was run on the motion corrected, downsampled image stacks using pre-807 configured parameters for dendrites and axons. All labeled ROIs were used for analysis (Fig. 3). 808 809 F1 score for model comparison 810 The output of d-NMF, CNMF, and suite2p were compared to manually labeled data in a pixel-811 wise fashion. The metric used was the F1 score, calculated as 812 813 This measure is a combination of recall (True Positive Rate, TP/(TP+FN)) and precision (Positive 815 Predictive Value, TP/(TP+FP)), as human-labeled datasets may have incomplete labeling, which 816 would yield poor estimates of false positives.

818
Fitness Trace 819 To generate ground truth data to evaluate the Fitness Trace (Fig. 4), 90 ROIs generated from d-820 NMF (30 apical dendrites, 30 soma, and 30 basal dendrites, when possible) were chosen at 821 random from 9 FOVs. Their activity traces were re-estimated by taking a weighted average 822 across the ROI, which was then de-trended using a rolling minimum of 45 frames, then z-823 scored. Putative transients were identified as contiguous time bins where the z-scored trace 824 exceeded 2 standard deviations. These event detection criteria was meant to be permissive, as 825 it was not performed on the de-mixed data, so false transients would likely be present. This was 826 by design, to contain examples of valid and invalid transients. Transients were then manually 827 reviewed and classified as valid or invalid using a previously published graphical user 828 interface 37 .

830
The image stack was de-trended and transformed into ΔF/F values, then z-scored across time.

831
For a given ROI, a square bounding box was defined around all non-zero values with a padding 832 of 3 pixels in each dimension. The Fitness Trace for a given ROI at time t was defined as the 833 correlation coefficient between the z-scored image stack at time t and the ROI, calculated only 834 using the pixels in that bounding box.

836
Transient Detection Methods 837 To describe the different detection methods (Fig. 4, Extended Data Fig. 4), it is useful to define 838 two types of activity traces. The first is the activity trace arrived at through d-NMF, that is, a 839 signal designed to be de-mixed from overlapping sources; we will denote this as De-mixed ΔF/F. 840 The second is a "simple" signal obtained through taking a weighted average across all of the 841 pixels of an ROI. We will denote this as Simple ΔF/F. ROIs were manually classified as apical dendrites, soma, or basal dendrites based on position in 854 the field of view. ROIs encompassing both a soma and one or more dendritic branches were 855 classified as soma. In all other instances we treat apical dendrites, soma, and basal dendrites as 856 independent populations, without assigning parent soma to individual dendrites. Thus, when 857 we refer to apical dendrites or basal dendrites, we are referring to dendrites which are some 858 distance from their parent cell body, which is often not visible in the same field of view.

860
We used two measures to quantify the degree of tuning of each rate map. Maps were divided 861 into 40 bins as described above. Spatial information content 88 was computed as 862 Where P i is the normalized occupancy in the ith bin, such that the sum across all P i = 1, and R i is 863 the normalized value of the rate map in the ith bin, such that the sum across all R i = 1. 864 Sparsity was computed as 865 Where R i is the value of the rate map in the ith bin.

868
Each ROI served as its own control to determine statistical significance of spatial tuning. For 869 each ROI, the transient times were shifted by n time bins, for n = -250 to 250. The information 870 content at a shift of 0 was then normalized to the 99 th percentile of information content at non-871 0 shift. A normalized information content value greater than 1 indicated significant tuning. 872 Additionally, a minimum of 4 transients must have been detected for an ROI to be called 873 significantly tuned.

875
Tuning Parametrization 876 Tuning curves were parametrized as a mixture of Von Mises functions, as described 877 previously 89 . Rate maps were circularly smoothed with a sigma of 15 cm (3 bins), and the curve 878 was fit by increasing the number of Von Mises functions until the residuals were below 25% of 879 the maximum value of the original rate map. Place field width was defined as the width of a 880 fitted component at 50% of the component's amplitude (i.e. full width at half max). For ROIs 881 with multiple place fields, their width was defined as the mean width of all place fields. 882 883 884 Quantifying Overlap 885 We defined two ROIs as overlapping if they shared at least 8 pixels and had a correlation 886 coefficient between their activity traces less than 0.75. This avoids inflating the overlap amount 887 from ROIs that are part of the same dendritic tree.

889
Stability 890 Within-day stability (Fig. 6) was assessed by splitting the data from Familiar Day 1 in half. Rate 891 maps were re-estimated for each ROI using only the data in each session half. Tuning curve (TC) 892 correlation was defined as the correlation coefficient between the tuning curves for each half.

893
This measure was computed for each ROI and then averaged to obtain a single value for the 894 session. To be included in the calculation for TC correlation, an ROI must have had a minimum 895 of 4 calcium transients and be significantly tuned across the entire session. Population vector 896 (PV) correlation was computed by constructing the rate estimates of all ROIs for a given 897 position bin in the first half and taking the correlation coefficient between that vector and the 898 corresponding vector in the second half. This yields a value for each position bin, which was 899 then averaged to obtain a single value for a session. To be included in the calculation for PV 900 correlation, an ROI must have had a minimum of 4 calcium transients, with no criteria for 901 significant tuning.

903
To quantify across-day stability (Fig. 6, Extended Data Figs. 9-10), image stacks from Familiar 904 Days 1 and 2 were stitched into a single image stack and motion corrected together. d-NMF was 905 then run on the combined image stack to obtain ROIs that spanned both sessions. Stability 906 analysis was then performed as for within-day estimates. To be included in across-day stability 907 analyses, an ROI must have had at least 4 detected calcium transients in both days and 908 significantly tuned in at least one day. The activity criterion ensured that correlation estimates 909 were not skewed by ROIs that were not visible on one of the recording days.

911
Structural Stability 912 To verify that differences in tuning stability were not due to the inability to track ROIs across 913 days, we additionally quantified the structural stability of ROIs (Extended Data Fig. 10). For each 914 ROI we took the average of frames at which significant calcium transients were detected for 915 Familiar Days 1 and 2. The structural stability was defined as the correlation coefficient 916 between these averages, using only data within the modified boundary of the ROI. The 917 modified boundary of an ROI was defined by first thresholding the ROI above 20% of its 918 maximum pixel intensity, then dilating the resulting mask by one pixel. 919 920 Population Vector Decoding 921 To test the impact of increased across-day stability in apical dendrites, we performed 922 population vector decoding in a manner similar to that described previously 50 (Fig. 7). A 923 separate decoder was constructed for each FOV for each of apical dendrites, soma, and basal 924 dendrites. Template tuning curves for each ROI were constructed as above. For within-day 925 decoding, data from the first half of Familiar Day 1 was used to define the template. For across-926 day decoding, all data from Familiar Day 1 was used to define the template. As above for 927 defining Population Vector correlations, only ROIs that were active in both days and 928 significantly tuned in at least one day were included for decoding.

930
Time-varying rate vectors for each ROI were constructed using data from either the second half 931 of Familiar Day 1 (within-day decoding) or all data from Familiar Day 2 (across-day decoding), 932 using 667 ms bins smoothed with a Gaussian smoothing kernel with σ = 3.3 seconds. For each 933 time point in the decoded portion of the data, the decoded position was the position 934 corresponding to the highest correlation with the template matrix. Any time points with zero 935 activity across all relevant ROIs had undefined correlation with any position. The most recently 936 decoded position was copied in to those undefined frames. 937 938