Quasi-Periodic Patterns of Neural Activity improve Classification of Alzheimer’s Disease in Mice

Resting state (rs)fMRI allows measurement of brain functional connectivity and has identified default mode (DMN) and task positive (TPN) network disruptions as promising biomarkers for Alzheimer’s disease (AD). Quasi-periodic patterns (QPPs) of neural activity describe recurring spatiotemporal patterns that display DMN with TPN anti-correlation. We reasoned that QPPs could provide new insights into AD network dysfunction and improve disease diagnosis. We therefore used rsfMRI to investigate QPPs in old TG2576 mice, a model of amyloidosis, and age-matched controls. Multiple QPPs were determined and compared across groups. Using linear regression, we removed their contribution from the functional scans and assessed how they reflected functional connectivity. Lastly, we used elastic net regression to determine if QPPs improved disease classification. We present three prominent findings: (1) Compared to controls, TG2576 mice were marked by opposing neural dynamics in which DMN areas were anti-correlated and displayed diminished anti-correlation with the TPN. (2) QPPs reflected lowered DMN functional connectivity in TG2576 mice and revealed significantly decreased DMN-TPN anti-correlations. (3) QPP-derived measures significantly improved classification compared to conventional functional connectivity measures. Altogether, our findings provide insight into the neural dynamics of aberrant network connectivity in AD and indicate that QPPs might serve as a translational diagnostic tool.


Animals -MRI handling procedures
During MRI handling procedures, animals were anesthetized with 2 % isoflurane. Animals were positioned in the magnet by securing their heads with ear bars and fixing incisors over a bite bar.
Ophthalmic ointment was applied to the eyes and a rectal temperature probe was used to monitor animal body core temperature. Body temperature was kept stable at 37 °C via a hot air supply (MR-compatible Small Animal Heating System, SA Instruments, Inc.). A pressure sensitive pad was used to monitor breathing rate and a fiber-optic pulse oximeter, positioned over the tail, to asses heart rate and O 2 saturation (MR-compatible Small Animal Monitoring and Gating system, SA Instruments, Inc.).
Animals then received a 0.3 mg/kg bolus injection of medetomidine (Domitor, Pfizer, Karlsruhe, Germany), after which isoflurane was lowered to 0.4 %. Starting at 15 min post-bolus, a subcutaneous catheter allowed continuous infusion of 0.6 mg/kg/h medetomidine. Functional resting state scans were acquired 30 min post-bolus, lasting 20 min. Great care was taken to keep procedures and conditions identical across animals, with preparatory handling never exceeding 10 min. The choice for a relatively higher anesthesia regimen was based on prior knowledge at our lab, and observations made during pilot tests, that animals would not remain sedated with low anesthesia levels. These observations can be explained given the known differences in neuronal activity and behavior across different mouse strains, and given the increased anxiety levels in multiple mouse models of Alzheimer's disease [1][2][3] .

Identifying quasi-periodic pattern temporal lengths -low contrast QPPs
We used a data-driven processing tool that allows one to determine the most representative window length for a specific type of spatiotemporal pattern. We term this strategy fractional average correlation (FA) and used it to identify the window length of the global-signal like lowcontrast QPPs. This was done independently in both WT and TG animals, indicating for each group an ideal window length at 6 s (Supplementary Fig. S2).
A set of global-signal like low-contrast QPPs was chosen so that each QPP in the set was of a different window length, ranging from 4.5 s to 15 s (at 3 s no global-signal like QPPs were clearly observed). A total of eight QPPs was thus chosen. Each individual QPP in this set was then subdivided into all possible consecutive fractions of a fixed length, specified by the smallest window size investigated. In our analysis the smallest window size was 6 TRs (3 s), meaning that a QPP of e.g. 18 TRs was divided into 13 fractions of 6 consecutive images, each of them shifted by 1 TR (0.5 s). Each of the individual fractions from a QPP at a respective window size were then treated as a reference, and the maximal cross-correlation was calculated with respect to a non-fractioned 'target' QPP at another window size. The average of resultant cross-correlation values from this comparison indicates how many fractions the reference and target QPP had in common, i.e. the FA. FA thus allows one to determine a measure of QPP similarity, which is independent of their spatiotemporal length. The FA-value was then calculated for all combinations of window sizes. All FA-values were represented in an n x n matrix, where n indicates the number of QPPs and each column indicates the FA-values for a respective reference QPP with each respective (target) QPP in the set. The FA-matrix was finally averaged across its columns to obtain the overall FA for each QPP in the set.
This approach is particularly useful because target QPPs at smaller window sizes than the reference QPP under investigation intrinsically have less fractions in common, given that they only represent a subpart of the larger QPP. Comparing larger reference QPPs with non-matching subparts in short target QPPs forces a decrease in the FA-value. In contrast, long target QPPs likely contain the full spatiotemporal pattern present in the reference QPP and will therefore cause a high FA. When FA-values are averaged across columns (i.e. all reference QPP FAvalues with a respective target QPP), this will lead to low average FA-values for QPPs shorter than then ideal window length and to high average FA-values for QPPs longer than the ideal window length. The tipping point of increasing FA, before a plateau is reached, reflects the optimal window size (Supplementary Fig. S2). A gradual decrease in the plateau is the result of lower QPP integrity at longer window sizes, given that these QPPs displayed lower occurrence rates and therefore less image frames were averaged in order to establish the QPP. Elastic net regression is a regularization technique for the estimation of generalized linear models. It operates by combining LASSO regularisation (L1) and ridge regression (L2) penalty terms. These penalty terms constrain the size of estimated predictor coefficients, which can be reduced to zero 5 . This effectively allows only the most relevant predictors to be part of the regression model. Thus, for each model, two hyper parameters need to be set: λ, which determines the two penalty sizes, and α, which provides the relative weight of both penalty terms and therefore determines the balance towards either LASSO or ridge optimization.

