OSCR

Local and global patterns support medical imaging as a biomarker of ageing.

Code ↔ Paper

6 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 6 matches
  1. [1] § Results › Lifestyle and Environment ↔ source/analysis/perform_phewas.py, lines 154–213 · score 0.96 · family history, sun exposure, mental health, sexual factors, social support, general health
  2. [2] § Results › Chronic Diseases ↔ source/analysis/disease_analysis.py, lines 51–78 · score 0.92 · ischemic heart disease, liver disease subjects, COPD subjects, MS subjects, scoliosis subjects, hypertension subjects
  3. [3] § Methods › Lifestyle and Environment ↔ source/analysis/perform_phewas.py, lines 320–399 · score 0.63 · rank normalisation, imaging phenotypes, variables, confounding, PheWAS, correlate
  4. [4] § Methods › Deep Learning Models ↔ source/cnn/model/resnet.py, lines 185–193 · score 0.61 · resnet18, ImageNet, pre trained, models
  5. [5] § Methods › Survival Analysis ↔ source/analysis/survival.py, lines 1–44 · score 0.60 · Cox proportional hazards, survival, death, female, brain, body
  6. [6] § Methods › Dataset › UK Biobank ↔ source/cnn/train_utils.py, lines 272–295 · score 0.57 · largest connected component, voxels, mask, training

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

