Mapping the computational similarity of individual neurons within large-scale ensemble recordings using the SIMNETS analysis framework

Abstract

The expansion of large-scale neural recording capabilities has provided new opportunities to examine multi-scale cortical network activity at single neuron resolution. At the same time, the growing scale and complexity of these datasets introduce new conceptual and technical challenges beyond what can be addressed using traditional analysis techniques. Here, we present the Similarity Networks (SIMNETS) analysis framework: an efficient and scalable pipeline designed to embed simultaneously recorded neurons into low dimensional maps according to the intrinsic relationship between their spike trains, making it possible to identify and visualize groups of neurons performing similar computations. The critical innovation is the use of pairwise spike train similarity (SSIM) matrices to capture the intrinsic relationship between the spike trains emitted by a neuron at different points in time (i.e., different experimental conditions), reflecting how the neuron responds to time-varying internal and external drives and making it possible to easily compare the information processing properties across neuronal populations. We use three publicly available neural population test datasets from the visual, motor, and hippocampal CA1 brain regions to validate the SIMNETS framework and demonstrate how it can be used to identify putative subnetworks (i.e., clusters of neurons with similar computational properties). Our analysis pipeline includes a novel statistical test designed to evaluate the likelihood of detecting spurious neuron clusters to validate network structure results. The SIMNETS framework provides a way to rapidly examine the computational structure of neuronal networks at multiple scales based on the intrinsic structure of single unit spike trains.

1 Introduction

The neural processes underpinning complex sensory, cognitive, and behavioral phenomena engage complex activation patterns across brain-wide networks (Fox et al., 2005; Greicius et al., 2003; Raichle et al., 2001; Sporns, 2013; Thomas Yeo et al., 2011). Within large-scale networks, smaller clusters of computationally interrelated neurons have been proposed to embody fundamental input–output processing modules (‘subnetworks’) responsible for perceptual integration, dexterous motor control, and memory storage or retrieval (Briggman et al., 2006; Harris, 2005; Truccolo et al., 2010; Yuste, 2015; Zagha et al., 2015). While recent technological advances have made it possible to record from ever larger neuronal populations at single neuron resolution (Buzsáki, 2004; Chung et al., 2022; Donoghue, 2002; Maynard et al., 1997; Paulk et al., 2022; Steinmetz et al., 2021; Steinmetz et al., 2018), progress in understanding the computational organization of large networks has been hindered because of the lack of appropriate conceptual and technical frameworks that are mathematically principled; amenable to a range of statistical tools, flexible, scalable; and fast enough to work on large datasets (Brown et al., 2004; Gao and Ganguli, 2015; Paninski and Cunningham, 2018; Stevenson and Kording, 2011; Zador et al., 2022). Objectively identifying and characterizing subnetwork structure within and across brain areas would greatly simplify the process of tracking information flow within cortical circuits, modeling multi-scale neural dynamics, and understanding the general principles of neural information processing (Alivisatos et al., 2013; Kohn et al., 2020; Urai et al., 2022; Yuste, 2015).

How can computational subnetworks be identified? In cortical areas where the tuning properties of individual neurons are well-established, such as the primary visual (V1) and motor (M1) regions, the most straightforward approach to assessing computational inter-relationships would be to calculate tuning parameter similarities across neurons (Berman et al., 1987; Georgopoulos et al., 1982; Livingstone et al., 1996). However, simple parametric tuning models often fail to capture the temporal complexities of single neuron outputs and can perform poorly or even break down entirely under more ecologically or ethologically relevant experimental conditions (David et al., 2004; Hosman et al., 2021; Olshausen and Field, 2005; Pospisil and Bair, 2021). Moreover, classic single-neuron tuning function estimation techniques can miss or “average away” computationally-relevant features of individual spike trains (Churchland and Shenoy, 2007; Cunningham and Yu, 2014; Malik et al., 2011).

A variety of useful and important tools have been developed to study coordinated neuronal activity based on the underlying premise that similar spike patterns shared by a pair of neurons imply similar information processing properties (Abeles and Gat, 2001; Aertsen and Gerstein, 1985; Amarasingham et al., 2012; Gerstein et al., 1985; Gerstein and Michalski, 1981; Grün et al., 2002a; Kiani et al., 2015; Lopes-dos-Santos et al., 2013). These methods group neurons according to the similarity of spiking statistics (e.g., coincident spiking or correlated spike rate fluctuation). However, the mechanistic origin and computational significance of these correlations have proven more complex than initially theorized and are not fully understood (Brody, 1999; Cole et al., 2016; de la Rocha et al., 2007; Friston, 2011; Roudi et al., 2015; Smith and Kohn, 2008; Soudry et al., 2013). This approach has limited the utility and explanatory power of inter-neuronal correlations as a means for establishing the computational relatedness of neurons within a common feature space. These interpretational issues are further confounded by the significant statistical and computational challenges that come with implementing these methods at scale (Amarasingham et al., 2012; Brody, 1999; Cohen and Kohn, 2011; de la Rocha et al., 2007; Soudry et al., 2013).

Here, we present SIMNETS, an unsupervised relational analysis framework designed to generate low-dimensional neuron maps that quantify and support a multi-scale view of the computational similarity relations among individual neurons. For the purposes of our analysis, we define computation as the mapping of a set of inputs to a set of outputs according to a given set of rules (Agüera y Arcas et al., 2003; Kanerva, 2009). Under these conditions, the computational equivalence of any two neural systems should not be sought solely in the actual pattern of their spiking outputs, but in the relations of their outputs to one another within each system (Kanerva, 2009; Shepard and Chipman, 1970). Based on this principle, SIMNETS compares neurons using the intrinsic relational structure of their firing patterns, represented by a spike train similarity (SSIM) matrix that captures the relative changes in activity across a set of predetermined time windows. In this way, each SSIM matrix can be considered to represent the “output space” of a neuron, which serves as a computational “fingerprint” within the context of a given dataset (Vargas-Irwin et al., 2015a). By employing the single neuron SSIM matrix representation of individual neuron output spaces, we can efficiently quantify the computational relationships among all neuron pairs within a population of concurrently recorded single units using relatively straightforward techniques, such as vector correlation or regression. As we will demonstrate, this type of inter-neuronal 2nd-order spike train comparison—the Computational Similarity (CS) score—becomes an effective means for identifying neurons that share computationally similar output spaces, even in cases where diverse encoding strategies are employed. For instance, consider two neurons responding to a specific subset of trials – one through elevated firing rates and the other through specific spike timing sequences. Despite distinct coding mechanisms, an appropriate similarity metric can capture the underlying congruency in trial-dependent response modulation. This approach enables the assessment of more nuanced computational relationships among biophysically diverse neurons. In short, neurons that have the same pattern of change in their responses across trials (i.e., intrinsics spike train similarity structure) are considered to be performing similar computations. For very large neuron populations, we can simplify the population-level interpretation of the set of CS scores by projecting the high-dimensional CS matrix into a lower-dimensional coordinate space, such that each point in this latent space represents a neuron and the distance between them corresponds to their degree of computational similarity. The geometric structure of the low-dimensional neuron embedding can greatly facilitate the identification of discrete clusters and gradients of functionally related neurons. As we will demonstrate, this general framework can efficiently scale to large numbers of trials and neurons, making the method particularly suited for examining large-scale neural recording datasets that were collected under ecologically or ethologically relevant experimental conditions.

2 Methods

The flow chart presented in Figure 1 provides a general overview of the main steps used for analyzing neural ensemble data using the SIMNETS framework:

• Selecting spike trains. Time series data representing the activity of N simultaneously recorded neurons (e.g., N spike time sequences) are split into S equal-duration time segments corresponding to experimentally relevant periods of interest.

• Generating a Spike train Similarity (SSIM) matrix for each neuron. Use an appropriate metric (e.g., VP edit-distance; Victor and Purpura, 1996) to calculate the intrinsic pairwise similarities among each neuron’s set of S spike trains, resulting in a set of N SxS SSIM matrices.

• Calculating the NxN Computational Similarity (CS) matrix. Calculate the pairwise similarities (e.g., Pearson’s correlation r) among all pairs of single neuron SSIM matrices. The resulting set of CS scores is represented as a NxN population CS matrix, where each entry depicts the similarities among the output space of a neuron pair and depicts the overall pattern of relationships among the neurons for the experiment.

• Visualizing results using dimensionality reduction (DR). An appropriate DR method (t-distributed stochastic neighbor embedding [t-SNE], multidimensional scaling [MDS], etc.) can be applied to the SSIM matrices generated in step 2 to visualize the single neuron SSIM matrix as SSIM maps (example provided in Figure 2C) where each point represents an individual spike train, and the distance between them represents their relative similarity (Maaten and Hinton, 2008). This step makes it easier to interpret the neuron’s relationship to task variables (e.g., stimulus inputs, movement types, outlier trials) and ultimately their computational role within the context of higher-level network structures. Applying DR to the CS matrix generated in step 3 will produce a low-dimensional CS map, with each point corresponding to an individual neuron, and the distance between them corresponds to their relative computational similarity (i.e., the similarity between their respective SSIM matrices).

Figure 1

Figure 2

The CS map is the main output of the algorithm, capturing the computational organization of the ensemble in a latent space representation and serving as a starting point for further analysis aimed at identifying groups of neurons with similar computational properties. The following sections and figures describe specific steps (Figures 2, 3) and subsequent analysis (Figure 4) in more detail.

Figure 3

Figure 4

Amarasingham et al., 2012; Schemper, 1984 for further discussion). The silhouette values of the original data (without shuffling) are overlaid in black. The highest silhouette value (highlighted in red) is used to select the optimal number of clusters. Clusters corresponding to the highest (Ksc = 4) are shown in the inset.

2.1 Representing single neuron output spaces using SSIM matrices

One of the key innovations of the SIMNETS framework consists of using single neuron SSIM matrices as a “computational fingerprint” to represent the intrinsic structure of the output space of individual neurons. This is accomplished by performing within-neuron measures of pairwise spike train similarity across all spike trains for each neuron, resulting in an SxS single neuron SSIM matrix for each of the N neurons (Figure 1 step 2; Figure 2B). Multiple approaches to estimate the similarity between pairs of spike trains have been proposed (Paiva et al., 2010) and could be used for this step. Here, we use the spike train metric proposed by Victor and Purpura (Houghton and Victor, 2010; Victor and Purpura, 1997; Victor and Purpura, 1996; See Box 1 for description).

