This script loads in predicted and measured blood-type proportions and investigates how they changed during the intervention.


Setup

Load required packages

library(tidyverse)
library(SummarizedExperiment)
library(lmerTest)
library(lme4)

Add to data

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')

Cell changes

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

Forest plot

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

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')

Session Info

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())