如何用R的ggplot2绘制Log Pearson III分布？轴标签与拟合线问题

问题描述 #

初始问题 #

• 已设置x轴断点为c(0.99, 0.9, 0.5, 0.2, 0.1, 0.02, 0.01)，但标签未显示0.99、0.02和0.01
• 如何在图中绘制拟合的Log Pearson III分布

后续问题 #

改用scale_x_log10替代scale_x_continuous后，Log Pearson III的拟合线无法正确显示，需修复该问题

解决方案 #

1. x轴标签不显示问题修复 #

初始代码使用scales::probability_trans("norm", lower.tail = FALSE)进行轴转换，该转换会将概率值映射到正态分位数空间，导致设置的断点不在转换后的轴范围内，因此标签无法显示。改用对数反转轴（scale_x_log10(trans = "reverse")）可直接保留设置的断点，同时符合水文频率分析的常用可视化逻辑。

2. 拟合线显示问题修复 #

拟合线需基于拟合参数，通过qpearsonIII计算对应概率的分位数，再转换回原始流量尺度（取10的幂，因拟合对象是对数转换后的数据）。使用scale_x_log10后，x轴输入为原始概率值的对数反转结果，而qpearsonIII需要原始概率值，因此需在geom_function中对输入x值进行反转对数转换，还原为原始概率后再计算分位数。

修正后的完整代码 #

library(fitdistrplus)
library(PearsonDS)
library(ggplot2)
library(e1071)

low_flows <- data.frame(Year = 1981:2010, Flow = c(0.03357143, 0.01328571, 0.02285714, 0.02657143, 0.04957143, 0.04085714, 0.16900000, 0.01057143, 
                                                   0.04128571, 0.10871429, 0.08771429, 0.09585714, 0.22057143, 0.11571428, 0.08300000, 0.11257143, 
                                                   0.13614286, 0.07742857, 0.09785714, 0.04728571, 0.04300000, 0.08385714, 0.02828571, 0.07271429, 
                                                   0.21428571, 0.10142857, 0.04400000, 0.10928571, 0.12471429, 0.31500000))

# 对数转换流量数据
Log_mydata <- log10(low_flows$Flow)

# 计算初始参数
m <- mean(Log_mydata)
v <- stats::var(Log_mydata)
g <- e1071::skewness(Log_mydata, type = 2)
my_shape <- (2/g)^2
my_scale <- sqrt(v)/sqrt(my_shape) * sign(g)
my_location <- m - my_scale * my_shape
start <- list(shape = my_shape, location = my_location, scale = my_scale)

# 定义Pearson III分布的密度、分布、分位数函数
dPIII <- function(x, shape, location, scale) PearsonDS::dpearsonIII(x, 
                                                                     shape, location, scale, log = FALSE)
pPIII <- function(q, shape, location, scale) PearsonDS::ppearsonIII(q, 
                                                                     shape, location, scale, lower.tail = TRUE, log.p = FALSE)
qPIII <- function(p, shape, location, scale) PearsonDS::qpearsonIII(p, 
                                                                     shape, location, scale, lower.tail = TRUE, log.p = FALSE)

# 拟合Log Pearson III分布
fit_lp3 <- fitdist(Log_mydata, distr = "PIII", start = start)

# 准备绘图数据：排序流量并计算经验概率
plotdata <- low_flows[order(low_flows$Flow), ]
plotdata$prob <- seq_len(nrow(plotdata))/(nrow(plotdata) + 1)

prob_scale_points <- c(0.99, 0.9, 0.5, 0.2, 0.1, 0.02, 0.01)

# 绘制图形：修复拟合线和x轴标签问题
ggplot(data = plotdata, aes(x = prob, y = Flow)) + 
  geom_point(size = 2, color="blue") + 
  # 关键修复：将对数反转后的x值还原为原始概率值，再计算分位数
  geom_function(fun = ~ 10^(qpearsonIII(10^(.x), params = fit_lp3$estimate)),
                color = "black", linewidth = 1.5, alpha = 0.7) + 
  scale_x_log10(name = "概率", 
                trans = "reverse", 
                breaks = prob_scale_points, 
                limits = rev(range(prob_scale_points)), 
                sec.axis = dup_axis(name = "重现期", 
                                    labels = function(x) {ifelse(1/x < 2, round(1/x, 2), round(1/x, 0))})) + 
  scale_y_log10(name = "流量", breaks = seq(0.05, 0.4, 0.05)) + 
  theme_bw()  + 
  theme(panel.border = element_rect(colour = "black", fill = NA, linewidth = 1), 
        panel.grid = element_line(linewidth = 0.2), 
        axis.title = element_text(size = 12), axis.text = element_text(size = 10), 
        axis.title.x.top = element_text(size = 12), 
        legend.text = element_text(size = 10), legend.title = element_text(size = 10))

关键修复说明 #

• x轴标签显示：使用scale_x_log10(trans = "reverse")直接对概率轴进行对数反转，无需正态转换，确保设置的所有断点prob_scale_points正常显示。
• 拟合线修复：scale_x_log10(trans = "reverse")会将x轴输入值转换为log10(1/prob)，因此在geom_function中通过10^(.x)将输入x值还原为原始概率，再传入qpearsonIII计算对应分位数，最后用10^()转换回原始流量尺度。

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