R语言 在ggplot中,如何填充两条正态曲线之间的区域

htzpubme  于 2023-02-27  发布在  其他
关注(0)|答案(2)|浏览(191)

我有两条正态曲线,我想填充两条曲线之间的右侧区域,因此左侧曲线低于y限制,右侧曲线高于y限制。为了绘制曲线,我使用stat_function(),因此ggplot绘制曲线时没有在aes()中定义y列。我已经绘制了曲线和X轴之间的填充区域。但我需要两条曲线之间的区域,而使用NA清空左侧曲线的技巧似乎并不像我预期的那样有效。
生成图的代码在一个函数中,因为我需要绘制几个不同的正态曲线对。
我该怎么做呢?

library(ggplot2)
library(ggthemes)

graf_normal <- function(Xmedia1, Xdt1, Xmedia2, Xdt2) {

  Xmin1 <- Xmedia1-4*Xdt1
  Xmax1 <- Xmedia1+4*Xdt1

  Xmin2 <- Xmedia2-4*Xdt2
  Xmax2 <- Xmedia2+4*Xdt2

  Ymax1 <- max(dnorm(Xmedia1, Xmedia1, Xdt1))
  Ymax2 <- max(dnorm(Xmedia2, Xmedia2, Xdt2))

  Xmin <- min(Xmin1, Xmin2)
  Xmax <- max(Xmax1, Xmax2)

  ggplot(data.frame(X = c(Xmin, Xmax)), aes(x = X)) +
  geom_hline(yintercept = 0, colour = "grey", linewidth = 1) +
  stat_function(fun = dnorm, 
              args = c(Xmedia1, Xdt1), 
              linewidth = 1, 
              colour = "grey") +
  stat_function(fun = dnorm, 
              args = c(Xmedia2, Xdt2), 
              linewidth = 1, 
              colour = "black") +
  geom_segment(aes(x = Xmedia1, y = 0, xend = Xmedia1, yend = Ymax1), 
             linetype = "dashed", 
             linewidth = 0, 
             colour = "grey") +
  geom_segment(aes(x = Xmedia2, y = 0, xend = Xmedia2, yend = Ymax2), 
             linetype = "dashed", 
             linewidth = 0, 
             colour = "black") +
  ####################################################################
  stat_function(fun = dnorm,
            args = c(Xmedia2, Xdt2),
            xlim = c(Xmedia2+1.5*Xdt2,Xmax2),
            geom = "area",
            fill = "red",
            alpha = 0.5) +
  stat_function(fun = dnorm,
              args = c(Xmedia1, Xdt1),
              xlim = c(Xmedia1,Xmax1),
              geom = "area",
              fill = NA,
              alpha = 0.01) +

  ##################################################################

  theme(
    line = element_blank(),
    axis.line.y = element_blank(),
    axis.text.x = element_blank(),
    axis.text.y = element_blank(),
    axis. Ticks = element_blank(),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    legend. Position = "none",
    panel. Grid = element_blank(),
    panel. Background = element_rect(fill = "lightgray", colour = NA),
 ) +
 xlim(c(Xmin, Xmax)) 
}

g1 <- graf_normal(250, 7, 253, 7)

g1

我得到的两条曲线图如下:

谢谢你,

    • EDIT:**使用@stephan的代码和数据过滤,我已经能够做到这一点,使用geom_ribbon()更容易:

酷的方式区分重叠区!完整的代码:

library(ggplot2)

graf_normal <- function(Xmedia1, Xdt1, Xmedia2, Xdt2, n = 1000) {
  x1 <- Xmedia1 + 4 * Xdt1 * seq(-1, 1, length. Out = n)
  x2 <- Xmedia2 + 4 * Xdt2 * seq(-1, 1, length. Out = n)

  dat <- data. Frame(
    x = union(x1, x2)
  )
  dat$y1 <- dnorm(dat$x, Xmedia1, Xdt1)
  dat$y2 <- dnorm(dat$x, Xmedia2, Xdt2)

  Ymax1 <- dnorm(Xmedia1, Xmedia1, Xdt1)
  Ymax2 <- dnorm(Xmedia2, Xmedia2, Xdt2)

  ggplot(dat, aes(x)) +
    geom_hline(yintercept = 0, colour = "grey", linewidth = 1) +
    geom_ribbon(
      data =  subset(dat, x >= Xmedia2 + 1.5 * Xdt2),
      aes(ymin = y1, ymax = y2),
      fill = "red", alpha = 0.8
    ) +
    geom_ribbon(
      data =  subset(dat, (x <= Xmedia2 + 1.5 * Xdt2) & (y2 > y1)),
      aes(ymin = y1, ymax = y2),
      fill = "red", alpha = 0.2
    ) +
    geom_ribbon(
      data = subset(dat, x <= Xmedia1 - 1.5 * Xdt2),
      aes(ymin = y1, ymax = y2),
      fill = "blue", alpha = 0.8
    ) +
    geom_ribbon(
      data = subset(dat, (x <= Xmedia2 ) & (y1 > y2)),
      aes(ymin = y1, ymax = y2),
      fill = "blue", alpha = 0.2
    ) +
    annotate(
      geom = "segment",
      x = c(Xmedia1, Xmedia2), y = 0,
      xend = c(Xmedia1, Xmedia2), yend = c(Ymax1, Ymax2),
      linetype = "dashed",
      linewidth = 1,
      colour = c("grey", "black")
    ) +
    geom_line(aes(y = y1), linewidth = 1, colour = "grey") +
    geom_line(aes(y = y2), linewidth = 1, colour = "black") +
    theme(
      line = element_blank(),
      axis.line.y = element_blank(),
      axis. Text = element_blank(),
      axis. Ticks = element_blank(),
      axis. Title = element_blank(),
      legend. Position = "none",
      panel. Grid = element_blank(),
      panel. Background = element_rect(fill = "lightgray", colour = NA),
    )
}

