Last updated: 2021-09-24

Checks: 7 0

Knit directory: femNATCD_MethSeq/

This reproducible R Markdown analysis was created with workflowr (version 1.6.2). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.


Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.

Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.

The command set.seed(20210128) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.

Great job! Recording the operating system, R version, and package versions is critical for reproducibility.

Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.

Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.

Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.

The results in this page were generated with repository version 56bb68f. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.

Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:


Ignored files:
    Ignored:    .Rhistory
    Ignored:    .Rproj.user/
    Ignored:    analysis/.Rhistory
    Ignored:    code/.Rhistory
    Ignored:    data/Epicounts.csv
    Ignored:    data/Epimeta.csv
    Ignored:    data/Epitpm.csv
    Ignored:    data/KangUnivers.txt
    Ignored:    data/Kang_DataPreprocessing.RData
    Ignored:    data/Kang_dataset_genesMod_version2.txt
    Ignored:    data/PatMeta.csv
    Ignored:    data/ProcessedData.RData
    Ignored:    data/RTrawdata/
    Ignored:    data/SNPCommonFilt.csv
    Ignored:    data/femNAT_PC20.txt
    Ignored:    output/BrainMod_Enrichemnt.pdf
    Ignored:    output/Brain_Module_Heatmap.pdf
    Ignored:    output/DMR_Results.csv
    Ignored:    output/GOres.xlsx
    Ignored:    output/LME_GOplot.pdf
    Ignored:    output/LME_Results.csv
    Ignored:    output/LME_Results_Sig.csv
    Ignored:    output/LME_tophit.svg
    Ignored:    output/ProcessedData.RData
    Ignored:    output/RNAvsMETplots.pdf
    Ignored:    output/Regions_GOplot.pdf
    Ignored:    output/ResultsgroupComp.txt
    Ignored:    output/SEM_summary_groupEpi_M15.txt
    Ignored:    output/SEM_summary_groupEpi_M2.txt
    Ignored:    output/SEM_summary_groupEpi_M_all.txt
    Ignored:    output/SEM_summary_groupEpi_TopHit.txt
    Ignored:    output/SEM_summary_groupEpi_all.txt
    Ignored:    output/SEMplot_Epi_M15.html
    Ignored:    output/SEMplot_Epi_M15.png
    Ignored:    output/SEMplot_Epi_M15_files/
    Ignored:    output/SEMplot_Epi_M2.html
    Ignored:    output/SEMplot_Epi_M2.png
    Ignored:    output/SEMplot_Epi_M2_files/
    Ignored:    output/SEMplot_Epi_M_all.html
    Ignored:    output/SEMplot_Epi_M_all.png
    Ignored:    output/SEMplot_Epi_M_all_files/
    Ignored:    output/SEMplot_Epi_TopHit.html
    Ignored:    output/SEMplot_Epi_TopHit.png
    Ignored:    output/SEMplot_Epi_TopHit_files/
    Ignored:    output/SEMplot_Epi_all.html
    Ignored:    output/SEMplot_Epi_all.png
    Ignored:    output/SEMplot_Epi_all_files/
    Ignored:    output/barplots.pdf
    Ignored:    output/circos_DMR_tags.svg
    Ignored:    output/circos_LME_tags.svg
    Ignored:    output/clinFact.RData
    Ignored:    output/dds_filt_analyzed.RData
    Ignored:    output/designh0.RData
    Ignored:    output/designh1.RData
    Ignored:    output/envFact.RData
    Ignored:    output/functional_Enrichemnt.pdf
    Ignored:    output/gostres.pdf
    Ignored:    output/modelFact.RData
    Ignored:    output/resdmr.RData
    Ignored:    output/resultsdmr_table.RData
    Ignored:    output/table1_filtered.Rmd
    Ignored:    output/table1_filtered.docx
    Ignored:    output/table1_unfiltered.Rmd
    Ignored:    output/table1_unfiltered.docx
    Ignored:    setup_built.R

Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.


