O.2 Report - Review

Focus on OXT by Gender

1 Read and Merge

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Read and Merge
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
wd <- "E:/Cinetic idei noi/EXPERIMENTE OGL Frontiers (O.2 & O.0.3 & O.0.2)"
setwd(wd)

DataTrust <- rio::import(file.path(wd, "O.2 REZULTATE\\O.2 Date PrelucrareSPSS cu NEO si STAI", "O.2 Date pt Trust BUN cu NEO si STAIY.xlsx"))
DataDG <- rio::import(file.path(wd, "O.2 REZULTATE\\O.2 Date PrelucrareSPSS cu NEO si STAI", "O.2 Date pt DG BUN cu NEO si STAIY.xlsx"))
DataVAS <- rio::import(file.path(wd, "O.2 REZULTATE\\O.2 VAS,IOS", "O.2 Date PrelucrareSPSS.xlsx"))
DataBIO <- rio::import(file.path(wd, "O.2 BIO", "O.2 Ox si Cortizol.xlsx"))
  
Data_merge1 <- merge(DataVAS, DataTrust)  
Data_merge2 <- merge(Data_merge1, DataDG) 
Data_merge3 <- merge(Data_merge2, DataBIO)

Data <- Data_merge3

test_names <- unique(unlist(lapply(list(DataTrust, DataDG, DataVAS, DataBIO), names)))
merge_names <- names(Data)

if(identical(merge_names[order(merge_names)], test_names[order(test_names)])){    # the order matters in identical()
  cat("**Merge was succesful**")
  rm("Data_merge1", "Data_merge2", "Data_merge3", "DataBIO", "DataDG", "DataTrust", "DataVAS", "test_names", "merge_names")
}else cat("**Merge unsuccesful**") 

Merge was succesful

# Gender Dataframe
Data_Gen <- rio::import(file.path(wd, "Gen varsta O03 O02 O2.xlsx"), which = "O.2")

Data_Gen_merged <- 
  Data %>%
  tidyr::separate(Indicativ,  c("ID_tag", "ID", "study_tag"), "\\s+") %>%                               # split on white space
  select(-c("ID_tag", "study_tag")) %>%
  mutate(ID = as.numeric(as.character(ID))) %>%
  dplyr::left_join(., Data_Gen, by = c("ID")) %>%
  select(1:7, Gen, Varsta, everything())

Data_Gen_merged_long <-
  Data_Gen_merged %>%
  dplyr::select(ID, Exercitiul, Gen, starts_with("Vas"), starts_with("Ox")) %>%
  dplyr::rename(VasS_pre_Poz = VasS_prePoz,
                VasB_pre_Poz = VasB_prePoz,
                VasS_post_Poz = VasS_postPoz,
                VasB_post_Poz = VasB_postPoz,
                VasS_pre_Neg = VasS_preNeg,
                VasB_pre_Neg = VasB_preNeg,
                VasS_post_Neg = VasS_postNeg,
                VasB_post_Neg = VasB_postNeg) %>%
  pivot_longer(-c(ID, Exercitiul, Gen), names_sep = "_", names_to = c("Vars", "PrePost", "Cond"), values_to = "Vals") %>%
  rstatix::convert_as_factor(ID, Exercitiul, Gen, Vars, PrePost, Cond)



Data_fem <- 
  Data_Gen_merged_long %>%
  filter(Gen == "F")

Data_masc <- 
  Data_Gen_merged_long %>%
  filter(Gen == "M")

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Derive new variables
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Data$D_VasS_Poz <- Data[, "VasS_postPoz"] - Data[, "VasS_prePoz"] 
Data$D_VasS_Neg <- Data[, "VasS_postNeg"] - Data[, "VasS_preNeg"]
Data$D_VasB_Poz <- Data[, "VasB_postPoz"] - Data[, "VasB_prePoz"] 
Data$D_VasB_Neg <- Data[, "VasB_postNeg"] - Data[, "VasB_preNeg"]
Data$D_IOS_Poz <- Data[, "IOS_postPoz"] - Data[, "IOS_prePoz"] 
Data$D_IOS_Neg <- Data[, "IOS_postNeg"] - Data[, "IOS_preNeg"]