Python · 612 lines · 24 KB · no license · 2 matches

  1. # %%
  2. import numpy as np
  3. import pandas as pd
  4. import ast
  5. import os
  6. from tqdm import tqdm
  7. import seaborn as sns
  8. import scipy.stats
  9. import math
  10. import csv
  11. import matplotlib.pyplot as plt
  13. import utils
  14. import paths
  15. #%%
  16. ############################################################################################################################################################################
  17. #### This script is used to generate the Manhattan plot and the correlation plot between non-imaging phenotypes and imaging phenotypes ###
  18. #### The script is adapted from Bai et al. 2020, "A population-based phenome-wide association study of cardiac and aortic structure and function", Nature medicine. ###
  19. ############################################################################################################################################################################
  21. #%%
  22. def normalise(x):
  23. return (x - np.mean(x)) / np.std(x)
  26. def rank_normalise(x):
  27. # Rank-based inverse normal transform
  28. # Please refer to the function inormal() in the FSLNets package
  30. # Get the rank of the values in x
  31. ri = np.argsort(np.argsort(x))
  33. # Correct for the ranks of repeated values
  34. # argsort assign different ranks for these values
  35. # We fill them with the same value
  36. u, inv_idx = np.unique(x, return_inverse=True)
  37. sii = np.sort(inv_idx)
  38. repeated_idx = np.unique(sii[np.diff(np.append(sii, 1)) == 0])
  39. for i in repeated_idx:
  40. ri[inv_idx == i] = np.mean(ri[inv_idx == i])
  42. # Perform inverse normal transform
  43. # ri + 1 so that the rank starts from 1, to be consistent with Karla's Matlab code
  44. # p squashes the rank into the range of [0, 1]
  45. # erfinv can generate a distribution from 2 * p - 1 with 0 mean and 1 standard deviation
  46. N = len(x)
  47. ri = ri + 1
  48. c = 3.0 / 8
  49. p = (ri - c) / (N - 2 * c + 1)
  50. y = math.sqrt(2) * scipy.special.erfinv(2 * p - 1)
  51. return y
  53. def p_adjust_fdr(p):
  54. # FDR correction for multiple testing
  55. # The implementation provides consistent results as R p.adjust function.
  56. # Use a new np array, otherwise p will be modified.
  57. p2 = np.zeros(p.shape, dtype=np.float32)
  58. idx = np.argsort(p)
  59. n = len(p)
  60. p2[idx] = (p[idx] * n) / np.arange(1, n + 1)
  61. p2[p2 > 1] = 1
  62. return p2
  65. def fdr_threshold(p, q):
  66. # p: vector of p-values
  67. # q: false discovery rate level
  68. #
  69. # return values
  70. # pID: p-value threshold based on independence or positive dependence
  71. # pN: nonparametric p-value threshold
  72. #
  73. # This function takes a vector of p-values and a FDR rate.
  74. # It returns two p-value thresholds, one based on an assumption of
  75. # independence or positive dependence, and one that makes no assumptions
  76. # about how the tests are correlated. For imaging data, an assumption of
  77. # positive dependence is reasonable, so it should be OK to use the first
  78. # (more sensitive) threshold.
  79. #
  80. # This implementation provides consistent results as the FDR function at
  81. # https://warwick.ac.uk/fac/sci/statistics/staff/academic-research/nichols/software/fdr
  82. # https://warwick.ac.uk/fac/sci/statistics/staff/academic-research/nichols/software/fdr/FDR.m
  83. p2 = p[~np.isnan(p)]
  84. p2 = np.sort(p2)
  85. n = len(p2)
  86. I = np.arange(1, n + 1)
  87. cVID = 1
  88. cVN = np.sum(1 / I)
  90. idx = np.nonzero(p2 <= ((I * q) / (n * cVID)))[0]
  91. pID = p2[np.max(idx)] if len(idx) >= 1 else 0
  93. idx = np.nonzero(p2 <= ((I * q) / (n * cVN)))[0]
  94. pN = p2[np.max(idx)] if len(idx) >= 1 else 0
  95. return pID, pN
  96. #%%
  97. ####################
  98. ### Loading Data ###
  99. ####################
  101. predictions_whole_body = pd.read_csv(paths.path_to_predictions_whole_body)
  102. predictions_brain = pd.read_csv(paths.path_to_predictions_brain)
  103. predictions_lungs = pd.read_csv(paths.organ_paths["lungs"])
  104. predictions_heart = pd.read_csv(paths.organ_paths["heart"])
  105. predictions_spine = pd.read_csv(paths.organ_paths["spine"])
  106. predictions_liver = pd.read_csv(paths.organ_paths["liver"])
  107. predictions_muscle = pd.read_csv(paths.organ_paths["muscle"])
  108. predictions_intestine = pd.read_csv(paths.organ_paths["intestine"])
  110. predictions_whole_body.rename(columns={"age_gap":"age_gap_whole_body"}, inplace=True)
  111. predictions_whole_body.rename(columns={"mean_prediction": "whole_body_pred"}, inplace=True)
  112. predictions_brain.rename(columns={"age_gap":"age_gap_brain"}, inplace=True)
  113. predictions_brain.rename(columns={"mean_prediction": "brain_pred"}, inplace=True)
  114. predictions_brain = predictions_brain[["eid", "age_gap_brain", "brain_pred"]]
  115. predictions_lungs.rename(columns={"age_gap":"age_gap_lungs"}, inplace=True)
  116. predictions_lungs.rename(columns={"mean_prediction": "lungs_pred"}, inplace=True)
  117. predictions_lungs = predictions_lungs[["eid", "age_gap_lungs", "lungs_pred"]]
  118. predictions_heart.rename(columns={"age_gap":"age_gap_heart"}, inplace=True)
  119. predictions_heart.rename(columns={"mean_prediction": "heart_pred"}, inplace=True)
  120. predictions_heart = predictions_heart[["eid", "age_gap_heart", "heart_pred"]]
  121. predictions_spine.rename(columns={"age_gap":"age_gap_spine"}, inplace=True)
  122. predictions_spine.rename(columns={"mean_prediction": "spine_pred"}, inplace=True)
  123. predictions_spine = predictions_spine[["eid", "age_gap_spine", "spine_pred"]]
  124. predictions_liver.rename(columns={"age_gap":"age_gap_liver"}, inplace=True)
  125. predictions_liver.rename(columns={"mean_prediction": "liver_pred"}, inplace=True)
  126. predictions_liver = predictions_liver[["eid", "age_gap_liver", "liver_pred"]]
  127. predictions_muscle.rename(columns={"age_gap":"age_gap_muscle"}, inplace=True)
  128. predictions_muscle.rename(columns={"mean_prediction": "muscle_pred"}, inplace=True)
  129. predictions_muscle = predictions_muscle[["eid", "age_gap_muscle", "muscle_pred"]]
  130. predictions_intestine.rename(columns={"age_gap":"age_gap_intestine"}, inplace=True)
  131. predictions_intestine.rename(columns={"mean_prediction": "intestine_pred"}, inplace=True)
  132. predictions_intestine = predictions_intestine[["eid", "age_gap_intestine", "intestine_pred"]]
  134. predictions = pd.merge(predictions_whole_body, predictions_brain, on="eid", how="left")
  135. predictions = pd.merge(predictions, predictions_lungs, on="eid", how="left")
  136. predictions = pd.merge(predictions, predictions_heart, on="eid", how="left")
  137. predictions = pd.merge(predictions, predictions_spine, on="eid", how="left")
  138. predictions = pd.merge(predictions, predictions_liver, on="eid", how="left")
  139. predictions = pd.merge(predictions, predictions_muscle, on="eid", how="left")
  140. predictions = pd.merge(predictions, predictions_intestine, on="eid", how="left")
  142. predictions = predictions[["eid", "age_gap_whole_body",
  143. "age_gap_brain",
  144. "whole_body_pred",
  145. "brain_pred",
  146. "age_gap_lungs", "lungs_pred",
  147. "age_gap_heart", "heart_pred",
  148. "age_gap_spine", "spine_pred",
  149. "age_gap_liver", "liver_pred",
  150. "age_gap_muscle", "muscle_pred",
  151. "age_gap_intestine", "intestine_pred"
  152. ]]
  154. #%%
  155. ################################################
  156. ### Merging data with non-imaging phenoytpes ###
  157. ################################################
  158. merged_data = pd.merge(predictions, pd.read_csv(paths.pain_path), on="eid", how="left")
  159. merged_data = pd.merge(merged_data, pd.read_csv(paths.sleep_path), on="eid", how="left")
  160. merged_data = pd.merge(merged_data, pd.read_csv(paths.sexual_factors_path), on="eid", how="left")
  161. merged_data = pd.merge(merged_data, pd.read_csv(paths.alcohol_path), on="eid", how="left")
  162. merged_data = pd.merge(merged_data, pd.read_csv(paths.diet_path), on="eid", how="left")
  163. merged_data = pd.merge(merged_data, pd.read_csv(paths.employment_path), on="eid", how="left")
  164. merged_data = pd.merge(merged_data, pd.read_csv(paths.family_history_path), on="eid", how="left")
  165. merged_data = pd.merge(merged_data, pd.read_csv(paths.female_sex_specific_path), on="eid", how="left")
  166. merged_data = pd.merge(merged_data, pd.read_csv(paths.general_health_path), on="eid", how="left")
  167. merged_data = pd.merge(merged_data, pd.read_csv(paths.male_sex_specific_path), on="eid", how="left")
  168. merged_data = pd.merge(merged_data, pd.read_csv(paths.medical_conditions_path), on="eid", how="left")
  169. merged_data = pd.merge(merged_data, pd.read_csv(paths.mental_health_path), on="eid", how="left")
  170. merged_data = pd.merge(merged_data, pd.read_csv(paths.smoking_path), on="eid", how="left")
  171. merged_data = pd.merge(merged_data, pd.read_csv(paths.social_support_path), on="eid", how="left")
  172. merged_data = pd.merge(merged_data, pd.read_csv(paths.sun_exposure_path), on="eid", how="left")
  173. merged_data = pd.merge(merged_data, pd.read_csv(paths.telomeres_path), on="eid", how="left")
  174. merged_data = pd.merge(merged_data, pd.read_csv(paths.month_of_birth_path), on="eid", how="left")
  175. merged_data = pd.merge(merged_data, pd.read_csv(paths.deaths_path, sep="\t"), on="eid", how="left")
  176. merged_data = pd.merge(merged_data, pd.read_csv(paths.body_composition_path), on="eid", how="left")
  177. met = pd.read_csv(paths.path_to_MET_features)
  178. merged_data = pd.merge(merged_data, met, on="eid", how="left")
  180. fiftythree = utils.get_feature("53-2.0")
  182. merged_data = pd.merge(merged_data, fiftythree, on="eid", how="left")
  184. loneliness = utils.get_feature("2020-2.0")
  185. merged_data = pd.merge(merged_data, loneliness, on="eid", how="left")
  187. year_of_birth = utils.get_feature("34-0.0")
  188. merged_data = pd.merge(merged_data, year_of_birth, on="eid", how="left")
  189. merged_data = merged_data.dropna(subset=["age_gap_lungs", "lungs_pred", "age_gap_heart", "heart_pred",
  190. "age_gap_spine", "spine_pred","age_gap_liver", "liver_pred",
  191. "age_gap_muscle", "muscle_pred",
  192. "age_gap_intestine", "intestine_pred",
  193. ])
  194. df = merged_data
  195. df.dropna(subset=["34-0.0", "52-0.0", "21002-2.0"], inplace=True)
  196. predictions = predictions[predictions["eid"].isin(df["eid"])]
  197. # %%
  198. # # # # # # # # # # # # # # # # # # # #
  199. # Step 3: confounding factors
  200. # # # # # # # # # # # # # # # # # # # #
  201. sex = df['31-0.0'].values
  202. # Age provided by UK Biobank (21003-2.0) seems to be floored, i.e. with half a year error.
  203. # To get more accurate age values, we calculate age by date.
  204. age = np.zeros(len(df))
  205. for i in range(len(df)):
  206. # Calculate age
  207. d1 = datetime.date(int(df.iloc[i]['34-0.0']), int(df.iloc[i]['52-0.0']), 15)
  208. s = df.iloc[i]['53-2.0']
  209. d2 = datetime.date(int(s[:4]), int(s[5:7]), int(s[8:]))
  210. age[i] = np.round((d2 - d1).days / 365.25, 1)
  211. weight = df['21002-2.0'].values
  212. bmi = df['21001-2.0'].values
  213. height = np.round(np.sqrt(weight / bmi) * 100)
  215. # %%
  216. valid_idx = ~np.isnan(age) & ~np.isnan(sex) & ~np.isnan(weight) & ~np.isnan(height)
  217. df = df[valid_idx]
  218. df_idp = predictions[valid_idx]
  219. sex = sex[valid_idx]
  220. age = age[valid_idx]
  221. sex_age = sex * age # Sex and age interaction
  222. weight = weight[valid_idx]
  223. height = height[valid_idx]
  224. conf = np.stack((sex, age, sex_age, weight, height), axis=1)
  225. df_conf = pd.DataFrame(conf, index=df.index, columns=['Sex', 'Age', 'Sex * Age', 'Weight', 'Height'])
  226. # %%
  227. # Remove confounding factors from df
  228. df = df.drop(columns=[
  229. '34-0.0',
  230. '52-0.0',
  231. '21002-2.0',
  232. 'whole_body_pred',
  233. 'brain_pred',
  234. 'lungs_pred',
  235. 'heart_pred',
  236. 'spine_pred',
  237. 'liver_pred',
  238. 'muscle_pred',
  239. 'intestine_pred',
  240. 'age_gap_brain',
  241. 'age_gap_whole_body',
  242. 'age_gap_lungs',
  243. 'age_gap_heart',
  244. 'age_gap_spine',
  245. 'age_gap_liver',
  246. "age_gap_muscle",
  247. "age_gap_intestine",
  248. 'date_of_death',
  249. 'ins_index',
  250. 'dsource',
  251. 'source',
  252. "eid",
  253. ])
  254. # %%
  255. # # # # # # # # # # # # # # # # # # # #
  256. # Step 4: clean, de-confound and normalise data
  257. # # # # # # # # # # # # # # # # # # # #
  258. # This part of the code was adapted from Karla Miller's Matlab code at
  259. # https://www.fmrib.ox.ac.uk/ukbiobank/gwaspaper/index.html
  260. # Step 4.1: cleaning
  261. n_subj, n_col = df.shape
  262. bad_vars = []
  263. for i in range(n_col):
  264. val = df.iloc[:, i]
  266. # Discard columns which not numbers
  267. if not np.issubdtype(df.dtypes[i], np.number):
  268. bad_vars += [i]
  269. continue
  271. # Assume negative values are invalid, set them to NaN
  272. # There are also a lot of empty values, which are already NaN
  273. val[val < 0] = np.nan
  274. df.iloc[:, i] = val
  276. # Valid indices
  277. valid_idx = ~np.isnan(val)
  279. # Discard columns with more than 90% missing data
  280. if np.sum(valid_idx) < (0.1 * n_subj):
  281. bad_vars += [i]
  282. continue
  284. # Discard columns with over 95% elements with the exactly same value
  285. val_unique, counts = np.unique(val[valid_idx], return_counts=True)
  286. if np.max(counts) >= (0.95 * np.sum(valid_idx)):
  287. bad_vars += [i]
  288. continue
  289. #%%
  290. merged_data["2644-2.0"]
  291. # merged_data.drop(columns=[["2644-2.0","20161-2.0"]], inplace=True)
  292. merged_data.drop(columns=["20161-2.0"], inplace=True)
  293. # %%
  294. for i in range(n_col):
  295. for j in range(i + 1, n_col):
  296. if i in bad_vars or j in bad_vars:
  297. continue
  298. # Discard columns with very high correlation
  299. val_i = df.iloc[:, i]
  300. val_j = df.iloc[:, j]
  301. valid_idx = ~np.isnan(val_i) & ~np.isnan(val_j)
  302. if np.sum(valid_idx) == 0:
  303. continue
  304. # print(val_i[valid_idx], val_j[valid_idx])
  305. cc, _ = scipy.stats.pearsonr(val_i[valid_idx], val_j[valid_idx])
  307. if cc > 0.9999:
  308. # Keep the column with more valid elements
  309. if np.sum(~np.isnan(val_i)) > np.sum(~np.isnan(val_j)):
  310. bad_vars += [j]
  311. else:
  312. bad_vars += [i]
  313. # %%
  314. # The cleaned data
  315. bad_vars = np.unique(bad_vars)
  316. keep_vars = sorted(list(set(np.arange(n_col)) - set(bad_vars)))
  317. df = df.iloc[:, keep_vars]
  318. print('{0} columns kept after data cleaning.'.format(df.shape[1]))
  320. # Step 4.2: normalise confounding factors
  321. conf = (conf - np.mean(conf, axis=0)) / np.std(conf, axis=0)
  322. #%%
  323. # Step 4.3: normalise non imaging phenotypes
  324. df_cont = pd.read_csv(paths.path_to_continuous_features, index_col=0)
  325. #%%
  326. n_col = df.shape[1]
  327. for i in range(n_col):
  328. val = df.iloc[:, i]
  329. valid_idx = ~np.isnan(val)
  330. x = val[valid_idx]
  332. # Field ID
  333. field_id = int(df.columns[i][1].split('-')[0])
  334. is_continuous = df_cont.loc[field_id]['continuous']
  336. if is_continuous:
  337. # If it is a continuous variable, perform standard normalisation.
  338. x = normalise(x)
  339. else:
  340. # If we are not sure whether it is a continuous or categorical variable, convert it
  341. # into a continuous variable using rank-based inverse normal transform.``
  342. x = rank_normalise(x)
  343. df.iloc[:, i][valid_idx] = x# The cleaned data
  344. #%%
  345. # Step 4.4: de-confound and normalise IDPs
  346. try:
  347. df_idp = df_idp.drop(columns=["eid"])
  348. df_idp = df_idp.drop(columns=["split"])
  349. except:
  350. print("Columns already dropped")
  352. n_row = conf.shape[1]
  353. n_col = df_idp.shape[1]
  354. beta = np.zeros((n_row, n_col))
  355. for i in range(n_col):
  356. val = df_idp.iloc[:, i]
  357. print(val)
  358. valid_idx = ~np.isnan(val)
  359. x = val[valid_idx]
  360. beta[:, i] = np.dot(np.linalg.pinv(conf[valid_idx]), x)
  361. x = x - np.dot(conf[valid_idx], beta[:, i])
  362. x = normalise(x)
  363. df_idp.iloc[:, i][valid_idx] = x
  365. df_idp = df_idp[[
  366. "age_gap_brain",
  367. "age_gap_whole_body",
  368. "age_gap_lungs",
  369. "age_gap_heart",
  370. "age_gap_spine",
  371. "age_gap_liver",
  372. "age_gap_muscle",
  373. "age_gap_intestine"
  374. ]]
  375. # %%
  376. # # # # # # # # # # # # # # # # # # # # #
  377. # # Step 5: uni-variate correlation studies
  378. # # # # # # # # # # # # # # # # # # # # #
  379. M = df_idp.shape[1]
  380. N = df.shape[1]
  381. corr = np.zeros((M, N))
  382. corr_p = np.zeros((M, N))
  383. # features = []
  384. for i in range(M):
  385. features = []
  387. for j in range(N):
  388. # Remove NaNs
  389. x = df_idp.iloc[:, i]
  390. y = df.iloc[:, j]
  391. valid_idx = ~np.isnan(x) & ~np.isnan(y)
  392. x = x[valid_idx]
  393. y = y[valid_idx]
  395. # Pearson correlation
  396. cc, p_val = scipy.stats.pearsonr(x, y)
  397. corr[i, j] = cc
  398. corr_p[i, j] = p_val
  399. features.append(df.columns[j])
  401. # For p-value of 0, assign it with the mininal positive floating value
  402. # so that we can calculate the logarithm for the Manhattan plot
  403. corr_p[corr_p == 0] = np.finfo(np.float64).tiny
  405. # Logarithm
  406. log_corr_p = - np.log10(corr_p)
  408. # Save the tables
  409. df_corr = pd.DataFrame(corr, index=df_idp.columns, columns=df.columns)
  410. df_p = pd.DataFrame(corr_p, index=df_idp.columns, columns=df.columns)
  411. df_log_p = pd.DataFrame(log_corr_p, index=df_idp.columns, columns=df.columns)
  413. corr = df_corr.values
  414. corr_p = df_p.values
  415. log_corr_p = df_log_p.values
  417. # Bonferroni correction
  418. M, N = corr.shape
  419. p_bonf = 0.05 / (M * N)
  421. # FDR correction
  422. p_fdr, _ = fdr_threshold(corr_p.flatten(), 0.05)
  424. # Number of phenotypes that is significantly associated with at least one of the IDPs
  425. print('p_bonf = {0}'.format(p_bonf))
  426. print('p_fdr = {0}'.format(p_fdr))
  427. print('Number of correlations reaching Bonferroni threshold = {0}'.format(np.sum(corr_p < p_bonf)))
  428. print('Number of correlations reaching FDR threshold = {0}'.format(np.sum(corr_p < p_fdr)))
  429. print('Number of phenotypes reaching Bonferroni threshold = {0}'.format(np.sum(np.sum(corr_p < p_bonf, axis=0) > 0)))
  430. print('Number of phenotypes reaching FDR threshold = {0}'.format(np.sum(np.sum(corr_p < p_fdr, axis=0) > 0)))
  431. # %%
  432. # # # # # # # # # # # # # # # # # # # #
  433. # Step 6: Manhattan plot
  434. # # # # # # # # # # # # # # # # # # # #
  435. # The category for each column
  436. # They should be in ascending order. Otherwise, we need to sort the category IDs.
  438. pain_columns = pd.read_csv(paths.pain_path).columns
  439. sleep_columns = pd.read_csv(paths.sleep_path).columns
  440. sexual_columns = pd.read_csv(paths.sexual_factors_path).columns
  441. alcohol_columns = pd.read_csv(paths.alcohol_path).columns
  442. diet_columns = pd.read_csv(paths.diet_path).columns
  443. employment_columns = pd.read_csv(paths.employment_path).columns
  444. family_history_columns = pd.read_csv(paths.family_history_path).columns
  445. smoking_columns = pd.read_csv(paths.smoking_path).columns
  446. female_specific_columns = pd.read_csv(paths.female_sex_specific_path).columns
  447. general_health_columns = pd.read_csv(paths.general_health_path).columns
  448. male_specific_columns = pd.read_csv(paths.male_sex_specific_path).columns
  449. medical_conditions_columns = pd.read_csv(paths.medical_conditions_path).columns
  450. mental_health_columns = pd.read_csv(paths.mental_health_path).columns
  451. social_support_columns = pd.read_csv(paths.social_support_path).columns
  452. sun_exposure_columns = pd.read_csv(paths.sun_exposure_path).columns
  453. telomeres_columns = pd.read_csv(paths.telomeres_path).columns
  454. body_composition_columns = pd.read_csv(paths.body_composition_path).columns
  455. met = pd.read_csv(paths.path_to_MET_features).columns
  457. def get_category(field_id):
  459. if field_id in pain_columns:
  460. return "pain"
  461. elif field_id in sleep_columns:
  462. return "sleep"
  463. elif field_id in sexual_columns:
  464. return "sexual"
  465. elif field_id in alcohol_columns:
  466. return "alcohol"
  467. elif field_id in diet_columns:
  468. return "diet"
  469. elif field_id in employment_columns:
  470. return "employment"
  471. elif field_id in family_history_columns:
  472. return "family_history"
  473. elif field_id in smoking_columns:
  474. return "smoking"
  475. elif field_id in female_specific_columns:
  476. return "female_specific"
  477. elif field_id in general_health_columns:
  478. return "general_health"
  479. elif field_id in male_specific_columns:
  480. return "male_specific"
  481. elif field_id in medical_conditions_columns:
  482. return "medical_conditions"
  483. elif field_id in mental_health_columns:
  484. return "mental_health"
  485. elif field_id in social_support_columns or field_id=="2020-2.0" or field_id=="2020-2.0" or field_id=="6160-2.0_y" or field_id=="6160-0.0":
  486. return "social_support"
  487. elif field_id in sun_exposure_columns:
  488. return "sun_exposure"
  489. elif field_id in telomeres_columns:
  490. return "telomeres"
  491. elif field_id in body_composition_columns:
  492. return "body_composition"
  493. elif field_id in met:
  494. return "physical activity"
  495. else:
  496. return "general"
  498. # The Manhattan plot
  499. table = []
  500. for i in range(M):
  501. if df_idp.columns[i] == 'age_gap_whole_body':
  502. s = 'whole_body'
  503. elif df_idp.columns[i] == 'age_gap_brain':
  504. s = 'brain'
  505. elif df_idp.columns[i] == 'age_gap_lungs':
  506. s = 'lungs'
  507. elif df_idp.columns[i] == 'age_gap_heart':
  508. s = 'heart'
  509. elif df_idp.columns[i] == 'age_gap_spine':
  510. s = 'spine'
  511. elif df_idp.columns[i] == 'age_gap_liver':
  512. s = 'liver'
  513. elif df_idp.columns[i] == 'age_gap_muscle':
  514. s = 'muscle'
  515. elif df_idp.columns[i] == 'age_gap_intestine':
  516. s = 'intestine'
  518. for j in range(N):
  519. line = [j, log_corr_p[i, j], corr[i, j], abs(corr[i, j]), s, features[j]]
  520. table += [line]
  522. df_table = pd.DataFrame(table, columns=['x', 'log_p', 'corr', 'Correlation', 'Modality', "Feature"])
  523. #%%
  524. category = np.array([get_category(field_id) for field_id in tqdm(df_table["Feature"])])
  525. df_table["Category"] = category
  526. df_table = df_table.sort_values(by=["Category","Feature"])
  527. #%%
  528. df_table = df_table.sort_values(by=["Category","Feature", "Modality"])
  529. #%%
  530. df_table["x2"] = np.repeat(np.arange(1, N+1), M)
  531. #%%
  532. df_table = df_table.sample(frac=1)
  533. #%%
  534. df_table.reset_index(drop=True, inplace=True)
  535. #%%
  536. pallette = {
  537. 'heart': '#D81B60',
  538. 'muscle': "#1E88E5",
  539. 'whole_body': "#FFBF00",
  540. 'liver': '#004D40',
  541. 'brain': '#979E41',
  542. 'spine': '#51D9B8',
  543. 'lungs': '#AA58B3',
  544. 'intestine': '#C56800',
  545. }
  546. #%%
  547. plt.figure(figsize=(16, 8))
  548. ax = sns.scatterplot(x='x2', y='log_p', hue='Modality', size='Correlation',
  549. sizes=(30, 200), size_norm=mpl.colors.Normalize(vmin=0, vmax=0.3),
  550. data=df_table, alpha=0.8, palette=pallette)
  551. plt.plot([0, N], [-math.log10(p_bonf), -math.log10(p_bonf)], 'k--', linewidth=1, alpha=0.9)
  552. plt.text(-1, -math.log10(p_bonf), 'Bonf', horizontalalignment='right', fontsize=15)
  553. handles, labels = ax.get_legend_handles_labels()
  554. # labels[-1] = '0.3'
  555. ax.legend(handles, labels)
  557. xticks = []
  558. xticklabels = []
  559. category_list = df_table["Category"].values
  560. unique_category = np.unique(category_list)
  561. print(unique_category)
  562. for c in unique_category:
  563. x = np.max(df_table[df_table["Category"]==c]["x2"]) + 0.5
  564. plt.plot([x, x], [0, 310], 'k-.', linewidth=0.5, alpha=0.2)
  565. xticks += [np.mean(df_table[df_table["Category"]==c]["x2"])]
  566. xticklabels += [c]
  568. plt.xlim(0, N)
  569. plt.ylim(0, 110)
  570. plt.xticks(xticks, xticklabels, fontsize=15, rotation=90)
  571. plt.yticks(fontsize=15)
  572. plt.xlabel('Non-imaging phenotypes', fontsize=15, fontweight='bold')
  573. plt.ylabel(r'$-\log_{10}(p)$', fontsize=15, fontweight='bold')
  574. fig = plt.gcf()
  575. fig.set_size_inches(16, 8)
  576. plt.tight_layout()
  577. plt.savefig('../../outputs/figures/manhattan_phewas.pdf', bbox_inches='tight')
  578. plt.show()
  579. # %%
  580. plt.figure(figsize=(16, 8))
  582. ax = sns.scatterplot(x='x2', y='corr', hue='Modality', size='Correlation',
  583. sizes=(30, 100), size_norm=mpl.colors.Normalize(vmin=0, vmax=0.2),
  584. data=df_table, alpha=0.8, palette=palette)
  585. plt.hlines(0, xmin=0, xmax=N*2, linewidth=1, alpha=0.8, color="gray", linestyle='--')
  586. handles, labels = ax.get_legend_handles_labels()
  587. # labels[-1] = '0.25'
  588. ax.legend(handles, labels)
  590. xticks = []
  591. xticklabels = []
  592. category_list = df_table["Category"].values
  593. unique_category = np.unique(category_list)
  594. print(unique_category)
  595. for c in unique_category:
  597. x = np.max(df_table[df_table["Category"]==c]["x2"]) + 0.5
  598. plt.plot([x, x], [-0.15, 0.21], 'k-.', linewidth=0.5)
  599. xticks += [np.mean(df_table[df_table["Category"]==c]["x2"])]
  600. xticklabels += [c]
  602. plt.xlim(0, N)
  603. plt.xticks(xticks, xticklabels, fontsize=15, rotation=90)
  604. plt.yticks(fontsize=15)
  605. plt.xlabel('Non-imaging phenotypes', fontsize=15, fontweight='bold')
  606. plt.ylabel('Correlation', fontsize=15, fontweight='bold')
  608. fig = plt.gcf()
  609. fig.set_size_inches(16, 8)
  610. plt.tight_layout()
  611. plt.savefig('../../outputs/figures/correlations_phewas.pdf', bbox_inches='tight')
  612. plt.show()