In the current work, each single neuron SSIM matrix represents the relationship (i.e., similarity) between all spike trains for a single neuron for the specific time windows selected for all simultaneously recorded neurons. Figure 2. presents three examples of SSIM matrices derived from the spike trains of simulated neuron responses (Linear Non-Linear Poisson functional model; Meyer et al., 2016). Entries equal to zero (black) in the matrices represent identical spike trains (zero VP cost), while higher values represent increasing differences between the spike trains. For ease of viewing, the rows/columns of the matrices have been sorted according to simulated experimental conditions. Note that a specific ordering of the spike trains is not required for the subsequent comparisons (and is only included for illustrative purposes). It can be useful to graphically illustrate the relationship between individual spike trains across the experiment captured by the single neuron SSIM matrices using dimensionality reduction (Figure 2B). Here, we use PCA-initialized tSNE (Maaten and Hinton, 2008): a combined unsupervised dimensionality reduction technique that aims to preserve the local neighborhood structure, via t-SNE, as well as information related to the global shape, via PCA-initialization (Kobak and Linderman, 2021; Lee et al., 2015). The non-random PCA-initialization also ensures reproducibility across iterations (Kobak and Linderman, 2021). This mapping step is not necessary for the following steps of the SIMNETS analysis pipeline, but these maps can be useful to visualize and compare the features coded in the population. For example, in Figure 2B, the single neuron SSIM maps illustrate the separation of orientation selectivity in simulated neuron n1, the spatial frequency selectivity in n2, and the combined spatial and orientation selectivity of n3 (See Supplementary Figure S1; supplementary material for additional discussion and examples). Changes in the shape of the SSIM output space map may also indicate the influence of other variables which are not explicitly manipulated.

2.2 Computational similarity map: comparing single neuron output spaces

Each SSIM matrix characterizes the intrinsic geometry of the output space for a specific neuron. The next step of the algorithm compares these “computational fingerprints” across the full set of recorded neurons. This second-order comparison is made by correlating each pair of SSIM matrices (Figure 1, step 3; Figure 3). Here we use Pearson’s correlation (Figure 3A). Each correlation value is entered into a symmetrical (N x N) matrix (Figure 3B) representing the relationships across neurons, which we term the Computational Similarity (CS) matrix. Applying dimensionality reduction (here, PCA-initialized t-SNE; Maaten and Hinton, 2008) to the CS matrix creates a CS map (Figure 1, step 4; Figure 3C). Note that, while each colored point in a single neuron SSIM map represents an individual spike train (Figure 2C), each point in the CS map represents a neuron (Figure 3C). In the population CS maps, the distance between neurons is a relative measure of the similarity of computational fingerprints (the SSIM matrix). Therefore, neurons performing similar computations will tend to aggregate in clusters within the population CS map of neurons.

2.3 Identifying potential subnetworks: unsupervised cluster detection & validation procedure

After neurons are embedded in the CS map, the task of finding groups of neurons representing potential neuronal subnetworks can be addressed using a wide variety of clustering algorithms. Here, we use the k-means clustering algorithm, an efficient centroid-based clustering method, however, other density-based, distribution-based, or hierarchical-clustering methods would also be appropriate. The k-means algorithm works by iteratively assigning neurons into the pre-specified number of k clusters until the optimal solution is reached. However, the challenge is to ensure that the clustering structures identified by the k-means algorithm reflect genuine clusters within the CS neuron map and not false discoveries. To address this issue, we combined the k-means algorithm with two different cluster validation techniques, a silhouette graphical analysis and our novel shuffle-based statistical test, which together help avoid false cluster discovery (Figure 4).

Silhouette analysis is a validation technique used to determine the most likely number of clusters in a dataset (Figure 4; Rousseeuw, 1987). A silhouette value represents the ratio of the between to within-cluster distances: a point within an ideal cluster will be close to members of the same cluster and far from points assigned to different clusters, resulting in a high silhouette value. Choosing the cluster number with the highest average silhouette value—the Silhouette Coefficient (SC) score—maximizes the separation between potential clusters (Figure 4C, red point). We use a novel shuffled-based resampling procedure to assess the statistical significance of the SC score to avoid false cluster discovery. The test relies on creating 1000’s of surrogate datasets through a shuffling or randomization procedure that is applied to the rows/columns of the single neuron SSIM matrices We use the matrix permutation technique used in the Mantel test, which was specifically designed to randomize pairwise distance matrices. This is accomplished by permuting both rows and columns together in order to preserve matrix symmetry. In this way, it is possible to produce a new distance matrix for each neuron where the relationship between the trials (i.e., which trials are most similar) is randomized. This procedure effectively destroys any dependencies among the neuron’s single neuron SSIM matrices and, consequently, the associated clusters with the CS neuron map (Figure 4B). The empirical SC score is only considered statistically meaningful if it falls outside a distribution of surrogate SC values calculated from the 1000’s of surrogate datasets created in this manner (Figure 4C, see blue dots; for further insights on conditional resampling; see Amarasingham et al., 2012; Schemper, 1984). This unsupervised cluster detection and validation procedure supports the identification of statistically meaningful clusters of computationally similar neurons within the CS neuron maps.

3 Results

We apply the SIMNETS analysis framework to four different datasets. First, we apply the algorithm to a population of simulated neurons for which spike patterns were constructed to simulate three subnetworks encoding information using firing rates, precise spike timing, or a combination of the two (Figure 5). Next, we apply SIMNETS analysis to three publicly available in-vivo experimental population recordings datasets from (1) macaque primary visual cortex (Kohn and Smith, 2016; Teeters and Sommer, 2009), (2) macaque primary motor cortex (Rao and Donoghue, 2014), and (3) rat CA1 hippocampal region field (Pastalkova et al., 2015; Wang et al., 2015) to illustrate the types of results the SIMNETS framework produces and how they may be leveraged for further exploratory analyses and hypothesis testing (Figures 6–11; for algorithm inputs/outputs for each dataset see: Supplementary Table S1). Finally, we highlight some of the positive attributes of the SIMNETS analysis framework and Toolbox that are particularly advantageous when analyzing large-scale neural recording data (Figure 12).

Figure 5

Figure 6

Figure 7

Figure 8

Figure 9

Figure 10

Figure 11

Figure 12

3.1 Simulated neuron population

We evaluated the performance of the SIMNETS algorithm using simulated neural data with known ground truth. This allowed us to verify that the algorithm clusters neurons with similar informational content even when they use different encoding schemes. We simulated a population of 180 neurons (N) composed of three computationally distinct subnetworks, SN1, SN2, and SN3 (Figure 5A). We generated an output structure for these subnetworks across three simulated experimental conditions (C1, C2, and C3). Each subnetwork emitted a similar “baseline pattern” for two of the conditions and a different pattern for a unique “preferred” condition (C1 for SN1, C2 for SN2, and C3 for SN3). Neurons from a common subnetwork were designed to encode similar information (i.e., respond differentially to one specific condition) using different encoding formats; for example, 20 neurons in SN1 encoded condition A through a change in spike rates (rate-coding), another set of 20 neurons in the same subnetwork responded with a change in the timing of their spikes (temporal-coding neurons), while the remaining neurons responded with changes in both spike rates and the spike time variations (mixed coding; e.g., Figure 5B). The other two subnetworks were also generated to include examples of rate coding, temporal coding, and mixed coding artificial neurons in the same proportions. We simulated 30, one-second-long spike trains for each of the 180 artificial neurons, which included 10 repetitions of each condition (S = 30 spike trains per neuron; See Supplementary Table S1 for analysis summary).

SIMNETS can be used to quantify spike train similarity and explore single unit encoding properties across an experiment without assuming what coding schemes are relevant or what variables are generating spike train differences. It is also possible to examine the temporal accuracy of single unit spiking. To illustrate this point, the SIMNETS algorithm was applied using three different temporal accuracy values for the VP spike train similarity metric: a pure rate code (q = 0, i.e., no added cost for shifting spikes in time), 100 ms (q = 10), 5 ms (q = 200). The q parameter determines how sensitive the similarity metric is to differences in the timing of individual spikes, by setting an upper limit on how far spikes can be shifted in time (See Figure 5C, for three different single neuron SSIM matrices calculated using different q values). With a setting of q = 0, the simulated neurons operating with a rate-based encoding scheme and mixed coding neurons are grouped into three functionally distinct clusters in the CS map (Figure 5D, 2nd column), while temporal coding neurons form a single cluster at the center of the map (Figure 5D, 3^rdcolumn). As the value of q increases, the VP algorithm “cost” function becomes sensitive to differences in spike timing in addition to the total number of spikes (Victor and Purpura, 1996). SIMNETS correctly groups all neurons into three distinct clusters that reflect the ground-truth functional ensemble assignments (Figure 5F). At the highest q value (Figure 5H), the CS map shows sub-groupings within each of the three ground-truth subnetworks that correspond to neuron’s encoding formats (Figure 5H, 3^rd column); however, the optimal number of clusters remains in agreement with the ground-truth functional subnetwork assignments (K_sc = 3). By specifying a higher partition value for the k-means clustering step of the algorithm (e.g., k = 9), the sub-groupings within the detected clusters are defined by the coding properties of the neurons. For a demonstration of the interaction between the cluster number and the SIMNETS hyperparameters, perplexity and q, see Supplementary Figure S2.

To evaluate whether SIMNETS could successfully identify subnetworks, we compared the distribution of Euclidean distance between individual neurons in the SIMNETS CS map for neuron pairs within and between the artificially generated (‘ground truth’) subnetworks (Figures 5E,G,I, histogram). The within-subnetwork similarity was significantly higher than between-subnetwork values in all cases (Mann–Whitney p < 0.001). For q values > 0, there was no overlap between the two distributions, indicating the complete separation of the artificially generated subnetworks. Our results demonstrate that the SIMNETS algorithm can accurately separate neurons according to their computational properties (the latent variables upon which they operate), even if they employ different coding schemes to represent information using spike rates or temporal patterns.

3.2 Application to empirical multi-neuronal recordings

We next applied SIMNETS to three distinct datasets of simultaneously-recorded, multiple single-unit recordings from three different areas in two species: (1) nonhuman primate (NHP) primary visual cortex (V1) during visual stimulation (Kohn and Smith, 2016; Smith and Kohn, 2008), (2) NHP motor cortex (M1) neural recordings during a center-out reaching task (Rao and Donoghue, 2014), and (3) rat CA1 hippocampal region recordings during a left–right alternation maze task (Pastalkova et al., 2015; Pastalkova et al., 2008). These datasets were chosen because the functional roles of the neurons in these areas have been reasonably well-established and the information content of the single neuron and ensemble spiking patterns extensively described (Gur and Snodderly, 2007; Livingstone and Hubel, 1984; Rao and Donoghue, 2014; Smith and Kohn, 2008; Truccolo et al., 2008).

We note that these results are meant to illustrate the capabilities and utilities of the SIMNETS analysis pipeline and are not intended to represent a full characterization of the role of any identified putative subnetworks. For a full analysis and interpretation, many more datasets (sessions/subjects) would be required, as would various parameter sweeps to evaluate computational processes.

3.3 Macaque V1 neural population recordings during visual stimulation

We analyzed a dataset of 112 simultaneously recorded V1 neurons (recorded using a 96-channel “Utah” electrode array) during the presentation of drifting sinusoidal gratings in an anesthetized Macaque monkey (Kohn and Smith, 2016; Figure 6A). We extracted 1 s of spiking data from the first 30 repetitions of each stimulus (S = 360), starting 0.28 s after stimulus onset (see Figure 6B for raster plots and single neuron SSIM matrices for two example neurons) to capture the neural stimulus response period (eliminating response delay). SIMNETS was applied to these NxS spike trains generated with a fixed temporal accuracy setting of 1/q = 50 ms. The resulting NxN CS neuron matrix is projected into a lower-dimensional CS space upon which further statistical and clustering analyses are performed (Figures 6D–F). An unsupervised k-means cluster analysis identified 3 clusters of related single neurons in the V1 CS map (Figures 6D,G; SC = 0.78, peak max average silhouette value), indicating that computational features did not fall along a functional continuum but instead were organized into three functionally distinct putative subnetworks.

To validate the performance of the SIMNETS algorithm, we examined the relationship between the structure of the CS map and the known tuning properties of V1 neurons. We used classic parametric tuning models to estimate the receptive field properties of all neurons and labeled the neurons in the CS neuron map according to their “preferred” stimulus orientation (Figure 6E) and tuning depth (Figure 6F). Each neuron’s preferred orientation was calculated by fitting a Gaussian distribution to the stimulus-dependent firing rates and finding the orientation that maximizes the function over the range of orientations angles (θ = [0,180]) across both left and right drift directions. Tuning depth was calculated as the normalized difference between the peak and trough of the tuning function (See Supplementary methods). Again, we emphasize that the SIMNETS procedure is an entirely unsupervised approach and does not rely on any knowledge of the experimentally-induced neuron tuning functions to generate the CS map or to identify putative CS neuron clusters. We simply labeled each neuron with values obtained from its relationship to an experimental variable (here, grating direction) after the neurons were plotted in the SIMNETS CS map. A pairwise circular-linear correlation (r_cl) statistical analysis revealed a strong positive correlation between neurons’ distance in the low-dimensional CS neuron map and the differences in their preferred orientation parameters (Pearson, r_cl = 0.89). This relationship confirmed that neurons with similar preferred orientation were generally mapped near one another in the CS map (Figure 6E). The exception to this trend was observed for neurons with weaker tuning strengths (Figure 6F, brown points), which aggregated at the center of the CS neuron map. Overall, the distribution of the neurons’ preferred orientations corresponded well with the three observed clusters in the CS neuron space (Figure 6H). Collectively, these results confirm that the neurons’ sensitivity to the different stimulus orientations is a major contributor to the structure of V1 computational architecture. The CS map also reveals the presence of three computationally distinct groups of neurons within the context of this experimental paradigm.

The overall agreement between known V1 tuning features and the CS map generated using SIMNETS suggests that the algorithm will be able to identify computational architectures for neural systems where the tuning functions are not known a priori.

The DR techniques used to visualize the computational geometry representations of the single neuron SSIM maps (introduced in Figure 2) can be expanded to represent the combined output space of multiple neurons simultaneously, such as the CS neuron clusters (i.e., putative subnetworks) or the full neuron population. For multi-neuron SSIM maps, each point corresponds to the collective spiking pattern generated by all neurons within an ensemble (e.g., CS neuron cluster or full population) on a single trial, and the distance between points corresponds to the difference (spike train distance) in the neurons collective spike pattern across time. Here, we use this approach to visualize the relationship between activity patterns for the full neuron ensemble (Figure 7B), the three CS neuron clusters (Figure 7C), and example single neurons from each of the three CS neuron clusters (Figure 7D). Standard single neuron plots (parametric tuning functions and rasters) were also included for comparison.

3.4 Macaque M1 neural population recordings center-out reaching task

As a second example, we applied the SIMNETS algorithm to a dataset of 103 M1 neurons recorded using a 96-channel electrode array in a macaque performing a planar, 8-direction arm reaching task (see Methods for more details). We extracted 1-s long spike trains (S = 114) from each neuron during all trials where the monkey successfully reached the cued target, starting 0.1 s before movement onset (Figure 8A; see Supplementary Table S1 for algorithm inputs and parameters). SIMNETS was applied to these NxS spike trains generated with a fixed temporal accuracy setting of 1/q = 50 ms. SIMNETS analysis revealed that the neurons were organized into three neuron clusters in the M1 CS map (Figures 8D,G, SC score = 0.71, M1 peak average silhouette).

To examine the relationship between these clusters and functional properties of the neurons, we characterized each neuron’s preferred reach direction (PD) angle and tuning depth (Figure 8F) in a similar manner to the V1 dataset. Each neuron’s PD was estimated by fitting a von Mises distribution to the firing rates as a function of direction (Supplementary Figure S3; Mardia and Zemroch, 1975) A circular-linear correlation (rcl) analysis between the difference in PD and CS map distance between pairs of neurons revealed a significant positive relationship (Pearson, rcl = 0.92; p = 0.001), confirming that neurons with similar preferred directions tended to cluster together in the CS map, while neurons with dissimilar preferred directions were spatially distant (Figures 8D,E, right plot). After plotting the distribution of preferred reach direction, we again observed that the neuron tuning properties were similar across each of the three detected CS neuron clusters (Figure 8H).

SIMNETS analysis suggests that the computational properties of M1 neurons are not uniformly distributed along a functional continuum, as evidenced by the distinct CS neuron clusters identified in the population CS space. The separation between the CS neuron clusters is most apparent for CS neuron cluster 2 and cluster 3 (Figure 8D, orange and yellow CS clusters), which mirrors the discontinuity in the distribution of preferred reach angles (Figure 8E). These results agree with previous findings, which support the hypothesis that the biomechanical constraints of the limb are reflected in an uneven distribution of preferred directions among motor cortical neurons (Lillicrap and Scott, 2013).

As with the V1 dataset, we used a multi-scale latent space analysis to visualize the relationship between spiking patterns across different scales within the CS map (Figure 9A), from the full ensemble (Figure 9B) to CS neuron clusters indicative of potential subnetworks (Figure 9C) and individual neurons functional properties (Figure 9D). Standard single neuron plots (parametric tuning functions and raster plots) are again included for comparison (Figure 9C). Each M1 neuron had a unique single neuron SSIM map that captured several computationally relevant features of the neurons’ spike train outputs (Figure 9C; SSIM matrices/maps for two example neurons). In contrast to the regular conical and circular parabolic shapes of the V1 single neuron SSIM spaces (Figure 7C), the M1 SSIM maps displayed more irregular and heterogenous global geometries (Figure 9D).

3.5 SIMNETS subnetworks in a population of rat hippocampal CA1 neurons

Lastly, we assessed the CS structure of a population of hippocampal neurons using SIMNETS, a region that is both functionally and structurally quite different from the primary visual and motor neocortical areas examined above. We applied SIMNETS to a publicly available dataset of rat CA1 hippocampal neurons recorded using multi-site silicon probes in a rat performing a left/right-alternation navigation task in a “figure-8” maze (Pastalkova et al., 2015; Figure 10A). All neurons with firing rates greater than 5 Hz were included in the analysis (N = 80/106 recorded neurons). The rat performed 17 correct trials (T = 17, trials), taking an average of 4.3 s to reach the reward location at the end of the left or right arm of the maze. Spike trains used for the SIMNETS algorithm were obtained by extracting 1-s of spiking activity for each 80 mm increment of distance traveled by the rat (Figure 10B, left, example single neuron raster plots with S spike trains; see Supplementary methods). This resulted in a set of S = 442 spike trains (17 trials x 26 time-epochs) per neuron.

SIMNETS was applied to the set of SxN neurons, using the VP temporal sensitivity setting of q = 35 (Figures 10B,C; see Supplementary Table S1 for all algorithm inputs/outputs). After generating the low-dimensional CA1 population CS neuron map, we characterized the organizing properties of the neuronal relationships using the unsupervised k-mean cluster detection/validation procedure. This procedure identified Ksc = 5 statistically significant CS neuron clusters within the CA1 CS neuron map (SC = 0.70, peak max average silhouette value; Figure 10D), suggesting that the neurons are organized into five computationally distinct clusters of neurons. To verify this, we used two different classic functional characterization approaches to assign functional ids to the neurons within the SIMNETS CS map and assess their functional organization within the context of the CS neuron map (Supplementary Figure 10F). Specifically, we categorized the CA1 neurons as having place cell-like activity (n = 60, place-cells, PCs) or lacking place cell-like activity (n = 20, non-place-cells, non-PCs) using an information-theoretic measure (Skaggs et al., 1992) and a spatial firing rate tuning functions (Figure 10B, example single neuron raster plots; see Supplementary methods for analysis details). As is reported in other studies, some of the neurons categorized as place-cells exhibited complex firing responses, such as multiple peak place-dependent firing responses (Supplementary Figure S3D, show example single neuron response functions). An examination of the physiological properties of the identified k-means neuron clusters revealed that one cluster was almost entirely (95%) composed of non-PCs (Figures 10E,F, K-means cluster c3, purple), while the other four k-means CS neuron clusters were entirely or almost entirely (>92%) composed of PC’s (Figures 10E,F). This indicated one of at least three possibilities: (1) the “non-PC CS cluster,” cluster c3 corresponded to a cluster of untuned, “uninformative” neurons, (2) that the individual neurons computed a common task-related signal that was unrelated to place-coding, or (3) that place-coding signals were only apparent at the level of the neuronal cluster rather than the individual neurons (Figure 10E, cluster c3). While a full exploration of these different possibilities is beyond the scope of this work, we do provide a demonstration of how the SIMNETS CS neuron maps may be used to visualize the relationship between spiking patterns across scales for hypothesis formation and testing (Figure 11).

Once again, we used a multi-scale latent space analysis to visualize the relationship between spiking patterns across different scales within the CS map (Figure 11A), from the full ensemble (Figure 11B) to CS neuron clusters indicative of potential subnetworks (Figure 11C) and individual neurons functional properties (Figure 11D). As expected, the spatial layout of the task-space was reflected in the topology of the population-level SSIMS map (Figure 11B), while each cluster-level SSIMS maps captured the place-dependent modulation of the clusters spiking patterns as the rat traversed a specific region of the maze (Figure 11C). Surprisingly, the “non-PC neuron cluster,” (Figure 11C, cluster c3) did not appear to reflect a cluster of non-responsive or “uninformative” neurons; instead, it had a cylindrical shape that resulted from at least two dominant patterns of spike pattern modulation: a within-trial periodic spiking pattern in the x-y plane (Figure 11B, cluster c3, x-y viewing angle) and a location or distance-dependent variance along the z-axis (See Supplementary Figure S3C, for x-z viewing angle). This cylindrical shape was also observed in the single neuron SSIM maps of the neurons within cluster c3 (see Figure 11C, purple box/stars, neurons c3-n7 and c3-n20), however, the distance-dependent x-z variance was less apparent at the single neuron level. Overall, this suggests that the activity of the non-PC neuron cluster, c3, displays dynamics that may reflect a different task variable that is in some way correlated with the rat’s performance on the task (e.g., head direction, speed profile along maze arms) or, alternatively, the intrinsic circuit dynamics (Pastalkova et al., 2008; Richard et al., 2013; Wang et al., 2015). Latent spaces displaying 3-dimensional ring topology reflecting animal position and running direction encoding have been recently described in CA1 (Esparza et al., 2025). Further investigation of this phenomenon is outside of the scope of the present work but is an example of how SIMNETS may be used to quickly uncover different types of functional relationships among neural representations that might not be otherwise evident when using classic population and single unit analyses. These results highlight the advantages of applying SIMNETS to neural recordings where tuning properties or other coded variables are not readily apparent or unknown a priori.

3.6 Computational efficiency and analysis run-time

The SIMNETS framework can rapidly analyze large numbers of neurons using a standard personal computer (e.g., 6-Core Intel Core i7, 64GB of RAM). With the current Matlab® implementation (see included SIMNETS software package), the simulated neuron population (N = 180, S = 30) and Macaque M1 neuron population (N = 103, S = 360) were processed in less than 2 s. Larger datasets, such as the Macaque V1 neuron population (N = 112; S = 360) and rat hippocampal CA1 neuron population (N = 80, S = 442), take approximately 10 s to process. This is significantly faster than a typical cross-correlation-based functional connectivity analysis. For example, a dataset of 1,000 neurons can take an hour to processes with a standard cross-correlation analysis (Figure 12, blue line; utilzing Matlab® xcorr function), whereas our implementation of the SIMNETS algorithm takes only about a minute to process the same dataset (Figure 12, black line). The SIMNETS analysis is faster because it involves fewer spike train comparison operations (NxSxS) than the cross-correlations analysis (N2xSxS, across all time-lags), and the number of spike train comparisons per added neuron has a favorable linear scaling. If parallel or distributed computing resources are used to calculate the single neuron SSIM matrices, the execution time can be reduced even further to seconds or milliseconds. This would be especially useful for datasets with larger trial numbers (e.g., S > 1,000), which are more computationally intensive than datasets with larger neuron numbers.

4 Discussion

Advances in multi-electrode recording technology have made it possible to record or image the activity of thousands of individual neurons simultaneously (Chung et al., 2022; Hsieh and Zheng, 2019; Paulk et al., 2022; Steinmetz et al., 2021; Vogt, 2019). The shift in focus from studying single neurons to population levels responses has brought about a new era of discovery, innovation, and collaboration (Buzsáki, 2004; Cunningham and Yu, 2014; Hurwitz et al., 2021; Kohn et al., 2020; Paninski and Cunningham, 2018; Paulk et al., 2022; Saxena and Cunningham, 2019; Urai et al., 2022); highlighting the need for new theoretical and analytical frameworks that can bridge the “representational gap” between neuron-level and population-level computation.

With this problem in mind, we devised SIMNETS, a novel relational analysis framework. SIMNETS was designed to generate low-dimensional maps of computational similarity across populations of simultaneously recorded neurons. The critical difference between SIMNETS and traditional functional connectivity measures is that we assess similarities among the neurons’ intrinsic latent space geometries, rather than the similarities of the neurons’ response time course (e.g., correlation between firing rates). We emphasize that in contrast to neuronal functional connectivity analyses that measure time series correlations among the firing rates vectors or spiking patterns between different neurons (Abeles and Gat, 2001; Aertsen and Gerstein, 1985; Gerstein and Michalski, 1981; Grün et al., 2002b; Humphries, 2011; Kiani et al., 2015; Stuart et al., 2002) SIMNETS only directly compares spike trains generated by the same neuron, i.e., self-similarity. We performed spike train comparisons using the VP metric, a point process edit distance measure that can be tested over a variety of temporal accuracy settings. By quantifying the similarity of a neuron’s spike train outputs across all trials (a within-neuron 1^st-order comparison), we generate a geometric representation of each neuron’s spike train output space, which we use as a “computational fingerprint” which is captured in its spike train similarity matrix (Figure 2). This approach takes into account changes in spiking related to all the variables in the experiment (controlled, uncontrolled, and noise) in an unsupervised manner, that is, without the use of trial labels or neuron functional labels (Kohn and Smith, 2016; Pastalkova et al., 2015; Rao and Donoghue, 2014). SIMNETS projects the neuron population into a common coordinate space that captures the relationships between the neuron’s computational fingerprints. The distance between neurons in the CS space provides a quantifiable measure of computational similarity, and CS neuron clusters reveal potential subnetworks with distinct computational properties. We used a dimensionality reduction method that aims to preserve both local and global structure (e.g., PCA initialized t-SNE; Maaten and Hinton, 2008; Kobak and Linderman, 2021; Lee et al., 2015). The structure of the CS maps can provide additional insight into the type of information that is being computed potential and subnetwork organization. SIMNETS could help resolve long debated issues as to whether classes of neurons that represent a unique set of features belong to a discrete computing module, or alternatively, whether there exists a functional gradient of computational properties.

We see the novelty of the SIMNETS framework as combining similarity metrics, dimensionality reduction and clustering to organize neurons into functional groups. The specific combination of algorithms chosen for the implementation of SIMNETS presented here were selected because they have been thoroughly validated (on an individual basis) through previous work. However, it is important to note that each of these steps can be accomplished using a variety of methods and could potentially benefit from recent developments in each of these areas. Other analysis pipelines could also benefit from employing SIMNETS as a preprocessing step. For example, applying DR techniques to examine the low dimensional latent spaces associated with clusters in the SIMNETS CS map provides another way to probe the potential computational properties of putative subnetworks (Figures 7, 9, 11). This strategy could take advantage of recent advances in DR methods that incorporate tensor analysis (slice TCA, Pellegrino et al., 2024), ensemble dynamics (NoMAD, Karpowicz et al., 2025), or explicitly include task variables to enhance low dimensional representations (DAD, Dyer et al., 2017; CEBRA, Schneider et al., 2022).

Finally, the combination of a short processing time (< 1 min per 1,000 neurons), a computational complexity that scales near-linearly with the size of the neuron population (Figure 12), and ease of parallelization, makes SIMNETS an efficient tool for exploring very large-scale neuron populations with readily available hardware.

4.1 Algorithm validation

We tested and validated SIMNETS using simulated data with known ground truth as well as three different experimental neural recording test datasets to demonstrate the general applicability of the analysis framework across a range of coding schemes, neural systems (visual, motor, cognitive), behaviors (awake or anesthetized), and different species (rat and macaque). This included a demonstration of a novel shuffle-based statistical procedure based on the Mantel test for identifying non-spurious CS neuron clusters within a low-dimensional population CS neuron map.

Our analysis of simulated data with ground-truth subnetworks and known single neuron properties (i.e., encoded content and encoding timescale) demonstrates the current implementation of the SIMNETS algorithm can capture a neuron’s local and global computational relationships, even under conditions where computationally similar neurons relied on very different encoding schemes (e.g., rate coding vs. fine temporal coding schemes). This agnosticism can be beneficial when analyzing neurons that generate highly heterogeneous spike train outputs because of differences in their internal machinery (e.g., inhibitory vs. excitatory neurons) yet still reflect task-relevant information (Churchland and Shenoy, 2007; Diba et al., 2014; Martina et al., 1998). We demonstrated how the VP metric temporal accuracy parameter (q) could be used to interrogate the effect of temporal resolution on the computational configuration of the neurons.

The application of SIMNETS to three publicly available multi-neuron recording datasets validated the capabilities of the method in revealing known dominant feature representation in V1, M1, and CA1 (orientation, direction, place, respectively) without imposing stimulus or movement-driven parametric tuning models a priori. Of note, the results show that the computational properties of neurons are not uniformly distributed, but instead are organized into apparent subnetworks (i.e., CS clusters in the low-dimensional CS map) that emphasize certain collections of feature representations. Our results also suggest that SIMNETS may be able to detect distinct subnetwork structures hypothesized to support ensemble place-coding or complex feature conjunctions (Buzsáki, 2019; Diba et al., 2014). Although it was beyond the scope of this report to demonstrate the functional or computational significance of the detected putative subnetworks, our results suggest that the detected clusters are physiologically meaningful; our subsequent unsupervised cluster validation statistical procedure further provided support to these findings. Finally, our choice of datasets allowed us to demonstrate that this method generalizes well to neural recordings from a variety of brain regions, including sensory areas (Kohn and Smith, 2016), motor areas (Rao and Donoghue, 2014), and memory/cognitive areas (Pastalkova et al., 2008); different species, including rat and non-human primate; and across different recording technologies including laminar silicone probes (Buzsáki, 2004) and fixed-dept multi-electrode arrays (Maynard et al., 1997).

Collectively, our results demonstrate that single neuron SSIM matrices provide a simple, scalable, and powerful format for characterizing the computational properties of neurons, as well as identifying computational similarities and differences between them. In addition, examining the arrangement of neurons in the CS maps may suggest novel coding schemes that are yet unexplored by experimenters, generating new hypotheses for future experiments.

4.2 Comparison to existing approaches

The SIMNETS analysis framework builds upon rich theoretical and mathematical literature spanning multiple domains and disciplines. Geometric and metric space representational models of similarity have a long history of application in the field of psychology where they have been used to model the perceptual relationships between sensory stimuli as latent or low-dimensional perceptual metric-space (Forstmann, 2011; Goldstone et al., 1991; Hyman, 1974; Shepard and Chipman, 1970; Shepard, 1987; Shepard, 1964; Thurstone, 2017; Tversky, 1977; Zenker and Gärdenfors, 2015). This approach has recently been adapted to study the intrinsic structure of neural representational spaces in visual, motor, and cognitive brain regions (for reviews, see; Edelman et al., 1998; Haxby, 2012; Kriegeskorte et al., 2008; Kriegeskorte and Kievit, 2013; Lehky et al., 2013). Several mathematical variations of this framework have been developed to study the information content of macroscale fMRI BOLD signals (Diedrichsen and Kriegeskorte, 2017; Haxby et al., 2001; Kriegeskorte and Bandettini, 2007), and both population-level and neuron-level spiking patterns (Houghton and Victor, 2010; Kiani et al., 2007; Vargas-Irwin et al., 2015a, 2015b; Victor and Purpura, 1996). To the best of our knowledge, this is the first work that uses a measure of 2nd-order point-process similarity analysis to generate latent space embedding of a collection of individual neurons.

The concept of a low-dimensional embedding that captures the functional relationship between spiking neurons was introduced in the seminal papers by (Baker and Gerstein, 2000; Gerstein and Michalski, 1981) describing the use of “Gravitational Clustering” (Wright, 1977): a neuron clustering and visualization tool for identifying groups of neurons with synchronous spiking patterns. This method is based on an analogy of the physics of the gravitational forces governing the dynamics and interactions of macroscopic particles. It treats the N neurons as N particles moving within an N-dimensional space, where charges that influence the attractive and repulsive interactions between particles are dictated by the temporal dynamics of pairwise synchronous spiking activity between neurons. The result is a visualization of particle clusters (and their trajectories) that represent dynamically evolving assemblies of synchronously active neurons. More recent work by Kiani et al. (2015) involves the application of the dimensionality reduction techniques (i.e., multidimensional scaling) to pairwise measures of between-neuron spike rate covariations to detect natural functional modules in a population of pre-arcuate and motor cortex neurons (Kiani et al., 2015). Although the goal of this method is like that of SIMNETS—to the extent that it makes use of single-trial information to group neurons in an unsupervised manner—this approach is applied to binned spike rates and appears to capture the covariation information carried in the absolute firing rates of the neurons’ responses across trials (Kiani et al., 2015), rather than the information carried in the neurons’ intrinsic spike train relational geometry. However, our analysis of the synthetic neurons using temporal codes, as well as rat CA1 neurons demonstrate that SIMNETS can detect a broader range of meaningful information processing motifs that reflect both condition-independent and condition-dependent information carried in the fine temporal structure of spike train outputs.

Yang et al. (2019) identified functional groupings of artificial neurons in a neural network by applying t-SNE to activation-based measures of single neuron selectivity across 20 different tasks (i.e., task variance; Yang et al., 2019). The method is similar to the SIMNETS method in that it projects neurons onto a CS map according to the similarity of measures of output variance. This method requires knowledge of trial labels within tasks to group neurons according to the similarities of their changes in selectivity across tasks. By contrast, SIMNETS does not require a priori knowledge about the similarities or differences of the information encoded on certain trials (i.e., trial labels), only that the trials are recorded simultaneously. This feature of SIMNETS is particularly valuable when trying to identify groups of functionally similar neurons in experiments with awake and freely moving animals, where the assumption of repeatable perceptual, cognitive, or behavioral states across trials is not always possible.

Several previous studies have used pairwise comparisons between the spike trains of different neurons to identify putative subnetworks or cell assemblies (Aertsen et al., 1989; De Blasi et al., 2019; Gerstein et al., 1985; Grün et al., 2002a; Humphries, 2011; Roudi et al., 2015; Singer and Gray, 1995). As with the gravitational clustering method, these studies have operated under the working hypothesis that the detection of time-series spike time or rate correlations is a signature of a potential functional link (Cole et al., 2016; Roudi et al., 2015; Singer and Gray, 1995; Ts’o et al., 1986). By contrast, SIMNETS identifies neurons with similar informational content even if they exhibit heterogeneous spike train outputs or utilize different encoding timescales (e.g., rate vs. precise spike timing). For the simulated neuron population, SIMNETS successfully identified clusters of the simulated neurons according to their ground-truth functional subnetworks and by increasing the sensitivity of the VP metric, the CS map was able to reveal further sub-groupings within the three identified functional CS clusters that highlighted subtle differences their output space geometries that corresponded to the different encoding timescale. This feature of SIMNETS could be particularly useful for determining if neurons encoded information across heterogenous timescales, e.g., fast vs. slow-timescale inhibitory neurons, or if they compute across a hierarchy of timescales (Cusinato et al., 2022; Fortenberry et al., 2012; Gorchetchnikov and Grossberg, 2007; Honey et al., 2012). Additionally, this feature of SIMNETS enabled us to cluster physiologically and functionally distinct groups of neurons in the CA1 dataset (i.e., place-cells vs. non-place cells), which were not obvious when VP metric was insensitive to spike times (i.e., VP metric q = 0; data not shown).

Capturing the computational similarity between pairs of neurons in a more general way could be accomplished using Information theoretic approaches (Shannon, 1997), which could include estimating the shared or mutual information between the output spike trains of different neurons. Other approaches focus on the asymmetry of the predictive power between variables at different lags, resulting in “directed” estimates of functional connectivity such as Granger Causality, Transfer Entropy, or the Directed Transfer Function (Granger, 1983; Ito et al., 2011; Ursino et al., 2020). These strategies are based on estimating joint probability distributions across the activity patterns of pairs of neurons. A similar approach can be applied to relationships between multiple neurons using Generalized Linear Models (Chen et al., 2009; Truccolo et al., 2005). The number of possible activity patterns is very large, so this type of calculation can be challenging even when relatively large amounts of data are available. Using large amounts of data presents additional problems for this approach because it assumes that relationships between neurons remain relatively constant across the time span used to fit the models. Hence, every additional neuron requires exponentially more data and computing power and analyzing hundreds to thousands of simultaneously recorded neurons becomes intractable.

4.3 Mitigating experimenter bias

A considerable amount of research in systems neuroscience has focused on identifying new classes of neurons based on their information-processing properties. The standard approach for many of these experiments involves recording single unit activity while a certain experimental variable of interest is manipulated (for example, providing systematic changes in stimulus features, or eliciting different behavioral responses; Meyer et al., 2016). Standard statistical tests (ANOVA, etc.) are then used to determine if each neuron displays significant changes in firing rate across the experimental conditions. The percentage of significant neurons is usually reported as a functionally distinct “class” of neurons sensitive to the variable of interest. It is common to exclude neurons that do not reach statistical significance or cannot be fit using a predetermined tuning model from further analysis. This approach is prone to both selection and confirmation bias, and ultimately produces “classes” of neurons identified based on arbitrary statistical thresholds imposed on what are likely continuous distributions of properties. The SIMNETS analysis framework is an unsupervised approach to determine if neurons are organized across a functional continuum or are organized into statistically separable functional classes, thereby mitigating the experimenter’s bias inherent in parametric neural discovery methods. The SIMNETS approach is different in that it seeks to organize the neurons into a low dimensional computational similarity (CS) map solely based on the intrinsic structure of their firing rates, without making any a priori assumptions about their tuning properties. This makes it possible to detect clusters of neurons with similar informational content and can also reveal gradients of functional properties if they exist. This is not to say that experimental design does not play an important role: the informational content of neuronal spiking is still evaluated within the specific context of the recording. However, SIMNETS does not target specific dependent variables, but instead tries to map the full population of neurons onto a single unified CS map. It is therefore possible to identify ‘orphan’ clusters of neurons that display variations in their spiking related to uncontrolled variables. While assigning a function to an orphan cluster related to an uncontrolled variable may be challenging, it highlights the potential impact of additional variables in the neuronal population, and presents the opportunity to design additional experiments and potentially identify new classes of neurons. In addition to providing a principled way to determine if a consistent organization of information processing modules can be found across sessions and subjects, we believe that the ability to intuitively visualize relationships within networks of neurons will provide a unique perspective leading to new data-driven hypotheses and experimental refinement.

4.4 Limitations and considerations

Several important limitations of SIMNETS are worth noting. First, estimates of similarity using spike train metrics require that the time windows of interest be of equal length, making it difficult to compare neural responses with different time courses. This weakness is common to all trial-averaging models commonly used in the literature that we are aware of. It is important to note that the optimal length of the time window will depend on the specific experimental context and should be adjusted on a case-by-case basis, either based on previous work (as done here) or empirically. Using spike train similarity metrics makes it possible to use relatively long time windows without losing information encoded at fine temporal resolution. This feature makes it easier to adopt time windows that fully encompass behavioral events of interest.

Second, although the SIMNETS framework does not require a priori assumptions about the variables potentially encoded by neural activity, experimental design and data selection will still have a direct effect on the results obtained. For example, a set of neurons identified as a functional subnetwork could separate into smaller groups with different computational properties when additional task conditions are added to the analysis. Thus, the functional properties identified using SIMNETS are only valid within the context of the data examined and may not necessarily extrapolate to different experimental conditions. Ultimately, the computational space of the population is only as rich as the experimental design allows (Gao and Ganguli, 2015).

Third, SIMNETS does not provide a means to assess sources of variance. For example, subsets of neurons may exhibit systematic changes in spike train outputs through an experimental session, reward experience, timing, attention, or many other variables which can be revealed by the neural data; however, it will be up to the experimenter to identify their source. Nevertheless, SIMNETS can guide hypothesis generation for experiments that identify these influences.

Fourth, although peak silhouette values are indicative of the optimal number of clusters according to clustering statistics, they should be taken as a guideline and not an absolute measure of the number of subnetworks. The present study focused on a set of diverse recordings from various sources to highlight the performance of the SIMNETS algorithm for different types of data. An in-depth study of potential subnetwork structure using SIMNETS should include many similar datasets, ideally sampling multiple subjects and sessions to properly assess the number of clusters present across many CS maps before drawing conclusions regarding network architecture.

Fifth, it is possible that neuronal subnetworks are constantly re-arranged depending on changing external or internal signals, e.g., ethological demands, task demands, or task stages, or attentional levels. For example, a neuron could potentially exhibit rapid changes in its computational/functional interrelationships with other neurons during processing epochs that are shorter than the selected analysis time window. The current version of the SIMNETS algorithm was not designed to distinguish between rapidly changing sub-network memberships. However, applying the SIMNETS algorithm multiple times over different epochs can help with determining if different network configurations are engaged at different times.

References

• 1

  AbelesM.GatI. (2001). Detecting precise firing sequences in experimental data. J. Neurosci. Methods107, 141–154. doi: 10.1016/s0165-0270(01)00364-8

  • CrossRef
  • Google Scholar
  • View reference in article

• 2

  AertsenA. M.GersteinG. L. (1985). Evaluation of neuronal connectivity: sensitivity of cross-correlation. Brain Res.340, 341–354.

  • Google Scholar
  • View reference in article

• 3

  AertsenA. M.GersteinG. L.HabibM. K.PalmG. (1989). Dynamics of neuronal firing correlation: modulation of effective connectivity. J. Neurophysiol.61, 900–917.

  • Google Scholar
  • View reference in article

• 4

  Agüera y ArcasB.FairhallA. L.BialekW. (2003). Computation in a single neuron: Hodgkin and Huxley revisited. Neural Comput.15, 1715–1749. doi: 10.1162/08997660360675017

  • CrossRef
  • Google Scholar
  • View reference in article

• 5

  AlivisatosA. P.ChunM.ChurchG. M.DeisserothK.DonoghueJ. P.GreenspanR. J.et al. (2013). The brain activity map. Science339, 1284–1285. doi: 10.1126/science.1236939

  • CrossRef
  • Google Scholar
  • View reference in article

• 6

  AmarasinghamA.HarrisonM. T.HatsopoulosN. G.GemanS. (2012). Conditional modeling and the jitter method of spike resampling. J. Neurophysiol.107, 517–531. doi: 10.1152/jn.00633.2011

  • CrossRef
  • Google Scholar
  • View reference in article

• 7

  AronovD.VictorJ. D. (2004). Non-Euclidean properties of spike train metric spaces. Phys. Rev. E Stat. Nonlin. Soft Matter Phys.69:061905.

  • Google Scholar

• 8

  BakerS. N.GersteinG. L. (2000). Improvements to the sensitivity of gravitational clustering for multiple neuron recordings. Neural Comput.12, 2597–2620. doi: 10.1162/089976600300014863

  • CrossRef
  • Google Scholar
  • View reference in article

• 9

  BermanN. E.WilkesM. E.PayneB. R. (1987). Organization of orientation and direction selectivity in areas 17 and 18 of cat cerebral cortex. J. Neurophysiol.58, 676–699.

  • Google Scholar
  • View reference in article

• 10

  BriggmanK. L.AbarbanelH.DiKristanW. B. (2006). From crawling to cognition: analyzing the dynamical interactions among populations of neurons. Curr. Opin. Neurobiol.6, 135–44. doi: 10.1016/j.conb.2006.03.014

  • CrossRef
  • Google Scholar
  • View reference in article

