# Endfoot white paper
The data used in the endfeet project are from two “sources”. The first group of data is the data from the sleep project with GCaMP6f, the second group of data is about vessel diameters with recordings of the endfeet via GLT-GFP and lumen via texas red.

The data used directly from the sleep project without edits are:
* Arteriole, venous and capillary ROIs.
* The EMG and the normalized ECoG at 512 Hz for GCaMP.
* The episodes and state definitions. In the case with the GLT the StateProcessor class from the sleep project was directly used to define the wake states.

# Methods: Code

## Camera

ROIs were placed on the snout and wheel on the infrared surveillance camera to detect whisking and locomotion. A trace of the average absolute difference between each frame was created for each ROI. This was calculated by 1) converting each frame to a grayscale image and 2) then subtracting each frame by the previous frame. 3) Subsequently the absolute value of each frame was calculated and 4) the average value within each ROI was calculated.

As there always is some level of noise or movement in the camera the traces displayed baseline level of activity. In an attempt to normalize the traces we estimated the baseline level of activity in each trace and subtracted it. The baseline was estimated by rounding the traces to 3 decimal places and calculate the mode. In trials where the mouse was mostly running the baseline could erroneously be estimated to the level representing locomotion, to avoid this values above 5 were excluded from the baseline estimate. The activity after subtracting the baseline seemed to be normalized enough that further normalizing of the variance was not necessary.

## Electrophysiology

The EMG and ECoG are originally saved in csv files in the trial folders at 10k Hz. The csv file contains both of the channels, but which channel contains EMG and ECoG differs between experiment days. The channel assignment is written in the Excel sheets that also contain the manual sleep scoring. Because of this the Excel sheets are parsed before the electrophysiology is loaded.

With the sleep project data the previously extracted EMG and normalized ECoG at 512 Hz is used. For the GLT dataset the electrophysiology is loaded from the original files and resampled to 512 Hz. The ECoG from the sleep project was normalized to the average total power in the 0.5–30 Hz frequency range during NREM sleep per mouse and average power was calculated using the MATLAB bandpower() function. This was done by dividing the traces by the square root of the bandpower. The normalization in the sleep project can be found at dogbane.processing_other.normalized_eeg.

When we plot the EMG we usually subtract the median value of the trial and then use the absolute value. When plotting ECoG total power we first filter the normalized ECoG traces with a 0.5 to 30 Hz zero phase two-pass bandpass filter of order 2 using the butter and the filtfilt functions. We then apply the hilbert transform and use the square absolute value to estimate the instantaneous power. The traces are lastely smoothed with a 0.2 second gaussian kernel.

## Sleep-wake state scoring

The wake states were defined with the same algorithm as in the sleep project. The following steps are done to define wake states: **1**) Locomotion is defined as periods above 4 a.u. on the wheel camera trace. **2**) Gaps in locomotion periods less than 5 seconds apart are defined as locomotion. **3**) Whisking is defined as periods above 1.5 a.u. on the whisking camera trace. **4**) Gaps in whisking periods less than 5 seconds apart are defined as whisking. **5**) Locomotion periods shorter than 5 seconds are removed. These periods will then rather be defined as whisking. **6**) Whisking is redefined as periods that have whisking but do not have locomotion. This is to prevent defining a period as both whisking and locomotion. **7**) Whisking periods that are shorter than 3 seconds are defined as twitching. **8**) Twitching is ignored. This makes it possible to have short whisking in quiet wakefulness. **9**) Quiet wakefulness is defined as periods that do not have locomotion, whisking or twitching. **10**) The last 3 seconds of each quiet wakefulness episode is removed. **11**) Quiet wakefulness periods less than 5 seconds is removed. **12**) The 10 first seconds of each recording is ignored. Because this is done after step 11 there are some quiet episodes that are shorter than 5 seconds.

Sleep states were identified from filtered ECoG (0.5–30 Hz) and EMG signals (100–1000 Hz) based on standard criteria for rodent sleep17,34 (Supplementary Fig. 1b). NREM sleep was defined as high amplitude delta (0.5–4 Hz) ECoG activity and low EMG activity; IS was defined as an increase in theta (5–9 Hz) and sigma (9–16Hz) ECoG activity, and a concomitant decrease in delta ECoG activity; REM sleep was defined as low amplitude theta ECoG activity with theta/delta ratio >0.5 and low EMG activity. Wakefulness states were identified using the infrared-sensitive surveillance camera video by drawing ROIs over the running wheel and mouse snout (Supplementary Fig. 1a). The signal in the wheel and snout ROIs was quantified by calculating the mean absolute pixel difference between consecutive frames in the respective ROIs. Voluntary locomotion was identified as signals above a threshold in the wheel ROI. Spontaneous whisking was defined in the snout ROI. Quiet wakefulness was defined as wakefulness with no signal above threshold in both ROIs.

