diff --git a/boxplot_stat_all.R b/boxplot_stat_all.R index 5288fb6..a71969a 100755 --- a/boxplot_stat_all.R +++ b/boxplot_stat_all.R @@ -28,7 +28,7 @@ fp_npa = length(unique(lf_fp_npa$mosaic)); fp_npa cat("\nPlotting boxplots with stats for:", my_sample_npa , "\n========================================================\n") -plots_npa = doMyPlots(lf_fp_npa) +plots_npa = doMyPlotsStats(lf_fp_npa) npa_plot = ggpubr::ggarrange(plotlist = plots_npa , align = "hv" , ncol = 7 @@ -94,7 +94,7 @@ fp_serum = length(unique(lf_fp_serum$mosaic)); fp_serum cat("\nPlotting boxplots with stats for:", my_sample_serum , "\n========================================================\n") -plots_serum = doMyPlots(lf_fp_serum) +plots_serum = doMyPlotsStats(lf_fp_serum) serum_plot = ggpubr::ggarrange(plotlist = plots_serum , align = "hv" , ncol = 7 diff --git a/boxplot_stat_function.R b/boxplot_stat_function.R index 070adf9..e8f0e1c 100755 --- a/boxplot_stat_function.R +++ b/boxplot_stat_function.R @@ -5,7 +5,7 @@ getwd() ############################################################ # TASK: boxplots with stats ############################################################ -doMyPlots <- function(df) { +doMyPlotsStats <- function(df) { mediators = levels(as.factor(df$mediator)) plots <- list() @@ -61,4 +61,5 @@ doMyPlots <- function(df) { return(plots) } -############################################################### \ No newline at end of file +############################################################### + diff --git a/boxplot_stat_function_test.R b/boxplot_stat_function_test.R index f9a86fd..f75f2d3 100755 --- a/boxplot_stat_function_test.R +++ b/boxplot_stat_function_test.R @@ -22,8 +22,16 @@ pdf(output_boxplot_stats, width=22, height=16) #med_names = c("eotaxin3", "il12p70", "itac", "il13") #lf_test = lf_fp_npa[lf_fp_npa$mediator%in%med_names,] -#baz=cowplot::plot_grid(plotlist=foo, align = 'hv', ncol=2, nrow=2) -foo_annot = annotate_figure(foo +baz=cowplot::plot_grid(plotlist=foo, align = 'hv', ncol=2, nrow=2) +baz + +sam_plot = ggpubr::ggarrange(plotlist = foo + , align = "hv" + , ncol = 7 + , nrow = 5 + , common.legend = T) +sam_plot +foo_annot = annotate_figure(sam_plot , top = text_grob(my_sample_name , color = "blue" , face = "bold" @@ -41,3 +49,80 @@ foo_annot = annotate_figure(foo foo_annot dev.off() ################################################################## +# ind components of function +my_sample_name = "SAM" +to_remove = c("ifna2a" + , "il29") + +lf_fp_sam2 = lf_fp_sam[!lf_fp_sam$mediator%in%to_remove,] +lf_fp_sam3 = lf_fp_sam[lf_fp_sam$mediator%in%to_remove,] + +df_test = lf_fp_sam3 +df_test = lf_fp_sam2 +table(df_test$mediator) + +p = ggplot(df_test)+ geom_boxplot(aes(x = timepoint + , y = value + , color = obesity + #, palette = c("#00BFC4", "#F8766D") +))+ + scale_colour_manual(values=c("#00BFC4", "#F8766D")) + + theme(axis.text.x = element_text(size = 15) + , axis.text.y = element_text(size = 15 + , angle = 0 + , hjust = 1 + , vjust = 0) + , axis.title.x = element_blank() + , axis.title.y = element_blank() + , legend.position = "none" + , plot.subtitle = element_text(size = 20, hjust = 0.5) + , plot.title = element_text(size = 20, hjust = 0.5)) + + labs(title = "TEST") +p + +#-------- +# stats +#--------- +to_remove = c("ifna2a" + , "il29") +lf_fp_sam3 = lf_fp_sam[lf_fp_sam$mediator%in%to_remove,] + +lf_fp_sam4 = lf_fp_sam3[lf_fp_sam3$timepoint == "t3",] +table(lf_fp_sam4$timepoint) + +wilcox_test(value~obesity, data = lf_fp_sam4) + +table(lf_fp_sam3$mediator, lf_fp_sam3$timepoint) +time = c("t3") +lf_fp_sam3 = lf_fp_sam3[!lf_fp_sam3$timepoint%in%time,] +table(lf_fp_sam3$mediator, lf_fp_sam3$timepoint) + +med = c("il29") +time2 = c("t2") + +lf_fp_sam3 = lf_fp_sam3[!lf_fp_sam3$mediator%in%med && lf_fp_sam3$timepoint%in%time2,] +table(lf_fp_sam3$mediator, lf_fp_sam3$timepoint) + +df_test = lf_fp_sam3 + + +table(df_test$mediator) + + +stat_df <- df_test %>% + group_by(timepoint, mediator) %>% + wilcox_test(value ~ obesity, paired = F) %>% + add_significance("p") +stat_df$p_format = round(stat_df$p, digits = 3) + + +stat_df <- stat_df %>% + add_xy_position(x = "timepoint", dodge = 0.8) +p2 = p + stat_pvalue_manual(stat_df + #, y.position = max_y + , label = "{p_format} {p.signif}" + , hide.ns = T + , tip.length = 0)+ + scale_y_continuous(expand = expansion(mult = c(0.05, 0.25))) + +p2 diff --git a/boxplot_stat_log.R b/boxplot_stat_log.R new file mode 100644 index 0000000..0953113 --- /dev/null +++ b/boxplot_stat_log.R @@ -0,0 +1,153 @@ +#!/usr/bin/Rscript +getwd() +setwd("~/git/mosaic_2020/") +getwd() +############################################################### +#============= +# Input +#============= +source("boxplot_stat_function.R") +source("plot_data.R") + +#============= +# Output +#============= +outfile_bp = paste0("boxplots_stats_npa_serum_LOG", ".pdf") +output_boxplot_stats = paste0(outdir_plots, outfile_bp); output_boxplot_stats +pdf(output_boxplot_stats, width=22, height=16) +############################################################### +#=============================== +# data assignment for plots +#=============================== +#---------------------------------------- +#------- +# NPA +#------- +# vl data juggling with 0! +vl_temp_fix = lf_fp_npa[lf_fp_npa$mediator=="vl_pfu_ul",] +vl_temp_fix$value[1:20] +vl_temp_fix$value[vl_temp_fix$value == 0] <- 1 +vl_temp_fix$value[1:20] + +orig_nrow = nrow(lf_fp_npa); orig_nrow; table(lf_fp_npa$mediator) +lf_fp_npa = lf_fp_npa[!lf_fp_npa$mediator=="vl_pfu_ul",] +temp_nrow = nrow(lf_fp_npa); temp_nrow; table(lf_fp_npa$mediator) + +# add vl_temp_fix back +lf_fp_npa = rbind(lf_fp_npa, vl_temp_fix) +revised_nrow = nrow(lf_fp_npa); revised_nrow + +orig_nrow == revised_nrow + +head(lf_fp_npa$value) + +# log all values and reassign +lf_fp_npa$value = log10(lf_fp_npa$value) +head(lf_fp_npa$value) +vl_temp_fix2 = lf_fp_npa[lf_fp_npa$mediator=="vl_pfu_ul",] +vl_temp_fix2$value[1:20] + +#------------------------------------------ +my_sample_npa = "NPA" +fp_npa = length(unique(lf_fp_npa$mosaic)); fp_npa + +cat("\nPlotting boxplots with stats for:", my_sample_npa + , "\n========================================================\n") + +plots_npa = doMyPlotsStats(lf_fp_npa) +npa_plot = ggpubr::ggarrange(plotlist = plots_npa + , align = "hv" + , ncol = 7 + , nrow = 5 + , common.legend = T) +#npa_plot +npa_plot_annot = annotate_figure(npa_plot + , top = text_grob(my_sample_npa + , color = "blue" + , face = "bold" + , size = 14) + , bottom = text_grob(paste0("Mosaic data\nFlu positive adults (n=", fp_npa, ")") + , color = "blue" + , hjust = 1 + , x = 0.98 + , face = "italic" + , size = 10) + , left = text_grob("Levels (pg/ml)" + , color = "black" + , rot = 90 + , size = 18)) +npa_plot_annot +#dev.off() +#------------- +# SAM +#------------- +#my_sample_sam = "SAM" +#fp_sam = length(unique(lf_fp_sam$mosaic)); fp_sam + +#cat("\nPlotting boxplots with stats for:", my_sample_sam +# , "\n========================================================\n") + +#plots_sam = doMyPlots(lf_fp_sam) +#sam_plot = ggpubr::ggarrange(plotlist = plots_sam +# , align = "hv" +# , ncol = 7 +# , nrow = 5 +# , common.legend = T) +#sam_plot +#sam_plot_annot = annotate_figure(sam_plot +# , top = text_grob(my_sample_sam +# , color = "blue" +# , face = "bold" +# , size = 14) +# , bottom = text_grob(paste0("Mosaic data\nFlu positive adults (n=", fp_sam, ")") +# , color = "blue" +# , hjust = 1 +# , x = 0.98 +# , face = "italic" +# , size = 10) +# , left = text_grob("Levels (pg/ml)" +# , color = "black" +# , rot = 90 +# , size = 18)) +#sam_plot_annot +#dev.off() +#------------------------------------------ +#------------- +# SERUM +#------------- +head(lf_fp_serum$value) +# reassign +lf_fp_serum$value = log10(lf_fp_serum$value) +head(lf_fp_serum$value) + +my_sample_serum = "serum" +fp_serum = length(unique(lf_fp_serum$mosaic)); fp_serum + +cat("\nPlotting boxplots with stats for:", my_sample_serum + , "\n========================================================\n") + +plots_serum = doMyPlotsStats(lf_fp_serum) +serum_plot = ggpubr::ggarrange(plotlist = plots_serum + , align = "hv" + , ncol = 7 + , nrow = 5 + , common.legend = T) +#serum_plot +serum_plot_annot = annotate_figure(serum_plot + , top = text_grob(my_sample_serum + , color = "blue" + , face = "bold" + , size = 14) + , bottom = text_grob(paste0("Mosaic data\nFlu positive adults (n=", fp_serum, ")") + , color = "blue" + , hjust = 1 + , x = 0.98 + , face = "italic" + , size = 10) + , left = text_grob("Levels (pg/ml)" + , color = "black" + , rot = 90 + , size = 18)) +serum_plot_annot +dev.off() +############################################################### diff --git a/boxplot_with_stats_test.R b/boxplot_with_stats_test.R index da73758..f513cc2 100755 --- a/boxplot_with_stats_test.R +++ b/boxplot_with_stats_test.R @@ -9,7 +9,7 @@ getwd() #============= # Input #============= -source("boxplot_linear.R") +source("plot_data.R") ####################################################### med_names = c("eotaxin3", "il12p70", "itac", "il13") lf_test = lf_fp_npa[lf_fp_npa$mediator%in%med_names,] diff --git a/plot_data.R b/plot_data.R index 5d34538..a1150ca 100755 --- a/plot_data.R +++ b/plot_data.R @@ -21,7 +21,7 @@ metadata_all$mosaic[metadata_all$adult==1 & metadata_all$age<=18] #----------- # npa #----------- -#wf_fp_npa = npa_wf[npa_wf$flustat == 1,] +wf_fp_npa = npa_wf[npa_wf$flustat == 1,] lf_fp_npa = npa_lf[npa_lf$flustat == 1,] lf_fp_npa$timepoint = paste0("t", lf_fp_npa$timepoint) lf_fp_npa$timepoint = as.factor(lf_fp_npa$timepoint) @@ -35,7 +35,7 @@ table(lf_fp_npa$mediator) #----------- # sam #----------- -#wf_fp_sam = samm_wf[samm_wf$flustat == 1,] +wf_fp_sam = sam_wf[sam_wf$flustat == 1,] lf_fp_sam = sam_lf[sam_lf$flustat == 1,] lf_fp_sam$timepoint = paste0("t", lf_fp_sam$timepoint) lf_fp_sam$timepoint = as.factor(lf_fp_sam$timepoint) @@ -49,7 +49,7 @@ table(lf_fp_sam$mediator) #----------- # serum #----------- -#wf_fp_serum = serum_wf[serum_wf$flustat == 1,] +wf_fp_serum = serum_wf[serum_wf$flustat == 1,] lf_fp_serum = serum_lf[serum_lf$flustat == 1,] lf_fp_serum$timepoint = paste0("t", lf_fp_serum$timepoint) lf_fp_serum$timepoint = as.factor(lf_fp_serum$timepoint)