Skip to content
Permalink
main
Switch branches/tags

Name already in use

A tag already exists with the provided branch name. Many Git commands accept both tag and branch names, so creating this branch may cause unexpected behavior. Are you sure you want to create this branch?
Go to file
 
 
Cannot retrieve contributors at this time
# Analysis script for the paper:
# "DNA methylation patterns of FKBP5 regulatory regions in brain and blood of
# humanized mice and humans"
# ==============================================================================
# Author: Natan Yusupov, natan_yusupov@psych.mpg.de
# Description:
# The script analyzes DNAm data in three functional intronic glucocorticoid-responsive
# elements (introns 2, 5 and 7) quantified by targeted bisulfite sequencing (HAM-TBS)
# in three tissues/brain regions (blood, prefrontal cortex and hippocampus) of mice
# carrying the human FK506-binding protein 5 (FKBP5) gene of FKBP5 at baseline,
# in cases of differing genotype (rs1360780 single nucleotide polymorphism), and
# following application of dexamethasone. Further analyses compare DNAm patterns
# in the humanized mouse to those in human peripheral blood and human postmortem
# brain prefrontal cortex.
# Summary of steps performed in the script "analysis_script.R":
# 1. Load libraries and data
# 2. Prepare data for analysis
# 3. Mean and SD in humanized mouse tissues (Supplementary Table 5)
# 4. Spearman correlation coefficients between CpGs of humanized mouse tissue
# (Supplementary Table 7 and Supplementary Figure 2)
# 5. Comparison DNAm in human blood at baseline study 1 and 2 (Supplementary Table 8)
# 6. Mean and SD in human and humanized mouse blood (Supplementary Table 9)
# 7. Delta mean of DNAm in blood between humanized mouse and humans (Supplementary Table 10)
# 8. Mean and SD in human and humanized mouse PFC (Supplementary Table 11)
# 9. Delta mean of DNAm in PFC between humanized mouse and humans (Supplementary Table 12)
# 10. Delta mean of DNAm in PFC and blood between humanized mouse and humans
# aged 20-29 years (Supplementary Table 13)
# 11. Multivariate linear regressions in humanized mouse (Supplementary Tables
# 14, 15, 16, and parts of Supplementary Figures 5, 8, 9)
# A. Blood time point 0
# B. Blood time point 4
# C. Blood time point 24 hours
# D. PFC time point 0
# E. PFC time point 4 hours
# F. PFC time point 24 hours
# G. HIP time point 0
# H. HIP time point 4 hours
# I. HIP time point 24 hours
# 12. Interindividual variability of DNAm in humanized mouse tissues (Supplementary Figure 1)
# 13. Correlation matrix of CpGs in the humanized mouse tissues (Supplementary Figure 3)
# 14. Genotype effects on DNAm in humanized mouse and human blood (Supplementary Figure 5)
# 15. Genotype effects on DNAm in humanized mouse and human PFC (Supplementary Figure 6)
# 16. Dexamethasone effects on DNAm in humanized mouse and human blood (Supplementary Figure 7)
# 17. Nominal interactions effects of dexamethasone with rs1360780 on DNAm
# in humanized mouse hippocampus (Supplementary Figure 9)
# 18. DNAm patterns in different tissues of humanized mouse (Figure 2)
# 19. DNAm patterns at baseline in blood and PFC of humanized mouse and humans (Figure 3)
# 20. Age dependent DNAm patterns in blood and PFC between humanized mouse and human (Figure 4)
# 21. Dexamethasone effects on DNAm per CpG in humanized mouse and human blood (Figure 5)
# 22. Session information
# ==============================================================================
# 1. Load libraries and data
library(tidyverse)
library(readxl)
library(writexl)
library(broom)
library(ggpubr)
library(Hmisc)
library(trio)
library(corrplot)
library(janitor)
library(reshape)
intron_2 <- readRDS("data/cpg_locations_intron2_vector.RDS")
intron_5 <- readRDS("data/cpg_locations_intron5_vector.RDS")
intron_7 <- readRDS("data/cpg_locations_intron7_vector.RDS")
cpg_locations <- readRDS("data/cpg_locations_vector.RDS")
path <- "/Users/natan_yusupov/Desktop/github_repos/Fkbp5_DNAm_HAMTBS_humanized_mouse_humans/"
human_blood_study1_df <- read.csv(paste0(path, "data/human_blood_study1_df.csv"),
sep = ",", row.names = 1, check.names = FALSE)
human_blood_study1_df <- human_blood_study1_df %>%
mutate(sex = as.factor(sex), rs1360780_T = as.factor(rs1360780_T), agebin = as.factor(agebin),
sample = paste0(sample, "_1"))
human_blood_study2_df <- read.csv(paste0(path, "data/human_blood_study2_df.csv"),
sep = ",", row.names = 1, check.names = FALSE)
human_blood_study2_df <- human_blood_study2_df %>%
mutate(Dex__0_baseline__1_3h__2_24h_ = as.factor(Dex__0_baseline__1_3h__2_24h_),
rs1360780_T = as.factor(rs1360780_T), sex = as.factor(sex),
dex_dichot = as.factor(dex_dichot),
sample = paste0(sample, "_2"))
human_postmortem_brain_study3_df <- read.csv(paste0(path, "data/human_postmortem_brain_study3_df.csv"),
sep = ",", row.names = 1, check.names = FALSE)
human_postmortem_brain_study3_df <- human_postmortem_brain_study3_df %>%
mutate(agebin = as.factor(agebin), sex = as.factor(sex),
rs1360780_T = as.factor(rs1360780_T), rs1360780_anyT = as.factor(rs1360780_anyT),
sample = paste0(sample, "_3"))
humanized_mouse_df <- read.csv(paste0(path, "data/humanized_mouse_df.csv"),
sep = ",", row.names = 1, check.names = FALSE)
humanized_mouse_df <- humanized_mouse_df %>%
mutate(Genotype = as.factor(Genotype), Timegroup = as.factor(Timegroup),
Group = as.factor(Group), Dissecter = as.factor(Dissecter), Tissue = as.factor(Tissue))
results_wiechmann_df <- read.csv(paste0(path, "data/results_wiechmann_df.csv"),
sep = ",", row.names = 1, check.names = FALSE)
# 2. Prepare data for analysis
methylation_human_blood_study1 <- human_blood_study1_df %>%
column_to_rownames(var = "sample") %>%
dplyr::select(contains(cpg_locations))
methylation_human_blood_study2_baseline <- human_blood_study2_df %>%
filter(dex_dichot == "0") # only baseline
methylation_humouse_blood <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with("_corrected")) %>%
filter(Tissue == "Blood", Group %in% c("No Treatment", "V")) %>% # only untreated
dplyr::select(sample, contains(cpg_locations)) %>%
column_to_rownames("sample") %>%
mutate(study = "Humanized Mouse")
covariates_humouse <- humanized_mouse_df %>%
dplyr::select(sample, Group, Genotype, Tissue)
methylation_blood_combined <- rbind(methylation_human_blood_study1, methylation_humouse_blood)
methylation_blood_combined_long <- methylation_blood_combined %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
methylation_human_postmortem_brain <- human_postmortem_brain_study3_df %>%
dplyr::select(sample, ends_with("_corrected")) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, contains(cpg_locations)) %>%
mutate(study = "Human") %>%
column_to_rownames("sample")
methylation_humouse_pfc <- humanized_mouse_df %>%
filter(Tissue == "PFC", Group %in% c("No Treatment", "V")) %>% # only untreated
dplyr::select(sample, ends_with("_corrected")) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, contains(cpg_locations)) %>%
column_to_rownames("sample") %>%
mutate(study = "Humanized Mouse")
methylation_pfc_combined <- rbind(methylation_human_postmortem_brain, methylation_humouse_pfc)
methylation_pfc_combined_long <- methylation_pfc_combined %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
# 3. Mean and SD in humanized mouse tissues (Supplementary Table 5)
methylation_humouse_all_tissues <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
filter(Group %in% c("No Treatment", "V")) %>% # only untreated
dplyr::select(sample, Tissue, contains(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected"))
methylation_humouse_all_tissues_long <- methylation_humouse_all_tissues %>%
pivot_longer(cols = -c(sample, Tissue), names_to = "CpG", values_to = "methylation")
methylation_humouse_all_tissues_mean_sd <- methylation_humouse_all_tissues %>%
group_by(Tissue) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -Tissue, names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation)
# 4. Spearman correlation coefficients between CpGs of humanized mouse tissue
# (Supplementary Table 7 and Supplementary Figure 2)
humouse_tissue_cpg_level <- methylation_humouse_all_tissues %>%
mutate(sample = str_sub(sample, 1, 6)) %>%
filter(Tissue == "Blood") %>%
dplyr::select(-Tissue) %>%
rename_at(.vars = vars(`35558386`:`35608022`), .funs = ~ paste0(., "_Blood")) %>%
left_join(methylation_humouse_all_tissues %>% mutate(sample = str_sub(sample, 1, 6)) %>%
filter(Tissue == "PFC") %>% dplyr::select(-Tissue) %>%
rename_at(.vars = vars(`35558386`:`35608022`), .funs = ~ paste0(., "_PFC"))) %>%
left_join(methylation_humouse_all_tissues %>% mutate(sample = str_sub(sample, 1, 6)) %>%
filter(Tissue == "HIP") %>% dplyr::select(-Tissue) %>%
rename_at(.vars = vars(`35558386`:`35608022`), .funs = ~ paste0(., "_HIP")))
corr_brain_blood_humouse <- humouse_tissue_cpg_level %>%
tibble::remove_rownames() %>%
column_to_rownames("sample") %>%
as.matrix() %>%
rcorr(type = "spearman")
corr_matrix_brain_blood_humouse <- corr_brain_blood_humouse$r
corr_matrix_brain_blood_humouse_df <- corr_matrix_brain_blood_humouse %>%
as.data.frame() %>%
rownames_to_column("CpG")
corr_brain_blood_humouse_pvalue <- corr_brain_blood_humouse$P
plot_corr_brain_blood_humouse <- corrplot(corr_brain_blood_humouse$r, method = "color",
p.mat = corr_brain_blood_humouse_pvalue,
insig = "pch", pch = "ns", sig.level = 0.05,
tl.cex = 0.4, pch.cex = 0.3, tl.col = "darkblue",
tl.srt = 45, mar=c(0,0,2,0))
cor(humouse_tissue_cpg_level$`35570224_Blood`, humouse_tissue_cpg_level$`35570224_PFC`, use="complete.obs", method = "spearman")
cor(humouse_tissue_cpg_level$`35570224_Blood`, humouse_tissue_cpg_level$`35570224_HIP`, use="complete.obs", method = "spearman")
cor(humouse_tissue_cpg_level$`35570224_PFC`, humouse_tissue_cpg_level$`35570224_HIP`, use="complete.obs", method = "spearman")
# 5. Comparison DNAm in human blood at baseline study 1 and 2 (Supplementary Table 8)
outcome_variables <- methylation_human_blood_study1 %>%
dplyr::select(starts_with("35")) %>%
colnames()
human_blood_study1_df_subset <- human_blood_study1_df %>%
dplyr::select(-rs1360780_T, -agebin) %>%
na.omit()
regressors <- c("age", "sex", "CD8T", "CD4T", "NK", "Bcell", "Mono", "Gran")
lm_tidy <- function() {
results <- list()
for (i in outcome_variables) {
f <- as.formula(paste0("`", i, "`", "~", paste(regressors, collapse = "+")))
tidy <- augment(lm(formula = f, data = human_blood_study1_df_subset, na.action = na.omit))
tidy <- tidy %>% dplyr::select(.resid)
colnames(tidy) <- as.character(i)
results[[i]] <- tidy
}
results
}
reisduals <- lm_tidy()
methylation_human_blood_study1_regressed_out <- as.data.frame(do.call(cbind, reisduals))
rownames(methylation_human_blood_study1_regressed_out) <- human_blood_study1_df_subset$sample
methylation_human_blood_study1_regressed_out <- methylation_human_blood_study1_regressed_out %>%
mutate(study = as.factor("Human Study 1"))
outcome_variables <- methylation_human_blood_study2_baseline %>%
dplyr::select(starts_with("35")) %>%
colnames()
regressors <- c("age", "sex", "Neutrophil", "Monocyte", "Bcell", "NKcell", "Tcell")
lm_tidy <- function() {
results <- list()
for (i in outcome_variables) {
f <- as.formula(paste0("`", i, "`", "~", paste(regressors, collapse = "+")))
tidy <- augment(lm(formula = f, data = methylation_human_blood_study2_baseline, na.action = na.omit))
tidy <- tidy %>% dplyr::select(.resid)
colnames(tidy) <- as.character(i)
results[[i]] <- tidy
}
results
}
reisduals <- lm_tidy()
methylation_human_blood_study2_regressed_out <- as.data.frame(do.call(cbind, reisduals))
rownames(methylation_human_blood_study2_regressed_out) <- methylation_human_blood_study2_baseline$sample
methylation_human_blood_study2_regressed_out <- methylation_human_blood_study2_regressed_out %>%
mutate(study = as.factor("Human Study 2"))
methylation_human_blood_study1 <- methylation_human_blood_study1 %>%
mutate(study = as.factor("Human Study 1"))
methylation_human_blood_study2_baseline <- methylation_human_blood_study2_baseline %>%
mutate(study = as.factor("Human Study 2"))
methylation_human_blood_2_studies <- methylation_human_blood_study2_baseline %>%
dplyr::select(sample, starts_with("35"), study) %>%
bind_rows(methylation_human_blood_study1 %>% rownames_to_column("sample")) %>%
relocate(study, .after = last_col())
methylation_human_blood_2_studies_long <- methylation_human_blood_2_studies %>%
dplyr::select(!starts_with("3560")) %>% # remove intron 2 since missing in study 2
pivot_longer(cols = -c(sample, study), names_to = "CpG", values_to = "methylation")
methylation_human_blood_2_studies_regressed_out <- methylation_human_blood_study2_regressed_out %>%
bind_rows(methylation_human_blood_study1_regressed_out) %>%
relocate(study, .after = last_col())
methylation_human_blood_2_studies_regressed_out_long <- methylation_human_blood_2_studies_regressed_out %>%
rownames_to_column("sample") %>%
dplyr::select(!starts_with("3560")) %>% # remove intron 2 since missing in study 2
pivot_longer(cols = -c(sample, study), names_to = "CpG", values_to = "methylation")
t_test_human_blood_studies <- compare_means(methylation ~ study,
data = methylation_human_blood_2_studies_long,
group.by = "CpG", method = "t.test") %>%
mutate(p.adj = p*16,
p.adj.signif = case_when(p < 0.05/16 & p >= 0.01/16 ~ "*",
p < 0.01/16 & p >= 0.001/16 ~ "**",
p < 0.001/16 & p >= 0.0001/16 ~ "***",
p < 0.0001/16 ~ "****",
p > 0.05/16 ~ "ns",
TRUE ~ "error"))
t_test_human_blood_studies_regressed_out <- compare_means(methylation ~ study,
data = methylation_human_blood_2_studies_regressed_out_long,
group.by = "CpG", method = "t.test")
# 6. Mean and SD in human and humanized mouse blood (Supplementary Table 9)
methylation_human_humouse_blood_mean_sd <- methylation_blood_combined %>%
group_by(study) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -study, names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation)
methylationn_human_blood_study2_mean_sd <- results_wiechmann_df %>%
filter(CpG %in% cpg_locations) %>%
dplyr::select(CpG, mean = "mean_baseline", sd = "sd_baseline") %>%
mutate(study = "Human Study 2")
methylation_blood_baseline_three_sets_df <- rbind(methylationn_human_blood_study2_mean_sd,
methylation_human_humouse_blood_mean_sd)
# 7. Delta mean of DNAm in blood between humanized mouse and humans (Supplementary Table 10)
delta_baseline_blood_three_sets_df <- methylation_blood_baseline_three_sets_df %>%
dplyr::select(-sd) %>%
pivot_wider(names_from = study, values_from = mean) %>%
dplyr::select(CpG, `Humanized Mouse`, `Human Study 1`, `Human Study 2`) %>%
mutate(`Delta DNAm Human Study 1` = `Humanized Mouse` - `Human Study 1`,
`Delta DNAm Human Study 2` = `Humanized Mouse` - `Human Study 2`)
# 8. Mean and SD in human and humanized mouse PFC (Supplementary Table 11)
methylation_human_postmortem_humouse_mean_sd <- methylation_pfc_combined %>%
group_by(study) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -study, names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation)
# 9. Delta mean of DNAm in PFC between humanized mouse and humans (Supplementary Table 12)
delta_baseline_pfc_human_humouse_df <- methylation_human_postmortem_humouse_mean_sd %>%
dplyr::select(-sd) %>%
pivot_wider(names_from = study, values_from = mean) %>%
dplyr::select(CpG, `Humanized Mouse`, Human) %>%
mutate(`Delta DNAm` = `Humanized Mouse` - Human)
# 10. Delta mean of DNAm in PFC and blood between humanized mouse and humans
# aged 20-29 years (Supplementary Table 13)
methylation_human_blood_agebins <- methylation_human_blood_study1 %>%
rownames_to_column("sample") %>%
left_join(human_blood_study1_df %>% dplyr::select(sample, agebin)) %>%
mutate(study = "Human")
methylation_blood_combined_agebins <- bind_rows(methylation_human_blood_agebins,
methylation_humouse_blood %>% rownames_to_column("sample"))
methylation_blood_combined_agebins_long <- methylation_blood_combined_agebins %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
methylation_human_humouse_blood_mean_sd_agebins <- methylation_blood_combined_agebins %>%
group_by(study, agebin) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -c(study, agebin), names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation) %>%
mutate(agebin = case_when(!is.na(agebin) ~ paste0("Human ", agebin),
TRUE ~ "Humanized Mouse"))
methylation_human_postmortem_brain_agebins <- methylation_human_postmortem_brain %>%
rownames_to_column("sample") %>%
left_join(human_postmortem_brain_study3_df %>% dplyr::select(sample, agebin), by = "sample")
methylation_pfc_combined_agebins <- bind_rows(methylation_human_postmortem_brain_agebins,
methylation_humouse_pfc %>% rownames_to_column("sample"))
methylation_pfc_combined_agebins_long <- methylation_pfc_combined_agebins %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
methylation_human_humouse_pfc_mean_sd_agebins <- methylation_pfc_combined_agebins %>%
group_by(study, agebin) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -c(study, agebin), names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation) %>%
mutate(agebin = case_when(!is.na(agebin) ~ paste0("Human ", agebin),
TRUE ~ "Humanized Mouse"))
delta_baseline_blood_human_humouse_agebins_df <- methylation_human_humouse_blood_mean_sd_agebins %>%
ungroup() %>%
dplyr::select(-c(sd, study)) %>%
filter(agebin %in% c("Human 20-29", "Humanized Mouse")) %>%
pivot_wider(names_from = agebin, values_from = mean) %>%
mutate(`Delta DNAm` = `Humanized Mouse` - `Human 20-29`)
delta_baseline_pfc_human_humouse_agebins_df <- methylation_human_humouse_pfc_mean_sd_agebins %>%
ungroup() %>%
dplyr::select(-c(sd, study)) %>%
filter(agebin %in% c("Human 20-29", "Humanized Mouse")) %>%
pivot_wider(names_from = agebin, values_from = mean) %>%
mutate(`Delta DNAm` = `Humanized Mouse` - `Human 20-29`)
# 11. Multivariate linear regressions in humanized mouse (Supplementary Tables 14, 15, 16,
# and parts of Supplementary Figures 5, 8, 9)
cpgs_remain_blood <- iqr_cpg_tissue %>% # filter out CpGs with IQR<1
filter(Tissue == "Blood") %>%
filter(IQR > 1) %>%
pull(CpG)
cpgs_remain_pfc <- iqr_cpg_tissue %>%
filter(Tissue == "PFC") %>%
filter(IQR > 1) %>%
pull(CpG)
cpgs_remain_hip <- iqr_cpg_tissue %>%
filter(Tissue == "HIP") %>%
filter(IQR > 1) %>%
pull(CpG)
covariates_mlr_humouse <- humanized_mouse_df %>%
dplyr::select(sample, Genotype, Timegroup, Group, Dissecter, column, Tissue)
methylation_humouse_mval <- humanized_mouse_df %>%
dplyr::select(sample, Tissue, ends_with("_mval")) %>%
rename_with(~str_remove(., "_mval"))
reduced_df_wide_blood <- methylation_humouse_mval %>%
filter(Tissue == "Blood") %>%
dplyr::select(sample, all_of(cpgs_remain_blood)) %>%
left_join(covariates_mlr_humouse, by = "sample")
reduced_df_wide_pfc <- methylation_humouse_mval %>%
filter(Tissue == "PFC") %>%
dplyr::select(sample, all_of(cpgs_remain_pfc)) %>%
left_join(covariates_mlr_humouse, by = "sample")
reduced_df_wide_hip <- methylation_humouse_mval %>%
filter(Tissue == "HIP") %>%
dplyr::select(sample, all_of(cpgs_remain_hip)) %>%
left_join(covariates_mlr_humouse, by = "sample")
reduced_df_mval_all_tissues <- dplyr::bind_rows(reduced_df_wide_blood, reduced_df_wide_pfc, reduced_df_wide_hip)
reduced_df_mval_all_tissues <- reduced_df_mval_all_tissues %>%
mutate(Group = fct_relevel(Group, c("No Treatment", "V", "DEX")),
Timegroup = fct_relevel(Timegroup, c("t0", "t4", "t24")),
Tissue = fct_relevel(Tissue, c("Blood", "PFC", "HIP")))
# A. Blood time point 0
outcome <- reduced_df_wide_blood %>%
dplyr::select(starts_with("35")) %>%
colnames() %>% sort()
predictor <- "Genotype"
lm_summary <- function(tissue, timegroup){ # function for MLR in different data subsets
data_sub <- reduced_df_mval_all_tissues %>% filter(Tissue == tissue, Timegroup %in% timegroup)
results <- list()
for(i in outcome){
f <- as.formula(paste0("`", i, "`", "~", paste(predictor, collapse = "*")))
tidy <- tidy(lm(formula = f, data = data_sub, na.action = na.omit))
results[[i]] <- tidy
}
results
}
blood_0 <- lm_summary("Blood", "t0")
blood_0_df <- purrr::map_dfr(blood_0, ~ .x, .id = "CpG") %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_blood_0 <- blood_0_df %>%
filter(term != "(Intercept)") %>%
ggplot(aes(x = estimate, y = CpG, fill = p.value < 0.05,
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1, "*", "")))) +
geom_col(width = 0.2) +
theme_bw() +
geom_text(vjust = 0.5, hjust = -1, size = 3) +
scale_fill_manual(values= c("grey", "black")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "Genotype Effect in Blood") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# B. Blood time point 4
predictor <- c("Genotype", "Group")
blood_4 <- lm_summary("Blood", "t4")
blood_4_df <- purrr::map_dfr(blood_4, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_blood_4 <- blood_4_df %>%
filter(!term %in% c("GenotypeRiA", "(Intercept)")) %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = ifelse(p.value < 0.05, 1, 0.2),
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = -1, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse(limit = c(0,-1)) +
labs(x = "Effect Size", y = "CpG", title = "Dexamethasone Effect in Blood (4 Hours)") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# C. Blood time point 24 hours
blood_24 <- lm_summary("Blood", "t24")
blood_24_df <- purrr::map_dfr(blood_24, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_blood_24 <- blood_24_df %>%
filter(!term %in% c("GenotypeRiA", "(Intercept)")) %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = ifelse(p.value < 0.05, 1, 0.2),
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = -1, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "After 24h - Blood") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# D. PFC time point 0
outcome <- reduced_df_wide_pfc %>%
dplyr::select(starts_with("35")) %>%
colnames() %>% sort()
predictor <- "Genotype"
covariates <- c("Dissecter", "column")
lm_summary <- function(tissue, timegroup){
data_sub <- final_df_humouse %>% filter(Tissue == tissue, Timegroup == timegroup)
results <- list()
for(i in outcome){
f <- as.formula(paste0("`", i, "`", "~", paste(covariates, collapse = "+"), "+",
paste(predictor, collapse = "*")))
tidy <- tidy(lm(formula = f, data = data_sub, na.action = na.omit))
results[[i]] <- tidy
}
results
}
pfc_0 <- lm_summary("PFC", "t0")
pfc_0_df <- purrr::map_dfr(pfc_0, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_pfc_0 <- pfc_0_df %>%
ggplot(aes(x = estimate, y = CpG, fill = p.value < 0.05,
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(values= c("grey", "black")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "Genotype Effect in PFC") +
theme(plot.title = element_text(size = 12, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 8),
axis.text.x = element_text(size = 8, hjust = 1, angle = 90),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
annotate("rect", ymin = -Inf, ymax = "35569751", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35578891", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# E. PFC time point 4 hours
predictor <- c("Genotype", "Group")
pfc_4 <- lm_summary("PFC", "t4")
pfc_4_df <- purrr::map_dfr(pfc_4, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_pfc_4 <- pfc_4_df %>%
filter(term != "GenotypeRiA") %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = ifelse(p.value < 0.05, 1, 0.2),
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "After 4h - PFC") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# F. PFC time point 24 hours
pfc_24 <- lm_summary("PFC", "t24")
pfc_24_df <- purrr::map_dfr(pfc_24, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_pfc_24 <- pfc_24_df %>%
filter(term != "GenotypeRiA") %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = ifelse(p.value < 0.05, 1, 0.2),
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "After 24h - PFC") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
pfc_0_df <- purrr::map_dfr(pfc_0, ~ .x, .id = "CpG") %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
pfc_4_df <- purrr::map_dfr(pfc_4, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
pfc_24_df <- purrr::map_dfr(pfc_24, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
# G. HIP time point 0
outcome <- reduced_df_wide_hip %>%
dplyr::select(starts_with("35")) %>%
colnames() %>% sort()
predictor <- "Genotype"
hip_0 <- lm_summary("HIP", "t0")
hip_0_df <- purrr::map_dfr(hip_0, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_hip_0 <- hip_0_df %>%
ggplot(aes(x = estimate, y = CpG, fill = p.value < 0.05,
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(values= c("grey", "black")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "Genotype Effect in HIP") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# H. HIP time point 4 hours
predictor <- c("Genotype", "Group")
hip_4 <- lm_summary("HIP", "t4")
hip_4_df <- purrr::map_dfr(hip_4, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_hip_4 <- hip_4_df %>%
filter(term != "GenotypeRiA") %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = p.value < 0.05,
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "Dexamethasone Effect in HIP (4 Hours)") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
# I. HIP time point 24 hours
hip_24 <- lm_summary("HIP", "t24")
hip_24_df <- purrr::map_dfr(hip_24, ~ .x, .id = "CpG") %>%
filter(!term %in% c("(Intercept)", "DissecterLi"), !str_detect(term, "column")) %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
plot_hip_24 <- hip_24_df %>%
filter(term != "GenotypeRiA") %>%
ggplot(aes(x = estimate, y = CpG, fill = term, alpha = p.value < 0.05,
label = ifelse(q.value < 0.05, "**", ifelse(q.value < 0.1 & p.value < 0.05, "*", "")))) +
geom_col(position = "dodge", width = 0.4) +
geom_text(vjust = 0.5, hjust = 1.5, size = 3) +
theme_bw() +
scale_fill_manual(labels = c("GenotypeRiA:DEX", "DEX"),
values = c("goldenrod1", "darkgreen")) +
scale_x_reverse() +
labs(x = "Effect Size", y = "CpG", title = "Dexamethasone Effect in HIP (24 Hours)") +
theme(plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12, hjust = 1, angle = 90),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
guides(color=guide_legend("alpha"), alpha = FALSE) +
annotate("rect", ymin = -Inf, ymax = "35558721", xmin = -Inf, xmax = Inf, fill = "lightgreen", alpha = 0.2) +
annotate("rect", ymin = "35569751", ymax = "35578891", xmin = -Inf, xmax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", ymin = "35607856", ymax = Inf, xmin = -Inf, xmax = Inf, fill = "red4", alpha = 0.2)
hip_0_df <- purrr::map_dfr(hip_0, ~ .x, .id = "CpG") %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
hip_4_df <- purrr::map_dfr(hip_4, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
hip_24_df <- purrr::map_dfr(hip_24, ~ .x, .id = "CpG") %>%
group_by(term) %>%
mutate(q.value = p.adjust(p.value, method = "fdr"))
# 12. Interindividual variability of DNAm in humanized mouse tissues (Supplementary Figure 1)
iqr_cpg_tissue <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, Tissue, contains(cpg_locations)) %>%
group_by(Tissue) %>%
summarise(across(-c(sample), IQR, na.rm = TRUE)) %>%
pivot_longer(cols = -Tissue, names_to = "CpG", values_to = "IQR") %>%
arrange(desc(IQR)) %>%
mutate(functional_region = as.factor(case_when(
CpG %in% intron_7 ~ "Intron 7",
CpG %in% intron_5 ~ "Intron 5",
CpG %in% intron_2 ~ "Intron 2",
TRUE ~ "error")))
mean_iqr_cpg_tissue <- iqr_cpg_tissue %>%
group_by(Tissue, functional_region) %>%
summarise(Mean_IQR = mean(IQR))
iqr_cpg_cutoff <- iqr_cpg_tissue %>%
mutate(IQR = ifelse(CpG == "35607969" & Tissue == "Blood", 0.4, IQR)) %>%
# change for visualization, otherwise rounds up and gets into wrong bin by ggplot
ggplot(aes(x = IQR)) +
geom_histogram(boundary=0.997, bins = 30) +
geom_vline(color = "darkred", xintercept = 1) +
theme_bw() +
scale_x_continuous(breaks = seq(0, 21, by = 1)) +
theme(plot.title = element_text(size = 12, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 10),
axis.text.y = element_text(size = 10, hjust = 1),
axis.text.x = element_text(size = 10, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
facet_wrap(~Tissue, nrow = 3)
iqr_cpg_by_region_blood <- iqr_cpg_tissue %>%
filter(Tissue == "Blood") %>%
ggplot(aes(x = reorder(CpG, -IQR), y = IQR, fill = functional_region)) +
geom_col() +
geom_hline(color = "darkred", yintercept = 1) +
labs(y = "IQR", x = "CpG", title = "IQR of DNAm Levels Across CpGs - Blood") +
theme_bw() +
theme(plot.title = element_text(size = 12, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 10),
axis.text.x = element_text(size = 8, hjust = 1, angle = 45),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
scale_y_continuous(breaks = seq(0, 21, by = 1), limits = c(0,21)) +
scale_fill_manual("legend", values = alpha(c("Intron 7" = "lightgreen",
"Intron 5" = "royalblue4",
"Intron 2" = "red4"), 0.6))
iqr_cpg_by_region_pfc <- iqr_cpg_tissue %>%
filter(Tissue == "PFC") %>%
ggplot(aes(x = reorder(CpG, -IQR), y = IQR, fill = functional_region)) +
geom_col() +
geom_hline(color = "darkred", yintercept = 1) +
labs(y = "IQR", x = "CpG", title = "IQR of DNAm Levels Across CpGs - PFC") +
theme_bw() +
theme(plot.title = element_text(size = 12, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 10),
axis.text.x = element_text(size = 8, hjust = 1, angle = 45),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
scale_y_continuous(breaks = seq(0, 21, by = 1), limits = c(0,21)) +
scale_fill_manual("legend", values = alpha(c("Intron 7" = "lightgreen",
"Intron 5" = "royalblue4",
"Intron 2" = "red4"), 0.6))
iqr_cpg_by_region_hip <- iqr_cpg_tissue %>%
filter(Tissue == "HIP") %>%
ggplot(aes(x = reorder(CpG, -IQR), y = IQR, fill = functional_region)) +
geom_col() +
geom_hline(color = "darkred", yintercept = 1) +
labs(y = "IQR", x = "CpG", title = "IQR of DNAm Levels Across CpGs - HIP") +
theme_bw() +
theme(plot.title = element_text(size = 12, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 10),
axis.text.x = element_text(size = 8, hjust = 1, angle = 45),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "bottom") +
scale_y_continuous(breaks = seq(0, 21, by = 1), limits = c(0,21)) +
scale_fill_manual("legend", values = alpha(c("Intron 7" = "lightgreen",
"Intron 5" = "royalblue4",
"Intron 2" = "red4"), 0.6))
# 13. Correlation matrix of CpGs in the humanized mouse tissues (Supplementary Figure 3)
methylation_humouse_all_tissues <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with("_corrected")) %>%
dplyr::select(sample, Tissue,contains(cpg_locations))
methylation_humouse_all_tissues_long <- methylation_humouse_all_tissues %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
reduced_df_long_blood <- methylation_humouse_all_tissues_long %>%
filter(Tissue == "Blood", CpG %in% cpgs_remain_blood) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
reduced_df_long_pfc <- methylation_humouse_all_tissues_long %>%
filter(Tissue == "PFC", CpG %in% cpgs_remain_pfc) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
reduced_df_long_hip <- methylation_humouse_all_tissues_long %>%
filter(Tissue == "HIP", CpG %in% cpgs_remain_hip) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
methylation_humouse_all_tissues_corrected <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, Tissue, contains(cpg_locations))
methylation_humouse_all_tissues_corrected_long <- methylation_humouse_all_tissues_corrected %>%
pivot_longer(cols = "35558386":"35608022", names_to = "CpG", values_to = "methylation")
reduced_df_long_corrected_blood <- methylation_humouse_all_tissues_corrected_long %>%
filter(Tissue == "Blood", CpG %in% cpgs_remain_blood) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
reduced_df_long_corrected_pfc <- methylation_humouse_all_tissues_corrected_long %>%
filter(Tissue == "PFC", CpG %in% cpgs_remain_pfc) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
reduced_df_long_corrected_hip <- methylation_humouse_all_tissues_corrected_long %>%
filter(Tissue == "HIP", CpG %in% cpgs_remain_hip) %>%
left_join(covariates_humouse %>% dplyr::select(-Tissue), by = "sample")
r <- reduced_df_long_blood %>%
pivot_wider(id_cols = sample, names_from = CpG, values_from = methylation) %>%
column_to_rownames("sample") %>%
as.matrix() %>%
rcorr(type = "pearson")
r_pvalue <- r$P
p_corr <- corrplot(r$r, method = "color", type = "upper", p.mat = r_pvalue, insig = "pch",
pch = "ns", sig.level = 0.05, tl.cex = 0.8, pch.cex = 0.5, tl.col = "darkblue",
tl.srt = 45, title = "Blood", mar=c(0,0,2,0))
r <- reduced_df_long_corrected_pfc %>%
pivot_wider(id_cols = sample, names_from = CpG, values_from = methylation) %>%
column_to_rownames("sample") %>%
as.matrix() %>%
rcorr(type = "pearson")
r_pvalue <- r$P
p_corr <- corrplot(r$r, method = "color", type = "upper", p.mat = r_pvalue, insig = "pch",
pch = "ns", sig.level = 0.05, tl.cex = 0.8, pch.cex = 0.5, tl.col = "darkblue",
tl.srt = 45, title = "Prefrontal Cortex", mar=c(0,0,2,0))
r <- reduced_df_long_corrected_hip %>%
pivot_wider(id_cols = sample, names_from = CpG, values_from = methylation) %>%
column_to_rownames("sample") %>%
as.matrix() %>%
rcorr(type = "pearson")
r_pvalue <- r$P
p_corr <- corrplot(r$r, method = "color", type = "upper", p.mat = r_pvalue, insig = "pch",
pch = "ns", sig.level = 0.05, tl.cex = 0.8, pch.cex = 0.5, tl.col = "darkblue",
tl.srt = 45, title = "Hippocampus", mar=c(0,0,2,0))
# 14. Genotype effects on DNAm in humanized mouse and human blood (Supplementary Figure 5)
delta_genotype_humouse_blood <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
filter(Tissue == "Blood", Group %in% c("No Treatment", "V")) %>% # only untreated
dplyr::select(sample, Genotype, contains(c(intron_7, intron_5))) %>%
group_by(Genotype) %>%
summarise(across("35558386":"35578891",
.fns = ~ mean(.x, na.rm = TRUE))) %>%
t() %>%
as.data.frame() %>%
dplyr::slice(-1) %>%
dplyr::rename(`C/G` = "V1", `A/T` = "V2") %>%
mutate(`C/G` = as.numeric(`C/G`), `A/T` = as.numeric(`A/T`),
delta_genotype_humouse = `A/T` - `C/G`) %>%
rownames_to_column("CpG") %>%
dplyr::select(-`A/T`, -`C/G`) %>%
dplyr::rename(`Humanized Mouse` = "delta_genotype_humouse") %>%
t() %>% row_to_names(row_number = 1)
delta_genotype_human_blood <- methylation_human_blood_study2_baseline %>%
dplyr::select(sample, Genotype = rs1360780_T, contains(cpg_locations)) %>%
filter(Genotype %in% c(0, 2)) %>%
group_by(Genotype) %>%
summarise(across("35558386":"35578891",
.fns = ~ mean(.x, na.rm = TRUE))) %>%
t() %>%
as.data.frame() %>%
dplyr::slice(-1) %>%
dplyr::rename(`C/G` = "V1", `A/T` = "V2") %>%
mutate(`C/G` = as.numeric(`C/G`), `A/T` = as.numeric(`A/T`),
delta_genotype_human = `A/T` - `C/G`) %>%
rownames_to_column("CpG") %>%
dplyr::select(-`A/T`, -`C/G`) %>%
dplyr::rename(Human = "delta_genotype_human") %>%
t() %>% row_to_names(row_number = 1)
delta_genotype_human_humouse_blood <- rbind(delta_genotype_humouse_blood, delta_genotype_human_blood)
heatmap_genotype_blood <- delta_genotype_human_humouse_blood %>%
melt() %>% # reshape data as needed for heatmap
mutate(X2 = as.character(X2), value = as.numeric(value)) %>%
dplyr::rename(Percentage = "value") %>%
mutate(Percentage = round(Percentage, digits = 1)) %>%
ggplot(aes(x = X2, y = X1)) +
geom_tile(color = "black", aes(fill = Percentage)) +
scale_fill_gradient2(midpoint = 0, limits = c(-10, 10), space = "Lab", name="DNA methylation [%]") +
geom_text(aes(label = Percentage), size = 5, colour = "black", angle = 90) +
labs(title = Delta ~ "DNA Methylation of Risk Allele Homozygosity - Blood") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, size = 14),
axis.text.y = element_text(angle = 90, hjust = 0.5, size = 20),
axis.title = element_blank(),
plot.title = element_text(face = "bold", size = 24),
legend.title = element_text(face = "bold", size = 16),
legend.text= element_text(size=16),
legend.position = "bottom") +
scale_x_discrete(expand=c(0,0)) + scale_y_discrete(expand=c(0,0))
# 15. Genotype effects on DNAm in humanized mouse and human PFC (Supplementary Figure 6)
delta_genotype_humouse_pfc <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
filter(Tissue == "PFC", Group %in% c("No Treatment", "V")) %>% # only untreated
dplyr::select(sample, Timegroup, Genotype, contains(cpg_locations)) %>%
group_by(Genotype) %>%
summarise(across("35558386":"35608022",
.fns = ~ mean(.x, na.rm = TRUE))) %>%
t() %>% as.data.frame() %>% dplyr::slice(-1) %>%
dplyr::rename(`C/G` = "V1", `A/T` = "V2") %>%
mutate(`C/G` = as.numeric(`C/G`), `A/T` = as.numeric(`A/T`),
delta_genotype_humouse = `A/T` - `C/G`) %>%
rownames_to_column("CpG") %>%
dplyr::select(-`A/T`, -`C/G`) %>%
dplyr::rename(`Humanized Mouse` = "delta_genotype_humouse") %>%
t() %>% row_to_names(row_number = 1)
delta_genotype_human_pfc <- human_postmortem_brain_study3_df %>%
dplyr::select(-ends_with("_mval")) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, Genotype = rs1360780_T, contains(cpg_locations)) %>%
filter(Genotype %in% c("CC", "TT")) %>%
group_by(Genotype) %>%
summarise(across("35558386":"35608022",
.fns = ~ mean(.x, na.rm = TRUE))) %>%
t() %>%
as.data.frame() %>%
dplyr::slice(-1) %>%
dplyr::rename(`C/G` = "V1", `A/T` = "V2") %>%
mutate(`C/G` = as.numeric(`C/G`), `A/T` = as.numeric(`A/T`),
delta_genotype_human = `A/T` - `C/G`) %>%
rownames_to_column("CpG") %>%
dplyr::select(-`A/T`, -`C/G`) %>%
dplyr::rename(Human = "delta_genotype_human") %>%
t() %>% row_to_names(row_number = 1)
delta_genotype_human_humouse_pfc <- rbind(delta_genotype_humouse_pfc, delta_genotype_human_pfc)
heatmap_genotype_pfc <- delta_genotype_human_humouse_pfc %>%
melt() %>%
mutate(X2 = as.character(X2), value = as.numeric(value)) %>%
dplyr::rename(Percentage = "value") %>%
mutate(Percentage = round(Percentage, digits = 1)) %>%
ggplot(aes(x = X2, y = X1)) +
geom_tile(color = "black", aes(fill = Percentage)) +
geom_text(aes(label = Percentage), size = 5, colour = "black", angle = 90) +
scale_fill_gradient2(midpoint = 0, limits = c(-10, 10), space = "Lab", name="DNA methylation [%]") +
labs(title = Delta ~ "DNA Methylation of Risk Allele Homozygosity - PFC") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, size = 14),
axis.text.y = element_text(angle = 90, hjust = 0.5, size = 20),
axis.title = element_blank(),
plot.title = element_text(face = "bold", size = 24),
legend.title = element_text(face = "bold", size = 16),
legend.text= element_text(size=16),
legend.position = "bottom") +
scale_x_discrete(expand=c(0,0)) + scale_y_discrete(expand=c(0,0))
# 16. Dexamethasone effects on DNAm in humanized mouse and human blood (Supplementary Figure 7)
delta_dex_human_humouse_blood <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
filter(Tissue == "Blood", Timegroup == c("t0", "t4"),
Group %in% c("No Treatment", "DEX")) %>%
dplyr::select(sample, Timegroup, Group, contains(c(intron_7, intron_5))) %>%
group_by(Timegroup, Group) %>%
summarise(across("35558386":"35578891",
.fns = ~ mean(.x, na.rm = TRUE))) %>%
t() %>%
as.data.frame() %>%
dplyr::slice(-1,-2) %>%
dplyr::rename(t0 = "V1", t4 = "V2") %>%
mutate(t0 = as.numeric(t0), t4 = as.numeric(t4), delta_dex_humouse = t4 - t0) %>%
rownames_to_column("CpG") %>%
dplyr::select(-t0, -t4) %>%
left_join(results_wiechmann_df %>% dplyr::select(CpG, delta_dex_human) %>%
mutate(CpG = as.character(CpG))) %>%
dplyr::rename(`Humanized Mouse` = "delta_dex_humouse", Human = "delta_dex_human") %>%
t() %>% row_to_names(row_number = 1)
heatmap_dex_blood <- delta_dex_human_humouse_blood %>%
melt() %>%
mutate(X2 = as.character(X2), value = as.numeric(value)) %>%
dplyr::rename(Percentage = "value") %>%
mutate(Percentage = round(Percentage, digits = 1)) %>%
ggplot(aes(x = X2, y = X1)) +
geom_tile(color = "black", aes(fill = Percentage)) +
scale_fill_gradient2(midpoint = 0, limits = c(-26, 26), space = "Lab", name="DNA methylation [%]") +
geom_text(aes(label = Percentage), size = 5, colour = "black", angle = 90) +
labs(title = Delta ~ "DNA Methylation After Dexamethasone") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, size = 14),
axis.text.y = element_text(angle = 90, hjust = 0.5, size = 20),
axis.title = element_blank(),
plot.title = element_text(face = "bold", size = 24),
legend.title = element_text(face = "bold", size = 16),
legend.text= element_text(size=10),
legend.position = "bottom") +
scale_x_discrete(expand=c(0,0)) + scale_y_discrete(expand=c(0,0))
# 17. Nominal interactions effects of dexamethasone with rs1360780 on DNAm
# in humanized mouse hippocampus (Supplementary Figure 9)
methylation_humouse_all_tissues_complete <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, Tissue, Genotype, Group, Timegroup, contains(cpg_locations)) %>%
mutate(Group = fct_relevel(Group, c("No Treatment", "V", "DEX")),
Timegroup = fct_relevel(Timegroup, c("t0", "t4", "t24")),
Tissue = fct_relevel(Tissue, c("Blood", "PFC", "HIP")))
intron_5_35570224_effects <- ggplot(methylation_humouse_all_tissues_complete, aes(x = Group, y = `35570224`)) +
geom_boxplot(aes(fill = Genotype), outlier.size = 1) +
geom_point(aes(fill = Genotype), position = position_jitterdodge(jitter.height = 0, jitter.width = 0.2), size = 1) +
facet_grid(Tissue ~ Timegroup, scales = "free") +
scale_fill_brewer(palette = "Paired") +
labs(title = "Intron 5 - 35570224", y = "DNA Methylation [%]") +
theme_bw() +
theme(legend.position = "bottom",
plot.title = element_text(size = 16, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text.x = element_text(size = 12),
axis.text.y = element_text(size = 12, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 1,
linetype = "solid", colour = "darkgray")) +
stat_compare_means(aes(group = Genotype), method = "t.test", label = "p.signif") +
scale_y_continuous(expand = expansion(mult = c(0.05, 0.15)))
# 18. DNAm patterns in different tissues of humanized mouse (Figure 2)
intron_7_humouse_all_tissues <- methylation_humouse_all_tissues_mean_sd %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = Tissue, group = Tissue)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("red", "black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 7") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_humouse_all_tissues
intron_5_humouse_all_tissues <- methylation_humouse_all_tissues_mean_sd %>%
filter(CpG %in% intron_5) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = Tissue, group = Tissue)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("red", "black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 5") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = "35570224", ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", xmin = "35578739", xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
intron_5_humouse_all_tissues
intron_2_humouse_all_tissues <- methylation_humouse_all_tissues_mean_sd %>%
filter(CpG %in% intron_2) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = Tissue, group = Tissue)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("red", "black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 2") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "red4", alpha = 0.2)
intron_2_humouse_all_tissues
# 19. DNAm patterns at baseline in blood and PFC of humanized mouse and humans (Figure 3)
intron_7_blood_human_humouse <- methylation_blood_baseline_three_sets_df %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "darkgrey", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 7") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_blood_human_humouse
intron_5_blood_human_humouse <- methylation_blood_baseline_three_sets_df %>%
filter(CpG %in% intron_5) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "darkgrey", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 5") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = "35570224", ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", xmin = "35578739", xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
intron_5_blood_human_humouse
intron_2_blood_human_humouse <- methylation_blood_baseline_three_sets_df %>%
filter(CpG %in% intron_2) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 2") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "red4", alpha = 0.2)
intron_2_blood_human_humouse
intron_7_pfc_humouse_human <- methylation_human_postmortem_humouse_mean_sd %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 7") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_pfc_humouse_human
intron_5_pfc_humouse_human <- methylation_human_postmortem_humouse_mean_sd %>%
filter(CpG %in% intron_5) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 5") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = "35570224", ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", xmin = "35578739", xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
intron_5_pfc_humouse_human
intron_2_pfc_humouse_human <- methylation_human_postmortem_humouse_mean_sd %>%
filter(CpG %in% intron_2) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = study, group = study)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_manual(values=c("black", "dodgerblue3")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2) +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 2") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),legend.text = element_text(size=15),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "red4", alpha = 0.2)
intron_2_pfc_humouse_human
# 20. Age dependent DNAm patterns in blood and PFC between humanized mouse and human (Figure 4)
intron_7_blood_humouse_human_agebins <- methylation_human_humouse_blood_mean_sd_agebins %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = agebin, group = agebin)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_brewer(palette = "Dark2") +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 7") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=10),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_blood_humouse_human_agebins
intron_5_blood_humouse_human_agebins <- methylation_human_humouse_blood_mean_sd_agebins %>%
filter(CpG %in% intron_5) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = agebin, group = agebin)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_brewer(palette = "Dark2") +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 5") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=10),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = "35570224", ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", xmin = "35578739", xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
intron_5_blood_humouse_human_agebins
intron_7_pfc_humouse_human_agebins <- methylation_human_humouse_pfc_mean_sd_agebins %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = agebin, group = agebin)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_brewer(palette = "Dark2") +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 7") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=10),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_pfc_humouse_human_agebins
intron_5_pfc_humouse_human_agebins <- methylation_human_humouse_pfc_mean_sd_agebins %>%
filter(CpG %in% intron_5) %>%
ggplot(aes(x = as.factor(CpG), y = mean, col = agebin, group = agebin)) +
geom_point(alpha = 0.5) +
geom_line(size = 0.7) +
ylim(0,100) +
scale_color_brewer(palette = "Dark2") +
labs(x = "CpG", y = "DNA Methylation [%]", title = "Intron 5") +
theme_bw() +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_text(size = 16),
axis.title.y = element_text(size = 16),
axis.text.x = element_text(size = 14, hjust = 1, angle = 45),
axis.text.y = element_text(size = 14, hjust = 1),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(), legend.text = element_text(size=10),
legend.position = "bottom") +
guides(color = guide_legend(override.aes = list(size = 3, alpha = 0.5))) +
annotate("rect", xmin = -Inf, xmax = "35570224", ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2) +
annotate("rect", xmin = "35578739", xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
intron_5_pfc_humouse_human_agebins
# 21. Dexamethasone effects on DNAm per CpG in humanized mouse and human blood (Figure 5)
course_human_df <- human_blood_study2_df %>%
pivot_longer(cols = `35558386`:`35578891`, names_to = "CpG", values_to = "value") %>%
mutate(Group = case_when(dex_dichot == "1" ~ "Dexamethasone",
dex_dichot == "0" ~ "Non",
TRUE ~ "error"),
Timegroup = case_when(Dex__0_baseline__1_3h__2_24h_ == "0" ~ "t0",
Dex__0_baseline__1_3h__2_24h_ == "1" ~ "t3",
Dex__0_baseline__1_3h__2_24h_ == "2" ~ "t24",
TRUE ~ "error")) %>%
dplyr::select(CpG, Timegroup, Group, value) %>%
filter(CpG %in% c(intron_7, intron_5)) %>%
group_by(Timegroup, Group, CpG) %>%
summarise(mean = mean(value), sd = sd(value)) %>%
mutate(Timegroup = factor(Timegroup, levels = c("t0", "t3", "t24")),
Group = as.factor(Group))
intron_5_blood_human_part1 <- course_human_df %>%
filter(CpG %in% c("35569751", "35569757", "35569777", "35569896", "35569922", "35570224")) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
scale_color_manual(values=c("black", "darkblue")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(limits = c(0,14), breaks = seq(0, 14, by = 2)) +
theme_bw() +
facet_wrap(~CpG) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=1)
intron_5_blood_human_part2 <- course_human_df %>%
filter(CpG %in% c("35578739", "35578830", "35578891")) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG, ncol = 2) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=8)
intron_7_blood_human <- course_human_df %>%
filter(CpG %in% intron_7) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
course_mouse_df <- humanized_mouse_df %>%
dplyr::select(-ends_with("_mval"), -ends_with(cpg_locations)) %>%
rename_with(~str_remove(., "_corrected")) %>%
dplyr::select(sample, Tissue, Timegroup, Group, contains(cpg_locations)) %>%
group_by(Tissue, Timegroup, Group) %>%
summarise(across("35558386":"35608022",
.fns = list(mean = ~ mean(.x, na.rm = TRUE), sd = ~ sd(.x, na.rm = TRUE)),
.names = "{.fn}_{.col}")) %>%
pivot_longer(cols = -c(Tissue, Timegroup, Group), names_to = "which" , values_to = "methylation") %>%
separate(which, into = c("which","CpG"), sep = "_") %>%
spread(key=which, value=methylation) %>%
mutate(Group = as.factor(ifelse(Group %in% c("V", "No Treatment"), "Non", "Dexamethasone")),
Timegroup = fct_relevel(Timegroup, c("t0", "t4", "t24")))
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=2)
intron_2_blood_mouse <- course_mouse_df %>%
filter(Tissue == "Blood", CpG %in% intron_2) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group, col = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
scale_color_manual(values=c("black", "darkred")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.y = element_text(size = 16, face = "bold"),
axis.title.x = element_blank(),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "red4", alpha = 0.2)
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=4)
intron_5_blood_mouse_part1 <- course_mouse_df %>%
filter(Tissue == "Blood",
CpG %in% c("35569751", "35569757", "35569777", "35569896", "35569922", "35570224")) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group, col = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
scale_color_manual(values=c("black", "darkblue")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=1)
intron_5_blood_mouse_part2 <- course_mouse_df %>%
filter(Tissue == "Blood", CpG %in% c("35578739", "35578830", "35578891")) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group, col = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
scale_color_manual(values=c("black", "darkblue")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG, ncol = 2) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "royalblue4", alpha = 0.2)
f <- function(y) seq(floor(min(y)), ceiling(max(y)),by=2)
intron_7_blood_mouse <- course_mouse_df %>%
filter(Tissue == "Blood", CpG %in% intron_7) %>%
ggplot(aes(x = Timegroup, y = mean, group = Group, col = Group)) +
geom_point(alpha = 0.5, position = position_dodge(width = 0.2)) +
scale_color_manual(values=c("black", "darkgreen")) +
geom_errorbar(aes(ymin = mean-sd, ymax = mean+sd), width = .2, position = position_dodge(width = 0.2)) +
labs(y = "Mean DNA Methylation [%]") +
scale_y_continuous(breaks=f) +
theme_bw() +
facet_wrap(~CpG) +
theme(plot.title = element_text(size = 20, face = "bold", hjust = 0.5),
axis.title.x = element_blank(),
axis.title.y = element_text(size = 16, face = "bold"),
axis.text.x = element_text(size = 14),
axis.text.y = element_text(size = 12),
legend.background = element_rect(fill = "lightgray", size = 0.5,
linetype = "solid", colour = "darkgray"),
legend.title = element_blank(),
legend.position = "none") +
annotate("rect", xmin = -Inf, xmax = Inf, ymin = -Inf, ymax = Inf, fill = "lightgreen", alpha = 0.2)
intron_7_blood_mouse
# 22. Session information
utils:::print.sessionInfo(sessionInfo()[-8])
# End of script