This script saves dot plots of the selected CpG, gene, and trait in muscle for Figure 4B.
Load packages
library(tidyverse)
library(ggrepel)
library(GenomicRanges)
library(ggpubr)
library(DNAmArray)
library(SummarizedExperiment)
Load top CpG results and methData
load('../GOTO_Data/GOTO_methData-filtered.Rdata')
load('../GOTO_Data/GOTO_results-top-muscle-adj.Rdata')
Save only adipose samples
methData <- methData[ , methData$tissue == 'fat']
Convert betas to data frame
beta_rows <- methData$Basename
betas <- as.data.frame(t(assay(methData)))
rownames(betas) <- beta_rows
Save selected CpGs
betas <- betas %>% dplyr::select(all_of(c("cg02649849",
"cg14434922",
"cg14504259")))
dim(betas)
## [1] 178 3
cpg_list <- colnames(betas)
Load data
load("../GOTO_Data/GOTO_targets-filtered.Rdata")
Filter adipose
targets <- targets %>%
filter(tissue == 'fat')
Save covariates
targets <- targets %>%
dplyr::select(IOP2_ID, timepoint, age, sex, smoke,
plate, array_row, Basename,
bmi, wc, perc_body_fat, fasting_insulin,
leptin, adiponectin, interleukin_6,
hdl_size_fasted, hdl_c, systolic_bp) %>%
mutate(ID = paste0(IOP2_ID, '_', timepoint))
ID_list <- targets %>% dplyr::select(ID)
Bind data frames
limma_base <- cbind(targets,
betas)
dim(limma_base)
## [1] 178 22
Make box plot
plot <- limma_base %>%
ggplot(aes(x = timepoint, y = cg14504259,
fill = timepoint)) +
geom_line(aes(group = IOP2_ID), color = '#DDDDDD', linewidth=1) +
geom_dotplot(binaxis = "y", binwidth = 0.02,
stackdir = "center",
dotsize = 0.9, stroke = 2) +
scale_fill_manual(values=c('#EE6677', '#66CCEE'))+
ggtitle('') +
ylab('') + xlab('') +
scale_x_discrete(expand = c(0,0)) +
theme(
panel.background = element_rect(fill = 'white',
color = 'black'),
legend.position = 'none',
plot.title =element_text(hjust=0.5),
text=element_text(size=16))
print(plot)
Save
png('../GOTO_Data/Figures/Figure_4B-cg14504259.png')
print(plot)
dev.off()
## png
## 2
Make box plot
plot <- limma_base %>%
ggplot(aes(x = timepoint, y = cg14434922,
fill = timepoint)) +
geom_line(aes(group = IOP2_ID), color = '#DDDDDD', linewidth=1) +
geom_dotplot(binaxis = "y", binwidth = 0.015,
stackdir = "center",
dotsize = 0.8, stroke = 2) +
scale_fill_manual(values=c('#EE6677', '#66CCEE'))+
ggtitle('') +
ylab('') + xlab('') +
scale_x_discrete(expand = c(0,0)) +
theme(
panel.background = element_rect(fill = 'white',
color = 'black'),
legend.position = 'none',
plot.title =element_text(hjust=0.5),
text=element_text(size=16))
print(plot)
Save
png('../GOTO_Data/Figures/Figure_4B-cg14434922.png')
print(plot)
dev.off()
## png
## 2
Make box plot
plot <- limma_base %>%
ggplot(aes(x = timepoint, y = cg02649849,
fill = timepoint)) +
geom_line(aes(group = IOP2_ID), color = '#DDDDDD', linewidth=1) +
geom_dotplot(binaxis = "y", binwidth = 0.015,
stackdir = "center",
dotsize = 0.8, stroke = 2) +
scale_fill_manual(values=c('#EE6677', '#66CCEE'))+
ggtitle('') +
ylab('') + xlab('') +
scale_x_discrete(expand = c(0,0)) +
theme(
panel.background = element_rect(fill = 'white',
color = 'black'),
legend.position = 'none',
plot.title =element_text(hjust=0.5),
text=element_text(size=16))
print(plot)
Save
png('../GOTO_Data/Figures/Figure_4B-cg02649849.png')
print(plot)
dev.off()
## png
## 2
ID list
ID_list <- targets %>% mutate(
ID = paste0(IOP2_ID, '_', timepoint)
) %>% dplyr::select(ID)
Load RNAseq data
load('../GOTO_Data/DEGs/expData_adipose.RData')
goto_exp <- expData %>%
dplyr::select(IOP2_ID, timepoint, final_plate,
starts_with('ENS')) %>%
mutate(
ID = str_c(IOP2_ID, '_', timepoint)
)
exp_df <- left_join(ID_list, goto_exp, by = 'ID')
dim(exp_df)
## [1] 178 41655
Gene symbols
ens2gene <- cinaR::grch37
name <- "ENSG00000064218"
df <- exp_df %>% filter(!is.na(timepoint))
df$IOP2_ID <- as.factor(df$IOP2_ID)
df$timepoint <- as.factor(df$timepoint)
Dot plot
plot <-df %>%
ggplot(aes(x = timepoint, y = ENSG00000064218,
fill = timepoint)) +
geom_line(aes(group = IOP2_ID), color = '#DDDDDD', linewidth=1) +
geom_dotplot(binaxis = "y", binwidth = 0.5,
stackdir = "center",
dotsize = 0.8, stroke = 2) +
scale_fill_manual(values=c('#EE6677', '#66CCEE'))+
ggtitle('') +
ylab('') + xlab('') +
scale_x_discrete(expand = c(0,0)) +
theme(
panel.background = element_rect(fill = 'white',
color = 'black'),
legend.position = 'none',
plot.title =element_text(hjust=0.5),
text=element_text(size=16))
print(plot)
Save
png('../GOTO_Data/Figures/Figure_4B-DMRT3.png')
print(plot)
dev.off()
## png
## 2
plot <- targets %>%
ggplot(aes(x = timepoint, y = perc_body_fat,
fill = timepoint)) +
geom_line(aes(group = IOP2_ID), color = '#DDDDDD', linewidth=1) +
geom_dotplot(binaxis = "y", binwidth = 2,
stackdir = "center",
dotsize = 0.8, stroke = 2) +
scale_fill_manual(values=c('#EE6677', '#66CCEE'))+
ggtitle('') +
ylab('') + xlab('') +
scale_x_discrete(expand = c(0,0)) +
theme(
panel.background = element_rect(fill = 'white',
color = 'black'),
legend.position = 'none',
plot.title =element_text(hjust=0.5),
text=element_text(size=16))
print(plot)
Save
png('../GOTO_Data/Figures/Figure_4B-bodyfat.png')
print(plot)
dev.off()
## png
## 2
sessionInfo()
## R version 4.2.2 (2022-10-31)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Rocky Linux 8.10 (Green Obsidian)
##
## Matrix products: default
## BLAS/LAPACK: /usr/lib64/libopenblas-r0.3.15.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel grid stats4 stats graphics grDevices utils
## [8] datasets methods base
##
## other attached packages:
## [1] DNAmArray_2.0.0
## [2] pls_2.8-2
## [3] FDb.InfiniumMethylation.hg19_2.2.0
## [4] org.Hs.eg.db_3.14.0
## [5] TxDb.Hsapiens.UCSC.hg19.knownGene_3.2.2
## [6] minfi_1.40.0
## [7] bumphunter_1.36.0
## [8] locfit_1.5-9.8
## [9] iterators_1.0.14
## [10] foreach_1.5.2
## [11] Biostrings_2.62.0
## [12] XVector_0.34.0
## [13] ggpubr_0.4.0
## [14] ggrepel_0.9.1
## [15] circlize_0.4.15
## [16] ComplexHeatmap_2.14.0
## [17] RColorBrewer_1.1-3
## [18] pheatmap_1.0.12
## [19] lmerTest_3.1-3
## [20] lme4_1.1-30
## [21] Matrix_1.5-4.1
## [22] SummarizedExperiment_1.24.0
## [23] MatrixGenerics_1.10.0
## [24] matrixStats_1.0.0
## [25] AnnotationHub_3.2.2
## [26] BiocFileCache_2.2.1
## [27] dbplyr_2.2.1
## [28] GenomicFeatures_1.46.5
## [29] AnnotationDbi_1.56.2
## [30] Biobase_2.58.0
## [31] GenomicRanges_1.46.1
## [32] GenomeInfoDb_1.34.9
## [33] IRanges_2.32.0
## [34] S4Vectors_0.36.2
## [35] BiocGenerics_0.44.0
## [36] forcats_0.5.2
## [37] stringr_1.5.0
## [38] dplyr_1.1.3
## [39] purrr_0.3.4
## [40] readr_2.1.2
## [41] tidyr_1.2.1
## [42] tibble_3.2.1
## [43] ggplot2_3.4.3
## [44] tidyverse_1.3.2
## [45] cinaR_0.2.3
## [46] edgeR_3.40.2
## [47] limma_3.54.2
## [48] rmarkdown_2.16
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.3 tidyselect_1.2.0
## [3] RSQLite_2.2.17 BiocParallel_1.32.6
## [5] munsell_0.5.0 preprocessCore_1.60.2
## [7] codetools_0.2-19 withr_2.5.0
## [9] colorspace_2.1-0 filelock_1.0.2
## [11] highr_0.10 knitr_1.43
## [13] rstudioapi_0.14 ggsignif_0.6.3
## [15] labeling_0.4.2 GenomeInfoDbData_1.2.9
## [17] farver_2.1.1 bit64_4.0.5
## [19] rhdf5_2.42.1 vctrs_0.6.3
## [21] generics_0.1.3 xfun_0.39
## [23] timechange_0.2.0 R6_2.5.1
## [25] doParallel_1.0.17 illuminaio_0.40.0
## [27] clue_0.3-64 reshape_0.8.9
## [29] rhdf5filters_1.10.1 bitops_1.0-7
## [31] cachem_1.0.8 DelayedArray_0.24.0
## [33] assertthat_0.2.1 promises_1.2.0.1
## [35] BiocIO_1.8.0 scales_1.2.1
## [37] vroom_1.5.7 googlesheets4_1.0.1
## [39] gtable_0.3.3 Cairo_1.6-0
## [41] rlang_1.1.1 genefilter_1.76.0
## [43] GlobalOptions_0.1.2 splines_4.2.2
## [45] rtracklayer_1.54.0 rstatix_0.7.0
## [47] GEOquery_2.62.2 gargle_1.5.0
## [49] htm2txt_2.2.2 broom_1.0.1
## [51] BiocManager_1.30.21 yaml_2.3.7
## [53] reshape2_1.4.4 abind_1.4-5
## [55] modelr_0.1.9 backports_1.4.1
## [57] httpuv_1.6.11 tools_4.2.2
## [59] nor1mix_1.3-0 ellipsis_0.3.2
## [61] jquerylib_0.1.4 siggenes_1.68.0
## [63] Rcpp_1.0.10 plyr_1.8.8
## [65] sparseMatrixStats_1.10.0 progress_1.2.2
## [67] zlibbioc_1.44.0 RCurl_1.98-1.12
## [69] prettyunits_1.1.1 openssl_2.0.6
## [71] GetoptLong_1.0.5 haven_2.5.1
## [73] cluster_2.1.4 fs_1.6.2
## [75] magrittr_2.0.3 data.table_1.14.8
## [77] magick_2.7.4 reprex_2.0.2
## [79] googledrive_2.0.0 hms_1.1.2
## [81] mime_0.12 evaluate_0.21
## [83] xtable_1.8-4 XML_3.99-0.14
## [85] mclust_6.0.0 readxl_1.4.1
## [87] shape_1.4.6 compiler_4.2.2
## [89] biomaRt_2.50.3 crayon_1.5.2
## [91] minqa_1.2.5 htmltools_0.5.5
## [93] later_1.3.1 tzdb_0.4.0
## [95] lubridate_1.9.2 DBI_1.1.3
## [97] MASS_7.3-60 rappdirs_0.3.3
## [99] boot_1.3-28.1 car_3.1-0
## [101] cli_3.6.1 quadprog_1.5-8
## [103] pkgconfig_2.0.3 GenomicAlignments_1.30.0
## [105] numDeriv_2016.8-1.1 xml2_1.3.4
## [107] annotate_1.72.0 bslib_0.5.0
## [109] rngtools_1.5.2 multtest_2.50.0
## [111] beanplot_1.3.1 rvest_1.0.3
## [113] scrime_1.3.5 doRNG_1.8.6
## [115] digest_0.6.31 base64_2.0.1
## [117] cellranger_1.1.0 DelayedMatrixStats_1.16.0
## [119] restfulr_0.0.15 curl_5.0.1
## [121] shiny_1.7.2 Rsamtools_2.10.0
## [123] rjson_0.2.21 nloptr_2.0.3
## [125] lifecycle_1.0.3 nlme_3.1-162
## [127] jsonlite_1.8.5 Rhdf5lib_1.20.0
## [129] carData_3.0-5 askpass_1.1
## [131] fansi_1.0.4 pillar_1.9.0
## [133] lattice_0.21-8 survival_3.5-5
## [135] KEGGREST_1.34.0 fastmap_1.1.1
## [137] httr_1.4.6 interactiveDisplayBase_1.32.0
## [139] glue_1.6.2 png_0.1-8
## [141] BiocVersion_3.16.0 bit_4.0.5
## [143] stringi_1.7.12 sass_0.4.6
## [145] HDF5Array_1.22.1 blob_1.2.4
## [147] memoise_2.0.1
Clear
rm(list=ls())