【问题标题】:r : ecdf over histogramr : 直方图上的 ecdf
【发布时间】:2015-05-31 02:15:42
【问题描述】:

在 R 中,使用 ecdf 我可以绘制经验累积分布函数

plot(ecdf(mydata))

使用hist,我可以绘制数据的直方图

hist(mydata)

如何在同一个图中绘制直方图和 ecdf?

编辑

我尝试做类似的东西

https://mathematica.stackexchange.com/questions/18723/how-do-i-overlay-a-histogram-with-a-plot-of-cdf

【问题讨论】:

  • 请务必查看我参考的问题和答案,了解为什么 ggplot2 不支持这种视觉效果。也就是说,在基础 R 中是可能的。但肯定需要不止一条线。

标签: r plot histogram cdf


【解决方案1】:

也有点晚了,这是另一个解决方案,它用第二个 y 轴扩展了 @Christoph 的解决方案。

par(mar = c(5,5,2,5))
set.seed(15)
dt <- rnorm(500, 50, 10)
h <- hist(
  dt,
  breaks = seq(0, 100, 1),
  xlim = c(0,100))

par(new = T)

ec <- ecdf(dt)
plot(x = h$mids, y=ec(h$mids)*max(h$counts), col = rgb(0,0,0,alpha=0), axes=F, xlab=NA, ylab=NA)
lines(x = h$mids, y=ec(h$mids)*max(h$counts), col ='red')
axis(4, at=seq(from = 0, to = max(h$counts), length.out = 11), labels=seq(0, 1, 0.1), col = 'red', col.axis = 'red')
mtext(side = 4, line = 3, 'Cumulative Density', col = 'red')

诀窍如下:不要在绘图中添加线,而是在顶部绘制另一个绘图,这就是我们需要par(new = T) 的原因。然后您必须稍后添加 y 轴(否则它将绘制在左侧的 y 轴上)。

学分去here(@tim_yates 回答)和there

