DOI: 10.7291/d1n96w
title: Data-driven segmentation of cortical calcium dynamics

description: Demixing signals in transcranial videos of neuronal calcium flux across
 the cerebral hemispheres is a key step before mapping features of cortical
 organization. Here we demonstrate that independent component analysis can
 optimally recover neural signal content in widefield recordings of
 neuronal cortical calcium dynamics captured at a minimum sampling rate of
 1.5×106 pixels per one-hundred millisecond frame for seventeen minutes
 with a magnification ratio of 1:1. We show that a set of spatial and
 temporal metrics obtained from the components can be used to build a
 random forest classifier, which separates neural activity and artifact
 components automatically at human performance. Using this data, we
 establish functional segmentation of the mouse cortex to provide a map of
 ~115 domains per hemisphere, in which extracted time courses maximally
 represent the underlying signal in each recording. Domain maps revealed
 substantial regional motifs, with higher order cortical regions presenting
 large, eccentric domains compared with smaller, more circular ones in
 primary sensory areas. This workflow of data-driven video decomposition
 and machine classification of signal sources can greatly enhance high
 quality mapping of complex cerebral dynamics.
description: <strong>Mice </strong>
 All animal studies were conducted in
 accordance with the University of California, Santa Cruz Office of Animal
 Research Oversight and Institutional Animal Care and Use Committee
 protocols. P21-22 Snap25 GCaMP6s transgenic mice (JAX: 025111), Cx3cr1 GFP
 (JAX: 005582), and Aldh1 GFP (MGI: 3843271) were maintained on a C57/Bl6
 background in UCSCs mouse facilities. To identify Snap25 GCaMP expressing
 mice, a single common forward primer (5’-CCC AGT TGA GAT TGG AAA GTG-3’)
 was used in conjunction with either transgene specific reverse primer
 (5’-ACT TCG CAC AGG ATC CAA GA-3’; 230 band size) or control reverse
 primer (5’-CTG GTT TTG TTG GAA TCA GC-3’; 498 band size). The expression
 of this transgene resulted in pan-neuronal expression of GCaMP6s
 throughout the nervous system. To identify GFP expressing mice a forward
 (5’-CCT ACG GCG TGC AGT GCT TCA GC-3’) and reverse (5’-CGG CGA GCT GCA CGC
 TGC GTC CTC-3’; 400 band size) PCR amplification was used to identify
 which animals had the GFP transgene. At the end of each recording session,
 the animal was either euthanized or perfused and the brain
 dissected. <strong>Surgical
 procedure</strong> All mice
 were anesthetized with isoflurane (2.5% in pure oxygen) for the procedure.
 Body temperature was maintained at 30C for the duration of the surgery and
 recovery using a feedback-regulated heating pad. Lidocaine (1%) was
 applied subcutaneously in the scalp, followed by careful removal of skin
 above the skull. Ophthalmic ointment was used to protect the eyes during
 the surgery. The cranium was attached to two head bars using
 cyanoacrylate, one across the occipital bone of the skull and the other on
 the lateral parietal bone. After the surgery was complete, mice were
 transferred to a rotating disk for the duration of the recording which
 started 1 hour after the last isoflurane was administered. At the end of
 the recording session, the animal was euthanized and the brain
 dissected. <strong>Recording
 calcium dynamics</strong> In-vivo wide-field fluorescence recordings were
 collected in a minimally invasive manner. Imaging through the skull by
 single-photon excitation light from two blue LED lights (470 nm; Thorlabs
 M470L3) produces a green fluorescent signal that is collected through
 coupled 50mm Nikon lenses (f5.6 / f1.2, optical magnification ∼ 1x) into a
 scientific CMOS camera (PCO Edge 5.5MP, 6.5µm pixel resolution). The top
 lens in the tandem lens system was used to focus on the cortical surface,
 thereby lowering the magnification slightly; anatomical representation for
 each pixel corresponded to 6.9±0.2µm (min: 6.7µm, max: 7.2µm). Excitation
 light was filtered with a 480/30 nm bandpass (Chroma Technology AT480/30x)
 and the emission signal was filtered with 520/36 nm bandpass (Edmund
 Optics 67-044). Data collection was performed in a dark, quiet room with
 minimal changes in ambient light or sound, thus the brain activity
 recorded was at resting state without direct stimulation. Raw data was
 written directly as a set of 16-bit multi-image TIFF files. Concurrent
 body camera recording was taken using infrared (IR) raspberry pi cameras
 (Raspberry Pi NoIR V2), recorded at 30 fps synced to the CMOS.
 The total amount of data recorded for each
 animal was generally at least 40 min and the amount of time in between
 video segments was less than 1 minute. All analyzed data consisted of two
 sequential full spatial resolution recordings concatenated together giving
 20 min of data sampled at 10 frames per second for each video
 decomposition. All video segments
 consisted of a set of continuously collected images at 10 or 20 frames per
 second for 10 minutes. The total amount of recorded data for each animal
 was generally at least 40 min and the amount of time in between video
 segments was less than 1 minute. When more than 10 minutes of video data
 were used for single decompositions, multiple videos were concatenated
 together. All analyses were conducted on data recorded at 10Hz, except
 exploration of effects of resolution on data quality. Spatial resolution analyses were performed on a
 single 10 minute recording at 10Hz. Spatial down sampling was conducted by
 taking the mean between groups of pixels in a dxd square, where d is the
 integer down sampling factor. Temporal resolution analyses were performed
 on a single 10 minute recording at 20Hz. This data was temporally down
 sampled by taking the mean between an integer number of subsequent
 frames. <strong>ICA
 decomposition and saving</strong> ICA was performed using FastICA, implemented
 through python’s <em>sklearn</em> decomposition. The ICA
 decomposition was applied to the spatially flattened (xy,t) 2-D
 representation of the video data under the cortical ROI mask. The mean
 time series is pre-subtracted from the array before SVD decomposition or
 ICA decomposition, since ICA cannot separate sources with a mean signal
 effect. The filtered, unfiltered mean, ICA components, mixing matrix, and
 associated metadata are all saved. Data is stored and saved in this
 flattened format for storage optimization. Components are locally
 spatially reconstructed for visualization in the GUI. Requesting the full number of components resulted
 in extremely lengthy processing times. To reduce the processing time, the
 data was preprocessed through Singular Value Decomposition (SVD)
 whitening, and noise components were cropped. To ensure that no signal was
 lost, and there were ample dimensions left for ICA separation, the
 inflection point between SVD signal and noise floor was identified, and
 SVD components were reduced to include components equal to 5 times the SVD
 signal to noise cutoff value. This cutoff multiplier can be adjusted while
 ICA projecting. After calculating
 and sorting the ICA results, excessive noise components are removed from
 the dataset for compression. The cutoff was determined by identifying the
 inflection point in the lag-1 autocorrelation distribution with a
 two-peaked KDE fit. Components were saved such that 75% of the components
 saved were signal or artifact, and 25% of the components saved were noise
 associated. If not enough noise components were returned by the ICA
 decomposition, there is a risk that signals were not sufficiently unmixed,
 so the ICA decomposition was repeated with a higher SVD cutoff until
 enough additional noise components were included. This was not necessary
 for any of our experiments in this paper, but can be useful in ensuring
 signal capture in recordings of smaller areas, such as cranial
 windows. For resolution down
 sampling analysis, the number of components requested was kept constant
 for all resolution values to avoid confounding the results. The cutoff for
 the highest resolution tested was used for all lower resolution
 decompositions. ICA returns
 components that are unsorted and often spatiotemporally inverted.
 Components were first sorted by their time series variance. The spatial
 histogram of sIC values across each component can be visualized as a
 single-tailed gaussian distribution centered around 0, where the tail
 represents the spatial domain of each component. The two edges of the
 distribution are first identified. The boundary closer to 0 is taken as
 the edge of the central noise distribution, and that boundary is used to
 define the dynamic threshold on the opposite side of 0. Any values outside
 of this noise distribution and is part of the wider tail are included in
 the binarized domain of the component. If the tail was negative, the
 component was inverted spatially and temporally by multiplying each by a
 factor of -1. In this way, components were all identified as positive
 effectors for visualization purposes. Movie rebuilding is not affected by
 this process.
description: All code used in this paper is available at
 github.com/ackmanlab/pyseas or as a package (pySEAS: python Signal
 Extraction and Segmentation) in the python package index (pip install
 seas). In pySEAS, we have built a graphical user
 interface (GUI) to explore the data (see Documentation or Github for
 example scripts). We have also built hdf5manager
 package associated with pySEAS that uses h5py packages, however the user
 can use their own HDF5 reader if interested.