Classification analysis -Additional details
In elastic net regression, two risks are represented by: 1) over-fitting models with too many predictors, and 2) overestimating classification accuracy by selecting only those hyper parameters that produce the best results. To avoid these biases, we employed a nested 10-fold cross-validation approach 6 . The inner loop was used to determine the hyper parameters with a grid search approach, while the outer loop was used to evaluate the resultant regression models.
More specifically, this means that in the outer loop, 90% of the data was used as a training set and 10% as a validation set, which was repeated 10 times to ensure that each subject was part of the validation set once. Each training set of the outer loop was then respectively used to perform a grid search on a set of α values and λ values. For each combination of hyper parameters tested on the inner loop, a 10-fold cross-validation was performed to determine the optimal values, which were used to construct the regression model that was evaluated in the outer loop.
To evaluate classification performance of the resultant models, we constructed receiver operating characteristic (ROC) curves and calculated the area under the curve (AUC). The entire crossvalidation procedure described above was repeated 10 times for each model, in order to increase the reliability of the cross-validation error and to determine the mean AUC values 6

Independent component analysis (ICA) and regression
Group ICA was performed on WT and TG groups separately (Supplementary Fig. S3), and on both groups simultaneously (Supplementary Fig. S10), using the GIFT toolbox (v4.0b) 8 . This allowed us to determine RSNs for the respective analyses. For the separate WT and TG analyses, the number of independent components was set to five for both groups. If a higher number was used this led to splitting of bilateral components in the WT group. The presented RSNs of the TG group (Supplementary Fig. S3) also appeared split when three or four components were used, indicating that the observation of unilateral components in TG was not solely the result of running the analysis with five components. For the combined group analysis, the number of components was set to 4, based on the same heuristics as described above. The ICA was run on variance-normalized data, filtered between 0.01-0.2Hz, using the Infomax algorithm with no auto-filling of data reduction values. Stability analysis was performed using the ICASSO algorithm, rerunning the ICA 15 times with a minimal cluster size of 5 and maximal of 10. Using the built-in GIFT functionality, SPM12 (Statistical parametric mapping) was used to obtain contrast images per subject for all components. The latter were used in one sample T-tests for all components to obtain significant group level RSNs (p < 0.05, FDR corrected, threshold 4 voxels).
To evaluate the contribution of RSNs to FC in both groups (Supplementary Fig. S10), the GIFT toolbox was used to regress group-level spatial maps (RSNs) into individual subject level data. WT and TG QPPs were compared across groups, using pSTC (cfr. Fig. 1). Panels display the group mean ratio of each QPP's subject-average occurrences in the reference group over its subject-average occurrences in the target group. QPPs were hereby split into two significance groups and in each they were sorted by increasing ratios: non-significant group differences (twosample T-test, p>0.05) are indicated in grey, while significant ones are color-coded in correspondence with the reference group. The ideology of this procedure is to determine QPPs that are the most dissimilar between groups with regard to their occurrence rate. A) At 3 s window lengths, WT QPPs displayed higher occurrence ratios versus TG, and vice versa. WT QPP ratios versus TG were more often significant compared to TG QPP ratios versus WT. B) At 6 s window lengths, WT QPPs displayed many significant and higher occurrence ratios versus TG. TG QPPs displayed mostly lower and non-significant ratios versus WT. A small set of highcontrast QPPs, similar to those observed in (A), were removed from the 6 s QPPs, based on a 95 % confidence interval of QPP (cfr. Supplementary Fig. 1). C-E) For further investigation, one QPP was selected from each group comparison, from the significant highest occurrence rate ratios (A-B, black arrows). A total of three QPPs were thus selected. Note their nearly identical spatiotemporal similarity to the QPPs presented in Fig. 2. F) Subject mean occurrence rates for all QPPs, based on STC and pSTC (two-sample T-test, **p<0.01, ***p<0.001, error flags show standard error). G-H). Similar as in Fig. 5, the effect of QPP regression was evaluated for DMNand DMN-TPN-like FC. Note that the results are nearly identical. Overall this figure confirms that both selection criteria presented in the current study determine consistent results.