【讨论】:

    【解决方案2】:

    有两种方法可以解决这个问题。一种是忽略不同的比例并在直方图中使用相对频率。这导致更难阅读直方图。第二种方法是改变一个或另一个元素的比例。

    我怀疑this question 很快就会引起你的兴趣,尤其是@hadley 的回答。

    ggplot2 单尺度

    这是ggplot2 中的解决方案。我不确定您是否会对结果感到满意,因为 CDF 和直方图(计数或相对)在完全不同的视觉尺度上。请注意,此解决方案的数据位于名为 mydata 的数据框中,所需变量位于 x 中。

    library(ggplot2)
    set.seed(27272)
    mydata <- data.frame(x=  rexp(333, rate=4) + rnorm(333))
    
     ggplot(mydata, aes(x)) + 
         stat_ecdf(color="red") + 
         geom_bar(aes(y = (..count..)/sum(..count..))) 
    

    base R 多尺度

    在这里,我将重新调整经验 CDF,使其最大值不是 1,而是具有最高相对频率的 bin。

    h  <- hist(mydata$x, freq=F)
    ec <- ecdf(mydata$x)
    lines(x = knots(ec), 
        y=(1:length(mydata$x))/length(mydata$x) * max(h$density), 
        col ='red')
    

    【讨论】:

      【解决方案3】:

      您可以尝试使用第二个轴的 ggplot 方法

      set.seed(15)
      a <- rnorm(500, 50, 10)
      
      # calculate ecdf with binsize 30
      binsize=30
      df <- tibble(x=seq(min(a), max(a), diff(range(a))/binsize)) %>% 
              bind_cols(Ecdf=with(.,ecdf(a)(x))) %>% 
              mutate(Ecdf_scaled=Ecdf*max(a))
      # plot
      ggplot() + 
        geom_histogram(aes(a), bins = binsize) +
        geom_line(data = df, aes(x=x, y=Ecdf_scaled), color=2, size = 2) + 
        scale_y_continuous(name = "Density",sec.axis = sec_axis(trans = ~./max(a), name = "Ecdf"))
      

      编辑

      由于缩放错误,我添加了第二个解决方案,提前计算所有内容:

      binsize=30
      a_range= floor(range(a)) +c(0,1)
      
      b <- seq(a_range[1], a_range[2], round(diff(a_range)/binsize)) %>% floor() 
      
      
      df_hist <- tibble(a) %>% 
        mutate(gr = cut(a,b, labels = floor(b[-1]), include.lowest = T, right = T)) %>% 
        count(gr) %>% 
        mutate(gr = as.character(gr) %>% as.numeric()) 
      
      # calculate ecdf with binsize 30
      df <- tibble(x=b) %>% 
        bind_cols(Ecdf=with(.,ecdf(a)(x))) %>% 
        mutate(Ecdf_scaled=Ecdf*max(df_hist$n))
        
      ggplot(df_hist, aes(gr, n)) + 
         geom_col(width = 2, color = "white") + 
         geom_line(data = df, aes(x=x, y=Ecdf*max(df_hist$n)), color=2, size = 2) +
         scale_y_continuous(name = "Density",sec.axis = sec_axis(trans = ~./max(df_hist$n), name = "Ecdf"))
      

      【讨论】:

      • ECDF线的归一化不正确(我也不知道怎么做更好)。最好根据最高的直方图条进行归一化。
      • 感谢您的快速响应和更新。与此同时,我想出了一个我自己的答案,你可能想验证它的正确性(我不是 R 的专家!)。
      【解决方案4】:

      正如已经指出的那样,这是有问题的,因为您要合并的图具有如此不同的 y 尺度。你可以试试

      set.seed(15)
      mydata<-runif(50)
      hist(mydata, freq=F)
      lines(ecdf(mydata))
      

      得到

      【讨论】:

        【解决方案5】:

        虽然有点晚了......另一个使用预设垃圾箱的版本:

        set.seed(15)
        dt <- rnorm(500, 50, 10)
        h <- hist(
            dt,
            breaks = seq(0, 100, 1),
            xlim = c(0,100))
            ec <- ecdf(dt)
            lines(x = h$mids, y=ec(h$mids)*max(h$counts), col ='red')
            lines(x = c(0,100), y=c(1,1)*max(h$counts), col ='red', lty = 3) # indicates 100%
            lines(x = c(which.min(abs(ec(h$mids) - 0.9)), which.min(abs(ec(h$mids) - 0.9))), # indicates where 90% is reached
                  y = c(0, max(h$counts)), col ='black', lty = 3)
        

        (只有第二个 y 轴还没有工作……)

        【讨论】:

          【解决方案6】:

          除了之前的答案,我想让 ggplot 进行繁琐的计算(与 @Roman's solution 相比,它已根据我的要求进行了更新),即计算并绘制直方图 计算并叠加 ECDF。我想出了以下(伪代码):

          # 1. Prepare the plot
          plot <- ggplot() + geom_hist(...)
          
          # 2. Get the max value of Y axis as calculated in the previous step
          maxPlotY <- max(ggplot_build(plot)$data[[1]]$y)
          
          # 3. Overlay scaled ECDF and add secondary axis
          plot +
            stat_ecdf(aes(y=..y..*maxPlotY)) +
            scale_y_continuous(name = "Density", sec.axis = sec_axis(trans = ~./maxPlotY, name = "ECDF"))
          

          这样您就不需要事先计算所有内容并将结果提供给ggpplot。放轻松,让它为你做一切!

          【讨论】:

            猜你喜欢
            • 2016-06-02
            • 2022-11-15
            • 1970-01-01
            • 1970-01-01
            • 1970-01-01
            • 1970-01-01
            • 2015-02-09
            • 1970-01-01
            • 2013-06-17
            相关资源
            最近更新 更多