graf_normal(250, 7, 253, 7)

然而,代码并不适用于所有曲线,工作!:

graf_normal(250, 7, 253, 3)

svmlkihl

svmlkihl1#

填充正态曲线之间的区域的一个选项是使用ggh4x::stat_difference,然而这需要手动计算密度值,并且通过geom_line而不是依赖于stat_function()来绘制:

library(ggplot2)
library(ggh4x)
graf_normal <- function(Xmedia1, Xdt1, Xmedia2, Xdt2, n = 101) {
  x1 <- Xmedia1 + 4 * Xdt1 * seq(-1, 1, length.out = n)
  x2 <- Xmedia2 + 4 * Xdt2 * seq(-1, 1, length.out = n)
  
  dat <- data.frame(
    x = union(x1, x2)
  )
  dat$y1 <- dnorm(dat$x, Xmedia1, Xdt1)
  dat$y2 <- dnorm(dat$x, Xmedia2, Xdt2)

  Ymax1 <- dnorm(Xmedia1, Xmedia1, Xdt1)
  Ymax2 <- dnorm(Xmedia2, Xmedia2, Xdt2)

  ggplot(dat, aes(x)) +
    geom_hline(yintercept = 0, colour = "grey", linewidth = 1) +
    ggh4x::stat_difference(
      data = ~ subset(.x, x >= Xmedia2 + 1.5 * Xdt2),
      aes(ymin = y1, ymax = y2)
    ) +
    annotate(
      geom = "segment",
      x = c(Xmedia1, Xmedia2), y = 0,
      xend = c(Xmedia1, Xmedia2), yend = c(Ymax1, Ymax2),
      linetype = "dashed",
      linewidth = 1,
      colour = c("grey", "black")
    ) +
    geom_line(aes(y = y1), linewidth = 1, colour = "grey") +
    geom_line(aes(y = y2), linewidth = 1, colour = "black") +
    scale_fill_manual(values = c(scales::alpha("red", .5), "transparent")) +
    theme(
      line = element_blank(),
      axis.line.y = element_blank(),
      axis.text = element_blank(),
      axis.ticks = element_blank(),
      axis.title = element_blank(),
      legend.position = "none",
      panel.grid = element_blank(),
      panel.background = element_rect(fill = "lightgray", colour = NA),
    )
}

graf_normal(250, 7, 253, 7)

编辑实际上stat_differnce并不是你的情况所需要的。考虑的太复杂了。正如@JuanRiera在他的评论中提到的,我们可以用geom_ribbon来填充这个区域:

library(ggplot2)

graf_normal <- function(Xmedia1, Xdt1, Xmedia2, Xdt2, n = 101) {
  x1 <- Xmedia1 + 4 * Xdt1 * seq(-1, 1, length.out = n)
  x2 <- Xmedia2 + 4 * Xdt2 * seq(-1, 1, length.out = n)

  dat <- data.frame(
    x = union(x1, x2)
  )
  dat$y1 <- dnorm(dat$x, Xmedia1, Xdt1)
  dat$y2 <- dnorm(dat$x, Xmedia2, Xdt2)

  Ymax1 <- dnorm(Xmedia1, Xmedia1, Xdt1)
  Ymax2 <- dnorm(Xmedia2, Xmedia2, Xdt2)

  ggplot(dat, aes(x)) +
    geom_hline(yintercept = 0, colour = "grey", linewidth = 1) +
    geom_ribbon(aes(ymin = y1, ymax = y2),
      data = subset(dat, x >= Xmedia2 + 1.5 * Xdt2),
      fill = "red", alpha = 0.5
    ) +
    annotate(
      geom = "segment",
      x = c(Xmedia1, Xmedia2), y = 0,
      xend = c(Xmedia1, Xmedia2), yend = c(Ymax1, Ymax2),
      linetype = "dashed",
      linewidth = 1,
      colour = c("grey", "black")
    ) +
    geom_line(aes(y = y1), linewidth = 1, colour = "grey") +
    geom_line(aes(y = y2), linewidth = 1, colour = "black") +
    theme(
      line = element_blank(),
      axis.line.y = element_blank(),
      axis.text = element_blank(),
      axis.ticks = element_blank(),
      axis.title = element_blank(),
      legend.position = "none",
      panel.grid = element_blank(),
      panel.background = element_rect(fill = "lightgray", colour = NA),
    )
}