• 11

  BrodyC. D. (1999). Correlations without synchrony. Neural Comput.11, 1537–1551.

  • Google Scholar
  • View reference in article

• 12

  BrownE. N.KassR. E.MitraP. P. (2004). Multiple neural spike train data analysis: state-of-the-art and future challenges. Nat. Neurosci.7, 456–461. doi: 10.1038/nn1228

  • CrossRef
  • Google Scholar
  • View reference in article

• 13

  BuzsákiG. (2004). Large-scale recording of neuronal ensembles. Nat. Neurosci.7, 446–451. doi: 10.1038/nn1233

  • CrossRef
  • Google Scholar
  • View reference in article

• 14

  BuzsákiG. (2019). “Internally Organized Cell Assembly Trajectories,” In The Brain from Inside Out. New York: Oxford University Press, p. 165–98.

  • Google Scholar
  • View reference in article

• 15

  ChenZ.PutrinoD. F.BaD. E.GhoshS.BarbieriR.BrownE. N. (2009). A regularized point process generalized linear model for assessing the functional connectivity in the cat motor cortex. Annu Int Conf IEEE Eng Med Biol Soc2009, 5006–5009. doi: 10.1109/IEMBS.2009.5334610

  • CrossRef
  • Google Scholar
  • View reference in article

• 16

  ChungJ. E.SellersK. K.LeonardM. K.GwilliamsL.XuD.DoughertyM. E.et al. (2022). High-density single-unit human cortical recordings using the Neuropixels probe. Neuron110, 2409–2421.e3. doi: 10.1016/j.neuron.2022.05.007

  • CrossRef
  • Google Scholar
  • View reference in article

• 17

  ChurchlandM. M.ShenoyK. V. (2007). Temporal complexity and heterogeneity of single-neuron activity in premotor and motor cortex. J. Neurophysiol.97, 4235–4257. doi: 10.1152/jn.00095.2007

  • CrossRef
  • Google Scholar
  • View reference in article

• 18

  CohenM. R.KohnA. (2011). Measuring and interpreting neuronal correlations. Nat. Neurosci.14, 811–819. doi: 10.1038/nn.2842

  • CrossRef
  • Google Scholar
  • View reference in article

• 19

  ColeM. W.YangG. J.MurrayJ. D.RepovšG.AnticevicA. (2016). Functional connectivity change as shared signal dynamics. J. Neurosci. Methods259, 22–39. doi: 10.1016/j.jneumeth.2015.11.011

  • CrossRef
  • Google Scholar
  • View reference in article

• 20

  CunninghamJ. P.YuB. M. (2014). Dimensionality reduction for large-scale neural recordings. Nat. Neurosci.17, 1500–1509. doi: 10.1038/nn.3776

  • CrossRef
  • Google Scholar
  • View reference in article

• 21

  CusinatoR.AlnesS. L.van MarenE.BoccalaroI.LedergerberD.AdamantidisA.et al. (2022). Intrinsic neural timescales in the temporal lobe support an auditory processing hierarchy. The Journal of Neuroscience: The Official Journal of the Society for Neuroscience. 43, 3696–3707. doi: 10.1101/2022.09.27.509695

  • CrossRef
  • Google Scholar
  • View reference in article

• 22

  DavidS. V.VinjeW. E.GallantJ. L. (2004). Natural stimulus statistics alter the receptive field structure of v1 neurons. J. Neurosci.24, 6991–7006. doi: 10.1523/JNEUROSCI.1422-04.2004

  • CrossRef
  • Google Scholar
  • View reference in article

• 23

  De BlasiS.CibaM.BahmerA.ThielemannC. (2019). Total spiking probability edges: a cross-correlation based method for effective connectivity estimation of cortical spiking neurons. J. Neurosci. Methods312, 169–181. doi: 10.1016/j.jneumeth.2018.11.013

  • CrossRef
  • Google Scholar
  • View reference in article

• 24

  de la RochaJ.DoironB.Shea-BrownE.JosićK.ReyesA. (2007). Correlation between neural spike trains increases with firing rate. Nature448, 802–806. doi: 10.1038/nature06028

  • CrossRef
  • Google Scholar
  • View reference in article

• 25

  DibaK.AmarasinghamA.MizusekiK.BuzsákiG. (2014). Millisecond timescale synchrony among hippocampal neurons. J. Neurosci.34, 14984–14994. doi: 10.1523/JNEUROSCI.1091-14.2014

  • CrossRef
  • Google Scholar
  • View reference in article

• 26

  DiedrichsenJ.KriegeskorteN. (2017). Representational models: a common framework for understanding encoding, pattern-component, and representational-similarity analysis. PLoS Comput. Biol.13:e1005508. doi: 10.1371/journal.pcbi.1005508

  • CrossRef
  • Google Scholar
  • View reference in article

• 27

  DonoghueJ. P. (2002). Connecting cortex to machines: recent advances in brain interfaces. Nat. Neurosci.5, 1085–1088. doi: 10.1038/nn947

  • CrossRef
  • Google Scholar
  • View reference in article

• 28

  DyerE. L.AzarM. G.PerichM. G.FernandesH. L.NaufelS.MillerL. E.et al. (2017). A cryptography-based approach for movement decoding. Nat. Biomed. Eng.1, 967–976. doi: 10.1038/s41551-017-0169-7

  • CrossRef
  • Google Scholar
  • View reference in article

• 29

  EdelmanS.Grill-SpectorK.KushnirT.MalachR. (1998). Toward direct visualization of the internal shape representation space by fMRI. Psychobiology26, 309–321.

  • Google Scholar
  • View reference in article

• 30

  EsparzaJ.QuintanillaJ. P.CidE.MedeirosA. C.GallegoJ. A.de la PridaL. M. (2025). Cell-type-specific manifold analysis discloses independent geometric transformations in the hippocampal spatial code. Neuron. 113, 1098–109.

  • Google Scholar

• 31

  FernandezJ. M.FarellB. (2009). Is perceptual space inherently non-Euclidean?J. Math. Psychol.53, 86–91.

  • Google Scholar

• 32

  ForstmannB. U. (2011). Reciprocal relations between cognitive neuroscience and formal cognitive models: opposites attract? Trends Cogn. Trends Cogn. Sci.15, 272–279. doi: 10.1016/j.tics.2011.04.002

  • CrossRef
  • Google Scholar
  • View reference in article

• 33

  FortenberryB.GorchetchnikovA.GrossbergS. (2012). Learned integration of visual, vestibular, and motor cues in multiple brain regions computes head direction during visually guided navigation. Hippocampus22, 2219–2237. doi: 10.1002/hipo.22040

  • CrossRef
  • Google Scholar
  • View reference in article

• 34

  FoxM. D.SnyderA. Z.VincentJ. L.CorbettaM.Van EssenD. C.RaichleM. E. (2005). The human brain is intrinsically organized into dynamic, anticorrelated functional networks. Proc. Natl. Acad. Sci. USA102, 9673–9678. doi: 10.1073/pnas.0504136102

  • CrossRef
  • Google Scholar
  • View reference in article

• 35

  FristonK. J. (2011). Functional and effective connectivity: a review. Brain Connect.1, 13–36. doi: 10.1089/brain.2011.0008

  • CrossRef
  • Google Scholar
  • View reference in article

• 36

  GaoP.GanguliS. (2015). On simplicity and complexity in the brave new world of large-scale neuroscience. Curr. Opin. Neurobiol.32, 148–155. doi: 10.1016/j.conb.2015.04.003

  • CrossRef
  • Google Scholar
  • View reference in article

• 37

  GeorgopoulosA. P.KalaskaJ. F.CaminitiR.MasseyJ. T. (1982). On the relations between the direction of two-dimensional arm movements and cell discharge in primate motor cortex. J. Neurosci.2, 1527–1537.

  • Google Scholar
  • View reference in article

• 38

  GersteinG. L.MichalskiA. (1981). “Firing synchrony in a neural group: putative sensory code,” in Neural communication and control. Pergamon, p. 93–102.

  • Google Scholar
  • View reference in article

• 39

  GersteinG. L.PerkelD. H.DayhoffJ. E. (1985). Cooperative firing activity in simultaneously recorded populations of neurons: detection and measurement. J. Neurosci.5, 881–889.

  • Google Scholar
  • View reference in article

• 40

  GoldstoneR. L.MedinD. L.GentnerD. (1991). Relational similarity and the nonindependence of features in similarity judgments. Cogn. Psychol.23, 222–262.

  • Google Scholar
  • View reference in article

• 41

  GorchetchnikovA.GrossbergS. (2007). Space, time and learning in the hippocampus: how fine spatial and temporal scales are expanded into population codes for behavioral control. Neural Netw.20, 182–193. doi: 10.1016/j.neunet.2006.11.007

  • CrossRef
  • Google Scholar
  • View reference in article

• 42

  GrangerC. W. J. (1983). “Generating Mechanisms, Models, and Causality,” in Advances in Econometrics. Cambridge: Cambridge University Press, p. 237–254.

  • Google Scholar
  • View reference in article

• 43

  GreiciusM. D.KrasnowB.ReissA. L.MenonV. (2003). Functional connectivity in the resting brain: a network analysis of the default mode hypothesis. Proc. Natl. Acad. Sci. USA100, 253–258. doi: 10.1073/pnas.0135058100

  • CrossRef
  • Google Scholar
  • View reference in article

• 44

  GrünS.DiesmannM.AertsenA. (2002a). Unitary events in multiple single-neuron spiking activity: I. Detection and significance. Neural Comput.14, 43–80. doi: 10.1162/089976602753284455

  • CrossRef
  • Google Scholar
  • View reference in article

• 45

  GrünS.DiesmannM.AertsenA. (2002b). Unitary events in multiple single-neuron spiking activity: II. Nonstationary data. Neural Comput.14, 81–119. doi: 10.1162/089976602753284464

  • CrossRef
  • Google Scholar
  • View reference in article

• 46

  GurM.SnodderlyD. M. (2007). Direction selectivity in V1 of alert monkeys: evidence for parallel pathways for motion processing. J. Physiol.585, 383–400. doi: 10.1113/jphysiol.2007.143040

  • CrossRef
  • Google Scholar

• 47

  HarrisK. D. (2005). Opinion: neural signatures of cell assembly organization. Nat. Rev. Neurosci.6, 399–407. doi: 10.1038/nrn1669

  • CrossRef
  • Google Scholar

• 48

  HaxbyJ. V. (2012). Multivariate pattern analysis of fMRI: the early beginnings. NeuroImage62, 852–855. doi: 10.1016/j.neuroimage.2012.03.016

  • CrossRef
  • Google Scholar

• 49

  HaxbyJ. V.GobbiniM. I.FureyM. L.IshaiA.SchoutenJ. L.PietriniP. (2001). Distributed and overlapping representations of faces and objects in ventral temporal cortex. Science293, 2425–2430. doi: 10.1126/science.1063736

  • CrossRef
  • Google Scholar

