2

我想对我拥有的一些数据点执行对数 Pearson III 拟合。但是，每次我尝试它时，我都会收到错误消息，我真的不知道该怎么办。我也许应该补充一点，我几天前才开始使用 R，所以，我不是这方面的专家。

重要的代码部分，没有导入内容的部分等等是这样的：

pIIIpars<-list(shape=1, location=1, scale=1) 

dPIII<-function(x, shape, location, scale) PearsonDS::dpearsonIII(x, shape=1, location=1, scale=1, params=pIIIpars, log=FALSE)

pPIII<-function(q, shape, location, scale) PearsonDS::ppearsonIII(q, shape=1, location=1, scale=1, params=pIIIpars, lower.tail = TRUE, log.p = FALSE)

qPIII<-function(p, shape, location, scale) PearsonDS::qpearsonIII(p, shape=1, location=1, scale=1, params=pIIIpars, lower.tail = TRUE, log.p = FALSE)

fitPIII<-fitdistrplus::fitdist(flowdata3$OEP, distr="PIII", method="mle", start=list("shape"=5000, "location"=5000, "scale"=5000))

summary(fitPIII)

plot(fitPIII)

我正在使用 PearsonDS 包来定义 Log Pearson III 分布，并使用 fitdistrplus 来进行拟合。

我总是得到的错误信息是这样的：

[1] "Error in optim(par = vstart, fn = fnobj, fix.arg = fix.arg, obs = data,  : \n  function cannot be evaluated at initial parameters\n"
attr(,"class")
[1] "try-error"
attr(,"condition")
<simpleError in optim(par = vstart, fn = fnobj, fix.arg = fix.arg, obs = data,     ddistnam = ddistname, hessian = TRUE, method = meth, lower = lower,     upper = upper, ...): function cannot be evaluated at initial parameters>
Error in fitdistrplus::fitdist(flowdata3$OEP, distr = "PIII", method = "mle",  : 
  the function mle failed to estimate the parameters, 
                with the error code 100

我确实理解错误消息，它只是；如果这不是传递初始值的正确方法，那是什么？有人有想法吗？

干杯，罗伯特

4

2 回答 2

6

以下示例遵循 Kite (2004) 并匹配他的结果。

# Annual maximum discharge data for the St Mary River at Stillwater Nova Scotia (Kite, 2004)
# PIII fit to the logs of the discharges

StMary <- c(565,294,303,569,232,405,228,232,394,238,524,368,464,411,368,487,394,
            337,385,351,518,365,515,280,289,255,334,456,479,334,394,348,428,337,
            311,453,328,564,527,510,371,824,292,345,442,360,371,544,552,651,190,
            202,405,583,725,232,974,456,289,348,564,479,303,603,514,377,318,342,
            593,378,255,292)

LStMary <- log(StMary)

m <- mean(LStMary)
v <- var(LStMary)
s <- sd(LStMary)
g <- e1071::skewness(LStMary, type=1)

# Correct the sample skew for bias using the recommendation of 
# Bobee, B. and R. Robitaille (1977). "The use of the Pearson Type 3 and Log Pearson Type 3 distributions revisited." 
# Water Resources Reseach 13(2): 427-443, as used by Kite

n <- length(StMary)
g <- g*(sqrt(n*(n-1))/(n-2))*(1+8.5/n)

# We will use method of moment estimates as starting values for the MLE search

my.shape <- (2/g)^2
my.scale <- sqrt(v)/sqrt(my.shape)*sign(g) # modified as recommended by Carl Schwarz
my.location <- m-sqrt(v * my.shape)

my.param <- list(shape=my.shape, scale=my.scale, location=my.location)


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)

fitdistrplus::fitdist(LStMary, distr="PIII", method="mle", start=my.param)

另请注意，MLE 估计可能并不总是适用。参见风筝 (2004, p119)。他引用了 Matalas 和 Wallis (1973) 的话，他们指出如果样本偏差很小，那么解决方案是不可能的。您可以在矩估计器的方法中看到这一点，因为随着 g 变为零，形状参数会爆炸。

Kite, GW (2004) 水文学中的频率和风险分析。水资源出版物

北卡罗来纳州马塔拉斯和 JR 沃利斯 (1973)。“尤里卡！它符合 Pearson Type 3 Distribution。” 水资源研究 9（2）：281-289。

于 2013-03-22T11:56:27.137 回答
2

如果数据正偏斜，即长右尾，则上述代码将起作用。但是，如果要将数据拟合到负偏态分布，则需要将 my.scale 乘以偏态系数的符号。

请参阅第 24 页的公式 (8)、(9) 和 (10)

确定洪水流量和频率的指南 公告 17

可在 https://acwi.gov/hydrology/Frequency/b17c/获得

所以修改my.scale为 

my.scale <- sqrt(v)/sqrt(my.shape)*sign(g)

如果不修改尺度，初始值往往不会导致拟合函数收敛。

于 2017-01-19T05:44:46.933 回答