These are the previous versions of the repository in which changes were made to the R Markdown (analysis/02_1_Main_LME.Rmd) and HTML (docs/02_1_Main_LME.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.

File Version Author Date Message
html 43333cc achiocch 2021-09-24 adds new build
html b497cc9 achiocch 2021-09-17 Build site.
html ccbc9e4 achiocch 2021-08-10 Build site.
html 70cd649 achiocch 2021-08-06 Build site.
html 2a53a87 achiocch 2021-08-06 Build site.
html e3a9ae3 achiocch 2021-08-04 Build site.
html f6bbdc0 achiocch 2021-08-04 Build site.
html 710904a achiocch 2021-08-03 Build site.
Rmd 16112c3 achiocch 2021-08-03 wflow_publish(c(“analysis/", "code/”, “docs/”))
html 16112c3 achiocch 2021-08-03 wflow_publish(c(“analysis/", "code/”, “docs/”))
html cde8384 achiocch 2021-08-03 Build site.
html d761be4 achiocch 2021-07-31 Build site.
html b452d2f achiocch 2021-07-30 wflow_publish(c(“analysis/", "code/”, "docs/*"))
Rmd 1a9f36f achiocch 2021-07-30 reviewed analysis
html 2734c4e achiocch 2021-05-08 Build site.
html a847823 achiocch 2021-05-08 wflow_publish(c(“analysis/", "code/”, "docs/*"))
html 9cc52f7 achiocch 2021-05-08 Build site.
html 158d0b4 achiocch 2021-05-08 wflow_publish(c(“analysis/", "code/”, "docs/*"))
html 0f262d1 achiocch 2021-05-07 Build site.
html 5167b90 achiocch 2021-05-07 wflow_publish(c(“analysis/", "code/”, "docs/*"))
html 05aac7f achiocch 2021-04-23 Build site.
html f5c5265 achiocch 2021-04-19 Build site.
html 17f1eec achiocch 2021-04-10 Build site.
Rmd 4dee231 achiocch 2021-04-10 wflow_publish(c(“analysis/", "code/”, "docs/*"))
html 91de221 achiocch 2021-04-05 Build site.
html 4ea1bba achiocch 2021-02-25 Build site.
Rmd 6c21638 achiocch 2021-02-25 wflow_publish(c(“analysis/", "code/”, "docs/*"), update = F)
html 6c21638 achiocch 2021-02-25 wflow_publish(c(“analysis/", "code/”, "docs/*"), update = F)

## load Data  ####
load(paste0(Home,"/output/ProcessedData.RData"))

define model

# test SVAs for significant difference between groups

svas = names(colData(dds_filt))[str_detect(names(colData(dds_filt)), "^V")]

testframe = as.data.frame(colData(dds_filt)) %>% dplyr::select(group, all_of(svas))

svas_sig = c()

for (s in svas){
  tmp = t.test(testframe[,s]~testframe[,"group"])$p.value
  if(tmp <= 0.05){
    svas_sig=c(svas_sig,s) 
  }
}

#add significant sva to model

sva = paste0(svas_sig, collapse = "+")

dds_filt$group <- relevel(dds_filt$group, ref = "CTRL")
tmpdes=gsub(" group", paste0(sva,"+group"), as.character(design(dds_filt))[2])
design(dds_filt) = formula(paste0("~",tmpdes))

designmatrix = model.matrix(design(dds_filt), colData(dds_filt))

model formula: Epi ~ 0 + Age + int_dis + medication + contraceptives + cigday_1 + V8 + group

analyze

if (reanalyze){
  dds_filt = DESeq(dds_filt)
  savefile=TRUE
} else 
{
  if(file.exists(paste0(Home,"/output/dds_filt_analyzed.RData"))){
    load(paste0(Home,"/output/dds_filt_analyzed.RData"))
    savefile=FALSE
    
  }else {
    print("no preprocessed file found data is reanalyzed")
    dds_filt = DESeq(dds_filt)
    savefile=TRUE
  }
}

res = results(dds_filt)
summary(res)

out of 216102 with nonzero total read count
adjusted p-value < 0.1
LFC > 0 (up)       : 0, 0%
LFC < 0 (down)     : 1, 0.00046%
outliers [1]       : 0, 0%
low counts [2]     : 0, 0%
(mean count < 4)
[1] see 'cooksCutoff' argument of ?results
[2] see 'independentFiltering' argument of ?results
res <- res %>% cbind(as.data.frame(rowRanges(dds_filt)))
res$ID = paste0("tag_",1:nrow(res), "_",res$gene)
rownames(res)=res$ID
res$CHR = as.numeric(res$seqnames)

restop = as.data.frame(res[order(res$pvalue, decreasing = F), ])
restop = restop[restop$pvalue<=10E-3,]

if(savefile){
  save(dds_filt, file=paste0(Home,"/output/dds_filt_analyzed.RData"))}

sumarize results linear model

topN = 10
DESeq2::plotDispEsts(dds_filt)

lp <- -log10(res$pvalue)
ord <- order(lp, decreasing = TRUE)[1:topN]

qqman::qq(res$pvalue, main="uncorrected p-value")

qqman::manhattan(res, chr="CHR", bp="start", p="pvalue", snp="ID", col = c("#366697", "#BB8300"))

ressum= res[ord,]
display_tab(ressum)
baseMean log2FoldChange lfcSE stat pvalue padj seqnames start end width strand raw_gene gene distance feature cpg tf_binding samples_per_target mean std baseMean.1 baseVar allZero dispGeneEst dispGeneIter dispFit dispersion dispIter dispOutlier dispMAP Age int_disno int_disyes medicationyes contraceptivesyes cigday_1 V8 groupCD SE_Age SE_int_disno SE_int_disyes SE_medicationyes SE_contraceptivesyes SE_cigday_1 SE_V8 SE_groupCD WaldStatistic_Age WaldStatistic_int_disno WaldStatistic_int_disyes WaldStatistic_medicationyes WaldStatistic_contraceptivesyes WaldStatistic_cigday_1 WaldStatistic_V8 WaldStatistic_groupCD WaldPvalue_Age WaldPvalue_int_disno WaldPvalue_int_disyes WaldPvalue_medicationyes WaldPvalue_contraceptivesyes WaldPvalue_cigday_1 WaldPvalue_V8 WaldPvalue_groupCD betaConv betaIter deviance maxCooks ID CHR
tag_134778_MIR4500HG,SLITRK5 13.294 -1.128 0.201 -5.600 0 0.005 chr13 88324376 88324381 6
MIR4500HG,SLITRK5 MIR4500HG,SLITRK5 0 upstream cpg Pol2//GR 100 15.706 12.631 13.294 61.760 FALSE 0.177 11 0.171 0.176 11 FALSE 0.176 -0.045 4.741 5.031 0.167 -0.126 0.031 0.706 -1.128 0.054 0.854 0.884 0.214 0.192 0.017 0.879 0.201 -0.839 5.551 5.689 0.781 -0.660 1.763 0.803 -5.600 0.401 0.000 0.000 0.435 0.509 0.078 0.422 0 TRUE 6 645.105 NA tag_134778_MIR4500HG,SLITRK5 13
tag_204782_PDGFB 14.798 0.792 0.174 4.558 0 0.308 chr22 39629982 39629987 6
PDGFB PDGFB 0 intronic NoCpg 100 16.873 10.585 14.798 67.579 FALSE 0.131 11 0.154 0.135 11 FALSE 0.135 0.118 1.825 1.468 -0.032 -0.272 -0.009 0.172 0.792 0.048 0.767 0.792 0.184 0.169 0.015 0.810 0.174 2.448 2.380 1.854 -0.176 -1.610 -0.606 0.212 4.558 0.014 0.017 0.064 0.860 0.107 0.545 0.832 0 TRUE 6 650.103 NA tag_204782_PDGFB 22
tag_35466_MIR4792 5.316 -1.543 0.339 -4.552 0 0.308 chr3 24561792 24561797 6
MIR4792 MIR4792 0 downstream NoCpg Pol2//CCNT2//ZBTB7A_(SC-34508)//GABP//HDAC2_(SC-6296)//TAF1//Pol2-4H8//NRSF//GATA-2//ELF1_(SC-631)//GATA2_(CG2-96)//TAL1_(SC-12984) 86 5.882 5.642 5.316 24.526 FALSE 0.494 11 0.420 0.481 12 FALSE 0.481 -0.058 3.348 3.980 0.265 0.228 0.076 1.156 -1.543 0.090 1.417 1.466 0.352 0.311 0.028 1.430 0.339 -0.647 2.363 2.715 0.752 0.731 2.711 0.808 -4.552 0.517 0.018 0.007 0.452 0.465 0.007 0.419 0 TRUE 6 525.501 NA tag_35466_MIR4792 3
tag_120467_PVRL1 36.541 -0.468 0.103 -4.538 0 0.308 chr11 119544539 119544544 6
PVRL1 PVRL1 0 intronic NoCpg Pol2-4H8//NFKB//BATF//YY1_(C-20)//CTCF//TCF12//CHD2_(N-1250)//RXRA//Pol2//SMC3_(ab9263) 102 42.863 26.350 36.541 118.993 FALSE 0.038 7 0.065 0.042 11 FALSE 0.042 -0.030 5.753 5.685 0.303 0.117 0.016 -0.529 -0.468 0.028 0.444 0.459 0.109 0.099 0.009 0.461 0.103 -1.066 12.948 12.392 2.774 1.178 1.746 -1.148 -4.538 0.287 0.000 0.000 0.006 0.239 0.081 0.251 0 TRUE 4 725.887 NA tag_120467_PVRL1 11
tag_137745_KLHDC2 18.591 -0.686 0.154 -4.447 0 0.321 chr14 50234989 50234994 6
KLHDC2 KLHDC2 0 UTR5 cpg Pol2//TAF1//TBP//CCNT2//Sin3Ak-20//HMGN3//Pol2-4H8//HEY1//Oct-2 102 22.373 17.109 18.591 61.296 FALSE 0.093 8 0.123 0.098 11 FALSE 0.098 -0.028 4.755 4.974 0.334 -0.145 0.039 0.049 -0.686 0.042 0.659 0.682 0.162 0.147 0.013 0.679 0.154 -0.677 7.211 7.297 2.066 -0.987 2.979 0.072 -4.447 0.498 0.000 0.000 0.039 0.324 0.003 0.943 0 TRUE 4 667.003 NA tag_137745_KLHDC2 14
tag_197061_BMP7 9.825 1.043 0.235 4.442 0 0.321 chr20 55812238 55812243 6
BMP7 BMP7 0 intronic NoCpg 96 10.961 8.664 9.825 93.487 FALSE 0.272 12 0.229 0.264 12 FALSE 0.264 0.011 2.921 2.234 -0.190 -0.073 -0.031 0.077 1.043 0.065 1.034 1.069 0.250 0.229 0.021 1.102 0.235 0.165 2.825 2.091 -0.761 -0.318 -1.521 0.070 4.442 0.869 0.005 0.037 0.447 0.751 0.128 0.944 0 TRUE 6 614.616 NA tag_197061_BMP7 20
tag_68218_BEND3 8.151 -0.968 0.226 -4.280 0 0.479 chr6 107436782 107436787 6
BEND3 BEND3 0 upstream cpg E2F1//Pol2//HA-E2F1//AP-2gamma//Egr-1//E2F6_(H-50)//Pol2-4H8 98 9.745 7.996 8.151 24.042 FALSE 0.179 8 0.275 0.195 11 FALSE 0.195 -0.060 4.147 4.160 0.356 0.056 0.048 -1.633 -0.968 0.060 0.955 0.987 0.239 0.214 0.019 0.943 0.226 -0.985 4.342 4.217 1.494 0.262 2.480 -1.733 -4.280 0.324 0.000 0.000 0.135 0.793 0.013 0.083 0 TRUE 5 563.376 NA tag_68218_BEND3 6
tag_10516_C1orf52 6.161 -1.091 0.255 -4.278 0 0.479 chr1 85722202 85722207 6
C1orf52 C1orf52 0 intronic NoCpg 92 7.520 6.511 6.161 19.477 FALSE 0.210 10 0.363 0.239 9 FALSE 0.239 -0.088 4.098 4.510 0.592 -0.160 0.064 -0.449 -1.091 0.068 1.071 1.107 0.263 0.241 0.021 1.092 0.255 -1.299 3.827 4.075 2.253 -0.667 3.001 -0.411 -4.278 0.194 0.000 0.000 0.024 0.505 0.003 0.681 0 TRUE 6 529.701 NA tag_10516_C1orf52 1
tag_15433_METTL13 9.501 -0.858 0.202 -4.251 0 0.479 chr1 171751562 171751567 6
METTL13 METTL13 0 intronic NoCpg 100 11.255 9.282 9.501 30.942 FALSE 0.137 11 0.237 0.154 11 FALSE 0.154 -0.019 3.838 3.822 0.235 0.085 0.001 -0.378 -0.858 0.054 0.859 0.888 0.217 0.194 0.018 0.901 0.202 -0.346 4.471 4.307 1.084 0.436 0.046 -0.420 -4.251 0.730 0.000 0.000 0.278 0.663 0.964 0.674 0 TRUE 4 575.449 NA tag_15433_METTL13 1
tag_55541_ANKRD34B 7.256 1.048 0.247 4.240 0 0.479 chr5 79855772 79855777 6
ANKRD34B ANKRD34B 0 exonic NoCpg 92 8.422 6.364 7.256 28.876 FALSE 0.262 8 0.309 0.271 11 FALSE 0.271 0.154 0.071 0.263 -0.411 0.037 -0.076 -2.467 1.048 0.069 1.100 1.136 0.265 0.239 0.022 1.105 0.247 2.223 0.065 0.232 -1.551 0.153 -3.434 -2.233 4.240 0.026 0.948 0.817 0.121 0.878 0.001 0.026 0 TRUE 6 562.754 NA tag_55541_ANKRD34B 5
write.csv(res, file = paste0(Home,"/output/LME_Results.csv"), row.names = T)
write.csv2(restop, file = paste0(Home,"/output/LME_Results_Sig.csv"), row.names = T)
EnhancedVolcano(res,
                lab = rownames(res),
                x = 'log2FoldChange',
                y = 'pvalue',
                pCutoff = 10e-3,
                FCcutoff = 1)


sessionInfo()
R version 4.0.3 (2020-10-10)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19042)