Data$D_Sam1_Poz <- Data[, "Sam1_postPoz"] - Data[, "Sam1_prePoz"] 
Data$D_Sam1_Neg <- Data[, "Sam1_postNeg"] - Data[, "Sam1_preNeg"]
Data$D_Sam2_Poz <- Data[, "Sam2_postPoz"] - Data[, "Sam2_prePoz"] 
Data$D_Sam2_Neg <- Data[, "Sam2_postNeg"] - Data[, "Sam2_preNeg"]
Data$D_Sam3_Poz <- Data[, "Sam3_postPoz"] - Data[, "Sam3_prePoz"] 
Data$D_Sam3_Neg <- Data[, "Sam3_postNeg"] - Data[, "Sam3_preNeg"]

Data$D_DG_Poz <- Data[, "DG_postPozTot"] - Data[, "DG_prePozTot"] 
Data$D_DG_Neg <- Data[, "DG_postNegTot"] - Data[, "DG_preNegTot"]

Data$D_TrustMin_Poz <- Data[, "TrustMinPozPost"] - Data[, "TrustMinPozPre"] 
Data$D_TrustMin_Neg <- Data[, "TrustMinNegPost"] - Data[, "TrustMinNegPre"]
Data$D_TrustTot_Poz <- Data[, "TrustTotPozPost"] - Data[, "TrustTotPozPre"] 
Data$D_TrustTot_Neg <- Data[, "TrustTotNegPost"] - Data[, "TrustTotNegPre"]

Data$D_Cort_Poz <- Data[, "Cort_post_Poz"] - Data[, "Cort_pre_Poz"] 
Data$D_Cort_Neg <- Data[, "Cort_post_Neg"] - Data[, "Cort_pre_Neg"]
Data$D_Ox_Poz <- Data[, "Ox_post_Poz"] - Data[, "Ox_pre_Poz"] 
Data$D_Ox_Neg <- Data[, "Ox_post_Neg"] - Data[, "Ox_pre_Neg"]

–>

1.1 Dataset - Long Format

Data_Gen_merged_long %>%
    DT::datatable(
      extensions = 'Buttons',
      options = list(pageLength = 10,
                     scrollX='500px',
                     dom = 'Bfrtip',
                     buttons = c('excel', "csv")))

1.2 Dataset - Wide Format

Data_Gen_merged %>%
  dplyr::select(-2) %>%
    DT::datatable(
      extensions = 'Buttons',
      options = list(pageLength = 10,
                     scrollX='500px',
                     dom = 'Bfrtip',
                     buttons = c('excel', "csv")))

2 Define functions

# Define Function for Two-way rmANOVA
# library(tidyverse)
# library(ggpubr)
# library(rstatix)
# library(rlang)