## Sleep-wake state transition scoring
Sleep-to-wake transitions were determined as previously described. For the transition from NREM and IS to wakefulness, onset of wakefulness was determined by the first 16 sign of ECoG desynchronization (activation) . During the transitions from REM to wakefulness, end of REM was identified by the interruption of sustained theta waves and the onset of desynchronized ECoG.

## Vasomotion recording
Vessel diameters were measured using a custom tool which used the GLT1-eGFP fluorescence channel to measure the endfoot tube diameter and Texas Red-labeled dextran channel to measure the lumen diameter (Supplementary video 1). The diameter of the endfoot tube was measured along a manually drawn line crossing the center of the vessel. Moving from the center towards the periphery, width of the vessel was defined as the distance between the first index of half maximum fluorescence at either side of the center of the vessel. The vessel lumen diameter was similarly measured at the downstroke at half maximum fluorescence moving from center towards to periphery. For each trial we created traces of the endfoot tube and vessel diameters in micrometers. To reduce noise, we performed a 1 s window averaging of the imaging time series. Distance was sampled at 2 Hz. The size of the PVS was calculated by subtracting the lumen diameter from the endfoot tube diameter. The vessel diameter was ignored at timepoints where the vessel tool was unable to estimate the diameter. Where required, the missing values were replaced by linear interpolation. Vessel traces were inspected and artifacts were manually removed. The missing values were replaced with linear interpolation.

## Vasomotion calculation
Arterioles and venules were identified by their morphology and direction of blood flow. Vessel baseline for Figs. 1 and 2, and Supplementary Figs. 2 and 3 was defined from periods of quiet wakefulness. Separately for the endfoot tube, lumen and perivascular trace the baseline diameters were defined as the median of all samples within the baseline periods. The baseline diameters were subtracted from the continuous diameter traces. The vessel diameter change (Fig. 1, Supplementay Fig. 2) in each wake-sleep state was then calculated by taking the median diameter within each episode.

## Vasomotion oscillation power
For Fig 2. vasomotion spectrograms were created with short-time Fourier transform using the MATLAB function spectrogram() with a window of 30 seconds, overlap equal to the window minus 1 sample and number of discreet Fourier transform samples equal to the size of the window. Within each sleep-wake state episode the average power in the 0.1–0.3 Hz range was calculated with the MATLAB bandpower() function on the diameter traces. The 0.1–0.3 Hz frequency range was chosen to best illustrate the slow oscillations we observed in the spectrograms. To create instantaneous power traces the vessel traces were first filtered the traces using a 0.1–0.3 Hz zero phase two-pass bandpass filter of order 2 using the butter() and the filtfilt() MATLAB functions. Subsequently we applied the Hilbert transform and calculated the square absolute value.

## Ca2+ event detection
We manually defined ROIs over endfoot compartments abutting arterioles (A-EF), capillaries (C-EF) and venules (V-EF). Arterioles and venules were identified by their morphology and direction of blood flow, whereas capillaries were defined as vessels smaller than 8 μm in diameter. A-EF and V-EF ROIs were manually moved with a custom made tool to follow the vasomotion of the vessels. We then used our event-based ROA algorithm16 to detect Ca2+ events and calculate their frequency, density, size and duration across the sleep-states.

## Ca2+ signaling and vasomotion at sleep-wake state transitions
DF/F0 traces from endfoot ROIs and vasomotion traces were aligned to the start of episodes. As there were small differences in the sampling frequency from the two-photon microscope the vasomotion traces were resampled to exactly 2Hz and the ROI DF/F0 traces to 30 Hz. For ROI traces, average fluorescence in the period -10 to -5 s was used to define F0. Additionally the ROI traces were smoothed with a kernel on the form t ⋅ exp(-t / τ) (with τ = 0.25 s) that rises quickly, and tapers off with an exponential decay, that smoothed the traces without skewing the activity to earlier time points. The average dilation between -10 to -5 s was used as a reference and subtracted from each transition.

