OSCR

Peripheral blood microarray-based transcriptomic and epigenetic analyses identify immune, inflammation, and metabolic dysregulation in Alzheimer's disease.

Code ↔ Paper

16 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 16 matches
  1. [1] § Methods › DNA methylation profiling ↔ R/DNAm-preprocessing.R, lines 107–150 · score 0.85 · bisulfite conversion rates, rmSNPandCH, cross hybridizing, sex chromosomes, genomic, technical
  2. [2] § Methods › Functional enrichment and visualization ↔ R/differential_analyses.R, lines 723–852 · score 0.81 · normalized enrichment score, MSigDB, fold change, NES, BH FDR, fgsea
  3. [3] § Methods › Functional enrichment and visualization ↔ R/MOFA.R, lines 361–502 · score 0.81 · normalized enrichment score, MSigDB, NES, ranked, fgsea, ORA
  4. [4] § Methods › Amyloid PET and MRI-based atrophy measures ↔ R/phenotype.R, lines 538–597 · score 0.74 · cubic age, meta ROI cortical, cortical thickness, Hippocampal volumes, ICV, ADNI
  5. [5] § Methods › Gene expression profiling ↔ R/phenotype.R, lines 599–637 · score 0.69 · Affymetrix Human Genome, U219, chip, Gene expression, ADNI, probes
  6. [6] § Methods › Network and multi-omics analyses ↔ R/WGCNA.R, lines 845–940 · score 0.69 · principal component, module eigengene, limma, linear, fit, matrices
  7. [7] § Methods › Survival analysis for MCI to AD progression ↔ R/WGCNA.R, lines 1418–1563 · score 0.60 · module assignments, original modules, MEs, AD progression, predictors, covariates
  8. [8] § Methods › DNA methylation profiling ↔ R/phenotype.R, lines 1131–1173 · score 0.58 · ADNI Genetics Core, DNA Methylation, QC
  9. [9] § Results › Co-methylation network analysis reveals immune-related dysregulation in AD ↔ R/wgcna2cytoscape.R, lines 151–274 · score 0.56 · co methylation network, CpGs, Hub gene, WGCNA, modules
  10. [10] § Results › Gene co-expression network analysis (WGCNA) identifies robust immune networks ↔ R/WGCNA.R, lines 608–633 · score 0.55 · S13a, S13c, WGCNA, women, gene
  11. [11] § Methods › Network and multi-omics analyses ↔ R/MOFA.R, lines 35–89 · score 0.55 · MOFA model, limma, trained, variance, linear, matched
  12. [12] § Methods › DNA methylation profiling ↔ R/phenotype.R, lines 1393–1471 · score 0.55 · immune cell, DNA methylation, EpiDISH, blood, profiling
  13. [13] § Results › Participant characteristics ↔ R/tables.R, lines 18–85 · score 0.55 · amyloid PET positive, hippocampal volume, S1a, carriers, variables, plasma
  14. [14] § Results › Gene co-expression network analysis (WGCNA) identifies robust immune networks ↔ R/WGCNA.R, lines 446–485 · score 0.54 · module membership, gene co expression, gene modules, MM, WGCNA, network
  15. [15] § Methods › Sensitivity analyses ↔ R/tables.R, lines 18–85 · score 0.54 · amyloid positive, carrier status, amyloid negative, amyloid PET, MCI, CU
  16. [16] § Results › AD-related gene modules are linked to AD endophenotypes ↔ R/phenotype.R, lines 236–281 · score 0.53 · mini mental state, examination, cognitive, MMSE

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