perform_phewas.py at commit a772251, no license · at the source

Overview

Authors: Tamara T. Mueller1, Sophie Starck1, Rozafë Llalloshi1, Georgios Kaissis1, Alexander Ziller1, Robert Graf1,2, Christopher Schlett3, Steffen Ringhof3, Fabian Bamberg3, Mark Wielpütz4, Henry Völzke4, Michael Leitzmann5, Thoralf Niendorf6, Thomas Keil7,8,9, Lilian Krist7, Tobias Pischon6, André Karch10, Klaus Berger10, Jan Kirschke2, Daniel Rueckert1,11, Rickmer Braren12,13,14
ORCID iDs: Lilian Krist, Tobias Pischon, Jan Kirschke
14 affiliations
  1. Lab for AI in Medicine and Healthcare, Technical University of Munich,Munich, Germany
  2. Department of Diagnostic and Interventional Neuroradiology, School of Medicine, TUM University Hospital,Munich, Germany
  3. Department of Diagnostic and Interventional Radiology, Medical Center, University of Freiburg, Faculty of Medicine,Freiburg, Germany
  4. Institut für Community Medicine, University Medicine Greifswald,Greifswald, Germany
  5. Institute for Epidemiology and Preventive Medicine, University of Regensburg,Regensburg, Germany
  6. Berlin Ultrahigh Field Facility (B.U.F.F.), Max Delbrück Center for Molecular Medicine, Helmholtz Association,Berlin, Germany
  7. Institute of Social Medicine, Epidemiology and Health Economics, Charité - Universitätsmedizin Berlin,Berlin, Germany
  8. Institute of Clinical Epidemiology and Biometry, University of Würzburg,Würzburg, Germany
  9. State Institute of Health I, Bavarian Health and Food Safety Authority,Erlangen, Germany
  10. Institute of Epidemiology and Social Medicine, University of Münster,Münster, Germany
  11. Department of Computing, Imperial College London,London, United Kingdom
  12. Institute for Diagnostic and Interventional Radiology, Klinikum rechts der Isar, Technical University of Munich,Munich, Germany
  13. German Cancer Consortium (DKTK), Munich partner site,Heidelberg, Germany
  14. Department of Diagnostic and Interventional Radiology and Nuclear Medicine, University Medical Center Hamburg-Eppendorf,Hamburg, Germany