tw_rmANOVA_func <- 
  function(data, id_var, cond_var, time_var, value_var, 
           assum_check = TRUE, posthoc_sig_interac = FALSE, posthoc_ns_interac = FALSE,
           p_adjust_method = "bonferroni"){
  
  # input dataframe needs to have columns names diffrent from "variable" and "value" because it collides with rstatix::shapiro_test
    
  id_var_enq <- rlang::enquo(id_var)  
  cond_var_enq <- rlang::enquo(cond_var)
  cond_var_name <- rlang::as_name(cond_var_enq)
  time_var_enq <- rlang::enquo(time_var)
  time_var_name <- rlang::as_name(time_var_enq)
  value_var_enq <- rlang::enquo(value_var)
  value_var_name <- rlang::as_name(value_var_enq)
    
  # Assumptions
  if(assum_check){
    cat("\n Outliers \n")
    data %>%
      dplyr::group_by(!!cond_var_enq, !!time_var_enq) %>%
      rstatix::identify_outliers(!!value_var_enq) %>%                                  # outliers (needs to be 0)
      print()
    
    cat("\n Normality assumption (p>.05) \n")
    data %>%
      dplyr::group_by(!!cond_var_enq, !!time_var_enq) %>%
      rstatix::shapiro_test(!!value_var_enq) %>%                                        # normality assumption (p>.05)
      print()
    
    qq_plot <- 
      ggpubr::ggqqplot(data = data, value_var_name, ggtheme = theme_bw(), title = "QQ Plot") +
        ggplot2::facet_grid(vars(!!time_var_enq), vars(!!cond_var_enq), labeller = "label_both")    # QQ plot
  }
  
  # Two-way rmANOVA - check for interaction (ex. F(2, 22) = 30.4, p < 0.0001)
  cat("\n Two-way rmANOVA \n")
  res_aov <- anova_test(                         # automatically does sphericity Mauchly’s test
    data = data, dv = !!value_var_enq, wid = !!id_var_enq,                             
    within = c(!!cond_var_enq, !!time_var_enq)
  )
  get_anova_table(res_aov) %>%  # ges: Greenhouse-Geisser sphericity correction is automatically applied to factors violating the sphericity assumption  
  print()
  
  #- Procedure for a significant two-way interaction -
  if(posthoc_sig_interac){
    cat("\n Effect of treatment at each time point \n")
    one_way <- 
      data %>%
      group_by(!!time_var_enq) %>%
      anova_test(dv = !!value_var_enq, wid = !!id_var_enq, within = !!cond_var_enq) %>%
      get_anova_table() %>%
      adjust_pvalue(method = "bonferroni")
    one_way %>% print()
    
    cat("\n Pairwise comparisons between treatment groups \n")
    pwc <-
      data %>%
      group_by(!!time_var_enq) %>%
      pairwise_t_test(
        as.formula(paste0(value_var_name, " ~ ", cond_var_name)),
        paired = TRUE,
        p.adjust.method = p_adjust_method
      )
    pwc %>% print()

    cat("\n Effect of time at each level of treatment - One-way ANOVA \n")
    one_way2 <-
      data %>%
      group_by(!!cond_var_enq) %>%
      anova_test(dv = !!value_var_enq, wid = !!id_var_enq, within = !!time_var_enq) %>%
      get_anova_table() %>%
      adjust_pvalue(method = p_adjust_method)
    one_way2  %>% print()
    cat("\n Pairwise comparisons between time points \n")
    pwc2 <-
      data %>%
      group_by(!!cond_var_enq) %>%
      pairwise_t_test(
        as.formula(paste0(value_var_name, " ~ ", time_var_name)),     # paste formula, not quosure
        paired = TRUE,
        p.adjust.method = p_adjust_method
      )
    pwc2  %>% print()
  }
  
  #- Procedure for non-significant two-way interaction- 
  # If the interaction is not significant, you need to interpret the main effects for each of the two variables: treatment and time.
  if(posthoc_ns_interac){
    cat("\n Comparisons for treatment variable \n")
    pwc_cond <-
      data %>%
      pairwise_t_test(
        as.formula(paste0(value_var_name, " ~ ", cond_var_name)),     # paste formula, not quosure             
        paired = TRUE,
        p.adjust.method = p_adjust_method
      )
    pwc_cond %>% print()
    cat("\n Comparisons for time variable \n")
    pwc_time <-
      data %>% 
      pairwise_t_test(
        as.formula(paste0(value_var_name, " ~ ", time_var_name)),     # paste formula, not quosure
        paired = TRUE,
        p.adjust.method = p_adjust_method
      )
    pwc_time %>% print()
  }
  
  # Visualization
  bx_plot <- 
    ggboxplot(data, x = time_var_name, y = value_var_name,
    color = cond_var_name, palette = "jco")
  pwc <- 
    pwc %>% 
    add_xy_position(x = time_var_name)
  bx_plot <- 
    bx_plot + 
    stat_pvalue_manual(pwc, tip.length = 0, hide.ns = TRUE) +
    labs(
      subtitle = get_test_label(res_aov, detailed = TRUE),
      caption = get_pwc_label(pwc)
    )
  
  if(assum_check){ 
    list(qq_plot, bx_plot)
  }else{
    bx_plot
  } 
  
}

