Validate blunted microglia immune response using FACS sorted microglia
Differential gene expression for design: Group + Trt + batch + Group:Trt
dds <- DESeqDataSetFromTximport(txi.rsem, colData = metadata, design = ~ group + trt + batch + group:trt)
samp <- colData(dds)$ID != "S16_MS6" & colData(dds)$ID != "S21_M41" # low read counts
# filter low expressed genes
keep <- rowSums(counts(dds)) >= 75
table(keep)
## keep
## FALSE TRUE
## 37036 18451
dds <- dds[keep,samp]
# add in gene symbols
genes <- gsub("\\..*","",rownames(dds))
# symbols <- mapIds(org.Mm.eg.db, keys=genes, column="SYMBOL", keytype="ENSEMBL", multiVals="first")
library(AnnotationHub)
ah <- AnnotationHub()
edb <- ah[["AH53222"]]
symbol <- mapIds(edb, keys=genes, column="SYMBOL", keytype="GENEID", multiVals="first")
mcols(dds) <- cbind(mcols(dds), symbol)
dds <- DESeq(dds)
LvSinCON <- (results(dds, name="trt_LPS_vs_SAL", alpha=0.05))
LvSinMIA <- (results(dds, list( c("trt_LPS_vs_SAL","groupMIA.trtLPS")), alpha=0.05))
LvSinCON$gene <- symbol
LvSinMIA$gene <- symbol
LPS Response in CON vs MIA
LPS.Response <- merge(as.data.frame(LvSinCON[,c(7,1,2,3,4,5,6)]),as.data.frame(LvSinMIA[,c(2,3,4,5,6)]), by="row.names", all.x=TRUE)
rownames(LPS.Response) <- LPS.Response$Row.names
colnames(LPS.Response) <- c("gene.id","gene.name","baseMean","log2FC.CON","lfcSE.CON","stat.CON","pval.CON","padj.CON","log2FC.MIA","lfcSE.MIA","stat.MIA","pval.MIA","padj.MIA")
#Assign the Group variable for coloring 4-way plot
for (i in c(1:nrow(LPS.Response))){
if (LPS.Response$padj.CON[i]<0.05 & (!is.na(LPS.Response$padj.CON[i]))){
LPS.Response$CONSig[i]='YES'
}
else {
LPS.Response$CONSig[i]='NO'
}
}
for (i in c(1:nrow(LPS.Response))){
if (LPS.Response$padj.MIA[i]<0.05 & (!is.na(LPS.Response$padj.MIA[i]))){
LPS.Response$MIASig[i]='YES'
}
else {
LPS.Response$MIASig[i]='NO'
}
}
for (i in c(1:nrow(LPS.Response))){
if (LPS.Response$CONSig[i]=='YES' & LPS.Response$MIASig[i]=='NO'){
LPS.Response$Group[i]=1 }
else if (LPS.Response$CONSig[i]=='NO' & LPS.Response$MIASig[i]=='YES'){
LPS.Response$Group[i]=2}
else if (LPS.Response$CONSig[i]=='YES' & LPS.Response$MIASig[i]=='YES'){
LPS.Response$Group[i]=3}
else {
LPS.Response$Group[i]=4}
}
LPS.Response.final<-subset(LPS.Response, LPS.Response$Group != 4)
ggplotRegression <- function (fit) {
require(ggplot2)
ggplot(fit$model, aes_string(x = names(fit$model)[2], y = names(fit$model)[1])) +
geom_point(color = "grey50") +
stat_smooth(method = "lm", col = "blue", size=2, se = TRUE) +
labs(subtitle = paste(" Slope =",signif(fit$coef[[2]], 3),
" P =",signif(summary(fit)$coef[2,4], 3),
"Adj R2 = ",signif(summary(fit)$adj.r.squared, 3),
"Intercept =",signif(fit$coef[[1]],3 )
))
}
ggplotRegression(lm(LPS.Response.final$log2FC.MIA ~ LPS.Response.final$log2FC.CON, data=LPS.Response.final)) + geom_abline(slope=1) + labs(x="CON logFC (LPS/SAL)", y="MIA logFC (LPS/SAL)") + theme_cowplot() + xlim(-30,35) + ylim(-30,35)
Statistical Concordance of MIA and CON LPS Response
library(epiR)
epi.ccc(LPS.Response.final$log2FC.CON, LPS.Response.final$log2FC.MIA, ci = "z-transform", conf.level = 0.95)$rho.c
## est lower upper
## 1 0.749 0.738 0.759
sessionInfo()
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] epiR_2.0.26 survival_3.2-11
## [3] ensembldb_2.14.1 AnnotationFilter_1.14.0
## [5] GenomicFeatures_1.42.3 AnnotationHub_2.22.1
## [7] BiocFileCache_1.14.0 dbplyr_2.1.1
## [9] org.Mm.eg.db_3.12.0 AnnotationDbi_1.52.0
## [11] DESeq2_1.30.1 SummarizedExperiment_1.20.0
## [13] Biobase_2.50.0 MatrixGenerics_1.2.1
## [15] matrixStats_0.59.0 GenomicRanges_1.42.0
## [17] GenomeInfoDb_1.26.7 IRanges_2.24.1
## [19] S4Vectors_0.28.1 BiocGenerics_0.36.1
## [21] forcats_0.5.1 stringr_1.4.0
## [23] dplyr_1.0.7 purrr_0.3.4
## [25] readr_1.4.0 tidyr_1.1.3
## [27] tibble_3.1.2 tidyverse_1.3.1
## [29] cowplot_1.1.1 ggplot2_3.3.5
##
## loaded via a namespace (and not attached):
## [1] colorspace_2.0-2 ellipsis_0.3.2
## [3] XVector_0.30.0 fs_1.5.0
## [5] rstudioapi_0.13 farver_2.1.0
## [7] bit64_4.0.5 interactiveDisplayBase_1.28.0
## [9] fansi_0.5.0 lubridate_1.7.10
## [11] xml2_1.3.2 splines_4.0.3
## [13] cachem_1.0.5 geneplotter_1.68.0
## [15] knitr_1.33 jsonlite_1.7.2
## [17] Rsamtools_2.6.0 broom_0.7.8
## [19] annotate_1.68.0 shiny_1.6.0
## [21] BiocManager_1.30.16 compiler_4.0.3
## [23] httr_1.4.2 backports_1.2.1
## [25] lazyeval_0.2.2 assertthat_0.2.1
## [27] Matrix_1.3-4 fastmap_1.1.0
## [29] cli_3.0.0 later_1.2.0
## [31] prettyunits_1.1.1 htmltools_0.5.1.1
## [33] tools_4.0.3 gtable_0.3.0
## [35] glue_1.4.2 GenomeInfoDbData_1.2.4
## [37] rappdirs_0.3.3 Rcpp_1.0.7
## [39] cellranger_1.1.0 jquerylib_0.1.4
## [41] Biostrings_2.58.0 vctrs_0.3.8
## [43] nlme_3.1-152 rtracklayer_1.50.0
## [45] xfun_0.24 rvest_1.0.0
## [47] mime_0.11 lifecycle_1.0.0
## [49] XML_3.99-0.6 zlibbioc_1.36.0
## [51] scales_1.1.1 ProtGenerics_1.22.0
## [53] hms_1.1.0 promises_1.2.0.1
## [55] RColorBrewer_1.1-2 yaml_2.2.1
## [57] curl_4.3.2 memoise_2.0.0
## [59] pander_0.6.4 sass_0.4.0
## [61] biomaRt_2.46.3 stringi_1.6.2
## [63] RSQLite_2.2.7 highr_0.9
## [65] BiocVersion_3.12.0 genefilter_1.72.1
## [67] BiocParallel_1.24.1 rlang_0.4.11
## [69] pkgconfig_2.0.3 bitops_1.0-7
## [71] evaluate_0.14 lattice_0.20-44
## [73] labeling_0.4.2 GenomicAlignments_1.26.0
## [75] bit_4.0.4 tidyselect_1.1.1
## [77] magrittr_2.0.1 R6_2.5.0
## [79] generics_0.1.0 DelayedArray_0.16.3
## [81] DBI_1.1.1 mgcv_1.8-36
## [83] pillar_1.6.1 haven_2.4.1
## [85] withr_2.4.2 RCurl_1.98-1.3
## [87] modelr_0.1.8 crayon_1.4.1
## [89] utf8_1.2.1 rmarkdown_2.9
## [91] progress_1.2.2 locfit_1.5-9.4
## [93] grid_4.0.3 readxl_1.3.1
## [95] blob_1.2.1 reprex_2.0.0
## [97] digest_0.6.27 xtable_1.8-4
## [99] httpuv_1.6.1 openssl_1.4.4
## [101] munsell_0.5.0 BiasedUrn_1.07
## [103] bslib_0.2.5.1 askpass_1.1