Figure 1a,b: Confounder analysis in drug-host-microbiome associations
Stacked bar charts show variance explained (R squared) by explanatory variable group and type of measurement.
Violin plot representing confounder analysis of features associated with T2D group compared to control.
Vertical axis shows Cliff’s delta, with density plots showing all features significant before drug deconfounding (two-sided MWU FDR < 0.1); number of features in each violin are listed.
All “naïve associations” (yellow violin) are separated into confounded or ambiguously/confidently deconfounded (blue, purple and red violins). About 45% of 2845 raw associations are reducible to medication effects (blue violin). Green violins show breakdown of significant drug confounders by drug category.
Required file in the input_data folder:
- Supplementary_Table_5_2019-09-13434.xlsx
Supplementary Table 5 Multivariate breakdown of variance. The table includes multivariate breakdown of variance for each feature space by each predictor, a summary of all interaction terms in the models, and the results of the confounder analysis for all patient groups analyzed in the MetaCardis cohort.
Load necessary libraries.
library (ggplot2)
library (readxl)Plot bar plot of explained variance per feature type
Read input data from Supplementary Table 5.
fileFolder = './input_data/'
data_rsq <- read_excel(paste0(fileFolder,
"Supplementary_Table_5_2019-09-13434.xlsx"),
sheet = "Summarized variance explained")Select feature space to plot.
# Species abundance, KEGG module abundance, CARD Antibiotic resistance gene abundance, Host lipoproteins, host serum metabolites
select_space <- c ("MGS", "KEGG module", "CARD ARO", "Lipoproteins", "Serum")
subdata <- subset (data_rsq,`Feature space` %in% select_space)
subdata = subdata[subdata$Category != 'Intersection',]
subdata = subdata[subdata$Category != 'Stool form',]
# sort types of features
subdata$`Feature space` <- factor(subdata$`Feature space`,
levels = select_space)Assemble bar plot.
p <- ggplot (subdata[order(subdata$Category), ], aes (x = `Feature space`, y = `Total adjusted R squared from significant factors`)) +
theme_classic () +
geom_bar (aes (fill = Category), stat = "identity", color = "black") +
xlab ("") + ylab ("Variance explained") + scale_fill_brewer (type = "qual", palette = "Paired")
print (p)
#Uncomment to print to file
#pdf (file = paste0 ("fig1a_bar_rsq.pdf"))
#print (p)
#dev.off ()Plot violin plots
Define accessory function.
f <- function (y) c (label = length (y), y = median (y))Read input data from Supplementary Table 5.
fileFolder = './input_data/'
data <- read_excel(paste0(fileFolder,
"Supplementary_Table_5_2019-09-13434.xlsx"),
sheet = "Confounder analysis")Prepare data for plotting.
data$StatusD <- paste (data$Status, as.character (sign (data$`Effect size`)))
# Get unique feature spaces
spaces <- unique(data$`Feature space`)
data <- subset (data,`Feature space` %in% spaces)
# Get unique contrasts (patient groups)
contrasts <- unique(data$Contrast)
#Get unique confounder status
stati <- unique(data$Status)
data <- subset (data, Status %in% stati)Plot violins for T2D group.
Change group index or loop over all contrasts to plot for a different patient group.
# Group 4 stands for T2D group.
for (i in c(4)) { #for (i in 1:7)
subdata <- subset (data, Contrast == contrasts [i])
p <- ggplot (subdata, aes (x = as.factor (Status), y = `Effect size`)) + scale_x_discrete(limits = rev(levels(as.factor (data$Status)))) +
geom_violin (aes (fill = Status), scale = "count", draw_quantiles = c (0.25, 0.5, 0.75), size = 0.5, alpha = 0.3, width = 30, position = position_dodge ()) +
geom_hline (yintercept = 0) +
xlab ("") +
ylab ("Effect size (Cliff's delta)") +
ggtitle(contrasts[i]) +
stat_summary (aes (col = subdata$Status, y = - 1), fun.data = f, geom = "text", vjust = 0, size = 4, angle = 45) +
theme (axis.text.x = element_text (angle = -90, hjust = 0, size = 12),
axis.text.y = element_text (size = 16),
axis.title.y = element_text (size = 16), legend.position = "none") +
scale_fill_manual (values = c ("#00FF00", "#00FF00", "#00FF00", "#00FF00", "#00FF00", "#00FF00", "#0000FF", "#FF00FF", "#FF0000", "#DDCC00")) +
scale_color_manual (values = c ("#009900", "#009900", "#009900", "#009900", "#009900", "#009900", "#0000FF", "#990099", "#990000", "#998800"))
#Uncomment to print to file
#pdf (file = paste0 ("fig1b_flame_", i, ".pdf"), width = 20, height = 6)
print (p)
#dev.off ()
}
sessionInfo()## R version 3.6.1 (2019-07-05)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 17763)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=English_United States.1252
## [2] LC_CTYPE=English_United States.1252
## [3] LC_MONETARY=English_United States.1252
## [4] LC_NUMERIC=C
## [5] LC_TIME=English_United States.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] readxl_1.3.1 ggplot2_3.3.3
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.6 highr_0.9 RColorBrewer_1.1-2 cellranger_1.1.0
## [5] bslib_0.2.4 compiler_3.6.1 pillar_1.6.0 jquerylib_0.1.4
## [9] tools_3.6.1 digest_0.6.27 jsonlite_1.7.2 evaluate_0.14
## [13] lifecycle_1.0.0 tibble_3.1.0 gtable_0.3.0 pkgconfig_2.0.3
## [17] rlang_0.4.10 DBI_1.1.1 yaml_2.2.0 xfun_0.22
## [21] withr_2.4.2 stringr_1.4.0 dplyr_1.0.5 knitr_1.33
## [25] generics_0.1.0 sass_0.3.1 vctrs_0.3.6 grid_3.6.1
## [29] tidyselect_1.1.0 glue_1.4.2 R6_2.5.0 fansi_0.4.2
## [33] rmarkdown_2.7 farver_2.1.0 purrr_0.3.4 magrittr_2.0.1
## [37] scales_1.1.1 ellipsis_0.3.1 htmltools_0.5.1.1 assertthat_0.2.1
## [41] colorspace_2.0-0 labeling_0.4.2 utf8_1.2.1 stringi_1.5.3
## [45] munsell_0.5.0 crayon_1.4.1