Application of Time-Aware PC algorithm to compute Causal Functional Connectivity in Alzheimer's Disease from fMRI data

Functional Connectivity between brain regions is known to be altered in Alzheimer's disease, and promises to be a biomarker for early diagnosis of the disease. While several approaches for functional connectivity obtain an un-directed network representing stochastic associations (correlations) between brain regions, association does not necessarily imply causation. In contrast, Causal Functional Connectivity is more informative, providing a directed network representing causal relationships between brain regions. In this paper, we obtained the causal functional connectome for the whole brain from recordings of resting-state functional magnetic resonance imaging (rs-fMRI) for subjects from three clinical groups: cognitively normal, mild cognitive impairment, and Alzheimer's disease. We applied the recently developed Time-aware PC (TPC) algorithm to infer the causal functional connectome for the whole brain. TPC supports model-free estimation of whole brain causal functional connectivity based on directed graphical modeling in a time series setting. We then perform an exploratory analysis to identify the causal brain connections between brain regions which have altered strengths between pairs of subject groups, and over the three subject groups, based on edge-wise p-values from statistical tests. We used the altered causal brain connections thus obtained to compile a comprehensive list of brain regions impacted by Alzheimer's disease according to the current data set. The brain regions thus identified are found to be in agreement with literature on brain regions impacted by Alzheimer's disease, published by researchers from clinical/medical institutions.


