Filling areas under two density curves in ggplot

The issue here is lag, it will shift the values along a vector, but it will trim one end off and fill the other with NA, e.g.

library(dplyr)
x <- 1:10
x
#>  [1]  1  2  3  4  5  6  7  8  9 10
lag(x, 2)
#>  [1] NA NA  1  2  3  4  5  6  7  8

What you need to do instead is generate new density values for a different normal distribution where the x values are shifted (or equivalently the mean has been shifted), below I've shifted it one unit to the right.

Not sure what packages you are using in your code so I have replaced the data generation with dnorm from base, which returns the density values for a given vector of x values.

library(dplyr)
library(ggplot2)

data <- data.frame(x = seq(-5, 5, length.out = 100)) %>%
  mutate(
    y = dnorm(x),
    e = dnorm(x - 1)
  )

ggplot(data) +
  geom_line(aes(x=x, y=y)) +
  geom_line(aes(x=x, y=e)) +
  geom_vline(aes(xintercept=0)) +
  geom_ribbon(data = subset(data, x<0), aes(x=x, ymin=0, ymax=y), fill = "blue", alpha = .3) +
  geom_ribbon(data = subset(data, e<y & x<0), aes(x=x, ymin=0, ymax=e), fill = "light blue", alpha = .5) +
  geom_ribbon(data = subset(data, y>e & x>0), aes(x=x, ymin=e, ymax=y), fill = "pink", alpha = .3) +
  geom_ribbon(data = subset(data, e<y & x>0), aes(x=x, ymin=0, ymax=e), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, y<=e & x>0), aes(x=x, ymin=0, ymax=y), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, e>y & x>0), aes(x=x, ymin=y, ymax=e), fill = "red", alpha = .3)

For your code to create data is might be something like this

library(dplyr)
library(ggplot2)

data <- distribution_normal(n = 100, mean = 0, sd = 1) %>%
  density() %>%
  as.data.frame()

shifted <- distribution_normal(n = 100, mean = 1, sd = 1) %>%
  density() %>%
  as.data.frame()

data <- data %>%
  mutate(e = shifted$y)

ggplot(data) +
  geom_line(aes(x=x, y=y)) +
  geom_line(aes(x=x, y=e)) +
  geom_vline(aes(xintercept=0)) +
  geom_ribbon(data = subset(data, x<0), aes(x=x, ymin=0, ymax=y), fill = "blue", alpha = .3) +
  geom_ribbon(data = subset(data, e<y & x<0), aes(x=x, ymin=0, ymax=e), fill = "light blue", alpha = .5) +
  geom_ribbon(data = subset(data, y>e & x>0), aes(x=x, ymin=e, ymax=y), fill = "pink", alpha = .3) +
  geom_ribbon(data = subset(data, e<y & x>0), aes(x=x, ymin=0, ymax=e), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, y<=e & x>0), aes(x=x, ymin=0, ymax=y), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, e>y & x>0), aes(x=x, ymin=y, ymax=e), fill = "red", alpha = .3)

Can't say for certain that this will work for you as I don't know the packages you've used

Its more noticeable in my code that yours, but there are some vertical gaps at the intercepts/boundaries of the ribbons. This is due to there not being a value exactly at these points, so to ensure they don't appear would be to calculate the x values and density values exactly at those points, which will likely be complicated at the density intercepts, as you'd need to calculate the exact x value at which it happens and then the corresponding density value. The brute force alternative is to plot a lot more values than needed, so that it is beyond notice (i.e. change the default n value in density for your way or make the length.out value in mine a lot larger e.g. 20000)

Here is a base R approach based off your original code.

library(bayestestR)
data <- distribution_normal(n = 100, mean = 0, sd = 1) %>%
  density() %>%
  as.data.frame()

original_length <- nrow(data)
step_size <- diff(data[1:2,1])
data <- rbind(data, data.frame(x = (step_size * 1:100) + max(data$x), y = 0))
data$e <- 0
data$e[seq(100,original_length+99)] <- data$y[seq(1,original_length)]

Thanks to the contributors above for their solutions. My solution is based off Olly Beagley's solution, and also shows how I got rid of the weird white spaces in the geom_ribbon.

I used the actual formula for a normal distribution instead of generating values using dnorm(x) or distribution_normal():

sigma <- 1
mean <- 0

sigma_2 <- 1
mean_2 <- 1

data <- data.frame(x = seq(-10,10,.001))

data <- data %>%
  mutate(y = 1/(sigma*sqrt(2*pi))*exp(-((x-mean)^2/(2*sigma^2))),
         y2 = 1/(sigma_2*sqrt(2*pi))*exp(-((x-mean_2)^2/(2*sigma_2^2)))
         )

ggplot(data) +
  geom_vline(aes(xintercept=0), size=1) +
  geom_line(aes(x=x, y=y)) +
  geom_line(aes(x=x, y=y2)) +
  geom_ribbon(data=subset(data, y>y2 & x <  0), aes(x=x, ymin=y2, ymax=y), fill = "blue", alpha = .3) +
  geom_ribbon(data = subset(data, y2<y & x< 0), aes(x=x, ymin=0, ymax=y2), fill = "light blue", alpha = .5) +
  geom_ribbon(data = subset(data, y>y2 & x> 0), aes(x=x, ymin=y2, ymax=y), fill = "pink", alpha = .3) +
  geom_ribbon(data = subset(data, y2<y & x>=0), aes(x=x, ymin=0, ymax=y2), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, y<=y2 & x>0), aes(x=x, ymin=0, ymax=y), fill = "purple", alpha = .3) +
  geom_ribbon(data = subset(data, y2>y & x>=x_intercept), aes(x=x, ymin=y, ymax=y2), fill = "red", alpha = .3) +
  xlim(-5,5)

This generated a y-value on each curve for each value of x: