Recovering directional brain networks under temporal undersampling: an application to schizophrenia

Technical Summary: Recovering Directional Brain Networks under Temporal Undersampling

Problem Statement
Functional connectivity research in schizophrenia has established a robust "thalamocortical signature": reduced coupling between the thalamus and prefrontal cortex, and increased coupling between the thalamus and sensory/motor cortices. However, standard functional connectivity measures are undirected and derived from Blood-Oxygenation-Level-Dependent (BOLD) fMRI signals. A fundamental limitation of fMRI is temporal undersampling: the BOLD signal is sampled at intervals (TR) of seconds, which is one to two orders of magnitude slower than the underlying neural dynamics. This mismatch causes causal graphs estimated at the measured timescale to be low-pass filtered images of the true process. Consequently, undersampling can hide true edges, fabricate spurious ones, and, most critically for directional claims, reverse the apparent orientation of causal relationships. Existing methods that treat the measured BOLD timescale as the causal timescale implicitly inherit these distortions, making directed interpretations potentially unreliable.

Methodology: RnR (Reason & Refine)
The authors introduce RnR, a causal discovery framework designed to estimate directed causal graphs from fMRI data while explicitly accounting for unknown temporal undersampling rates. RnR builds upon the Rate-Agnostic Structure Learning (RASL) framework and improves upon the solver-based RASL (sRASL) with five key innovations tailored for noisy, real-world fMRI data:

1. Equivalence Class Recovery: Instead of returning a single best-fit graph, RnR recovers an equivalence class of plausible directed graphs $(G, u)$, where $G$ is a causal-timescale graph and $u$ is the sampling rate that would produce the observed measurement graph $H$. It returns all candidates within a relative cost tolerance ($\Delta = 1.9 c_{min}$) of the optimal solution.
2. Biological Plausibility Constraints: The method imposes a target edge density (10–30%) to ensure solutions remain within biologically realistic bounds, penalizing candidates that deviate from this range.
3. Hierarchical Cost Optimization: The cost function is split into three components minimized in a fixed order: density error, bidirected-structure error (representing hidden common causes), and directed-orientation error. This moves from low-information to high-information constraints.
4. Evidence-Weighted Recovery: Edge weights are derived from the underlying first-order estimator (PCMCI in this study). The solver assigns higher costs to dropping strongly supported links and lower costs to dropping weak ones, allowing the cost function to select the most probable pre-image (e.g., distinguishing a bidirected edge from a mutual 2-cycle).
5. Meta-Solver Capability: RnR acts as a meta-solver on top of existing estimators (like PCMCI), offering the option to replace mutual 2-cycles (the footprint of confounders in first-order estimators) with bidirected edges, thereby folding undersampling awareness into the upstream output.

The workflow involves estimating a measurement graph $H$ using PCMCI, converting it to a compressed causal graph, and solving the inverse problem using Answer Set Programming (ASP) to find the set of $(G, u)$ pairs.

Key Contributions

• Directional Recovery: The method successfully recovers directed causal graphs from resting-state fMRI, moving beyond undirected correlation to infer causal direction.
• Stability Classification: It classifies edge orientations based on their stability across inferred sampling rates. An edge is labeled "undersampling-robust" if its direction is invariant across all voting rates, and "undersampling-fragile" if the direction flips decisively with the rate.
• Resolution-Dependent Undersampling: The paper demonstrates through simulation that coarse spatial resolution (aggregating neurons into regions/voxels) artificially lowers the perceived undersampling rate. This explains why single-timescale methods have historically succeeded in recovering undirected signatures but predicts that they will fail as spatial resolution increases toward the neural scale.

Results
The method was applied to the multi-site FBIRN cohort (311 participants: 151 healthy controls, 160 schizophrenia patients) using 10 spatially independent components covering seven functional domains.

1. Synthetic Validation: In simulations with known ground truth, RnR significantly outperformed five single-timescale baselines (PCMCI, FASK, GIMME, MVAR, MVGC). RnR improved F1 scores by an average of 64% and maintained high recall for group differences even as the sampling interval lengthened, whereas baselines suffered systematic bias and performance degradation.
2. Schizophrenia Connectome: RnR recovered 14 FDR-significant directed edges differentiating schizophrenia patients from controls, compared to $\le 3$ for any other method.
   • Primary Finding: The method reproduced the field's most replicated finding—the sensory-to-visual hyperconnectivity hub—but with directional specificity. It identified a robust, unidirectional edge from the postcentral gyrus (somatosensory) component to the primary visual cortex (calcarine).
   • Thalamocortical Orientation: The analysis revealed that the elevated thalamo-sensorimotor coupling is oriented from the cortex (postcentral gyrus) to the thalamus, rather than thalamus-to-cortex, adding a directional layer to the canonical signature.
   • Robustness: At the $N=10$ node resolution, 11 of the 14 significant edges were classified as "undersampling-robust," meaning their direction was invariant across the inferred rates ($u \in \{2, 3, 4, 5\}$). No edges were classified as decisively fragile.
3. Comparison with GCM: When compared to Granger Causality Mapping (GCM) on a six-component network, RnR identified a hierarchical structure driven by the Default Mode Network (VMPFC $\to$ rFIC), whereas GCM emphasized bidirectional temporal dynamics. This suggests RnR captures structural precedence that temporal-only methods miss.

Significance and Claims
The paper claims that treating temporal undersampling as a fundamental property of the measurement, rather than an artifact to be ignored, enables directed interpretations grounded in the data-generating process.

• Reframing Reproducibility: The authors argue that the success of two decades of single-timescale connectivity analysis is largely due to the "shielding" effect of coarse spatial resolution, which suppresses the perceived undersampling rate. As spatial resolution improves (e.g., moving from regions to finer parcellations), the perceived rate will rise, and single-timescale methods will increasingly fail to recover correct directions.
• Modest Interpretation: The authors explicitly state that their directed claims are conditional on modeling assumptions (first-order dynamics, stationarity, causal sufficiency at the component level). They do not claim to have discovered monosynaptic projections; rather, they have identified relations between spatially distributed independent components.
• Novelty vs. Stability: The method reframes "novel" directed edges not as unsupported findings, but as structures that undirected measures cannot resolve. It also reframes apparent sign disagreements in the literature as expected consequences of temporal subsampling (direction reversals) rather than failures of replication, provided the edges are not robust.
• Future Requirements: The work posits that as the field moves toward higher-resolution effective connectivity, explicit modeling of the sampling rate (as done by RnR) will become essential, as the "shielding" provided by spatial aggregation will no longer suffice.