Introduction
Alzheimer's disease (AD) is the most common age-related progressive neurodegenerative disorder.It typically begins with a preclinical phase and advances through mild cognitive impairment (MCI) to clinically significant AD, which is a form of dementia (Querfurth and LaFerla, 2010).Despite significant efforts to identify biomarkers for AD, it still relies on clinical diagnosis, and early and accurate prediction of the disease remains limited (Laske et al., 2015;Li et al., 2019).Abnormal resting-state functional connectivity (FC) between brain regions has been observed as early as two decades before brain atrophy and the emergence of AD symptoms (Ashraf et al., 2015;Nakamura et al., 2017).Therefore, resting-state FC can potentially determine the relative risk of developing AD (Brier et al., 2014;Sheline and Raichle, 2013).
Resting-state functional magnetic resonance imaging (rs-fMRI) records the bloodoxygen-level-dependent (BOLD) signal from different brain regions while individuals are awake and not engaged in any specific task.The BOLD signal is popularly used to infer functional connectivity between brain regions partly due to the advantage that BOLD signal provides high spatial resolution (Yamasaki et al., 2012;Sporns, 2013;Liu et al., 2015;Xue et al., 2019).
Functional connectivity refers to the stochastic relationship between brain regions with respect to their activity over time.Popularly, functional connectivity involves measuring statistical association between signals from different brain regions.The statistical association measures are either pairwise associations between pairs of brain regions such as Pearson's correlation, or multivariate i.e. incorporating multi-regional interactions such as undirected graphical models (Biswas and Shlizerman, 2022a).Detailed technical explanations of functional connectivity in fMRI can be found in Chen et al. (2017); Keilholz et al. (2017); Scarapicchia et al. (2018).The findings from studies using FC (Wang et al., 2007;Kim et al., 2016), and meta-analyses (Jacobs et al., 2013;Li et al., 2015;Badhwar et al., 2017) indicate a decrease in connectivity in several brain regions in relation to Alzheimer's disease (AD), such as the hippocampus and posterior cingulate cortex.These regions play a role in memory and attentional processing.On the other hand, some studies have found an increase in connectivity within brain regions in early stages of AD and MCI (Gour et al., 2014;Bozzali et al., 2015;Hillary and Grafman, 2017).This is a well known phenomenon, where increase in FC between certain brain regions occurs when the communication between other brain regions is impaired.Such hyperconnectivity has been interpreted as a compensatory mechanism where alternative paths within the brain's network are recruited (Hillary and Grafman, 2017;Marek and Dosenbach, 2022;Oldham and Fornito, 2019).
In contrast to associative FC, causal FC represents functional connectivity between brain regions more informatively by a directed graph, with nodes as the brain regions, directed edges between nodes indicating causal relationships between the brain regions, and weights of the directed edges quantifying the strength of the corresponding causal relationship (Spirtes et al., 2000).However, functional connectomics studies in general, and in relation to fMRI from Alzheimer's disease in particular, have predominantly used associative measures of FC (Reid et al., 2019).There are a few studies focusing on alterations in CFC in relation to Alzheimer's disease (Rytsar et al., 2011;Khatri et al., 2021), however this area is largely unexplored.This is partly due to the lack of methods that can infer the CFC in a desirable manner as explained next.
Several properties are desirable in the context of causal modeling of functional connectivity (Biswas and Shlizerman, 2022a;Smith et al., 2011).Specifically, the CFC should represent causality while free of limiting assumptions such as linearity of interactions.In addition, since the activity of brain regions are related over time, such temporal relationships should be incorporated in defining causal relationships in neural activity.The estimation of CFC should be computationally feasible for the whole brain functional connectivity, instead of limiting to a smaller brain network.It is also desirable to capture beyond-pairwise multi-regional cause and effect interactions between brain regions.Furthermore, since the BOLD signal occurs and is sampled at a temporal resolution that is far slower than the neuronal activity, thereby causal effects often appear as contemporaneous (Granger, 1969;Smith et al., 2011).Therefore, the causal model in fMRI data should support contemporaneous interactions between brain regions.
Among the methods for finding CFC, Dynamic Causal Model (DCM) requires a mechanistic biological model and compares different model hypotheses based on evidence from data, and is unsuitable for estimating the CFC of the whole brain (Friston et al., 2003;Smith et al., 2011).On the other hand, Granger Causality typically assumes a vector auto-regressive linear model for activity of brain regions over time, and it tells whether a regions's past is predictive of another's future (Granger, 2001).Furthermore, GC does not include contemporaneous interactions.This is a drawback since fMRI data often consists of contemporaneous interactions (Smith et al., 2011).In contrast, Directed Graphical Modeling (DGM) has the advantage that it does not require the specification of a parametric equation of the neural activity over time, it is predictive of the consequence of interventions, and supports estimation of whole brain CFC.Furthermore, the approach inherently goes beyond pairwise interactions to include multiregional interactions between brain regions, along with estimation of the cause and effect of such interactions.The Time-aware PC (TPC) algorithm is a recent method for computing the CFC based on DGM in a time series setting (Biswas and Shlizerman, 2022b).In addition, TPC also incorporates contemporaneous interactions among brain regions.A detailed comparative analysis of approaches to find causal functional connectivity is provided in Biswas and Shlizerman (2022a,b).With the development of methodologies such as Time-aware PC, it would be possible to infer the whole brain CFC with the aforementioned desirable properties.
In this paper, we apply the TPC algorithm to infer the causal functional connectivity between brain regions from resting-state fMRI data.By applying the algorithm to the fMRI of subjects, we estimate the subject-specific CFC for all subjects in the dataset.It is noteworthy that different subjects are in different clinical categories: Cognitively Normal (CN), Mild Cognitive Impairment (MCI), and Alzheimer's Disease (AD).To identify which causal connections differ between brain regions across pairs of clinical categories, we utilize Welch's t-test, comparing the weights of causal functional connections of subjects for a pair of clinical categories.This analysis reveals the p-values of causal links between brain regions to exhibit differences between subjects in distinct clinical categories, such as with cognitively normal vs. Alzheimer's disease.Additionally, we employ a Kruskal-Wallis H-test, which is a non-parametric version of ANOVA test, to find p-values of causal links to exhibit differences across subjects of the three clinical categories.These links provide insights into the causal connectivity connections that are relevant in exhibiting differences among the three clinical categories.We then compile a comprehensive list of brain regions impacted by Alzheimer's disease based on the altered causal links obtained from the current dataset.Notably, the obtained brain regions are consistent with existing literature, with each such publication being a report from a team involving a clinical setting and at least one medical expert, thereby validating the approach.
Table 1 includes a summary of the participants' demographic and medical information.
2.2.Image Acquisition.The acquisition of fMRI images was performed using a Siemens Magnetom Allegra scanner.The fMRI images were obtained using an echo planar imaging (EPI) sequence at a field strength of 3.0 Tesla, with a repetition time (TR) of 2.08 seconds, an echo time (TE) of 30 milliseconds, and a flip angle of 70 degrees.
The matrix size was 64 × 64 pixels, there were 32 axial slices parallel to AC-PC plane, in plane resolution was 3 × 3 mm 2 , slice thickness was 2.5 mm.Resting scans lasted for 7 mins and 20 secs for a total of 220 volumes during which subjects were instructed to keep their eyes closed, to not think of anything in particular and to refrain from falling asleep.
2.3.fMRI Preprocessing.The fMRI pre-processing steps were carried out using the CONN toolbox version 21a, which utilizes the Statistical Parametric Mapping (SPM12), both of which are MATLAB-based cross-platform software (Nieto-Castanon, 2021; Friston et al., 1994).We used the default preprocessing pipeline in CONN, consisting of the following steps in order: functional realignment and unwarp (subject motion estimation and correction), functional centring to (0,0,0) coordinates (translation), slice-time correction with interleaved slice order, outlier identification using Artifact Detection and Removal Tool (ART), segmentation into gray matter, white matter and cerebrospinal fluid tissue, and direct normalization into standard Montreal Neurological Institute (MNI) brain space, and lastly, smoothing using spatial convolution with a Gaussian kernel of 8mm full width half maximum.This pipeline was followed by detrending, and bandpass filtering (0.001-0.1 Hz) to remove low-frequency scanner drift and physiological noise of the fMRI images.The first four time points have been filtered out to remove any artifact.
For the extraction of Regions-Of-Interest (ROIs), the automated anatomical labeling (AAL) atlas was utilized on the preprocessed rs-fMRI dataset (Tzourio-Mazoyer et al., 2002).The list of all regions in AAL atlas is provided in Appendix A along with their abbreviated, short, and full region names.This specific parcellation method has been demonstrated to be optimal for studying the functional connectivity between brain regions (Arslan et al., 2018).The voxels within each ROI were averaged, resulting in a time series for each ROI.
2.4.Inference of causal functional connectivity: Time-aware PC algorithm.The Time-aware PC (TPC) Algorithm finds causal functional connectivity between brain regions from time series based on Directed Graphical Models (DGM) (Spirtes et al., 2000;Pearl, 2009;Biswas and Shlizerman, 2022a,b;Biswas and Mukherjee, 2022).While traditional DGM is applicable to static data, TPC extends the applicability of DGM to CFC inference in time series by firstly implementing the Directed Markov Property (DMP) to model causal spatial and temporal interactions in the time series by an unrolled Directed Acyclic Graph (DAG) of the time series.The unrolled DAG consists of nodes (v, t), for region of interest v and time t, and edge (v 1 , t 1 ) → (v 2 , t 2 ) reflecting causal interaction from the BOLD signal in region v 1 at time t 1 to the BOLD signal in region v 2 at time t 2 .The estimation of the unrolled DAG is carried out by first transforming the time series into sequential variables with a maximum time delay of interaction τ , and then applying the Peter-Clark (PC) algorithm to infer the unrolled DAG based on the sequential variables (Kalisch and Bühlman, 2007).TPC then rolls the DAG back to obtain the CFC graph between the regions of interest (see Figure 1) (Biswas and Shlizerman, 2022b).We consider τ = 1 for our analyses, which would include interactions of the BOLD signal between regions of interest with a maximum time delay of 2.08 s, the TR of the fMRI acquisition.
The CFC outcome of this methodology is interpretable in the following manner: An edge from region i → j in the CFC estimate represents significant causal interaction from brain region i at preceding times to region j at following times.The model and the approach are non-parametric, meaning that it does not require the specification of a parametric dynamical equation for neural activity.The method captures beyondpairwise multivariate interactions between brain regions.It also supports estimation of the CFC for the whole brain in a computationally feasible manner.It also allows for time delays in interactions between the brain units as well as the presence of feedback-loops.Furthermore, it has been shown that if the neural activity obeys an arbitrary dynamical process, the model outcome of TPC is consistent with respect to the causal relationships implied by the dynamical process and is predictive of counterfactual queries such as ablation or modulation (Biswas and Shlizerman, 2022b).