## Statistical analysis
Statistical analyses were performed in MATLAB using the fitlme() for linear mixed effects regression. As fixed effects, we included the states REM, IS, NREM, quiet wakefulness, whisking and locomotion. The effect of mouse identity on state was set as a random effect. In analyses of average power in the 0.1–0.3 Hz range (Fig. 2) the samples were log transformed before running the analysis to accommodate for a skewed distribution. Episodes shorter than 20 s were removed. Outliers in average power and diameter change were excluded from the analysis. For ROA data (Fig. 3 and supplementary Fig. 3) the interaction between state and ROI group (capillary, arteriole, venule) was included in the statistical models. For the ROA data analyses the samples within identical state and ROI group were averaged per trial and log transformed performing the statistical models.

## Vessel diameter

### Vessel tool
Vessel diameters were measured using Knut’s vessel tool which uses the first channel to measure the endfoot tube diameter and the second channel to measure the lumen diameter. The diameter of the endfoot tube is measured at the upstroke at half maximum fluorescence and the lumen is measured at the downstroke at half maximum fluorescence. The sampling rate was set to 2 Hz. To get a good image of endfoot and lumen 30 images (~1 second) are averaged together for each sample. This means the effective sampling rate is not truly 2 Hz because by averaging 30 frames the true sampling rate would be 1 Hz. The vessel diameter was set to not-a-number (NaN) at timepoints the vessel tool was unable to estimate the diameter. In calculations that required valid samples the missing values were replaced with linear interpolation.

For each trial we created traces of the endfoot tube and vessel diameters in micrometers. For trials that had both lumen and endfoot measurements another trace of the paravascular space was estimated by subtracting the endfoot tube diameter from the lumen diameter.

### Vessel trace correction
GLT trials containing vessel traces with artefacts were labeled with either `ignore` or `artifact`. Trials labeled with `ignore` were completely ignored from the analysis. Trials labeled with `artifact` were inspected and samples deemed outliers were manually replaced by NaN values using a custom built tool.

In transitions and estimations of power the missing values were replaced with linear interpolation.

**Note:** In GCaMP the ignoring and artefact correction of the vessel traces does not function as expected and needs further work.

### Vasomotion

In each trial we manually defined periods that would serve as a baseline measurement we could compare the other measurements with. Separately for the endfoot tube, lumen and perivascular trace the baseline dilation was defined as the median of all samples within the baseline periods. The baseline distance was then subtracted from each of the traces. The vasomotion in each wake-sleep state was then calculated by taking the median within each episode.

### Power estimation
Vasomotion spectrograms were calculated with short-time Fourier transform using the `spectrogram` Matlab function with a window of 30 seconds, overlap equal to the window minus 1 sample and number of DFT samples equal to the window.

The average power in the 0.1-0.3 Hz range was calculated with the Matlab `bandpow` function on the vessel traces within each sleep-wake state episode. The 0.1-0.3 Hz frequency range was chosen from looking at the spectrograms.

To estimate instantaneous power the vessel traces were first filtered with a 0.1 to 0.3 Hz zero phase two-pass bandpass filter of order 2 using the `butter` and the `filtfilt` Matlab functions. We then apply the Hilbert transform to obtain the analytical signals and took the square absolute value.

## ROIs

ROI traces were calculated from averaging the pixels inside each ROI. ROIs was moved with a custom made tool to follow the dilation of the vessels. The baseline (F0) was defined as the mode of the traces and used to create ΔF/F0 traces. The traces were rounded to the 3rd decimal point before calculating the baseline, without any rounding the mode can give wrong results. The traces were not resampled to 30Hz before events were detected as they were in the sleep project. This is because I noticed the number of detected events were different if the traces were resampled and I thought it would be better to use the original traces.

ΔF/F0 from each ROI was filtered using a gaussian filter (σ = 0.25s). The smoothed trace (gaussian filter, σ = 120 s) was subtracted from the ΔF/F0 to remove drift. Noise traces were approximated by subtracting the smoothed trace from the raw ΔF/F0 trace. Ca2+ events were detected as increases in ΔF/F0 larger than 2.5 times of the GCaMP6f noise trace. The duration of the signal was defined as the time points where the ΔF/F0 crossed 0 and the amplitude was the peak of the ΔF/F0 within the duration. Events starting before 1 second were ignored.

The ROI events are detected on the filtered version of the traces, these traces are saved alongside the events. When we plot ROI traces we usually use the filtered traces.