# ex. - run on long format
# tw_rmANOVA_func(data = selfesteem2, id_var = id, cond_var = treatment, time_var = time, value_var = scores, 
#                 posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

## Wilcoxon signed rank test on paired samples

# library(tidyverse)
# library(ggpubr)
# library(rstatix)
# library(rlang)

wilcoxon_paired_func <- function(data, id_var, time_var, value_var, 
                                 assum_check = TRUE){
  
  id_var_enq <- rlang::enquo(id_var)  
  id_var_name <- rlang::as_name(id_var_enq)
  time_var_enq <- rlang::enquo(time_var)
  time_var_name <- rlang::as_name(time_var_enq)
  value_var_enq <- rlang::enquo(value_var)
  value_var_name <- rlang::as_name(value_var_enq)
  
  cat("\n Summary statistics \n")
  data %>%
    group_by(!!time_var_enq) %>%
    get_summary_stats(!!value_var_enq, type = c("common")) %>%
    print() 
  
  if(assum_check){  # The test assumes that differences between paired samples should be distributed symmetrically around the median.
  data2 <- 
    data %>% 
    group_by(!!id_var_enq) %>% 
    mutate(differences = !!value_var_enq - lag(!!value_var_enq))

  hist_dif_plot <- gghistogram(data2, x = "differences", y = "..density..", 
                               fill = "steelblue", add_density = TRUE)
  }
  
  cat("\n Wilcoxon signed rank test on paired samples \n")
  stat_test <- 
    data  %>%
    wilcox_test(as.formula(paste0(value_var_name, " ~ ", time_var_name)), 
                paired = TRUE) %>%
    add_significance()
  stat_test %>% print()
  
  cat("\n Effect size \n")
  ef_size <-
    data %>%
    wilcox_effsize(as.formula(paste0(value_var_name, " ~ ", time_var_name)), 
                   paired = TRUE)
  ef_size %>% print()
  
  # Visualization
  bx_plot <- ggpaired(data, x = time_var_name, y = value_var_name, id = id_var_name,
                      ylab = value_var_name, xlab = time_var_name)
  
  stat_test <- 
    stat_test %>% 
    add_xy_position(x = time_var_name)
  
  bx_plot <-
    bx_plot + 
      stat_pvalue_manual(stat_test, tip.length = 0) +
      labs(subtitle = get_test_label(stat_test, detailed = TRUE))
  
  if(assum_check){ 
    list(hist_dif_plot, bx_plot)
  }else{
    bx_plot
  }  
  
}

3 rmANOVA OXT

cat("## All subjects")

## All subjects

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = OXT, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd      F        p p<.05   ges
1         Cond   1  28  0.537 0.470000       0.008
2      PrePost   1  28 17.388 0.000266     * 0.063
3 Cond:PrePost   1  28  0.220 0.642000       0.001

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Females")

## Females

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = OXT, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd      F     p p<.05   ges
1         Cond   1  20  1.054 0.317       0.022
2      PrePost   1  20 10.660 0.004     * 0.055
3 Cond:PrePost   1  20  0.188 0.669       0.001

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Males")

## Males

Data_masc %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = OXT, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05      ges
1         Cond   1   7 0.031 0.865       0.002000
2      PrePost   1   7 6.360 0.040     * 0.087000
3 Cond:PrePost   1   7 0.031 0.866       0.000557

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

4 rmANOVA Stres

cat("## All subjects")