Institutions: Technical University of Munich (Germany); LMU Klinikum (Germany); University of Freiburg (Germany); Universitätsmedizin Greifswald (Germany); University of Regensburg (Germany); Max Delbrück Center (Germany); Charité - Universitätsmedizin Berlin (Germany); University of Würzburg (Germany); Bayerisches Landesamt für Gesundheit und Lebensmittelsicherheit (Germany); University of Münster (Germany); Imperial College London (United Kingdom); Deutsches Konsortium für Translationale Krebsforschung (Germany); University Medical Center Hamburg-Eppendorf (Germany)
Journal: Communications medicine, volume 6, issue 1, article 335
Dates: received 2 June 2025; accepted 3 June 2026; published online 13 June 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1038/s43856-026-01722-3 · PMID 42288701 · PMCID PMC13264634 · OpenAlex W4410519636
Open access: gold, a free copy (OpenAlex)
Preprint: doi:10.1101/2025.05.18.654705
Status: code verified
Categories: human (organism), clinical / translational (subfield)
Methods: Connectivity, Statistics, Machine learning, fMRI & imaging
Keywords: Biomarkers
Topic: Health, Environment, Cognitive Aging (Health, Toxicology and Mutagenesis, Environmental Science), according to OpenAlex
Funding: European Research Council (884622)
Citations: not cited yet (Europe PMC); 54 references in the paper

