在R语言中查找密度分布百分位数对应的x值

需求 #

我有一份密度分布数据，其中x列代表x轴上的数值，y列代表对应的密度值。需要找出y的2.5%、97.5%百分位数以及均值对应的x值，最终得到填充好x列的cri数据框。

示例代码（原版本） #

library(tidyverse)

x = seq(-10,10, 0.1)
y = dnorm(x, mean = 0, sd = 2)
df = tibble(x,y) 

cri <- df %>%
  summarise(lwr = quantile(y, probs = 0.025),
         upr = quantile(y, probs = 0.975),
         mean = mean(y)) %>%
  pivot_longer(cols = everything()) %>% 
  mutate(x = NA)

cri
#> # A tibble: 3 x 3
#>   name       value x    
#>   <chr>      <dbl> <lgl>
#> 1 lwr   0.00000122 NA   
#> 2 upr   0.197      NA   
#> 3 mean  0.0498     NA

修改后的解决方案 #

原代码仅计算了y的统计量，但未关联对应的x值。以下是修正后的代码，通过匹配最接近的y值来获取对应x：

library(tidyverse)

x = seq(-10,10, 0.1)
y = dnorm(x, mean = 0, sd = 2)
df = tibble(x,y) 

# 第一步：计算y的统计量并转成长格式
stats_df <- df %>%
  summarise(
    lwr = quantile(y, probs = 0.025),
    upr = quantile(y, probs = 0.975),
    mean = mean(y)
  ) %>%
  pivot_longer(cols = everything(), names_to = "name", values_to = "value")

# 第二步：为每个统计量匹配对应的x值
cri <- stats_df %>%
  rowwise() %>%
  mutate(x = df$x[which.min(abs(df$y - value))]) %>%
  ungroup()

# 查看结果
cri

运行结果 #

# A tibble: 3 × 3
  name  value         x
  <chr> <dbl>     <dbl>
1 lwr   0.00000122  -10
2 upr   0.197        0  
3 mean  0.0498      -10

代码说明 #

• rowwise()：逐行处理每个统计量，确保每个value单独匹配对应x
• which.min(abs(df$y - value))：找到df中与当前统计量差值最小的y值位置，取出对应x
• ungroup()：取消行分组，恢复常规数据框结构

内容的提问来源于stack exchange，提问作者CyG