• 50

  HoneyC. J.ThesenT.DonnerT. H.SilbertL. J.CarlsonC. E.DevinskyO.et al. (2012). Slow cortical dynamics and the accumulation of information over long timescales. Neuron76:668. doi: 10.1016/j.neuron.2012.10.024

  • CrossRef
  • Google Scholar

• 51

  HosmanT.HynesJ. B.SaabJ.WilcoxenK. G.BuchbinderB. R.SchmanskyN.et al. (2021). Auditory cues reveal intended movement information in middle frontal gyrus neuronal ensemble activity of a person with tetraplegia. Sci. Rep.11:98. doi: 10.1038/s41598-020-77616-8

  • CrossRef
  • Google Scholar

• 52

  HoughtonC.VictorJ. (2010). “Measuring representational distances--the spike-train metrics approach,” in Visual population codes--toward a common multivariate framework for cell recording and functional imaging. eds. KriegeskorteN.KreimanG. (Cambridge, MA: MIT Press), 391–416.

  • Google Scholar

• 53

  HsiehF.ZhengJ. (2019). Unraveling pattern-based mechanics defining self-organized recurrent behaviors in a complex system: a zebrafish’s calcium brain-wide imaging example. Front. Appl. Math. Stat.5:13. doi: 10.3389/fams.2019.00013

  • CrossRef
  • Google Scholar

• 54

  HumphriesM. D. (2011). Spike-train communities: finding groups of similar spike trains. J. Neurosci.31, 2321–2336. doi: 10.1523/JNEUROSCI.2853-10.2011

  • CrossRef
  • Google Scholar

• 55

  HurwitzC.KudryashovaN.OnkenA.HennigM. H. (2021). Building population models for large-scale neural recordings: opportunities and pitfalls. Curr. Opin. Neurobiol.70, 64–73. doi: 10.1016/j.conb.2021.07.003

  • CrossRef
  • Google Scholar

• 56

  HymanR. (1974). Stimuli and performance: the processing of information and structure. Science186, 730–731.

  • Google Scholar

• 57

  ItoS.HansenM. E.HeilandR.LumsdaineA.LitkeA. M.BeggsJ. M. (2011). Extending transfer entropy improves identification of effective connectivity in a spiking cortical network model. PLoS One6:e27431. doi: 10.1371/journal.pone.0027431

  • CrossRef
  • Google Scholar

• 58

  JefferyK. J.DonnettJ. G.BurgessN.O’KeefeJ. M. (1997). Directional control of hippocampal place fields. Exp. Brain Res.117, 131–142.

  • Google Scholar

• 59

  KanervaP. (2009). Hyperdimensional computing: an introduction to computing in distributed representation with high-dimensional random vectors. Cogn. Comput.1, 139–159. doi: 10.1007/s12559-009-9009-8

  • CrossRef
  • Google Scholar

• 60

  KarpowiczB. M.AliY. H.WimalasenaL. N.SedlerA. R.KeshtkaranM. R.BodkinK.et al. (2025). Stabilizing brain-computer interfaces through alignment of latent dynamics. Nat. Commun.16:4662. doi: 10.1038/s41467-025-59652-y

  • CrossRef
  • Google Scholar

• 61

  KianiR.CuevaC. J.ReppasJ. B.PeixotoD.RyuS. I.NewsomeW. T. (2015). Natural grouping of neural responses reveals spatially segregated clusters in prearcuate cortex. Neuron85, 1359–1373. doi: 10.1016/j.neuron.2015.02.014

  • CrossRef
  • Google Scholar

• 62

  KianiR.EstekyH.MirpourK.TanakaK. (2007). Object category structure in response patterns of neuronal population in monkey inferior temporal cortex. J. Neurophysiol.97, 4296–4309. doi: 10.1152/jn.00024.2007

  • CrossRef
  • Google Scholar

• 63

  KobakD.LindermanG. C. (2021). Initialization is critical for preserving global data structure in both t-SNE and UMAP. Nat. Biotechnol.39, 156–157. doi: 10.1038/s41587-020-00809-z

  • CrossRef
  • Google Scholar

• 64

  KohnA.JasperA. I.SemedoJ. D.GokcenE.MachensC. K.YuB. M. (2020). Principles of corticocortical communication: proposed schemes and design considerations. Trends Neurosci.43, 725–737. doi: 10.1016/j.tins.2020.07.001

  • CrossRef
  • Google Scholar

• 65

  KohnA.SmithM. A. (2016). Utah array extracellular recordings of spontaneous and visually evoked activity from anesthetized macaque primary visual cortex (V1). CRCNS org. doi: 10.6080/K0NC5Z4X

  • CrossRef
  • Google Scholar

• 66

  KriegeskorteN.BandettiniP. (2007). Analyzing for information, not activation, to exploit high-resolution fMRI. NeuroImage38, 649–662. doi: 10.1016/j.neuroimage.2007.02.022

  • CrossRef
  • Google Scholar

• 67

  KriegeskorteN.KievitR. A. (2013). Representational geometry: integrating cognition, computation, and the brain. Trends Cogn. Sci.17, 401–412. doi: 10.1016/j.tics.2013.06.007

  • CrossRef
  • Google Scholar

• 68

  KriegeskorteN.MurM.BandettiniP. (2008). Representational similarity analysis - connecting the branches of systems neuroscience. Front. Syst. Neurosci.2:4. doi: 10.3389/neuro.06.004.2008

  • CrossRef
  • Google Scholar

• 69

  LeeJ. A.Peluffo-OrdóñezD. H.VerleysenM. (2015). Multi-scale similarities in stochastic neighbour embedding: reducing dimensionality while preserving both local and global structure. Neurocomputing169, 246–261. doi: 10.1016/j.neucom.2014.12.095

  • CrossRef
  • Google Scholar

• 70

  LegendreP. (2000). Comparison of permutation methods for the partial correlation and partial mantel tests. J. Stat. Comput. Simul.67, 37–73. doi: 10.1080/00949650008812035

  • CrossRef
  • Google Scholar

• 71

  LehkyS. R.SerenoM. E.SerenoA. B. (2013). Population coding and the labeling problem: extrinsic versus intrinsic representations. Neural Comput.25, 2235–2264. doi: 10.1162/NECO_a_00486

  • CrossRef
  • Google Scholar

• 72

  LillicrapT. P.ScottS. H. (2013). Preference distributions of primary motor cortex neurons reflect control solutions optimized for limb biomechanics. Neuron77, 168–179. doi: 10.1016/j.neuron.2012.10.041

  • CrossRef
  • Google Scholar

• 73

  LivingstoneM. S.FreemanD. C.HubelD. H. (1996). Visual responses in V1 of freely viewing monkeys. Cold Spring Harb. Symp. Quant. Biol.61, 27–37.

  • Google Scholar

• 74

  LivingstoneM. S.HubelD. H. (1984). Specificity of intrinsic connections in primate primary visual cortex. J. Neurosci.4, 2830–2835. doi: 10.1523/jneurosci.04-11-02830.1984

  • CrossRef
  • Google Scholar

• 75

  Lopes-dos-SantosV.RibeiroS.TortA. B. L. (2013). Detecting cell assemblies in large neuronal populations. J. Neurosci. Methods220, 149–166. doi: 10.1016/j.jneumeth.2013.04.010

  • CrossRef
  • Google Scholar

• 76

  MaatenL. V. D.HintonG. (2008). Visualizing data using t-SNE. J. Mach. Learn. Res.9, 2579–2605.

  • Google Scholar

• 77

  MalikW. Q.TruccoloW.BrownE. N.HochbergL. R. (2011). Efficient decoding with steady-state Kalman filter in neural interface systems. IEEE Trans. Neural Syst. Rehabil. Eng.19, 25–34. doi: 10.1109/TNSRE.2010.2092443

  • CrossRef
  • Google Scholar

• 78

  MantelN. (1968). Simultaneous confidence intervals and experimental design with normal correlation. Biometrics24, 434–437.

  • Google Scholar

• 79

  MardiaK. V.ZemrochP. J. (1975). Algorithm AS 86: the Von Mises distribution function. J. R. Stat. Soc.: Ser. C: Appl. Stat.24:268.

  • Google Scholar

• 80

  MartinaM.SchultzJ. H.EhmkeH.MonyerH.JonasP. (1998). Functional and molecular differences between voltage-gated K channels of fast-spiking interneurons and pyramidal neurons of rat hippocampus. J. Neurosci.18, 8111–8125.

  • Google Scholar

• 81

  MaynardE. M.NordhausenC. T.NormannR. A. (1997). The Utah intracortical electrode Array: a recording structure for potential brain-computer interfaces. Electroencephalogr. Clin. Neurophysiol.102, 228–239.

  • Google Scholar

• 82

  MazurekM.KagerM.Van HooserS. D. (2014). Robust quantification of orientation selectivity and direction selectivity. Front. Neural Circuits8:92.

  • Google Scholar

• 83

  MeyerA. F.WilliamsonR. S.LindenJ. F.SahaniM. (2016). Models of neuronal stimulus-response functions: elaboration, estimation, and evaluation. Front. Syst. Neurosci.10:109. doi: 10.3389/fnsys.2016.00109

  • CrossRef
  • Google Scholar

• 84

  OlshausenB. A.FieldD. J. (2005). How close are we to understanding v1?Neural Comput.17, 1665–1699. doi: 10.1162/0899766054026639

  • CrossRef
  • Google Scholar

• 85

  PaivaA. R.ParkI.PríncipeJ. C. (2010). A comparison of binless spike train measures. Neural Comput. Appl.19, 405–419.

  • Google Scholar

• 86

  PaninskiL.CunninghamJ. P. (2018). Neural data science: accelerating the experiment-analysis-theory cycle in large-scale neuroscience. Curr. Opin. Neurobiol.50, 232–241. doi: 10.1016/j.conb.2018.04.007

  • CrossRef
  • Google Scholar

• 87

  PastalkovaE.ItskovV.AmarasinghamA.BuzsákiG. (2008). Internally generated cell assembly sequences in the rat hippocampus. Science321, 1322–1327. doi: 10.1126/science.1159775

  • CrossRef
  • Google Scholar

• 88

  PastalkovaE.WangY.MizusekiK.BuzsákiG. (2015). Simultaneous extracellular recordings from left and right hippocampal areas CA1 and right entorhinal cortex from a rat performing a left/right alternation task and other behaviors. Available at: CRCNS.org. https://crcns.org/data-sets/hc/hc-5

  • Google Scholar

• 89

  PaulkA. C.KfirY.KhannaA. R.MustrophM. L.TrautmannE. M.SoperD. J.et al. (2022). Large-scale neural recordings with single neuron resolution using Neuropixels probes in human cortex. Nat. Neurosci.25, 252–263. doi: 10.1038/s41593-021-00997-0

  • CrossRef
  • Google Scholar

• 90

  PellegrinoA.SteinH.Cayco-GajicN. A. (2024). Dimensionality reduction beyond neural subspaces with slice tensor component analysis. Nat. Neurosci.27, 1199–1210. doi: 10.1038/s41593-024-01626-2

  • CrossRef
  • Google Scholar