## All subjects

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasS, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05       ges
1         Cond   1  29 0.004 0.949       0.0000565
2      PrePost   1  29 3.940 0.057       0.0060000
3 Cond:PrePost   1  29 4.259 0.048     * 0.0090000

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Females")

## Females

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasS, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05      ges
1         Cond   1  21 0.014 0.907       0.000288
2      PrePost   1  21 1.407 0.249       0.003000
3 Cond:PrePost   1  21 5.206 0.033     * 0.017000

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Males")

## Males

Data_masc %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasS, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05      ges
1         Cond   1   7 0.198 0.670       0.010000
2      PrePost   1   7 4.505 0.071       0.028000
3 Cond:PrePost   1   7 0.022 0.885       0.000121

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

5 rmANOVA WB

cat("## All subjects")

## All subjects

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasB, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05   ges
1         Cond   1  29 1.566 0.221       0.008
2      PrePost   1  29 2.072 0.161       0.003
3 Cond:PrePost   1  29 5.086 0.032     * 0.019

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Females")

## Females

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasB, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05   ges
1         Cond   1  21 0.607 0.445       0.005
2      PrePost   1  21 0.850 0.367       0.003
3 Cond:PrePost   1  21 4.814 0.040     * 0.033

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

cat("## Males")

## Males

Data_masc %>%
  select(ID, Vars, PrePost, Cond, Vals) %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%                     # for complete cases
  ungroup() %>%
  tw_rmANOVA_func(id_var = ID, cond_var = Cond, time_var = PrePost, value_var = VasB, 
                  posthoc_sig_interac = TRUE, posthoc_ns_interac = TRUE)

 Outliers 

 Normality assumption (p>.05) 

 Two-way rmANOVA 
ANOVA Table (type III tests)

        Effect DFn DFd     F     p p<.05   ges
1         Cond   1   7 1.735 0.229       0.021
2      PrePost   1   7 5.885 0.046     * 0.007
3 Cond:PrePost   1   7 0.414 0.540       0.001

 Effect of treatment at each time point 

 Pairwise comparisons between treatment groups 

 Effect of time at each level of treatment - One-way ANOVA 

 Pairwise comparisons between time points 

 Comparisons for treatment variable 

 Comparisons for time variable 

[[1]]

[[2]]

6 Wilcoxon paired on Order of Conditions - OXT

cat("## Whole Sample - Poz first, Neg second")

## Whole Sample - Poz first, Neg second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = OXT, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Whole Sample - Neg first, Poz second")

## Whole Sample - Neg first, Poz second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = OXT, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Poz first, Neg second")

## Females - Poz first, Neg second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = OXT, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Neg first, Poz second")

## Females - Neg first, Poz second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "Ox") %>%
  dplyr::rename(OXT = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = OXT, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

7 Wilcoxon paired on Order of Conditions - Stres

cat("## Whole Sample - Poz first, Neg second")

## Whole Sample - Poz first, Neg second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasS, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Whole Sample - Neg first, Poz second")

## Whole Sample - Neg first, Poz second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasS, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Poz first, Neg second")

## Females - Poz first, Neg second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasS, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Neg first, Poz second")

## Females - Neg first, Poz second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "VasS") %>%
  dplyr::rename(VasS = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasS, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

8 Wilcoxon paired on Order of Conditions - WB

cat("## Whole Sample - Poz first, Neg second")

## Whole Sample - Poz first, Neg second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasB, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Whole Sample - Neg first, Poz second")

## Whole Sample - Neg first, Poz second

Data_Gen_merged_long %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasB, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Poz first, Neg second")

## Females - Poz first, Neg second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Pozitiv - Negativ") %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasB, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

cat("## Females - Neg first, Poz second")

## Females - Neg first, Poz second

Data_fem %>%
  select(ID, Vars, PrePost, Cond, Vals, Exercitiul) %>%
  filter(Exercitiul == "Negativ - Pozitiv") %>%
  filter(Vars == "VasB") %>%
  dplyr::rename(VasB = Vals) %>%
  drop_na() %>%
  group_by(ID) %>% 
  filter(n() > 3) %>%
  ungroup() %>%
  wilcoxon_paired_func(id_var = ID, time_var = PrePost, value_var = VasB, assum_check = TRUE)

 Summary statistics 

 Wilcoxon signed rank test on paired samples 

 Effect size 

[[1]]

[[2]]

9 Session Info

R version 3.6.1 (2019-07-05)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 8.1 x64 (build 9600)

Matrix products: default

locale:
[1] LC_COLLATE=Romanian_Romania.1250  LC_CTYPE=Romanian_Romania.1250    LC_MONETARY=Romanian_Romania.1250 LC_NUMERIC=C                     
[5] LC_TIME=Romanian_Romania.1250    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] rlang_0.4.6        rstatix_0.5.0      rio_0.5.16         summarytools_0.8.8 DT_0.5             ggpubr_0.2.5       magrittr_1.5       broom_0.5.6       
 [9] papaja_0.1.0.9842  psych_1.9.12.31    forcats_0.5.0      stringr_1.4.0      dplyr_0.8.5        purrr_0.3.3        readr_1.3.1        tidyr_1.0.2       
[17] tibble_3.0.0       ggplot2_3.3.0      tidyverse_1.3.0    pacman_0.5.1      

loaded via a namespace (and not attached):
 [1] nlme_3.1-140       bitops_1.0-6       matrixStats_0.54.0 fs_1.4.1           lubridate_1.7.4    httr_1.4.1         ggsci_2.9          tools_3.6.1       
 [9] backports_1.1.6    R6_2.4.1           DBI_1.0.0          colorspace_1.4-1   withr_2.1.2        tidyselect_1.0.0   mnormt_1.5-6       curl_4.3          
[17] compiler_3.6.1     cli_2.0.2          rvest_0.3.5        xml2_1.3.1         sandwich_2.5-0     labeling_0.3       scales_1.1.0       mvtnorm_1.1-0     
[25] digest_0.6.25      foreign_0.8-71     pkgconfig_2.0.3    htmltools_0.4.0    dbplyr_1.4.3       htmlwidgets_1.5.1  readxl_1.3.1       rstudioapi_0.11   
[33] pryr_0.1.4         shiny_1.2.0        farver_2.0.3       generics_0.0.2     zoo_1.8-4          jsonlite_1.6.1     crosstalk_1.0.0    zip_1.0.0         
[41] car_3.0-7          RCurl_1.95-4.11    modeltools_0.2-22  rapportools_1.0    Matrix_1.2-17      Rcpp_1.0.4.6       munsell_0.5.0      fansi_0.4.1       
[49] abind_1.4-5        lifecycle_0.2.0    multcomp_1.4-8     stringi_1.4.6      yaml_2.2.1         carData_3.0-2      MASS_7.3-51.4      plyr_1.8.6        
[57] grid_3.6.1         parallel_3.6.1     promises_1.0.1     crayon_1.3.4       lattice_0.20-38    splines_3.6.1      haven_2.2.0        pander_0.6.3      
[65] hms_0.5.3          knitr_1.28         pillar_1.4.3       ggsignif_0.4.0     stats4_3.6.1       codetools_0.2-16   reprex_0.3.0       glue_1.4.0        
[73] data.table_1.12.8  modelr_0.1.6       vctrs_0.2.4        httpuv_1.4.5       cellranger_1.1.0   gtable_0.3.0       assertthat_0.2.1   xfun_0.13         
[81] openxlsx_4.1.0     coin_1.2-2         mime_0.9           xtable_1.8-4       later_0.7.5        survival_2.44-1.1  TH.data_1.0-9      ellipsis_0.3.0    

 

---

A work by Claudiu Papasteri

claudiu.papasteri@gmail.com