tanu/mosaic_2020
1
0
Fork 0
Code Issues Pull requests Releases Wiki Activity

updated un[aired stats with n_obs_complete and removed mediator vitd from sam and npa

Browse source
This commit is contained in:
Tanushree Tunstall 2020-10-30 13:04:45 +00:00
parent 1b5d353060
commit a2dbecd03d
5 changed files with 165 additions and 65 deletions
Download patch file Download diff file Expand all files Collapse all files

17
boxplot.R
View file

@ -30,12 +30,22 @@ lf_fp_npa$timepoint = paste0("t", lf_fp_npa$timepoint)
lf_fp_npa$timepoint = as.factor(lf_fp_npa$timepoint) lf_fp_npa$timepoint = as.factor(lf_fp_npa$timepoint)
lf_fp_npa$obesity = as.factor(lf_fp_npa$obesity) lf_fp_npa$obesity = as.factor(lf_fp_npa$obesity)
table(lf_fp_npa$mediator)
head(lf_fp_npa$value[lf_fp_npa$mediator == "vitd"])
lf_fp_npa = lf_fp_npa[!lf_fp_npa$mediator == "vitd",]
#wf_fp_sam = samm_wf[samm_wf$flustat == 1,] #wf_fp_sam = samm_wf[samm_wf$flustat == 1,]
lf_fp_sam = sam_lf[sam_lf$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 = paste0("t", lf_fp_sam$timepoint)
lf_fp_sam$timepoint = as.factor(lf_fp_sam$timepoint) lf_fp_sam$timepoint = as.factor(lf_fp_sam$timepoint)
lf_fp_sam$obesity = as.factor(lf_fp_sam$obesity) lf_fp_sam$obesity = as.factor(lf_fp_sam$obesity)
table(lf_fp_sam$mediator)
head(lf_fp_sam$value[lf_fp_sam$mediator == "vitd"])
lf_fp_sam = lf_fp_sam[!lf_fp_sam$mediator == "vitd",]
#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 = serum_lf[serum_lf$flustat == 1,]
lf_fp_serum$timepoint = paste0("t", lf_fp_serum$timepoint) lf_fp_serum$timepoint = paste0("t", lf_fp_serum$timepoint)
@ -52,11 +62,6 @@ pdf(output_boxplot, width = 20, height = 15)
if (is.factor(lf_fp_npa$timepoint) && is.factor(lf_fp_npa$timepoint)){ if (is.factor(lf_fp_npa$timepoint) && is.factor(lf_fp_npa$timepoint)){
cat ("PASS: required groups are factors") cat ("PASS: required groups are factors")
} }
table(lf_fp_npa$mediator)
head(lf_fp_npa$value[lf_fp_npa$mediator == "vitd"])
lf_fp_npa = lf_fp_npa[!lf_fp_npa$mediator == "vitd",]
#------------------------------------------ #------------------------------------------
title_npa_linear = "Boxplot: NPA (Linear)" title_npa_linear = "Boxplot: NPA (Linear)"
#----------------------------------------- #-----------------------------------------
@ -114,8 +119,6 @@ bxp_npa_log
if (is.factor(lf_fp_sam$timepoint) && is.factor(lf_fp_sam$timepoint)){ if (is.factor(lf_fp_sam$timepoint) && is.factor(lf_fp_sam$timepoint)){
cat ("PASS: required groups are factors") cat ("PASS: required groups are factors")
} }
table(lf_fp_sam$mediator)
lf_fp_sam = lf_fp_sam[!lf_fp_sam$mediator == "vitd",]
#------------------------------------------ #------------------------------------------
title_sam_linear = "Boxplot: sam (Linear)" title_sam_linear = "Boxplot: sam (Linear)"

71
flu_stats_unpaired_npa.R
View file

@ -30,7 +30,7 @@ flu_stats_time_unpaired = paste0(outdir_stats, outfile_name)
wf = npa_wf[npa_wf$flustat == 1,] wf = npa_wf[npa_wf$flustat == 1,]
lf = npa_lf[npa_lf$flustat == 1,] lf = npa_lf[npa_lf$flustat == 1,]
lf$timepoint = paste0("t", lf$timepoint) lf$timepoint = paste0("t", lf$timepoint)
lf = lf[!lf$mediator == "vitd",]
######################################################################## ########################################################################
# clear variables # clear variables
rm(sam_lf, sam_wf rm(sam_lf, sam_wf
@ -43,12 +43,12 @@ rm(pivot_cols)
# sanity checks # sanity checks
table(lf$timepoint) table(lf$timepoint)
if (table(lf$flustat) == table(npa_lf$flustat)[[2]]){ #if (table(lf$flustat) == table(npa_lf$flustat)[[2]]){
cat("Analysing Flu positive patients for:", my_sample_type) # cat("Analysing Flu positive patients for:", my_sample_type)
}else{ #}else{
cat("FAIL: problem with subsetting data for:", my_sample_type) # cat("FAIL: problem with subsetting data for:", my_sample_type)
quit() # quit()
} #}
######################################################################## ########################################################################
# Unpaired stats at each timepoint b/w groups: wilcoxon UNpaired analysis # Unpaired stats at each timepoint b/w groups: wilcoxon UNpaired analysis
# with correction # with correction
@ -58,6 +58,7 @@ my_adjust_method = "BH"
#============== #==============
# unpaired: t1 # unpaired: t1
#==============
lf_t1 = lf[lf$timepoint == "t1",] lf_t1 = lf[lf$timepoint == "t1",]
sum(is.na(lf_t1$value)) sum(is.na(lf_t1$value))
@ -80,10 +81,20 @@ stats_un_t1$timepoint = "t1"
stats_un_t1 = as.data.frame(stats_un_t1) stats_un_t1 = as.data.frame(stats_un_t1)
class(stats_un_t1) class(stats_un_t1)
# calculate n_obs for each mediator #----------------------------------------
n_t1 = data.frame(table(lf_t1_comp$mediator)) # calculate n_obs for each mediator: t1
colnames(n_t1) = c("mediator", "n_obs") #----------------------------------------
n_t1$mediator = as.character(n_t1$mediator) #n_t1 = data.frame(table(lf_t1_comp$mediator))
n_t1_all = data.frame(table(lf_t1$mediator))
colnames(n_t1_all) = c("mediator", "n_obs")
n_t1_all$mediator = as.character(n_t1_all$mediator)
n_t1_comp = data.frame(table(lf_t1_comp$mediator))
colnames(n_t1_comp) = c("mediator", "n_obs_complete")
n_t1_comp$mediator = as.character(n_t1_comp$mediator)
merge_cols = intersect(names(n_t1_all), names(n_t1_comp)); merge_cols
n_t1= merge(n_t1_all, n_t1_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -130,10 +141,20 @@ stats_un_t2$timepoint = "t2"
stats_un_t2 = as.data.frame(stats_un_t2) stats_un_t2 = as.data.frame(stats_un_t2)
class(stats_un_t2) class(stats_un_t2)
# calculate n_obs for each mediator #----------------------------------------
n_t2 = data.frame(table(lf_t2_comp$mediator)) # calculate n_obs for each mediator: t2
colnames(n_t2) = c("mediator", "n_obs") #----------------------------------------
n_t2$mediator = as.character(n_t2$mediator) #n_t2 = data.frame(table(lf_t2_comp$mediator))
n_t2_all = data.frame(table(lf_t2$mediator))
colnames(n_t2_all) = c("mediator", "n_obs")
n_t2_all$mediator = as.character(n_t2_all$mediator)
n_t2_comp = data.frame(table(lf_t2_comp$mediator))
colnames(n_t2_comp) = c("mediator", "n_obs_complete")
n_t2_comp$mediator = as.character(n_t2_comp$mediator)
merge_cols = intersect(names(n_t2_all), names(n_t2_comp)); merge_cols
n_t2= merge(n_t2_all, n_t2_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -180,10 +201,20 @@ stats_un_t3$timepoint = "t3"
stats_un_t3 = as.data.frame(stats_un_t3) stats_un_t3 = as.data.frame(stats_un_t3)
class(stats_un_t3) class(stats_un_t3)
# calculate n_obs for each mediator #----------------------------------------
n_t3 = data.frame(table(lf_t3_comp$mediator)) # calculate n_obs for each mediator: t3
colnames(n_t3) = c("mediator", "n_obs") #----------------------------------------
n_t3$mediator = as.character(n_t3$mediator) #n_t3 = data.frame(table(lf_t3_comp$mediator))
n_t3_all = data.frame(table(lf_t3$mediator))
colnames(n_t3_all) = c("mediator", "n_obs")
n_t3_all$mediator = as.character(n_t3_all$mediator)
n_t3_comp = data.frame(table(lf_t3_comp$mediator))
colnames(n_t3_comp) = c("mediator", "n_obs_complete")
n_t3_comp$mediator = as.character(n_t3_comp$mediator)
merge_cols = intersect(names(n_t3_all), names(n_t3_comp)); merge_cols
n_t3= merge(n_t3_all, n_t3_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -291,6 +322,7 @@ my_col_order2 = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"
@ -320,6 +352,7 @@ colnames(combined_unpaired_stats_f) = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"

73
flu_stats_unpaired_sam.R
View file

@ -30,6 +30,7 @@ flu_stats_time_unpaired = paste0(outdir_stats, outfile_name)
wf = sam_wf[sam_wf$flustat == 1,] wf = sam_wf[sam_wf$flustat == 1,]
lf = sam_lf[sam_lf$flustat == 1,] lf = sam_lf[sam_lf$flustat == 1,]
lf$timepoint = paste0("t", lf$timepoint) lf$timepoint = paste0("t", lf$timepoint)
lf = lf[!lf$mediator == "vitd",]
######################################################################## ########################################################################
# clear variables # clear variables
rm(npa_lf, npa_wf rm(npa_lf, npa_wf
@ -43,12 +44,12 @@ rm(pivot_cols)
table(lf$timepoint) table(lf$timepoint)
length(unique(lf$mosaic)) length(unique(lf$mosaic))
if (table(lf$flustat) == table(sam_lf$flustat)[[2]]){ #if (table(lf$flustat) == table(sam_lf$flustat)[[2]]){
cat("Analysing Flu positive patients for:", my_sample_type) # cat("Analysing Flu positive patients for:", my_sample_type)
}else{ #}else{
cat("FAIL: problem with subsetting data for:", my_sample_type) # cat("FAIL: problem with subsetting data for:", my_sample_type)
quit() # quit()
} #}
######################################################################## ########################################################################
# Unpaired stats at each timepoint b/w groups: wilcoxon UNpaired analysis # Unpaired stats at each timepoint b/w groups: wilcoxon UNpaired analysis
@ -70,8 +71,8 @@ foo = lf_t1[which(is.na(lf_t1$value)),]
lf_t1_comp = lf_t1[-which(is.na(lf_t1$value)),] lf_t1_comp = lf_t1[-which(is.na(lf_t1$value)),]
stats_un_t1 = compare_means(value~obesity stats_un_t1 = compare_means(value~obesity
, group.by = "mediator" , group.by = "mediator"
#, data = lf_t1 , data = lf_t1
, data = lf_t1_comp #, data = lf_t1_comp
, paired = FALSE , paired = FALSE
, p.adjust.method = my_adjust_method) , p.adjust.method = my_adjust_method)
@ -82,10 +83,20 @@ stats_un_t1$timepoint = "t1"
stats_un_t1 = as.data.frame(stats_un_t1) stats_un_t1 = as.data.frame(stats_un_t1)
class(stats_un_t1) class(stats_un_t1)
# calculate n_obs for each mediator #----------------------------------------
n_t1 = data.frame(table(lf_t1_comp$mediator)) # calculate n_obs for each mediator: t1
colnames(n_t1) = c("mediator", "n_obs") #----------------------------------------
n_t1$mediator = as.character(n_t1$mediator) #n_t1 = data.frame(table(lf_t1_comp$mediator))
n_t1_all = data.frame(table(lf_t1$mediator))
colnames(n_t1_all) = c("mediator", "n_obs")
n_t1_all$mediator = as.character(n_t1_all$mediator)
n_t1_comp = data.frame(table(lf_t1_comp$mediator))
colnames(n_t1_comp) = c("mediator", "n_obs_complete")
n_t1_comp$mediator = as.character(n_t1_comp$mediator)
merge_cols = intersect(names(n_t1_all), names(n_t1_comp)); merge_cols
n_t1= merge(n_t1_all, n_t1_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -131,10 +142,20 @@ stats_un_t2$timepoint = "t2"
stats_un_t2 = as.data.frame(stats_un_t2) stats_un_t2 = as.data.frame(stats_un_t2)
class(stats_un_t2) class(stats_un_t2)
# calculate n_obs for each mediator #----------------------------------------
n_t2 = data.frame(table(lf_t2_comp$mediator)) # calculate n_obs for each mediator: t2
colnames(n_t2) = c("mediator", "n_obs") #----------------------------------------
n_t2$mediator = as.character(n_t2$mediator) #n_t2 = data.frame(table(lf_t2_comp$mediator))
n_t2_all = data.frame(table(lf_t2$mediator))
colnames(n_t2_all) = c("mediator", "n_obs")
n_t2_all$mediator = as.character(n_t2_all$mediator)
n_t2_comp = data.frame(table(lf_t2_comp$mediator))
colnames(n_t2_comp) = c("mediator", "n_obs_complete")
n_t2_comp$mediator = as.character(n_t2_comp$mediator)
merge_cols = intersect(names(n_t2_all), names(n_t2_comp)); merge_cols
n_t2= merge(n_t2_all, n_t2_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -180,10 +201,20 @@ stats_un_t3$timepoint = "t3"
stats_un_t3 = as.data.frame(stats_un_t3) stats_un_t3 = as.data.frame(stats_un_t3)
class(stats_un_t3) class(stats_un_t3)
# calculate n_obs for each mediator #----------------------------------------
n_t3 = data.frame(table(lf_t3_comp$mediator)) # calculate n_obs for each mediator: t3
colnames(n_t3) = c("mediator", "n_obs") #----------------------------------------
n_t3$mediator = as.character(n_t3$mediator) #n_t3 = data.frame(table(lf_t3_comp$mediator))
n_t3_all = data.frame(table(lf_t3$mediator))
colnames(n_t3_all) = c("mediator", "n_obs")
n_t3_all$mediator = as.character(n_t3_all$mediator)
n_t3_comp = data.frame(table(lf_t3_comp$mediator))
colnames(n_t3_comp) = c("mediator", "n_obs_complete")
n_t3_comp$mediator = as.character(n_t3_comp$mediator)
merge_cols = intersect(names(n_t3_all), names(n_t3_comp)); merge_cols
n_t3= merge(n_t3_all, n_t3_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -294,6 +325,7 @@ my_col_order2 = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"
@ -323,6 +355,7 @@ colnames(combined_unpaired_stats_f) = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"

67
flu_stats_unpaired_serum.R
View file

@ -67,8 +67,8 @@ ci = which(is.na(lf_t1$value))
lf_t1_comp = lf_t1[-which(is.na(lf_t1$value)),] lf_t1_comp = lf_t1[-which(is.na(lf_t1$value)),]
stats_un_t1 = compare_means(value~obesity stats_un_t1 = compare_means(value~obesity
, group.by = "mediator" , group.by = "mediator"
#, data = lf_t1 , data = lf_t1
, data = lf_t1_comp #, data = lf_t1_comp
, paired = FALSE , paired = FALSE
, p.adjust.method = my_adjust_method) , p.adjust.method = my_adjust_method)
@ -79,10 +79,20 @@ stats_un_t1$timepoint = "t1"
stats_un_t1 = as.data.frame(stats_un_t1) stats_un_t1 = as.data.frame(stats_un_t1)
class(stats_un_t1) class(stats_un_t1)
# calculate n_obs for each mediator #----------------------------------------
n_t1 = data.frame(table(lf_t1_comp$mediator)) # calculate n_obs for each mediator: t1
colnames(n_t1) = c("mediator", "n_obs") #----------------------------------------
n_t1$mediator = as.character(n_t1$mediator) #n_t1 = data.frame(table(lf_t1_comp$mediator))
n_t1_all = data.frame(table(lf_t1$mediator))
colnames(n_t1_all) = c("mediator", "n_obs")
n_t1_all$mediator = as.character(n_t1_all$mediator)
n_t1_comp = data.frame(table(lf_t1_comp$mediator))
colnames(n_t1_comp) = c("mediator", "n_obs_complete")
n_t1_comp$mediator = as.character(n_t1_comp$mediator)
merge_cols = intersect(names(n_t1_all), names(n_t1_comp)); merge_cols
n_t1= merge(n_t1_all, n_t1_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -119,8 +129,8 @@ lf_t2_comp = lf_t2[-which(is.na(lf_t2$value)),]
stats_un_t2 = compare_means(value~obesity stats_un_t2 = compare_means(value~obesity
, group.by = "mediator" , group.by = "mediator"
#, data = lf_t2 , data = lf_t2
, data = lf_t2_comp #, data = lf_t2_comp
, paired = FALSE , paired = FALSE
, p.adjust.method = my_adjust_method) , p.adjust.method = my_adjust_method)
# add timepoint and convert to df # add timepoint and convert to df
@ -128,10 +138,20 @@ stats_un_t2$timepoint = "t2"
stats_un_t2 = as.data.frame(stats_un_t2) stats_un_t2 = as.data.frame(stats_un_t2)
class(stats_un_t2) class(stats_un_t2)
# calculate n_obs for each mediator #----------------------------------------
n_t2 = data.frame(table(lf_t2_comp$mediator)) # calculate n_obs for each mediator: t2
colnames(n_t2) = c("mediator", "n_obs") #----------------------------------------
n_t2$mediator = as.character(n_t2$mediator) #n_t2 = data.frame(table(lf_t2_comp$mediator))
n_t2_all = data.frame(table(lf_t2$mediator))
colnames(n_t2_all) = c("mediator", "n_obs")
n_t2_all$mediator = as.character(n_t2_all$mediator)
n_t2_comp = data.frame(table(lf_t2_comp$mediator))
colnames(n_t2_comp) = c("mediator", "n_obs_complete")
n_t2_comp$mediator = as.character(n_t2_comp$mediator)
merge_cols = intersect(names(n_t2_all), names(n_t2_comp)); merge_cols
n_t2= merge(n_t2_all, n_t2_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -168,8 +188,8 @@ lf_t3_comp = lf_t3[-which(is.na(lf_t3$value)),]
stats_un_t3 = compare_means(value~obesity stats_un_t3 = compare_means(value~obesity
, group.by = "mediator" , group.by = "mediator"
#, data = lf_t3 , data = lf_t3
, data = lf_t3_comp #, data = lf_t3_comp
, paired = FALSE , paired = FALSE
, p.adjust.method = my_adjust_method) , p.adjust.method = my_adjust_method)
@ -178,11 +198,20 @@ stats_un_t3$timepoint = "t3"
stats_un_t3 = as.data.frame(stats_un_t3) stats_un_t3 = as.data.frame(stats_un_t3)
class(stats_un_t3) class(stats_un_t3)
#----------------------------------------
# calculate n_obs for each mediator: t3
#----------------------------------------
#n_t3 = data.frame(table(lf_t3_comp$mediator))
n_t3_all = data.frame(table(lf_t3$mediator))
colnames(n_t3_all) = c("mediator", "n_obs")
n_t3_all$mediator = as.character(n_t3_all$mediator)
# calculate n_obs for each mediator n_t3_comp = data.frame(table(lf_t3_comp$mediator))
n_t3 = data.frame(table(lf_t3_comp$mediator)) colnames(n_t3_comp) = c("mediator", "n_obs_complete")
colnames(n_t3) = c("mediator", "n_obs") n_t3_comp$mediator = as.character(n_t3_comp$mediator)
n_t3$mediator = as.character(n_t3$mediator)
merge_cols = intersect(names(n_t3_all), names(n_t3_comp)); merge_cols
n_t3 = merge(n_t3_all, n_t3_comp, by = merge_cols, all = T)
#================================== #==================================
# Merge: merge stats + n_obs df # Merge: merge stats + n_obs df
@ -289,6 +318,7 @@ my_col_order2 = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"
@ -318,6 +348,7 @@ colnames(combined_unpaired_stats_f) = c("mediator"
, "timepoint" , "timepoint"
, "sample_type" , "sample_type"
, "n_obs" , "n_obs"
, "n_obs_complete"
, "group1" , "group1"
, "group2" , "group2"
, "method" , "method"

2
summary_stats_time_group_npa.R
View file

@ -31,7 +31,7 @@ wf = npa_wf[npa_wf$flustat == 1,]
lf = npa_lf[npa_lf$flustat == 1,] lf = npa_lf[npa_lf$flustat == 1,]
lf$timepoint = paste0("t", lf$timepoint) lf$timepoint = paste0("t", lf$timepoint)
lf = lf[!lf$mediator == "vitd",] lf = lf[!lf$mediator == "vitd",]
########################################################################
#======================================================= #=======================================================
# summary stats by timepoint and obesity: each mediator # summary stats by timepoint and obesity: each mediator