Abstract

Background: Understanding human ageing across multiple organs is essential for characterising individual health trajectories and identifying abnormal ageing processes. Multi-organ imaging provides an opportunity to quantify biological ageing beyond chronological age. The aim of this study is to assess organ-specific and whole-body ageing patterns and their associations with disease and lifestyle factors.

Methods: In this large-scale study, we evaluate biological ageing patterns using 70,000 MRI scans from the UK Biobank and the German National Cohort. We employ 3D ResNet-18 models to predict chronological age from various body regions (brain, heart, liver, spine, lungs, muscle, and intestine) and the whole body. From these predictions, we derive “age gaps” relative to a strictly healthy reference cohort, which enables the identification of accelerated ageing patterns. We then evaluate associations with chronic diseases and lifestyle factors, and a virtual ageing framework was developed to explore counterfactual scenarios by substituting anatomical regions across subjects, quantifying local impacts on global biological age.

Results: Here we show significant associations between detected accelerated ageing and specific chronic diseases, including multiple sclerosis and chronic obstructive pulmonary disease, as well as lifestyle factors such as smoking and physical activity. Virtual substitution of anatomical regions demonstrates that local substitutions can influence global ageing patterns.

Conclusions: This study demonstrates that multi-organ imaging enables the detection of abnormal ageing patterns at both local and global levels. The presented framework provides a foundation for improved risk stratification and supports the development of personalised approaches to health assessment and disease prevention.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repository

