使用 fitdist 函数(fitdistrplus 包)估计尺度和形状参数的尺度
Scaling for estimating scale and shape parameter with fitdist function (fitdistrplus package)
如标题中所述,我在 R(fitdistrplus
包)中的 fitdist
函数遇到缩放问题。
请看下面的代码:
# Initialize arrays for storing result
fit_store_scale <- rep(NA, 3)
fit_store_shape <- rep(NA, 3)
# load data
data1 <- c(7.616593e-05, 5.313253e-05, 1.604328e-04, 6.482365e-05,
4.217499e-05, 6.759114e-05, 3.531301e-05, 1.934228e-05,
6.263665e-05, 8.796205e-06)
data2 <- c(7.616593e-06, 5.313253e-06, 1.604328e-05, 6.482365e-06,
4.217499e-06, 6.759114e-06, 3.531301e-06, 1.934228e-06,
6.263665e-06, 8.796205e-07)
data3 <- c(7.616593e-07, 5.313253e-07, 1.604328e-06, 6.482365e-07,
4.217499e-07, 6.759114e-07, 3.531301e-07, 1.934228e-07,
6.263665e-07, 8.796205e-08)
# form data frame
data <- data.frame(data1, data2, data3)
# set scaling factor
scaling <- 1 #works without warnings and errors at:
#10000 (data1), 100000 (data2) or
#1000000 (data3)
# store scale and shape parameter of data1, data2 and data3 in Array
for(i in 1:3)
{
fit.w1 <- fitdist(data[[i]]*scaling,"weibull", method = "mle")
fit_store_scale[i] <- fit.w1$estimate[[2]]*1/scaling
#1/scaling is needed for correcting scale parameter
fit_store_shape[i] <- fit.w1$estimate[[1]]
}
我有三个数据向量,它们存储在一个数据框中。现在我想使用 fitdist
函数为每一列数据(data1
、data2
和 data3
)分别估计比例和形状参数,最后将它们存储在 fit_store_scale
和 fit_store_shape
分别。
这里的问题是 fitdist
函数在没有适当的比例因子的情况下无法工作,并且 data1
、data2
和 data3
需要不同的因子。我正在寻找一种解决方案来自动为每列数据确定最佳比例因子,从而使 fitdist
函数最终起作用。
解决这个问题的一种方法是通过按 10^j
:
缩放来继续尝试拟合分布
for(i in 1:3)
{
j <- 0
while(inherits(try(fitdist(data[[i]] * 10^j, "weibull", method = "mle"), silent = TRUE), "try-error"))
{
j <- j + 1
}
cat("\nFor data[[", i, "]], used j =", j, "\n\n")
fit.w1 <- fitdist(data[[i]] * 10^j, "weibull", method = "mle")
fit_store_scale[i] <- fit.w1$estimate[[2]] * 1/10^j
#1/scaling is needed for correcting scale parameter
fit_store_shape[i] <- fit.w1$estimate[[1]]
}
# For data[[ 1 ]], used j = 2
# For data[[ 2 ]], used j = 3
# For data[[ 3 ]], used j = 4
# > fit_store_scale
# [1] 6.590503e-05 6.590503e-06 6.590503e-07
# > fit_store_shape
# [1] 1.56613 1.56613 1.56613
也就是说,对于 data[[1]]
,我们成功地使用了 j = 2
(按 10^2 == 100
缩放),对于 data[[2]]
,我们使用了 j = 3 == 10^3 == 1,000
,并且对于 data[[3]]
,我们使用 j = 4 == 10^4 == 10,000
.
最终,这将找到 10 的最小幂来缩放数据并实现拟合。有关此 approach/theme.
的变体,请参阅 ?fitdist
下的示例 #14
如果您不是完全执着于 fitdist
,您可以使用更稳健的东西——以下使用对数尺度上的参数拟合 Weibull,并使用 Nelder-Mead 而不是基于梯度的方法。
拟合这些数据似乎没有任何问题。
dd <- data.frame(data1,data2,data3)
library("bbmle")
fx <- function(x) {
m1 <- mle2(y~dweibull(shape=exp(logshape),scale=exp(logscale)),
data=data.frame(y=x),start=list(logshape=0,logscale=0),
method="Nelder-Mead")
exp(coef(m1))
}
t(sapply(dd,fx)) ## not quite the output format you asked for,
## but easy enough to convert.
## logshape logscale
## data1 1.565941 6.589057e-05
## data2 1.565941 6.589054e-06
## data3 1.565941 6.589055e-07
对于您具有标准分布 (d*()
) 函数的任何分布,此方法应该相当有效。
如标题中所述,我在 R(fitdistrplus
包)中的 fitdist
函数遇到缩放问题。
请看下面的代码:
# Initialize arrays for storing result
fit_store_scale <- rep(NA, 3)
fit_store_shape <- rep(NA, 3)
# load data
data1 <- c(7.616593e-05, 5.313253e-05, 1.604328e-04, 6.482365e-05,
4.217499e-05, 6.759114e-05, 3.531301e-05, 1.934228e-05,
6.263665e-05, 8.796205e-06)
data2 <- c(7.616593e-06, 5.313253e-06, 1.604328e-05, 6.482365e-06,
4.217499e-06, 6.759114e-06, 3.531301e-06, 1.934228e-06,
6.263665e-06, 8.796205e-07)
data3 <- c(7.616593e-07, 5.313253e-07, 1.604328e-06, 6.482365e-07,
4.217499e-07, 6.759114e-07, 3.531301e-07, 1.934228e-07,
6.263665e-07, 8.796205e-08)
# form data frame
data <- data.frame(data1, data2, data3)
# set scaling factor
scaling <- 1 #works without warnings and errors at:
#10000 (data1), 100000 (data2) or
#1000000 (data3)
# store scale and shape parameter of data1, data2 and data3 in Array
for(i in 1:3)
{
fit.w1 <- fitdist(data[[i]]*scaling,"weibull", method = "mle")
fit_store_scale[i] <- fit.w1$estimate[[2]]*1/scaling
#1/scaling is needed for correcting scale parameter
fit_store_shape[i] <- fit.w1$estimate[[1]]
}
我有三个数据向量,它们存储在一个数据框中。现在我想使用 fitdist
函数为每一列数据(data1
、data2
和 data3
)分别估计比例和形状参数,最后将它们存储在 fit_store_scale
和 fit_store_shape
分别。
这里的问题是 fitdist
函数在没有适当的比例因子的情况下无法工作,并且 data1
、data2
和 data3
需要不同的因子。我正在寻找一种解决方案来自动为每列数据确定最佳比例因子,从而使 fitdist
函数最终起作用。
解决这个问题的一种方法是通过按 10^j
:
for(i in 1:3)
{
j <- 0
while(inherits(try(fitdist(data[[i]] * 10^j, "weibull", method = "mle"), silent = TRUE), "try-error"))
{
j <- j + 1
}
cat("\nFor data[[", i, "]], used j =", j, "\n\n")
fit.w1 <- fitdist(data[[i]] * 10^j, "weibull", method = "mle")
fit_store_scale[i] <- fit.w1$estimate[[2]] * 1/10^j
#1/scaling is needed for correcting scale parameter
fit_store_shape[i] <- fit.w1$estimate[[1]]
}
# For data[[ 1 ]], used j = 2
# For data[[ 2 ]], used j = 3
# For data[[ 3 ]], used j = 4
# > fit_store_scale
# [1] 6.590503e-05 6.590503e-06 6.590503e-07
# > fit_store_shape
# [1] 1.56613 1.56613 1.56613
也就是说,对于 data[[1]]
,我们成功地使用了 j = 2
(按 10^2 == 100
缩放),对于 data[[2]]
,我们使用了 j = 3 == 10^3 == 1,000
,并且对于 data[[3]]
,我们使用 j = 4 == 10^4 == 10,000
.
最终,这将找到 10 的最小幂来缩放数据并实现拟合。有关此 approach/theme.
的变体,请参阅?fitdist
下的示例 #14
如果您不是完全执着于 fitdist
,您可以使用更稳健的东西——以下使用对数尺度上的参数拟合 Weibull,并使用 Nelder-Mead 而不是基于梯度的方法。
拟合这些数据似乎没有任何问题。
dd <- data.frame(data1,data2,data3)
library("bbmle")
fx <- function(x) {
m1 <- mle2(y~dweibull(shape=exp(logshape),scale=exp(logscale)),
data=data.frame(y=x),start=list(logshape=0,logscale=0),
method="Nelder-Mead")
exp(coef(m1))
}
t(sapply(dd,fx)) ## not quite the output format you asked for,
## but easy enough to convert.
## logshape logscale
## data1 1.565941 6.589057e-05
## data2 1.565941 6.589054e-06
## data3 1.565941 6.589055e-07
对于您具有标准分布 (d*()
) 函数的任何分布,此方法应该相当有效。