knitr::opts_chunk$set(comment = NA,
  fig.path = '01_figures_aim4_q05_q08_univar/',
  fig.width = 20,
  fig.height = 12,
  results = 'markup',
  tidy = FALSE,
  message = FALSE,
  warning = FALSE,
  echo = FALSE)
alpha <- 0.05

r packages

library(devtools)
library(prettydoc)
# avoid 'objects' in workspace
rm(list = ls())

setwd(Sys.getenv('HOME'))

load_all('~/git/chmitools')

load data

Main bullet points in load data section. chmi.phen() load data by group and aim. filter data by t_igg, antigen_ag and t_point_ag/t2_point_ag columns.

# load data
dat <- chmi.phen(
  study_chmi = 'antibodies',
  aim_ab_chmi = 'aim_4')

dat_fold <- chmi.phen(
  study_chmi = 'antibodies',
  aim_ab_chmi = 'aim_4',
  fold_change = TRUE)

# list to perform different function
  n_igg <- sort(unique(dat$t_igg))
  l_t2point_ag_1 <- c('C-1', 'D7')
  l_t2point_ag_2 <- c('D7','D11_D13')
  l_group <- c('t_igg', 'antigen_ag', 't2_point_ag')
  l_varsy <- c('ppp_time_ori', 'ppp_tbs_qpcr')
  l_varsx <- c('log10_mfi', 'ratio_mfi')
  l_dbl_beta <- c('beta', 'ci_lower', 'ci_upper')
  l_pval <- c('raw_pval', 'adjust_pval')
  l_showed <- c(l_group, l_dbl_or, l_pval)

# subsets for questions
dat_1 <- dat %>%
  filter(t2_point_ag %in% l_t2point_ag_1) %>%
  filter(pcr_pos_pid == '1') %>%
  filter(!is.na(log10_mfi)) %>%
  droplevels() ## Q5 & Q7

dat_2 <- dat_fold %>%
    filter(t2_point_ag %in% l_t2point_ag_2) %>%
  filter(tbs_pos_qpcr == '1') %>%
    filter(!is.na(ratio_mfi)) %>%
    droplevels() %>%
  mutate(semi_t2_point_ag = factor(paste0(immune_status_ag, '_', t2_point_ag))) %>%
  filter(semi_t2_point_ag  != 'naive_D11_D13') ## Q6 & Q8

tables: odd ratios & p_values in glm models

Odd ratios, lower & upper confidence interval, and adjusted p_values. chmi.stat.glm_or() obtain all parameters for all groups. 1st subset (dat_1): variable (pcr_pos_pid), t2_point_ag ('C-1', 'D7', 'D11_D13'). 2nd subset (dat_2): variable (pcr_pos_pid), t2_point_ag ('D7 / C-1'). 3rd subset (dat_1): variable (tbs_pos_qpcr), t2_point_ag ('C-1', 'D7'). 4th subset (dat_2): variable (tbs_pos_qpcr), t2_point_ag ('D7 / C-1'). we are used p.adjust() function to adjust p-values by Benjamini & Hochberg method.

# glm models for `l_varsy` & `log10_mfi`
  tab_1 <- chmi.stat.lm_beta(
        phen = dat_1,
        l_group = l_group,
        l_varsy = l_varsy,
        l_varsx = l_varsx[1],
        arg_multivar = 'univariate',
        arg_info_model = NULL,
      arg_info_beta = NULL,
        format_pval = TRUE,
        padj_method = 'BH')

# glm models for `l_varsy` & `ratio_mfi`
  tab_2 <- chmi.stat.lm_beta(
        phen = dat_2,
        l_group = l_group,
        l_varsy = l_varsy,
        l_varsx = l_varsx[2],
        arg_multivar = 'univariate',
      arg_info_model = NULL,
      arg_info_beta = NULL,
        format_pval = TRUE,
        padj_method = 'BH')


# join `tab_1` & `tab_2`
    tab <- bind_rows(tab_1, tab_2) %>%
        mutate_at(vars(c('t2_point_ag')),  ~ fct_relevel(., c('C-1', 'D7', 'D11_D13'))) %>%
        mutate_at(vars(c('vars_x', 'adjust_text', 'adjust_signif')), as.factor)

tables for q01 & q03

# 'tab_1' datatable
  datatable(tab_1 %>% select(l_showed),
    class = 'cell-border stripe',
    filter = 'top',
    extensions = c('Buttons', 'FixedColumns'),
      options = list(dom = 'Bfrtip',
        buttons = c('csv', 'copy'),
        pageLength = 10,
        scrollX = TRUE,
        text = 'Download')) %>%
    formatRound(l_dbl_or, 3) %>%
    formatRound(l_pval, 5)

 ```

## tables for `q02` & `q04`

```r

# 'tab_2' datatable
  datatable(tab_2 %>% select(l_showed),
    class = 'cell-border stripe',
    filter = 'top',
    extensions = c('Buttons', 'FixedColumns'),
      options = list(dom = 'Bfrtip',
        buttons = c('csv', 'copy'),
        pageLength = 10,
        scrollX = TRUE,
        text = 'Download')) %>%
    formatRound(l_dbl_or, 3) %>%
    formatRound(l_pval, 5)

plots for odd_ratio, raw_pval & adjust_pval

# heatmaps plots
h_beta <- chmi.stat_plot.beta_heatmap(
    phen = tab,
    arg_plot = 'beta')

h_raw <- chmi.stat_plot.beta_heatmap(
    phen = tab,
    arg_plot = 'raw_pval')

h_adj <- chmi.stat_plot.beta_heatmap(
    phen = tab,
    arg_plot = 'adjust_pval')

# pointrange plots
p_odd <- chmi.stat_plot.pointrange(
    phen = tab,
    arg_igg = n_igg)

heatmaps for odd ratios

# heatmap
print(h_odd)

heatmap for raw p_values

# loop
print(h_raw)

heatmap for adjusted p_values

# loop
print(h_adj)

point_range for odd ratios

# point_range
print(p_odd)

save in r_data folder

# save
save.image(file = '~/git/chmitools/projects/04_aim4/r_data/01_aim4_q1_q4_univar.RData')


mvazquezs/chmitools documentation built on May 1, 2020, 2:06 a.m.