OSCR

Stable individual differences dominate adult brain volume variation until later life.

Code ↔ Paper

2 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 2 matches
  1. [1] § Results › Dynamical aging norm curves for brain structures ↔ scripts/create_mcmc_figures.ipynb, lines 42–48 · score 0.74 · brain stem volume, supratentorial volumes, lateral ventricles, Gray matter volume, hippocampal
  2. [2] § Methods › Data ↔ scripts/create_mcmc_figures.ipynb, lines 42–48 · score 0.64 · supratentorial volume, lateral ventricles, hippocampal volume

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

Jupyter notebook · 766 lines · 33 KB · no license · 2 matches

  1. # %%
  2. import torch
  3. import pyro
  4. import numpy as np
  5. import matplotlib.pyplot as plt
  6. import matplotlib as mpl
  7. import matplotlib.ticker as ticker
  10. import os
  11. import sys
  13. module_path = os.path.abspath(os.path.join('..'))
  14. if module_path not in sys.path:
  15. sys.path.append(module_path)
  16. from src.ghmm_mcmc import GaussianHMM as ghmm
  18. #import ghmm_mcmc.GaussianHMM as ghmm
  20. #/ess/p23/home/p23-edvardgr/Code/gaussian-hmm-mcmc-main/src/ghmm_mcmc/FastGaussianHMM.py
  21. from pyro.infer import MCMC, NUTS, Predictive
  22. from pathlib import Path
  24. import warnings
  25. warnings.filterwarnings("ignore")
  28. # %%
  29. #plt.style.use('ggplot')
  30. plt.style.use('seaborn-v0_8-whitegrid')
  31. torch.set_default_dtype(torch.float64)
  33. torch.autograd.set_grad_enabled(mode=False)
  35. #mpl.rcParams['font.family'] = 'sans-serif'
  36. #mpl.rcParams['font.sans-serif'] = ['DejaVu Sans']
  38. out_dir = Path('/tsd/p23/home/p23-edvardgr/file-export') / 'edvard_export'
  42. # %%
  43. pretty_names = {"TotalHippocampus": "Total Hippocampus Volume",
  44. "TotalGrayVol": "Gray Matter Volume",
  45. "TotalWhiteVol": "White Matter Volume",
  46. "SupraTentorialVolNotVent": "Supratentorial Volume",
  47. "TotalLateralVentricle": "Total Lateral Ventricles Volume",
  48. "Brain-Stem": "Brain Stem Volume"}
  50. # %%
  51. def plt_samplesd_line(y, ax, x=None, plot_method='samples', subsample_rate=10, **args):
  52. """
  53. Args:
  54. x - tensor (samples, time)
  55. type - string, 'CI': 95%, and 50% confidence interval with median or 'samples': plot the lines directly
  56. args - plt args
  57. """
  58. with torch.no_grad():
  59. y=y[::subsample_rate]
  60. n_timesteps = y.shape[1]
  61. n_samples = y.shape[0]
  62. if x is None:
  63. x = np.arange(n_timesteps)
  65. args_ci = {}
  66. for key, value in args.items():
  67. if not key in ['label', 'linestyle', 'ls']:
  68. args_ci[key] = value
  70. if plot_method == 'CI':
  71. q = torch.tensor([0.025, 0.25, 0.5, 0.75, 0.975])
  72. quantiles = torch.quantile(y, q, dim=0, keepdim=False)
  73. # Plot first the 95% CI
  74. ax.fill_between(x, quantiles[0], quantiles[4], alpha=0.25, lw=0, **args_ci)
  75. ax.fill_between(x, quantiles[1], quantiles[3], alpha=0.25, lw=0, **args_ci)
  76. ax.plot(x, quantiles[2], **args)
  77. return ax
  78. elif plot_method == 'samples':
  79. n_max = 1000
  80. if n_samples > n_max:
  81. perm = torch.randperm(n_samples)
  82. idx = perm[:n_max]
  83. y = y[idx]
  84. alpha = max(1 / y.shape[0], 0.003)
  85. ax.plot(x, y.T, alpha=alpha, **args_ci)
  86. return ax
  88. # %%
  89. age_interval = 1.0
  91. n_controlpoints=9
  92. n_states = 1
  94. #post_fix = '_years_final_run_uio_ucam' #'_years_final_run' #'_long_icv' # '_strong_priors_smooth_spline2'
  95. post_fix = '_years_final_run'
  96. #features=["prior", "TotalHippocampus", "TotalGrayVol", "TotalWhiteVol", "SupraTentorialVolNotVent", "TotalLateralVentricle", "Brain-Stem"]
  97. features=["TotalHippocampus", "TotalGrayVol", "TotalWhiteVol", "SupraTentorialVolNotVent", "TotalLateralVentricle", "Brain-Stem"]
  98. #features = ["TotalHippocampus", "TotalGrayVol"]
  99. #features=["TotalHippocampus", "TotalGrayVol"]
  100. #features=["TotalHippocampus", "TotalGrayVol", "SupraTentorialVolNotVent", "Brain-Stem"]
  101. #features = features[:2]
  103. #dataset_path = Path('/tsd/p23/home/p23-edvardgr/Storage/p23-edvardgr') /'Didac_Data/processed_data' / f'feature_age_interval_{age_interval}_uio_ucam.npy'
  104. dataset_path = Path('/tsd/p23/home/p23-edvardgr/Storage/p23-edvardgr') /'Didac_Data/processed_data' / f'feature_age_interval_{age_interval}.npy'
  106. dataset = np.load(dataset_path, allow_pickle=True).item()
  108. device = torch.device('cpu')
  110. plot_method = 'CI'
  114. # %%
  115. # Plot the std of the dataset as function of age
  116. fig, axs = plt.subplots(2, 3, figsize=(7, 4))
  117. axs = axs.flatten()
  118. for ax, feature in zip(axs, features):
  119. data = dataset['features'][feature]
  120. #sd = np.nanstd(np.log(data), axis=0) #/np.nanmean(data, axis=0)
  121. sd = np.nanstd(data, axis=0) #/np.nanmean(data, axis=0)
  122. #sd = np.nanmean(data, axis=0)
  123. ax.plot(sd)
  124. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  125. ages_normed = (ages - 18)/age_interval
  126. ax.set_xticks(ages_normed)
  127. ax.set_xticklabels(ages)
  128. ax.set_ylabel("mm\u00B3")
  129. ax.set_xlabel("Age")
  130. ax.set_xlim(((18 - 18)/age_interval, (90 - 18)/age_interval))
  131. ax.set_title(feature)
  132. ax.yaxis.set_major_locator(ticker.MaxNLocator(nbins=4))
  134. fig.suptitle("Standard Deviation of Brain Structure Volumes over the Adult Lifespan")
  136. fig.tight_layout()
  138. fig.savefig(out_dir / 'std_volumes.pdf', dpi=300, bbox_inches = 'tight')
  141. # %%
  142. fig, axs = plt.subplots(2, len(features)//2, figsize=(7, 4))
  143. axs = axs.flatten()
  145. for feature_idx, feature in enumerate(features):
  147. feature_name = feature
  148. if feature == "prior":
  149. feature_name = 'TotalHippocampus'
  151. print(feature, feature_name)
  152. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  153. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  154. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  157. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  159. n_timesteps = dataset['features'][feature_name].shape[-1]
  162. priors = {}
  163. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  164. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  165. priors['z0_corr_concentration'] = 1000.0
  166. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([0.2, 0.01])}},
  167. 'icv': {'type': 'r', 'dim': 1, 'priors': {'loc': torch.tensor([[0.], [0.]]), 'scale': torch.tensor([[1.], [0.01]])}} }
  168. n_sites = dataset['covariates']['site'].shape[-1]
  169. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  171. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  172. n_states, device=device, priors=priors,
  173. z0_mu_covariates=z0_mu_covariates,
  174. y_bias_covariates=y_bias_covariates)
  176. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  177. # Subsample the mcmc_samples
  178. new_samples = {}
  180. for key, value in mcmc_samples.items():
  181. new_samples[key] = value[::10]
  183. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  186. #for age_idx, start_age in enumerate([18, 35, 50, 60, 70]):
  187. # start_age_adjusted = int((start_age-18)/age_interval)
  188. # nac, nac_R = model.get_auto_covariance_coeff(t0 = start_age_adjusted, unconstrained_params=new_samples)
  189. # nac = nac[..., 0, 1].squeeze()
  190. # #nac = nac_R.squeeze()
  191. # plt_samplesd_line(nac[:, start_age_adjusted:], axs[feature_idx], x=np.arange(start_age_adjusted, n_timesteps), plot_method=plot_method, subsample_rate=1, color=colors[age_idx], label=f"{start_age}")
  193. corr_coef = model.get_correlation_coeff(unconstrained_params=new_samples).squeeze()
  194. plt_samplesd_line(corr_coef, axs[feature_idx], plot_method=plot_method, subsample_rate=1, color='black', label="Volume")
  196. Delta_Ts = [1, 2, 4, 8, 16]
  198. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  199. for plot_idx, Delta_T in enumerate(Delta_Ts):
  200. corr_dy_v = model.get_dy_v_covariance(Delta_T=Delta_T, unconstrained_params=new_samples).squeeze()
  201. if Delta_T == 1:
  202. label = f"{Delta_T} Yr Interval"
  203. else:
  204. label = f"{Delta_T} Yrs Interval"
  205. plt_samplesd_line(corr_dy_v, axs[feature_idx], x=np.arange(Delta_T, n_timesteps),
  206. plot_method=plot_method, subsample_rate=1, color=colors[plot_idx], label=label, ls='dotted')
  208. axs[feature_idx].set_ylim((0.8, 1.05))
  209. axs[feature_idx].set_title(f'{pretty_names[feature]}')
  211. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  212. ages_normed = (ages - 18)/age_interval
  213. axs[feature_idx].set_xticks(ages_normed)
  214. axs[feature_idx].set_xticklabels(ages)
  215. axs[feature_idx].set_xlim(((18 - 18)/age_interval, (90 - 18)/age_interval))
  217. y_ticks = np.linspace(0, 1, 6)
  218. axs[feature_idx].set_yticks(y_ticks)
  219. if feature_idx//3 == 1:
  220. axs[feature_idx].set_xlabel('Age')
  221. if feature_idx%3 == 0:
  222. axs[feature_idx].set_ylabel(r'$\rho$')
  223. # Assume all subplots share the same legend entries;
  224. # Get handles and labels from one of the axes.
  225. handles, labels = axs[0].get_legend_handles_labels()
  227. # Create a legend on the figure (centered below all subplots)
  228. fig.legend(handles, labels, loc='lower center', ncol=len(labels)//2,
  229. bbox_to_anchor=(0.5, -0.10))
  231. # Adjust the layout to provide space for the legend
  232. fig.subplots_adjust(bottom=0.25)
  233. fig.suptitle('Correlation of Estimates of Change Rate with Latent Change Rate')
  235. fig.tight_layout()
  237. fig.savefig(out_dir / 'correlation_vol_changerate.pdf', dpi=300, bbox_inches = 'tight')
  239. # %%
  240. fig, axs = plt.subplots(2, len(features)//2, figsize=(7, 4))
  241. axs = axs.flatten()
  243. for feature_idx, feature in enumerate(features):
  245. feature_name = feature
  246. if feature == "prior":
  247. feature_name = 'TotalHippocampus'
  249. print(feature, feature_name)
  250. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  251. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  252. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  255. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  257. n_timesteps = dataset['features'][feature_name].shape[-1]
  260. priors = {}
  261. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  262. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  263. priors['z0_corr_concentration'] = 1000.0
  264. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([0.2, 0.01])}},
  265. 'icv': {'type': 'r', 'dim': 1, 'priors': {'loc': torch.tensor([[0.], [0.]]), 'scale': torch.tensor([[1.], [0.01]])}} }
  266. n_sites = dataset['covariates']['site'].shape[-1]
  267. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  269. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  270. n_states, device=device, priors=priors,
  271. z0_mu_covariates=z0_mu_covariates,
  272. y_bias_covariates=y_bias_covariates)
  274. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  275. # Subsample the mcmc_samples
  276. new_samples = {}
  278. for key, value in mcmc_samples.items():
  279. new_samples[key] = value[::10]
  281. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  284. start_age_adjusted = 0
  285. nac, nac_R = model.get_auto_covariance_coeff(t0 = start_age_adjusted, unconstrained_params=new_samples)
  286. r2 = nac[..., 0, 0].squeeze()**2
  287. #nac = nac_R.squeeze()
  288. plt_samplesd_line(r2[:, start_age_adjusted:], axs[feature_idx],
  289. x=np.arange(start_age_adjusted, n_timesteps),
  290. plot_method=plot_method, subsample_rate=1,
  291. color=colors[0], label=f"Latent volume")
  293. r2_R = nac_R.squeeze()**2
  294. plt_samplesd_line(r2_R[:, start_age_adjusted:], axs[feature_idx],
  295. x=np.arange(start_age_adjusted, n_timesteps),
  296. plot_method=plot_method, subsample_rate=1,
  297. color=colors[1], label=f"Measured volume")
  299. print_ages = np.array([30, 40, 50, 60, 70, 80, 90])
  300. print(feature)
  301. print('R2 latent:')
  302. for print_age in print_ages:
  303. print(f"{print_age}: {torch.mean(r2[:, int((print_age-18)/age_interval)]).item():.4f} - ", end='')
  304. print('R2 measured:')
  305. for print_age in print_ages:
  306. print(f"{print_age}: {torch.mean(r2_R[:, int((print_age-18)/age_interval)]).item():.4f} - ", end='')
  308. axs[feature_idx].set_ylim((0.8, 1.05))
  309. axs[feature_idx].set_title(f'{pretty_names[feature]}')
  311. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  312. ages_normed = (ages - 18)/age_interval
  313. axs[feature_idx].set_xticks(ages_normed)
  314. axs[feature_idx].set_xticklabels(ages)
  315. axs[feature_idx].set_xlim(((18 - 18)/age_interval, (90 - 18)/age_interval))
  317. y_ticks = np.linspace(0, 1, 6)
  318. axs[feature_idx].set_yticks(y_ticks)
  320. if feature_idx//3 == 1:
  321. axs[feature_idx].set_xlabel('Age')
  322. if feature_idx%3 == 0:
  323. axs[feature_idx].set_ylabel('R\u00B2')
  325. # Assume all subplots share the same legend entries;
  326. # Get handles and labels from one of the axes.
  327. handles, labels = axs[0].get_legend_handles_labels()
  329. # Create a legend on the figure (centered below all subplots)
  330. fig.legend(handles, labels, loc='lower center', ncol=len(labels), bbox_to_anchor=(0.5, -0.05))
  332. # Adjust the layout to provide space for the legend
  333. fig.subplots_adjust(bottom=0.25)
  334. fig.suptitle('Explained Variance from Volume Levels at Age 18')
  336. fig.tight_layout()
  339. fig.savefig(out_dir / 'explained_variance.pdf', dpi=300, bbox_inches = 'tight')
  341. # %%
  342. fig, axs = plt.subplots(2, len(features)//2, figsize=(7, 4))
  343. axs = axs.flatten()
  345. for feature_idx, feature in enumerate(features):
  347. feature_name = feature
  348. if feature == "prior":
  349. feature_name = 'TotalHippocampus'
  351. print(feature, feature_name)
  352. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  353. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  354. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  357. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  359. n_timesteps = dataset['features'][feature_name].shape[-1]
  362. priors = {}
  363. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  364. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  365. priors['z0_corr_concentration'] = 1000.0
  366. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([0.2, 0.01])}},
  367. 'icv': {'type': 'r', 'dim': 1, 'priors': {'loc': torch.tensor([[0.], [0.]]), 'scale': torch.tensor([[1.], [0.01]])}} }
  368. n_sites = dataset['covariates']['site'].shape[-1]
  369. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  371. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  372. n_states, device=device, priors=priors,
  373. z0_mu_covariates=z0_mu_covariates,
  374. y_bias_covariates=y_bias_covariates)
  376. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  377. # Subsample the mcmc_samples
  378. new_samples = {}
  380. for key, value in mcmc_samples.items():
  381. new_samples[key] = value[::10]
  383. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  385. for age_idx, start_age in enumerate([18, 35, 50, 60, 70]):
  386. start_age_adjusted = int((start_age-18)/age_interval)
  387. nac, nac_R = model.get_auto_covariance_coeff(t0 = start_age_adjusted, unconstrained_params=new_samples)
  388. nac = nac[..., 1, 1].squeeze()
  389. plt_samplesd_line(nac[:, start_age_adjusted:], axs[feature_idx], x=np.arange(start_age_adjusted, n_timesteps), plot_method=plot_method, subsample_rate=1, color=colors[age_idx], label=f"{start_age}")
  392. axs[feature_idx].set_ylim((0.8, 1.05))
  393. axs[feature_idx].set_title(f'{pretty_names[feature]}')
  395. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  396. ages_normed = (ages - 18)/age_interval
  397. axs[feature_idx].set_xticks(ages_normed)
  398. axs[feature_idx].set_xticklabels(ages)
  399. axs[feature_idx].set_xlim(((18 - 18)/age_interval, (90 - 18)/age_interval))
  401. y_ticks = np.linspace(0, 1, 6)
  402. axs[feature_idx].set_yticks(y_ticks)
  403. if feature_idx//3 == 1:
  404. axs[feature_idx].set_xlabel('Age')
  405. if feature_idx%3 == 0:
  406. axs[feature_idx].set_ylabel(r"$\rho_{AC}$")
  407. # Assume all subplots share the same legend entries;
  408. # Get handles and labels from one of the axes.
  409. handles, labels = axs[0].get_legend_handles_labels()
  411. # Create a legend on the figure (centered below all subplots)
  412. fig.legend(handles, labels, loc='lower center', ncol=len(labels), bbox_to_anchor=(0.5, -0.05))
  414. # Adjust the layout to provide space for the legend
  415. fig.subplots_adjust(bottom=0.25)
  416. fig.suptitle('Autocorrelation of Latent Change Rate')
  418. fig.tight_layout()
  421. fig.savefig(out_dir / 'correlation_changerate.pdf', dpi=300, bbox_inches = 'tight')
  423. # %%
  424. for key in mcmc_samples.keys():
  425. print(f"'{key}': '',")
  427. # %%
  428. #fig, axs_features = plt.subplots(3, len(features), figsize=(len(features)/6*20, 5), sharex=True)
  429. fig = plt.figure(figsize=(8, 20))
  430. rows = 1
  431. subfigs = fig.subfigures(rows, int(np.ceil(len(features)/rows))).flatten()
  433. pretty_coef_names = {'R_sigma': r'$\sigma_R$',
  434. 'b_t_coeffs': r'$b(t)$',
  435. 'sigma_p_coeffs': r'$\sigma_q(t)$',
  436. 'y_bias.cat_coeffs.site.params': r'$\beta(site)$',
  437. 'z0_mu.bias': r'$\beta_{sv0}$',
  438. 'z0_mu.cat_coeffs.sex.params': r'$\beta_{sv0,sex}$',
  439. 'z0_mu.regressors.icv': r'$\beta_{sv0,ICV}$',
  440. 'z0_sigma': r'$\sigma_{sv0}$'}
  442. for feature_idx, feature in enumerate(features):
  444. feature_name = feature
  445. if feature == "prior":
  446. feature_name = 'TotalHippocampus'
  448. print(feature, feature_name)
  450. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  451. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  452. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  455. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  457. n_timesteps = dataset['features'][feature_name].shape[-1]
  460. priors = {}
  461. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  462. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  463. priors['z0_corr_concentration'] = 1000.0
  464. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([1., 0.01])}}}
  465. n_sites = dataset['covariates']['site'].shape[-1]
  466. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  468. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  469. n_states, device=device, priors=priors,
  470. z0_mu_covariates=z0_mu_covariates,
  471. y_bias_covariates=y_bias_covariates)
  473. model.norm_mean = dataset['norm_means'][feature_name]
  474. model.norm_std = dataset['norm_std'][feature_name]
  477. if feature == 'prior':
  478. mcmc_samples = Predictive(model, {}, num_samples=100, parallel=True)(y[[0]])
  479. else:
  480. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  482. n_params = len(list(pretty_coef_names.keys()))
  483. n_params = 0
  484. for idx, (param_name, pretty_coef_name) in enumerate(pretty_coef_names.items()):
  485. params = mcmc_samples[param_name]
  486. params = params.reshape(params.shape[0], -1).cpu().detach().numpy()
  487. n_params += params.shape[1]
  488. print(n_params)
  489. axs = subfigs[feature_idx].subplots(n_params, 1, sharex=True)
  490. #subfigs[feature_idx].subplots_adjust(hspace=0.4, wspace=0.4)
  491. subfigs[feature_idx].suptitle(pretty_names[feature], y=0.91, wrap=True)
  492. param_idx = 0
  493. for idx, (param_name, pretty_coef_name) in enumerate(pretty_coef_names.items()):
  494. params = mcmc_samples[param_name]
  495. params = params.reshape(params.shape[0], -1).cpu().detach().numpy()
  496. axs[param_idx].text(0, 1, pretty_coef_name, verticalalignment='bottom', horizontalalignment='right', transform=axs[param_idx].transAxes, rotation=90, rotation_mode='anchor')
  497. for sub_param_idx in range(params.shape[1]):
  498. axs[param_idx].plot(params[:, sub_param_idx], lw=0.5)
  499. axs[param_idx].axis('off')
  500. param_idx += 1
  501. axs[param_idx-1].plot([0, 1], [0, 0], color='black', transform=axs[param_idx-1].transAxes)
  503. #fig.tight_layout()
  504. #plt.show()
  506. fig.savefig(out_dir / 'trace_plots.pdf', dpi=300, bbox_inches = 'tight')
  508. # %%
  509. #fig, axs_features = plt.subplots(3, len(features), figsize=(len(features)/6*20, 5), sharex=True)
  510. fig = plt.figure(figsize=(7, 7), dpi=300)
  511. rows = 2
  512. cols = int(np.ceil(len(features)/rows))
  513. subfigs = fig.subfigures(rows, cols).flatten()
  515. for feature_idx, feature in enumerate(features):
  517. feature_name = feature
  518. if feature == "prior":
  519. feature_name = 'TotalHippocampus'
  521. print(feature, feature_name)
  523. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  524. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  525. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  528. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  530. n_timesteps = dataset['features'][feature_name].shape[-1]
  532. priors = {}
  533. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  534. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  535. priors['z0_corr_concentration'] = 1000.0
  536. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([0.2, 0.01])}},
  537. 'icv': {'type': 'r', 'dim': 1, 'priors': {'loc': torch.tensor([[0.], [0.]]), 'scale': torch.tensor([[1.], [0.01]])}} }
  539. n_sites = dataset['covariates']['site'].shape[-1]
  540. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  542. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  543. n_states, device=device, priors=priors,
  544. z0_mu_covariates=z0_mu_covariates,
  545. y_bias_covariates=y_bias_covariates)
  547. model.norm_mean = dataset['norm_means'][feature_name]
  548. model.norm_std = dataset['norm_std'][feature_name]
  551. if feature == 'prior':
  552. mcmc_samples = Predictive(model, {}, num_samples=100, parallel=True)(y[[0]])
  553. else:
  554. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  556. mu_t_g_t, Sigma_t_g_t, y_mu, y_Sigma, b_t, Q_t, _ = model.get_posterior_distribution(unconstrained_params=mcmc_samples)
  558. print(f'mu_t_g_t: {mu_t_g_t.shape}')
  559. print(f'Sigma_t_g_t: {Sigma_t_g_t.shape}')
  560. print(f'b_t: {b_t.shape}')
  561. print(f'Q_t: {Q_t.shape}')
  563. age_range = 100-18
  564. if not ('prior' in features):
  565. mu_t_g_t[..., 0] = mu_t_g_t[..., 0]*model.norm_std + model.norm_mean
  566. mu_t_g_t[..., 1] = mu_t_g_t[..., 1]*model.norm_std/(age_range)
  567. Sigma_t_g_t[..., 0, 0] = Sigma_t_g_t[..., 0, 0]*model.norm_std**2
  568. Sigma_t_g_t[..., 1, 1] = Sigma_t_g_t[..., 1, 1]*model.norm_std**2/(age_range**2)
  570. b_t = b_t*model.norm_std/(age_range**2)
  571. Q_t = Q_t*model.norm_std**2/(age_range**4)
  573. #axs = axs_features[:, feature_idx]
  574. axs = subfigs[feature_idx].subplots(2, 1)
  575. subfigs[feature_idx].subplots_adjust(hspace=0.4, wspace=0.2)
  576. subfigs[feature_idx].suptitle(pretty_names[feature])
  577. ax_s = axs[0]
  578. ax_v = axs[1]
  579. #ax_a = axs[2]
  581. if not ('prior' in features):
  582. for i in range(y.shape[0]):
  583. # Get the indices where the observations are not missing
  584. nan_idx = np.isnan(y[i, :, 0].detach().cpu().numpy())
  585. t = np.arange(y.shape[1])[~nan_idx]
  586. ax_s.plot(t, y[i, ~nan_idx, 0].detach().cpu().numpy(), '.-', alpha=0.02, markersize=2, color='black', rasterized=True)
  587. alpha = 0.003
  588. for idx in range(mu_t_g_t.shape[-1]):
  589. mean = mu_t_g_t[:, ..., idx].detach().squeeze().cpu()
  590. std = torch.sqrt(Sigma_t_g_t[:, ..., idx, idx].detach().squeeze().cpu())
  591. if idx == 0:
  592. Sigma_icv = (mcmc_samples['z0_mu.regressors.icv'].squeeze()[..., 0]**2 * 1.0)[:,None]
  593. std = torch.sqrt(std**2 + Sigma_icv)
  594. #axs[idx+1].plot(mean, label=f"State {idx}", color='blue', alpha=alpha)
  595. plt_samplesd_line(mean, axs[idx], plot_method=plot_method, color='blue')
  596. plt_samplesd_line(mean + 2*std, axs[idx], plot_method=plot_method, color='red')
  597. plt_samplesd_line(mean - 2*std, axs[idx], plot_method=plot_method, color='red')
  599. mean = mu_t_g_t[:, ..., :int((90-18)/age_interval), idx].detach().squeeze().cpu()
  600. std = torch.sqrt(Sigma_t_g_t[:, ..., :int((90-18)/age_interval), idx, idx].detach().squeeze().cpu())
  601. min_value = torch.min(torch.mean(mean, dim=0) - 2*std)
  602. max_value = torch.max(torch.mean(mean, dim=0) + 2*std)
  603. margin = (max_value - min_value)*0.075
  604. y_lim = (min_value - margin, max_value + margin)
  605. axs[idx].set_ylim(y_lim)
  607. # Plot the posterior of the acceleration
  608. #mean = b_t[:, ..., 0].detach().squeeze().cpu()
  609. #std = torch.sqrt(Q_t[:, ..., 0].detach().squeeze().cpu())
  611. #mean = b_t[:, ..., :int((90-18)/age_interval), 0].detach().squeeze().cpu()
  612. #std = torch.sqrt(Q_t[:, ..., :int((90-18)/age_interval), 0].detach().squeeze().cpu())
  616. #min_value = torch.min(torch.mean(mean, dim=0) - 2*std)
  617. #max_value = torch.max(torch.mean(mean, dim=0) + 2*std)
  618. #margin = (max_value - min_value)*0.075
  619. #y_lim = (min_value - margin, max_value + margin)
  620. #ax_a.set_ylim(y_lim)
  622. #plt_samplesd_line(mean, ax_a, plot_method=plot_method, color='blue')
  623. #plt_samplesd_line(mean + 2*std, ax_a, plot_method=plot_method, color='red')
  624. #plt_samplesd_line(mean - 2*std, ax_a, plot_method=plot_method, color='red')
  626. #ax_s.set_title('Volume')
  627. ax_v.set_title("Change Rate")
  628. #ax_a.set_title("Second derviative")
  630. if feature_idx%cols==0:
  631. ax_s.set_ylabel("mm\u00B3 ")
  632. ax_v.set_ylabel("mm\u00B3/year")
  633. #ax_a.set_ylabel("mm\u00B3/year\u00B2")
  635. ax_v.set_xlabel("Age")
  637. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  638. ages_normed = (ages - 18)/age_interval
  640. ax_s.set_xticks(ages_normed)
  641. ax_s.set_xticklabels(ages)
  642. ax_v.set_xticks(ages_normed)
  643. ax_v.set_xticklabels(ages)
  645. ax_s.ticklabel_format(style='sci', useOffset=False, axis='y', scilimits=(3, 3))
  646. ax_v.ticklabel_format(style='sci', useOffset=False, axis='y', scilimits=(3, 3))
  647. #ax_a.ticklabel_format(style='sci', useOffset=False, axis='y', scilimits=(3, 3))
  649. ax_s.yaxis.set_major_locator(ticker.MaxNLocator(nbins=6))
  650. ax_v.yaxis.set_major_locator(ticker.MaxNLocator(nbins=6))
  651. #ax_a.yaxis.set_major_locator(ticker.MaxNLocator(nbins=6))
  653. ax_s.set_xlim((0, (90-18)/age_interval))
  654. ax_v.set_xlim((0, (90-18)/age_interval))
  656. #fig.tight_layout()
  657. #plt.show()
  659. fig.savefig(out_dir / 'norm_curves.pdf', dpi=300, bbox_inches = 'tight')
  661. # %%
  662. fig, axs = plt.subplots(2, len(features)//2, figsize=(7, 4))
  663. axs = axs.flatten()
  665. for feature_idx, feature in enumerate(features):
  667. feature_name = feature
  668. if feature == "prior":
  669. feature_name = 'TotalHippocampus'
  671. print(feature, feature_name)
  672. results_dir = Path('/tsd/p23/home/p23-edvardgr/Code/gaussian-hmm-job-scripts/model_fitter_test/Results') / f'{feature_name}_age_int_{age_interval:.1f}{post_fix}'
  673. mcmc_samples_path = results_dir / 'mcmc_samples.pht'
  674. #meta_data_dict = torch.load(results_dir / 'meta_dict.pht', map_location=torch.device('cpu'), weights_only=False)
  677. y = torch.tensor(dataset['features'][feature_name][: ,:, np.newaxis], device=device, dtype=torch.float64)
  679. n_timesteps = dataset['features'][feature_name].shape[-1]
  682. priors = {}
  683. priors['R_sigma'] = {'loc': dataset['measurment_error'][feature_name]/dataset['norm_std'][feature_name],
  684. 'scale': 0.5*dataset['measurment_error'][feature_name] /dataset['norm_std'][feature_name]}
  685. priors['z0_corr_concentration'] = 1000.0
  686. z0_mu_covariates={'sex': {'type': 'c', 'dim': 2, 'priors': {'loc': torch.tensor([0., 0.]), 'scale': torch.tensor([0.2, 0.01])}},
  687. 'icv': {'type': 'r', 'dim': 1, 'priors': {'loc': torch.tensor([[0.], [0.]]), 'scale': torch.tensor([[1.], [0.01]])}} }
  688. n_sites = dataset['covariates']['site'].shape[-1]
  689. y_bias_covariates={'site': {'type': 'c', 'dim': n_sites, 'priors': {'loc': torch.tensor([0.]), 'scale': torch.tensor([0.1])}}}
  693. mcmc_samples = torch.load(mcmc_samples_path, map_location=torch.device('cpu'), weights_only=False)
  694. # Subsample the mcmc_samples
  695. new_samples = {}
  697. for key, value in mcmc_samples.items():
  698. new_samples[key] = value[::300]
  700. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  702. model = ghmm.GaussianStateModel(n_timesteps, n_controlpoints,
  703. n_states, device=device, priors=priors,
  704. z0_mu_covariates=z0_mu_covariates,
  705. y_bias_covariates=y_bias_covariates)
  707. #corr_coef = model.get_correlation_coeff(unconstrained_params=new_samples).squeeze()
  708. #plt_samplesd_line(corr_coef, axs[feature_idx], plot_method=plot_method, subsample_rate=1, color='black', label="$NAC_{s_k v_k}$")
  710. n_samples = [1, 5, 11]
  712. linestyles = ['dotted', (0, (1, 1)), (0, (0.5, 0.5))]
  713. age_sample_intervals = [1, 2, 4]
  714. colors = plt.rcParams['axes.prop_cycle'].by_key()['color']
  715. for age_int_idx, age_sample_interval in enumerate(age_sample_intervals):
  716. for plot_idx, n_sample in enumerate(n_samples):
  717. new_n_timesteps = int(n_timesteps*n_sample/age_sample_interval)
  718. model = ghmm.GaussianStateModel(new_n_timesteps, n_controlpoints,
  719. n_states, device=device, priors=priors,
  720. z0_mu_covariates=z0_mu_covariates,
  721. y_bias_covariates=y_bias_covariates)
  723. corr_dy_v = model.get_alpha_v_covariance(n=n_sample, unconstrained_params=new_samples).squeeze()
  724. if age_sample_interval == 1:
  725. label = f"{age_sample_interval} Yr Interval - {n_sample+1} Samples"
  726. else:
  727. label = f"{age_sample_interval} Yrs Interval - {n_sample+1} Samples"
  728. #age80_normalized = np.round((80-18)/(age_interval*age_sample_interval)*n_sample-n_sample)
  729. #age80_corr = torch.mean(corr_dy_v[:, age80_normalized]).item()
  730. #axs[feature_idx].plot((age80_normalized + n_sample)/n_sample*age_sample_interval , age80_corr, 'x', color='black')
  731. #print(f"Interval {age_sample_interval} years, N samples {n_sample +1}, corr age 80 {age80_corr}")
  734. plt_samplesd_line(corr_dy_v, axs[feature_idx], x=np.arange(n_sample+1, new_n_timesteps)/n_sample*age_sample_interval,
  735. plot_method=plot_method, subsample_rate=1, color=colors[age_int_idx], label=label, ls=linestyles[plot_idx])
  737. #axs[feature_idx].set_ylim((0.8, 1.05))
  738. axs[feature_idx].set_title(f'{pretty_names[feature]}')
  740. ages = np.array([20, 30, 40, 50, 60, 70, 80, 90])
  741. ages_normed = (ages - 18)/age_interval
  742. axs[feature_idx].set_xticks(ages_normed)
  743. axs[feature_idx].set_xticklabels(ages)
  744. axs[feature_idx].set_xlim(((18 - 18)/age_interval, (90 - 18)/age_interval))
  746. y_ticks = np.linspace(0, 1, 6)
  747. axs[feature_idx].set_yticks(y_ticks)
  748. if feature_idx//3 == 1:
  749. axs[feature_idx].set_xlabel('Age')
  750. if feature_idx%3 == 0:
  751. axs[feature_idx].set_ylabel(r'\rho')
  753. # Assume all subplots share the same legend entries;
  754. # Get handles and labels from one of the axes.
  755. handles, labels = axs[0].get_legend_handles_labels()
  757. # Create a legend on the figure (centered below all subplots)
  758. fig.legend(handles, labels, loc='lower center', ncol=3, bbox_to_anchor=(0.5, -0.15))
  760. # Adjust the layout to provide space for the legend
  761. fig.subplots_adjust(bottom=0.25)
  762. fig.suptitle('Correlation of Estimates of Change Rate with Latent Change Rate Using Two or More Time Points')
  764. fig.tight_layout()
  766. fig.savefig(out_dir / 'num_sessions_changerate.pdf', dpi=300, bbox_inches = 'tight')

create_mcmc_figures.ipynb at commit 49d09d4, no license · at the source

Overview

Authors: Edvard OS Grødem1,2, Didac Vidal-Pineiro1, Øystein Sørensen1, David Bartrés-Faz3,4,5, Andreas M Brandmaier6,7,8, Gabriele Cattaneo4,9, Pablo F Garrido1, Richard N Henson10, Simone Kühn8,11,12, Ulman Lindenberger6,8, Bradley J MacIntosh2, Lars Nyberg1,13,14,15, Alvaro Pascual-Leone16,17, Stephen M Smith18, Cristina Solé-Padullés3,5, Javier Solana-Sánchez4,9, Leiv Otto Watne19,20, Kristine B Walhovd1,2, Atle Bjørnerud2, Anders M Fjell1,2
ORCID iDs: Edvard OS Grødem, Øystein Sørensen, David Bartrés-Faz, Andreas M Brandmaier, Gabriele Cattaneo, Pablo F Garrido, Richard N Henson, Simone Kühn, Ulman Lindenberger, Bradley J MacIntosh, Lars Nyberg, Alvaro Pascual-Leone, Stephen M Smith, Cristina Solé-Padullés, Javier Solana-Sánchez, Kristine B Walhovd, Atle Bjørnerud, Anders M Fjell
20 affiliations
  1. Center for Lifespan Changes in Brain and Cognition (LCBC), Department of Psychology, University of Oslo, Oslo, Norway
  2. Computational Radiology and Artificial Intelligence (CRAI), Department of Radiology and Nuclear Medicine, Oslo University Hospital, Oslo, Norway
  3. Department of Medicine, Faculty of Medicine and Health Sciences, Institute of Neurosciences, University of Barcelona, Barcelona, Spain
  4. Institut Guttmann, Institut Universitari de Neurorehabilitació adscrit a la UAB, Badalona, Spain
  5. The August Pi i Sunyer Biomedical Research Institute (IDIBAPS), Hospital Clinic of Barcelona, Barcelona, Spain
  6. Center for Lifespan Psychology, Max Planck Institute for Human Development, Berlin, Germany
  7. Department of Psychology, MSB Medical School Berlin, Berlin, Germany
  8. Max Planck UCL Centre for Computational Psychiatry and Ageing Research, Berlin, Germany, and London, United Kingdom
  9. Fundació Institut d’Investigació en Ciències de la Salut Germans Trias i Pujol, Badalona, Spain
  10. MRC Cognition and Brain Sciences Unit, Department of Psychiatry, University of Cambridge, Cambridge, United Kingdom
  11. Center for Environmental Neuroscience, Max Planck Institute for Human Development, Berlin, Germany
  12. Department of Psychiatry and Psychotherapy, University Medical Center Hamburg-Eppendorf, Hamburg, Germany
  13. Umeå Center for Functional Brain Imaging, Umeå University, Umeå, Sweden
  14. Department of Medical and Translational Biology, Umeå University, Umeå, Sweden
  15. Department of Diagnostics and Intervention, Umeå University, Umeå, Sweden
  16. Hinda and Arthur Marcus Institute for Aging Research, Deanna and Sidney Wolk Center for Memory Health, Hebrew SeniorLife, Boston, MA, United States
  17. Department of Neurology, Harvard Medical School, Boston, MA, United States
  18. Oxford Centre for Functional MRI of the Brain (FMRIB/WIN), Oxford University, Oxford, United Kingdom
  19. Oslo Delirium Research Group, Institute of Clinical Medicine, Campus Ahus, University of Oslo, Oslo, Norway
  20. Department of Geriatric Medicine, Akershus University Hospital, Oslo, Norway
Institutions: Oslo University Hospital (Norway); University of Oslo (Norway); Universitat Autònoma de Barcelona (Spain); Institut Guttmann (Spain); Hospital Clínic de Barcelona (Spain); Consorci Institut D'Investigacions Biomediques August Pi I Sunyer (Spain); Universitat de Barcelona (Spain); MSB Medical School Berlin (Germany); Max Planck Institute for Human Development (Germany); University College London (United Kingdom); Institut d'Investigació en Ciències de la Salut Germans Trias i Pujol (Spain); MRC Cognition and Brain Sciences Unit (United Kingdom); University of Cambridge (United Kingdom); University Medical Center Hamburg-Eppendorf (Germany); Umeå University (Sweden); Harvard University (United States); Hebrew SeniorLife (United States); University of Oxford (United Kingdom); Wellcome Centre for Integrative Neuroimaging (United Kingdom); Akershus University Hospital (Norway)
Journal: Imaging neuroscience (Cambridge, Mass.), volume 4, article IMAG.a.1242
Dates: received 12 August 2025; accepted 13 April 2026; published online 26 May 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1162/imag.a.1242 · PMID 42212223 · PMCID PMC13214574 · OpenAlex W7158467521
Open access: diamond, a free copy (OpenAlex)
Status: code verified
Categories: structural MRI / diffusion (modality), human (organism), computational (subfield)
Methods: Smoothing, state filtering, decompositions
Keywords: stochastic dynamical system, brain aging, magnetic resonance imaging, methods, MRI, brain age, BAG
Topic: Functional Brain Connectivity Studies (Cognitive Neuroscience, Neuroscience), according to OpenAlex
Funding: European Research Council (283634, 313440, 725025)
Citations: not cited yet (Europe PMC); 106 references in the paper

Abstract

Individual differences in the volumes of brain structures are often linked to various conditions, including Alzheimer’s disease, schizophrenia, and overall brain health. However, it remains unclear to what extent these differences reflect individual levels present from young adulthood or diverging aging trajectories from later ages. In this study, we analyze the aging dynamics of the volumes of six brain structures based on magnetic resonance imaging (MRI) scans from a large cross-cohort longitudinal sample of cognitively healthy adults (n = 8,311 with 18,520 MRIs, ages from 18 to 97 years). From general assumptions about structural brain dynamics and measurement noise, a stochastic dynamical model was fitted to the data to estimate both the variability and persistence of structural changes across adulthood. Using this model, we calculated how much of the variance of volumetric differences between individuals can be attributed to stable levels from young adulthood versus systematic changes at older ages, as well as the theoretical sensitivity of longitudinal studies to detect individual differences in change. The findings were as follows: (1) Before age 60 years, inter-individual differences in neuroanatomical volumes almost exclusively reflect stable differences between individuals, while the influence from systematic differences in rate-of-change increases thereafter: up to 50% of the variation being due to differences in change at 80 years. In contrast, ventricular volume reflects differences in change from early adulthood. (2) Current brain-age models are unlikely to be sensitive to detect differences in aging trajectories. (3) Imaging studies have low reliability in detecting inter-individual brain changes before age 60 years. After 60 years, the study reliability increases sharply with longer intervals between scans and more modestly with additional intermediate observations. In conclusion, our results reinforce the view that it is critical to distinguish stable early adulthood levels from systematic differences in change when studying adult brain aging.

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 2 matches between paragraphs and lines of code.

EdvardGrodem/brain-trajectories

License: none: the authors keep all their rights
State: the link answers, verified on 28 September 2026
Evidence: files inventoried
Commit: 49d09d46b844859a100428ba04f9ef801b3bd210, 11 August 2025
Languages: Jupyter (5), Python (4)
Size: 11 files, 9 scripts
Software Heritage: not archived
Found in: “Data and Code Availability”
Holds: README, environment (setup.py), 5 notebooks
Not found: license file, CITATION.cff, tests, continuous integration, documentation
Tools: Matplotlib (8 files), NumPy (8 files), Pyro (6 files), PyTorch (6 files), seaborn (3 files), pandas (2 files), SciPy (2 files)
Availability: 1 check, the latest on 28 September 2026: the link answers
  • 28 September 2026: the link answers
10 files
At the source: github.com/EdvardGrodem/brain-trajectories

The paper's code and data availability statement is in the Data section.

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;
  • 9 scripts, each with its path and the digest of its content;
  • 2 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 and Code Availability

The raw data were gathered from many different datasets. Different agreements are required for each dataset. Most datasets are openly available with prespecified data usage agreements. For some datasets, such as UKB, fees may apply. Requests for LCBC, UB, and COGNORM should be submitted to the corresponding principal investigator. The code to fit the model and generate the figures is available at https://github.com/EdvardGrodem/brain-trajectories.

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, 28 September 2026: the first record

Recorded: type, language, journal, volume, pages, dates, 20 authors, 7 keywords, 1 funder, 103 references.

Cite

This paper

Grødem, E. O., Vidal-Pineiro, D., Sørensen, Ø., Bartrés-Faz, D., Brandmaier, A. M., Cattaneo, G., Garrido, P. F., Henson, R. N., Kühn, S., Lindenberger, U., MacIntosh, B. J., Nyberg, L., Pascual-Leone, A., Smith, S. M., Solé-Padullés, C., Solana-Sánchez, J., Watne, L. O., Walhovd, K. B., Bjørnerud, A., & Fjell, A. M. (2026). Stable individual differences dominate adult brain volume variation until later life. Imaging neuroscience (Cambridge, Mass.), 4, IMAG.a.1242. https://doi.org/10.1162/imag.a.1242

BibTeX

@article{grdem2026stable,
author = {Grødem, Edvard OS and Vidal-Pineiro, Didac and Sørensen, Øystein and Bartrés-Faz, David and Brandmaier, Andreas M and Cattaneo, Gabriele and Garrido, Pablo F and Henson, Richard N and Kühn, Simone and Lindenberger, Ulman and MacIntosh, Bradley J and Nyberg, Lars and Pascual-Leone, Alvaro and Smith, Stephen M and Solé-Padullés, Cristina and Solana-Sánchez, Javier and Watne, Leiv Otto and Walhovd, Kristine B and Bjørnerud, Atle and Fjell, Anders M},
title = {{Stable individual differences dominate adult brain volume variation until later life}},
journal = {Imaging neuroscience (Cambridge, Mass.)},
year = {2026},
month = may,
volume = {4},
pages = {IMAG.a.1242},
publisher = {MIT Press},
issn = {2837-6056},
doi = {10.1162/imag.a.1242},
url = {https://doi.org/10.1162/imag.a.1242},
pmid = {42212223},
pmcid = {PMC13214574}
}

RIS

TY - JOUR
AU - Grødem, Edvard OS
AU - Vidal-Pineiro, Didac
AU - Sørensen, Øystein
AU - Bartrés-Faz, David
AU - Brandmaier, Andreas M
AU - Cattaneo, Gabriele
AU - Garrido, Pablo F
AU - Henson, Richard N
AU - Kühn, Simone
AU - Lindenberger, Ulman
AU - MacIntosh, Bradley J
AU - Nyberg, Lars
AU - Pascual-Leone, Alvaro
AU - Smith, Stephen M
AU - Solé-Padullés, Cristina
AU - Solana-Sánchez, Javier
AU - Watne, Leiv Otto
AU - Walhovd, Kristine B
AU - Bjørnerud, Atle
AU - Fjell, Anders M
TI - Stable individual differences dominate adult brain volume variation until later life
T2 - Imaging neuroscience (Cambridge, Mass.)
J2 - Imaging Neurosci (Camb)
PY - 2026
DA - 2026/05/26
VL - 4
SP - IMAG.a.1242
SN - 2837-6056
PB - MIT Press
DO - 10.1162/imag.a.1242
UR - https://doi.org/10.1162/imag.a.1242
LA - en
ER -

CSL-JSON

{
"id": "10.1162/imag.a.1242",
"type": "article-journal",
"title": "Stable individual differences dominate adult brain volume variation until later life",
"container-title": "Imaging neuroscience (Cambridge, Mass.)",
"author": [
{
"family": "Grødem",
"given": "Edvard OS"
},
{
"family": "Vidal-Pineiro",
"given": "Didac"
},
{
"family": "Sørensen",
"given": "Øystein"
},
{
"family": "Bartrés-Faz",
"given": "David"
},
{
"family": "Brandmaier",
"given": "Andreas M"
},
{
"family": "Cattaneo",
"given": "Gabriele"
},
{
"family": "Garrido",
"given": "Pablo F"
},
{
"family": "Henson",
"given": "Richard N"
},
{
"family": "Kühn",
"given": "Simone"
},
{
"family": "Lindenberger",
"given": "Ulman"
},
{
"family": "MacIntosh",
"given": "Bradley J"
},
{
"family": "Nyberg",
"given": "Lars"
},
{
"family": "Pascual-Leone",
"given": "Alvaro"
},
{
"family": "Smith",
"given": "Stephen M"
},
{
"family": "Solé-Padullés",
"given": "Cristina"
},
{
"family": "Solana-Sánchez",
"given": "Javier"
},
{
"family": "Watne",
"given": "Leiv Otto"
},
{
"family": "Walhovd",
"given": "Kristine B"
},
{
"family": "Bjørnerud",
"given": "Atle"
},
{
"family": "Fjell",
"given": "Anders M"
}
],
"container-title-short": "Imaging Neurosci (Camb)",
"volume": "4",
"page": "IMAG.a.1242",
"DOI": "10.1162/imag.a.1242",
"PMID": "42212223",
"PMCID": "PMC13214574",
"ISSN": "2837-6056",
"publisher": "MIT Press",
"URL": "https://doi.org/10.1162/imag.a.1242",
"language": "en",
"issued": {
"date-parts": [
[
2026,
5,
26
]
]
}
}

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.1038/s41593-026-02363-4 [code]
Cortical thickness changes precede high levels of amyloid by at least 7 years.
Journal: Nature neuroscience
In common: structural MRI / diffusion, 2 references, 4 authors
[2] doi:10.1162/imag.a.1306 [code]
Beneficial effects of foreign language learning and aerobic exercise on dentate gyrus volume and mnemonic discrimination in healthy older adults: Results from a randomized controlled trial.
Journal: Imaging neuroscience (Cambridge, Mass.)
In common: 1 reference, 3 authors
[3] doi:10.1162/imag.a.1196 [code]
Failure to detect entorhinal grid-like signals in a passive navigation human fMRI study.
Journal: Imaging neuroscience (Cambridge, Mass.)
In common: 3 authors
[4] doi:10.1038/s41467-026-73668-y [code]
Convergent and divergent brain-cognition development in early adolescence.
Journal: Nature communications
In common: PyTorch, seaborn, pandas, 3 other tools, 5 references
[5] doi:10.1038/s41467-026-71271-9 [code]
Exposome-wide patterns predict brain health in aging.
Journal: Nature communications
In common: seaborn, pandas, SciPy, 2 other tools, 5 references
[6] 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, structural MRI / diffusion, 4 references
[7] doi:10.1038/s41467-026-74957-2 [code]
Joint trajectories of brain atrophy, white matter hyperintensities and cognition quantify brain maintenance.
Journal: Nature communications
In common: structural MRI / diffusion, 6 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: PyTorch, seaborn, pandas, 3 other tools, structural MRI / diffusion, 3 references
[9] doi:10.64898/2026.08.13.26360304 [code]
Lifespan brain structural variation reveals shared organization across mental health conditions
Journal: medRxiv (preprint)
In common: seaborn, pandas, SciPy, 2 other tools, author Bradley J. MacIntosh
[10] doi:10.7554/elife.104053 [code]
Brain cognition gaps reveal associations with dopamine and factors related to brain health through artificial intelligence prediction of functional connectome.
Journal: eLife
In common: pandas, SciPy, NumPy, 1 reference, author Lars Nyberg

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.