alunats
alunats

Reputation: 89

Problem with reference to variable in R function

I'd like to write a function to do ANOVA in batch, but there's a problem I can't solve. The problem is with the variable reference. The whole code is written correctly, because everything is calculated "on foot", but when I try to insert this code into the function, it doesn't work. Can someone take a look and point out where the problem is?

library(tidyverse)
library(ggpubr)
library(rstatix)
library(ggprism)

df <- data.frame(
  stringsAsFactors = FALSE,
  id = c(1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5),
  treat = c("o","o","o","o","o","j","j","j","j","j","z","z","z","z","z","w","w","w","w","w"),
  vo2 = c("47.48","42.74","45.23","51.65","49.11","51.00","43.82","49.88","54.61","52.20","51.31",
          "47.56","50.69","54.88","55.01","51.89","46.10","50.98","53.62","52.77"))

df$vo2 <- as.numeric(df$vo2)
df$treat <- factor(df$treat)

Everything works fine in the code below...

# Summary
group_by(df, treat) %>% 
  summarise(
    N = n(),
    Mean = mean(vo2, na.rm = TRUE),
    Sd = sd(vo2, na.rm = TRUE))

# ANOVA
res.aov <- anova_test(dv = vo2, wid = id, within = treat, data = df)
get_anova_table(res.aov, correction = c("auto"))

# Pairwise comparisons
pwc <- df %>%
  pairwise_t_test(vo2 ~ treat, paired = TRUE, conf.level = 0.95,
                  detailed = TRUE, p.adjust.method = "bonferroni")
pwc

pwc <- pwc %>% add_xy_position(x = "treat")

ggplot(df, aes(x = treat, y = vo2)) +
  stat_boxplot(aes(x = treat, y = vo2, color = treat), geom = 'errorbar', coef=1.5, width=0.4, linetype = 1) +
  geom_boxplot(aes(x = treat, y = vo2, color = treat, fill = treat)) +
  geom_jitter(aes(x = treat, y = vo2, color = treat, fill = treat), width = 0.2) +
  stat_summary(fun = mean, geom = "point", shape = 0, size = 2, color = "black", stroke = 1) +
  #xlab(deparse(substitute(x))) +   ylab(deparse(substitute(y))) +
  theme_prism(base_size = 14) + scale_x_discrete(guide = "prism_bracket") +
  scale_fill_prism(palette = "floral") +   scale_colour_prism(palette = "floral") +
  scale_y_continuous(expand = expansion(mult = c(0.02, 0.05))) +
  stat_pvalue_manual(pwc, tip.length = 0, hide.ns = TRUE) +
  labs(
    subtitle = get_test_label(res.aov, detailed = TRUE),
    caption = get_pwc_label(pwc)) +
  theme(legend.position = "NULL")


The result of this code is the chart below enter image description here

But the function doesn't work... The fix for "deparse(substitute(x)" only partially solved the problem, and I get a warning:

1: In mean.default(~"vo2", na.rm = TRUE): argument is not numeric or logical: returning an
2: In var(if (is.vector(x) || is.factor(x)) x else as.double(x), na.rm = na.rm) : NAs introduced by coercion.

Below is the full code for this function:

my_function <- function(df, x, y) { 
  
  x <- deparse(substitute(x))
  y <- deparse(substitute(y))
  
  formula <- as.formula(paste0(y, "~", x))
  
  # Summary
a <- group_by(df, {{x}}) %>% 
    summarise(
      N = n(),
      Mean = mean({{y}}, na.rm = TRUE),
      Sd = sd({{y}}, na.rm = TRUE)
    )
  # ANOVA
  res.aov <<- anova_test(dv = {{y}}, wid = id, within = {{x}}, data = df)
b <- get_anova_table(res.aov, correction = c("auto"))
  
  # Pairwise comparisons
pwc <-  pairwise_t_test(data = df, formula, paired = TRUE, conf.level = 0.95,
                          detailed = TRUE, p.adjust.method = "bonferroni")
  
pwc2 <<- pwc %>% add_xy_position(x = {{x}})
  
  # Plot
d <- ggplot(df, aes(x = {{x}}, y = {{y}})) +
       stat_boxplot(aes(x = {{x}}, y = {{y}}, color = {{x}}), geom = 'errorbar', coef=1.5, width=0.4, linetype = 1) +
       geom_boxplot(aes(x = {{x}}, y = {{y}}, color = {{x}}, fill = {{x}})) +
       geom_jitter(aes(x = {{x}}, y = {{y}}, color = {{x}}, fill = {{x}}), width = 0.2) +
    stat_summary(fun = mean, geom = "point", shape = 0, size = 2, color = "black", stroke = 1) +
    stat_pvalue_manual(pwc2, tip.length = 0, hide.ns = TRUE) +
      xlab(deparse(substitute(x))) + ylab(deparse(substitute(y))) +
      scale_x_discrete(guide = "prism_bracket") +
    theme_prism(base_size = 14) +
    theme(legend.position = "NULL")
  
  #ggsave(paste0(deparse(substitute(x)), "_",
  # deparse(substitute(y)), ".png"), width=160, height=90, units="mm", dpi=600)
  
   output <- list(a,b,pwc2,d)
  return(output)
}

my_function(df, treat, vo2)

I have a huge request for tips on how to solve this.

Upvotes: 0

Views: 165

Answers (1)

Edward
Edward

Reputation: 18493

The problem lies in the formula. Try adding the following code:

  x1 <- deparse(substitute(x))
  y1 <- deparse(substitute(y))
  
  formula <- as.formula(paste0(y1, "~", x1))

  # Pairwise comparisons
  pwc <- pairwise_t_test(data=df, formula, paired = TRUE, conf.level = 0.95,
                    detailed = TRUE, p.adjust.method = "bonferroni")

Upvotes: 1

Related Questions