# Methods

**Animals.** Animal care and use was performed in compliance with the Yale IACUC, U. S. Department of Health and Human Services and Institution guidelines. Neonatal Ai38 floxed GCaMP3 reporter mice (JAX no. 014538)58 crossed with Rx-Cre[^Swindell:2006] or SNAP25-GCaMP6 (Ai103) transgenic mice aged 2–13 days (P2-P13) after birth (P0) were used.

**Transgenic mice generation.** The Snap25 locus was targeted with homologous recombination to insert LSL-F2A-GFP at the endogenous stop codon, using a targeting vector containing the following components: 5’ arm – F3 – last ~300bp of intron 7 – exon 8 up to the endogenous stop codon – loxP – stop codons – PGK polyA – loxP – F2A-EGFP – WPRE – bGH polyA – AttB – pPGK – neomycin-resistant gene – PGKpA – F5 – mRNA splice acceptor – domain 2 from the hygromycin-resistant gene – SV40 polyA – AttP – 3’ arm. Targeting constructs were generated using a combination of molecular cloning, gene synthesis (GenScript, Piscataway, US) and Red/ET recombineering (Gene Bridges, Heidelberg, DE). The 129S6B6F1 ES cell line, G4, was used for the gene targeting. Correctly targeted neomycin resistant clones were identified by PCR, and then confirmed by Southern blot.  One ES clone that gave high percentage chimeras following blastocyst injection was used in a Flp recombinase-mediated cassette exchange (RMCE) transfection to switch the LSL-F2A-GFP expression unit for a T2A-GCaMP6s expression unit.  The replacement vector included: F3 – last ~250bp of intron 7 – exon 8 up to the endogenous stop codon – T2A-GCaMP6s – WPRE – bGH polyA – AttB – pPGK – domain 1 from the hygromycin-resistant gene – mRNA splice donor – F5. Following co-transfection with a pCAG-Flpe plasmid, hygromycin resistant colonies were screened for correct 5’ and 3’ junctions and for lack of the original vector sequence. Correct clones were used in blastocyst injections to obtain germline transmission. Resulting mice were crossed to the Rosa26-PhiC31 mice (JAX Stock # 007743) to delete the pPGK-hygro selection marker cassette (in between the AttB and AttP sites), and then backcrossed to C57BL/6J mice and maintained in C57BL/6J congenic background. 

**Surgical procedure for in vivo imaging.** Mice aged P2-P13 were deeply anesthetized with isoflurane (2.5%) in oxygen and then placed on a heating pad set to 35ºC via a isothermic temperature monitor (NPI TC–20, ALA Scientific). Local anesthesia was produced by subcutaneous injection (0.05 ml) of 1% Xylocaine (10 mg/ml lidocaine/0.01 mg/ml epinephrine, AstraZeneca) under the scalp. After removal of the scalp, steel head posts were fixed to the exposed skull using cyanoacrylate glue. A 1 hr recovery period in the dark under continuously delivered medical oxygen with isoflurane at 0% was allowed after surgical preparation and the mouse was surrounded by a cotton ball nest. This recovery period was the typical minimum time required for spontaneous waves of activity to develop in the visual system after the cessation of deep anesthesia[^Ackman:2012]. A red LED (Radioshack) and photodiode sampled at 25 kHz with a Power1401 (Cambridge Electronic Design) were positioned to monitor respiratory rate and limb/body movements.

**Wide field calcium imaging.** A 16 bit CMOS camera (pco.edge, PCO) coupled to a Zeiss AxioZoom V16 microscope with 1X macro objective was used to image transcranial calcium dynamics. Epifluorescent illumination was provided by a DC stabilized Hg2+ light source (X-Cite, EXFO) through a EGFP filter cube set (Zeiss) with the minimum illumination intensity that gave detectable calcium signals using a exposure of 200 msec. Image frames corresponding to a field of view of 6 x 8 mm, 10 x 12 mm, or 20 x 24 mm were acquired at a rate of 5 or 10 frames per second (200 ms or 100 ms frame period). Each recording consisted of a single, continuously acquired movie for 10 min.

**Calcium signal detection.** Image processing and calcium signal detection was performed using custom software (available at [http://github.com/ackman678/wholeBrainDX](http://github.com/ackman678/wholeBrainDX)) written in MATLAB (The Mathworks, Natick, MA). The mean pixel intensity at each pixel location, F0 was subtracted and normalized to each frame, *Ft* of the movie to form a dF/F array: *A = (Ft - F0)/F0*. A background estimate was calculated and subtracted from every frame of A with a top hat filter using a disk shaped structuring element with radius of 620 µm. Each frame was smoothed with a Gaussian having a standard deviation of 56 µm and a signal intensity threshold, T was computed using Otsu’s method on the set of histogram of pixel intensities in A that corresponded to the set of pixels at the 99th percentile from the Sobel gradient transformation of A. Calcium domain signals were automatically segmented as contiguously connected components in space and time from the binary mask array, *A > T*. Components located outside the cortical hemisphere boundaries or having an area < 50 pixels or a duration of 1 frame were ignored.

**Statistical analysis.** Data sets were analyzed and plotted using custom routines implemented in MATLAB (The Mathworks, Natick, MA) and in R (The R Project for Statistical Computing, http://www.r-project.org) with the ggplot2 plotting library (http://ggplot2.org). Distribution means were compared using two-sample Student’s t-tests unless otherwise noted (Wilcoxon Rank Sum test was used for small, non-normally distributed data sets) or using ANOVA followed by Tukey’s HSD post-hoc test when analyzing the effects of multiple grouping factors (p < 0.05 set as significance). Values are reported as means with the standard error of the mean or medians with the median absolute deviation unless otherwise noted. All boxplots report the median, the 25th and 75th percentiles as lower and upper box hinges (1st and 3rd quartiles), the data range as lower and upper whiskers (lowest and highest data values within 1.5 * IQR of the lower and upper hinge, where IQR is the difference between the 3rd and 1st quartiles), and outliers as individual gray points (data values beyond the whisker range).  

**Calcium domain analysis.** The mean width in the medial-lateral and height in the rostral-caudal dimensions of the bounding box fitted to each segmented calcium domain signal was taken to be the domain diameter. The number of contiguous frames (bounding box depth) for each segmented calcium domain was taken to be the domain duration. The mean and maximum pixel intensities from A within each domain were taken as the mean and maximum domain amplitudes. After functional parcellation of the cortical hemipsheres (**Fig. 1**; [Supplementary Fig. 1](SupplementaryFig1.pnf)), domains were assigned areal membership by intersection of the domain centroid with a cortical area’s pixel mask. The number of individual domains per recording within a hemisphere or cortical area was taken to be domain frequency. 

**Wave motion analysis.** Optical flow was computed using the Lucas-Kanade method on the binary movie array from all the segmented calcium domain masks for a recording. The motion magnitude, *R* was the velocity vector sum for each domain. The wave motion index for each domain was calculated as *R^2^/D^2^*, where *D* was the domain diameter.

**Motor movement analysis.** A binary movie array from all the segmented calcium domain masks for a recording was intersected with masks representing different cortical area parcellations. The total number of active pixels per frame expressed as a fraction of the total number of possible pixels per frame for each cortical area gave active pixel fraction time courses for each cortical area in each recording. Movement signals acquired with the photodiode were bandpass filtered using an 8-order, 1-20 Hz pass band elliptic filter, and then rectified and downsampled to the movie frame rate to give a movement time course signal that corresponded to displacements of the limbs and body excluding those from respiration. Cross-correlations between combinations of cortical active fraction time courses and the movement signal time course were computed, and the Pearson’s correlation coefficient at the zeroth time lag was obtained.

**Hemisphere correlation analysis.** A binary movie array from all the segmented calcium domain masks for a recording was intersected with masks representing the cortical hemispheres. The total number of active pixels per frame expressed as a fraction of the total number of possible pixels per frame for each cortical area gave active pixel fraction time courses for each hemisphere in each recording. The mean, normalized medial-lateral and anterior-posterior positions for all active pixels within each hemisphere during coactive frames gave spatial center of mass timecourses for each recording. Pearson’s correlation coefficient was calculated between the active pixel fraction time courses and the activity center of mass time courses to give the temporal and spatial correlation for each movie.

**Functional connectivity analysis.** A binary movie array from all the segmented calcium domain masks for a recording was intersected with masks representing different cortical area parcellations. The total number of active pixels per frame expressed as a fraction of the total number of possible pixels per frame for each cortical area gave active pixel fraction time courses for each cortical area in each recording. Correlation matrices were calculated for each recording by computing pairwise Pearson’s correlation coefficients, *r*, from the matrix containing the cortical active pixel fraction time courses. The mean correlation matrix for each age group was computed and then the binarized correlation matrix at *r* > 0.15 (the maximum *r* value where the 3 largest communities were connected in the graphs at all age groups) was used to form an adjacency matrix with each node representing a cortical area and each edge representing an association between a pair of nodes at weight, *r*.

**Network analysis.** Graph theoretical analyses were performed using the igraph network analysis software library (http://igraph.org). Community structure was detected within each functional association matrix using a greedy optimization algorithm that maximizes the graph modularity score to perform hierarchical clustering[^Newman:2004][^Clauset:2004], where the modularity score measures the fraction of edges within modules for a graph partition compared with that of a randomized equivalent network. Network graphs were plotted using an anatomical layout or using a force-directed graph layout[^Fruchterman:1991] with nodes colored by module membership and edges connecting nodes reflecting the edge weight, *r*. Node degree was the number of connections that link a vertex to the rest of the network. The average path length, *L* of a graph was the mean of the shortest paths (fewest number of edges) between all pairs of nodes. The random average path length, *L~r~* was the mean of the shortest paths in a set of 1000 equivalent random networks that had the same degree sequence as the original graph. The local clustering coefficient was the ratio of the triangles connected to the node and the triples centered on the node, measuring the probability that two neighbors of a node are also connected. The global clustering coefficient, *C* was the ratio of the triangles and connected triples in the graph. The random global clustering coefficient, *C~r~* was the mean of the clustering coefficients in a set of 1000 equivalent random networks that had the same degree sequence as the original graph. The small-world index was calculated as the ratio of the normalized clustering coefficient (*C/C~r~*) and the normalized path length (*L/L~r~*), where a small-world index > 1 indicates a small-world network organization[^Humphries:2008][^Heuvel:2014]. Node strength was the column sums in the weighted adjacency matrix. Betweenness centrality scores corresponded to the fraction of all shortest paths that pass through a node[^Brandes:2001]. Eigenvector centrality scores were the values of the first eigenvector of the association matrix[^Bonacich:2007][^Lohmann:2010], reflecting for each node the sum of direct and indirect connections of every length in a network.