R · 1,471 lines · 65 KB · no license · 5 matches

  1. library(SummarizedExperiment)
  2. library(tidyverse)
  3. library(lme4)
  4. library(EpiDISH)
  6. # Define Paths & Functions ####
  7. base_dir <- ".../data"
  8. proc_dir <- ".../processed"
  9. idat_dir <- ".../IDAT"
  11. addZeros <- function(x){
  12. sprintf("%04d", as.integer(x))
  13. }
  15. # Define function to match tables at same time or within 1-year of gene expression profiling
  16. gene_within1yr <- function(df1, df2_Edate, df2_VISCODE = FALSE) {
  17. df1 <- df1
  19. # If table contains both VISCODE and EXAMDATE, use both to determine closest match
  20. if (df2_VISCODE){
  21. df1 <- df1 %>%
  22. mutate(date_diff = abs(lubridate::time_length(difftime(GeneX_Edate, {{ df2_Edate }}), "years"))) %>%
  23. mutate(Edate_match = ifelse(date_diff <= 1, 1, 0),
  24. VISCODE_match = ifelse(VISCODE == VISCODE2, 1, 0)) %>%
  25. group_by(RID) %>%
  26. mutate(min_match_date = min(date_diff)) %>%
  27. # Same time = VISCODE_match; within 1-year = Edate_match
  28. mutate(match_priority = case_when(
  29. VISCODE_match == 1 ~ 1,
  30. Edate_match == 1 ~ 2,
  31. TRUE ~ 3
  32. )) %>%
  33. # Only keep matches
  34. filter(match_priority <= 2) %>%
  35. # If multiple matches, keep min date diff
  36. filter(min_match_date == date_diff) %>%
  37. ungroup() %>%
  38. distinct(RID, VISCODE, .keep_all = T) %>%
  39. dplyr::select(-c(Edate_match, VISCODE_match, match_priority, date_diff,
  40. min_match_date, VISCODE2))
  41. } else{
  42. # If table is missing VISCODE, use EXAMDATE only
  43. df1 <- df1 %>%
  44. mutate(date_diff = abs(lubridate::time_length(difftime(GeneX_Edate, {{ df2_Edate }}), "years"))) %>%
  45. mutate(Edate_match = ifelse(date_diff <= 1, 1, 0)) %>%
  46. group_by(RID) %>%
  47. mutate(min_match_date = min(date_diff)) %>%
  48. filter(Edate_match == 1) %>%
  49. filter(min_match_date == date_diff) %>%
  50. ungroup() %>%
  51. distinct(RID, VISCODE, .keep_all = T) %>%
  52. dplyr::select(-c(Edate_match, date_diff, min_match_date))
  53. }
  54. }
  56. # Define function to match tables at same time or within 1-year of DNA methylation profiling
  57. dnam_within1yr <- function(df1, df2_Edate, df2_VISCODE = FALSE) {
  58. df1 <- df1
  60. # If table contains both VISCODE and EXAMDATE, use both to determine closest match
  61. if (df2_VISCODE){
  62. df1 <- df1 %>%
  63. mutate(date_diff = abs(lubridate::time_length(difftime(DNAm_Edate, {{ df2_Edate }}), "years"))) %>%
  64. mutate(Edate_match = ifelse(date_diff <= 1, 1, 0),
  65. VISCODE_match = ifelse(VISCODE == VISCODE2, 1, 0)) %>%
  66. group_by(RID, DNAm_Edate) %>%
  67. mutate(min_match_date = min(date_diff)) %>%
  68. # Same time = VISCODE_match; within 1-year = Edate_match
  69. mutate(match_priority = case_when(
  70. VISCODE_match == 1 ~ 1,
  71. Edate_match == 1 ~ 2,
  72. TRUE ~ 3
  73. )) %>%
  74. # Only keep matches
  75. filter(match_priority <= 2) %>%
  76. # If multiple matches, keep min date diff
  77. filter(min_match_date == date_diff) %>%
  78. ungroup() %>%
  79. distinct(RID, barcodes, .keep_all = T) %>%
  80. dplyr::select(-c(Edate_match, VISCODE_match, match_priority, date_diff,
  81. min_match_date, VISCODE2))
  82. } else{
  83. # If table is missing VISCODE, use EXAMDATE only
  84. df1 <- df1 %>%
  85. mutate(date_diff = abs(lubridate::time_length(difftime(DNAm_Edate, {{ df2_Edate }}), "years"))) %>%
  86. mutate(Edate_match = ifelse(date_diff <= 1, 1, 0)) %>%
  87. group_by(RID, DNAm_Edate) %>%
  88. mutate(min_match_date = min(date_diff)) %>%
  89. filter(Edate_match == 1) %>%
  90. filter(min_match_date == date_diff) %>%
  91. ungroup() %>%
  92. distinct(RID, barcodes, .keep_all = T) %>%
  93. dplyr::select(-c(Edate_match, date_diff, min_match_date))
  94. }
  95. }
  97. #______________________________________Prepare ADNI Tables Prior to Merging With -Omics Tables______________________________________#####
  98. message("Now reading in and preparing ADNI tables prior to merging with -omics metadata...")
  99. ## Demographics Table ####
  100. # Read in and prepare demographics table
  101. demog <- read.csv(file.path(base_dir, "All_Subjects_PTDEMOG_28Oct2024.csv"), stringsAsFactors = FALSE)
  102. adni_merge <- read.csv(file.path(base_dir, "ADNIMERGE.csv"), stringsAsFactors = FALSE)
  103. adni_merge_unique <- adni_merge %>% distinct_at(vars(RID), .keep_all = TRUE)
  104. demog_unique <- demog %>% distinct_at(vars(RID), .keep_all = TRUE)
  106. demog_merged <- left_join(adni_merge_unique, demog_unique %>%
  107. dplyr::select(RID, PTDOB, PTDOBYY), by = "RID")
  109. demog_merged <- demog_merged %>%
  110. dplyr::select(RID, ORIGPROT, AGE, PTDOB, PTGENDER,
  111. PTEDUCAT, PTRACCAT, PTETHCAT) %>%
  112. dplyr::rename(AGE_bl = AGE,
  113. PTDOBMMYY = PTDOB)
  114. demog_merged <- demog_merged %>%
  115. mutate(day = 1,
  116. PTDOB = as.character(paste(day, PTDOBMMYY, sep = "/")))
  117. demog_merged$PTDOB <- as.Date(demog_merged$PTDOB, format = "%d/%m/%Y")
  118. demog_merged <- demog_merged %>%
  119. dplyr::select(-c(PTDOBMMYY, day))
  120. demog_merged <- relocate(demog_merged, "PTDOB", .before = "AGE_bl")
  121. demog_merged$RID <- addZeros(demog_merged$RID)
  123. ## Diagnosis Table ####
  124. # Read in and prepare diagnosis table
  125. dx <- read.csv(file.path(base_dir, "All_Subjects_DXSUM_22Nov2024.csv"), stringsAsFactors = FALSE)
  126. dx <- dx %>%
  127. dplyr::select(RID, DIAGNOSIS, EXAMDATE, VISCODE2) %>%
  128. dplyr::rename(DX_Edate = EXAMDATE,
  129. VISCODE = VISCODE2) %>%
  130. mutate(DXSUM = case_when(
  131. (DIAGNOSIS == 1) ~ "CU",
  132. (DIAGNOSIS == 2) ~ "MCI",
  133. (DIAGNOSIS == 3) ~ "Dementia"),
  134. VISCODE = case_when(
  135. VISCODE == "sc" ~ "bl",
  136. TRUE ~ VISCODE)
  137. ) %>% filter(!is.na(DXSUM))
  138. dx$RID <- addZeros(dx$RID)
  139. dx <- dx %>%
  140. arrange(RID, DX_Edate) %>%
  141. group_by(RID) %>%
  142. mutate(DX_last = DXSUM[[length(DXSUM)]]) %>%
  143. ungroup() %>%
  144. mutate(DX_Edate = as.Date(DX_Edate, "%Y-%m-%d")) %>%
  145. dplyr::select(RID, VISCODE, DX_Edate, DXSUM, DX_last)
  147. ## Registry Table ####
  148. # Read in and prepare registry table
  149. registry <- read.csv(file.path(base_dir, "REGISTRY_28Jul2023.csv"), stringsAsFactor = FALSE)
  150. registry$RID <- addZeros(registry$RID)
  151. registry <- registry %>% dplyr::select(RID, VISCODE, VISCODE2)
  152. registry$VISCODE <- ifelse(registry$VISCODE == "v03", "bl", registry$VISCODE)
  153. registry$VISCODE <- ifelse(registry$VISCODE == "v04", "m03", registry$VISCODE)
  154. registry$VISCODE <- ifelse(registry$VISCODE == "v05", "m06", registry$VISCODE)
  155. registry <- registry %>% distinct_at(vars(RID, VISCODE), .keep_all = TRUE)
  157. ## Neurofilament Light Chain (NfL) Table ####
  158. # Read in and prepare NfL table
  159. nfl_data <- read.csv(file.path(base_dir, "ADNI_BLENNOWPLASMANFLLONG_10_03_18_22Apr2025.csv"), stringsAsFactor = FALSE, check.names = FALSE)
  161. # Recode with correct VISCODE2
  162. nfl_data$VISCODE2[which(nfl_data$RID=="1097"&nfl_data$EXAMDATE=="2013-01-31")]<-"m72"
  164. nfl_data <- nfl_data %>%
  165. group_by(RID,EXAMDATE) %>%
  166. mutate(PLASMA_NFL = mean(PLASMA_NFL)) %>%
  167. dplyr::distinct_at(vars(RID,EXAMDATE), .keep_all = TRUE) %>%
  168. ungroup()
  170. nfl_data <- nfl_data %>%
  171. mutate(VISCODE2 = case_when(
  172. VISCODE2 == "sc" ~ "bl",
  173. TRUE ~ VISCODE2
  174. )) %>%
  175. dplyr::select(RID, VISCODE2, EXAMDATE, PLASMA_NFL) %>%
  176. dplyr::rename(NfL_Edate = EXAMDATE, NfL_B = PLASMA_NFL)
  177. nfl_data$NfL_Edate <- as.Date(nfl_data$NfL_Edate, "%Y-%m-%d")
  178. nfl_data$RID <- addZeros(nfl_data$RID)
  180. ## Estimated CSF AB/TAU-Positivity Onset Table ####
  181. # Read in and prepare estimated CSF AB/TAU positivity onset table
  182. csf_age <- read.csv(file.path(base_dir, "estimated_csf_ages.csv"), stringsAsFactor = FALSE, check.names = FALSE)
  183. colnames(csf_age)[2] <- "AB_AGE"
  184. colnames(csf_age)[3] <- "TAU_AGE"
  185. csf_age$AB_AGE <- ifelse(csf_age$AB_AGE == "NaN", NA, csf_age$AB_AGE)
  186. csf_age$TAU_AGE <- ifelse(csf_age$TAU_AGE == "NaN", NA, csf_age$TAU_AGE)
  187. csf_age$RID <- addZeros(csf_age$RID)
  188. csf_age <- csf_age[!duplicated(csf_age$RID),]
  190. ## APOE Genotypes Table ####
  191. # Read in and prepare APOE genotypes table
  192. apoeres <- read.csv(file.path(base_dir, "APOERES_24Aug2023.csv"), stringsAsFactor = FALSE, check.names = FALSE)
  193. apoeres$APOE4 <- as.numeric(apply(apoeres[, c('APGEN1', 'APGEN2')]==4, 1, sum)>=1)
  194. apoeres$APOE2 <- as.numeric(apply(apoeres[, c('APGEN1', 'APGEN2')]==2, 1, sum)>=1)
  195. apoeres$APOE3 <- as.numeric(apply(apoeres[, c('APGEN1', 'APGEN2')]==3, 1, sum)>=1)
  197. apoeres$APOE4_n <- apply(apoeres[, c('APGEN1', 'APGEN2')], 1, function(row) {
  198. if (any(row == 4)) {
  199. if (all(row == 4)) {
  200. return(2)
  201. } else {
  202. return(1)
  203. }
  204. } else {
  205. return(0)
  206. }
  207. })
  209. apoeres$APOE2_n <- apply(apoeres[, c('APGEN1', 'APGEN2')], 1, function(row) {
  210. if (any(row == 2)) {
  211. if (all(row == 2)) {
  212. return(2)
  213. } else {
  214. return(1)
  215. }
  216. } else {
  217. return(0)
  218. }
  219. })
  221. apoeres$APOE3_n <- apply(apoeres[, c('APGEN1', 'APGEN2')], 1, function(row) {
  222. if (any(row == 3)) {
  223. if (all(row == 3)) {
  224. return(2)
  225. } else {
  226. return(1)
  227. }
  228. } else {
  229. return(0)
  230. }
  231. })
  233. apoeres_df <- apoeres[c("RID", "APOE4", "APOE2", "APOE3", "APOE4_n", "APOE2_n", "APOE3_n")]
  234. apoeres_df$RID <- addZeros(apoeres_df$RID)
  236. ## Cognitive Performance Tables ####
  237. # Read in and prepare cognitive performance tables
  238. # Clinical Dementia Rating (CDR)
  239. cdr <- read.csv(file.path(base_dir, "All_Subjects_CDR_22Nov2024.csv"), stringsAsFactors = FALSE)
  240. cdr <- cdr[cdr$CDGLOBAL>=0,]
  241. cdr$CDRSB <- rowSums(cdr[,c('CDMEMORY', 'CDORIENT', 'CDJUDGE', 'CDCOMMUN', 'CDHOME', 'CDCARE')])
  242. cdr <- cdr %>%
  243. mutate(VISCODE2 = case_when(
  244. VISCODE2 == "sc" ~ "bl",
  245. TRUE ~ VISCODE2
  246. )) %>%
  247. dplyr::select(RID, VISCODE2, VISDATE, CDRSB) %>%
  248. dplyr::rename(VISCODE = VISCODE2)
  249. cdr$RID <- addZeros(cdr$RID)
  250. cdr$CDR_Edate <- as.Date(cdr$VISDATE, "%Y-%m-%d")
  252. # 13-Item Alzheimer's Disease Assessment Scale–Cognitive Subscale (ADAS13)
  253. adas <- read.csv(file.path(base_dir, "All_Subjects_ADAS_ADNIGO23_22Nov2024.csv"), stringsAsFactors = FALSE)
  254. adas <- adas %>%
  255. dplyr::select(RID, VISCODE2, TOTAL13) %>%
  256. dplyr::rename(VISCODE = VISCODE2) %>%
  257. dplyr::rename(ADAS13 = TOTAL13) %>%
  258. group_by(RID, VISCODE) %>%
  259. mutate(across(.cols = dplyr::where(is.numeric), .fns=~mean(.,na.rm = TRUE))) %>%
  260. distinct_at(vars(RID), .keep_all = TRUE)
  261. adas$RID <- addZeros(adas$RID)
  262. adas$ADAS13 <- ifelse(adas$ADAS13 == "NaN", NA, adas$ADAS13)
  264. # Mini-Mental State Examination (MMSE)
  265. mmse <- read.csv(file.path(base_dir, "All_Subjects_MMSE_22Nov2024.csv"), stringsAsFactors = FALSE)
  266. mmse <- mmse %>%
  267. dplyr::select(RID, VISCODE2, VISDATE, MMSCORE) %>%
  268. dplyr::rename(MMSE = MMSCORE) %>%
  269. mutate(VISCODE2 = case_when(
  270. VISCODE2 == "sc" ~ "bl",
  271. TRUE ~ VISCODE2
  272. )) %>%
  273. group_by(RID, VISCODE2) %>%
  274. mutate(across(.cols = dplyr::where(is.numeric), .fns=~mean(., na.rm = TRUE))) %>%
  275. distinct_at(vars(RID), .keep_all = TRUE) %>%
  276. ungroup()
  278. mmse$MMSE_Edate <- as.Date(mmse$VISDATE, "%Y-%m-%d")
  279. mmse <- mmse %>% dplyr::select(-VISDATE)
  280. mmse$MMSE <- ifelse(mmse$MMSE == "NaN", NA, mmse$MMSE)
  281. mmse$RID <- addZeros(mmse$RID)
  283. ## CSF Biomarkers Table ####
  284. # Read in and prepare CSF biomarkers table
  285. upenn_csf <- read.csv(file.path(base_dir, "All_Subjects_UPENNBIOMK_ROCHE_ELECSYS_23Nov2024.csv"), stringsAsFactors = FALSE)
  286. upenn_csf <- upenn_csf %>%
  287. arrange(desc(RUNDATE)) %>%
  288. distinct_at(vars(RID, VISCODE2), .keep_all = TRUE) #keep most recent replicate
  289. upenn_csf <- upenn_csf %>%
  290. dplyr::rename(ABETA40_csf = ABETA40,
  291. ABETA42_csf = ABETA42,
  292. tTAU_csf = TAU,
  293. pTAU181_csf = PTAU)
  294. upenn_csf <- upenn_csf %>%
  295. mutate(ABETA42_csf = round(ABETA42_csf),
  296. ptau_pos_csf = pTAU181_csf > 24,
  297. amyloid_pos_csf = ABETA42_csf < 980,
  298. ptau_ab_ratio_csf = (as.numeric(pTAU181_csf) / as.numeric(ABETA42_csf)),
  299. ad_pathology_pos_csf = ptau_ab_ratio_csf > 0.025)
  300. upenn_csf <- upenn_csf %>%
  301. mutate(across(.cols=c(ABETA42_csf, ABETA40_csf, tTAU_csf,pTAU181_csf, ptau_ab_ratio_csf), .fns=as.numeric))
  302. upenn_csf <- upenn_csf %>%
  303. mutate(across(.cols=c(ptau_pos_csf, amyloid_pos_csf, ad_pathology_pos_csf), .fns=as.numeric))
  305. # reassign records with VISCODE2 = "UNK" to correct VISCODE2
  306. upenn_csf$VISCODE2[upenn_csf$RID==89 & upenn_csf$VISCODE2=="UNK"]<-"m162"
  307. upenn_csf$VISCODE2[upenn_csf$RID==232 & upenn_csf$VISCODE2=="UNK"]<-"m36"
  308. upenn_csf$VISCODE2[upenn_csf$RID==726 & upenn_csf$VISCODE2=="UNK"]<-"m24"
  309. upenn_csf$VISCODE2[upenn_csf$RID==790 & upenn_csf$VISCODE2=="UNK"]<-"m18"
  310. upenn_csf$VISCODE2[upenn_csf$RID==1200 & upenn_csf$VISCODE2=="UNK"]<-"m36"
  311. upenn_csf$VISCODE2[upenn_csf$RID==2373 & upenn_csf$VISCODE2=="UNK"]<-"m120"
  312. upenn_csf$VISCODE2[upenn_csf$RID==4216 & upenn_csf$VISCODE2=="UNK"]<-"m102"
  313. upenn_csf$VISCODE2[upenn_csf$RID==4393 & upenn_csf$VISCODE2=="UNK"]<-"m96"
  314. upenn_csf$VISCODE2[upenn_csf$RID==4643 & upenn_csf$VISCODE2=="UNK"]<-"m72"
  315. upenn_csf$VISCODE2[upenn_csf$RID==5140 & upenn_csf$VISCODE2=="UNK"]<-"m84"
  316. upenn_csf$VISCODE2[upenn_csf$RID==6161 & upenn_csf$VISCODE2=="UNK"]<-"bl"
  317. upenn_csf$VISCODE2[upenn_csf$RID==6880 & upenn_csf$VISCODE2=="UNK"]<-"bl"
  318. upenn_csf$VISCODE2[upenn_csf$RID==6906 & upenn_csf$VISCODE2=="UNK"]<-"bl"
  320. upenn_csf <- upenn_csf %>%
  321. dplyr::select(RID, VISCODE2, EXAMDATE, ABETA40_csf, ABETA42_csf, tTAU_csf, pTAU181_csf,
  322. ptau_pos_csf, amyloid_pos_csf, ptau_ab_ratio_csf, ad_pathology_pos_csf) %>%
  323. dplyr::rename(VISCODE = VISCODE2,
  324. CSF_Edate = EXAMDATE)
  325. upenn_csf$CSF_Edate <- as.Date(upenn_csf$CSF_Edate, "%Y-%m-%d")
  326. upenn_csf$RID <- addZeros(upenn_csf$RID)
  328. # Identify the first/last negative/positive dates of CSF amyloid
  329. last_negative_date_csf <- upenn_csf %>%
  330. group_by(RID) %>%
  331. filter(amyloid_pos_csf == 0) %>%
  332. mutate(last_AB_neg_date = max(CSF_Edate)) %>%
  333. dplyr::select(RID, last_AB_neg_date) %>%
  334. ungroup() %>%
  335. distinct()
  337. first_positive_date_csf <- upenn_csf %>%
  338. group_by(RID) %>%
  339. filter(amyloid_pos_csf == 1) %>%
  340. mutate(first_AB_pos_date = min(CSF_Edate)) %>%
  341. dplyr::select(RID, first_AB_pos_date) %>%
  342. ungroup() %>%
  343. distinct()
  345. last_negative_date2_csf <- upenn_csf %>%
  346. group_by(RID) %>%
  347. filter(ptau_pos_csf == 0) %>%
  348. mutate(last_TAU_neg_date2 = max(CSF_Edate)) %>%
  349. dplyr::select(RID, last_TAU_neg_date2) %>%
  350. ungroup() %>%
  351. distinct()
  353. first_positive_date2_csf <- upenn_csf %>%
  354. group_by(RID) %>%
  355. filter(ptau_pos_csf == 1) %>%
  356. mutate(first_TAU_pos_date2 = min(CSF_Edate)) %>%
  357. dplyr::select(RID, first_TAU_pos_date2) %>%
  358. ungroup() %>%
  359. distinct()
  361. last_negative_date3_csf <- upenn_csf %>%
  362. group_by(RID) %>%
  363. filter(ad_pathology_pos_csf == 0) %>%
  364. mutate(last_AD_neg_date3 = max(CSF_Edate)) %>%
  365. dplyr::select(RID, last_AD_neg_date3) %>%
  366. ungroup() %>%
  367. distinct()
  369. first_positive_date3_csf <- upenn_csf %>%
  370. group_by(RID) %>%
  371. filter(ad_pathology_pos_csf == 1) %>%
  372. mutate(first_AD_pos_date3 = min(CSF_Edate)) %>%
  373. dplyr::select(RID, first_AD_pos_date3) %>%
  374. ungroup() %>%
  375. distinct()
  377. ## Amyloid PET Table ####
  378. # Read in and prepare amyloid PET table
  379. ab_pet <- read.csv(file.path(base_dir, "All_Subjects_UCBERKELEY_AMY_6MM_23Nov2024.csv"), stringsAsFactors = FALSE)
  380. ab_pet$RID <- addZeros(ab_pet$RID)
  381. ab_pet <- ab_pet %>%
  382. dplyr::rename(PET_Edate = SCANDATE,
  383. amyloid_pos_pet_cross = AMYLOID_STATUS,
  384. amyloid_pos_pet_long = AMYLOID_STATUS_COMPOSITE_REF) %>%
  385. mutate(amyloid_pos_pet_centiloid = case_when(CENTILOIDS > 22 ~ 1,
  386. TRUE ~ 0)) %>%
  387. dplyr::select(RID, PET_Edate, IMAGE_RESOLUTION, TRACER, SUMMARY_SUVR, CENTILOIDS, amyloid_pos_pet_cross,
  388. amyloid_pos_pet_long, amyloid_pos_pet_centiloid)
  389. ab_pet$PET_Edate <- as.Date(ab_pet$PET_Edate, "%Y-%m-%d")
  391. # Identify the first/last negative/positive dates of amyloid PET
  392. last_negative_date_pet <- ab_pet %>%
  393. group_by(RID) %>%
  394. filter(amyloid_pos_pet_cross == 0) %>%
  395. mutate(last_AB_neg_date = max(PET_Edate)) %>%
  396. dplyr::select(RID, last_AB_neg_date) %>%
  397. ungroup() %>%
  398. distinct()
  400. first_positive_date_pet <- ab_pet %>%
  401. group_by(RID) %>%
  402. filter(amyloid_pos_pet_cross == 1) %>%
  403. mutate(first_AB_pos_date = min(PET_Edate)) %>%
  404. dplyr::select(RID, first_AB_pos_date) %>%
  405. ungroup() %>%
  406. distinct()
  408. last_negative_date2_pet <- ab_pet %>%
  409. group_by(RID) %>%
  410. filter(amyloid_pos_pet_long == 0) %>%
  411. mutate(last_AB_neg_date2 = max(PET_Edate)) %>%
  412. dplyr::select(RID, last_AB_neg_date2) %>%
  413. ungroup() %>%
  414. distinct()
  416. first_positive_date2_pet <- ab_pet %>%
  417. group_by(RID) %>%
  418. filter(amyloid_pos_pet_long == 1) %>%
  419. mutate(first_AB_pos_date2 = min(PET_Edate)) %>%
  420. dplyr::select(RID, first_AB_pos_date2) %>%
  421. ungroup() %>%
  422. distinct()
  424. last_negative_date3_pet <- ab_pet %>%
  425. group_by(RID) %>%
  426. filter(amyloid_pos_pet_centiloid == 0) %>%
  427. mutate(last_AB_neg_date3 = max(PET_Edate)) %>%
  428. dplyr::select(RID, last_AB_neg_date3) %>%
  429. ungroup() %>%
  430. distinct()
  432. first_positive_date3_pet <- ab_pet %>%
  433. group_by(RID) %>%
  434. filter(amyloid_pos_pet_centiloid == 1) %>%
  435. mutate(first_AB_pos_date3 = min(PET_Edate)) %>%
  436. dplyr::select(RID, first_AB_pos_date3) %>%
  437. ungroup() %>%
  438. distinct()
  440. ## FreeSurfer MRI Table ####
  441. # Get amyloid negative RIDs that were consistently negative for every scan they had
  442. ab_negs <- ab_pet %>%
  443. group_by(RID) %>%
  444. filter(amyloid_pos_pet_cross == 0) %>%
  445. mutate(first_a_neg_date_pet = min(PET_Edate)) %>%
  446. dplyr::select(RID, first_a_neg_date_pet) %>%
  447. ungroup() %>%
  448. distinct()
  450. # # Read in and prepare cross-sectional FreeSurfer MRI table
  451. fs_cs <- read.csv(file.path(base_dir, "UCSF_FS_CS_5.1_2019.csv"), stringsAsFactors = FALSE)
  452. fs_cs <- fs_cs %>%
  453. filter(OVERALLQC == "Pass") %>%
  454. dplyr::select(-VISCODE) %>%
  455. dplyr::rename(ICV = ST10CV,
  456. VISCODE = VISCODE2) %>%
  457. group_by(RID, VISCODE) %>%
  458. mutate(n = n()) %>%
  459. filter(
  460. # Keep all single-scan visits
  461. n == 1 |
  462. # For multiple scans: keep Non-Accelerated if available
  463. (n > 1 & IMAGETYPE == "Non-Accelerated T1") |
  464. # If no Non-Accelerated, then keep Accelerated
  465. (n > 1 & !any(IMAGETYPE == "Non-Accelerated T1") & IMAGETYPE == "Accelerated T1")
  466. ) %>%
  467. ungroup() %>%
  468. dplyr::select(-n)
  470. fs_cs$MRIcs_Edate <- as.Date(fs_cs$EXAMDATE, "%Y-%m-%d")
  471. fs_cs$RID <- addZeros(fs_cs$RID)
  472. fs_cs <- fs_cs %>% dplyr::select(-EXAMDATE)
  473. fs_cs <- relocate(fs_cs, "MRIcs_Edate", .after = "VISCODE")
  475. # Get demographics and closest diagnosis information
  476. fs_cs <- merge(fs_cs, demog_merged, by = "RID", all.x = TRUE) %>%
  477. mutate(age = round(lubridate::time_length(difftime(MRIcs_Edate, PTDOB), "years"), digits = 1))
  478. fs_cs <- merge(fs_cs, dx %>%
  479. dplyr::select(-VISCODE) %>%
  480. filter(!is.na(DX_Edate)), by = "RID") %>%
  481. mutate(date_diff = abs(lubridate::time_length(difftime(MRIcs_Edate, DX_Edate), "years"))) %>%
  482. group_by(RID, VISCODE) %>%
  483. mutate(min_match_date = min(date_diff)) %>%
  484. ungroup() %>%
  485. filter(min_match_date == date_diff) %>%
  486. distinct()
  488. # For CU only, get earliest amyloid negative scan per RID
  489. fs_cs_temp <- fs_cs %>%
  490. filter(DXSUM == "CU")
  491. fs_cs_temp <- fs_cs_temp %>%
  492. left_join(ab_negs, by = "RID") %>%
  493. mutate(diff_time = abs(lubridate::time_length(difftime(MRIcs_Edate, first_a_neg_date_pet), "years"))) %>%
  494. group_by(RID) %>%
  495. filter(diff_time == min(diff_time)) %>%
  496. ungroup()
  497. fs_cs_temp <- fs_cs_temp %>%
  498. group_by(RID) %>%
  499. filter(row_number() == 1) %>%
  500. ungroup()
  502. # ICV adjustment for ROI volumes
  503. # Create list of volumes that need to be ICV-adjusted
  504. volumes_cs <- fs_cs_temp %>%
  505. dplyr::select(ICV | contains("CV") | contains("SV"))
  506. volumes_cs <- volumes_cs %>%
  507. dplyr::select(-ICV)
  509. # Perform ICV adjustment across all ROI volumes
  510. results_cs <- c()
  512. for (roi_num in 1:length(names(volumes_cs))){
  513. current_roi <- names(volumes_cs)[roi_num]
  514. results_cs <- append(results_cs, current_roi)
  515. }
  517. for (i in results_cs){
  518. fit <- lm(fs_cs_temp[[i]] ~ ICV, data = fs_cs_temp, na.action = na.exclude)
  519. fs_cs[paste0(i,'_ICV')] <- fs_cs[paste0(i)] - (fs_cs$ICV - mean(fs_cs_temp$ICV, na.rm = T)) * fit$coefficients[2]
  520. }
  522. fs_cs <- fs_cs %>%
  523. dplyr::select(-names(volumes_cs))
  525. # Cubic age adjustment for ROI volumes (captures normative aging)
  526. # Obtain list of variables that need to be age-adjusted
  527. # Add "_ICV" suffix to CV and SV in temp to match
  528. colnames(fs_cs_temp) <- ifelse(grepl("CV|SV", colnames(fs_cs_temp)),
  529. paste0(colnames(fs_cs_temp), "_ICV"),
  530. colnames(fs_cs_temp))
  531. fs_cs_temp <- fs_cs_temp %>% dplyr::rename(ICV = ICV_ICV)
  533. features_cs <- fs_cs_temp %>%
  534. dplyr::select(ICV | contains("CV") | contains("SV"))
  535. features_cs <- features_cs %>%
  536. dplyr::select(-ICV)
  538. # Perform age adjustment across all ROI volumes
  539. results_cs <- c()
  541. for (roi_num in 1:length(names(features_cs))){
  542. current_roi <- names(features_cs)[roi_num]
  543. results_cs <- append(results_cs, current_roi)
  544. }
  546. # Adjust for cubic age
  547. for (i in results_cs){
  548. fit <- lm(fs_cs_temp[[i]] ~ age + poly(age, 3, raw = TRUE)[,"3"], data = fs_cs_temp, na.action = na.exclude)
  549. fs_cs[paste0(i,'_age')] <- fs_cs[paste0(i)] - (fs_cs$age - mean(fs_cs_temp$age, na.rm = T)) * fit$coefficients[2]
  550. }
  552. # Compute global measure of hippocampal volume after ICV and age adjustments
  553. fs_cs <- fs_cs %>%
  554. dplyr::select(-names(features_cs)) %>%
  555. dplyr::rename(HVa_L = ST29SV_ICV_age,
  556. HVa_R = ST88SV_ICV_age) %>%
  557. mutate(HVa = HVa_L + HVa_R)
  559. # Compute meta-ROI cortical thickness
  560. fs_cs <- fs_cs %>%
  561. mutate(SA_total_L = ST24SA + ST32SA + ST40SA + ST26SA,
  562. SA_total_R = ST83SA + ST91SA + ST99SA + ST85SA,
  563. meta_ROI_L = (ST24TA*ST24SA + ST32TA*ST32SA + ST40TA*ST40SA + ST26TA*ST26SA) / SA_total_L,
  564. meta_ROI_R = (ST83TA*ST83SA + ST91TA*ST91SA + ST99TA*ST99SA + ST85TA*ST85SA) / SA_total_R,
  565. meta_ROI = (meta_ROI_L * SA_total_L + meta_ROI_R * SA_total_R) / (SA_total_L + SA_total_R))
  567. fs_cs <- fs_cs %>%
  568. dplyr::select(RID, VISCODE, STATUS, IMAGETYPE, MRIcs_Edate, HVa_L, HVa_R, HVa, meta_ROI_L, meta_ROI_R, meta_ROI)
  570. #______________________________________Merge Prepared ADNI Tables With Gene Expression Metadata______________________________________#####
  571. ## Read in and Preprocess Raw Gene Expression Data ####
  572. message("Now reading in and preparing raw gene expression data...")
  573. raw_exp <- read.csv(file.path(base_dir, "ADNI_Gene_Expression_Profile.csv"), stringsAsFactors = FALSE)
  575. # Isolate and minimally format gene expression-level data
  576. # Extract RIDs from PTID
  577. id <- as.data.frame(t(raw_exp[raw_exp$Phase == "SubjectID", ]))
  578. id$RID <- gsub("^.*_(\\d{4})$", "\\1", id[,1])
  579. # Extract expression-level data
  580. exp <- raw_exp
  581. colnames(exp) <- exp[8,]
  582. exp <- exp[9:nrow(exp),]
  583. # Update column names to RIDs
  584. colnames(exp)[4:ncol(exp)] <- id$RID[4:nrow(id)]
  585. # Combine gene symbol and probe set ID
  586. exp$GeneSet <- paste(exp$Symbol, exp$ProbeSet, sep = " / ")
  587. # Restrict to the expression data
  588. row.names(exp) <- NULL
  589. exp <- exp[, names(exp) != ""]
  590. exp <- exp %>%
  591. tibble::column_to_rownames("GeneSet") %>%
  592. dplyr::select(-c(ProbeSet, LocusLink, Symbol, ))
  593. # Ensure expression data is numeric
  594. exp <- exp %>%
  595. mutate(dplyr::across(everything(), ~ as.numeric(as.character(.x))))
  596. # Sort by RID
  597. exp <- exp[, order(names(exp))]
  599. ### Update Gene Symbols Using Affymetrix & MSigDB Annotation Files ####
  600. # Ensures that we maximize overlap with MSigDB's gene sets for functional enrichment analyses
  601. # There are gene symbols that were unintentionally converted to dates (e.g., 1-Sep = SEPTIN 1)
  602. message("Now updating the gene symbols in the gene expression matrix using Affymetrix and MSigDB annotation files...")
  604. # Read in and prepare Affymetrix's annotation file
  605. anno <- data.table::fread(file.path(base_dir, "AFFY-U219_GPL13667-15572.txt"), header = T, sep = "\t")
  606. anno <- anno %>%
  607. mutate(`Gene Symbol` = ifelse(`Gene Symbol` == "---" | `Gene Symbol` == "", NA, `Gene Symbol`),
  608. `Gene Symbol` = ifelse(`Gene Symbol` == "Schr2q13", "SEPTIN10", `Gene Symbol`),
  609. `Gene Symbol` = gsub("///", "||", `Gene Symbol`))
  611. # Read in and prepare a modified version of Affymetrix's annotation file tailored to MSigDB
  612. chip <- read.table(file.path(base_dir, "Affymetrix_Human_Genome_U219_MSigDB.v.7.4_custom.chip"), sep = "\t",
  613. comment.char = "", quote = "", stringsAsFactors = FALSE,
  614. fill = TRUE, header = T)
  615. chip_collapsed <- chip %>%
  616. group_by(Probe.Set.ID) %>%
  617. summarise(Gene.Symbol.New = paste(unique(Gene.Symbol), collapse = " || "), .groups = "drop")
  619. # Further refine the reformatted gene expression matrix
  620. exp$GeneSet <- row.names(exp)
  621. exp <- exp %>%
  622. mutate(Symbol = sub(" /.*", "", GeneSet),
  623. ProbeSet = sub(".*/ ", "", GeneSet)) %>%
  624. relocate(Symbol, .after = GeneSet) %>%
  625. relocate(ProbeSet, .after = Symbol)
  627. # Update gene symbols with the following priority as needed:
  628. # (1): MSigDB version, (2) Affymetrix version, and (3) ADNI version
  629. exp_final <- left_join(exp, anno %>%
  630. dplyr::select(c(ID, `Gene Symbol`)) %>%
  631. rename(ProbeSet = ID,
  632. Affy_Symbol = `Gene Symbol`), by = "ProbeSet")
  633. exp_final <- left_join(exp_final, chip_collapsed %>%
  634. dplyr::select(c(Probe.Set.ID, Gene.Symbol.New)) %>%
  635. rename(ProbeSet = Probe.Set.ID,
  636. MSigDB_Symbol = Gene.Symbol.New), by = "ProbeSet") %>%
  637. relocate(ProbeSet, .before = Symbol)
  639. # Apply the naming priority when updating the gene symbols
  640. exp_final <- exp_final %>%
  641. mutate(
  642. MSigDB_Symbol = case_when(
  643. is.na(MSigDB_Symbol) & !is.na(Affy_Symbol) ~ Affy_Symbol,
  644. is.na(MSigDB_Symbol) & is.na(Affy_Symbol) ~ Symbol,
  645. TRUE ~ MSigDB_Symbol
  646. ),
  647. # Remove all asterisks
  648. MSigDB_Symbol = gsub("\\*", "", MSigDB_Symbol),
  649. # Manually correct symbol for SEPTIN which was still in date format (e.g., 1-SEP)
  650. # Symbols for other months (Dec & Mar) were fixed using annotation files
  651. MSigDB_Symbol = ifelse(grepl("^[0-9]+-Sep$", MSigDB_Symbol, ignore.case=T),
  652. paste0("SEPTIN", sub("-Sep$", "", MSigDB_Symbol, ignore.case=T)),
  653. MSigDB_Symbol),
  654. Combined = paste(MSigDB_Symbol, ProbeSet, sep = " / ")
  655. ) %>%
  656. column_to_rownames('Combined') %>%
  657. dplyr::select(-c(GeneSet, ProbeSet, Symbol, Affy_Symbol, MSigDB_Symbol))
  659. message("Now saving the final gene expression matrix for analysis...")
  660. # Save the updated gene expression matrix for analysis
  661. write.csv(exp_final, file.path(base_dir, "gene_expression_final.csv"), row.names = T)
  663. # Clear up workspace
  664. rm(exp, exp_final, id, chip, chip_collapsed, anno)
  666. ## Create and Prepare Gene Expression Metadata Table ####
  667. # Extract gene expression metadata
  668. exp_pheno <- raw_exp[1:7,]
  669. exp_pheno <- exp_pheno %>%
  670. dplyr::select(where(~ !all(. == "")))
  671. exp_pheno <- rbind(colnames(exp_pheno), exp_pheno)
  672. exp_pheno <- as.data.frame(t(exp_pheno))
  673. colnames(exp_pheno) <- exp_pheno[1,]
  674. exp_pheno <- exp_pheno[-1, ]
  675. exp_pheno <- exp_pheno %>%
  676. dplyr::rename(VISCODE = Visit,
  677. PTID = SubjectID,
  678. GeneX_YearDrawn = YearofCollection,
  679. AffyPlate = `Affy Plate`)
  680. exp_pheno$RID <- gsub("^.*_(\\d{4})$", "\\1", exp_pheno$PTID)
  681. exp_pheno <- relocate(exp_pheno, "RID", .before = "Phase")
  682. # Fix one case of RIN written as "[5.1]" and ensure RIN is numeric
  683. exp_pheno$RIN[exp_pheno$RID == "0307"] <- "5.1"
  684. exp_pheno$RIN <- as.numeric(exp_pheno$RIN)
  685. # Sort by RID
  686. exp_pheno <- exp_pheno[order(exp_pheno$RID), ]
  687. rownames(exp_pheno) <- NULL
  688. exp_pheno$index <- seq.int(nrow(exp_pheno))
  690. ### Combine ADNI Tables with Gene Expression Metadata ####
  691. # Main demographics
  692. exp_pheno_merged <- left_join(demog_merged, exp_pheno, by = "RID")
  693. exp_pheno_merged <- exp_pheno_merged[!is.na(exp_pheno_merged$GeneX_YearDrawn),]
  695. # Update VISCODES using registry
  696. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v03", "bl", exp_pheno_merged$VISCODE)
  697. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v04", "m03", exp_pheno_merged$VISCODE)
  698. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v05", "m06", exp_pheno_merged$VISCODE)
  700. exp_pheno_merged <- left_join(exp_pheno_merged, registry, by = c("RID","VISCODE"))
  701. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v06", exp_pheno_merged$VISCODE2, exp_pheno_merged$VISCODE)
  702. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v11", exp_pheno_merged$VISCODE2, exp_pheno_merged$VISCODE)
  703. # One case of v02, so use bl instead since it is only 11 days from v02
  704. exp_pheno_merged$VISCODE <- ifelse(exp_pheno_merged$VISCODE == "v02", "bl", exp_pheno_merged$VISCODE)
  705. exp_pheno_merged <- exp_pheno_merged %>%
  706. dplyr::select(-c(VISCODE2))
  708. # Approximate EXAMDATE for gene expression profiling from ADNIMERGE to calculate age
  709. # This is necessary because the gene expression metadata only came with the year of collection
  710. adni_merge2 <- adni_merge %>% dplyr::select(RID, VISCODE, EXAMDATE)
  711. adni_merge2$RID <- as.character(adni_merge2$RID)
  712. adni_merge2$RID <- addZeros(adni_merge2$RID)
  713. exp_pheno_merged <- left_join(exp_pheno_merged, adni_merge2, by = c("RID","VISCODE"))
  714. exp_pheno_merged <- exp_pheno_merged %>%
  715. dplyr::rename(GeneX_Edate = EXAMDATE)
  716. exp_pheno_merged$GeneX_Edate <- as.Date(exp_pheno_merged$GeneX_Edate,"%Y-%m-%d")
  717. exp_pheno_merged$AGE <- round(time_length(difftime(exp_pheno_merged$GeneX_Edate,
  718. exp_pheno_merged$PTDOB),
  719. "years"), digits = 1)
  720. exp_pheno_merged <- relocate(exp_pheno_merged, "AGE", .before = "AGE_bl")
  721. exp_pheno_merged <- relocate(exp_pheno_merged, "index", .before = "RID")
  723. # APOE genotypes
  724. exp_pheno_merged <- left_join(exp_pheno_merged, apoeres_df, by = "RID")
  726. # CSF AB/TAU estimated age of positivity onset
  727. exp_pheno_merged <- left_join(exp_pheno_merged, csf_age, by = "RID")
  728. exp_pheno_merged <- exp_pheno_merged %>%
  729. mutate(AB_time = AGE - AB_AGE,
  730. TAU_time = AGE - TAU_AGE)
  732. # Diagnosis
  733. exp_pheno_merged1 <- merge(exp_pheno_merged, dx %>%
  734. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  736. exp_pheno_merged2 <- gene_within1yr(exp_pheno_merged1, DX_Edate, df2_VISCODE = TRUE)
  738. # Get back rest of the data
  739. unmatched <- anti_join(exp_pheno_merged, exp_pheno_merged2, by = c("RID", "GeneX_Edate"))
  740. exp_pheno_merged3 <- bind_rows(exp_pheno_merged2, unmatched)
  741. exp_pheno_merged3 <- exp_pheno_merged3 %>%
  742. distinct(index, .keep_all = TRUE)
  744. # Compute time2event and MCI progressor status
  745. mci <- left_join(dx, exp_pheno_merged3 %>%
  746. dplyr::select(RID, GeneX_Edate), by = "RID")
  748. mci <- mci %>%
  749. group_by(RID) %>%
  750. mutate(date_diff = abs(lubridate::time_length(difftime(GeneX_Edate, DX_Edate), "years"))) %>%
  751. group_by(RID, GeneX_Edate) %>%
  752. mutate(min_match_date = min(date_diff)) %>%
  753. ungroup()
  755. first <- mci %>%
  756. group_by(RID) %>%
  757. mutate(DX_bl = ifelse(min_match_date == date_diff, DXSUM, NA),
  758. bl_DX_date = min(DX_Edate[min_match_date == date_diff], na.rm = T),
  759. bl_DX_date = as.Date(bl_DX_date, "%Y-%m-%d")) %>%
  760. filter(!is.na(DX_bl)) %>%
  761. ungroup() %>%
  762. dplyr::select(RID, DX_bl, bl_DX_date)
  764. mci <- left_join(mci, first, by = "RID")
  766. mci <- mci %>%
  767. arrange(RID, DX_Edate) %>%
  768. group_by(RID) %>%
  769. mutate(last_DX_date = max(DX_Edate, na.rm = T),
  770. last_DX_date = as.Date(last_DX_date, "%Y-%m-%d"),
  771. window = DX_Edate >= GeneX_Edate & DX_Edate <= last_DX_date,
  772. window = ifelse(window == FALSE & bl_DX_date <= GeneX_Edate & min_match_date == date_diff, TRUE, window)) %>%
  773. filter(window == TRUE) %>%
  774. mutate(first_AD_date = min(DX_Edate[DXSUM == "Dementia" & window], na.rm = T),
  775. first_AD_date = as.Date(first_AD_date, "%Y-%m-%d"),
  776. last_MCI_date = max(DX_Edate[DXSUM == "MCI" & window], na.rm = T),
  777. last_MCI_date = as.Date(last_MCI_date, "%Y-%m-%d"),
  778. swap = DX_bl == "MCI" & window & first_AD_date <= last_MCI_date,
  779. AD_flag = any(DXSUM == "Dementia" & window),
  780. Progressor = case_when(
  781. DX_bl == "MCI" & AD_flag & swap ~ NA,
  782. DX_bl == "MCI" & AD_flag ~ 1 ,
  783. DX_bl == "MCI" & !AD_flag ~ 0,
  784. TRUE ~ NA
  785. )) %>%
  786. ungroup()
  788. mci[sapply(mci, is.infinite)] <- NA
  790. mci <- mci %>%
  791. group_by(RID) %>%
  792. mutate(time2event = ifelse(Progressor == 1,
  793. lubridate::time_length(difftime(first_AD_date, GeneX_Edate), "years"),
  794. lubridate::time_length(difftime(last_DX_date, GeneX_Edate), "years"))) %>%
  795. filter(time2event >= 0)
  797. exp_pheno_merged3 <- left_join(exp_pheno_merged3, mci %>%
  798. dplyr::select(RID, bl_DX_date, last_DX_date, first_AD_date, last_MCI_date, Progressor, time2event), by = "RID")
  799. exp_pheno_merged3 <- exp_pheno_merged3 %>%
  800. distinct(index, .keep_all = TRUE)
  802. exp_pheno_merged3 <- exp_pheno_merged3 %>%
  803. mutate(DXSUM2 = case_when(
  804. (DXSUM == "CU") ~ "CU",
  805. (DXSUM == "Dementia") ~ "Dementia",
  806. (Progressor == 0) ~ "Non-Progressor",
  807. (Progressor == 1) ~ "Progressor MCI",
  808. TRUE ~ DXSUM
  809. ))
  811. # Cognitive performance
  812. cdr <- left_join(cdr, exp_pheno_merged3 %>%
  813. dplyr::select(RID, GeneX_Edate), by = "RID")
  815. cdr <-cdr %>%
  816. arrange(RID, CDR_Edate) %>%
  817. group_by(RID) %>%
  818. mutate(CDR_time=lubridate::time_length(difftime(CDR_Edate, GeneX_Edate), "years")) %>%
  819. mutate(n = length(unique(CDR_Edate))) %>%
  820. ungroup()
  822. # Compute individualized slopes of CDRSB
  823. model <- lmer(CDRSB ~ CDR_time + (CDR_time | RID), subset = n > 1, data = cdr, na.action = na.omit)
  824. slopes <- ranef(model)$RID[2]
  825. colnames(slopes)[1] <- "CDRSB_slope"
  826. slopes$RID <- row.names(slopes)
  827. cdr <- left_join(cdr, slopes, by = "RID")
  829. cdr <- cdr %>%
  830. dplyr::select(-c(n, GeneX_Edate, CDR_time, VISDATE))
  831. cdr_slopes <- cdr %>%
  832. dplyr::select(RID, CDRSB_slope) %>%
  833. distinct()
  835. exp_pheno_merged4 <- merge(exp_pheno_merged3, cdr %>%
  836. dplyr::select(-CDRSB_slope) %>%
  837. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  839. exp_pheno_merged5 <- gene_within1yr(exp_pheno_merged4, CDR_Edate, df2_VISCODE = TRUE)
  841. # Get back rest of the data
  842. unmatched <- anti_join(exp_pheno_merged3, exp_pheno_merged5, by = c("RID", "GeneX_Edate"))
  843. exp_pheno_merged6 <- bind_rows(exp_pheno_merged5, unmatched)
  844. exp_pheno_merged6 <- exp_pheno_merged6[order(exp_pheno_merged6$RID), ]
  845. rownames(exp_pheno_merged6) <- NULL
  847. exp_pheno_merged6 <- left_join(exp_pheno_merged6, cdr_slopes, by = "RID")
  848. exp_pheno_merged6 <- left_join(exp_pheno_merged6, adas, by = c("RID", "VISCODE"))
  850. mmse <- left_join(mmse, exp_pheno_merged6 %>%
  851. dplyr::select(RID, GeneX_Edate), by = "RID")
  853. mmse <- mmse %>%
  854. arrange(RID, MMSE_Edate) %>%
  855. group_by(RID) %>%
  856. mutate(MMSE_time=lubridate::time_length(difftime(MMSE_Edate, GeneX_Edate), "years")) %>%
  857. mutate(n = length(unique(MMSE_Edate))) %>%
  858. ungroup()
  860. # Compute individualized slopes of MMSE
  861. model <- lmer(MMSE ~ MMSE_time + (MMSE_time | RID), subset = n > 1, data = mmse, na.action = na.omit)
  862. slopes <- ranef(model)$RID[2]
  863. colnames(slopes)[1] <- "MMSE_slope"
  864. slopes$RID <- row.names(slopes)
  865. mmse <- left_join(mmse, slopes, by = "RID")
  867. mmse <- mmse %>%
  868. dplyr::select(-c(n, GeneX_Edate, MMSE_time))
  870. mmse_slopes <- mmse %>%
  871. dplyr::select(RID, MMSE_slope) %>%
  872. distinct()
  874. exp_pheno_merged7 <- merge(exp_pheno_merged6, mmse %>%
  875. dplyr::select(-MMSE_slope), by = "RID", all.x = TRUE)
  877. exp_pheno_merged8 <- gene_within1yr(exp_pheno_merged7, MMSE_Edate, df2_VISCODE = TRUE)
  879. # Get back rest of the data
  880. unmatched <- anti_join(exp_pheno_merged6, exp_pheno_merged8, by = c("RID", "GeneX_Edate"))
  881. exp_pheno_merged9 <- bind_rows(exp_pheno_merged8, unmatched)
  882. exp_pheno_merged9 <- exp_pheno_merged9[order(exp_pheno_merged9$RID), ]
  883. rownames(exp_pheno_merged9) <- NULL
  885. exp_pheno_merged9 <- left_join(exp_pheno_merged9, mmse_slopes, by = "RID")
  887. # Plasma NfL
  888. nfl_data <- left_join(nfl_data, exp_pheno_merged9 %>%
  889. dplyr::select(RID, GeneX_Edate), by = "RID")
  891. nfl_data <- nfl_data %>%
  892. arrange(RID, NfL_Edate) %>%
  893. group_by(RID) %>%
  894. mutate(NfL_time=lubridate::time_length(difftime(NfL_Edate, GeneX_Edate), "years")) %>%
  895. mutate(n = length(unique(NfL_Edate))) %>%
  896. ungroup()
  898. # Compute individualized slopes of NfL
  899. model <- lmer(NfL_B ~ NfL_time + (NfL_time | RID), subset = n > 1, data = nfl_data, na.action = na.omit)
  900. slopes <- ranef(model)$RID[2]
  901. colnames(slopes)[1] <- "NfL_B_slope"
  902. slopes$RID <- row.names(slopes)
  903. nfl_data <- left_join(nfl_data, slopes, by = "RID")
  905. nfl_data <- nfl_data %>%
  906. dplyr::select(-c(n, GeneX_Edate, NfL_time))
  908. nfl_slopes <- nfl_data %>%
  909. dplyr::select(RID, NfL_B_slope) %>%
  910. distinct()
  912. exp_pheno_merged10 <- merge(exp_pheno_merged9, nfl_data %>%
  913. dplyr::select(-NfL_B_slope), by = "RID", all.x = TRUE)
  915. exp_pheno_merged11 <- gene_within1yr(exp_pheno_merged10, NfL_Edate, df2_VISCODE = TRUE)
  917. # Get back rest of the data
  918. unmatched <- anti_join(exp_pheno_merged9, exp_pheno_merged11, by = c("RID", "GeneX_Edate"))
  919. exp_pheno_merged12 <- bind_rows(exp_pheno_merged11, unmatched)
  920. exp_pheno_merged12 <- exp_pheno_merged12[order(exp_pheno_merged12$RID), ]
  921. rownames(exp_pheno_merged12) <- NULL
  923. exp_pheno_merged12 <- left_join(exp_pheno_merged12, nfl_slopes, by = "RID")
  925. # Amyloid PET
  926. ab_pet <- left_join(ab_pet, exp_pheno_merged12 %>%
  927. dplyr::select(RID, GeneX_Edate), by = "RID")
  929. ab_pet <- ab_pet %>%
  930. arrange(RID, PET_Edate) %>%
  931. group_by(RID) %>%
  932. mutate(PET_time=lubridate::time_length(difftime(PET_Edate, GeneX_Edate), "years")) %>%
  933. mutate(n = length(unique(PET_Edate))) %>%
  934. ungroup()
  936. # Compute individualized slopes of amyloid centiloids
  937. model <- lmer(CENTILOIDS ~ PET_time + (PET_time | RID), subset = n > 1, data = ab_pet, na.action = na.omit)
  938. slopes <- ranef(model)$RID[2]
  939. colnames(slopes)[1] <- "CENTILOIDS_slope"
  940. slopes$RID <- row.names(slopes)
  941. ab_pet <- left_join(ab_pet, slopes, by = "RID")
  943. ab_pet <- ab_pet %>%
  944. dplyr::select(-c(n, GeneX_Edate, PET_time))
  946. ab_slopes <- ab_pet %>%
  947. dplyr::select(RID, CENTILOIDS_slope) %>%
  948. distinct()
  950. exp_pheno_merged13 <- merge(exp_pheno_merged12, ab_pet %>%
  951. dplyr::select(-CENTILOIDS_slope), by = "RID", all.x = TRUE)
  953. exp_pheno_merged14 <- gene_within1yr(exp_pheno_merged13, PET_Edate, df2_VISCODE = FALSE)
  955. # Get back rest of the data
  956. unmatched <- anti_join(exp_pheno_merged12, exp_pheno_merged14, by = c("RID", "GeneX_Edate"))
  957. exp_pheno_merged15 <- bind_rows(exp_pheno_merged14, unmatched)
  958. exp_pheno_merged15 <- exp_pheno_merged15[order(exp_pheno_merged15$RID), ]
  959. rownames(exp_pheno_merged15) <- NULL
  961. exp_pheno_merged15 <- left_join(exp_pheno_merged15, ab_slopes, by = "RID")
  963. exp_pheno_merged15 <- merge(exp_pheno_merged15, last_negative_date_pet, by = "RID", all.x = TRUE)
  964. exp_pheno_merged15 <- merge(exp_pheno_merged15, first_positive_date_pet, by = "RID", all.x = TRUE)
  965. exp_pheno_merged15 <- merge(exp_pheno_merged15, last_negative_date2_pet, by = "RID", all.x = TRUE)
  966. exp_pheno_merged15 <- merge(exp_pheno_merged15, first_positive_date2_pet, by = "RID", all.x = TRUE)
  967. exp_pheno_merged15 <- merge(exp_pheno_merged15, last_negative_date3_pet, by = "RID", all.x = TRUE)
  968. exp_pheno_merged15 <- merge(exp_pheno_merged15, first_positive_date3_pet, by = "RID", all.x = TRUE)
  970. exp_pheno_merged15 <- exp_pheno_merged15[order(exp_pheno_merged15$RID), ]
  971. rownames(exp_pheno_merged15) <- NULL
  973. # Ensure amyloid PET positivity follows expected logic
  974. # Prioritize using available data, then apply expected logic where appropriate
  975. # (e.g., can progress from negative to positive; once positive stays positive)
  976. exp_pheno_merged16 <- exp_pheno_merged15 %>%
  977. group_by(RID) %>%
  978. mutate(amyloid_pos_pet_cross = case_when(
  979. !is.na(last_AB_neg_date) & !is.na(first_AB_pos_date) & is.na(amyloid_pos_pet_cross) & first_AB_pos_date < last_AB_neg_date ~ amyloid_pos_pet_cross,
  980. !is.na(last_AB_neg_date) & is.na(amyloid_pos_pet_cross) & GeneX_Edate <= last_AB_neg_date ~ 0,
  981. !is.na(first_AB_pos_date) & is.na(amyloid_pos_pet_cross) & GeneX_Edate >= first_AB_pos_date ~ 1,
  982. TRUE ~ amyloid_pos_pet_cross)) %>%
  983. ungroup()
  985. exp_pheno_merged16 <- exp_pheno_merged16 %>%
  986. group_by(RID) %>%
  987. mutate(amyloid_pos_pet_long = case_when(
  988. !is.na(last_AB_neg_date2) & !is.na(first_AB_pos_date2) & is.na(amyloid_pos_pet_long) & first_AB_pos_date2 < last_AB_neg_date2 ~ amyloid_pos_pet_long,
  989. !is.na(last_AB_neg_date2) & is.na(amyloid_pos_pet_long) & GeneX_Edate <= last_AB_neg_date2 ~ 0,
  990. !is.na(first_AB_pos_date2) & is.na(amyloid_pos_pet_long) & GeneX_Edate >= first_AB_pos_date2 ~ 1,
  991. TRUE ~ amyloid_pos_pet_long)) %>%
  992. ungroup()
  994. exp_pheno_merged16 <- exp_pheno_merged16 %>%
  995. group_by(RID) %>%
  996. mutate(amyloid_pos_pet_centiloid = case_when(
  997. !is.na(last_AB_neg_date3) & !is.na(first_AB_pos_date3) & is.na(amyloid_pos_pet_centiloid) & first_AB_pos_date3 < last_AB_neg_date3 ~ amyloid_pos_pet_centiloid,
  998. !is.na(last_AB_neg_date3) & is.na(amyloid_pos_pet_centiloid) & GeneX_Edate <= last_AB_neg_date3 ~ 0,
  999. !is.na(first_AB_pos_date3) & is.na(amyloid_pos_pet_centiloid) & GeneX_Edate >= first_AB_pos_date3 ~ 1,
  1000. TRUE ~ amyloid_pos_pet_centiloid)) %>%
  1001. ungroup() %>%
  1002. dplyr::select(-c("last_AB_neg_date", "first_AB_pos_date", "last_AB_neg_date2", "first_AB_pos_date2", "last_AB_neg_date3", "first_AB_pos_date3"))
  1004. # CSF biomarkers
  1005. upenn_csf <- left_join(upenn_csf, exp_pheno_merged16 %>%
  1006. dplyr::select(RID, GeneX_Edate), by = "RID")
  1008. upenn_csf <- upenn_csf %>%
  1009. arrange(RID, CSF_Edate) %>%
  1010. group_by(RID) %>%
  1011. mutate(CSF_time=lubridate::time_length(difftime(CSF_Edate, GeneX_Edate), "years")) %>%
  1012. mutate(n = length(unique(CSF_Edate))) %>%
  1013. ungroup()
  1015. # Compute individualized slopes of CSF p-tau181
  1016. model <- lmer(pTAU181_csf ~ CSF_time + (CSF_time | RID), subset = n > 1, data = upenn_csf, na.action = na.omit)
  1017. slopes <- ranef(model)$RID[2]
  1018. colnames(slopes)[1] <- "pTAU181_csf_slope"
  1019. slopes$RID <- row.names(slopes)
  1020. upenn_csf <- left_join(upenn_csf, slopes, by = "RID")
  1022. upenn_csf <- upenn_csf %>%
  1023. dplyr::select(-c(n, GeneX_Edate, CSF_time))
  1025. csf_slopes <- upenn_csf %>%
  1026. dplyr::select(RID, pTAU181_csf_slope) %>%
  1027. distinct()
  1029. exp_pheno_merged17 <- merge(exp_pheno_merged16, upenn_csf %>%
  1030. dplyr::select(-pTAU181_csf_slope) %>%
  1031. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  1033. exp_pheno_merged18 <- gene_within1yr(exp_pheno_merged17, CSF_Edate, df2_VISCODE = TRUE)
  1035. # Get back rest of the data
  1036. unmatched <- anti_join(exp_pheno_merged16, exp_pheno_merged18, by = c("RID", "GeneX_Edate"))
  1037. exp_pheno_merged19 <- bind_rows(exp_pheno_merged18, unmatched)
  1038. exp_pheno_merged19 <- exp_pheno_merged19[order(exp_pheno_merged19$RID), ]
  1039. rownames(exp_pheno_merged19) <- NULL
  1041. exp_pheno_merged19 <- left_join(exp_pheno_merged19, csf_slopes, by = "RID")
  1043. exp_pheno_merged19 <- merge(exp_pheno_merged19, last_negative_date_csf, by = "RID", all.x = TRUE)
  1044. exp_pheno_merged19 <- merge(exp_pheno_merged19, first_positive_date_csf, by = "RID", all.x = TRUE)
  1045. exp_pheno_merged19 <- merge(exp_pheno_merged19, last_negative_date2_csf, by = "RID", all.x = TRUE)
  1046. exp_pheno_merged19 <- merge(exp_pheno_merged19, first_positive_date2_csf, by = "RID", all.x = TRUE)
  1047. exp_pheno_merged19 <- merge(exp_pheno_merged19, last_negative_date3_csf, by = "RID", all.x = TRUE)
  1048. exp_pheno_merged19 <- merge(exp_pheno_merged19, first_positive_date3_csf, by = "RID", all.x = TRUE)
  1050. exp_pheno_merged19 <- exp_pheno_merged19[order(exp_pheno_merged19$RID), ]
  1051. rownames(exp_pheno_merged19) <- NULL
  1053. # Ensure CSF biomarker positivity follows expected logic
  1054. exp_pheno_merged20 <- exp_pheno_merged19 %>%
  1055. group_by(RID) %>%
  1056. mutate(amyloid_pos_csf = case_when(
  1057. !is.na(last_AB_neg_date) & !is.na(first_AB_pos_date) & is.na(amyloid_pos_csf) & first_AB_pos_date < last_AB_neg_date ~ amyloid_pos_csf,
  1058. !is.na(last_AB_neg_date) & is.na(amyloid_pos_csf) & GeneX_Edate <= last_AB_neg_date ~ 0,
  1059. !is.na(first_AB_pos_date) & is.na(amyloid_pos_csf) & GeneX_Edate >= first_AB_pos_date ~ 1,
  1060. TRUE ~ amyloid_pos_csf)) %>%
  1061. ungroup()
  1063. exp_pheno_merged20 <- exp_pheno_merged20 %>%
  1064. group_by(RID) %>%
  1065. mutate(ptau_pos_csf = case_when(
  1066. !is.na(last_TAU_neg_date2) & !is.na(first_TAU_pos_date2) & is.na(ptau_pos_csf) & first_TAU_pos_date2 < last_TAU_neg_date2 ~ ptau_pos_csf,
  1067. !is.na(last_TAU_neg_date2) & is.na(ptau_pos_csf) & GeneX_Edate <= last_TAU_neg_date2 ~ 0,
  1068. !is.na(first_TAU_pos_date2) & is.na(ptau_pos_csf) & GeneX_Edate >= first_TAU_pos_date2 ~ 1,
  1069. TRUE ~ ptau_pos_csf)) %>%
  1070. ungroup()
  1072. exp_pheno_merged20 <- exp_pheno_merged20 %>%
  1073. group_by(RID) %>%
  1074. mutate(ad_pathology_pos_csf = case_when(
  1075. !is.na(last_AD_neg_date3) & !is.na(first_AD_pos_date3) & is.na(ad_pathology_pos_csf) & first_AD_pos_date3 < last_AD_neg_date3 ~ ad_pathology_pos_csf,
  1076. !is.na(last_AD_neg_date3) & is.na(ad_pathology_pos_csf) & GeneX_Edate <= last_AD_neg_date3 ~ 0,
  1077. !is.na(first_AD_pos_date3) & is.na(ad_pathology_pos_csf) & GeneX_Edate >= first_AD_pos_date3 ~ 1,
  1078. TRUE ~ ad_pathology_pos_csf)) %>%
  1079. ungroup() %>%
  1080. dplyr::select(-c("last_AB_neg_date", "first_AB_pos_date", "last_TAU_neg_date2", "first_TAU_pos_date2", "last_AD_neg_date3", "first_AD_pos_date3", "index"))
  1082. # Hippocampal volume and meta-ROI cortical thickness
  1083. fs_cs <- left_join(fs_cs, exp_pheno_merged20 %>%
  1084. dplyr::select(RID, GeneX_Edate), by = "RID")
  1086. fs_cs <- fs_cs %>%
  1087. arrange(RID, MRIcs_Edate) %>%
  1088. group_by(RID) %>%
  1089. mutate(MRI_time=lubridate::time_length(difftime(MRIcs_Edate, GeneX_Edate), "years")) %>%
  1090. mutate(n = length(unique(MRIcs_Edate))) %>%
  1091. ungroup()
  1093. # Compute individualized slopes of hippocampal volume
  1094. model <- lmer(HVa ~ MRI_time + (MRI_time | RID), subset = n > 1, data = fs_cs, na.action = na.omit)
  1095. slopes <- ranef(model)$RID[2]
  1096. colnames(slopes)[1] <- "HVa_slope"
  1097. slopes$RID <- row.names(slopes)
  1098. fs_cs <- left_join(fs_cs, slopes, by = "RID")
  1100. # Compute individualized slopes of meta-ROI cortical thickness
  1101. model <- lmer(meta_ROI ~ MRI_time + (MRI_time | RID), subset = n > 1, data = fs_cs, na.action = na.omit)
  1102. slopes <- ranef(model)$RID[2]
  1103. colnames(slopes)[1] <- "meta_ROI_slope"
  1104. slopes$RID <- row.names(slopes)
  1105. fs_cs <- left_join(fs_cs, slopes, by = "RID")
  1107. fs_slopes <- fs_cs %>%
  1108. dplyr::select(RID, HVa_slope, meta_ROI_slope) %>%
  1109. distinct()
  1111. fs_cs <- fs_cs %>%
  1112. dplyr::select(-c(n, GeneX_Edate, MRI_time, HVa_slope, meta_ROI_slope))
  1114. exp_pheno_merged21 <- merge(exp_pheno_merged20, fs_cs %>%
  1115. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  1117. exp_pheno_merged22 <- gene_within1yr(exp_pheno_merged21, MRIcs_Edate, df2_VISCODE = TRUE)
  1119. # Get back rest of the data
  1120. unmatched <- anti_join(exp_pheno_merged20, exp_pheno_merged22, by = c("RID", "GeneX_Edate"))
  1121. exp_pheno_merged23 <- bind_rows(exp_pheno_merged22, unmatched)
  1122. exp_pheno_merged23 <- exp_pheno_merged23[order(exp_pheno_merged23$RID), ]
  1123. rownames(exp_pheno_merged23) <- NULL
  1125. exp_pheno_merged23 <- left_join(exp_pheno_merged23, fs_slopes, by = "RID")
  1127. message("Now saving final gene expression phenotype data for analysis...")
  1128. # Save final gene expression phenotype data for analysis
  1129. write.csv(exp_pheno_merged23, file.path(base_dir, "GeneX_pheno_final.csv"), row.names = F)
  1131. #______________________________________Merge Prepared ADNI Tables With DNA Methylation Metadata______________________________________#####
  1132. ## Create and Prepare DNA Methylation Metadata Table ####
  1133. message("Now reading in and preparing DNA methylation metadata table...")
  1134. dnam_pheno <- read.csv(file.path(base_dir, "Sample_Sheet_rv1.csv"), stringsAsFactor = FALSE)
  1135. # Exclude the 14 samples not in metadata (failed QC by ADNI Genetics Core)
  1136. dnam_pheno <- dnam_pheno[!is.na(dnam_pheno$RID), ]
  1137. dnam_pheno <- dnam_pheno %>%
  1138. dplyr::rename(DNAm_Edate = Edate,
  1139. DNAm_DateDrawn = DateDrawn) %>%
  1140. dplyr::select(-Basename)
  1141. dnam_pheno$RID <- addZeros(dnam_pheno$RID)
  1142. dnam_pheno$DNAm_Edate <- as.Date(dnam_pheno$DNAm_Edate, "%m/%d/%Y")
  1143. dnam_pheno$DNAm_DateDrawn <- as.Date(dnam_pheno$DNAm_DateDrawn, "%m/%d/%Y")
  1145. ### Combine ADNI Tables with DNA Methylation Metadata ####
  1146. # Main demographics
  1147. dnam_pheno_merged <- left_join(demog_merged, dnam_pheno, by = "RID")
  1148. dnam_pheno_merged <- dnam_pheno_merged[!is.na(dnam_pheno_merged$barcodes),]
  1149. dnam_pheno_merged$AGE <- round(time_length(difftime(dnam_pheno_merged$DNAm_DateDrawn,
  1150. dnam_pheno_merged$PTDOB),
  1151. "years"), digits = 1)
  1152. dnam_pheno_merged <- relocate(dnam_pheno_merged, "AGE", .before = "AGE_bl")
  1153. dnam_pheno_merged <- relocate(dnam_pheno_merged, "barcodes", .after = "RID")
  1155. # APOE genotypes
  1156. dnam_pheno_merged <- left_join(dnam_pheno_merged, apoeres_df, by = "RID")
  1158. # CSF AB/TAU estimated age of positivity onset
  1159. dnam_pheno_merged <- left_join(dnam_pheno_merged, csf_age, by = "RID")
  1160. dnam_pheno_merged <- dnam_pheno_merged %>%
  1161. mutate(AB_time = AGE - AB_AGE,
  1162. TAU_time = AGE - TAU_AGE)
  1164. # Diagnosis
  1165. dnam_pheno_merged1 <- merge(dnam_pheno_merged, dx, by = "RID", all.x = TRUE)
  1167. dnam_pheno_merged2 <- dnam_within1yr(dnam_pheno_merged1, DX_Edate, df2_VISCODE = FALSE)
  1169. # Get back rest of the data
  1170. unmatched <- anti_join(dnam_pheno_merged, dnam_pheno_merged2, by = c("RID", "DNAm_Edate"))
  1171. dnam_pheno_merged3 <- bind_rows(dnam_pheno_merged2, unmatched)
  1172. dnam_pheno_merged3 <- dnam_pheno_merged3 %>%
  1173. distinct(barcodes, .keep_all = TRUE)
  1175. # Compute time2event and MCI progressor status
  1176. mci <- left_join(dx, dnam_pheno_merged3 %>%
  1177. dplyr::select(RID, DNAm_Edate), by = "RID")
  1179. mci <- mci %>%
  1180. group_by(RID, DNAm_Edate) %>%
  1181. mutate(date_diff = abs(lubridate::time_length(difftime(DNAm_Edate, DX_Edate), "years")),
  1182. min_match_date = min(date_diff)) %>%
  1183. ungroup()
  1185. first <- mci %>%
  1186. arrange(RID, DNAm_Edate) %>%
  1187. group_by(RID) %>%
  1188. mutate(DX_bl = DXSUM[min_match_date == date_diff][1],
  1189. bl_DX_date = DX_Edate[min_match_date == date_diff][1],
  1190. bl_DX_date = as.Date(bl_DX_date, "%Y-%m-%d")) %>%
  1191. filter(!is.na(DX_bl)) %>%
  1192. ungroup() %>%
  1193. distinct(RID, .keep_all = T) %>%
  1194. dplyr::select(RID, DX_bl, bl_DX_date)
  1196. mci <- left_join(mci, first, by = "RID")
  1198. mci <- mci %>%
  1199. arrange(RID, DNAm_Edate, DX_Edate) %>%
  1200. group_by(RID) %>%
  1201. mutate(last_DX_date = max(DX_Edate, na.rm = T),
  1202. last_DX_date = as.Date(last_DX_date, "%Y-%m-%d"),
  1203. window = DX_Edate >= DNAm_Edate[[1]] & DX_Edate <= last_DX_date) %>%
  1204. filter(window == TRUE) %>%
  1205. mutate(first_AD_date = min(DX_Edate[DXSUM == "Dementia" & window], na.rm = T),
  1206. first_AD_date = as.Date(first_AD_date, "%Y-%m-%d"),
  1207. last_MCI_date = max(DX_Edate[DXSUM == "MCI" & window], na.rm = T),
  1208. last_MCI_date = as.Date(last_MCI_date, "%Y-%m-%d"),
  1209. swap = DX_bl == "MCI" & window & first_AD_date <= last_MCI_date,
  1210. AD_flag = any(DXSUM == "Dementia" & window),
  1211. Progressor = case_when(
  1212. DX_bl == "MCI" & AD_flag & swap ~ NA,
  1213. DX_bl == "MCI" & AD_flag ~ 1 ,
  1214. DX_bl == "MCI" & !AD_flag ~ 0,
  1215. TRUE ~ NA
  1216. )) %>%
  1217. ungroup()
  1219. mci[sapply(mci, is.infinite)] <- NA
  1221. mci <- mci %>%
  1222. group_by(RID) %>%
  1223. mutate(time2event = ifelse(Progressor == 1,
  1224. lubridate::time_length(difftime(first_AD_date, DNAm_Edate), "years"),
  1225. lubridate::time_length(difftime(last_DX_date, DNAm_Edate), "years"))) %>%
  1226. filter(time2event >= 0)
  1228. dnam_pheno_merged3 <- left_join(dnam_pheno_merged3, mci %>%
  1229. dplyr::select(RID, DX_bl, bl_DX_date, last_DX_date, first_AD_date, last_MCI_date, Progressor, time2event), by = "RID")
  1230. dnam_pheno_merged3 <- dnam_pheno_merged3 %>%
  1231. distinct(barcodes, .keep_all = TRUE)
  1233. dnam_pheno_merged3 <- dnam_pheno_merged3 %>%
  1234. mutate(DXSUM2 = case_when(
  1235. (Progressor == 0) ~ "Non-Progressor",
  1236. (Progressor == 1) ~ "Progressor MCI",
  1237. (DXSUM == "CU") ~ "CU",
  1238. (DXSUM == "Dementia") ~ "Dementia",
  1239. TRUE ~ DXSUM
  1240. ))
  1242. # Cognitive performance & individualized slopes of CDRSB and MMSE
  1243. dnam_pheno_merged4 <- merge(dnam_pheno_merged3, cdr %>%
  1244. dplyr::select(-CDRSB_slope) %>%
  1245. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  1247. dnam_pheno_merged5 <- dnam_within1yr(dnam_pheno_merged4, CDR_Edate, df2_VISCODE = TRUE)
  1249. # Get back rest of the data
  1250. unmatched <- anti_join(dnam_pheno_merged3, dnam_pheno_merged5, by = c("RID", "DNAm_Edate"))
  1251. dnam_pheno_merged6 <- bind_rows(dnam_pheno_merged5, unmatched)
  1252. dnam_pheno_merged6 <- dnam_pheno_merged6[order(dnam_pheno_merged6$RID, dnam_pheno_merged6$DNAm_Edate), ]
  1253. rownames(dnam_pheno_merged6) <- NULL
  1255. dnam_pheno_merged6 <- left_join(dnam_pheno_merged6, cdr_slopes, by = "RID")
  1257. dnam_pheno_merged6 <- left_join(dnam_pheno_merged6, adas, by = c("RID", "VISCODE"))
  1259. dnam_pheno_merged7 <- merge(dnam_pheno_merged6, mmse %>%
  1260. dplyr::select(-MMSE_slope), by = "RID", all.x = TRUE)
  1262. dnam_pheno_merged8 <- dnam_within1yr(dnam_pheno_merged7, MMSE_Edate, df2_VISCODE = TRUE)
  1264. # Get back rest of the data
  1265. unmatched <- anti_join(dnam_pheno_merged6, dnam_pheno_merged8, by = c("RID", "DNAm_Edate"))
  1266. dnam_pheno_merged9 <- bind_rows(dnam_pheno_merged8, unmatched)
  1267. dnam_pheno_merged9 <- dnam_pheno_merged9[order(dnam_pheno_merged9$RID, dnam_pheno_merged9$DNAm_Edate), ]
  1268. rownames(dnam_pheno_merged9) <- NULL
  1270. dnam_pheno_merged9 <- left_join(dnam_pheno_merged9, mmse_slopes, by = "RID")
  1272. # Plasma NfL & individualized slopes
  1273. dnam_pheno_merged10 <- merge(dnam_pheno_merged9, nfl_data %>%
  1274. dplyr::select(-NfL_B_slope), by = "RID", all.x = TRUE)
  1276. dnam_pheno_merged11 <- dnam_within1yr(dnam_pheno_merged10, NfL_Edate, df2_VISCODE = TRUE)
  1278. # Get back rest of the data
  1279. unmatched <- anti_join(dnam_pheno_merged9, dnam_pheno_merged11, by = c("RID", "DNAm_Edate"))
  1280. dnam_pheno_merged12 <- bind_rows(dnam_pheno_merged11, unmatched)
  1281. dnam_pheno_merged12 <- dnam_pheno_merged12[order(dnam_pheno_merged12$RID, dnam_pheno_merged12$DNAm_Edate), ]
  1282. rownames(dnam_pheno_merged12) <- NULL
  1284. dnam_pheno_merged12 <- left_join(dnam_pheno_merged12, nfl_slopes, by = "RID")
  1286. # Amyloid PET & individualized slopes of centiloids
  1287. dnam_pheno_merged13 <- merge(dnam_pheno_merged12, ab_pet %>%
  1288. dplyr::select(-CENTILOIDS_slope), by = "RID", all.x = TRUE)
  1290. dnam_pheno_merged14 <- dnam_within1yr(dnam_pheno_merged13, PET_Edate, df2_VISCODE = FALSE)
  1292. # Get back rest of the data
  1293. unmatched <- anti_join(dnam_pheno_merged12, dnam_pheno_merged14, by = c("RID", "DNAm_Edate"))
  1294. dnam_pheno_merged15 <- bind_rows(dnam_pheno_merged14, unmatched)
  1295. dnam_pheno_merged15 <- dnam_pheno_merged15[order(dnam_pheno_merged15$RID, dnam_pheno_merged15$DNAm_Edate), ]
  1296. rownames(dnam_pheno_merged15) <- NULL
  1298. dnam_pheno_merged15 <- left_join(dnam_pheno_merged15, ab_slopes, by = "RID")
  1300. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, last_negative_date_pet, by = "RID", all.x = TRUE)
  1301. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, first_positive_date_pet, by = "RID", all.x = TRUE)
  1302. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, last_negative_date2_pet, by = "RID", all.x = TRUE)
  1303. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, first_positive_date2_pet, by = "RID", all.x = TRUE)
  1304. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, last_negative_date3_pet, by = "RID", all.x = TRUE)
  1305. dnam_pheno_merged15 <- merge(dnam_pheno_merged15, first_positive_date3_pet, by = "RID", all.x = TRUE)
  1307. dnam_pheno_merged15 <- dnam_pheno_merged15[order(dnam_pheno_merged15$RID, dnam_pheno_merged15$DNAm_Edate), ]
  1308. rownames(dnam_pheno_merged15) <- NULL
  1310. # Ensure amyloid PET positivity follows expected logic
  1311. dnam_pheno_merged16 <- dnam_pheno_merged15 %>%
  1312. group_by(RID) %>%
  1313. mutate(amyloid_pos_pet_cross = case_when(
  1314. !is.na(last_AB_neg_date) & !is.na(first_AB_pos_date) & is.na(amyloid_pos_pet_cross) & first_AB_pos_date < last_AB_neg_date ~ amyloid_pos_pet_cross,
  1315. !is.na(last_AB_neg_date) & is.na(amyloid_pos_pet_cross) & DNAm_Edate <= last_AB_neg_date ~ 0,
  1316. !is.na(first_AB_pos_date) & is.na(amyloid_pos_pet_cross) & DNAm_Edate >= first_AB_pos_date ~ 1,
  1317. TRUE ~ amyloid_pos_pet_cross)) %>%
  1318. ungroup()
  1320. dnam_pheno_merged16 <- dnam_pheno_merged16 %>%
  1321. group_by(RID) %>%
  1322. mutate(amyloid_pos_pet_long = case_when(
  1323. !is.na(last_AB_neg_date2) & !is.na(first_AB_pos_date2) & is.na(amyloid_pos_pet_long) & first_AB_pos_date2 < last_AB_neg_date2 ~ amyloid_pos_pet_long,
  1324. !is.na(last_AB_neg_date2) & is.na(amyloid_pos_pet_long) & DNAm_Edate <= last_AB_neg_date2 ~ 0,
  1325. !is.na(first_AB_pos_date2) & is.na(amyloid_pos_pet_long) & DNAm_Edate >= first_AB_pos_date2 ~ 1,
  1326. TRUE ~ amyloid_pos_pet_long)) %>%
  1327. ungroup()
  1329. dnam_pheno_merged16 <- dnam_pheno_merged16 %>%
  1330. group_by(RID) %>%
  1331. mutate(amyloid_pos_pet_centiloid = case_when(
  1332. !is.na(last_AB_neg_date3) & !is.na(first_AB_pos_date3) & is.na(amyloid_pos_pet_centiloid) & first_AB_pos_date3 < last_AB_neg_date3 ~ amyloid_pos_pet_centiloid,
  1333. !is.na(last_AB_neg_date3) & is.na(amyloid_pos_pet_centiloid) & DNAm_Edate <= last_AB_neg_date3 ~ 0,
  1334. !is.na(first_AB_pos_date3) & is.na(amyloid_pos_pet_centiloid) & DNAm_Edate >= first_AB_pos_date3 ~ 1,
  1335. TRUE ~ amyloid_pos_pet_centiloid)) %>%
  1336. ungroup() %>%
  1337. dplyr::select(-c("last_AB_neg_date", "first_AB_pos_date", "last_AB_neg_date2", "first_AB_pos_date2", "last_AB_neg_date3", "first_AB_pos_date3"))
  1339. # CSF biomarkers & individualized slopes of p-tau181
  1340. dnam_pheno_merged17 <- merge(dnam_pheno_merged16, upenn_csf %>%
  1341. dplyr::select(-pTAU181_csf_slope) %>%
  1342. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  1344. dnam_pheno_merged18 <- dnam_within1yr(dnam_pheno_merged17, CSF_Edate, df2_VISCODE = TRUE)
  1346. # Get back rest of the data
  1347. unmatched <- anti_join(dnam_pheno_merged16, dnam_pheno_merged18, by = c("RID", "DNAm_Edate"))
  1348. dnam_pheno_merged19 <- bind_rows(dnam_pheno_merged18, unmatched)
  1349. dnam_pheno_merged19 <- dnam_pheno_merged19[order(dnam_pheno_merged19$RID, dnam_pheno_merged19$DNAm_Edate), ]
  1350. rownames(dnam_pheno_merged19) <- NULL
  1352. dnam_pheno_merged19 <- left_join(dnam_pheno_merged19, csf_slopes, by = "RID")
  1354. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, last_negative_date_csf, by = "RID", all.x = TRUE)
  1355. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, first_positive_date_csf, by = "RID", all.x = TRUE)
  1356. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, last_negative_date2_csf, by = "RID", all.x = TRUE)
  1357. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, first_positive_date2_csf, by = "RID", all.x = TRUE)
  1358. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, last_negative_date3_csf, by = "RID", all.x = TRUE)
  1359. dnam_pheno_merged19 <- merge(dnam_pheno_merged19, first_positive_date3_csf, by = "RID", all.x = TRUE)
  1361. dnam_pheno_merged19 <- dnam_pheno_merged19[order(dnam_pheno_merged19$RID, dnam_pheno_merged19$DNAm_Edate), ]
  1362. rownames(dnam_pheno_merged19) <- NULL
  1364. # Ensure CSF biomarker positivity follows expected logic
  1365. dnam_pheno_merged20 <- dnam_pheno_merged19 %>%
  1366. group_by(RID) %>%
  1367. mutate(amyloid_pos_csf = case_when(
  1368. !is.na(last_AB_neg_date) & !is.na(first_AB_pos_date) & is.na(amyloid_pos_csf) & first_AB_pos_date < last_AB_neg_date ~ amyloid_pos_csf,
  1369. !is.na(last_AB_neg_date) & is.na(amyloid_pos_csf) & DNAm_Edate <= last_AB_neg_date ~ 0,
  1370. !is.na(first_AB_pos_date) & is.na(amyloid_pos_csf) & DNAm_Edate >= first_AB_pos_date ~ 1,
  1371. TRUE ~ amyloid_pos_csf)) %>%
  1372. ungroup()
  1374. dnam_pheno_merged20 <- dnam_pheno_merged20 %>%
  1375. group_by(RID) %>%
  1376. mutate(ptau_pos_csf = case_when(
  1377. !is.na(last_TAU_neg_date2) & !is.na(first_TAU_pos_date2) & is.na(ptau_pos_csf) & first_TAU_pos_date2 < last_TAU_neg_date2 ~ ptau_pos_csf,
  1378. !is.na(last_TAU_neg_date2) & is.na(ptau_pos_csf) & DNAm_Edate <= last_TAU_neg_date2 ~ 0,
  1379. !is.na(first_TAU_pos_date2) & is.na(ptau_pos_csf) & DNAm_Edate >= first_TAU_pos_date2 ~ 1,
  1380. TRUE ~ ptau_pos_csf)) %>%
  1381. ungroup()
  1383. dnam_pheno_merged20 <- dnam_pheno_merged20 %>%
  1384. group_by(RID) %>%
  1385. mutate(ad_pathology_pos_csf = case_when(
  1386. !is.na(last_AD_neg_date3) & !is.na(first_AD_pos_date3) & is.na(ad_pathology_pos_csf) & first_AD_pos_date3 < last_AD_neg_date3 ~ ad_pathology_pos_csf,
  1387. !is.na(last_AD_neg_date3) & is.na(ad_pathology_pos_csf) & DNAm_Edate <= last_AD_neg_date3 ~ 0,
  1388. !is.na(first_AD_pos_date3) & is.na(ad_pathology_pos_csf) & DNAm_Edate >= first_AD_pos_date3 ~ 1,
  1389. TRUE ~ ad_pathology_pos_csf)) %>%
  1390. ungroup() %>%
  1391. dplyr::select(-c("last_AB_neg_date", "first_AB_pos_date", "last_TAU_neg_date2", "first_TAU_pos_date2", "last_AD_neg_date3", "first_AD_pos_date3"))
  1393. # Hippocampal volume and meta-ROI cortical thickness & individualized slopes
  1394. dnam_pheno_merged21 <- merge(dnam_pheno_merged20, fs_cs %>%
  1395. dplyr::rename(VISCODE2 = VISCODE), by = "RID", all.x = TRUE)
  1397. dnam_pheno_merged22 <- dnam_within1yr(dnam_pheno_merged21, MRIcs_Edate, df2_VISCODE = TRUE)
  1399. # Get back rest of the data
  1400. unmatched <- anti_join(dnam_pheno_merged20, dnam_pheno_merged22, by = c("RID", "DNAm_Edate"))
  1401. dnam_pheno_merged23 <- bind_rows(dnam_pheno_merged22, unmatched)
  1402. dnam_pheno_merged23 <- dnam_pheno_merged23[order(dnam_pheno_merged23$RID, dnam_pheno_merged23$DNAm_Edate), ]
  1403. rownames(dnam_pheno_merged23) <- NULL
  1405. dnam_pheno_merged23 <- left_join(dnam_pheno_merged23, fs_slopes, by = "RID")
  1407. # Estimate immune cell-type proportions
  1408. adni_se <- readRDS(file.path(proc_dir, "dasen_se_all_rv1.RDS"))
  1409. beta <- as.matrix(assays(adni_se)$DNAm)
  1410. rm(adni_se)
  1412. data(centDHSbloodDMC.m)
  1413. BloodFrac.m <- epidish(beta.m = beta, ref.m = centDHSbloodDMC.m, method = "RPC")$estF
  1414. cell_types <- as.data.frame(BloodFrac.m)
  1415. cell_types <- tibble::rownames_to_column(cell_types, "barcodes")
  1416. cell_types$Gran <- cell_types$Neutro + cell_types$Eosino
  1418. dnam_pheno_merged24 <- left_join(dnam_pheno_merged23, cell_types, by = "barcodes")
  1420. message("Now saving final (longitudinal) DNA methylation phenotype data...")
  1421. # Save final DNA methylation phenotype data for analysis
  1422. write.csv(dnam_pheno_merged24, file.path(base_dir, "DNAm_pheno_final.csv"), row.names = F)
  1424. ## Create Subset of DNA Methylation Phenotype Data Aligned to Time of Gene Expression Profiling ####
  1425. # Matched at same time-point or within 1-year of gene expression profiling
  1426. message("Now creating a cross-sectional subset of DNA methylation phenotype data based on time of gene expression profiling...")
  1427. # Read in phenotype data from both -omics that we just created
  1428. exp_pheno <- read.csv(file.path(base_dir, "GeneX_pheno_final.csv"), stringsAsFactor = FALSE, check.names = FALSE)
  1429. exp_pheno$RID <- addZeros(exp_pheno$RID)
  1430. dnam_pheno <- read.csv(file.path(base_dir, "DNAm_pheno_final.csv"), stringsAsFactor = FALSE)
  1431. dnam_pheno$RID <- addZeros(dnam_pheno$RID)
  1433. # Identify the cross-sectional matched time-point for DNA methylation based on time of gene expression profiling
  1434. merged <- merge(dnam_pheno %>%
  1435. dplyr::select(-c(DXSUM2, Progressor, time2event)) %>%
  1436. dplyr::rename(DNAm_DX_bl = DX_bl),
  1437. exp_pheno %>%
  1438. dplyr::select(RID, VISCODE, GeneX_Edate, DX_Edate, DXSUM, DXSUM2, Progressor, time2event) %>%
  1439. dplyr::rename(GeneX_VISCODE = VISCODE,
  1440. GeneX_DXSUM = DXSUM,
  1441. GeneX_DX_Edate = DX_Edate), by = "RID", all.x = TRUE) %>%
  1442. # Use time of gene expression profiling as the anchor
  1443. mutate(date_diff = abs(lubridate::time_length(difftime(DNAm_Edate, GeneX_Edate), "years"))) %>%
  1444. mutate(Edate_match = ifelse(date_diff <= 1, 1, 0),
  1445. VISCODE_match = ifelse(VISCODE == GeneX_VISCODE, 1, 0)) %>%
  1446. group_by(RID) %>%
  1447. mutate(min_match_date = min(date_diff)) %>%
  1448. # Same time = VISCODE_match; within 1-year = Edate_match
  1449. mutate(match_priority = case_when(
  1450. VISCODE_match == 1 ~ 1,
  1451. Edate_match == 1 ~ 2,
  1452. TRUE ~ 3
  1453. )) %>%
  1454. # Only keep matches
  1455. filter(match_priority <= 2) %>%
  1456. # If multiple matches, keep min date diff
  1457. filter(min_match_date == date_diff) %>%
  1458. # Replace any missing diagnosis, EXAMDATE, and VISCODE using the ones from gene expression pheno data
  1459. mutate(DXSUM = ifelse(is.na(DXSUM) & date_diff == 0 & !is.na(GeneX_DXSUM), GeneX_DXSUM, DXSUM),
  1460. VISCODE = ifelse(is.na(VISCODE) & date_diff == 0 & !is.na(GeneX_VISCODE), GeneX_VISCODE, VISCODE),
  1461. DX_Edate = ifelse(is.na(DX_Edate) & date_diff == 0 & !is.na(GeneX_DX_Edate), GeneX_DX_Edate, DX_Edate)
  1462. ) %>%
  1463. # Filter missing B, a convenient proxy for the replicates that were excluded during original DNA methylation preprocessing
  1464. filter(!is.na(B)) %>%
  1465. ungroup() %>%
  1466. dplyr::select(-c(Edate_match, VISCODE_match, match_priority, min_match_date, date_diff,
  1467. GeneX_DXSUM, GeneX_DX_Edate))
  1469. message("Now saving the final DNA methylation phenotype data for analysis...")
  1470. # Save the final DNA methylation phenotype data (same time or within 1-year of gene expression) for cross-sectional analysis
  1471. write.csv(merged, file.path(base_dir, "DNAm_pheno_final_shared.csv"), row.names = F)

