NetHGL Network-aware Hierarchical Graph Learning

Hierarchical Graph Learning Identifies Remission-Associated Network Organization in Treatment-Resistant Depression Following Stanford Neuromodulation Therapy

Four panels. A: Stanford Neuromodulation Therapy targets a left DLPFC site anticorrelated with the sgACC. B: scans per visit by clinical state in the active and sham arms. C: message passing among ROIs within networks, then among networks. D: the thalamus and five partner networks.
Overview. (A) SNT stimulates an individualized left DLPFC target functionally anticorrelated with the sgACC; participants received a 5-day course of active or sham SNT, 10 treatments per day. (B) Scans per visit by concurrent clinical state in each arm. (C) Message passing among ROIs within each functional network, then directed message passing among networks. (D) Attention attribution identifies a significant interaction between thalamus–partner balance and active SNT.

Abstract

Major depressive disorder is heterogeneous in symptoms and treatment response, which limits individualized care, particularly in treatment-resistant depression. Brain networks are consistently implicated, but their alterations are distributed across many regions.

We asked whether a graph model built on the brain's functional-network organization can classify each resting-state fMRI scan as remission or persistent depression at the time of scanning and identify network interactions associated with treatment. Using longitudinal data from a randomized, sham-controlled Stanford Neuromodulation Therapy (SNT) trial (NCT04739969), we developed NetHGL, a graph neural network that learns network-level representations from fMRI through intra-network message passing, hierarchical aggregation, and cross-network interactions.

NetHGL achieved higher mean performance than conventional machine-learning and graph-neural-network baselines, with an area under the receiver operating characteristic curve of 0.717 (95% CI, 0.665–0.770), balanced accuracy of 0.663 (0.614–0.712), and specificity of 0.700 (0.600–0.790) across 10 subject-level splits. Model interpretation revealed distributed remission-associated network differences, with thalamic interactions with cortical and striatal networks emerging as a prominent treatment-related component.

These findings suggest that NetHGL captures clinically relevant whole-brain states. They point to reorganization of these cortical–subcortical interactions as a candidate correlate of remission after active SNT.

  • 766resting-state fMRI scans
  • 42participants, 20 active SNT and 22 sham
  • 450 → 24ROIs pooled into networks
  • 0.717mean AUC over 10 balanced random splits

Contributions

Time-series-informed graph construction

ROI node features are summary statistics of each ROI's time series under within-ROI and whole-brain normalizations, accommodating scans of varying duration.

Hierarchical ROI-to-network modeling

The model learns community-aware ROI embeddings, pools them into network-level nodes with learned weights, and passes messages among networks. It retains ROI-level detail within networks while modeling higher-order interactions among them.

Network-level interpretability

Learned attention on network-level edges is analyzed to identify interactions associated with remission and treatment, generating hypotheses about pathways of treatment response.

Method

Network-aware Hierarchical Graph Learning

Each scan is parcellated with a combined atlas: 400 cortical regions from the Schaefer atlas, assigned to the 17 Yeo functional networks, plus 50 subcortical regions from the Tian Scale III atlas, grouped into 7 structures. That gives 450 ROIs in 24 networks. NetHGL passes messages among the ROIs of each network, pools every network into one node with learned weights, lets the 24 network nodes attend to one another, and combines them in a weighted ensemble that classifies the scan as depression or remission. The graph is the same for every scan; scans differ in their node features and in the attention the model assigns.

Two-panel architecture diagram. A, graph construction: fMRI and MADRS at each visit, BOLD fMRI, preprocessing and parcellation into 450 ROIs, node features, ROI graph, and the remitted or depressed label. B, network-aware hierarchical graph learning: within-network message passing with two GraphSAGE layers, ROI-to-network weighted pooling from 450 ROIs to 24 networks, network-level message passing with two GATv2 layers, and a network weighted ensemble that outputs remitted or depressed.
NetHGL architecture. (A) Graph construction. Each scan is parcellated into 450 ROIs, each ROI's time series is summarized as node features, and ROIs are connected within their functional network. The label is the concurrent clinical state from the MADRS. (B) Network-aware hierarchical graph learning. Two GraphSAGE layers pass messages within networks, learned weights pool the ROI embeddings into network nodes, two GATv2 layers pass messages among networks, and a weighted ensemble over networks predicts depression versus remission.

