#---------------------------------------------
# Statistic functions to plot in `CHMI` Study
#---------------------------------------------
### p_values & boxplot for 'log10_mfi'--------------------------------------------------------------
# '@export
chmi.stat_plot.pval_boxplot <- function(
phen, pval_phen,
l_igg,
var_x, var_y,
aes_color,
aes_geom = c('double', 'variable'))
## `variable` have two or more levels
## `double` only has one level
{
# enquos
var_x <- enquo(var_x)
var_y <- enquo(var_y)
aes_color <- enquo(aes_color)
# list
l_plots <- list()
# loop
for(i in seq_along(l_igg)) {
# plot
p_1 <- ggplot(
data = phen %>% filter(t_igg == l_igg[i]),
aes(x = !!var_x, y = !!var_y, color = !!aes_color)) +
geom_boxplot(aes(x = factor(!!var_x), y = !!var_y), outlier.shape = NA) +
geom_jitter(aes(x = factor(!!var_x), y = !!var_y), width = .2, height = 0, alpha = .25) +
facet_wrap(~ antigen_ag) +
scale_color_viridis(discrete = TRUE, begin = 0, end = .5, option = 'D') +
labs(title = paste0(l_igg[i]), x = paste0(quo_name(var_x)), y = paste0(quo_name(var_y)))
## argument `aes_geom`
if(aes_geom == 'variable') {
p_1 <- p_1 +
geom_text(
data = pval_phen %>% filter(t_igg == l_igg[i]),
aes(x = !!var_x, y = max(phen[, quo_name(var_y)] + 1), label = adjust_text),
family = 'Georgia', colour = 'black', size = 4)
} else {
p_1 <- p_1 +
geom_text(
data = pval_phen %>% filter(t_igg == l_igg[i]),
aes(x = 2, y = max(phen[, quo_name(var_y)] + 1), label = adjust_text),
family = 'Georgia', colour = 'black', size = 4)
}
# list
l_plots[[length(l_plots) + 1]] = p_1
}
# return
return(l_plots)
}
### boxplot & p_values for 'mfi' in 'fold changes'----------------------------------------
# '@export
chmi.stat_plot.pval_boxfold <- function(
phen, pval_phen,
l_group,
l_igg,
var_y, aes_color)
{
# enquos
var_y <- enquo(var_y)
aes_color <- enquo(aes_color)
# list
l_plots <- list()
# loop
for(i in seq_along(l_igg)) {
# plot
p_1 <- ggplot(
data = phen %>% filter(t_igg == l_igg[i]),
aes(x = t2_point_ag, y = !!var_y, color = !!aes_color)) +
geom_boxplot(outlier.shape = NA) +
geom_point(
position = position_jitterdodge(jitter.width = .2, jitter.height = 0),
alpha = .25) +
facet_wrap(~ antigen_ag) +
scale_color_viridis(discrete = TRUE, begin = 0, end = .5, option = 'D') +
geom_text(data = pval_phen %>% filter(t_igg == l_igg[i]),
aes(x = t2_point_ag, y = p_ymax, label = adjust_text),
family = 'Georgia', colour = 'black', size = 3) +
labs(title = paste0(l_igg[i]), x = 't2_point', y = paste0(quo_name(var_y)))
# list
l_plots[[length(l_plots) + 1]] = p_1
}
# return
return(l_plots)
}
### medians & geom_line for 'mfi'-----------------------------------------------------------
# '@export
chmi.stat_plot.median_line <- function(
phen,
l_igg,
var_x, var_y,
aes_color)
{
# enquos
var_x <- enquo(var_x)
var_y <- enquo(var_y)
aes_color <- enquo(aes_color)
# list
l_plots <- list()
# loop
for(i in seq_along(l_igg)) {
# plot
p_1 <- ggplot(
data = phen %>% filter(t_igg == l_igg[[i]]),
aes(x = !!var_x, y = !!var_y, color = !!aes_color, group = original_id)) +
geom_line(alpha = .2) +
stat_summary(fun.y = 'median',
aes(group = !!aes_color), geom = 'point', size = 2, alpha = .7) +
facet_wrap(~ antigen_ag) +
scale_color_viridis(discrete = T, begin = 0, end = .5, option = 'D')
# list
l_plots[[length(l_plots) + 1]] = p_1
}
# return
return(l_plots)
}
### or plots ---------------------------------------------------------------------------------------
# '@export
chmi.stat_plot.or_heatmap <- function(
phen,
l_varsy, l_varsx,
arg_igg = c('IgG', 'IgG1', 'IgG2', 'IgG3', 'IgG4', 'IgM'),
arg_plot = c('or', 'raw_pval', 'adjust_pval', 'or_raw', 'or_adj'),
stat = c('ggpubr'))
{
# stat
stat <- match.arg(stat, c('ggpubr'))
### variables to avoid & updates
l_vars <- c('null_deviance', 'df_null', 'log_lik', 'aic', 'bic', 'deviance', 'df_residual',
'r_squared', 'adj_r_squared', 'sigma', 'mod_statistic', 'mod_pval', 'df', 'deviance', 'df_residual',
'std_error', 't_value')
# fix round variables
phen[, c('or')] <- round(phen[, c('or')], 2)
phen[, c('raw_pval', 'adjust_pval')] <- round(phen[, c('raw_pval', 'adjust_pval')], 5)
### argument to filter `l_varsy`
if(!missing(l_varsy)) {
# filter 'phen'
phen <- phen %>%
filter_at(vars(vars_y), all_vars(. %in% l_varsy)) %>%
droplevels()
}
### argument to filter `l_varsx`
if(!missing(l_varsx)) {
# filter 'phen'
phen <- phen %>%
filter_at(vars(vars_x), all_vars(. %in% l_varsx)) %>%
droplevels()
}
### argument to filter `igg`
if(!missing(arg_igg)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(t_igg), all_vars(. %in% arg_igg)) %>%
droplevels()
}
# group_split
phen_1 <- phen %>%
select(-one_of(l_vars)) %>%
filter(complete.cases(.)) %>%
droplevels() %>%
group_by(vars_y, vars_x) %>%
group_split()
# list
l_or <- list()
l_raw <- list()
l_adj <- list()
l_orraw <- list()
l_oradj <- list()
### make all plots
# loop
for(i in seq_along(phen_1)) {
# plot 'or'
p_1 <- ggplot(
phen_1[[i]],
aes(x = t2_point_ag, y = antigen_ag, fill = or)) +
geom_tile() +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[i]] %>% select(or)),
max(phen_1[[i]] %>% select(or)), by = 1),
limits = c(
min(phen_1[[i]] %>% select(or)),
max(phen_1[[i]] %>% select(or)))) +
facet_wrap( ~ t_igg, nrow = 1, scales = 'free_x') +
labs(title = unique(paste0(
droplevels(phen_1[[i]][['vars_y']]), ' ~ ',
droplevels(phen_1[[i]][['vars_x']]))),
x = 't_point', y = 'antigen', fill = 'odd_ratio')
# plot 'raw_pval'
p_2 <- ggplot(
phen_1[[i]],
aes(x = t2_point_ag, y = antigen_ag, fill = raw_pval)) +
geom_tile() +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[i]] %>% select(raw_pval)),
max(phen_1[[i]] %>% select(raw_pval)), by = .5),
limits = c(
min(phen_1[[i]] %>% select(raw_pval)),
max(phen_1[[i]] %>% select(raw_pval)))) +
facet_wrap( ~ t_igg, nrow = 1, scales = 'free_x') +
labs(title = unique(paste0(
droplevels(phen_1[[i]][['vars_y']]), ' ~ ',
droplevels(phen_1[[i]][['vars_x']]))),
x = 't_point', y = 'antigen', fill = 'raw_pval')
# plot 'adj_pval'
p_3 <- ggplot(
phen_1[[i]],
aes(x = t2_point_ag, y = antigen_ag, fill = adjust_pval)) +
geom_tile() +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[i]] %>% select(adjust_pval)),
max(phen_1[[i]] %>% select(adjust_pval)), by = .5),
limits = c(
min(phen_1[[i]] %>% select(adjust_pval)),
max(phen_1[[i]] %>% select(adjust_pval)))) +
facet_wrap( ~ t_igg, nrow = 1, scales = 'free_x') +
labs(title = unique(paste0(
droplevels(phen_1[[i]][['vars_y']]), ' ~ ',
droplevels(phen_1[[i]][['vars_x']]))),
x = 't_point', y = 'antigen', fill = 'adj_pval')
# combine list
p_raw <- ggpubr::ggarrange(p_1, p_2,
ncol = 1, nrow = 2,
common.legend = FALSE,
legend = 'right',
labels = c('a. odd_ratio', 'b. raw_pvalues'))
p_adj <- ggpubr::ggarrange(p_1, p_3,
ncol = 1, nrow = 2,
common.legend = FALSE,
legend = 'right',
labels = c('a. odd_ratio', 'b. adjust_pvalues'))
# list
l_or[[length(l_or) + 1]] = p_1
l_raw[[length(l_raw) + 1]] = p_2
l_adj[[length(l_adj) + 1]] = p_3
l_orraw[[length(l_orraw) + 1]] = p_raw
l_oradj[[length(l_oradj) + 1]] = p_adj
}
### argument 'plots'
if(arg_plot == 'or') {
# return
return(l_or)
} else if(arg_plot == 'raw_pval') {
# return
return(l_raw)
} else if(arg_plot == 'adjust_pval') {
# return
return(l_adj)
} else if (arg_plot == 'or_raw') {
# return
return(l_orraw)
} else {
# return
return(l_oradj)
}
}
### beta plots ---------------------------------------------------------------------------------------
# '@export
chmi.stat_plot.heatmap <- function(
phen,
l_group,
l_varsy = c('NULL'),
l_varsx = c('NULL'),
arg_igg = c('IgG', 'IgG1', 'IgG2', 'IgG3', 'IgG4', 'IgM', 'NULL'),
arg_plot = c('beta', 'or'))
{
### list of variables to choose
l_vars_arg <- c(l_group, 'vars_y', 'vars_x', arg_plot)
### update 'data'
phen_1 <- phen %>%
select(one_of(l_vars_arg)) %>%
filter(complete.cases(.)) %>%
droplevels()
### argument to filter `l_varsy`, `l_vars_x`
# list
l_vars_1 <- c(l_varsy[2], l_varsx[1], arg_igg)
l_vars_2 <- c(l_varsy[1], l_varsx[1])
l_vars_3 <- c(l_varsy[1])
l_vars_4 <- NULL
### conditions
if(length(l_vars_1) != 0) {
phen_1 <- phen_1 %>% filter_at(vars(t_igg, vars_y, vars_x), all_vars(. %in% l_vars_1))
} else {
phen_1
}
### group_split
phen_1 <- phen_1
group_by_at(l_names[c(2:3)]) %>%
group_split()
# list
l_plots <- list()
### make all plots
# loop
for(i in seq_along(phen_1)) {
# plot 'or'
p_1 <- ggplot(
phen_1[[1]],
aes(x = t2_point_ag, y = antigen_ag, fill = beta)) +
geom_tile() +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[1]] %>% select(beta)),
max(phen_1[[1]] %>% select(beta)), by = 1),
limits = c(
min(phen_1[[1]] %>% select(beta)),
max(phen_1[[1]] %>% select(beta)))) +
facet_wrap( ~ t_igg, nrow = 1, scales = 'free_x') +
labs(title = unique(paste0(
droplevels(phen_1[[1]][['vars_y']]), ' ~ ',
droplevels(phen_1[[1]][['vars_x']]))),
x = 't_point', y = 'antigen', fill = paste0(arg_plot))
# list
l_plots[[length(l_plots) + 1]] = p_1
}
### return
return(l_plots)
}
### geom_point plots (for OR)-------------------------------------------------------------------------------
# '@export
chmi.stat_plot.pointrange_or <- function(
phen,
l_varsy, l_varsx,
arg_igg = c('IgG', 'IgG1', 'IgG2', 'IgG3', 'IgG4', 'IgM'),
stat = c('ggpubr'))
{
# stat
stat <- match.arg(stat, c('ggpubr'))
# variables to avoid
l_vars <- c('null_deviance', 'df_null', 'log_lik', 'aic', 'bic', 'deviance', 'df_residual',
'r_squared', 'adj_r_squared', 'sigma', 'mod_statistic', 'mod_pval', 'df', 'deviance', 'df_residual',
'std_error', 't_value')
# fix round variables
phen[, c('or', 'ci_lower', 'ci_upper')] <- round(phen[, c('or', 'ci_lower', 'ci_upper')], 2)
phen[, c('raw_pval', 'adjust_pval')] <- round(phen[, c('raw_pval', 'adjust_pval')], 5)
# argument to filter `l_varsy` | `l_varsx`
if(!missing(l_varsy)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(vars_y), all_vars(. %in% l_varsy)) %>%
droplevels()
}
if(!missing(l_varsx)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(vars_x), all_vars(. %in% l_varsx)) %>%
droplevels()
}
if(!missing(arg_igg)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(t_igg), all_vars(. %in% arg_igg)) %>%
droplevels()
}
# group_split
phen_1 <- phen %>%
select(-one_of(l_vars)) %>%
filter(complete.cases(.)) %>%
droplevels() %>%
group_by(vars_y, vars_x, t2_point_ag) %>%
group_split()
# list
l_plots <- list()
### make all plots
# loop
for(i in seq_along(phen_1)) {
# plot 'point_range'
p_1 <- ggplot(
phen_1[[i]],
aes(x = antigen_ag, y = or)) +
geom_linerange(aes(ymin = ci_lower, ymax = ci_upper, colour = or), size = 0.2) +
geom_pointrange(aes(ymin = ci_lower, ymax = ci_upper, colour = or), size = 0.2) +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[i]] %>% select(or)),
max(phen_1[[i]] %>% select(or)), by = 1),
limits = c(
min(phen_1[[i]] %>% select(or)),
max(phen_1[[i]] %>% select(or)))) +
labs(title = unique(paste0(
droplevels(phen_1[[i]][['vars_y']]), ' ~ ',
droplevels(phen_1[[i]][['vars_x']]), ' (',
droplevels(phen_1[[i]][['t2_point_ag']]), ')' )),
x = 'antigen', y = 'odd_ratio', fill = 'odd_ratio') +
geom_hline (yintercept = 0, linetype = 3) +
coord_flip() +
facet_wrap( ~ t_igg, ncol = length(levels(phen_1[[i]][['t_igg']])), scales = 'free_x') +
geom_text(
data = phen_1[[i]] %>% filter(raw_pval <= 0.05),
aes(label = paste(or)),
vjust = 0.7, hjust = -0.7, nudge_x = 0.4,
size = 3, family = 'Georgia', colour = 'black',
check_overlap = TRUE)
# list
l_plots[[length(l_plots) + 1]] = p_1
}
# return
return(l_plots)
}
### geom_point plots (for BETA)-------------------------------------------------------------------------------
# '@export
chmi.stat_plot.pointrange_beta <- function(
phen,
l_varsy, l_varsx,
arg_igg = c('IgG', 'IgG1', 'IgG2', 'IgG3', 'IgG4'),
stat = c('ggpubr'))
{
# stat
stat <- match.arg(stat, c('ggpubr'))
# variables to avoid
l_vars <- c('null_deviance', 'df_null', 'log_lik', 'aic', 'bic', 'deviance', 'df_residual',
'r_squared', 'adj_r_squared', 'sigma', 'mod_statistic', 'mod_pval', 'df', 'deviance', 'df_residual',
'std_error', 't_value')
# fix round variables
phen[, c('beta', 'ci_lower', 'ci_upper')] <- round(phen[, c('beta', 'ci_lower', 'ci_upper')], 2)
phen[, c('raw_pval', 'adjust_pval')] <- round(phen[, c('raw_pval', 'adjust_pval')], 5)
# argument to filter `l_varsy` | `l_varsx`
if(!missing(l_varsy)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(vars_y), all_vars(. %in% l_varsy)) %>%
droplevels()
}
if(!missing(l_varsx)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(vars_x), all_vars(. %in% l_varsx)) %>%
droplevels()
}
if(!missing(arg_igg)) {
# filter 'phen_1'
phen <- phen %>%
filter_at(vars(isotype), all_vars(. %in% arg_igg)) %>%
droplevels()
}
phen <- phen %>%
mutate(
beta_signif_raw = factor(ifelse(raw_pval <= 0.05, beta, '')),
beta_signif_adj = factor(ifelse(adjust_pval <= 0.05, paste0('*'), '')))
# group_split
phen_1 <- phen %>%
select(-analyte_block) %>%
select(-one_of(l_vars)) %>%
filter(complete.cases(.)) %>%
droplevels() %>%
group_by(subset, vars_x) %>%
group_split()
# list
l_plots <- list()
### make all plots
# loop
for(i in seq_along(phen_1)) {
# plot 'point_range'
p_1 <- ggplot(
phen_1[[1]],
aes(x = analyte_label, y = beta)) +
geom_linerange(aes(ymin = ci_lower, ymax = ci_upper), size = 1) +
geom_pointrange(aes(ymin = ci_lower, ymax = ci_upper), size = 0.6) +
scale_fill_gradientn(colours = c('#56B4E9', 'white', '#8B0000'),
breaks = seq(
min(phen_1[[1]] %>% select(beta)),
max(phen_1[[1]] %>% select(beta)), by = 0.5),
limits = c(
min(phen_1[[1]] %>% select(beta)),
max(phen_1[[1]] %>% select(beta)))) +
labs(title = paste0(phen_1[[1]][['subset']]), x = '', y = 'beta', fill = 'beta', size = 5) +
geom_hline (yintercept = 0, linetype = 'dotdash', size = 1.5, color = 'coral') +
#ylim(-0.8, 1) +
coord_flip() +
facet_wrap( ~ vars_x_label + isotype, ncol = 5, scales = 'free_x') +
theme(
plot.title = element_text(color = "black", size = 22, face = "bold", vjust = 6),
plot.margin = unit(c(4, 4, 4, 1.8), "cm"),
axis.title.x = element_text(color = "black", size = 16, face = "bold", vjust = -0.7),
axis.title.y = element_text(color = "black", size = 16, face = "bold"),
panel.border = element_rect(color = "black", fill = NA, size = 0.5),
panel.background = element_rect(fill = "white", colour = "white", size = 0.5, linetype = 'solid'),
panel.grid.major = element_line(size = 0.5, linetype = 'dashed', colour = "grey"),
panel.grid.minor = element_line(size = 0.25, linetype = 'dashed', colour = "grey"),
strip.text = element_text(color = "black", size = 14),
panel.spacing = unit(0.8, "lines"),
axis.text = element_text(size = 13)) +
geom_text(
data = phen_1[[1]],
aes(label = paste0(beta_signif_raw, beta_signif_adj)),
vjust = 0.3, hjust = 1.5, nudge_x = -0.4, nudge_y = 0.02,
size = 4, family = 'Arial', colour = 'black',
check_overlap = TRUE)
# list
l_plots[[length(l_plots) + 1]] = p_1
}
# return
return(l_plots)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.