Network Response

In this tutorial we will cover how to do statistical significance testing for a ‘network response’. This is the average state/mode time course epoched around events of interest. We test is the response is significantly different from zero.

Get network response

First we need to get the network response for each subject. This is the epoched state/mode time course averaged over trials, which would be a (subjects, states/modes, time) array. In this tutorial, we will simulate this.

import numpy as np

n_subjects = 20
n_time = 250
n_networks = 6

network_response = np.random.normal(size=(n_subjects, n_time, n_networks))

# Add non-zero network response for the first network for time points 30-80
network_response[:, 30:80, 0] += 1.5

# Add non-zero network response for the second network for time points 80-100
network_response[:, 80:100, 1] += 2

Statistical significance testing

osl-dynamics has the analysis.statistics.evoked_response_max_stat_perm function for doing GLM permutation stats testing see if a value is significantly different from zero. This function uses the maximum test statistic to control for multiple comparisons across the time points and networks.

from osl_dynamics.analysis import statistics

pvalues = statistics.evoked_response_max_stat_perm(
    network_response,
    n_perm=1000,
    n_jobs=4,
)
print(pvalues.shape)

QUEUEING TASKS | Running permutations:   0%|          | 0/1000 [00:00<?, ?it/s]
QUEUEING TASKS | Running permutations:   0%|          | 1/1000 [00:00<01:52,  8.87it/s]
QUEUEING TASKS | Running permutations: 100%|██████████| 1000/1000 [00:00<00:00, 7056.77it/s]

PROCESSING TASKS | Running permutations:   0%|          | 0/1000 [00:00<?, ?it/s]
PROCESSING TASKS | Running permutations:   2%|▏         | 20/1000 [00:00<00:04, 198.74it/s]
PROCESSING TASKS | Running permutations:   5%|▍         | 49/1000 [00:00<00:03, 251.23it/s]
PROCESSING TASKS | Running permutations:   8%|▊         | 77/1000 [00:00<00:03, 261.54it/s]
PROCESSING TASKS | Running permutations:  11%|█         | 107/1000 [00:00<00:03, 274.53it/s]
PROCESSING TASKS | Running permutations:  14%|█▎        | 137/1000 [00:00<00:03, 283.04it/s]
PROCESSING TASKS | Running permutations:  17%|█▋        | 168/1000 [00:00<00:02, 290.99it/s]
PROCESSING TASKS | Running permutations:  20%|█▉        | 198/1000 [00:00<00:02, 284.53it/s]
PROCESSING TASKS | Running permutations:  23%|██▎       | 227/1000 [00:00<00:02, 277.77it/s]
PROCESSING TASKS | Running permutations:  26%|██▌       | 256/1000 [00:00<00:02, 279.72it/s]
PROCESSING TASKS | Running permutations:  28%|██▊       | 285/1000 [00:01<00:02, 279.25it/s]
PROCESSING TASKS | Running permutations:  31%|███▏      | 314/1000 [00:01<00:02, 280.85it/s]
PROCESSING TASKS | Running permutations:  34%|███▍      | 345/1000 [00:01<00:02, 282.80it/s]
PROCESSING TASKS | Running permutations:  38%|███▊      | 375/1000 [00:01<00:02, 282.21it/s]
PROCESSING TASKS | Running permutations:  40%|████      | 404/1000 [00:01<00:02, 283.50it/s]
PROCESSING TASKS | Running permutations:  44%|████▎     | 436/1000 [00:01<00:01, 290.09it/s]
PROCESSING TASKS | Running permutations:  47%|████▋     | 466/1000 [00:01<00:01, 290.01it/s]
PROCESSING TASKS | Running permutations:  50%|████▉     | 496/1000 [00:01<00:01, 291.42it/s]
PROCESSING TASKS | Running permutations:  53%|█████▎    | 526/1000 [00:01<00:01, 283.20it/s]
PROCESSING TASKS | Running permutations:  56%|█████▌    | 555/1000 [00:01<00:01, 274.27it/s]
PROCESSING TASKS | Running permutations:  58%|█████▊    | 583/1000 [00:02<00:01, 271.29it/s]
PROCESSING TASKS | Running permutations:  62%|██████▏   | 615/1000 [00:02<00:01, 281.69it/s]
PROCESSING TASKS | Running permutations:  64%|██████▍   | 644/1000 [00:02<00:01, 283.52it/s]
PROCESSING TASKS | Running permutations:  67%|██████▋   | 674/1000 [00:02<00:01, 286.84it/s]
PROCESSING TASKS | Running permutations:  71%|███████   | 709/1000 [00:02<00:00, 301.33it/s]
PROCESSING TASKS | Running permutations:  74%|███████▍  | 740/1000 [00:02<00:00, 288.11it/s]
PROCESSING TASKS | Running permutations:  77%|███████▋  | 769/1000 [00:02<00:00, 278.88it/s]
PROCESSING TASKS | Running permutations:  80%|███████▉  | 798/1000 [00:02<00:00, 280.41it/s]
PROCESSING TASKS | Running permutations:  83%|████████▎ | 828/1000 [00:02<00:00, 281.49it/s]
PROCESSING TASKS | Running permutations:  86%|████████▌ | 861/1000 [00:03<00:00, 290.22it/s]
PROCESSING TASKS | Running permutations:  89%|████████▉ | 891/1000 [00:03<00:00, 290.57it/s]
PROCESSING TASKS | Running permutations:  92%|█████████▏| 923/1000 [00:03<00:00, 295.79it/s]
PROCESSING TASKS | Running permutations:  95%|█████████▌| 953/1000 [00:03<00:00, 291.94it/s]
PROCESSING TASKS | Running permutations:  98%|█████████▊| 984/1000 [00:03<00:00, 290.44it/s]
PROCESSING TASKS | Running permutations: 100%|██████████| 1000/1000 [00:03<00:00, 284.06it/s]

COLLECTING RESULTS | Running permutations:   0%|          | 0/1000 [00:00<?, ?it/s]
COLLECTING RESULTS | Running permutations: 100%|██████████| 1000/1000 [00:00<00:00, 406306.69it/s]
(250, 6)

Let’s plot the group average network response with the significant time points highlighted.

from osl_dynamics.utils import plotting

t = np.arange(n_time)
group_network_response = np.mean(network_response, axis=0)
fig, ax = plotting.plot_evoked_response(
    t,
    group_network_response,
    pvalues,
    x_label="Sample",
    y_label="Network Activation",
)

Note, this function also has a covariates argument that can be used to account for confounds.