Endogenous Sequential Cortical Activity Evoked by Visual Stimuli

Although the functional properties of individual neurons in primary visual cortex have been studied intensely, little is known about how neuronal groups could encode changing visual stimuli using temporal activity patterns. To explore this, we used in vivo two-photon calcium imaging to record the activity of neuronal populations in primary visual cortex of awake mice in the presence and absence of visual stimulation. Multidimensional analysis of the network activity allowed us to identify neuronal ensembles defined as groups of cells firing in synchrony. These synchronous groups of neurons were themselves activated in sequential temporal patterns, which repeated at much higher proportions than chance and were triggered by specific visual stimuli such as natural visual scenes. Interestingly, sequential patterns were also present in recordings of spontaneous activity without any sensory stimulation and were accompanied by precise firing sequences at the single-cell level. Moreover, intrinsic dynamics could be used to predict the occurrence of future neuronal ensembles. Our data demonstrate that visual stimuli recruit similar sequential patterns to the ones observed spontaneously, consistent with the hypothesis that already existing Hebbian cell assemblies firing in predefined temporal sequences could be the microcircuit substrate that encodes visual percepts changing in time.


Introduction
The coordinated activity of neurons firing in synchrony and organized in temporal sequences has been proposed as the generalized substrate for a wide variety of physiological computations and behaviors (Lorente de No, 1938;Hebb, 1949;Abeles, 1991;Seidemann et al., 1996;Baeg et al., 2003;Grinvald et al., 2003;Harris et al., 2003;Dragoi and Buzsáki, 2006;Fujisawa et al., 2008;Pastalkova et al., 2008;Buzsáki, 2010;Crowe et al., 2010;Dragoi and Tonegawa, 2011;Harvey et al., 2012). However, little is known about how sequential activity patterns are generated and whether they have a functional significance. For example, it has been proposed that specific groups of neurons firing in sequences could represent canonical microcircuit motifs (Hebbian cell assemblies) that can be recalled reliably (Hebb, 1949;Huyck, 2001;Lansner, 2009), allowing the formation of complex representations from basic multicellular functional units . Moreover, Hebb's cell assembly hypothesis postulates that neurons firing in sequential order will create activity cycles reverberating in the absence of external stimulation (Hebb, 1949). Accordingly, there is in vitro and in vivo evidence that the visual cortex can generate repeated patterns of spontaneous activity that resemble those observed under visual stimulation MacLean et al., 2005;MacLean et al., 2006;Miller et al., 2014). In addition, visual cortical neurons responding to the same orientation display similar spontaneous activity (Ch'ng and Reid, 2010) and have a higher probability to be connected among themselves (Ko et al., 2011), suggesting that synaptic plasticity, perhaps following Hebbian learning rules, could explain repeated temporal dynamics.
Using two-photon calcium imaging with single-cell resolution in vivo, we recently found that most cortical activity occurred in the form of coactive groups of cells defining neuronal ensembles (Miller et al., 2014). Although the same ensembles were repeatedly triggered by visual stimuli of the same orientation, these repetitions were never exact and neurons could promiscuously join different ensembles. Interestingly, the same ensembles that were visually evoked could be detected in the spontaneous activity, suggesting that cortical circuits may build a response to sensory stimuli using internal building blocks.
Although that recent work described the spatial structure of neuronal ensembles, it did not examine whether recurrent temporal patterns were also present. To further understand the temporal structure and properties of neuronal ensembles, we have now used independent analytical tools that capture spatiotempo-ral population dynamics  and applied it to the same dataset (Miller et al., 2014). As before, our approach allowed us to identify groups of neurons firing in synchrony that reliably represent specific visual scenes. However, in addition, we find that these ensembles of neurons from primary visual cortex encode natural stimuli with temporally structured network dynamics. Sequential activity patterns that emerge in the absence of visual stimulation are also recruited by recurrent natural scenes, demonstrating repeated temporally structured network activity in primary visual cortex of awake mice. Our results indicate that the representation of visual percepts could be implemented using endogenous spatiotemporal dynamics of neural circuits.

Note on present study
This study represents an independent analysis of experiments described previously (Miller et al., 2014). Although the raw data are mainly the same, the questions investigated and analytical tools are different, focusing on the detection and characterization of temporal patterns.

Animal surgery
All animal procedures were performed according to the guidelines of Columbia University Institutional Animal Care Guidelines and have been described previously (Miller et al., 2014). Briefly, experiments were performed on C57BL/6 male mice (n ϭ 4) or on parvalbumin-Cre (n ϭ 1) or somatostatin-Cre (n ϭ 2) ϫ LSL-tdTomato transgenic male mice obtained from The Jackson Laboratory at the age of postnatal day 40 -80. Mice were anesthetized with isoflurane (1-2%). A titanium head plate was attached to the skull using dental cement. The site for craniotomy was delimited using a dental drill over left V1. To facilitate head fixation, mice underwent training to maneuver on a spherical treadmill for 1-3 h for 2-3 d.

Functional multicell calcium imaging
On the imaging day, the part of the skull marked previously was removed while mice were anesthetized with isoflurane. Bulk loading was performed using the calcium indicator Oregon Green BAPTA-1 AM (Invitrogen) dissolved in 4 l of DMSO plus 20% Pluronic F-127 (Invitrogen) and diluted with 35 l of pipette solution (150 mM NaCl, 2.5 mM KCl, 10 mM HEPES, pH 7.4). Sulforhodamine 101 (Invitrogen) was added to the solution pipette to label astrocytes. Pipette solution was injected in layer 2/3 of V1 under visual control (10 p.s.i., 8 min per site) in 3 different locations. Cortical neurons were monitored with a two-photon microscope (Moveable Objective Microscope; Sutter Instruments) attached to a Ti: sapphire laser (Chameleon Vision II; Coherent) with 880 nm excitation light focused through a 20ϫ [0.95 numerical aperture (NA); Olympus] or 25ϫ (1.05 NA; Olympus) water-immersion objective. Scanning and image acquisition were controlled by the Sutter Instruments software Mscan. The optical field recorded was ϳ315 ϫ 315 m. Short videos (ϳ720 s) with a sample rate of 125-250 ms/frame were collected at time intervals of 5-10 min for up to 2 h.

