Insert mean line in density plot for multiple facets using geom_vline

Viewed 551

For the dataset below, I would like to use a density plot to compare distributions in the two sexes for the following variables: withdraw, affect and estress . The plots should also contain vertical lines depicting the mean of the variables for each sex, and drawn on the same figure with different facets.

My current code is as follows, however I have no idea how to use geom_vline to get the mean for each facet:

vars_2_plot <- names(estress.dat1)[!names(estress.dat1) %in% "sex"]
estress.dat1 %>%
  pivot_longer(cols = all_of(vars_2_plot),
    names_to = "Predictors", values_to = "Values"
    ) %>%
  ggplot( aes(x = Values, fill = factor(sex) )) +   
  geom_density(alpha = 0.3) +
  facet_wrap(~ Predictors, scales = "free") +
  theme_minimal() +
labs(fill = "sex", x = "Predictor values") 

My data:

structure(list(withdraw = c(3, 1, 3.66, 4.66, 4.33, 3, 1, 1, 
2, 4, 4.33, 1, 5, 1.66, 4, 1.33, 1.66, 1.66, 2, 4, 2, 3.66, 3.66, 
4, 1, 1.66, 4.66, 5, 2.33, 1.33, 1.66, 1.33, 4, 1.66, 3, 1, 2, 
2, 2.66, 3, 1, 1, 4.66, 1, 1, 2, 1.33, 1, 3.66, 2.33, 1, 2, 2.33, 
2.33, 2.66, 2.33, 3, 1.66, 1, 2.66, 5, 3.33, 1, 2.33, 1, 1.33, 
3.33, 5, 2, 1, 2.33, 3, 1.66, 1, 1, 2, 1, 1, 4.33, 2, 1, 3, 1, 
2.66, 3.33, 2.33, 1, 1, 3, 2.66, 1.33, 2, 1, 1, 1, 2, 3, 1.33, 
2.33, 4.66, 2, 2.66, 1, 1, 4, 4, 2, 3, 1, 3, 1, 2.33, 1, 1, 2, 
2.33, 3, 1, 3, 1, 4.33, 3.33, 1, 1.33, 3.33, 1.66, 2, 1, 1, 2.33, 
4, 3, 2.33, 4, 4, 1, 1.33, 1, 1, 1.66, 2, 2.66, 3.66, 3.33, 1, 
2.33, 2.33, 1, 1, 1, 4.66, 2.33, 3, 3.33, 4.33, 1.33, 1, 3.66, 
4.33, 1, 2.66, 6.33, 2.66, 1, 3.33, 4.33, 1, 1.33, 1, 1, 3, 1, 
4, 1.33, 1, 2, 1.66, 1, 4, 1.66, 1, 1, 1.66, 5, 2.33, 1.66, 3, 
2.66, 1, 3.33, 4, 3, 2.66, 1, 3.66, 3, 2, 4, 7, 2, 1, 3.66, 2.66, 
1, 1.66, 1, 1.33, 3, 4, 4.66, 4, 2.66, 1, 1, 1, 1.66, 2.33, 3, 
3, 1, 4, 2.33, 4.33, 2, 2.66, 2, 1, 1, 1.66, 1, 1, 1.33, 2.33, 
1, 4, 2, 1, 2, 3, 3.66, 3, 3.66, 4, 3.66, 1, 1, 1, 2, 1, 2.66, 
1, 4, 1, 2.33, 4.33, 1.66, 2.66, 3.33, 3, 1, 2.33, 4), affect = c(2.6, 
1, 2.4, 1.16, 1, 1.5, 1, 1.16, 1.33, 3, 3, 2, 1.83, 1.16, 1.16, 
1, 1.33, 1.5, 1.33, 3, 3.16, 1.5, 1.5, 1.83, 1, 1.33, 1.5, 1.33, 
2.16, 1.16, 2.16, 1.83, 1.4, 1, 1.33, 1, 2.33, 1.16, 1, 1, 1.16, 
1.5, 2.16, 1.5, 2.66, 1.5, 1.16, 1.83, 2.66, 1.33, 1, 1.33, 1.83, 
1, 1.5, 1, 1.33, 1.33, 1.33, 1.16, 3, 1, 1, 1.66, 1.5, 1, 1.83, 
3.66, 3, 1.66, 1.83, 1.5, 1.66, 1, 1.16, 2.16, 1.16, 1, 1.5, 
1, 1.5, 4.33, 1.16, 2.33, 3, 2.33, 1.16, 1, 1.33, 1, 1.33, 1, 
1.16, 1, 1, 1.66, 1.5, 1.66, 1.33, 3.16, 1, 1, 1, 1.33, 2, 2.16, 
1, 1, 1, 2.33, 1.16, 3.83, 1.66, 1.33, 1.66, 1, 1.16, 1.16, 1.33, 
2, 1.16, 1.33, 1, 1.16, 1.33, 1.16, 1, 1.83, 1.16, 1, 2.83, 1, 
1.66, 2.33, 1, 1, 2, 1, 1, 2, 2.83, 1.33, 2, 2, 1.16, 1.16, 4.16, 
1.33, 1.16, 1.5, 1, 2, 1.33, 1.66, 2.66, 1, 1, 1.5, 1.66, 1, 
1.5, 4, 1.16, 1, 2, 3.5, 1.33, 1.16, 1.33, 1.33, 1.16, 1, 1.5, 
1, 2, 2.66, 1.33, 1.16, 2.66, 5, 1.5, 1, 1.5, 1.5, 1.16, 1.5, 
1, 3.16, 1.83, 1.16, 1, 1.33, 3.5, 1, 2, 1.4, 1.66, 1.66, 5, 
2.16, 2, 1.5, 1.33, 1.66, 2.83, 1.16, 1, 1.16, 1.5, 1.5, 1.83, 
1.83, 1, 1, 1, 1.16, 1.83, 1.66, 2, 1, 1.66, 1.16, 2.5, 1.16, 
1, 2, 1, 1.16, 1.83, 1.16, 1.33, 1.16, 1.6, 1, 1.5, 1.33, 1, 
2.16, 1.16, 3.16, 2.33, 2.16, 2, 1.5, 1, 1.33, 1.83, 1.33, 1.16, 
1, 1, 1, 1, 1, 1.83, 1, 1.5, 2, 1.33, 1.16, 1.83, 1.5), estress = c(6, 
5, 5.5, 3, 4.5, 6, 5.5, 3, 5.5, 6, 5.5, 4, 3, 2.5, 3.5, 6, 4, 
6, 3.5, 4, 2.5, 4, 2.5, 4.5, 1, 4.5, 3, 2.5, 2.5, 6, 4, 6.5, 
4.5, 3.5, 2, 5, 6, 3, 3.5, 4.5, 3.5, 5, 6, 4.5, 5, 5, 3, 4, 6.5, 
4, 3, 5.5, 4, 5, 3, 3, 2.5, 3.5, 5.5, 5, 7, 5, 6, 5.5, 6, 5, 
5, 5, 7, 5.5, 6.5, 3, 6.5, 4.5, 3, 5, 4.5, 3.5, 4.5, 1, 4, 6.5, 
2, 6.5, 5.5, 5, 6, 5, 3, 4.5, 5.5, 5, 1.5, 3, 3, 4.5, 6, 4, 3, 
7, 4.5, 3.5, 3.5, 6, 4, 6, 2.5, 3, 3, 5, 5, 7, 6.5, 6.5, 2, 4.5, 
5.5, 4.5, 2.5, 3.5, 5.5, 4.5, 4, 5.5, 7, 5, 3, 2.5, 4, 5, 7, 
4, 5, 3.5, 6.5, 2, 6, 3.5, 1, 6.5, 7, 4, 5.5, 5.5, 5.5, 6.5, 
7, 4.5, 4.5, 6.5, 6, 5.5, 5, 6.5, 6.5, 5, 4.5, 4.5, 4.5, 7, 6, 
7, 4, 6, 7, 7, 4, 4, 5.5, 4.5, 2.5, 4, 5.5, 4, 4, 5.5, 4.5, 5.5, 
4.5, 6, 5, 5.5, 6, 4, 3.5, 5, 4.5, 3.5, 6, 3.5, 4, 4, 5.5, 7, 
5, 4.5, 3.5, 5, 1, 5.5, 6.5, 5.5, 3, 7, 7, 3, 5, 3, 3, 6, 2.5, 
7, 3, 1.5, 1.5, 5.5, 5.5, 6, 5.5, 6.5, 3, 6.5, 6.5, 4, 5, 7, 
4.5, 4, 4, 5.5, 7, 4, 4.5, 4, 3.5, 4, 4, 6, 5, 6, 6, 4.5, 5.5, 
4.5, 4.5, 6, 4.5, 2.5, 3.5, 1, 3.5, 3.5, 4.5, 4, 5, 3.5, 4, 4, 
3.5, 3, 5.5, 5.5), sex = structure(c(2L, 1L, 2L, 2L, 2L, 2L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 2L, 
2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 
2L, 2L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 
1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 
2L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 
1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 
1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 
2L, 1L, 2L, 2L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 1L, 
2L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 
1L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 1L, 
1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L
), .Label = c("0", "1"), class = "factor")), row.names = c(NA, 
-262L), class = "data.frame")
1 Answers

For the data= argument in geom_vline you need a data frame with a column identical to your facet variable, as well as the corresponding mean values which you may calculate using colmeans.

library(magrittr);library(tidyr);library(ggplot2)
vars_2_plot <- names(estress.dat1)[!names(estress.dat1) %in% "sex"]

estress.dat1 %>%
  pivot_longer(cols=all_of(vars_2_plot),
               names_to="Predictors", values_to="Values"
  ) %>%
  ggplot(aes(x=Values, fill=factor(sex) )) +   
  geom_density(alpha=0.3) +
  geom_vline(data=data.frame(means=colMeans(estress.dat1[-4]),
                             Predictors=colnames(estress.dat1[-4])), 
             aes(xintercept=means), color="red") +
  facet_wrap(~ Predictors, scales="free") +
  theme_minimal() +
  labs(fill="sex", x="Predictor values")

Result

enter image description here

Related