phenotype.R at commit 5874e91, no license · at the source

Overview

Authors: Brendan A. Mitchell1,2, Isabella Hausle1, Sara Smith3,4, Pamela Thropp5, Marina Sirota6,7, Duygu Tosun1
  1. Department of Radiology and Biomedical Imaging, University of California,San Francisco, CA USA
  2. UC Berkeley - UCSF Graduate Program in Bioengineering,Berkeley, CA USA
  3. Department of Medicine, Division of Rheumatology, University of California,San Francisco, CA USA
  4. CoLabs, University of California,San Francisco, CA USA
  5. Northern California Institute for Research and Education (NCIRE),San Francisco, CA USA
  6. Bakar Computational Health Sciences Institute, University of California,San Francisco, CA USA
  7. Department of Pediatrics, University of California,San Francisco, CA USA
Institutions: University of California, San Francisco (United States); University of California, Berkeley (United States); Northern California Institute for Research and Education (United States)
Journal: NPJ dementia, volume 2, issue 1, article 74
Dates: received 28 January 2026; accepted 9 May 2026; published online 31 August 2026; in print 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1038/s44400-026-00100-z · PMID 42687894 · PMCID PMC13529590 · OpenAlex W7204824944
Open access: hybrid, a free copy (OpenAlex)
Status: code verified
Categories: genetics / omics (modality), human (organism), Alzheimer's / dementia (population)
Methods: Statistics, Smoothing, state filtering, decompositions, fMRI & imaging, Preprocessing
Keywords: Biomarkers, Computational biology and bioinformatics, Diseases, Neurology, Neuroscience
Topic: Alzheimer's disease research and treatments (Physiology, Medicine), according to OpenAlex
Funding: NSF | National Science Foundation Graduate Research Fellowship Program (GRFP) (2445150); National Institute on Aging (R01AG060393, U19AG024904)
Citations: not cited yet (Europe PMC); 83 references in the paper

