This script loads in predicted and measured blood-type proportions and investigates how they changed during the intervention.
Load required packages
library(tidyverse)
library(SummarizedExperiment)
library(lmerTest)
library(lme4)
Load data
load("../GOTO_Data/GOTO_targets-filtered.Rdata")
load("../GOTO_Data/GOTO_methData-filtered.Rdata")
Cell names
cell_names <- c("cc_blood_meas_eos", "cc_blood_meas_baso", "cc_blood_meas_neut",
"cc_blood_meas_lymph", "cc_blood_meas_mono", "cc_blood_music_BNK",
"cc_blood_music_CD4T", "cc_blood_music_CD8T", "cc_blood_music_claM",
"cc_blood_music_CLP", "cc_blood_music_cMOP", "cc_blood_music_CMP",
"cc_blood_music_ery", "cc_blood_music_GMP", "cc_blood_music_hMDP",
"cc_blood_music_HSC", "cc_blood_music_immB", "cc_blood_music_interM",
"cc_blood_music_kineNK", "cc_blood_music_LMPP", "cc_blood_music_matureN",
"cc_blood_music_memB", "cc_blood_music_MEP", "cc_blood_music_metaN",
"cc_blood_music_MLP", "cc_blood_music_MPP", "cc_blood_music_myeN",
"cc_blood_music_naiB", "cc_blood_music_NKP", "cc_blood_music_nonM",
"cc_blood_music_plasma", "cc_blood_music_preB", "cc_blood_music_preM",
"cc_blood_music_proB", "cc_blood_music_proN", "cc_blood_music_regB",
"cc_blood_music_toxiNK", "cc_blood_idol_CD8T", "cc_blood_idol_CD4T",
"cc_blood_idol_NK", "cc_blood_idol_Bcell", "cc_blood_idol_Mono",
"cc_blood_idol_Neu", "cc_blood_ext_Bas", "cc_blood_ext_Bmem",
"cc_blood_ext_Bnv", "cc_blood_ext_CD4mem", "cc_blood_ext_CD4nv",
"cc_blood_ext_CD8mem", "cc_blood_ext_CD8nv", "cc_blood_ext_Eos",
"cc_blood_ext_Mono", "cc_blood_ext_Neu", "cc_blood_ext_NK",
"cc_blood_ext_Treg")
Keep only blood samples
targets <- targets %>% filter(tissue == 'fasted blood')
Run analysis
for(i in cell_names){
tryCatch({
base_mean <- mean((targets %>%
filter(timepoint == 0) %>% dplyr::select(all_of(i)))[,1], na.rm=TRUE)
base_sd <- sd((targets %>%
filter(timepoint == 0) %>% dplyr::select(all_of(i)))[,1], na.rm=TRUE)
fit <- lmer(get(i) ~ timepoint + age + sex + smoke + (1|IOP2_ID),
data = targets)
df <- data.frame(
cell = i,
mean = base_mean,
sd = base_sd,
es = summary(fit)$coefficients[2,1],
se = summary(fit)$coefficients[2,2],
p = summary(fit)$coefficients[2,5]
)}, error=function(e){})
if(i == cell_names[1]){
out <- df
} else {
out <- rbind(out, df)
}
}
Adjust p-values
out$padj <- p.adjust(out$p, method = 'fdr')
Inspect
out
## cell mean sd es se
## 1 cc_blood_meas_eos 3.572448980 2.50227239 -0.1282034359 0.144826480
## 2 cc_blood_meas_baso 0.529489796 0.30123927 0.0169492631 0.030226636
## 3 cc_blood_meas_neut 50.365306122 8.49938678 0.2145103307 0.665862174
## 4 cc_blood_meas_lymph 36.887857143 8.35822063 -0.0668066982 0.563027029
## 5 cc_blood_meas_mono 8.619795918 2.15162129 -0.0335627730 0.209167762
## 6 cc_blood_meas_mono 8.619795918 2.15162129 -0.0335627730 0.209167762
## 7 cc_blood_music_CD4T 0.292476827 1.99745302 1.0118318147 0.679145661
## 8 cc_blood_music_CD8T 18.804149558 8.74870228 -1.6043614568 0.948285622
## 9 cc_blood_music_claM 1.170128364 2.03774218 0.5045951652 0.329181020
## 10 cc_blood_music_claM 1.170128364 2.03774218 0.5045951652 0.329181020
## 11 cc_blood_music_cMOP 0.090823317 0.46365898 -0.0029254592 0.043178141
## 12 cc_blood_music_cMOP 0.090823317 0.46365898 -0.0029254592 0.043178141
## 13 cc_blood_music_ery 0.200073338 0.27429065 -0.0938340458 0.026784198
## 14 cc_blood_music_ery 0.200073338 0.27429065 -0.0938340458 0.026784198
## 15 cc_blood_music_hMDP 0.017040937 0.09766992 0.0119290478 0.023285655
## 16 cc_blood_music_hMDP 0.017040937 0.09766992 0.0119290478 0.023285655
## 17 cc_blood_music_immB 1.119572196 2.09667077 0.0305941866 0.115449371
## 18 cc_blood_music_interM 20.866926500 5.36707734 -0.5345511502 0.675672994
## 19 cc_blood_music_interM 20.866926500 5.36707734 -0.5345511502 0.675672994
## 20 cc_blood_music_interM 20.866926500 5.36707734 -0.5345511502 0.675672994
## 21 cc_blood_music_matureN 28.886669159 5.34590121 0.8121323638 0.606668338
## 22 cc_blood_music_matureN 28.886669159 5.34590121 0.8121323638 0.606668338
## 23 cc_blood_music_matureN 28.886669159 5.34590121 0.8121323638 0.606668338
## 24 cc_blood_music_metaN 24.484878804 4.91824037 0.0654190486 0.637049121
## 25 cc_blood_music_metaN 24.484878804 4.91824037 0.0654190486 0.637049121
## 26 cc_blood_music_metaN 24.484878804 4.91824037 0.0654190486 0.637049121
## 27 cc_blood_music_myeN 0.078717869 0.70846082 -0.0495652888 0.085420567
## 28 cc_blood_music_naiB 0.114359307 1.02303436 -0.0154298285 0.066192047
## 29 cc_blood_music_naiB 0.114359307 1.02303436 -0.0154298285 0.066192047
## 30 cc_blood_music_nonM 0.697287114 1.41953919 -0.2014784295 0.230293267
## 31 cc_blood_music_plasma 0.031678143 0.05776909 -0.0007911789 0.007833615
## 32 cc_blood_music_plasma 0.031678143 0.05776909 -0.0007911789 0.007833615
## 33 cc_blood_music_preM 3.133885498 1.60106385 0.0442845394 0.208315596
## 34 cc_blood_music_preM 3.133885498 1.60106385 0.0442845394 0.208315596
## 35 cc_blood_music_preM 3.133885498 1.60106385 0.0442845394 0.208315596
## 36 cc_blood_music_regB 0.001304611 0.01174150 -0.0010089701 0.001252167
## 37 cc_blood_music_toxiNK 0.010028458 0.02624329 0.0014794001 0.003205640
## 38 cc_blood_idol_CD8T 12.983266948 3.94190360 0.0374379301 0.209628651
## 39 cc_blood_idol_CD4T 15.735304887 5.71888076 0.0725776725 0.338227206
## 40 cc_blood_idol_NK 5.878544695 3.09030393 -0.2020364241 0.175332805
## 41 cc_blood_idol_Bcell 6.524220170 4.56760578 0.1303572515 0.128211288
## 42 cc_blood_idol_Mono 10.125862363 2.36821852 0.0569883077 0.193296831
## 43 cc_blood_idol_Neu 51.843173963 7.65722660 -0.1173757086 0.673422800
## 44 cc_blood_ext_Bas 0.678673469 0.62766732 0.0988558207 0.054777588
## 45 cc_blood_ext_Bmem 2.377653061 4.90761008 -0.0113838178 0.085383986
## 46 cc_blood_ext_Bnv 3.596938776 1.72662200 0.0958007570 0.081145307
## 47 cc_blood_ext_CD4mem 10.424183673 3.70244219 0.3126068658 0.268955124
## 48 cc_blood_ext_CD4nv 6.809489796 4.42993411 -0.1215789234 0.184747029
## 49 cc_blood_ext_CD8mem 7.269387755 4.78451687 -0.1022345282 0.215023779
## 50 cc_blood_ext_CD8nv 0.762755102 1.27056914 0.0829997276 0.070706112
## 51 cc_blood_ext_Eos 2.370510204 2.45403071 -0.0520323767 0.197745612
## 52 cc_blood_ext_Mono 7.353163265 2.01338695 0.0225512486 0.165156975
## 53 cc_blood_ext_Neu 47.721734694 8.77793614 -0.1896909439 0.747842700
## 54 cc_blood_ext_NK 5.590612245 2.43424571 -0.1422775222 0.154570089
## 55 cc_blood_ext_Treg 1.268265306 0.98248478 0.0268097219 0.073869179
## p padj
## 1 0.3781977983 0.94614913
## 2 0.5762625221 0.94614913
## 3 0.7480239208 0.94614913
## 4 0.9057904041 0.94614913
## 5 0.8728519972 0.94614913
## 6 0.8728519972 0.94614913
## 7 0.1382801416 0.94614913
## 8 0.0945378736 0.94614913
## 9 0.1292078445 0.94614913
## 10 0.1292078445 0.94614913
## 11 0.9461491339 0.94614913
## 12 0.9461491339 0.94614913
## 13 0.0007512441 0.02065921
## 14 0.0007512441 0.02065921
## 15 0.6091713988 0.94614913
## 16 0.6091713988 0.94614913
## 17 0.7916836185 0.94614913
## 18 0.4311935908 0.94614913
## 19 0.4311935908 0.94614913
## 20 0.4311935908 0.94614913
## 21 0.1844859109 0.94614913
## 22 0.1844859109 0.94614913
## 23 0.1844859109 0.94614913
## 24 0.9184653839 0.94614913
## 25 0.9184653839 0.94614913
## 26 0.9184653839 0.94614913
## 27 0.5625831908 0.94614913
## 28 0.8162691410 0.94614913
## 29 0.8162691410 0.94614913
## 30 0.3842150244 0.94614913
## 31 0.9198203376 0.94614913
## 32 0.9198203376 0.94614913
## 33 0.8322026063 0.94614913
## 34 0.8322026063 0.94614913
## 35 0.8322026063 0.94614913
## 36 0.4215968162 0.94614913
## 37 0.6459194618 0.94614913
## 38 0.8586234361 0.94614913
## 39 0.8305372041 0.94614913
## 40 0.2519791825 0.94614913
## 41 0.3116647888 0.94614913
## 42 0.7687550158 0.94614913
## 43 0.8619926537 0.94614913
## 44 0.0742044279 0.94614913
## 45 0.8941743062 0.94614913
## 46 0.2405860914 0.94614913
## 47 0.2479363533 0.94614913
## 48 0.5120012216 0.94614913
## 49 0.6354994641 0.94614913
## 50 0.2432740654 0.94614913
## 51 0.7930060229 0.94614913
## 52 0.8916720306 0.94614913
## 53 0.8002983680 0.94614913
## 54 0.3595779827 0.94614913
## 55 0.7174370050 0.94614913
Save MuSiC cells for figure
forest_df <- out %>%
filter(cell %in% c("cc_blood_music_BNK",
"cc_blood_music_CD4T", "cc_blood_music_CD8T", "cc_blood_music_claM",
"cc_blood_music_CLP", "cc_blood_music_cMOP", "cc_blood_music_CMP",
"cc_blood_music_ery", "cc_blood_music_GMP", "cc_blood_music_hMDP",
"cc_blood_music_HSC", "cc_blood_music_immB", "cc_blood_music_interM",
"cc_blood_music_kineNK", "cc_blood_music_LMPP", "cc_blood_music_matureN",
"cc_blood_music_memB", "cc_blood_music_MEP", "cc_blood_music_metaN",
"cc_blood_music_MLP", "cc_blood_music_MPP", "cc_blood_music_myeN",
"cc_blood_music_naiB", "cc_blood_music_NKP", "cc_blood_music_nonM",
"cc_blood_music_plasma", "cc_blood_music_preB", "cc_blood_music_preM",
"cc_blood_music_proB", "cc_blood_music_proN", "cc_blood_music_regB",
"cc_blood_music_toxiNK"))
Print plot
plot <- forest_df %>%
mutate(lower = es - (1.96 * se),
upper = es + (1.96 * se)) %>%
filter(mean > 0.05) %>%
ggplot(aes(x = es,
y = reorder(cell, mean),
xmin = lower,
xmax = upper)) +
geom_vline(xintercept=0, linewidth=1,
color='grey60', linetype='dashed') +
geom_errorbar(width=0.5,
linewidth=1,
position=position_dodge(width=0.9)) +
geom_point(size=5,
shape=21,
stroke=1.2,
position=position_dodge(width=0.9),
fill = "#D62839") +
xlab('Intervention effect') + ylab('') + xlim(c(-7,7)) +
theme(axis.text = element_text(size=14, color = '#373334'),
axis.title = element_text(size=16, hjust=0.5,
color = '#373334'),
text=element_text(size=14),
panel.background = element_rect(fill = 'white',
color='#373334'),
panel.grid.major = element_line(color = 'grey95'),
panel.grid.minor = element_line(color = 'grey95'),
plot.background = element_rect(fill = 'white'),
axis.ticks.x = element_line(size=1))
## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
print(plot)
Save plot
png(file="../GOTO_Data/Figures/Figure_5A.png")
print(plot)
dev.off()
## png
## 2
Save data
write_csv(out, file = '../GOTO_Data/Tables/ST16.csv')
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 stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] lmerTest_3.1-3
## [2] lme4_1.1-30
## [3] Matrix_1.5-4.1
## [4] IlluminaHumanMethylation450kmanifest_0.4.0
## [5] IlluminaHumanMethylationEPICmanifest_0.3.0
## [6] FlowSorted.BloodExtended.EPIC_1.1.1
## [7] FlowSorted.Blood.EPIC_1.12.1
## [8] IlluminaHumanMethylationEPICanno.ilm10b4.hg19_0.6.0
## [9] nlme_3.1-162
## [10] quadprog_1.5-8
## [11] genefilter_1.76.0
## [12] ExperimentHub_2.2.1
## [13] AnnotationHub_3.2.2
## [14] BiocFileCache_2.2.1
## [15] dbplyr_2.2.1
## [16] MuSiC_0.2.0
## [17] nnls_1.4
## [18] lubridate_1.9.2
## [19] DNAmArray_2.0.0
## [20] pls_2.8-2
## [21] FDb.InfiniumMethylation.hg19_2.2.0
## [22] org.Hs.eg.db_3.14.0
## [23] TxDb.Hsapiens.UCSC.hg19.knownGene_3.2.2
## [24] GenomicFeatures_1.46.5
## [25] AnnotationDbi_1.56.2
## [26] minfi_1.40.0
## [27] bumphunter_1.36.0
## [28] locfit_1.5-9.8
## [29] iterators_1.0.14
## [30] foreach_1.5.2
## [31] Biostrings_2.62.0
## [32] XVector_0.34.0
## [33] SummarizedExperiment_1.24.0
## [34] Biobase_2.58.0
## [35] MatrixGenerics_1.10.0
## [36] matrixStats_1.0.0
## [37] GenomicRanges_1.46.1
## [38] GenomeInfoDb_1.34.9
## [39] IRanges_2.32.0
## [40] S4Vectors_0.36.2
## [41] BiocGenerics_0.44.0
## [42] forcats_0.5.2
## [43] stringr_1.5.0
## [44] dplyr_1.1.3
## [45] purrr_0.3.4
## [46] readr_2.1.2
## [47] tidyr_1.2.1
## [48] tibble_3.2.1
## [49] ggplot2_3.4.3
## [50] tidyverse_1.3.2
## [51] cli_3.6.1
## [52] htmltools_0.5.5
## [53] rlang_1.1.1
## [54] rmarkdown_2.16
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.3 tidyselect_1.2.0
## [3] RSQLite_2.2.17 grid_4.2.2
## [5] BiocParallel_1.32.6 cinaR_0.2.3
## [7] munsell_0.5.0 codetools_0.2-19
## [9] preprocessCore_1.60.2 withr_2.5.0
## [11] colorspace_2.1-0 filelock_1.0.2
## [13] highr_0.10 knitr_1.43
## [15] rstudioapi_0.14 labeling_0.4.2
## [17] GenomeInfoDbData_1.2.9 MCMCpack_1.6-3
## [19] farver_2.1.1 bit64_4.0.5
## [21] rhdf5_2.42.1 coda_0.19-4
## [23] vctrs_0.6.3 generics_0.1.3
## [25] xfun_0.39 timechange_0.2.0
## [27] R6_2.5.1 illuminaio_0.40.0
## [29] bitops_1.0-7 rhdf5filters_1.10.1
## [31] cachem_1.0.8 reshape_0.8.9
## [33] DelayedArray_0.24.0 assertthat_0.2.1
## [35] vroom_1.5.7 promises_1.2.0.1
## [37] BiocIO_1.8.0 scales_1.2.1
## [39] googlesheets4_1.0.1 gtable_0.3.3
## [41] mcmc_0.9-7 MatrixModels_0.5-1
## [43] splines_4.2.2 rtracklayer_1.54.0
## [45] gargle_1.5.0 GEOquery_2.62.2
## [47] htm2txt_2.2.2 broom_1.0.1
## [49] BiocManager_1.30.21 yaml_2.3.7
## [51] reshape2_1.4.4 modelr_0.1.9
## [53] backports_1.4.1 httpuv_1.6.11
## [55] tools_4.2.2 nor1mix_1.3-0
## [57] ellipsis_0.3.2 jquerylib_0.1.4
## [59] RColorBrewer_1.1-3 siggenes_1.68.0
## [61] Rcpp_1.0.10 plyr_1.8.8
## [63] sparseMatrixStats_1.10.0 progress_1.2.2
## [65] zlibbioc_1.44.0 RCurl_1.98-1.12
## [67] prettyunits_1.1.1 openssl_2.0.6
## [69] haven_2.5.1 fs_1.6.2
## [71] magrittr_2.0.3 data.table_1.14.8
## [73] SparseM_1.81 reprex_2.0.2
## [75] googledrive_2.0.0 mime_0.12
## [77] hms_1.1.2 evaluate_0.21
## [79] xtable_1.8-4 XML_3.99-0.14
## [81] mclust_6.0.0 readxl_1.4.1
## [83] compiler_4.2.2 biomaRt_2.50.3
## [85] crayon_1.5.2 minqa_1.2.5
## [87] later_1.3.1 tzdb_0.4.0
## [89] DBI_1.1.3 MASS_7.3-60
## [91] rappdirs_0.3.3 boot_1.3-28.1
## [93] pkgconfig_2.0.3 GenomicAlignments_1.30.0
## [95] numDeriv_2016.8-1.1 xml2_1.3.4
## [97] annotate_1.72.0 bslib_0.5.0
## [99] rngtools_1.5.2 multtest_2.50.0
## [101] beanplot_1.3.1 rvest_1.0.3
## [103] doRNG_1.8.6 scrime_1.3.5
## [105] digest_0.6.31 base64_2.0.1
## [107] cellranger_1.1.0 edgeR_3.40.2
## [109] DelayedMatrixStats_1.16.0 restfulr_0.0.15
## [111] curl_5.0.1 shiny_1.7.2
## [113] Rsamtools_2.10.0 quantreg_5.94
## [115] nloptr_2.0.3 rjson_0.2.21
## [117] lifecycle_1.0.3 jsonlite_1.8.5
## [119] Rhdf5lib_1.20.0 askpass_1.1
## [121] limma_3.54.2 fansi_1.0.4
## [123] pillar_1.9.0 lattice_0.21-8
## [125] KEGGREST_1.34.0 fastmap_1.1.1
## [127] httr_1.4.6 survival_3.5-5
## [129] interactiveDisplayBase_1.32.0 glue_1.6.2
## [131] png_0.1-8 BiocVersion_3.16.0
## [133] bit_4.0.5 stringi_1.7.12
## [135] sass_0.4.6 HDF5Array_1.22.1
## [137] blob_1.2.4 memoise_2.0.1
Clear
rm(list=ls())