Results

Classification of concurrent clinical state

Relative to the seven baselines, NetHGL achieved the highest mean value of every evaluated metric across 10 balanced random participant-level splits. Its AUC advantage over every baseline was significant under the normal approximation (Holm p ≤ 0.048); the differences from GAT and BrainGNN remained significant under the t(33) sensitivity analysis.

Classification performance. Mean [95% participant-cluster bootstrap CI] across 10 matched participant-level splits. Positive class (+): depression. Best and second-best means are highlighted.
ModelACCBalanced ACCAUCSpecificityRecall (+)Precision (+)F1 (+)
SVM0.516[0.462–0.574]0.565[0.513–0.637]0.604[0.533–0.702]0.657[0.539–0.799]0.473[0.396–0.563]0.823[0.705–0.917]0.589[0.500–0.659]
Random Forest0.528[0.471–0.584]0.550[0.481–0.626]0.594[0.503–0.708]0.590[0.440–0.747]0.509[0.440–0.574]0.810[0.690–0.902]0.592[0.508–0.651]
CNN0.569[0.517–0.626]0.547[0.474–0.630]0.632[0.557–0.728]0.507[0.354–0.672]0.588[0.532–0.661]0.797[0.701–0.888]0.669[0.619–0.714]
GCN0.600[0.560–0.642]0.608[0.565–0.659]0.663[0.621–0.716]0.623[0.538–0.724]0.593[0.540–0.655]0.841[0.755–0.911]0.691[0.642–0.728]
GraphSAGE0.616[0.570–0.666]0.616[0.565–0.670]0.678[0.622–0.755]0.617[0.514–0.718]0.616[0.562–0.687]0.849[0.772–0.915]0.705[0.661–0.743]
GAT0.587[0.551–0.629]0.589[0.551–0.633]0.646[0.609–0.699]0.593[0.515–0.676]0.585[0.543–0.642]0.831[0.748–0.899]0.679[0.641–0.711]
BrainGNN0.599[0.560–0.643]0.621[0.593–0.656]0.650[0.614–0.700]0.663[0.611–0.731]0.580[0.529–0.637]0.860[0.780–0.920]0.682[0.634–0.721]
NetHGL (ours)0.643[0.598–0.692]0.663[0.614–0.712]0.717[0.665–0.770]0.700[0.600–0.790]0.626[0.571–0.689]0.877[0.806–0.933]0.727[0.683–0.764]

Ablations

Ablation study

All four ablations reduced mean performance across all metrics, and all four AUC reductions were significant (Holm p ≤ 0.0042). Removing network-level message passing produced the largest decline, followed by grouping the cortex into 7 instead of 17 networks.

Ablation study. Same splits, metrics and confidence intervals as above. ΔAUC is the difference from the full model.
ConfigurationAUCΔAUCBalanced ACCSpecificityRecall (+)F1 (+)
Cross-network ROI message passing0.668[0.616–0.725]−0.0490.629[0.586–0.675]0.653[0.557–0.753]0.604[0.560–0.661]0.702[0.661–0.734]
Mean ROI-to-network pooling0.662[0.613–0.718]−0.0550.608[0.561–0.658]0.633[0.529–0.720]0.583[0.532–0.654]0.683[0.639–0.724]
No network-level message passing0.547[0.511–0.593]−0.1700.544[0.512–0.584]0.557[0.484–0.646]0.531[0.481–0.573]0.634[0.572–0.670]
7 cortical networks0.629[0.573–0.696]−0.0880.600[0.554–0.655]0.623[0.523–0.740]0.577[0.522–0.648]0.678[0.635–0.714]
Full NetHGL (17 cortical networks)0.717[0.665–0.770]—0.663[0.614–0.712]0.700[0.600–0.790]0.626[0.571–0.689]0.727[0.683–0.764]

