Representational content determines the long-term dynamics of the CA1 spatial code


 Hippocampal subregion CA1 is thought to support episodic memory by reinstating a stable spatial code. Yet recent experiments have demonstrated that this code is largely unstable on a timescale of days, challenging its presumed function. While these dynamics may indeed reflect homogenous drift within the population, they may alternatively reflect distinct time-varying representational component(s) which coexists alongside other stable components. Here we adjudicate between these possibilities. To this end, we characterized the mouse CA1 spatial code over more than a month of daily experience in an extended geometric morph paradigm. We find that this code is governed by distinct representational components with different long-term dynamics, including stable components representing spatial geometry and prior experience. These components are mediated by separate neural ensembles with similar short-term spatial reliability and precision. Together, these results demonstrate that the long-term dynamics of the CA1 spatial code are defined by representational content, not homogenous drift.


Introduction
Hippocampal subregion CA1 represents a mixture of external and internal cuesincluding the shape of the environment 1 , visual landmarks 2,3 , objects 4 , task-relevant information 5 , and past experience [6][7][8][9] through changes to the spatial tuning properties of its principal cells 10 . This code is hypothesized to support spatial and episodic memory by later reinstating its content. However, longitudinal experiments in mice have challenged this view 11,12 , reporting that the vast majority of CA1 place cells change their spatial tuning properties more rapidly than expected (on the order of days), while only a subset of these cells (15-25%) provide a stable spatial representation on this timescale 12 .
These long-term representational dynamics may reflect two possible phenomena, each carrying its own theoretical implications. On one hand, these dynamics may reflect homogenous population-wide representational drift, as could result from gradual but persistent plasticity mechanisms throughout the hippocampal circuit 13,14 . In this case, drift would be independent of the content of the representation, with complete representational turnover on sufficiently long timescales. Such drift would thus pose a strong challenge to the persistent reinstatement hypothesis, and would at the very least motivate a revision to include compensatory downstream mechanisms which might keep pace with, filter out, or otherwise overcome this drift 14 . Alternatively, these dynamics may instead reflect time-varying representational component(s) alongside other stable components, as might arise from the convergence of both stable and time-varying inputs onto the population. In this case, we would expect time-varying dynamics to be restricted to distinct component(s) of the representation at the level of the population, and possibly at the level of individual neurons. Such component-specific dynamics would pose a minimal challenge to the persistent reinstatement hypothesis, as downstream readers tuned to a stable representational component would be blind to any variability in other orthogonal components.
To adjudicate between these possibilities, one must assay multiple representational components over long timescales and determine whether drift depends on, or is independent of, representational content. To this end, we characterized the mouse CA1 spatial code over more than a month of daily experience in an extended version of the classic geometric morph paradigm. We find that the CA1 spatial code is defined by distinct representational components mediated by separate neural ensembles with different long-term dynamics, including stable components reflecting the shape of space and prior experiences.

