Calculate standard error from correlation and add as error bars to barplot

I want to plot the correlation between several factors in my data set. If possible, I'd like to try to add error bars or whiskers to these plotted values. The example dataset of mtcars isn't the greatest, but it's usable. The actual data uses a time series along the x axis. I have spearman specified because that's the correlation used in my analysis, not because it's the correct choice in the mtcars data set. I've seen some other posts suggesting the use of cor.test and pulling the values from that, but I'm not sure how that would be applied back to the bar chart to use as error bars. Here's the code to create the basic bar chart below.

mtcarstest <- mtcars %>%
  group_by(cyl) %>%
  summarise(COR = cor(disp,hp, method = "spearman", use="complete.obs"))

ggplot(data = mtcarstest) +
  aes(x = cyl, y = COR) +
  geom_bar(stat = "identity") +
  theme_minimal()

enter image description here

1 answer

  • answered 2020-05-25 06:04 jay.sf

    The cor.test yields a list where actually everything is stored what you need. So just write a function that grabs the desired values. We can use by here, which yields a list, that we can rbind to get a matrix with perfect row names for plotting. The do.call is required for the rbind of the data frames of a list.

    res <- do.call(rbind, by(mtcars, mtcars$cyl, function(x) {
      rr <- with(x, cor.test(disp, hp, method="pearson"))
      return(c(rr$estimate, CI=rr$conf.int))
    }))
    # res
    #          cor        CI1       CI2
    # 4  0.4346051 -0.2235519 0.8205544
    # 6 -0.5136284 -0.9133933 0.3904543
    # 8  0.1182556 -0.4399266 0.6105282
    

    Note, that method="spearman" won't work with the mtcars data, so I used "pearson" here.

    To plot the data I recommend barplot which comes with R. We store the bar locations b <- and use them as x-coordinates for the arrows. For the y-coordinates we take the values from our matrix.

    b <- barplot(res[,1], ylim=range(res)*1.2, 
                 main="My Plot", xlab="cyl", ylab="Cor. disp ~ hp")
    arrows(b, res[,2], b, res[,3], code=3, angle=90, length=.1)
    abline(h=0)
    box()
    

    enter image description here