• 91

  PospisilD. A.BairW. (2021). Accounting for biases in the estimation of neuronal signal correlation. J. Neurosci.41, 5638–5651. doi: 10.1523/JNEUROSCI.2775-20.2021

  • CrossRef
  • Google Scholar

• 92

  RaichleM. E.MacLeodA. M.SnyderA. Z.PowersW. J.GusnardD. A.ShulmanG. L. (2001). A default mode of brain function. Proc. Natl. Acad. Sci. USA98, 676–682. doi: 10.1073/pnas.98.2.676

  • CrossRef
  • Google Scholar

• 93

  RaoN. G.DonoghueJ. P. (2014). Cue to action processing in motor cortex populations. J. Neurophysiol.111, 441–453. doi: 10.1152/jn.00274.2013

  • CrossRef
  • Google Scholar

• 94

  RichardG. R.TitizA.TylerA.HolmesG. L.ScottR. C.Lenck-SantiniP.-P. (2013). Speed modulation of hippocampal theta frequency correlates with spatial memory performance. Hippocampus23, 1269–1279. doi: 10.1002/hipo.22164

  • CrossRef
  • Google Scholar

• 95

  RoudiY.DunnB.HertzJ. (2015). Multi-neuronal activity and functional connectivity in cell assemblies. Curr. Opin. Neurobiol.32, 38–44. doi: 10.1016/j.conb.2014.10.011

  • CrossRef
  • Google Scholar

• 96

  RousseeuwP. J. (1987). Silhouettes: a graphical aid to the interpretation and validation of cluster analysis. J. Comput. Appl. Math.20, 53–65.

  • Google Scholar

• 97

  SaxenaS.CunninghamJ. P. (2019). Towards the neural population doctrine. Curr. Opin. Neurobiol.55, 103–111. doi: 10.1016/j.conb.2019.02.002

  • CrossRef
  • Google Scholar

• 98

  SchemperM. (1984). Analyses of associations with censored data by generalized Mantel and Breslow tests and generalized Kendall correlation coefficients. Biom. J.26, 309–318.

  • Google Scholar

• 99

  SchneiderS.LeeJ. H.MathisM. W. (2022). Learnable latent embeddings for joint behavioral and neural analysis. arXiv. Available at: http://arxiv.org/abs/2204.00673

  • Google Scholar

• 100

  ShannonC. E. (1997). The mathematical theory of communication. MD Comput.14, 306–317.

  • Google Scholar

• 101

  ShepardR. N. (1964). Attention and the metric structure of the stimulus space. J. Math. Psychol.1, 54–87.

  • Google Scholar

• 102

  ShepardR. N. (1987). Toward a universal law of generalization for psychological science. Science237, 1317–1323.

  • Google Scholar

• 103

  ShepardS.ChipmanR. (1970). Second-order isomorphism of internal representations: shapes of states. Cogn. Psychol.1, 1–17.

  • Google Scholar

• 104

  SingerW.GrayC. M. (1995). Neuronal synchrony: a versatile code for the definition of relations?Annu. Rev. Neurosci.18, 49–65.

  • Google Scholar

• 105

  SkaggsWEMcnaughtonBLGothardKMMarkusEJ (1992) “An information-theoretic approach to deciphering the hippocampal code”. Proceedings of the 5th International Conference on Neural Information Processing Systems

  • Google Scholar

• 106

  SmithM. A.KohnA. (2008). Spatial and temporal scales of neuronal correlation in primary visual cortex. J. Neurosci.28, 12591–12603. doi: 10.1523/JNEUROSCI.2929-08.2008

  • CrossRef
  • Google Scholar

• 107

  SoudryD.KeshriS.StinsonP.OhM.-H.IyengarG.PaninskiL. (2013). A shotgun sampling solution for the common input problem in neural connectivity inference. arXiv. Available at: http://arxiv.org/abs/1309.3724

  • Google Scholar

• 108

  SpornsO. (2013). Structure and function of complex brain networks. Dialogues Clin. Neurosci.15, 247–262. doi: 10.31887/dcns.2013.15.3/osporns

  • CrossRef
  • Google Scholar

• 109

  SteinmetzN. A.AydinC.LebedevaA.OkunM.PachitariuM.BauzaM.et al. (2021). Neuropixels 2.0: a miniaturized high-density probe for stable, long-term brain recordings. Science372:eabf4588. doi: 10.1126/science.abf4588

  • CrossRef
  • Google Scholar

• 110

  SteinmetzN. A.KochC.HarrisK. D.CarandiniM. (2018). Challenges and opportunities for large-scale electrophysiology with Neuropixels probes. Curr. Opin. Neurobiol.50, 92–100. doi: 10.1016/j.conb.2018.01.009

  • CrossRef
  • Google Scholar

• 111

  StevensonI. H.KordingK. P. (2011). How advances in neural recording affect data analysis. Nat. Neurosci.14, 139–142. doi: 10.1038/nn.2731

  • CrossRef
  • Google Scholar

• 112

  StuartL.WalterM.BorisyukR. (2002). Visualisation of synchronous firing in multidimensional spike trains. Bio Systems. 67, 265–79. doi: 10.1016/S0303-2647(02)000849

  • CrossRef
  • Google Scholar

• 113

  TeetersJ. L.SommerF. T. (2009). CRCNS.ORG: a repository of high-quality data sets and tools for computational neuroscience. BMC Neurosci. 10:S6. doi: 10.1186/1471-2202-10-S1-S6

  • CrossRef
  • Google Scholar

• 114

  Thomas YeoB. T.KrienenF. M.SepulcreJ.SabuncuM. R.LashkariD.HollinsheadM.et al. (2011). The organization of the human cerebral cortex estimated by intrinsic functional connectivity. J. Neurophysiol.106, 1125–1165. doi: 10.1152/jn.00338.2011

  • CrossRef
  • Google Scholar

• 115

  ThurstoneL. L. (2017). “A law of comparative judgment,” in Scaling. Milton Park, Abingdon, Oxfordshire, England: Routledge, p. 1–92.

  • Google Scholar

• 116

  TruccoloW.EdenU. T.FellowsM. R.DonoghueJ. P.BrownE. N. (2005). A point process framework for relating neural spiking activity to spiking history, neural ensemble, and extrinsic covariate effects. J. Neurophysiol.93, 1074–1089. doi: 10.1152/jn.00697.2004

  • CrossRef
  • Google Scholar

• 117

  TruccoloW.FriehsG. M.DonoghueJ. P.HochbergL. R. (2008). Primary motor cortex tuning to intended movement kinematics in humans with tetraplegia. J. Neurosci.28, 1163–1178. doi: 10.1523/JNEUROSCI.4415-07.2008

  • CrossRef
  • Google Scholar

• 118

  TruccoloW.HochbergL. R.DonoghueJ. P. (2010). Collective dynamics in human and monkey sensorimotor cortex: predicting single neuron spikes. Nat. Neurosci.13, 105–111. doi: 10.1038/nn.2455

  • CrossRef
  • Google Scholar

• 119

  Ts’oD. Y.GilbertC. D.WieselT. N. (1986). Relationships between horizontal interactions and functional architecture in cat striate cortex as revealed by cross-correlation analysis. J. Neurosci.6, 1160–1170.

  • Google Scholar

• 120

  TverskyA. (1977). Features of similarity. Psychol. Rev.84, 327–352. doi: 10.1037/0033-295X.84.4.327

  • CrossRef
  • Google Scholar

• 121

  UraiA. E.DoironB.LeiferA. M.ChurchlandA. K. (2022). Large-scale neural recordings call for new insights to link brain and behavior. Nat. Neurosci.25, 11–19. doi: 10.1038/s41593-021-00980-9

  • CrossRef
  • Google Scholar

• 122

  UrsinoM.RicciG.MagossoE. (2020). Transfer entropy as a measure of brain connectivity: a critical analysis with the help of neural mass models. Front. Comput. Neurosci.14:45. doi: 10.3389/fncom.2020.00045

  • CrossRef
  • Google Scholar

• 123

  Vargas-IrwinC. E.BrandmanD. M.ZimmermannJ. B.DonoghueJ. P.BlackM. J. (2015a). Spike train SIMilarity space (SSIMS): a framework for single neuron and ensemble data analysis. Neural Comput.27, 1–31. doi: 10.1162/NECO_a_00684

  • CrossRef
  • Google Scholar

• 124

  Vargas-IrwinC. E.FranquemontL.BlackM. J.DonoghueJ. P. (2015b). Linking objects to actions: encoding of target object and grasping strategy in primate ventral premotor cortex. J. Neurosci.35, 10888–10897. doi: 10.1523/JNEUROSCI.1574-15.2015

  • CrossRef
  • Google Scholar

• 125

  VictorJ. D. (2005). Spike Train Metrics. Curr. Opin. Neurobiol.15, 585–92.

  • Google Scholar

• 126

  VictorJ. D.PurpuraK. P. (1996). Nature and precision of temporal coding in visual cortex: a metric-space analysis. J. Neurophysiol.76, 1310–1326. doi: 10.1152/jn.1996.76.2.1310

  • CrossRef
  • Google Scholar

• 127

  VictorJ. D.PurpuraK. P. (1997). Metric-space analysis of spike trains: theory, algorithms, and application. Netw. Comput. Neural Syst.8, 127–64.

  • Google Scholar

• 128

  VogtN. (2019). Massively parallel intracellular recordings. Nat. Methods16:1079. doi: 10.1038/s41592-019-0644-z

  • CrossRef
  • Google Scholar

• 129

  WangY.RomaniS.LustigB.LeonardoA.PastalkovaE. (2015). Theta sequences are essential for internally generated hippocampal firing fields. Nat. Neurosci.18, 282–288. doi: 10.1038/nn.3904

  • CrossRef
  • Google Scholar

• 130

  WrightW. E. (1977). Gravitational clustering. Pattern Recogn.9, 151–66.

  • Google Scholar

• 131

  YangG. R.JoglekarM. R.SongH. F.NewsomeW. T.WangX.-J. (2019). Task representations in neural networks trained to perform many cognitive tasks. Nat. Neurosci.22, 297–306. doi: 10.1038/s41593-018-0310-2

  • CrossRef
  • Google Scholar

• 132

  YusteR. (2015). From the neuron doctrine to neural networks. Nat. Rev. Neurosci.16, 487–497. doi: 10.1038/nrn3962

  • CrossRef
  • Google Scholar

• 133

  ZadorA.RichardsB.ÖlveczkyB.EscolaS.BengioY.BoahenK.et al. (2022). Toward next-generation artificial intelligence: catalyzing the NeuroAI revolution. arXiv. Available at: http://arxiv.org/abs/2210.08340

  • Google Scholar

• 134

  ZaghaE.GeX.McCormickD. A. (2015). Competing neural ensembles in motor cortex gate goal-directed motor output. Neuron88, 565–577. doi: 10.1016/j.neuron.2015.09.044

  • CrossRef
  • Google Scholar

• 135

  ZenkerF.GärdenforsP. (eds.). (2015). Applications of conceptual spaces: The case for geometric knowledge representation. Cham, Switzerland: Springer.

  • Google Scholar