Grid cells have attracted broad attention because of their highly symmetric hexagonal firing patterns. Recently, research has shifted its focus from the global symmetry of grid cell activity to local distortions both in space and time, such as drifts in orientation, local defects of the hexagonal symmetry, and the decay and reappearance of grid patterns after changes in lighting condition. Here, we introduce a method that allows to visualize and quantify such local distortions, by assigning both a local grid score and a local orientation to each individual spike of a neuronal recording. The score is inspired by a standard measure from crystallography, which has been introduced to quantify local order in crystals. By averaging over spikes recorded within arbitrary regions or time periods, we can quantify local variations in symmetry and orientation of firing patterns in both space and time.
Neurons in the entorhinal cortex are often classified based on their spatial firing patterns. Grid cells, in particular, have attracted a lot of attention, because of their highly symmetric hexagonal activity [1,2]. To decide whether a given cell should be classified as a grid cell or a non-grid cell [3,3–5], researchers need a grid score: a number that quantifies the degree of hexagonal spatial symmetry in the firing pattern of the cell. Since the first reports of grid cells, the community has developed and refined a standard grid score [6–9], which relies on the following procedure. Spike locations are transformed into a rate map. The peaks in the spatial autocorrelogram of this rate map are used to obtain measures of the grid spacing and the grid orientation. The autocorrelogram of the rate map is then cropped, rotated and correlated with its unrotated copy. A grid score is finally obtained from the resulting correlation-vs-angle function at selected angles (for details, see S1 Appendix). This procedure assumes a consistent hexagonal pattern throughout the arena and results in a single grid score for the entire firing pattern. While this grid score has many benefits, it is inherently global and ignores local properties of firing patterns. Recent studies have shown, however, that the firing patterns of grid cells show interesting local effects, such as drifts in orientation [10] or local defects of the hexagonal symmetry [11,12]. These results call for new methods for quantifying the local symmetry and the local orientation of grids. Here, we introduce a method that assigns both a grid score and an orientation to each individual spike of a neuronal recording. This spike-based score is inspired by a standard measure from crystallography [13], which has been introduced to quantify local order in crystals. By averaging over spikes recorded within arbitrary regions or time periods, we can quantify local variations in symmetry and orientation of firing patterns in both space and time. A high spatial and temporal resolution can be achieved by averaging over multiple cells. To illustrate the potential of our method, we apply it to both artificially generated data and neuronal recordings from grid cells. In the following, we first introduce the method in detail. Comparing it to the conventional correlogram-based grid score [7], we show that both methods perform similarly when quantifying and classifying grid properties globally. We then illustrate the benefits of our spike-based score to analyze and visualize local properties of spike patterns in space and time.
The angular symmetry of an arrangement of N neighbors around a given atom in a crystal can be quantified by a complex number [13]: ψ (M) ≔ 1 N ∑ l = 1 N exp (i M ϕ l), where i is the imaginary unit and ϕl is the angle between an arbitrary reference vector and the vector pointing from the atom in question to its l-th neighbor [13]. In the following, we will use the horizontal axis from left to right as the reference vector (Fig 1a). M is a natural number that characterizes the symmetry of interest. For instance, M = 6 corresponds to a hexagonal arrangement of atoms and M = 4 to a quadratic arrangement. The computation of this score can be visualized as an addition of vectors in the complex plane (Fig 1a). The absolute value |ψ (M) | ranges from 0 to 1, and larger values indicate a more pronounced M-fold symmetry. If the N neighboring atoms lie on the corners of a regular M-polygon, the angle 1Marg (ψ (M) ) is the orientation of this polygon (Fig 1a). The function arg (x) maps to values between −π and π radians, so the resulting orientation is between −π/M and π/M. For hexagonal symmetry, this corresponds to values between −30 and 30 degrees. The ψ measure has been introduced to quantify symmetry in crystalline structures where neighborhood is clearly defined. For example, the nearest neighbors of a reference atom are all atoms within a distance of one lattice constant. For a given spike in a grid pattern, however, the nearest neighbors are most likely other spikes from the same grid field. To assess whether a given spike is part of a symmetric grid, we are hence interested not in neighboring spikes, but in spikes in neighboring grid fields, i. e. , spikes that are roughly one grid spacing away from the spike in question. To identify these relevant spikes, we define a neighborhood shell—a two dimensional annulus (Fig 1b). We then center this neighborhood shell at the location of the spike in question and consider all spikes within the shell as its neighbors. Assuming that the grid spacing does not change drastically across the environment, we use the same shell size for all spikes. To identify a suitable size of the shell, we use the histogram of pairwise distances between all spike locations. For periodic firing patterns, this histogram has pronounced peaks (Fig 1b). The first peak characterizes the size of the grid fields. The second peak corresponds to the grid spacing ℓ, and further peaks arise at multiples of the grid spacing. We therefore set the inner radius of the neighborhood shell to 5 6 ℓ and the outer radius to 7 6 ℓ (Fig 1b). To reduce the noise arising from finite numbers of spikes, we determine the peaks from a smoothed histogram, obtained by applying a Gaussian filter with a standard deviation of 1% of the maximal pairwise distance. The suggested shell width of 1 3 ℓ is motivated by the typical size of grid fields, but the suggested grid score is robust to changes of the shell width. Scaling the shell width with the grid spacing is commensurate with the experimental finding that the ratio of grid spacing and grid field size is constant [15]. The detection of the correct neighborhood shell as the second peak of the histogram of pairwise distances is a crucial step in the computation of the grid score. For noisy grid patterns, the first peak of the histogram—which indicates the size of the grid fields—is sometimes hard to detect, leading to a potential misclassification of the second peak as the first, and consequently a wrong neighborhood shell. In data sets where this is the case, we use the following remedy: We define a rough cutoff distance that is larger than the grid field size and smaller than the grid spacing, and use the first peak above the cutoff to determine the shell size. For the analyzed grid patterns, taking 15% of the side length of the arena for this cutoff worked reliably. In the dataset with very noisy and elliptic grids, peak detection can fail for many cells. In this case, we exploit the fact that nearby grid cells have a similar grid spacing [2] and determine the neighborhood shell from less noisy grid cells of the same module. Whenever we use either of these remedies, we highlight it in the main text. A large | ψ k (6) | value indicates that the neighbors of the k-th spike are arranged on the corners of a hexagon, but it does not require that all corners of the hexagon are occupied. If, for example, the firing fields of a cell are arranged on a line, this leads to high |ψ (6) | values although the structure does not resemble a grid (Fig 1c). This can lead to a wrongful classification of a cell as a grid cell, i. e. , to false positives. To avoid high local grid scores without an actual 6-fold symmetry, we compute | ψ k (M) | for a range of symmetries (M ∈ {2, …, 7}) and only assign a non-zero grid score to spike k if the hexagonal symmetry is the strongest: ψ^k≔{| ψk (6) |, if| ψk (6) |>| ψk (M) |∀M∈{2,3, 4,5, 7 }0, otherwise. We refer to ψ ^ k as the grid score of spike k. Comparing to symmetries with M > 7 did not have an effect on our results, so we restricted M to a maximum of 7. A global grid score Ψ and the global grid orientation Θ for a given cell is calculated by averaging over the respective properties of the individual spikes. A summary of the algorithm is shown in Fig 2. In the following, we compare the suggested local grid score to a conventional correlation-based grid score. Details on the calculation of this score, on the shuffling procedure for grid cell classification, and on the generation of artificial spike data can be found in S1 Appendix.
Grid scores are typically used for two purposes: Quantifying the hexagonality, i. e. , the gridness, of the spike map and classifying a cell as a grid cell. Given that the patterns of grid cells are rarely perfectly hexagonal, a grid score should classify a cell as a grid cell even if the hexagonal pattern is weakly distorted. For example, a grid score needs to be robust against small shifts of grid field locations or elliptic distortions like shearing. At the same time, a grid score should quantify such distortions by returning lower values for stronger distortions. So far we assumed that gridness is an arena-wide property of a firing pattern. We did not consider local defects. But recent experiments have shown that grid patterns are often locally distorted. Such distortions have been observed in both asymmetric [11,12] and quadratic environments [10,20], and are presumably caused by boundary effects. Moreover, it has been shown that grid patterns recorded in two arenas separated by a wall do not form a coherent grid, but merge into a more coherent pattern when the wall is removed [19]. We now show that our spike-based grid score is suitable to quantify and visualize such local properties. We first illustrate this for generated spike data and then for experimental recordings. To generate an example of a grid pattern with local defects, we first generated a perfect grid and then added random vectors to the field locations (cf. Fig 3a, see S1 Appendix for details) only in the eastern half of the arena (Fig 6a and 6b). To analyze gridness locally, we divide the arena into partitions and compute the Ψ score within each partition, i. e. , the mean over all ψ ^ k scores for all spikes k within the given partition (Fig 6c). Averaging over multiple cells highlights that the gridness decays along the west-east axis but stays stable along the south-north axis, and this trend is already visible for individual cells (Fig 6c and 6d). Note that due to the strong noise on the east side of the arena, we detect the neighborhood shell from the histogram of pairwise spike distances, using the first peak above a cutoff (see Methods). We now apply the same analysis to recent experimental data. Wernle et al. recorded grid cells from rats that explored the northern and southern half of a 2m × 2m box [19]. The two halves were initially separated by a wall. The cells fired in grid patterns independently in both halves of the arena (Fig 6e). Then the authors removed the wall and let the rat explore the entire arena. At the southern and northern boundaries, the grid fields typically stayed at the same locations as before wall removal. Close to the former location of the wall, the grid fields moved—and sometimes even merged—to form a more coherent grid pattern (Fig 6f). Quantifying the local increase of coherence requires a local grid score. We therefore use our spike-based grid score to show that the rearrangement of grid fields after wall removal indeed leads to a local increase of gridness. Given that the firing patterns of the grid cells are noisy and strongly elliptic, it is difficult to determine the peaks in the histograms of pairwise spike distances (Fig 6e and 6f). We thus set the central distance of the neighborhood shell manually to ℓ = 65cm—a value that we obtained by visual inspection of the histograms of pairwise spike distances of all 11 recorded cells, after wall removal. We divide the arena into three equal sized partitions. For south-north partitioning (Fig 6g), the grid score is typically lowest in the partition that contains the separating wall (Fig 6h). After wall removal, the grid score increases in all partitions, but the increase is most prominent in the partition that contained the separating wall (Fig 6h). In contrast, a west-east partitioning does not lead to a single partition with predominantly low grid scores. Moreover, the increase after wall removal is similar for all west-east partitions (Fig 6i and 6j). Our analysis thus quantifies the qualitative observation that two initially independent grids locally rearrange to form a more coherent grid. For less elliptic grid patterns, we assume our quantitative results to be even more pronounced. In summary, a spike-based grid score visually highlights local grid defects. Averaging spike scores in arbitrary partitions and potentially over many cells enables a systematic analysis of grid defects on a spatial scale smaller than the grid spacing. Early work on grid cells has shown that the firing patterns of anatomically nearby cells typically have a similar orientation [2], which tends to align with a cardinal axis of the environment [10,11]. While these findings suggest a constant grid orientation, experiments from the same group have shown that the orientation can drift within the arena or exhibit local distortions [10,20]. To quantify these experiments, a local measure of grid orientation is needed. In the following, we study local grid orientation using the angle θk of spike k (see the algorithm in Fig 2 for the definition) on generated (see S1 Appendix) and recorded firing patterns that were kindly provided by Stensola et al. [10]. When each spike is color-coded according to its grid orientation θk, the spikes of firing fields have the same color for a perfect hexagonal grid—except the spikes at the edge of a grid field (Fig 7a). In contrast, the color-coding shows a sudden switch if the orientation is changed abruptly (Fig 7b) or a color gradient if the grid orientation drifts continuously across the arena, both in generated data (Fig 7c) and in recordings (Fig 7d). If the drift in orientation is along the south-north axis of the box, the gradient follows this axis (Fig 7c and 7d). This can be highlighted by partitioning the arena and calculating the mean value of the orientation θk for all spikes k within each partition (Fig 7c and 7d). To highlight an orientation drift recorded in multiple cells, we average over the the orientation θk for spikes from several cells. For the recordings of Stensola et al. [10], this shows that the orientation increases monotonically from south to north, but varies less systematically from west to east, for a set of 5 cells recorded in the same rat (Fig 7e and 7f). Note that changes in grid orientation distort the correlogram, which can lead to low ρ scores (Fig 7b and 7c), even if the grid-like structure is clearly apparent. In contrast, the Ψ score is robust to distortions in orientation because of its locality. In summary, a spike-based orientation measure highlights local distortions in grid orientation. Averaging spike-based orientations, using arbitrary spatial partitionings, allows to further classify the distortion. Whether or not the spikes of a cell form a hexagonal spatial pattern can vary in time. For example, the activity of grid cells loses its periodicity when the light is turned off [21,22]. A decay of grid-like firing patterns can also be triggered by internal modifications, like the deactivation of hippocampal drive to the entorhinal cortex [23]. In the following, we show how the spike-based grid score can be used to study changes in grid cell activity with high temporal resolution. The standard procedure for computing a grid score in a given time interval is to take the spike locations from all spikes that occurred in this interval and compute the grid score for the resulting pattern. While this can be done with both the correlogram-based score and the spike-based score, it introduces constraints on the temporal resolution, because the temporal window must contain a sufficient coverage of the arena to obtain a meaningful score. This problem can be circumvented with a spike-based grid score, because for a given spike map, each spike k has a grid score ψ ^ k and a time tk. In other words, the grid score ψ ^ k at time tk quantifies the contribution of spike k to the grid pattern formed by all spikes. Just like the global grid score Ψ is calculated by averaging over the individual ψ ^ k values of all spikes, a temporally local grid score can be obtained by averaging the ψ ^ k values of the spikes in an arbitrary time interval. The resulting average quantifies how much the spikes in this interval contribute to the grid in the complete recording. To demonstrate that the spike-based grid score can characterize temporal changes in grid patterns, we consider generated data with a sudden transition from spatially homogeneous to grid cell activity (Fig 8a). Plotting the individual grid scores as a function of spike time highlights the transition to the grid pattern, although the grid scores ψ ^ k of the individual spikes are computed using the locations of all spikes (Fig 8a and 8b). The transition is emphasized by temporal smoothing (Fig 8c). Requirements for the smoothing filter are discussed later. Note that such an analysis first of all requires that a grid pattern can be detected in the full recording, so that a meaningful neighborhood shell can be identified. To facilitate the detection of the neighborhood shell in the presence of noise, we detect the first peak above a cutoff in the histogram of pairwise spike distances (see Methods). We now use this method to analyze spike patterns from cells in mice that forage in an arena in alternating trials of darkness and light [21]. For the recorded grid cells (Fig 8d), the filtered time evolution of grid scores ψ ^ k is strongly correlated with darkness and light—grid scores tend to be higher in light trials (Fig 8e, see S1 Appendix for details on the data preparation). The Pearson correlation coefficient between lighting condition and the smoothed grid scores (r = 0. 57 and r = 0. 53 for single cells and r = 0. 69 for the average over smoothed grid scores of two cells) indicate that only few simultaneously recorded cells are necessary to reach confident conclusions about temporal variations in the quality of grid patterns or contextual influences on grid cell activity. But how quickly does the grid pattern decay when the light is turned off? How rapidly does it reappear when the light is turned on again? Previously, these questions were answered with a temporal resolution of 10 seconds using map similarity [21] (see below for a discussion of map similarity). We now show that we can reach the same temporal resolution using the spike-based grid score. To this end, we create a reference spike map for each recorded grid cell in the data [21], which comprises the locations of all spikes that were fired during light trials. Note that these subsets of spikes show a much clearer grid pattern than that of spikes during darkness (Fig 8f). We then quantify how well the spikes of a given cell fit into the grid pattern of the cell’s reference spike map, by adding every spike k individually to the reference spike map and computing its ψ ^ k score as usual. We then average the ψ ^ k score of spikes that occur at similar times during light/dark trials. To this end, we divide each 2 minute trial into 12 × 10 second blocks. This division is done separately for light and dark trials. We then average the grid scores ψ ^ k over all spikes and over 73 recorded grid cells within each 10 second block. Plotting this average as a function of time shows that grids decay within 30 seconds after the light is turned off (Fig 8g). Grids need a similar time to reappear after the light is turned on (Fig 8g). In summary, the spike-based grid score readily contains temporal information which can be used to study distortions of the grid pattern in time. The time evolution of the grid score can be correlated with any other time varying variable to unravel its influence on the grid pattern, e. g. , light conditions, running speed or the concentration of neuromodulators.
We have shown spatially local grid scores in sub-partitions of the arena and temporally local grid scores in time intervals. The resulting values can be misleading if the spatial subinterval is too small or if the temporal subinterval is too short. For example, if we study the local grid score of a single cell in a partition that does not contain a grid field of that cell, the mean over all spike scores in that partition will be very low, suggesting a local defect. Such sampling biases could be mitigated by making the partition larger (or the interval longer, for temporal analyses), at the cost of resolution. For a single cell, the resolution is thus bounded by the grid spacing. However, a reduction of sampling biases without a loss of resolution can be obtained by averaging over multiple cells. It would be beneficial if these cells have different spatial phases, because this ensures the presence of a grid field in each partition and thus reduces sampling biases. Averaging over many cells thus enables a higher spatial resolution than the grid spacing and a higher temporal resolution than the time the rat needs to sample a sufficient number of grid fields in a single cell. The correlogram-based ρ score can be made more local by computing it in spatial subregions within the arena. A window can also be used to compute the sliding average of the ρ score in the entire arena [20], resulting in a spatial map of local grid scores. Selecting the size of the window is subject to a trade-off: while a small window is necessary for a high spatial resolution, the ρ score will fail to detect a grid if the window is too small—e. g. , smaller then the grid spacing. In contrast, the spike-based score can be computed in arbitrary spatial partitionings (Fig 7) by averaging over the spikes in each partition. The side length of a partition can even be smaller than the grid spacing, in particular if multiple cells are used to compute the score in each partition (see above). The flexibility to use partitionings of an arbitrary shape, e. g. , trapezoidal or triangular, can be important to study grid cells in environments with complex geometries [11,12]. We have shown that a spike-based grid score enables a straightforward analysis of the temporal stability of grid patterns and illustrated this on recorded data from mice foraging through an arena in light and darkness [21]. In our analysis of the decay and reappearance of grid cell activity (Fig 8), we reproduced the results shown in [21], where the authors used map similarity for their analysis. In the map similarity measure, a rate map is correlated with a reference map. As a reference map, the authors used recordings during times of light. In each 10 second block, they collected all spike locations from that block in all trials (separating dark and light trials) and used these spike locations to create a rate map. If the correlation of the resulting rate map with the reference map is high, the pattern is similar, indicating that the grid has not decayed yet, and vice versa. For grid patterns without defects in field locations and with similar firing rates for different grid fields, map similarity and our measure should lead to similar results. However, our measure focuses on how much a spike is part of the grid of the reference map. In contrast, map similarity only quantifies how much spikes are part of the rate pattern of the reference map, irrespective of whether it forms a grid. If for example, the reference map has a firing field that does not belong to the hexagonal pattern, spikes that are located in this firing field will lead to low values in our measure. In contrast, map similarity assigns high values to activity within any firing field of the reference map, be it a grid field or a defect. We have shown that the spike-based Ψ score is more robust to elliptic deformations of grid patterns than the correlogram-based ρ score (Fig 4c). This robustness can be crucial for the detection of grid cells, because many grid patterns are of elliptic shape. Quantifying the local improvement of a grid pattern in the scenario where two grids coalesce in a contiguous environment (Fig 6) was only possible because grid patterns could be detected, even though they were strongly elliptic. A correlogram-based grid score can also be used with elliptic grid patterns, but requires an explicit compensation of the ellipticity. To this end, the ellipticity is detected and removed from the correlogram before the grid score is computed, by compressing the correlogram along the principal axis of the elliptic distortion [24]. In principle, we could apply such a transformation to all spike locations and then compute the spike-based grid scores. We have shown, however, that such an additional step is not necessary as long as the ellipticity is not too strong. If a firing pattern is grid-like, but with different grid spacings in different parts of the arena, we cannot compute the neighborhood shell from all pairwise distances. Instead, the neighborhood would need to be determined individually for each spike, using the histogram of pairwise distances only between a given spike and all other spikes. Since we are not aware of local distortions to the grid spacing, we did not consider this. For grid fields with six neighboring fields, spikes in the center of the field receive the highest grid score and grid scores decay isotropically in all direction. In contrast, in grid fields at the boundaries, spikes with high ψ ^ k form a curved stripe (S1 Fig). The reason is the following: Being in the boundary grid field, moving towards the central neighbors perturbs the 6-fold symmetry stronger than moving tangentially to it (S1 Fig). Such a directional dependence is not present in grid fields that have neighboring grid fields on all sides. We focused on hexagonal patterns, i. e. , grids of 6-fold symmetry. Experiments suggest that there are also band-like firing patterns, i. e. , 2-fold symmetry [3]—though their existence has been questioned [25,26]. Moreover, boundary effects in arenas of complex shape could cause grid deformations that result in other symmetries. The suggested spike-based grid score can be readily applied to any M-fold symmetry. For example, to study quadratic symmetry, we would replace M = 6 by M = 4 and compare to the remaining symmetries, now including 6 (see the algorithm in S1 Appendix). Depending on the symmetry, the radius of the neighborhood shell needs to be chosen differently. For example, for quadratic symmetry, results are improved if the neighborhood shell is set to 2/3ℓ instead of ℓ, to only incorporate the four nearest neighboring fields (Fig 9d). For band-like symmetry, the shell radius should be ℓ/2, so that the next band is not included in the shell (Fig 9e). As the experiments on spatially tuned cells become more complex, new properties will be unveiled. A spike-based grid score could be used as a single unifying method to classify grid cells and to analyze their local and temporal characteristics. The presented method could further be useful for the analysis of grid patterns in other modalities [27,28] or in other contexts, e. g. , in the anatomical arrangement of grid cells [29] or retinal ganglion cells [30].
Grid cells are neurons in mammals whose activity depends on the location of the animal in a striking way: grid cells fire spikes at multiple locations that form a symmetric lattice of repeating hexagons. Recent experiments have shown that the symmetry of these patterns is not as stable as initially thought. For example, patterns show local distortions or drifts in orientation. Here, we propose a method to visualize and quantify such local distortions directly from the recorded spike locations. Our method is inspired by an approach used in crystallography and we show that it reliably detects distortions in grid patterns. Moreover, we demonstrate that it can be used to study the stability of grid patterns not only in space, but also in time. Our method enables researchers to analyze the spatial symmetry of neuronal activity in more detail and could thus contribute to the understanding of spatial representations in mammals.
ecology and environmental sciences
cell cycle and cell division
applied mathematics
cell processes
condensed matter physics
social sciences
light
neuroscience
electromagnetic radiation
simulation and modeling
algorithms
crystals
mathematics
materials science
vision
crystallography
research and analysis methods
solid state physics
animal cells
soil science
soil perturbation
physics
psychology
cellular neuroscience
cell biology
neurons
biology and life sciences
cellular types
physical sciences
materials
sensory perception
A local measure of symmetry and orientation for individual spikes of grid cells