graf_normal(250, 7, 253, 7)

xesrikrc

xesrikrc2#

一种方法是使用geom_polygon而不是stat_function
在函数中的ggplot()之前添加以下内容:

poly <- data.frame(xs = seq(Xmedia2+1.5*Xdt2, Xmax2, length.out = n)) |>
    transform(
      y1 = dnorm(xs, Xmedia1, Xdt1),
      y2 = dnorm(xs, Xmedia2, Xdt2)
    ) |>
    with(data.frame(X = c(xs, rev(xs)), Y = c(y1, rev(y2))))

然后将后面的两个stat_function调用替换为一个

geom_polygon(aes(X, Y), data = poly,
               fill = "red", alpha = 0.5) +

完整来源:

graf_normal <- function(Xmedia1, Xdt1, Xmedia2, Xdt2, n = 20) {

  Xmin1 <- Xmedia1-4*Xdt1
  Xmax1 <- Xmedia1+4*Xdt1

  Xmin2 <- Xmedia2-4*Xdt2
  Xmax2 <- Xmedia2+4*Xdt2

  Ymax1 <- max(dnorm(Xmedia1, Xmedia1, Xdt1))
  Ymax2 <- max(dnorm(Xmedia2, Xmedia2, Xdt2))

  Xmin <- min(Xmin1, Xmin2)
  Xmax <- max(Xmax1, Xmax2)

  poly <- data.frame(xs = seq(Xmedia2+1.5*Xdt2, Xmax2, length.out = n)) |>
    transform(
      y1 = dnorm(xs, Xmedia1, Xdt1),
      y2 = dnorm(xs, Xmedia2, Xdt2)
    ) |>
    with(data.frame(X = c(xs, rev(xs)), Y = c(y1, rev(y2))))

  ggplot(data.frame(X = c(Xmin, Xmax)), aes(x = X)) +
  geom_hline(yintercept = 0, colour = "grey", linewidth = 1) +
  geom_polygon(aes(X, Y), data = poly,
               fill = "red", alpha = 0.5) +
  stat_function(fun = dnorm, 
              args = c(Xmedia1, Xdt1), 
              linewidth = 1, 
              colour = "grey") +
  stat_function(fun = dnorm, 
              args = c(Xmedia2, Xdt2), 
              linewidth = 1, 
              colour = "black") +
  geom_segment(aes(x = Xmedia1, y = 0, xend = Xmedia1, yend = Ymax1), 
             linetype = "dashed", 
             linewidth = 0, 
             colour = "grey") +
  geom_segment(aes(x = Xmedia2, y = 0, xend = Xmedia2, yend = Ymax2), 
             linetype = "dashed", 
             linewidth = 0, 
             colour = "black") +
  ####################################################################
  # stat_function(fun = dnorm,
  #           args = c(Xmedia2, Xdt2),
  #           xlim = c(Xmedia2+1.5*Xdt2,Xmax2),
  #           geom = "area",
  #           fill = "red",
  #           alpha = 0.5) +
  # stat_function(fun = dnorm,
  #             args = c(Xmedia1, Xdt1),
  #             xlim = c(Xmedia1,Xmax1),
  #             geom = "area",
  #             fill = NA,
  #             alpha = 0.01) +
  ##################################################################
  theme(
    line = element_blank(),
    axis.line.y = element_blank(),
    axis.text.x = element_blank(),
    axis.text.y = element_blank(),
    axis.ticks = element_blank(),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    legend.position = "none",
    panel.grid = element_blank(),
    panel.background = element_rect(fill = "lightgray", colour = NA),
    ) +
    xlim(c(Xmin, Xmax)) 

}

注:我选择将geom_polygon移到绘图堆栈中更靠前的位置,这样dnorm线就会在红色区域的"顶部"。这是否重要取决于您的上下文和渲染引擎。

相关问题