Its files are read in the Code ↔ Paper reader above, with 6 matches between paragraphs and lines of code.

tamaramueller/age_biomarkers

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: a772251fc863ac522f143f918919b7659ec45cb3, 21 April 2026
Languages: Python (35)
Size: 47 files, 35 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README, environment (source/environment.yaml)
Not found: license file, CITATION.cff, tests, continuous integration, documentation
Tools: PyTorch (18 files), NumPy (16 files), pandas (11 files), NiBabel (7 files), SciPy (6 files), PyTorch Lightning (5 files), Matplotlib (5 files), seaborn (4 files), scikit-image (2 files), SimpleITK (1 file), statsmodels (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
36 files
At the source: github.com/tamaramueller/age_biomarkers/

Code availability

All files can be found in the corresponding source code repository at: https://github.com/tamaramueller/age_biomarkers/54.

Reproduced under the paper's license (CC BY), from the paper cited above.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 1 repository of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 35 scripts, each with its path and the digest of its content;
  • 6 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

No dataset and no data link were found in the paper.

Data availability

A were created with Biorender, detailed acknowledgement will be provided upon publication. Eligible researchers may access UK Biobank data on https://www.ukbiobank.ac.uk, upon registration. For this study, permission to access and analyse the UK Biobank data was approved under the application 87802. This project was conducted with data (Application No. NAKO-732) from the German National Cohort (NAKO) (https://www.nako.de). The NAKO is funded by the Federal Ministry of Education and Research (BMBF) [project funding reference numbers: 01ER1301A/B/C, 01ER1511D, 01ER1801A/B/C/D and 01ER2301A/B/C], federal states of Germany and the Helmholtz Association, the participating universities and the institutes of the Leibniz Association. We thank all participants who took part in the NAKO study and the staff of this research initiative. Figures 1 and 8A were created with Biorender, detailed acknowledgement will be provided upon publication. Figure 1 was created in BioRender. Mueller, T. (2025) https://BioRender.com/8ryioo7. Figure 8A was created in BioRender. Mueller, T. (2025) https://BioRender.com/8rm2rn2.

Reproduced under the paper's license (CC BY), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 21 authors, 1 keyword, 1 funder, 43 references.

Cite

This paper

Mueller, T. T., Starck, S., Llalloshi, R., Kaissis, G., Ziller, A., Graf, R., Schlett, C., Ringhof, S., Bamberg, F., Wielpütz, M., Völzke, H., Leitzmann, M., Niendorf, T., Keil, T., Krist, L., Pischon, T., Karch, A., Berger, K., Kirschke, J., . . . Braren, R. (2026). Local and global patterns support medical imaging as a biomarker of ageing. Communications medicine, 6(1), 335. https://doi.org/10.1038/s43856-026-01722-3

BibTeX

@article{mueller2026local,
author = {Mueller, Tamara T. and Starck, Sophie and Llalloshi, Rozafë and Kaissis, Georgios and Ziller, Alexander and Graf, Robert and Schlett, Christopher and Ringhof, Steffen and Bamberg, Fabian and Wielpütz, Mark and Völzke, Henry and Leitzmann, Michael and Niendorf, Thoralf and Keil, Thomas and Krist, Lilian and Pischon, Tobias and Karch, André and Berger, Klaus and Kirschke, Jan and Rueckert, Daniel and Braren, Rickmer},
title = {{Local and global patterns support medical imaging as a biomarker of ageing}},
journal = {Communications medicine},
year = {2026},
month = jun,
volume = {6},
number = {1},
pages = {335},
publisher = {Nature Publishing Group},
issn = {2730-664X},
doi = {10.1038/s43856-026-01722-3},
url = {https://doi.org/10.1038/s43856-026-01722-3},
pmid = {42288701},
pmcid = {PMC13264634}
}

RIS

TY - JOUR
AU - Mueller, Tamara T.
AU - Starck, Sophie
AU - Llalloshi, Rozafë
AU - Kaissis, Georgios
AU - Ziller, Alexander
AU - Graf, Robert
AU - Schlett, Christopher
AU - Ringhof, Steffen
AU - Bamberg, Fabian
AU - Wielpütz, Mark
AU - Völzke, Henry
AU - Leitzmann, Michael
AU - Niendorf, Thoralf
AU - Keil, Thomas
AU - Krist, Lilian
AU - Pischon, Tobias
AU - Karch, André
AU - Berger, Klaus
AU - Kirschke, Jan
AU - Rueckert, Daniel
AU - Braren, Rickmer
TI - Local and global patterns support medical imaging as a biomarker of ageing
T2 - Communications medicine
J2 - Commun Med (Lond)
PY - 2026
DA - 2026/06/13
VL - 6
IS - 1
SP - 335
SN - 2730-664X
PB - Nature Publishing Group
DO - 10.1038/s43856-026-01722-3
UR - https://doi.org/10.1038/s43856-026-01722-3
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s43856-026-01722-3",
"type": "article-journal",
"title": "Local and global patterns support medical imaging as a biomarker of ageing",
"container-title": "Communications medicine",
"author": [
{
"family": "Mueller",
"given": "Tamara T."
},
{
"family": "Starck",
"given": "Sophie"
},
{
"family": "Llalloshi",
"given": "Rozafë"
},
{
"family": "Kaissis",
"given": "Georgios"
},
{
"family": "Ziller",
"given": "Alexander"
},
{
"family": "Graf",
"given": "Robert"
},
{
"family": "Schlett",
"given": "Christopher"
},
{
"family": "Ringhof",
"given": "Steffen"
},
{
"family": "Bamberg",
"given": "Fabian"
},
{
"family": "Wielpütz",
"given": "Mark"
},
{
"family": "Völzke",
"given": "Henry"
},
{
"family": "Leitzmann",
"given": "Michael"
},
{
"family": "Niendorf",
"given": "Thoralf"
},
{
"family": "Keil",
"given": "Thomas"
},
{
"family": "Krist",
"given": "Lilian"
},
{
"family": "Pischon",
"given": "Tobias"
},
{
"family": "Karch",
"given": "André"
},
{
"family": "Berger",
"given": "Klaus"
},
{
"family": "Kirschke",
"given": "Jan"
},
{
"family": "Rueckert",
"given": "Daniel"
},
{
"family": "Braren",
"given": "Rickmer"
}
],
"container-title-short": "Commun Med (Lond)",
"volume": "6",
"issue": "1",
"page": "335",
"DOI": "10.1038/s43856-026-01722-3",
"PMID": "42288701",
"PMCID": "PMC13264634",
"ISSN": "2730-664X",
"publisher": "Nature Publishing Group",
"URL": "https://doi.org/10.1038/s43856-026-01722-3",
"language": "en",
"issued": {
"date-parts": [
[
2026,
6,
13
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1002/hbm.70425 [code]
Brain Age Estimation on T2-FLAIR Scans for Application to Multiple Sclerosis.
Journal: Human brain mapping
In common: PyTorch, pandas, SciPy, 1 other tool, clinical / translational, 8 references
[2] doi:10.1186/s40708-026-00316-y [code]
Generalizable and explainable deep learning for brain MRI: a multi-cohort evaluation of 3D architectures for age and sex prediction.
Journal: Brain informatics
In common: PyTorch, seaborn, pandas, 3 other tools, 6 references
[3] doi:10.1162/imag.a.1164 [code]
Bias and generalizability of brain age prediction models: A multi-cohort evaluation with anatomical and interpretability insights.
Journal: Imaging neuroscience (Cambridge, Mass.)
In common: NiBabel, statsmodels, PyTorch, 5 other tools, 4 references
[4] doi:10.1038/s41380-026-03691-4 [code]
Breaking the norm: population-scale deviations of brain structure in depression and anxiety.
Journal: Molecular psychiatry
In common: NiBabel, statsmodels, PyTorch, 5 other tools, clinical / translational, author Tobias Pischon
[5] doi:10.1038/s41398-026-04131-1 [code]
Multimodal phenotypic classification of generalized anxiety and panic using structural MRI data and psychosocial factors: machine learning results from the German National Cohort (NAKO) study.
Journal: Translational psychiatry
In common: clinical / translational, 1 reference, 2 authors
[6] doi:10.1007/s12021-026-09817-x [code]
Circle of Willis-Guided Localization for Simultaneous Detection and Classification of Large Vessel Occlusions in Brain CTA.
Journal: Neuroinformatics
In common: PyTorch Lightning, SimpleITK, scikit-image, 7 other tools
[7] doi:10.1093/radadv/umag025 [code]
OpenMAP-BrainAge: generalizable and interpretable brain age predictor from MRI.
Journal: Radiology advances
In common: SimpleITK, NiBabel, PyTorch, 4 other tools, 3 references
[8] doi:10.3389/frai.2026.1771088 [code]
Few-shot deployment of pretrained MRI transformers in brain imaging tasks.
Journal: Frontiers in artificial intelligence
In common: SimpleITK, scikit-image, NiBabel, 6 other tools, 1 reference
[9] doi:10.1038/s41398-026-04081-8 [code]
Functional system-specific brain aging across the Alzheimer's disease continuum.
Journal: Translational psychiatry
In common: NiBabel, statsmodels, seaborn, 4 other tools, clinical / translational, 3 references
[10] doi:10.1038/s42003-026-10205-z [code]
Source-space EEG alpha activity reveals brain age gaps due to neurodegeneration and disparity.
Journal: Communications biology
In common: statsmodels, seaborn, pandas, 3 other tools, 4 references

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.