Matrix products: default

locale:
[1] LC_COLLATE=German_Germany.1252  LC_CTYPE=German_Germany.1252   
[3] LC_MONETARY=German_Germany.1252 LC_NUMERIC=C                   
[5] LC_TIME=German_Germany.1252    

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets 
[8] methods   base     

other attached packages:
 [1] scales_1.1.1                compareGroups_4.4.6        
 [3] kableExtra_1.3.1            qqman_0.1.4                
 [5] EnhancedVolcano_1.8.0       ggrepel_0.9.1              
 [7] RCircos_1.2.1               doParallel_1.0.16          
 [9] iterators_1.0.13            foreach_1.5.1              
[11] DESeq2_1.30.0               SummarizedExperiment_1.20.0
[13] Biobase_2.50.0              MatrixGenerics_1.2.0       
[15] matrixStats_0.57.0          GenomicRanges_1.42.0       
[17] GenomeInfoDb_1.26.2         IRanges_2.24.1             
[19] S4Vectors_0.28.1            BiocGenerics_0.36.0        
[21] forcats_0.5.0               stringr_1.4.0              
[23] dplyr_1.0.2                 purrr_0.3.4                
[25] readr_1.4.0                 tidyr_1.1.2                
[27] tibble_3.0.4                ggplot2_3.3.3              
[29] tidyverse_1.3.0             RColorBrewer_1.1-2         
[31] knitr_1.30                  workflowr_1.6.2            