Results
We recorded daily from CA1 with calcium imaging (Fig. 1A) in mice freely exploring open environments for 32 days in an extended version of a geometric morph paradigm (5 mice, 160 sessions total; Fig. 1B; Fig. S1). In this paradigm, mice were first familiarized with two geometrically-distinct environments, and later tested in deformed (morphed) versions of these environments spanning the shapespace between the two familiar environments. Previous acute versions of this paradigm have demonstrated that a subset of CA1 cells exhibit experiencedependent attractor-like properties in their spatial tuning, with an abrupt transition in their preferred firing locations near the midpoint of the morph sequence, while other cells fired at geometrically similar locations across the morph sequence 8,9 . By employing more geometricallydissimilar familiar environments and extending this paradigm across a 32 day period, we can thus assay three potential components of the CA1 representationa representation of the spatial geometry (geometry), a representation of experience-dependent attractor dynamics (attractor), and a representation which varies with time (time).
We began by first characterizing these data within the morph sequence akin to previous acute reports. Following motion correction 15 , cells were segmented and calcium traces were extracted via constrained nonnegative matrix factorization 16,17 (Fig. S1). The likelihood of spiking events (henceforth firing rate) was inferred from deconvolution of the filtered calcium traces via a second-order autoregressive model 18 . Cell identity was tracked across recording sessions on the basis of extracted spatial footprints (n = 2745) 19 . We found that individual cells exhibited a variety of dynamics, including abrupt transitions in their preferred locations near the midpoint of the morph sequence, persistent firing at geometrically-similar locations across the morph sequence, and the appearance of punctate spatial tuning in a subset of environments (Fig. 1C). Across the population, the morph sequence produced sigmoidal transitions between the (partially-correlated) familiar environment rate maps ( Fig. 1D; Fig. S2). Fitting these transitions with sigmoidal functions allowed us to quantify dynamic properties across sequences within each mouse, revealing qualitatively distinct patterns within each mouse and a trend toward increased decorrelation of maps across sequences (Fig. S2).  Traditional within-morph sequence analyses only characterize the representational similarity between subsets of sessions. To take advantage of the full representational structure, we explored population-level representational similarity across all 32 days. To this end, we first computed the median rate map similarity across all tracked cells for each pairwise comparison of sessions ( Fig. 2A). We then reduced this representation to two dimensions via nonmetric multidimensional scaling, an unsupervised technique for visualizing high-dimensional structure in a digestible low-dimensional space that preserves the representational similarity between sessions as well as possible. The resultant embedding revealed two dimensions that strongly determined representational similarity in this paradigm: a spatial context component and a timevarying component (Fig. 2B). These components were pronounced in every mouse. Quantification of the embedded representations revealed that the contextual and time-varying components defined nearly orthogonal dimensions in this 2D subspace, with an absolute angular difference of 89.0° ± 5.8° (mean ± circular standard deviation) between context and time-varying dimensions (see Methods; all angles: [83.0°, 84.6°, 85.5°, 95.2°, 96.7°]; Rayleigh's test versus uniformity on the 0° to 180° range: p = 2.47e-3, z = 4.80). This reliable structure provides evidence that spatial context and time-varying components of the CA1 spatial code are representationally distinct at the level of the population. The observation of distinct population-level representational components could be the product of a distributed code across neurons that exhibit mixed-selectivity. Conversely, it could also be a product of separate neural ensembles coding for individual representational components. To adjudicate between these possibilities, for each cell with at least 28 comparisons we created a representational similarity matrix (RSM) of the pairwise rate map correlations between all tracked sessions (Fig. 3A) 20 . (Twenty-eight was chosen as the number of pairwise comparisons within a single bookend morph sequence.) Next, we constructed predicted RSMs on the basis of similarity in time, spatial geometry, and contextual group as defined by the attractor dynamics (see Methods). Using these matrices, we computed the amount of variance in each (Fisher-transformed) cell RSM explained by a model containing the three predicted RSMs as well as all possible interactions (total r 2 ). We then determined the amount of total explained variance attributable to each factor by quantifying the reduction in explained variance when selectively excluding that factor (and its interactions).
The results of this analysis revealed remarkable selectivity in representational content and long-term dynamics across simultaneously-recorded cells (Fig. S3). Across all cells with a full model r 2 of at least 0.1915 (90 th percentile of the shuffled distribution, n = 1316 of 2745, 47.9% of cells with sufficient comparisons, binomial test: p = ~0.0), explained variance was typically attributable to a single factor (Fig. 3B). Moreover, reorganization of RSMs for time-, attractor-, and geometry-loading cells revealed clear distinctions in their temporal and contextual dynamics ( Fig. 3C; Fig S4). To quantify these observations, we computed null distributions by randomly shuffling RSM values within each comparison across cells. The resulting null distributions thus preserve the population-level statistics while destroying any within-cell crosssession effects. Relative to these distributions, the first-ranked factor was attributed slightly but significantly more variance, while the second-and third-ranked factors were attributed significantly less variance than would be expected given the variability within the population alone ( Fig. 3D,E)a pattern which held across a range of inclusion criteria and grew more pronounced with more stringent criteria (Fig. S5). These results provide evidence that simultaneously-recorded neural ensembles are selectively tuned to distinct representational components, governed by differing long-term dynamics.
Are differently-loading cells distinguished by more than their long-term dynamics and representational content? To address this possibility, we grouped cells according to primary factor loading and compared a number of properties as follows. We first examined whether groups differed by the short-term reliability and precision of their spatial codes. For cells in all groups, full model r 2 was correlated with median within-session split-half correlation (SHC; Fig.  3F), indicating that a better model fit was associated with a more reliable spatial code in the short-term regardless of group. Median within-session SHCs of cells meeting full model r 2 criteria exceeded those of cells not meeting criteria and were similarly high in all groups (Fig.  3G), indicating comparable spatial reliability on a timescale of minutes. Moreover, despite the differences in long-term dynamics between groups (Fig. 3H), within-session SHCs remained high across all sessions regardless of group (Fig. 3I). To quantify this observation, for each cell we computed the least squares linear fit of across-session map similarity and within-session SHC as a function of temporal lag between sessions. Time-loading cells exhibited a modest but reliable decrease in across-session rate map correlations for each day separating compared sessions, which was significantly greater than that observed in attractor-loading or geometryloading cells (Fig. 3J). In contrast, time-loading cells exhibited significantly less of a decrease in within-session SHC as a function of lag, and differences between groups were relatively minor (Fig. 3J). Echoing these reliability results, the precision of the spatial code as measured by median spatial information content (SIC) was higher for cells meeting criteria and did not substantially differ between groups (Fig. 3K). Together, these results indicate that short-term spatial reliability and precision is predictive of selectively coding some representational content, but is not predictive of the specific content itself. Table S1. (A) Each cell RSM was fit with a general linear model consisting of time, attractor, and geometry components including interactions. Then the relative reduction in r 2 was computed when excluding each factor. Cell RSMs were Fisher-transformed before fitting. Note that model attractor and geometry RSMs tend to be correlated, limiting the amount of variance that can be attributed to one or the other. We investigated whether other coding properties may distinguish these cell groups. Some theories posit that the most excitable cells are more likely to be recruited into a representation 13 , and firing rate is known to correlate with spatial reliability and precision on short timescales 21 . In our data, cells meeting the full model r 2 criterion tended to have higher mean and peak firing rates than cells not meeting this criterion (Fig. 4A,B). Interestingly, attractor and geometry-loading cells exhibited higher firing rates than time-loading cells (Fig.  4A,B), possibly reflecting differences in the intrinsic properties and/or inputs specific to these ensembles. Biases in place field locations also differed by group (Fig. 4C). Time-loading cells tended to tile the navigable space more evenly with a bias toward corners, attractor cells tended to have fields near the center of the environment, and geometry-loading cells typically had fields near the deformed corners. Time-loading cells were tracked across more sessions than attractor-or geometry-loading cells (Fig. 4D), consistent with a lack of context-selectivity in these cells and indicating that variance attributed to time was not the result of cell dropout. Interestingly, all groups were tracked across fewer sessions than cells not meeting criteria, also consistent with the time-and context-dependence of grouped cells and suggesting that unclassified cells did not fail to reach criteria due to dropout. The proportion of each group within the population differentially predicted morph sequence properties: the proportions of geometryand attractor-loading cells were positively related to the decorrelation between familiar environment representations, while the proportion of time-loading cells was negatively related to this decorrelation (Fig. S2). Signal-to-noise differences between groups were minor (Fig. S5). No clear anatomical organization according to primary factor was observed (Fig. S6). Thus, ensembles can be partially distinguished by properties such as firing rate, field location, and the number of tracked sessions.  Table S1.

(A) Median mean firing rate (MFR) as a function of group for cells meeting criteria, and for cells not meeting criteria (unclassified). (B) Median peak firing rate (PFR) as a function of group and for cells not meeting criteria. (C)
Field density within the environment separately for each group. Field density was computed by first averaging the maps from all tracked sessions within each cell. Next, density at a given pixel was defined as the proportion of cells with activity exceeding 75% of their average map PFR at that pixel. Shuffled density distributions were computed by randomly permuting group assignment 1000 times. (D) Number of sessions in which a cell was tracked as a function of group and for cells not meeting criteria. **p < 0.01, ***p < 0.001

Discussion
Here we show that the CA1 spatial code is defined by distinct representational components with different long-term dynamics, including stable components representing spatial geometry and prior experience. Stable and time-varying components were distinct at the level of the population as well as individual cells. These components were not distinguishable by their shortterm spatial stability or precision, though they were partially differentiated by other characteristics such as firing rate and field location.
These findings have several important implications. First, by assaying representational content across extended experience, these results go beyond prior work 12,22 to establish that differing long-term dynamics governing the CA1 spatial code are mapped to distinct representational components, including time, spatial geometry, and experience-dependent attractor dynamics. Components representing the spatial context, such as those coding spatial geometry and prior experience in similar environments, are persistently reinstated across a timescale of weeks. This mix of persistent and time-varying representational dynamics is consistent with many theories of episodic memory and hippocampal function, such as indexing theory 23 and temporal context models 24 . In such theories, a time-varying component allows for the encoding and memory of the temporal relationship between different events, while the persistent context-dependent component allows for linking events occurring within the same spatial context across time.
These results also demonstrate that long-term changes in the CA1 spatial code do not reflect homogenous representational drift, which could result from uniform plasticity at the level of the population 13 . Rather, these changes likely reflect either the selective influence of plasticity on a distinct subpopulation within CA1 and/or the influence of time-varying input(s) mapped onto this subpopulation, such as those found in the lateral entorhinal cortex 25,26 and hippocampal subregion CA2 27 . Such time-varying inputs may include surrounding soundscapes, scentscapes, motivations, and intervening experience that cannot be completely controlled in this paradigm, as well as a representation of time per se.
Our finding that these distinct components coexist within the population code adds to a growing literature that challenges the traditional view of a representationally homogenous CA1 spatial code 28-30 , and in turn also challenges models which assume a single state space 31 or basis set 32 as the source of the CA1 spatial code. Notably, the short-term spatial reliability and precision of a cell was predictive of that cell coding some representational content but not the particular content itselfsuggesting that such properties alone are generally insufficient for inferring the representational content of a cell. Finally, these results provide strong evidence of the heterogenous and malleable nature of hippocampal remapping 33,34 , one that is linked to a dynamic recruitment and subsequent reinstatement of distinct neural ensembles 13

Subjects
Naive male mice (C57Bl/6, Charles River) were housed in pairs on a 12-hour light/dark cycle at 22°C and 40% humidity with food and water ad libitum. All experiments were carried out during the light portion of the light/dark cycle, and in accordance with McGill University and Douglas Hospital Research Centre Animal Use and Care Committee (protocol #20157725) and with Canadian Institutes of Health Research guidelines.

Surgeries
During all surgeries, mice were anesthetized via inhalation of a combination of oxygen and 5% Isoflurane before being transferred to the stereotaxic frame (David Kopf Instruments), where anesthesia was maintained via inhalation of oxygen and 0.5-2.5% Isoflurane for the duration of the surgery. Body temperature was maintained with a heating pad and eyes were hydrated with gel (Optixcare). Carprofen (10 ml kg -1 ) and saline (0.5 ml) were administered subcutaneously at the beginning of each surgery. Preparation for recordings involved three surgeries per mouse.
One to three weeks post-injection, either a 1.8mm (AKCA127, AKCA130, AKCA135) or 0.5mm (AKCA170, CAMCA130) diameter gradient refractive index (GRIN) lens (Go!Foton) was implanted above dorsal CA1 (Referenced to bregma: ML = 2.0 mm, AP = -2.1 mm; Referenced to brain surface: DV = -1.35 mm). Implantation of the 1.8mm diameter GRIN lens required aspiration of intervening cortical tissue, while no aspiration was required for implantation of the 0.5mm diameter GRIN lens. Results observed using 1.8-or 0.5-mm diameter GRIN lenses were similar. In addition to the GRIN lens, two stainless steel screws were threaded into the skull above the contralateral hippocampus and prefrontal cortex to stabilize the implant. Dental cement (C&B Metabond) was applied to secure the GRIN lens and anchor screws to the skull. A silicone adhesive (Kwik-Sil, World Precision Instruments) was applied to protect the top surface of the GRIN lens until the next surgery.
One to three weeks after lens implantation, an aluminum baseplate was affixed via dental cement (C&B Metabond) to the skull of the mouse, which would later secure the miniaturized fluorescent endoscope (miniscope) in place during recording. The miniscope/baseplate was mounted to a stereotaxic arm for lowering above the implanted GRIN lens until the field of view contained visible cell segments and dental cement was applied to affix the baseplate to the skull. A polyoxymethylene cap with a metal nut weighing ~3 g was affixed to the baseplate when the mice were not being recorded, to protect the baseplate and lens, as well as to simulate the weight of the miniscope.
After surgery, animals were continuously monitored until they recovered. For the initial three days after surgery mice were provided with a soft diet supplemented with Carprofen for pain management (MediGel CPF, ~5mg kg -1 each day). Familiarization to both environments (while recording in the room A-associated environment to monitor imaging quality and habituate the mouse to recording) began 3 to 7 days following the baseplate surgery.

Data acquisition
In vivo calcium videos were recorded with a UCLA miniscope (v3; miniscope.org) containing a monochrome CMOS imaging sensor (MT9V032C12STM, ON Semiconductor) connected to a custom data acquisition (DAQ) box (miniscope.org) with a lightweight, flexible coaxial cable. The DAQ was connected to a PC with a USB 3.0 SuperSpeed cable and controlled with Miniscope custom acquisition software (miniscope.org). The outgoing excitation LED was set to between 2-8% (~0.05-0.2 mW), depending on the mouse to maximize signal quality with the minimum possible excitation light to mitigate the risk of photobleaching. Gain was adjusted to match the dynamic range of the recorded video to the fluctuations of the calcium signal for each recording to avoid saturation. Behavioral video data were recorded by a webcam mounted above the environment. Behavioral video recording parameters were adjusted such that only the red LED on the CMOS of the miniscope was visible. The DAQ simultaneously acquired behavioral and cellular imaging streams at 30 Hz as uncompressed avi files and all recorded frames were timestamped for post-hoc alignment.
All recording environments were constructed of a grey Lego base and black Lego bricks (Lego, Inc) according to the dimensions specified in the main text and supplemental figures. All external walls had a height of 22 cm. During recording, the environment was dimly lit by a nearby computer screen, which could serve as directional cue. During familiarization and all unrecorded experience the environments were well-lit by overhead room lighting. All sessions were 20 min, and only one session was recorded per day to avoid photobleaching. The mouse was always placed in the same corner at the start of the session and was allowed to explore the environment for 15 to 30 s prior to the start of data acquisition. Following each recording the environment was cleaned with disinfectant (Prevail).

Data preprocessing
Calcium imaging data were preprocessed prior to analyses via a pipeline of open source MATLAB (MathWorks; version R2015a) functions to correct for motion artifacts 15 , segment cells and extract transients 16,17 , and infer the likelihood of spiking events via deconvolution of the filtered transient trace through a second-order autoregressive model 18 . The motion-corrected calcium imaging data were manually inspected to ensure that motion correction was effective and did not introduce additional artifacts. Following this preprocessing pipeline, the spatial footprints of all cells were manually verified to remove lens artifacts. Position data were inferred from the onboard miniscope red LED offline following recording using a custom written MATLAB (MathWorks) script and were manually corrected if needed. Cells were tracked across sessions on the basis of their spatial footprints and/or centroids 19 .

Data analysis
All analyses were conducted using the vector of likelihood of spiking events inferred via deconvolution, treating this vector as if it were the firing rate of the cell (henceforth firing rate) 18 . Similar results were observed when the likelihood of spiking was inferred as a binarized vector indicating periods when calcium transients were rising, instead of deconvolution. Only periods of movement for which the mouse speed exceeded 2 cm s -1 were included in analyses, though no qualitative differences were observed if this criterion was relaxed.
Rate maps were constructed by first binning the position data into pixels corresponding to a 2.5 cm x 2.5 cm grid of locations. Then the mean firing rate was computed for each pixel and then smoothed with a 4 cm standard deviation isometric Gaussian kernel. For all comparisons between rate maps, similarity was measured as the Pearson's correlation between corresponding pixels. Individual cell representational similarity matrices (RSMs) were thus computed as the Pearson's correlation between pairwise rate map comparisons for tracked sessions. Population-level RSMs were computed by taking the median across all cells for each pairwise comparison.
Within-sequence analyses summarized morph sequence dynamics with a transition plot, which captured the similarity of all 6 morph environment maps to the 2 familiar environment maps. To this end, rate map correlations between each environment and familiar rate maps at the beginning and end of the morph sequence were computed. Only comparisons between cells whose within-session split-half rate map correlation (SHC) exceed the 99 th percentile of a shuffled control for at least one of the compared sessions were included. The shuffled control was computed for each cell by randomly circularly shifting its firing rate vector relative to the position data by at least 30 sec and recomputing the SHC 1000 times to create the null distribution. The median across all pairwise cell comparisons was taken as the measure of central tendency for transition plots. To characterize these plots, the Fisher-transformed median values comparing the morph environments to each of the two familiar environments were both fit to minimize mean squared error with a 4-parameter sigmoid of the form where is the position of the current environment in the shapespace (arbitrarily chosen to range from 1 to 6), and parameters determine the shape of the sigmoid. The intersection of these two sigmoids was taken as the measure of the transition point between familiar maps. The maximum absolute difference between these two sigmoids within the sampled shapespace was taken as the measure of maximum map decorrelation. Sigmoidality was computed as the mean absolute difference between each sigmoid and a linear fit from endpoint to endpoint, with endpoints normalized to range from zero to one.
Two-dimensional embeddings of the population representation were computed and quantified as follows. First, the population RSM was computed by taking the median across cells for each comparison. Next, this RSM was computed to a distance matrix as 1-(RSM) with the diagonal set to a distance of 0. Finally, this distance matrix was reduced to a twodimensional embedding via Kruskal's nonmetric multidimensional scaling, implemented by the MATLAB function mdscale with the default parameterizations. To estimate the context dimension within this embedding, the average difference separating neighboring-in-time familiar environment recordings was computed. To estimate the time dimension within this embedding, the average difference between the first and last recordings of each familiar environment was computed.
Characterization of individual cell RSM variance was computed as follows. Only cells tracked across at least 28 pairwise comparisons (8 sessions) were included. Model RSMs were computed on the basis of temporal, geometric, and attractor group similarity, as schematized in Figure 3a. The time RSM was computed as the difference in session number between pairwise comparisons. The geometry RSM was the difference in environment shape between pairwise comparisons, assuming that consecutive environments mapped linear steps through the shapespace separating the two familiar environments. The attractor RSM was computed such that if a pairwise comparison was between two sessions on the same side of the withinsequence transition point then those environments were assigned a similarity of 1; otherwise, that comparison was assigned a similarity of 0. Models consisting of all factors and their interactions (full model) as well as models with each factor and their interactions removed were then fit to the Fisher-transform individual cell RSMs as described in the main text. Notably, geometry and attractor factors were often correlated (though the extent of their correlation differed depending on transition point dynamics and the particular subset of tracked sessions for each cell), often limiting the amount of explained variance that could be attributed to one or the other. The time RSMs tended to be uncorrelated with either the geometric or attractor RSMs.
Signal-to-noise (SNR) for each dF/F trace was computed as a measure of data quality prior to spike estimation via deconvolution as described in previously 35 . Briefly, because calcium transients around the baseline can only be positive, for each trace the variance of a normal noise distribution was estimated by scaling the standard deviation of trace values below baseline via a half-normal distribution, such that: where is the detrended filtered calcium trace and ( < 0) are all filtered trace values below baseline. Next, we z-score the trace such that ( ) = NOISE ⁄ . These z-scored values are next transformed into probabilities ( ) = (− ( )), where denotes the cumulative distribution function of the standard normal distribution. Then we compute the lowest probability event over a time window of 0.4 s ( = 12 frames), the duration around which a GCamp6F transient is typically at its maximum, such that: where min is the minimum value across all timepoints from = 1 to = the length of − + 1.
Finally, we define the average peak SNR as: where −1 is the quantile function of the standard normal distribution (logit function).
Spatial information content (SIC) was computed from the whole-session rate maps of each cell as described previously 36 via the equation: where is the rate map pixel index, is the probability of sampling pixel , is the mean firing rate at pixel , and ̅ is the mean firing rate across all pixels.

Histological validation of expression and recording targets
After experiments, animals were perfused to verify GRIN lens placement. Mice were deeply anesthetized and intracardially perfused with 4% paraformaldehyde in PBS. Brains were dissected and post-fixed with the same fixative. Coronal sections (50 μm) of the entire hippocampus were cut using a vibratome and sections were mounted directly on glass slides. Sections were split and half of all sections were stained for DAPI and mounted with Fluoromount-G (Southern Biotechnology) to localize GRIN lens placement and to evaluate virus expression. Due to the large imageable surface but restricted miniscope field of view (~0.5 mm x ~0.8 mm), we were unable to determine more specific localization of populations within the hippocampus for mice recorded with 1.8 mm lenses.

Statistics and reproducibility
All statistical tests are noted where the corresponding results are reported throughout the main text and supplement. All tests were uncorrected 2-tailed tests unless otherwise noted. Z-values for nonparametric Wilcoxon tests were not estimated or reported for comparisons with fewer than 15 datapoints. In all cases, the median served as the measure of central tendency unless otherwise noted. Box plots portray the minimum and maximum (whiskers), upper and lower quartiles (boxes), and median (cinch).

Code availability
All custom code written for reported analyses are publicly available at [insert Github link] or via request to the corresponding authors.

Data availability
The complete dataset for all experiments are publicly available at [insert DRYAD link] or via request to the corresponding authors.