Rolled CFC Between Regions
Unrolled DAG for Time Series . . .
Steps conveying the concept of the TPC algorithm to infer CFC from observed neural time series data: First the neural time series is transformed to form sequential samples with a maximum time delay of interaction, τ (here τ = 1).Then, Peter-Clark (PC) algorithm is applied on the sequential samples to obtain the unrolled DAG satisfying the Directed Markov Property.Finally the unrolled DAG is transformed to obtain the Rolled CFC between regions.
It is noteworthy that implementing the DMP on the unrolled DAG to model causal relationships over time enables contemporaneous interactions e.g. from region u to region v at time t (Biswas and Shlizerman, 2022b).Such contemporaneous interactions are represented by the edge (u, t) → (v, t) in the unrolled DAG, and presence of such an edge in the unrolled DAG would be reflected as an edge u → v in the Rolled CFC outcome.Such contemporaneous interactions are especially relevant in fMRI due to the relatively slow temporal resolution of the BOLD signal compared to the underlying neural activity (Smith et al., 2011).
2.5.Group-wise Comparisons of Estimated CFC.Using the subject-specific CFC estimated by TPC algorithm, we perform further statistical tests.We use Welch's t-test to obtain p-values of connections to exhibit greater or less weight in one disease stage compared to another disease stage (Yuen, 1974).We use the Kruskal-Wallis H-test, which is a non-parametric version of the ANOVA test, to obtain p-values of connections to exhibit unequal weight in either of the four disease stages (Kruskal and Wallis, 1952).