### TODO: Region-of-activity (ROA) algorithm

From the sleep paper:

ROA algorithm consists of the following steps (see Supplementary Fig. 3): (A) Preprocessing. Imaging data was corrected for movement artifacts and smoothed in the spatial domain (gaussian smoothing, σ = 2 pixels) producing the time series F. (B) Calculating ΔF/F0 time series. A baseline image (F0) was calculated by smoothing the pre-processed time series (F) in time (moving average filter, width 1.0 second), resulting in a lowpass filtered time series (FLP), and calculating the mode of each of the pixels over time. The pre-processed time series (F) was then subtracted and divided by the baseline image (ΔF/F0 = (F - F0) / F0), resulting in a ΔF/F0 time series (S). (C) Calculating a noise-based activity threshold and thresholding the data. A moving average filter (width, w = 1.0 second) was then applied to the ΔF/F0 time series to create a smoothed, lowpass filtered time series (SLP). A highpass filtered time series (SHP) was then created by subtracting the smoothed time series (SLP) from the original ΔF/F0 time series (S). The highpass filtering of S was done to estimate noise in our time series. To estimate the variance of the noise in S we assumed it could be approximated by variance of SHP. As SLP is a moving average filtered version of S we would then expect that the standard deviation of the noise in SLP is a factor of 1/〖√L〗_filtlower than S, where Lfilt is the length (number of frames) of the moving average filter. A standard deviation image (σ) (the noise approximation of SLP) was created by calculating the standard deviation of SHP and diving by 〖√L〗_filt. Artifacts from fluorescence drifts20 were removed by subtracting a 10 s moving average of SLP from SLP,, producing a final, bandpass filtered time series SBP. Voxels in SBP were then thresholded by the corresponding pixels in the σ image and multiplied by a common factor k (k = 5 was used in the analysis), resulting in a 3D matrix of active or inactive voxels. Connected components were then detected and the descriptive properties of the events were extracted. Adjacent active voxels in space and time were assigned to single ROAs. As vessel walls move within the field-of-view with vascular dilation and constriction, the regions over and immediately surrounding blood vessels were prone to artifacts and manually masked out before connecting voxels. For each ROA the starting time, maximal spatial extent (µm²), volume (µm² · s) and duration (s) was recorded.

ROA frequency was calculated by counting all ROAs with starting times within a given frame and subsequently dividing by the sampled area (total field-of-view minus vessel masks) and the time per frame. The percentage of active voxels in a particular state episode was calculated by dividing the number of active voxels by the total number of voxels in that episode, while excluding the ignored areas. 3D renderings of ROA activity (Fig. 2) were made by outlining the ROAs and plotting with the MATLAB function patch(). The time resolution was decimated for improved performance and visual representation.

### TODO: ROA in endfeet ROIs

The ROA algorithm assumes the structure under each pixel does not move. If any part of the FOV moves it is likely that false events are detected and real events go undetected. To detect ROAs in the endfeet, which move as the vessel dilates, we cropped an area around each ROI and adjusted for the movement so the center of the ROIs do not change in time. The area cropped around each ROI was a square with 10 pixel padding around the bounding box of the ROIs. For each ROI the ROA frequency and density trace was calculated.

### TODO: Transitions
Strict vs non strict.

### TODO: Statistical analysis

## Data structure

Processing of the data is done from the main.m script. The main script contains links to other scripts that opens a GUI with buttons that contain various processing and plotting steps. The GUI works by selecting the trials which should be processed then pressing the buttons. As the processing code changes often details are not explained here. Most of the processing uses the save_var and load_var functions to move data.

## Comments

* The moving ROI tool is implemented with a special ROI class not used anywhere else. Saving the ROIs with the roimanager_lite again will remove the moving ROI data.
* The bar plots we have made with the ROIs do not have very good residuals. By cleaning up the data by correcting artifacts or removing outliers we could probably improve it quite much.
* Using a higher sampling frequency with the vessel tool seemed to just make the artifacts higher “resolution”.

# Methods: Animals and procedures

## Animals
Male C57BL/6J (Janvier Labs) and Glt1-eGFP mice were housed on a 12-h light/dark cycle (lights on at 8 AM), 1–4 mice per cage. Each animal underwent surgery at the age of 8–10 weeks, followed by accommodation to being head-restrained and two-photon imaging (2–3 times per week) for up to 2 months. Adequate measures were taken to minimize pain and discomfort. Sample size was based on our previous published studies using similar techniques. No randomization or blinding was performed. All procedures were approved by the Norwegian Food Safety Authority (project number: 11983).