Visual stimulation
Visual stimuli were generated using MATLAB (The MathWorks) Psychophysics Toolbox and displayed on a LCD monitor positioned 15 cm from the right eye at 45°to the long axis of the animal. Spontaneous activity in layer 2/3 of visual cortex was recorded under complete darkness in an optically isolated laser room at the beginning of the experiments and sometimes in the middle of the experiments (with a monitor and room lights turned off). The imaging setup and the objective were completely enclosed with blackout fabric and a black electrical tape and there was no detectable light coming into the mouse eye during the spontaneous condition. Visual stimuli consisted of full-field drifting gratings or natural scenes. Sine or square-wave gratings (100% contrast, 0.035 cycles/degree, 2 cycles/s) drifting in 8 different directions randomly were presented for 5 s, followed by 5 s of mean luminescence. Natural scenes consisted of 10 distinct natural scenes displayed repeatedly for 30 s using the MATLAB Psychophysics Toolbox synchronized with image acquisition.

Image processing
Image processing was performed with ImageJ version 1.42q software and custom-made programs written in MATLAB as described previously (Mao et al., 2001;Cossart et al., 2003;. Acquired images were processed to correct motion artifacts (Dombeck et al., 2007;Evangelidis and Psarakis, 2008;Kaifosh et al., 2013). Active neurons were automatically identified, coordinates for each cell in a 2D plane were assigned, and the calcium transients were measured as a function of time. Calcium transients were computed as changes in fluorescence: (F i Ϫ F 0 )/F 0 , where F i denotes the fluorescence intensity at any frame and F 0 denotes the basal fluorescence of each neuron (Miller et al., 2014). Spike probability was inferred from calcium signals using a fast, non-negative deconvolution method (Vogelstein et al., 2010) and then a threshold of 3 SDs above noise was determined from spike probabilities of the entire population in each experiment. The spike onsets inferred from calcium signals were used to represent the activity of the network. We constructed an N ϫ F binary matrix, where N denotes the number of active neurons and F represents the total number of frames for each video. For the analysis, only calcium transients elicited by neurons were considered. Each row in the binary matrix represents the activity of one neuron. To visualize all of the active neurons, the binary matrix was plotted as a raster plot in which ones are represented by dots. The time histogram of the raster plot illustrates the overall behavior of the network (Mao et al., 2001;Cossart et al., 2003;.

Analytical tools
Vectorization of the network activity. To identify neuronal pools firing in synchrony, temporal vectors were constructed that represent the simultaneous activation of different neurons. We tested the significance of the peaks of synchrony that we observed against the null hypothesis that the synchronous firing of neuronal pools is given by a random process (Shmiel et al., 2005(Shmiel et al., , 2006. We generated 1000 shuffled raster plots (see "Reshuffling methods" section below) and compared the distribution of the random peaks against the peaks of synchrony observed in the real data (Mao et al., 2001;Cossart et al., 2003;. Peaks of synchrony are defined by all of the adjacent high-activity frames in a given time window . Only the peaks of synchrony with more than the cells expected by chance ( p Ͻ 0.01) were considered for further analysis. Synchronously active neuronal pools define population vectors representing the overall activity of the network. Such vectors can be used to describe the network activity as a function of time (Schreiber et al., 2003;Stopfer et al., 2003;Brown et al., 2005;Sasaki et al., 2007;. The population vectors representing the network activity define a multidimensional space in which the number of dimensions is given by the total number of active cells in high-activity frames. It has been shown that the temporal vectorization of the network activity allows the discrimination of similar patterns repeated at different times (Schreiber et al., 2003;Ikegaya et al., 2004;Brown et al., 2005;. Term frequency-inverse document frequency normalization. After the identification of the vectors with low probability to be random events ( p Ͻ 0.01), we normalized the [N ϫ T] binary matrix (where N ϭ number of active cells and T ϭ the number of significant vectors) that represents the activity of the network. The normalization was based on a classic term frequency (TF)-inverse document frequency (IDF) search algorithm widely used in the field of natural language processing to recognize words that are key to relevant documents (Lan et al., 2009;Islamaj Dogan and Lu, 2010;Garbarine et al., 2011). The main purpose of the normalization is to find the neurons that are more relevant to specific visual stimuli. The TF-IDF matrix elements are given by TF*IDF for each active neuron at each significant frame. TF measures active neurons normalized by the total number of active neurons for each frame. However, cells that are frequently active have little effect in defining neuronal ensembles even though they can be key players orchestrating the transitions between different groups . Therefore, to reduce the impact of highly active neurons in the definition of a neuronal ensemble, we computed IDF, which measures the number of times that a neuron appears in the total number of significant vectors. For this particular case, IDF was defined as follows: IDF(i) ϭ log(total vectors/vectors(i)), where vectors(i) denote all of the significant vectors where neuron i appears. TF-IDF captures how important a word is in a collection of documents. In the context of the present work, TF-IDF defines the importance of a specific neuron in a collection of significant vectors that represent neuronal ensembles.
Representation of multidimensional datasets in time. To construct similarity maps, all of the possible vector pairs taken from the normalized TF-IDF matrix were compared. The similarity index between a pair of vectors is defined by their normalized inner product (Schreiber et al., 2003;Sasaki et al., 2007; as follows: cos() ϭ A ⅐ B /͉͉A͉͉ ͉͉B͉͉, which represents the cosine of the angle between two vectors. If two population vectors in an "n" dimensional space point in the same direction, the similarity index has a value of 1. In neural network context, high similarity values between two different time points describe repeated groups of neurons firing in synchrony. The inner product between different elements of a multidimensional space normalized by the magnitude of the elements has been widely used in search engines because it represents a computationally efficient algorithm that captures the similarities between documents. The visualization of the angles between all of the possible vector pairs represents a multidimensional reduction from "n" dimensions mapped into a 3D space, where the first two dimensions denote the significant vectors as a function of time and the third dimension indicates the angle between them. Significant patterns can be identified by vectors denoted by repetitive structures in the similarity maps (Victor and Purpura, 1996;Levy et al., 2001;Schreiber et al., 2003;Morelli et al., 2006;Kreuz et al., 2007;Tiesinga et al., 2008). To determine whether similar patterns of activity repeated in a statistically significant manner, we shuffled the normalized TF-IDF matrix, preserving the dimensionality of population vectors, and compared the probability distribution of similarity coefficients from real data and shuffled data.
Identification of neuronal ensembles. To detect different groups from multidimensional datasets, the cosine similarity of all possible vector pairs was used. Defined structures in the similarity maps (independent of the distance function chosen) have been widely used to detect recurrent patterns (Victor and Purpura, 1996;Levy et al., 2001;Schreiber et al., 2003;Morelli et al., 2006;Kreuz et al., 2007;Tiesinga et al., 2008). However, in those approaches, the independent variable considered to perform the computations was given by the identity of single neurons. Similarity maps generated by vectors in which dimensionality is defined by individual trials give valuable information about neuronal correlations, but cannot capture directly the sequential activity patterns that define a network. The main difference with our approach is that we consider the number of active cells as the dimension of the vectors, so our similarity maps can reveal temporal patterns (Sasaki et al., 2007;Carrillo-Reid et al., 2011). Neuronal ensembles are defined by a subset of neurons that reocurr at different times. To highlight significant patterns from similarity maps, we transform them in a binary matrix M of size (T ϫ T ). The cutoff value represents the percentage of total cells that are coactive in two different frames above chance levels. A neuronal ensemble is defined by groups of neurons with synchronous, recurrent, and alternate activity (Hebb, 1949;Carrillo-Reid et al., 2011). To identify neuronal ensembles, we obtained the singular value decomposition (SVD) of the significant patterns matrix. SVD has been widely used in signal processing analysis because it allows the identification of components representing the main features of a given system. Formally, the SVD of a real matrix can be expressed as the multiplication of three factors, U, Α, and V T , where U and V T represent orthonormal bases and the elements of Α are the singular values. In our particular case, the matrix M, which contains the information of the significant patterns of activity, is a real symmetric matrix; so V ϭ U. Therefore, the SVD is given by M ϭ V Α V T . The number of singular values that defines the total number of neuronal ensembles was determined by the rate of singular values decaying above chance levels. Rapid decaying singular values can reproduce original data with high accuracy. SVD is a method to factorize a dataset into linearly independent components that in this case represent neuronal ensembles. Analogous results were obtained with similar approaches such as locally linear embedding and clustering analysis . However, the main advantage of using SVD is the unsupervised identification of neuronal ensembles without the necessity to define a specific number of clusters; this fact is crucial to identify the most representative cells during the course of an experiment, opening the possibility to manipulate ad libitum identified neurons with single-cell resolution.
Identification of closed cycles with graph theory. After the identification of the neuronal ensembles, we looked for closed trajectories in the sequential activity patterns that are a characteristic of Hebbian cell assemblies (Hebb, 1949;. Neuronal ensembles transitions were plotted as a function of time (Lee and Wilson, 2002;Carrillo-Reid et al., 2011). From these transitions, we constructed isomorphic directed graphs with the properties of Hamiltonian or Eulerian cycles. To define sequential patterns mathematically that occur above chance levels, we just considered Hamiltonian or Eulerian closed cycles. The number of closed cycles that occurred in spontaneous activity or visually evoked activity was significantly different compared with the number of closed cycles from shuffled data. A Hamiltonian cycle of a directed graph (digraph) is a closed trajectory that contains every vertex of the graph exactly once. An Euler cycle is a closed trajectory passing through every edge exactly once (Diestel, 2005). Digraphs, represented by adjacency matrices, allow the automatic identification of template sequences. Template sequences indicate repetitive transitions between neuronal ensembles and are useful to identify spatiotemporal patterns with high probability of recurrence. The use of graph theory applied to neuronal ensembles transitions revealed that complex cycles of activity can be composed by smaller cycles or "primitives" that could represent basic procedures (Hebb, 1949;Bienenstock and Geman, 1995;Hammer, 2003;Abeles et al., 2004;. Cross-correlation analysis. We performed the normalized crosscorrelation between all of the possible combinations of neurons that belong to the same neuronal ensemble and between all of the combinations of neurons belonging to neuronal ensembles firing in sequence repeatedly. Note that, because of the multidimensional nature of the network activity, it is extremely difficult and time consuming to define neuronal ensembles from an exhaustive cross-correlation of all possible neuron pairs Carrillo-Reid et al., 2011). Cells belonging to the same neuronal ensemble present high cross-correlation values at zero time lag, whereas neurons from two sequential neuronal ensembles have high cross-correlation values at specific time lags that reflect the temporal structure between two defined groups.
Granger causality. To determine whether the occurrence of a given neuronal ensemble "causes" the emergence of another neuronal ensemble, we performed a Granger causality test (Granger, 1969) between the activation time courses of all possible combinations of neuronal ensembles. Although causality per se is a philosophical concept, the Granger approach infers a practical version of causality by assessing the degree to which time course "b" can predict time course "a" better than time course "a" alone. We estimated the activation time course for each neuronal ensemble based on the similarity index function of the most representative cells. We then performed a multiple regression of each ensemble onto several time-shifted versions of itself to obtain R 2 and then calculated the change in R 2 after adding time-shifted versions of other ensemble activation time courses (termed herein "Granger coefficients"). Granger coefficients reflected R 2 change in the time course of ensemble "a" gained by considering the back-shifted time course of another ensemble "b" in addition to the back-shifted time course "a." To discard nonzero Granger coefficients that could be occurring by chance, we shuffled the time courses of each ensemble 1000 times, calculating the random distribution for each shuffle separately for each experimental condition.
Doublets analysis. To identify sequential pairs of neurons with fixed time intervals, we performed an exhaustive matching algorithm that searches for a given time interval in all of the possible combinations of active neurons (Ayzenshtat et al., 2010). We consider as events only doublets that are present above chance levels compared with the network activity shuffled in time or in space. Doublet analysis emphasizes the existence of pairs of neurons with structured temporal activity.
Reshuffling methods. We performed four different spatiotemporal shuffling paradigms of our datasets, each catered specifically to a point of interest to robustly assess statistical significance in the absence of obvious parametric test for population analysis. These approaches involved shuffling 1000 times the binary matrices representing the overall activity of the network and/or ensemble activation time courses as follows: (1) preserving the interburst intervals but shuffling their orders across time within cells (used to identify the maximum number of cells coactive by chance), (2) preserving the activity within time frames but shuffling the events across cells (used to establish the stability and robustness of ensembles while maintaining the dimensionality of the population vectors), (3) preserving the temporal relationships between cells but shuffling the nature and direction of their sequences by swapping spikes at the pairwise-level, and (4) randomly time shifting the activity within each cell (used for assessing the significance of Granger coefficients by holding autocorrelations constant within cells).

Detecting neuronal coactivations in mouse visual cortex
To measure neuronal responses to visual stimulation in neuronal populations, we performed in vivo two-photon calcium imaging from layer 2/3 of primary visual cortex from awake, freely moving mice (Fig. 1A) because calcium imaging allows the spatial and temporal study of neural networks with single-cell resolution (Yuste and Katz, 1991;Ohki et al., 2005;Grewe and Helmchen, 2009;Wallace and Kerr, 2010). For each experiment, the imaged field of view was ϳ315 ϫ 315 m, containing on average 80 Ϯ 25 active neurons ( Fig. 1B; n ϭ 7 mice). For analysis, we assigned 2D spatial regions of interest to each cell in the field of view and measured their calcium transients as a function of time (Fig. 1C). To identify action potential activity represented by calcium transients, we estimated the spike probability as described previously (Vogelstein et al., 2010) ( Fig. 1D) and considered for analysis only those signals above a threshold Ͼ3 times SDs from the background fluorescence noise (Mao et al., 2001;Cossart et al., 2003;. These suprathreshold signals were used to generate binary arrays that depict population activity as raster plots (Fig. 1E), where each row represents an active neuron and each column in the raster plot signals the overall activity of the network at a specific time.
In these raster plots, neuronal activity was widely distributed across the entire population of imaged cells (Fig. 1E). Moreover, the time histogram of population activity ( per 300 s, each epoch comprises 6 Ϯ 3 active frames; n ϭ 7 mice). Only peaks of synchronous activity with negligible probability of occurring by chance were considered for further analysis ( Fig. 1E; p Ͻ 0.01; see Materials and Methods).

Multidimensional vectors depict repeated activity patterns
Raster plots can be understood as multidimensional arrays where each point in time defines a population vector ( Fig. 2A) the dimensionality of which corresponds to the total number of active cells (Stopfer et al., 2003;Brown et al., 2005;Sasaki et al., 2007;. Because cells with high activity levels can bias the identification of groups of neurons responding to specific visual stimuli, we performed TF-IDF of the binary arrays (Fig. 2B). TF-IDF has been used widely in search engines to lower the weight of frequent terms. In the context of neural activity, word vectors represent the activity of each cell and document vectors represent a population vector (Fig. 2B).
To detect repetitions from population vectors above chance levels, we computed the similarity coefficients (Fig. 2C) of the normalized TF-IDF matrix from all of the possible combinations of vectors among different experimental conditions and compared them with the similarity coefficients obtained from shuffled data preserving the activity within time frames (Fig. 2D). Similarity coefficients above chance levels demonstrated that similar population vectors occurred in spontaneous and visually evoked activity. These experiments demonstrated that population vectors defining network activity in primary visual cortex repeat above chance levels, which is consistent with our past work (Miller et al., 2014).

Identification of neuronal groups activated by visual stimuli
Even though responses to visual stimuli of layer 2/3 neurons in primary visual cortex have been widely characterized in many experimental conditions (Victor and Purpura, 1996;Tiesinga et al., 2008;Sawinski et al., 2009;Kampa et al., 2011;Marshel et al., 2011;Espinosa and Stryker, 2012;Gavornik and Bear, 2014), the identification of neuronal groups during the course of an experiment in which coordinated activity represents specific visual percepts remains unexplored because optical or electrical techniques that provide largescale multineuronal data are still relatively recent.
To identify groups of neurons that respond to a given orientation repeatedly, we used drifting gratings at four different orientations and two opposite directions. Drifting gratings were presented randomly intermixed with interstimulus intervals of the same duration (Fig. 3A), as described previously (Miller et al., 2014). Visual responses to drifting gratings evoked peaks of synchronous activity reflecting the overall behavior of the network (Fig. 3A, bottom).
Because individual cells can respond to multiple orientations (Victor and Purpura, 1996;Mrsic-Flogel et al., 2007;Marshel et al., 2011;Espinosa and Stryker, 2012), it is difficult to identify a representative group of neurons that respond to specific visual stimuli (Kampa et al., 2011). To find neuronal ensembles responding reliably to a given visual stimulus, we analyzed the network activity as a whole using population vectors because vectorization allows rigorous and quantitative comparison of activity patterns at different times (Schreiber et al., 2003;Brown et al., 2005;Sasaki et al., 2007;. To detect repeated patterns of activity, we constructed similarity maps of the normalized inner product of all possible vector pairs. Similarity maps still reflect the temporal characteristics of the network activity and can be used to visualize similar population vectors occurring at specific time points (Sasaki et al., 2007;Carrillo-Reid et al., 2011). Significant patterns obtained from these similarity maps represent the existence of stable structures in the responses (Victor and Purpura, 1996;Levy et al., 2001;Schreiber et al., 2003;Morelli et al., 2006;Kreuz et al., 2007;Tiesinga et al., 2008) that occurred above chance levels (Fig. 2D) and can be used to identify the time when similar population vectors were present (Fig. 3B). The matrix of significant patterns can be factorized by SVD in the three factors, U, Α, and V T . The strongest singular values that define the total number of neuronal ensembles were determined by the rate of singular values decaying above chance levels (Fig. 3C). On average, for the experimental condition of drifting gratings, we identified 6 Ϯ 1 ensembles ( Fig. 3C; n ϭ 7  data. Each factor from SVD depicts repetitive features of a multidimensional array that in our case represents the time when a specific orientation was presented (Fig. 3D). Repeated coactive patterns were readily observable and corresponded to distinctive orientation gratings (Fig. 3E), which is consistent with our previous work (Miller et al., 2014). Although our analytical tools allowed us to define neuronal ensembles that distinguished specific orientations accurately, we could not separate between opposite directions at the population level, consistent with the fact that direction selectivity of individual neurons in V1 is lower than orientation selectivity (Marshel et al., 2011). These results confirmed the existence of neuronal ensembles in primary visual cortex that respond to a given orientation and further extend our previous work (Miller et al., 2014) identifying automatically the time when a given ensemble occurs.

Cells with diverse orientation selectivity compose neuronal ensembles
To study the relation between neuronal ensembles and orientation selectivity, a known functional property of the mammalian primary visual cortex (Hubel and Wiesel, 1959), population responses to visual stimuli were studied. We identified the most representative neurons of each neuronal ensemble (base elements or core neurons) using a classic search engine algorithm to recognize words that are key to documents ( Fig. 2; see Materials and Methods). In the context of population dynamics, each neuron represents a word and documents define the collection of high-dimensional vectors comprising each group. Calcium transients from neurons belonging to a given ensemble showed that neuronal ensembles could be composed of cells that just responded to a given orientation and of cells responding to multiple orientations (Fig. 3F ). Base elements (core neurons) consisted of groups of 14 Ϯ 5 cells (n ϭ 7 animals) and reproduced the overall network behavior (Fig. 3G).
Ensembles did not simply reflect spatial clusters (Fig. 3H ) because the centroid of cells belonging to each ensemble was located in the middle of the optical field (Fig. 3I, left) and the mean distance between all of the neurons from each ensemble covered most of the recorded area (Fig. 3I, right). This altogether confirmed a "salt and pepper' anatomical distribution of neurons (Espinosa and Stryker, 2012) despite highly tuned visual responses. Cells belonging to core ensembles had low overlap (Fig.  3J ), indicating the possibility to manipulate identified subnetworks with single-cell resolution to drive behavioral outcomes. Interestingly, even if core ensembles were very selective for a given orientation, they could still have neurons with low orientation selectivity (Fig. 3K ). Our analysis thus demonstrated that specific assortments of neurons reliably respond to a given visual stimuli as a group, a first step to studying the temporal characteristics of neuronal ensembles.

Natural scenes trigger sequences of neuronal ensembles
We then turned our attention to the temporal analysis of the dynamics between neuronal ensembles using natural scenes as visual stimulus. Sequential activity patterns between neuronal groups have been proposed as associated cells that can act as a temporally closed system after sensory stimulation has ceased, constituting the implementation of a mental representation (Hebb, 1949). To identify repetitive sequential patterns between neuronal ensembles, we created an artificial experimental paradigm in which we presented 10 different natural scenes for 30 s repeated in the same sequential order many times. This enabled us to be certain that the transitions observed were evoked by our visual stimuli and not by external factors independent of visual stimulation. The expectation was to find closed cycles of activity between neuronal ensembles that were precisely related to recurrent stimuli repeated many times. To describe sequential activity patterns from cortical networks in vivo, we applied graph theory techniques to the temporal transitions between neuronal ensembles identifying template pathways embedded in closed cycles ( Fig. 4A; see Materials and Methods). For this analysis, we just considered Hamiltonian or Eulerian closed cycles that are represented by unique adjacency matrixes ( Fig. 4B; see Materials and Methods) and allow the identification of recurrent transitions between neuronal ensembles Carrillo-Reid et al., 2011).
Closed cycles have unique mathematical properties that allow their identification from a sequence of events. It is important to highlight that each transition implies specific groups of neurons firing in synchrony and recurrently above chance levels determined from shuffled data. To investigate whether closed cycles between neuronal ensembles in our experiments occurred above chance levels, we used a conservative approach. We shuffled the transitions describing a given experiment and compared the total number of Hamiltonian or Eulerian closed cycles that appeared in real data with the total number of closed cycles that can emerge by chance in spontaneous and visually evoked conditions. Indeed, closed cycles describing our experiments were significantly different from closed cycles emerging by chance (Fig. 4C). Moreover, the same closed cycles that appeared in real data were never found in shuffled data. We also wondered what the probability was of finding a given sequence above chance levels for a specific number of vertices. To answer this question, we calculated the probability distribution of finding a random sequence in a graph as a function of the number of vertices (Fig. 4D). In graphs with Ͼ4 vertices, the probability of finding a given sequence was higher than chance levels ( p Ͻ 0.05). Note that, for our analysis, each vertex represents a neuronal ensemble. This analysis demonstrated that the closed cycles and sequential patterns identified in our experiments could not be explained by chance occurrence.
We already demonstrated that neural activity defined by highdimensional vectors allows the reliable representation of specific drifting gratings (Fig. 3); however, visual experience in the real world is composed of the successive representation of complex scenes with defined spatiotemporal properties. We therefore tested the hypothesis that specific neuronal ensembles encoded different natural scenes. We expected that natural scenes presented repeatedly generate recurrent transitions between a constant subset of temporally organized ensembles. The network responses captured by raster plots (Fig. 5A, top) showed that natural scenes evoked peaks of synchronous activity without apparent periodicity (Fig. 5A, bottom) even if the same scenes were repeated in exactly the same order. However, similarity maps representing the angle between these high-dimensional population vectors depicted stereotypical trajectories showing repetitive responses illustrated by parallel lines to the main diagonal (Fig.  5B). The strongest factors of these high-dimensional arrays demonstrated that natural scenes could also be encoded by neuronal ensembles (Fig. 5C). As in the case of specific orientations, neuronal ensembles can distinguish accurately between different natural scenes (67 Ϯ 14% of matching; n ϭ 7 animals). Graph theory applied to the temporal transitions between neuronal ensembles allowed us to identify highly repetitive sequences embedded in closed cycles (Figs. 5C, colored arrows), as we expected from a recurrent set of natural scenes. Calcium transients from the most representative neurons of each ensemble demonstrated the exis-tence of sequential activity between neuronal pairs (Fig. 5D, colored arrows). To further corroborate that neuronal ensembles and sequential activity patterns were not an artifact of our analytical methods, we computed the cross-correlation of all possible combinations of neurons between cells belonging to each neuronal ensemble (Fig. 5E, top left) and cells belonging to different neuronal ensembles firing in sequence (Fig. 5E, bottom left). Accordingly, the cross-correlation values between cells belonging to the same neuronal ensemble displayed maxima at zero time lag, whereas cells belonging to sequential patterns showed high crosscorrelation values at time lags determined by the interval between different natural scenes (maximum normalized correlation coefficients: 0.18 Ϯ 0.082 for neuronal ensembles; 0.12 Ϯ 0.023 for sequential patterns; n ϭ 7 mice; for clarity, just one sequence is highlighted). The presence of cross-correlation peaks at nonzero lag (Fig. 5E, right) indicated that, during natural scenes, neuronal groups are temporally tied to visual stimuli. Similar to drifting gratings, the spatial locations of neurons belonging to each neuronal ensemble (Fig. 5F ) defining natural scenes presented a widespread distribution (Fig.  5G). Conversely, as for drifting-grating conditions, there was a low overlap between core ensembles (Fig. 5H ). Interestingly, neurons composing a given core ensemble exhibited either high or low orientation selectivity (Fig. 5I ). These results demonstrate that specific neuronal ensembles respond recurrently to a given set of natural scenes and that multidimensional population vectors comprise the temporal characteristics of the overall network activity and thus could serve to encode the temporal structure of the stimulus.

Spontaneous neuronal coactivations show sequential activity patterns
Our previous results showed that different natural scenes changing in time could be encoded by a stereotyped succession of neuronal groups or ensembles. However, whether primary visual cortex can generate sequential activity patterns in the absence of repetitive visual stimuli remains unknown. If endogenous temporal cortical patterns are the substrate of dynamic visual experience, it is predicted that the spontaneous activity of neuronal populations should also generate synchronous, recurrent, and sequential activity, as has been recently demonstrated in auditory cortex (Luczak et al., 2013). To test this, we recorded the spontaneous activity of cortical neurons that were never exposed to our natural scenes stimuli and applied the same analytical tools that we used previously to detect sequential activity patterns in population activity.
Spontaneous activity was recorded with monitors and lights turned off. The imaging setup and the objective were completely enclosed with blackout fabric and black electrical tape so that no detectable light cues were coming into the mouse's eyes during spontaneous recordings. Consistent with previous work (Miller et al., 2014), raster plots representing the overall activity of the recorded neurons (Fig. 6A) showed that the number of synchrony peaks between spontaneous and evoked activity were not significantly different (Fig. 6A, bottom; 76 Ϯ 35 peaks in sponta-A B C D Figure 4. Identification of closed cycles of activity describing sequential activity patterns using graph theory. A, Cell assemblies described by Hebb. Each letter represents a group of neurons firing in synchrony. Arrows depict the transition between two neuronal ensembles. Numbers indicate the transition between different neuronal ensembles as a time function. Letters on bottom indicate the "sentence" describing all of the transitions. The graph depicting sequential patterns defined by Hebb can be transformed by graph theory in an isomorphic graph in which each vertex represents only one vertex of the original graph. Numbers on bottom indicate transitions between neuronal ensembles in function of time. B, Closed cycles of activity extracted from the isomorphic graph that depicts the transitions between neuronal ensembles. Each closed cycle is defined by a unique adjacency matrix. Note that some pathways can be part of different closed cycles. C, Total number of closed cycles that can be found during spontaneous ( p Ͻ 0.0001; Mann-Whitney test) and visually evoked activity ( p ϭ 0.0007; Mann-Whitney test) in real and shuffled data. Note that exactly the same closed cycles that appeared in real data were never found in shuffled data. D, Probability distribution of finding a specific sequence in a closed cycle as a function of the number of vertices. Note that, for closed cycles with Ͼ4 vertices, the probability of finding a specific sequence is above chance levels ( p Ͻ 0.05).
neous activity vs 75 Ϯ 23 peaks in natural scenes; n ϭ 7 mice). As in the case of natural scenes, similarity maps of population vectors showed repetitive structures suggesting the existence of neuronal ensembles firing in sequences (Fig. 6B). The strongest factors of our multidimensional datasets representing spontaneous activity of awake mice demonstrated recurrent neuronal ensembles ( Fig. 6C; 5 Ϯ 1 neuronal ensembles/300 s; n ϭ 7 mice). Graph theory techniques applied to the transitions between ensembles allowed us to identify recurrent sequential patters between neuronal ensembles generated endogenously (Figs. 6C,D, colored arrows). Intriguingly, calcium transients from pairs of neurons belonging to a given sequential pattern during spontaneous activity showed 300 Ϯ 85% more repetitive transitions (Fig. 6D, colored arrows) than the ones observed at the population level in the same time window. Note that, at the population level, many cells must fire synchronously to be considered part of a neuronal ensemble. Neurons belonging to the same ensemble showed high crosscorrelation peaks at zero time lag (Fig. 6E, top left), whereas cells belonging to a recurrent sequence showed high cross-correlation values at nonzero time lags (Fig. 6E, bottom left, only one sequence is highlighted). Conversely, cross-correlation peaks at nonzero time lags showed that spontaneous ensembles do not appear periodically (Fig. 6E, top right) and that sequential patterns between ensembles have fixed time lags (Fig. 6E, bottom  right). Similar to drifting gratings and natural scenes the spatial distribution of neurons belonging to each neuronal ensemble during spontaneous activity (Fig. 6F ) was widespread (Fig. 6G). In addition, neurons belonging to core ensembles showed even less overlap than all of the neurons associated with a given ensem- ble (Fig. 6H ) and, as it was for drifting gratings and natural scenes, the elements composing spontaneous ensembles displayed a wide range of orientation selectivity indexes (Fig. 6I). Our results demonstrate the existence of preferred transitions between two different neuronal groups during spontaneous activity.

Similarity between spontaneous and evoked sequential neuronal coactivations
Finally, we investigated whether the temporal patterns found in spontaneous activity were related to those under visual stimulation.
To do this, we first recorded spontaneous population activity in the absence of visual stimuli, making sure that no other visual cues were seen by the mice (see Materials and Methods), and then presented natural scenes repeated in the same order several times. We already demonstrated the existence of highly similar population vectors between spontaneous and visually evoked activity above chance levels ( Fig. 2). Because population vectors represent neuronal ensembles, our working hypothesis was that neuronal ensembles evoked by visual stimulation also appeared spontaneously in the absence of external stimuli. Indeed, the strongest factors resulting from the network activity demonstrated that 85 Ϯ 10% ( Fig. 7A; n ϭ 7 mice) of the ensembles that appeared before natural scenes exposure also recurred during visual stimulation, which was consistent with our previous work (Miller et al., 2014). The representation of neuronal ensemble transitions as directed graphs revealed that the same recurrent transitions between neuronal ensembles generated endogenously in the spontaneous activity were also evoked by visual stimulation (Figs. 7A,B, colored arrows; 50 Ϯ 15% of sequential activity patterns evoked by visual stimuli were present in spontaneous activity; n ϭ 7 mice). Calcium transients from the most representative elements of each neuronal ensemble illustrate that the same neurons had recur- rent transitions during spontaneous activity and visual stimuli (Fig. 7C, colored arrows).
To investigate whether recurrent transitions observed in both spontaneous and visually evoked activity occurred above chance levels, we shuffled the transitions of spontaneous activity and searched for the same recurrent transitions during natural scenes and in shuffled data. This analysis demonstrated that recurrent sequences considered in the present analysis were indeed above chance levels (Fig. 7D). Interestingly, the mean time interval of these transitions did not significantly differ between spontaneous and natural scenes conditions (Fig. 7E), suggesting that pairs of neurons belonging to a given sequential pattern could have fixed time intervals. These data showed that specific sequential activity patterns occurring in spontaneous activity above chance levels could be entrained by visual stimulation as if visual stimuli used temporal motifs already present in cortical microcircuits.

Precise firing sequences underlie sequential patterns
The fact that similar sequential patterns existed at different time suggests that spontaneous cortical activity could be represented by a highly organized temporal structure, such as precise firing sequences previously described in monkey cortex during behavioral tasks (Abeles et al., 1993). To test this hypothesis, we used a completely independent analysis focusing on repetition of neuronal pairs (doublets), defined as a temporal sequence between two cells with a fixed time interval (Fig. 8A). We only considered as doublets those that were repeated above chance levels after spatial and temporal shuffling (Fig. 8B). Exhaustive template matching of doublets from spontaneous activity of cortical cells (Fig. 8C; 20 Ϯ 10 doublets for 5 min intervals; n ϭ 5 animals) demonstrated multiple pathways between pairs of neurons reverberating over relative long time periods (one neuron can be part of 5 Ϯ 3 doublets); 73 Ϯ 8% of neurons that participate in doublets also take part in the peaks of synchronous activity that defined neuronal ensembles (Fig. 8C, bottom). Doublets analysis revealed the existence of an underlying temporal structure between neuronal pairs of endogenous generated sequential activity patterns. Interestingly, doublets appeared more frequently than sequential transitions between neuronal ensembles, suggesting that the synchronization of neurons that define an ensemble can be composed by multiple motifs from a wider area than the A B C D E Figure 7. Similarity between endogenous and evoked sequential activity patterns. A, Neuronal ensembles sorted in time representing spontaneous activity (left) can also be evoked by repetitive visual stimulation (right). Note the existence of stereotyped activity during natural stimuli conditions reflecting the same series of scenes presented many times. Note that sequential activity patterns elicited by natural stimuli were also present in spontaneous activity before exposure to visual stimuli. B, Graph theory techniques applied to the transitions between neuronal ensembles depicted sequential activity patterns (colored arrows) that appeared spontaneously (left) and in the presence of natural stimuli (right). Each colored arrow defines a pathway between two neuronal ensembles that is present in both conditions. C, Calcium transients of the most representative cells belonging to each neuronal ensemble in the absence (left) and the presence (right) of visual stimuli. Percentages indicate the conditional probability to go from one neuronal ensemble to another. Note that sequential patterns recruited by natural stimuli are also present in spontaneous activity. Visual stimulation triggers endogenous generated neuronal ensembles, indicating that sequential activity patterns recruited by sensory input are imprinted in the network. D, Total number of sequences that can be found in spontaneous activity that were recruited by visual stimulation compared with the number of sequences that can be found in shuffled data. Note that the recurrence between spontaneous and natural scenes is statistically significant compared with the shuffled condition ( p ϭ 0.0078; Wilcoxon matched-pairs signed-rank test). E, Interval of recurrent sequences in spontaneous activity compared with natural scenes is not significantly different ( p ϭ 0.1488; Mann-Whitney test).
optical field recorded. Our experiments demonstrated, by two completely independent analytical approaches, that temporally structured sequential activity generated spontaneously (in the absence of visual stimuli) can be repeated in the same order at different time points, showing the existence of precise firing sequences in these recordings (Abeles et al., 1993).

Prediction of future neuronal ensembles
If neuronal ensembles occur in repeated sequences, then the activation of some ensembles, either spontaneous or evoked, should predict the occurrence of the next ensemble within a given sequence. To test this hypothesis, we first calculated Grangercausal coefficients on estimated activation time courses for each neuronal ensemble pair transition ( Fig. 9A; see Materials and Methods). Granger causality was developed in economics, but has been also used in neuroscience as a straightforward computation to assess the potential for causal relationships between variables. We compared our estimated Granger coefficients (estimates of causal influence of one ensemble on the activation of another in a given sequence) to a bootstrapped random distribution of coefficients from time-shuffled data to identify ensemble transitions that are statistically significant and compared these coefficients for the same ensemble pairs between natural scenes and spontaneous conditions. Based on the max-imum interval between activated ensembles in spontaneous activity (Fig. 6E, bottom right), we included time lags of up to 15 s in our analysis. Many ensemble transitions (Fig. 9B, coefficients, black triangles) were statistically significant above the 95 th percentile level (Fig. 9B, shaded area) for both spontaneous and natural scenes conditions. Because natural scenes stimuli displayed temporal regularity (repeats) and likely evoked temporally structured activity, it is expected that the positive tail bootstrapped random distribution of coefficients is longer than in the spontaneous case, yet many of the Granger coefficients remain demonstrably significant (also as expected). From the scatterplot of coefficients (Fig.  9C), it is interesting that, whereas some significant coefficients (Fig. 9C, red circles) have high values only in the externally driven natural scenes condition, the opposite is not true for the spontaneous condition. Therefore, spontaneous sequential activity patterns that occur above chance level could be entrained by natural stimuli, suggesting an intrinsic temporal structure in the network.
We demonstrated that neuronal ensembles define sequential activity patterns. Because the activation of a specific group of cells is temporally tied to another set of neurons, visually evoked activity recruits a specific subset of sequential patterns that could also be observed spontaneously (Fig. 10). Therefore, sensory evoked activity could recruit preexistent temporal pathways from spontaneous activity.

Discussion
We used two-photon calcium imaging of layer 2/3 of primary visual cortex from freely moving mice combined with neuronal population analysis to demonstrate that: (1) analysis of network activity using high-dimensional population vectors allows the identification of groups of neurons (ensembles) that repeatedly respond to specific visual stimuli, (2) neuronal ensembles repeat in sequential patterns that can be observed in endogenous and visually evoked activity, (3) temporally precise firing sequences are present at the single-cell level, and (4) sequential patterns can be used to predict the occurrence of future neuronal ensembles.

Identification of neuronal ensembles in multineuronal recordings
Although the methodology of calcium imaging has allowed the study of multineuronal activity with single-cell resolution for more than two decades ( Dombeck et al., 2009;Harvey et al., 2009;Miller et al., 2014), the analytical tools with which to characterize network activity are still in development Churchland et al., 2012;Cunningham and Yu, 2014), even though they are as important as the optical tools with which to perform the recordings. In this study, we used linear algebra to extract the spatiotemporal features of network activity and demonstrate the existence of specific groups of cells in primary visual cortex with synchronous, recurrent, and alternating activity.
The main difference between our analytical tools and previous approaches searching for temporal properties of network activity resides in the representation of the overall activity as an array of multidimensional population vectors in which the main focus is on time points. Because each population vector is independent of the recording length our approach allows faster identification of recurrent groups of cells firing together. An alternative approach to finding ensembles could be to compute all possible correlations between all cells recorded. However, cross-correlation be- Figure 10. Diagram summarizing sequential activity patterns. Visual stimuli recruit preexistent sequential activity patterns. Each neuronal ensemble is defined by a group of cells firing in synchrony. Numbers indicate different neuronal ensembles (colored neurons denote active cells). The same groups of cells could define neuronal ensembles during spontaneous activity and visual stimulation. Arrows represent sequential transitions between different neuronal ensembles. Sequential activity patterns entrained by visual stimuli could drive behavioral outcomes, perform pattern completion, or transfer time information to other cortical areas. Red arrows denote sequential patterns recruited by visual stimulation; blue arrows indicate preexistent sequential patterns present in spontaneous activity. Sensory evoked activity can recruit preexistent temporal pathways found in spontaneous activity.
tween pairs of neurons depends on the recording length, making computation time longer.
The understanding of the spatiotemporal properties from network activity considered as multidimensional arrays involves the use of dimensionality reduction techniques (Sasaki et al., 2007;Carrillo-Reid et al., 2011) and a priori selection of a few components that preserve the features of interest (Churchland et al., 2012). Because we were interested in the nonsupervised identification of groups that can represent a series of visual stimuli, we used an alternative approach consisting of the visualization of the angles between all of the possible pairs of vectors (Sasaki et al., 2007;. The angles between all of the vectors define similarity maps, which represent a reduction of the dataset dimension but allow the visualization of the original temporal features. The factorization of the significant patterns obtained from similarity maps can then be used to extract the time when a specific group of neurons fired recurrently (Fig. 3). The main difference between SVD factorization and other dimensionality reduction algorithms is that the number of neuronal ensembles is given by the magnitude of the singular values, so it is not necessary to define a priori the number of groups. We demonstrated that the strongest factors taken from SVD represent groups of neurons firing recurrently to specific visual stimuli (Fig. 3).

Evoked neuronal ensembles recapitulate spontaneous ones
Groups of neurons responding to visual stimuli were also present spontaneously in the absence of visual stimulation. Because our recordings from spontaneous activity were performed in the absence of any visual cues or external factors (see Materials and Methods), the most parsimonious interpretation is that neuronal ensembles found in spontaneous activity represent groups of cells that were imprinted in cortical circuits during the animal's life. The fact that neurons with similar tuning properties have higher probability to be connected (Ko et al., 2011) supports this hypothesis. Indeed, spontaneous activity in auditory cortex also shapes the responses to sensory stimuli Harris et al., 2011;Luczak et al., 2013).
Interestingly, neuronal ensembles, responding to a specific visual stimulus, are spatiotemporally flexible. Therefore, neurons belonging to different ensembles overlapped (ϳ30%; Figs. 3J, 5H, 6H ). At the same time, the most representative neurons of each ensemble (core) have significantly less overlap (ϳ10%) with other cores (Figs. 3J, 5H, 6H ), indicating that a subset of neurons remains stable over long time periods. Recently, it has been shown that strongly coupled neurons remain stable during sensory stimuli or spontaneous activity, supporting the idea of a finite repertory of stable microcircuit motifs (Okun et al., 2015). The stability of specific groups of neurons opens the possibility of manipulating targeted populations using optogenetics with single-cell resolution. The relation between population coupling and cortical ensembles deserves further study because it may reveal general connectivity rules between groups of neurons related to physiological events. Interestingly, neurons with broad orientation tuning can be part of core ensembles (Figs. 3K, 5I, 6I ), suggesting that they could encode different aspects of visual stimuli, the spatiotemporal environment, or the context as a whole. Because in primary visual cortex, interneurons have broad orientation properties (Kerlin et al., 2010) and high population coupling (Okun et al., 2015), they could have a pivotal role in the orchestration of neuronal ensembles.
The widespread spatial distribution of neuronal ensembles could suggest an optimized cortical architecture to allow parallel processing of external events and fastest exchange of information between different representations. However, the lack of a distinct anatomical organization found in the present study could be a limitation of the 2D recording techniques. The development of real time holographic imaging will be useful for our understanding of cortical microcircuits in 3D (Nikolenko et al., 2008).

Hebbian cell assemblies in cortical microcircuits
Our experiments show that visual stimulation is able to recruit endogenous generated sequential patterns with temporally structured characteristics. The existence of closed cycles of activity with high recurrence was proposed by Hebb as a mechanism to capture, maintain, and consolidate relevant sensory inputs (Hebb, 1949). Our experiments are consistent with Hebb's hypothesis that cortical networks are indeed reverberant microcircuits composed of neuronal ensembles that are activated in closed pathways. Sequential activity patterns have been suggested as a mechanism underlying delay-based tasks (Abeles et al., 1993;Seidemann et al., 1996;Funahashi, 2006), spatial navigation (O'Keefe and Recce, 1993), complex motor actions (Hahnloser et al., 2002), and sensory perception (Luczak et al., 2007;Shusterman et al., 2011). The observation of recurrent sequences suggests that prolonged maintenance of neural activity could be mediated by the continuous replaying of sequential patterns, which is contrary to the view that prolonged neural activity represents tonic firing of individual cells.
Because our goal was the study of ensemble temporal properties, we generated recurrent sequential patterns artificially by repetitive stimuli (natural scenes). Graph theory allowed us to identify closed cycles of activity between ensemble transitions that occur above chance level (Fig. 4). A very small fraction of all possible sequences occurring above chance level is sufficient for successfully reconstruct visual scenes (Fig. 7). Further characterization of the origin, mechanisms, and function of endogenous neuronal ensembles in primary visual cortex that are not evoked by visual stimuli goes beyond the main goal of our study but deserves further investigation and the development of behavioral experimental paradigms.
Even though the activity patterns described in the present work are limited by the temporal resolution of two-photon imaging, the fact that we were able to identify groups of neurons responding reliably to repeated natural stimuli supports the idea of sequential patterns as being physiologically relevant events in primary visual cortex. In fact, the high redundancy of doublet pathways (Fig. 8) suggests dynamic networks with high flexibility that could allow the representation of infinite environments with a finite number of primitives. Therefore, at the mesoscale level, behavioral relevant timescales could be achieved by the composition of small spatiotemporal building blocks (Bienenstock, 1995). In addition, because sequential patterns have the capability to track long time delays, fixed time intervals between diverse ensembles could carry information about how much time passed but also about future events (Fig. 9). Consistent with this, sequential activity patterns have been proposed recently as being the underlying mechanism for episodic memory (Eichenbaum, 2013) that, in the case of visual cortex, can represent pattern completion (Gavornik and Bear, 2014) or stimulus processing (Ayzenshtat et al., 2010). In addition, abnormal neuronal firing patterns are associated with different neuropsychiatric disorders (Schnitzler and Gross, 2005;Uhlhaas and Singer, 2006). Therefore, during pathological states, dissonances in sequential activity between neuronal ensembles could generate abnormal rhythms as if neurons were not able to follow the external world. The ability to identify spatiotemporally structured groups of neurons that are related to physiological relevant events opens the possibility of restoring the normal flow of information using optogenetics with single-cell resolution.