R/bartest.R

Defines functions bartest

Documented in bartest

bartest <- function(measure, group, alternative = '', legend = TRUE) {
  require(ggplot2)
  require(compute.es)
  require(lsr)

  legend = legend

  group <- droplevels(as.factor(group))

  if(nlevels(as.factor(group)) == 2) {
  group <- as.factor(group)
  gg <- levels(as.factor(group))
  lev <- (as.numeric(as.factor(group)))
  sem1 = sd(na.exclude(measure[lev == 1]))/sqrt(length(na.exclude(measure[lev == 1])))
  mean1 = mean(na.exclude(measure[lev == 1]))
  sem2 = sd(na.exclude(measure[lev == 2]))/sqrt(length(na.exclude(measure[lev == 2])))
  mean2 = mean(na.exclude(measure[lev == 2]))

  df <- data.frame(
    group   = c(gg[[1]], gg[[2]]),
    measure = c(mean1, mean2),
    sem     = c(sem1, sem2)
  )


  Ns <- tapply(measure,group,length)
  means <- tapply(measure,group,mean, na.rm = T)
  sds <- tapply(measure,group, sd, na.rm = T)

  var.res <- var.test(measure ~ group)

  if(var.res$p.value < .05) {
    var.equal = F
  } else {var.equal = T}

  if(missing(alternative)) {
    t.res <- t.test(measure ~ group, var.equal = var.equal)
  }
  else {
    t.res <- t.test(measure ~ group, var.equal = var.equal, alternative = alternative)
  }

  am <- means[[1]]; bm <- means[[2]]
  asd <- sds[[1]]; bsd <- sds[[2]]
  aN <- Ns[[1]]; bN <- Ns[[2]]

  es <- compute.es::tes(t.res$statistic, aN, bN, verbose = F)

  space <- ''
  a <- paste(as.character(levels(group)[1]),' ',
             '(M=',round(am,2),', SD=', round(asd,2),')',
             ' - N=',aN,
             sep = '')
  b <- paste(as.character(levels(group)[2]),' ',
             '(M=',round(bm,2),', SD=', round(bsd,2),')',
             ' - N=',bN,
             sep = '')
  t <- paste('t-test: ','t(', round(t.res$parameter, 2),')=',round(t.res$statistic,2),', p=', round(t.res$p.value,4), sep = '')
  d <- paste('d=', es$d, sep = '')
  v <- paste('var: ','F(', var.res$parameter[1],',',var.res$parameter[2],')=',round(var.res$statistic,2),', p=', round(var.res$p.value,4), sep = '')

  text <- paste(space,a,b,v,t,d,space, sep = '\n')

  plot <- ggplot(df, aes(x=group, y=measure, fill=group)) +
    geom_bar(position=position_dodge(), stat="identity", fill=c("white", "Grey80"),
             colour="#000000", size=.4) +
    geom_errorbar(aes(ymin=measure-sem, ymax=measure+sem), width=.2, position=position_dodge(.8)) +
    guides(fill=F) +
    theme(panel.background = element_rect(fill = "white")) +
    labs(x = '', y = '', title = paste(substitute(measure))[[3]]) +
    {if(legend)geom_label(aes(x = 1.5, y = measure[[1]]/2, label = text), fill = "white")} +
    theme(plot.title = element_text(hjust = 0.5))

  print(paste(substitute(measure))[[3]])
  print(a)
  print(b)
  print(v)
  print(t)
  print(d)
  return(plot)
  }

  ##############################################################################.
  ##############################################################################.

  else {
    group <- as.factor(group)
    gg <- levels(as.factor(group))
    lev <- (as.numeric(as.factor(group)))
    sem1 = sd(na.exclude(measure[lev == 1]))/sqrt(length(na.exclude(measure[lev == 1])))
    mean1 = mean(na.exclude(measure[lev == 1]))
    sem2 = sd(na.exclude(measure[lev == 2]))/sqrt(length(na.exclude(measure[lev == 2])))
    mean2 = mean(na.exclude(measure[lev == 2]))
    sem3 = sd(na.exclude(measure[lev == 3]))/sqrt(length(na.exclude(measure[lev == 3])))
    mean3 = mean(na.exclude(measure[lev == 3]))

    df <- data.frame(
      group   = c(gg[[1]], gg[[2]], gg[[3]]),
      measure = c(mean1, mean2, mean3),
      sem     = c(sem1, sem2, sem3)
    )

    Ns <- tapply(measure,group,length)

    x <- aov(measure ~ group); xx <- summary(x); y <- TukeyHSD(x)

    AOV.res <- paste(paste('F',paste('(',paste(summary(x)[[1]]$Df[1],summary(x)[[1]]$Df[2], sep=','),
                                     ')', sep=''),
                           '=',
                           paste(round(summary(x)[[1]][["F value"]][[1]], digits = 2)),
                           sep=''),
                     paste('p=', (round(summary(x)[[1]][["Pr(>F)"]][[1]], digits = 5)), sep = ''),
                     sep = '; ')

    e <- lsr::etaSquared(x)
    eta <- paste('eta.sq', ' = ',round(e[1],2))

    p2.1 <- paste('p = ', round(y$`group`[1,4], 4), sep='')
    p3.1 <- paste('p = ', round(y$`group`[2,4], 4), sep='')
    p3.2 <- paste('p = ', round(y$`group`[3,4], 4), sep='')

    n2.1 <- paste(gg[[2]], '-', gg[[1]],'  ' ,sep = '')
    n3.1 <- paste(gg[[3]], '-', gg[[1]],'  ' ,sep = '')
    n3.2 <- paste(gg[[3]], '-', gg[[2]],'  ' ,sep = '')

    r2.1 <- paste(n2.1, p2.1)
    r3.1 <- paste(n3.1, p3.1)
    r3.2 <- paste(n3.2, p3.2)

    text.res <- paste(AOV.res, eta, r2.1, r3.1, r3.2, sep = '\n')

    mid <- mean(measure, na.rm = T)/2

    plot <- ggplot2::ggplot(df, aes(x=group, y=measure, fill=group)) +
      geom_bar(position=position_dodge(), stat="identity", fill=c("white", "Grey80", "Grey10"),
               colour="#000000", size=.4) +
      geom_errorbar(aes(ymin=measure-sem, ymax=measure+sem), width=.2, position=position_dodge(.8)) +
      guides(fill=F) +
      theme(panel.background = element_rect(fill = "white")) +
      labs(x = '', y = '',title = paste(substitute(measure))[[3]]) +
      {if(legend)geom_label(aes(x = 2, y = mid, label = text.res), fill = "white")} +
      theme(plot.title = element_text(hjust = 0.5))

    print(paste(substitute(measure))[[3]])
    print(Ns)
    print(AOV.res)
    print(eta)
    print(r2.1)
    print(r3.1)
    print(r3.2)


    return(plot)

  }
}
# bartest(ToothGrowth$len, ToothGrowth$supp, legend = F)
# bartest(ToothGrowth$len, ToothGrowth$supp, legend = T)
# bartest(ToothGrowth$len, ToothGrowth$dose, legend = T)
alemiani/explora documentation built on May 28, 2019, 4:54 p.m.