Findings

What the learned network interactions show

Network balance is a network's outgoing minus incoming attribution toward the remission output, computed on held-out test scans. The directions are those of the model's message passing; the physiological direction of influence remains open.

Remitted scans shift most cortical networks toward outgoing attribution

In a descriptive contrast, remitted scans showed a broad cortical shift toward greater outgoing relative to incoming attribution compared with depressed scans, particularly across control, default-mode, dorsal-attention, and salience/ventral-attention networks. Peripheral visual, somatomotor, and temporal parietal networks shifted toward incoming dominance.

Subcortical structures mostly moved the other way: the hippocampus, caudate, amygdala, and nucleus accumbens shifted toward lower outgoing dominance, while the thalamus and globus pallidus shifted toward outgoing dominance.

Largest contrasts (SD units): caudate −0.41, peripheral visual −0.36, temporal parietal −0.36, ContC +0.32.

Brain maps of the remitted-minus-depressed contrast of network balance. Cortical surfaces show most control, default-mode, dorsal-attention and salience networks in red, and peripheral visual, somatomotor A and temporal parietal networks in blue. Subcortical slices show the caudate, amygdala, hippocampus and nucleus accumbens in blue, and the globus pallidus and thalamus in light red.
Remitted minus depressed network balance, in SD units. Red: more outgoing attribution in remitted scans; blue: more incoming.

Discussion

Limitations and future work

Limitations

  1. Sample size. The sample was small relative to the model's complexity. Discrimination was moderate (AUC 0.717, 95% CI 0.665–0.770), and the AUC advantage over most baselines was not robust to the small-cluster t(33) correction, so performance does not yet support robust individual-level prediction or target selection.
  2. Single-hop attention. Network-level attention captures single-hop directed relationships; multi-hop propagation, recurrent interactions, and circuit motifs lie outside it.
  3. One cohort. The model was evaluated on a single SNT dataset (NCT04739969). External validation across sites, scanners, acquisition protocols, preprocessing pipelines, and diagnostic groups is needed.
  4. Time. Remission occurred mainly after treatment began, so remission and depression scans systematically differ in visit and time since treatment. The classifier was not adjusted for visit, and although the attribution models adjusted for days since initial treatment, the model may partly reflect time-related differences.
  5. Concurrent state. The labels describe clinical state at the time of scanning. The model classifies scans acquired during remission or depression; predicting future treatment response would have greater clinical utility.
  6. Repeated scans. The 766 scans come from 42 participants. Participant-level splits prevent leakage, and the attribution analyses used participant-clustered errors and participant-preserving permutations, but each test set held seven participants and only six sham-arm participants contributed remission scans, so sham-arm estimates and the treatment interactions have limited precision.
  7. Attribution. Attribution scores reflect the features the model uses to separate remission from persistent depression. The findings are hypothesis-generating and call for prospective and mechanistic validation before any use as a biomarker of remission.

Future work

  • Multimodal graphs. Integrate structural and functional connectivity, with anatomical pathways constraining functional interactions.
  • Prospective prediction. With larger longitudinal cohorts, move from concurrent-state classification to prediction of individual treatment response.
  • Change over treatment. Identify network changes that precede or track improvement.
  • Targeting. If the thalamic findings are confirmed, evaluate subject-specific target–thalamus connectivity as a soft prior, alongside target–sgACC connectivity, for ranking candidate cortical sites.

Cite

BibTeX

@article{nethgl2026,
  title  = {{NetHGL}: Hierarchical Graph Learning Identifies Remission-Associated
            Network Organization in Treatment-Resistant Depression Following
            Stanford Neuromodulation Therapy},
  year   = {2026},
  note   = {Preprint}
}