Interpretation of DESeq2 results with 3 factorial design with interaction terms.
2
0
Entering edit mode
s-rhussien • 0
@2a9652c8
Last seen 21 days ago
Egypt

Assalam Alykom,

I'm doing a 3 factorial DESeq2 analysis with interaction terms, I want to extract results to answer specific biological questions. Here is the code for the design:


dds<-DESeqDataSetFromMatrix(countData = CountData, colData = metaData, design = ~ age + gender + viralLoad + age:viralLoad + gender:viralLoad)
dds <- DESeq(dds)
resultsNames(dds)

# [1] "Intercept"               "age_old_vs_young"        "gender_F_vs_M"           "viralLoad_high_vs_N.A"  
# [5] "viralLoad_low_vs_N.A"    "viralLoad_medium_vs_N.A" "ageold.viralLoadhigh"    "ageold.viralLoadlow"    
# [9] "ageold.viralLoadmedium"  "genderF.viralLoadhigh"   "genderF.viralLoadlow"    "genderF.viralLoadmedium"

If I want to answer the following question: 1- For females, what are the DEGs in high viral load vs. control? Reading the example here, I drew the following design graph (reference levels are: age: young, viral load: N.A, gender: M): design graph I use the following coefficients: viralLoad_high_vs_N.A + genderF.viralLoadhigh. I use the following code to extract this data:

res <- results(dds, contrast=c(0,0,0,1,0,0,0,0,0,1,0,0), alpha=0.05)

But what if I want to answer: 2- For males, what are the DEGs in high viral load vs. control? How should I do this? as there is no interaction term such as: genderM.viralLoadhigh in the list of resultsNames(dds). Could I use something like: viralLoad_high_vs_N.A - gender_F_vs_M - genderF.viralLoadhigh ?

sessionInfo( )

R version 4.0.4 (2021-02-15)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19042)

Matrix products: default

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

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

other attached packages:
 [1] pheatmap_1.0.12             RColorBrewer_1.1-2          EnhancedVolcano_1.8.0      
 [4] ggrepel_0.9.1               forcats_0.5.1               stringr_1.4.0              
 [7] purrr_0.3.4                 readr_2.0.1                 tidyr_1.1.3                
[10] tibble_3.1.3                ggplot2_3.3.5               tidyverse_1.3.1            
[13] DESeq2_1.30.1               SummarizedExperiment_1.20.0 MatrixGenerics_1.2.1       
[16] matrixStats_0.60.0          dplyr_1.0.7                 GEOquery_2.58.0            
[19] Biobase_2.50.0              rtracklayer_1.49.5          GenomicRanges_1.42.0       
[22] GenomeInfoDb_1.26.7         Biostrings_2.58.0           XVector_0.30.0             
[25] IRanges_2.24.1              S4Vectors_0.28.1            AnnotationHub_2.22.1       
[28] BiocFileCache_1.14.0        dbplyr_2.1.1                BiocGenerics_0.36.1        

loaded via a namespace (and not attached):
  [1] ggbeeswarm_0.6.0              colorspace_2.0-2              ellipsis_0.3.2               
  [4] fs_1.5.0                      rstudioapi_0.13               bit64_4.0.5                  
  [7] interactiveDisplayBase_1.28.0 AnnotationDbi_1.52.0          fansi_0.5.0                  
 [10] lubridate_1.7.10              xml2_1.3.2                    splines_4.0.4                
 [13] extrafont_0.17                cachem_1.0.5                  geneplotter_1.68.0           
 [16] knitr_1.33                    jsonlite_1.7.2                Rsamtools_2.6.0              
 [19] Rttf2pt1_1.3.9                broom_0.7.9                   annotate_1.68.0              
 [22] shiny_1.6.0                   BiocManager_1.30.16           compiler_4.0.4               
 [25] httr_1.4.2                    backports_1.2.1               assertthat_0.2.1             
 [28] Matrix_1.3-4                  fastmap_1.1.0                 limma_3.46.0                 
 [31] cli_3.0.1                     later_1.3.0                   htmltools_0.5.1.1            
 [34] tools_4.0.4                   gtable_0.3.0                  glue_1.4.2                   
 [37] GenomeInfoDbData_1.2.4        maps_3.3.0                    rappdirs_0.3.3               
 [40] tinytex_0.33                  Rcpp_1.0.7                    cellranger_1.1.0             
 [43] vctrs_0.3.8                   ggalt_0.4.0                   extrafontdb_1.0              
 [46] xfun_0.23                     rvest_1.0.1                   mime_0.11                    
 [49] lifecycle_1.0.0               XML_3.99-0.6                  MASS_7.3-54                  
 [52] zlibbioc_1.36.0               scales_1.1.1                  vroom_1.5.4                  
 [55] hms_1.1.0                     promises_1.2.0.1              proj4_1.0-10.1               
 [58] yaml_2.2.1                    curl_4.3.2                    memoise_2.0.0                
 [61] ggrastr_0.2.3                 stringi_1.7.4                 RSQLite_2.2.7                
 [64] BiocVersion_3.12.0            genefilter_1.72.1             BiocParallel_1.24.1          
 [67] rlang_0.4.11                  pkgconfig_2.0.3               bitops_1.0-7                 
 [70] evaluate_0.14                 lattice_0.20-44               GenomicAlignments_1.26.0     
 [73] bit_4.0.4                     tidyselect_1.1.1              magrittr_2.0.1               
 [76] R6_2.5.1                      generics_0.1.0                DelayedArray_0.16.3          
 [79] DBI_1.1.1                     pillar_1.6.2                  haven_2.4.3                  
 [82] withr_2.4.2                   ash_1.0-15                    survival_3.2-11              
 [85] RCurl_1.98-1.3                modelr_0.1.8                  crayon_1.4.1                 
 [88] KernSmooth_2.23-20            utf8_1.2.2                    tzdb_0.1.2                   
 [91] rmarkdown_2.10                locfit_1.5-9.4                grid_4.0.4                   
 [94] readxl_1.3.1                  blob_1.2.2                    reprex_2.0.1                 
 [97] digest_0.6.27                 xtable_1.8-4                  httpuv_1.6.2                 
[100] munsell_0.5.0                 beeswarm_0.4.0                vipor_0.4.5
r DESeq DESeq2 • 256 views
ADD COMMENT
2
Entering edit mode
@mikelove
Last seen 22 hours ago
United States

Unfortunately, I only have time to answer software related questions on the support site. For statistical consultation about design and contrasts, I recommend meeting with a local statistician at your institute.

ADD COMMENT
0
Entering edit mode

Thanks very much

ADD REPLY
2
Entering edit mode
swbarnes2 ▴ 890
@swbarnes2-14086
Last seen 7 hours ago
San Diego

I think the simplest way to do your first question is to make a new column that joins sex and load, make your design sex_load + age, and then use contrasts to specify that you want to compare F_high to F_control.

See here. I think your questions fall in the category of "Can be done simpler without interactions"

http://bioconductor.org/packages/release/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#interactions

ADD COMMENT
0
Entering edit mode

You mean that if I have a factor with more than 2 levels, I can compare any 2 levels of that factor using contrast=c(factor_name, 1st_desired_level, 2nd_desired_level)?

ADD REPLY
0
Entering edit mode

Yes, that's just how the vignette says to use contrasts.

ADD REPLY
0
Entering edit mode

Thanks very much

ADD REPLY

Login before adding your answer.

Traffic: 323 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6