Abstract

Leveraging multi-omics to better understand the molecular signatures and pathways underlying Alzheimer’s disease (AD) pathogenesis is critical for early diagnosis and disease modifying interventions. We performed peripheral blood transcriptome (N = 669) and epigenome microarray analyses (N = 553) on non-Hispanic white participants from the Alzheimer’s Disease Neuroimaging Initiative (ADNI) to identify molecular signatures of AD. We identified specific transcripts (e.g., MAPK14, GM2A, CD177) and co-expression networks that were dysregulated in AD, marked by a strong influence of APOE ε4 genotype, and characterized by a consistent pattern of immune activation, inflammation, and metabolic suppression. Further, these peripheral signatures were linked to central AD pathology (amyloid PET, CSF p-tau181) and neurodegeneration (plasma NfL, regional atrophy), with two genes, MXD3 and NR4A1, identified as protective against progression from MCI to AD. Our work emphasizes the importance of APOE genotypes in AD pathophysiology and highlights potential targets for biomarker discovery and personalized therapeutic strategies.

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

cind/adni-mrna-dnam-analysis

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: 5874e916559bb984d4614f4299af1611056a7fa8, 19 June 2026
Languages: R (16)
Size: 17 files, 16 scripts
Software Heritage: not archived
Found in: “Code availability”
Holds: README
Not found: license file, CITATION.cff, environment file, tests, continuous integration, documentation
Tools: tidyverse (15 files), WGCNA (9 files), limma (7 files), broom (4 files), reshape2 (4 files), ComplexHeatmap (3 files), survival (3 files), clusterProfiler (2 files), data.table (2 files), ggpubr (1 file), lme4 (1 file), pheatmap (1 file), UMAP (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
17 files
At the source: github.com/cind/adni-mrna-dnam-analysis/

Code availability

The code used in this manuscript can be found at https://github.com/cind/adni-mrna-dnam-analysis/.

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

Data used in the preparation of this article were obtained from the ADNI database (adni.loni.usc.edu). As such, the investigators within the ADNI contributed to the design and implementation of ADNI and/or provided data but did not participate in the analysis or writing of this report. A complete listing of ADNI investigators can be found at: http://adni.loni.usc.edu/wp-content/uploads/how_to_apply/ADNI_Acknowledgement_List.pdf.

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 2, 28 September 2026

  • Publisher: n/a → Springer Science+Business Media

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 6 authors, 5 keywords, 2 funders, 79 references.

Cite

This paper

Mitchell, B. A., Hausle, I., Smith, S., Thropp, P., Sirota, M., & Tosun, D. (2026). Peripheral blood microarray-based transcriptomic and epigenetic analyses identify immune, inflammation, and metabolic dysregulation in Alzheimer's disease. NPJ dementia, 2(1), 74. https://doi.org/10.1038/s44400-026-00100-z

BibTeX

@article{mitchell2026peripheral,
author = {Mitchell, Brendan A. and Hausle, Isabella and Smith, Sara and Thropp, Pamela and Sirota, Marina and Tosun, Duygu},
title = {{Peripheral blood microarray-based transcriptomic and epigenetic analyses identify immune, inflammation, and metabolic dysregulation in Alzheimer's disease}},
journal = {NPJ dementia},
year = {2026},
month = aug,
volume = {2},
number = {1},
pages = {74},
publisher = {Springer Science+Business Media},
issn = {3005-1940},
doi = {10.1038/s44400-026-00100-z},
url = {https://doi.org/10.1038/s44400-026-00100-z},
pmid = {42687894},
pmcid = {PMC13529590}
}

RIS

TY - JOUR
AU - Mitchell, Brendan A.
AU - Hausle, Isabella
AU - Smith, Sara
AU - Thropp, Pamela
AU - Sirota, Marina
AU - Tosun, Duygu
TI - Peripheral blood microarray-based transcriptomic and epigenetic analyses identify immune, inflammation, and metabolic dysregulation in Alzheimer's disease
T2 - NPJ dementia
J2 - NPJ Dement
PY - 2026
DA - 2026/08/31
VL - 2
IS - 1
SP - 74
SN - 3005-1940
PB - Springer Science+Business Media
DO - 10.1038/s44400-026-00100-z
UR - https://doi.org/10.1038/s44400-026-00100-z
LA - en
ER -

CSL-JSON

{
"id": "10.1038/s44400-026-00100-z",
"type": "article-journal",
"title": "Peripheral blood microarray-based transcriptomic and epigenetic analyses identify immune, inflammation, and metabolic dysregulation in Alzheimer's disease",
"container-title": "NPJ dementia",
"author": [
{
"family": "Mitchell",
"given": "Brendan A."
},
{
"family": "Hausle",
"given": "Isabella"
},
{
"family": "Smith",
"given": "Sara"
},
{
"family": "Thropp",
"given": "Pamela"
},
{
"family": "Sirota",
"given": "Marina"
},
{
"family": "Tosun",
"given": "Duygu"
}
],
"container-title-short": "NPJ Dement",
"volume": "2",
"issue": "1",
"page": "74",
"DOI": "10.1038/s44400-026-00100-z",
"PMID": "42687894",
"PMCID": "PMC13529590",
"ISSN": "3005-1940",
"publisher": "Springer Science+Business Media",
"URL": "https://doi.org/10.1038/s44400-026-00100-z",
"language": "en",
"issued": {
"date-parts": [
[
2026,
8,
31
]
]
}
}

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.1016/j.xcrm.2026.102682 [code]
TET CpG sequence-context-specific DNA demethylation shapes progression of IDH-mutant gliomas.
Journal: Cell reports. Medicine
In common: survival, limma, UMAP, 9 other tools, genetics / omics, 2 references
[2] doi:10.1093/bioinformatics/btag592 [code]
Network-based stratification of allele-specific expression reveals patient subgroups in Huntington's disease.
Journal: Bioinformatics (Oxford, England)
In common: survival, WGCNA, limma, 10 other tools, genetics / omics
[3] doi:10.1038/s42003-026-10957-8 [code]
Brain defence by the extracellular matrix protein Cochlin.
Journal: Communications biology
In common: survival, WGCNA, limma, 10 other tools
[4] doi:10.1016/j.xcrm.2026.102766 [code]
A longitudinal single-cell and spatial multiomic atlas of pediatric high-grade glioma.
Journal: Cell reports. Medicine
In common: WGCNA, limma, UMAP, 8 other tools, genetics / omics, 2 references
[5] doi:10.1016/j.isci.2026.115657 [code]
Integration of machine learning to develop a disulfidptosis model for predicting glioma prognosis, immunotherapy response, and drug.
Journal: iScience
In common: survival, WGCNA, limma, 8 other tools
[6] doi:10.1038/s41593-026-02367-0 [code]
A reproducible three-dimensional model of human brain tissue to investigate physiological and disease-associated microglia phenotypes.
Journal: Nature neuroscience
In common: WGCNA, limma, clusterProfiler, 6 other tools, Alzheimer's / dementia, 2 references
[7] doi:10.3390/ijms27104466 [code]
Uncovering the Key Circuit FOSL2/FOS/EGR3/EGR1, Contributing to the Hyperexcitability of Excitatory Neurons in the Epileptic Temporal Cortex and Hippocampus.
Journal: International journal of molecular sciences
In common: limma, UMAP, clusterProfiler, 7 other tools, genetics / omics, 1 reference
[8] doi:10.1038/s41380-026-03629-w [code]
Maternal fasting during early gestation induces epigenetic alterations and schizophrenia-related phenotypes.
Journal: Molecular psychiatry
In common: WGCNA, limma, clusterProfiler, 4 other tools, genetics / omics, 3 references
[9] doi:10.1038/s41467-026-77170-3 [code]
DNA methylation profiling identifies long-range epigenetic silencing of clustered protocadherins as a key determinant of meningioma progression.
Journal: Nature communications
In common: survival, limma, broom, 6 other tools, genetics / omics
[10] doi:10.1038/s41467-026-73305-8 [code]
Comparative analysis of the cellular landscape in mammalian striatum.
Journal: Nature communications
In common: WGCNA, clusterProfiler, ComplexHeatmap, 5 other tools, genetics / omics, 2 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.