-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathrunner.py
145 lines (116 loc) · 4.24 KB
/
runner.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
import argparse
import csv
import requests
import warnings
warnings.filterwarnings("ignore")
warnings.simplefilter("ignore")
from pathlib import Path
import datalad.api as api
from datalad.api import Dataset, create, copy_file, remove, clone, get
import numpy as np
import pandas as pd
from scipy.spatial.distance import squareform
from tqdm.auto import tqdm
studies = ["CCNP", "BHRC", "NKI", "HBN", "PNC"]
git_repos = [f"https://github.com/ReproBrainChart/{s}_CPAC.git" for s in studies]
### Study Parameters ###
default_dict = dict(
task="rest",
run=None,
acq=None,
atlas="CC200",
space="MNI152NLin6ASym",
reg="36Parameter",
)
study_parameters = {study: default_dict.copy() for study in studies}
study_parameters["PNC"]["acq"] = "singleband"
study_parameters["NKI"]["acq"] = "645"
for s in ["CCNP", "BHRC", "HBN"]:
if s == "CCNP":
run = "01"
else:
run = "1"
study_parameters[s]["run"] = run
### Acquisition Times in ms ###
acquisition_times = dict(
CCNP=2.5,
BHRC=2.000,
NKI=0.645,
HBN=dict(SI=1.450, CBIC=0.800, RU=0.800, CUNY=0.8),
PNC=3.000,
)
def compute_dynamic_connectome(fpath, window_length=60, tr=None):
"""
Parameters
----------
window_length : int
Length of window to compute dynamic connectome in seconds. Default=45.
tr : int
TR in seconds
"""
assert tr is not None
df = pd.read_csv(fpath)
tr_per_window = int(np.round(window_length / tr))
corrs = []
for i in range(len(df) - tr_per_window):
corr = np.corrcoef(df[i : i + tr_per_window], rowvar=False)
idx = np.triu_indices_from(corr)
corrs.append(corr[idx])
return np.array(corrs)
def main(args):
out_path = Path(args.output)
# Grab the metadata files
for study in studies:
study_path = out_path / f"{study}"
if not study_path.exists():
study_path.mkdir(exist_ok=True)
url = f"https://raw.githubusercontent.com/ReproBrainChart/{study}_BIDS/main/study-{study}_desc-participants.tsv"
response = requests.get(url)
reader = csv.reader(response.text.splitlines(), skipinitialspace=True)
with open(study_path / f"{study}_desc-participants.tsv", "w") as f:
w = csv.writer(f)
w.writerows(reader)
# Datalad clone the datasets
api.clone(source=f"https://github.com/ReproBrainChart/{study}_CPAC.git", git_clone_opts=["-b", "complete-pass-0.1"], path=out_path / f"{study}_CPAC", )
for study, study_parameter in study_parameters.items():
# load metadata
df = pd.read_csv(out_path / f"{study}/{study}_desc-participants.tsv", delimiter="\t")
df = df[~df["p_factor_mcelroy_harmonized_all_samples"].isnull()]
print(f"Computing dynamic connectomes for {study}; Total files={len(df)}")
# Setup glob string
glob_str = "*".join(
[f"{k}-{v}" for k, v in study_parameter.items() if v is not None]
)
glob_str = "**/*" + glob_str + "*Mean_timeseries.1D"
p = out_path / f"{study}_CPAC/cpac_RBCv0"
files = []
# Loop over each row of metadata
for _, row in df.iterrows():
sub_path = p / f"sub-{row.participant_id}/ses-{row.session_id}"
# print
tmp = list(sub_path.glob(glob_str))
if len(tmp) == 0:
continue
else:
files += tmp
# Download all the files
api.get(files, dataset=out_path / f"{study}_CPAC/")
for file in files:
splits = file.name.split("_")
if study == "HBN":
site = splits[1][11:]
tr = acquisition_times[study][site]
else:
tr = acquisition_times[study]
# save file
participant_id = splits[0]
sub_path = out_path / f"{study}/{participant_id}"
sub_path.mkdir(parents=True, exist_ok=True)
out_fpath = sub_path / file.name
arr = compute_dynamic_connectome(file, tr=tr)
np.save(out_fpath, arr)
if __name__ == '__main__':
parser = argparse.ArgumentParser()
parser.add_argument("-o", "--output", help="output path", type=str, default="/output")
args = parser.parse_args()
main(args)