3.1.
Subject-specific Causal Functional Connectivity.Figure 2 shows the causal functional connectivity (CFC) estimated using TPC algorithm for an example subject in the cognitively normal group.In Figure 2-a, the CFC is represented in the form of a matrix, whose entry (i, j) indicates the presence of connectivity from region index i → j, and the value at entry (i, j) represents the weight of that causal connection.A positive value (blue) of the weight is indicative of excitatory influence whereas a negative 3.4.Brain Regions with altered connections.In Table 2, we list 9 brain regions (and 5 additional regions) that correspond to altered causal functional connections to or from them between subjects of CN, MCI and AD groups with edge-wise p-value less than 0.05 (less than 0.1).The subject-specific causal functional connectomes have been estimated using TPC algorithm (see Section 3.2).The brain regions found are in agreement with existing publications cited in Table 2-right column.

Discussion
In this study, we have obtained the causal functional connectivity of the whole brain from resting state fMRI time series.We used the recently developed Time-aware PC (TPC) algorithm based on directed graphical modeling in time series, to compute the causal functional connectivity.In the dataset, the subjects belonged to three clinical categories: cognitively normal (CN), mild cognitive impairment (MCI) and Alzheimer's disease (AD).We performed group-wise comparisons of the subject-specific causal functional connectivity to identify which causal functional connections are altered between pairs of subject groups.The altered causal functional connections between CN and AD were used to obtain brain regions involved with such altered connections in AD.This resulted in the identification of 12 brain regions where causal functional connections to or from those regions are altered in Alzheimer's disease with p-value less than 0.05, and 5 additional regions with p-value less than 0.1.
It is noteworthy that while several studies have concluded decreased connectivity in MCI and AD compared to CN (Jacobs et al., 2013;Li et al., 2015;Badhwar et al., 2017), prominent researchers have highlighted that MCI and early stages of AD can involve an increase in functional connectivity between brain regions.This increase occurs when the communication between specific brain regions is impaired.It has been interpreted as a compensatory mechanism where alternative paths within the brain's network are recruited (Hillary and Grafman, 2017;Marek and Dosenbach, 2022;Oldham and Fornito, 2019).This explains the presence of causal functional connections estimated by TPC algorithm, whose weight in AD is greater compared to CN in addition to connections with weight in AD less than that in CN (see Figure 3).The following are causal functional connections found by TPC, that have edge-wise p-value less than 0.05 for strength in AD greater than that in CN, and are lowest 5 in p-value: Lobule IV, V of cerebellar hemisphere Left → Lobule IV, V of cerebellar hemisphere Right; Superior frontal gyrus, dorsolateral Left → Superior frontal gyrus, medial Left; Middle occipital gyrus Left → Middle occipital gyrus Right; Middle temporal gyrus Left → Middle temporal gyrus Right; Heschl's gyrus Right → Superior temporal gyrus Right.
In the short term, the augmentation of functional connectivity along alternative pathways exhibits efficiency and adaptability of the brain.However, it is imperative to acknowledge the susceptibility of these densely interconnected hubs to beta-amyloid deposition, which can elicit secondary damage through metabolic stress, ultimately culminating in system breakdown (Hillary and Grafman, 2017).Consequently, the initial state of hyperconnectivity observed in neurodegenerative disorders may gradually transition into hypoconnectivity among the engaged pathways, thereby contributing to cognitive decline as the disease advances (Marek and Dosenbach, 2022).
Based on the causal functional connectome outcome alone, this study has been able to identify many brain regions related to Alzheimer's disease, which have been reported across more than 30 different studies, using different feature extraction methods and advanced imaging technologies.Therefore, this study demonstrates the promise of a causal functional connectivity approach based on directed graphical models in time series and estimated by TPC algorithm.
Although most of the regions linked to AD have been identified, certain regions such as the Thalamus ( Štepán-Buksakowska et al., 2014), that have been also linked with AD has not been identified by the methodology.This can be due to lack of significance given the low sample size or choice of significance level.However, the study demonstrates the potential in using the methodology to a larger dataset.Using a larger dataset, the causal functional connectomics methodology can be used for 1) computation of a subject's causal functional connectome from their fMRI data, 2) identification of specific connections as biomarkers for Alzheimer's disease, and 3) using biomarker connections for the early prognosis and diagnosis of Alzheimer's disease.
In this paper, we have demonstrated the following: (a) Application of the TPC algorithm to compute whole-brain CFC for each subject, (b) Interpretation of CFC in the context of AD using domain (neuropathological) knowledge, and (c) Exploratory analysis for edge-wise differences and corresponding brain regions with altered connectivity in subjects from pairs of clinical categories (CN, MCI and AD), and among the three clinical categories.The findings are consistent with published medical literature.In summary, our results show the promise of computing the whole-brain CFC from fMRI data using the TPC algorithm to gain prognostic and diagnostic insights.

Figure 2 .
Figure 2. CFC outcome of an example subject who is cognitively normal.The CFC is obtained by TPC algorithm.(a) The CFC is represented by a matrix, whose entry (i, j) represents the connection of region i → j. (b)The CFC is visualized with graph edges on the Frontal, Axial and Lateral brain maps (left to right).The nodes correspond to brain region centers, ranging from superficial (light gray) to deeper (darker gray) regions, in the AAL brain atlas.

Figure 3 .Figure 3 .Figure 4 .
Figure 3. Causal functional connections with p-values for altered weights in pariwise comparisons between groups less than 0.05 (green) and less than 0.1 (red), as obtained by Welch's t-test.(a1)-(a3).The connections are represented in matrix format with a non-zero entry in (i, j) corresponding to the edge i → j. (b1)-(b3).The connections are represented by graph edges on the brain schematics.

Table 1 .
Subject demographic information summary

Table 2 .
Brain regions corresponding to altered causal functional connections, whose strength is altered between CN, MCI and AD subject groups with p-values less than 0.05 (green) and additionally with p-values less than 0.1 (red) based on Kruskal-Wallis H-test.The causal functional connections are obtained by TPC algorithm.