## Cloning and virus production
Serotype 2/1 recombinant adeno-associated virus (rAAV) from plasmid construct pAA V-GF AP-GCaMP6f was generated as previously described and purified by A VB Sepharose affinity chromatography following titration with real-time PCR (rAAV titers about 1.0–6.0 × 10^-12 viral genomes/mL, TaqMan Assay, Applied Biosystems Inc.).

## Surgical procedures and intrinsic imaging
Mice were anesthetized with isoflurane. Two silver wires (200 μm, non-insulated, GoodFellow) were inserted epidurally into 2 burr holes overlying the right parietal hemisphere for ECoG recordings, and two stainless steel wires (50 μm thickness, insulated except 1 mm tip, GoodFellow) were implanted in the nuchal muscles for EMG recordings. The skull over the left hemisphere was thinned for intrinsic signal imaging, a custom-made titanium head-bar was glued to the skull and the implant sealed with a dental cement cap. After two days, representations of individual whiskers in the barrel cortex were mapped by intrinsic optical imaging. The brain region activated by single whisker deflection (10 Hz, 6 s) was identified by increased red light absorption. After two days, chronic window implantation and virus injection was performed as described previously22. A round craniotomy of 2.5 mm diameter was made over the barrel cortex using the intrinsic optical imaging map as a reference. The virus mixture (70 nL at 35 nL/min) was injected at two different locations in the barrel cortex at a depth of 400 μm. A window made of 2 circular coverslips of 2.5 and 3.5 mm was glued together by ultraviolet curing glue53, was then centered in the craniotomy and fastened by dental cement. Mice with implant complications were excluded from the study.

## Behavioral training
Mice were housed in an enriched environment with a freely spinning wheel in their home cages. One week before imaging, mice were habituated to be head-fixed on a freely spinning wheel under the two-photon microscope. Each mouse was trained head-fixed daily before the imaging for increasing durations ranging from 10 min on the first day to 70 min on the last. Mice that showed signs of stress and did not accommodate to head-restraint were not included in the study.

## Two-photon imaging
Four weeks after the surgery, mice were imaged under a laser scanning two-photon microscope (Ultima IV from Bruker/Prairie Technologies) as previously described22, with a Nikon 16 × 0.8 NA water-immersion objective (model CFI75 LWD 16XW). The fluorescence of GCaMP6f was excited at 999 nm with a Spectra-Physics InSight DeepSee laser, and emitted photons were detected with Peltier cooled photomultiplier tubes (model 7422PA-40 by Hamamatsu Photonics K.K.). Images (512 x 512 pixels) were acquired at 30 Hz in layer 2/3 of barrel cortex.

## Head-fixed sleep protocol
To assist sleep in a head-fixed position, we adjusted the running disc position to mimic the body’s natural position observed during unrestrained sleep25. We observed that locking the movement of the wheel once the mouse showed first signs of sleep, such as delta waves in ECoG and eyes closing, facilitated falling asleep. The imaging sessions of sleeping mice started at 9–10 AM (ZT 1–2), the beginning of the light phase, and lasted until 3–6 PM (ZT 7–10). The mice did not have access to food or water while sleeping under the microscope, however, in natural conditions mice feed almost exclusively during the dark phase, ZT 12–2426–28. First signs of drowsiness, such as high delta power in the ECoG signal and eyes closing, were observed 15–45 min after head-fixation, and typically mice spent 90–120 min head-fixed under the microscope before falling asleep. Mice that did not show any signs of sleep within the first 2 hours of head-fixation were removed from the microscope. The mice had an exact replica of the microscope running disc in their home cage, and in our hands this made a large difference in aiding the mice to fall asleep (i.e. initially we tried various types of stages like a spherical treadmill and a tube for immobilization with little success). Mice were not sleep deprived or manipulated in any other way before imaging to induce sleep.

## Behavior and electrophysiology recording
ECoG and EMG signals were recorded using a Multiclamp 700B amplifier with headstage CV-7B, and digitized by Digidata 1440 (both Molecular Devices). Mouse behavior was recorded by an infrared-sensitive surveillance camera. Data acquisition was synchronized by a custom-written LabVIEW software (National Instruments).