我正在对我的流量数据进行对数皮尔逊 III 分布拟合。拟合后,我想绘制观测值和拟合分布。
然而,这个数字并不是我所期望的:
1. x 轴标签没有显示 0.99 0.02 和 0.01,尽管我已将其设置breaks
为 c(0.99、0.9、0.5、0.2、0.1、0.02、0.01)。
2.如何在图中绘制拟合的 Log Pearson III 分布?我提供了一个示例图供参考,谢谢您的帮助。
library(fitdistrplus)
library(PearsonDS)
library(ggplot2)
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)
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)
fit_lp3 <- fitdist(Log_mydata, distr = "PIII", start = start)
plotdata = low_flows
plotdata <- plotdata[order(plotdata$Flow), ]
plotdata$prob <- seq_len(length(plotdata$Flow))/(length(plotdata$Flow) + 1)
prob_scale_points = c(0.99, 0.9, 0.5, 0.2, 0.1, 0.02, 0.01)
ggplot(data = plotdata, aes(x = prob, y = Flow)) +
geom_point(size = 2) +
xlab("Probability") +
scale_x_continuous(transform = scales::probability_trans("norm", lower.tail = FALSE),
breaks = prob_scale_points,
sec.axis = dup_axis(name = "Return Period",
labels = function(x) {ifelse(1/x < 2, round(1/x, 2), round(1/x, 0))})) +
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))
我的目标如下:
按照 Allan 的建议,Log Pearson III 可以正确绘制了!为了简化我的原始代码,我尝试使用scale_x_log10
而不是scale_x_continuous
。但是,Log Pearson III 的拟合线没有正确绘制。我想知道如何解决这个问题?
ggplot(data = plotdata, aes(x = prob, y = Flow)) +
geom_point(size = 2, color="blue") +
geom_function(fun = ~ 10^(qpearsonIII(.x, params = fit_lp3$estimate)),
color = "black", linewidth = 1.5, alpha = 0.7) +
scale_x_log10(name = "Probability",
trans = "reverse",
breaks = prob_scale_points,
limits = rev(range(prob_scale_points)),
sec.axis = dup_axis(name = "Return Period",
labels = function(x) {ifelse(1/x < 2, round(1/x, 2), round(1/x, 0))})) +
scale_y_log10(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))
这是图。我们可以看到拟合的线不正确。
您可以使用
geom_function
绘制带有拟合参数的输出qPearsonIII
,但您应该绘制记录的数据以匹配其输出。这意味着重新标记 y 轴以使其显示为对数刻度。为了按所需方式获取 x 轴,请包括
limits
以及breaks
。