loaded via a namespace (and not attached):
  [1] uuid_0.1-4             readxl_1.3.1           backports_1.2.0       
  [4] systemfonts_0.3.2      splines_4.0.3          BiocParallel_1.24.1   
  [7] digest_0.6.27          htmltools_0.5.1.1      fansi_0.4.1           
 [10] magrittr_2.0.1         Rsolnp_1.16            memoise_2.0.0         
 [13] annotate_1.68.0        modelr_0.1.8           extrafont_0.17        
 [16] officer_0.3.16         extrafontdb_1.0        colorspace_2.0-0      
 [19] blob_1.2.1             rvest_0.3.6            haven_2.3.1           
 [22] xfun_0.20              crayon_1.3.4           RCurl_1.98-1.2        
 [25] jsonlite_1.7.2         genefilter_1.72.0      survival_3.2-7        
 [28] glue_1.4.2             gtable_0.3.0           zlibbioc_1.36.0       
 [31] XVector_0.30.0         webshot_0.5.2          DelayedArray_0.16.0   
 [34] proj4_1.0-10.1         Rttf2pt1_1.3.8         maps_3.3.0            
 [37] DBI_1.1.1              Rcpp_1.0.5             viridisLite_0.3.0     
 [40] xtable_1.8-4           bit_4.0.4              truncnorm_1.0-8       
 [43] httr_1.4.2             calibrate_1.7.7        ellipsis_0.3.1        
 [46] mice_3.12.0            farver_2.0.3           pkgconfig_2.0.3       
 [49] XML_3.99-0.5           dbplyr_2.0.0           locfit_1.5-9.4        
 [52] labeling_0.4.2         tidyselect_1.1.0       rlang_0.4.10          
 [55] later_1.1.0.1          AnnotationDbi_1.52.0   munsell_0.5.0         
 [58] cellranger_1.1.0       tools_4.0.3            cachem_1.0.1          
 [61] cli_2.2.0              generics_0.1.0         RSQLite_2.2.2         
 [64] broom_0.7.3            evaluate_0.14          fastmap_1.1.0         
 [67] yaml_2.2.1             bit64_4.0.5            fs_1.5.0              
 [70] zip_2.1.1              whisker_0.4            ash_1.0-15            
 [73] ggrastr_0.2.1          xml2_1.3.2             compiler_4.0.3        
 [76] rstudioapi_0.13        beeswarm_0.2.3         reprex_1.0.0          
 [79] geneplotter_1.68.0     stringi_1.5.3          HardyWeinberg_1.7.1   
 [82] highr_0.8              ps_1.5.0               ggalt_0.4.0           
 [85] gdtools_0.2.3          lattice_0.20-41        Matrix_1.2-18         
 [88] vctrs_0.3.6            pillar_1.4.7           lifecycle_0.2.0       
 [91] data.table_1.13.6      flextable_0.6.2        bitops_1.0-6          
 [94] httpuv_1.5.5           R6_2.5.0               promises_1.1.1        
 [97] KernSmooth_2.23-18     vipor_0.4.5            writexl_1.3.1         
[100] codetools_0.2-18       MASS_7.3-53            assertthat_0.2.1      
[103] chron_2.3-56           rprojroot_2.0.2        withr_2.4.1           
[106] GenomeInfoDbData_1.2.4 hms_1.0.0              grid_4.0.3            
[109] rmarkdown_2.6          git2r_0.28.0           base64enc_0.1-3       
[112] lubridate_1.7.9.2      ggbeeswarm_0.6.0