将分布拟合到R中的给定频率值

Question

我的频率值随时间(x轴单位)而变化,如下图所示.在一些归一化之后,这些值可以被视为某些分布的密度函数的数据点.

问:假设这些频率点来自威布尔分布T,我如何将最佳威布尔密度函数拟合到这些点,以便从中推断出分布T参数？

sample <- c(7787,3056,2359,1759,1819,1189,1077,1080,985,622,648,518,
            611,1037,727,489,432,371,1125,69,595,624)

plot(1:length(sample), sample, type = "l")
points(1:length(sample), sample)

更新.为了防止被误解,我想补充一点解释.通过说我的频率值随时间变化(x轴单位)我的意思是我有数据说我有:

• 7787价值实现1
• 3056价值实现2
• 2359实现价值3 ......等

某种方式实现我的目标(我认为不正确)将创建一组这些实现:

# Loop to simulate values 
set.values <- c()
for(i in 1:length(sample)){
  set.values <<- c(set.values, rep(i, times = sample[i]))
}

hist(set.values)
lines(1:length(sample), sample)
points(1:length(sample), sample)

并使用fitdistr在set.values:

f2 <- fitdistr(set.values, 'weibull')
f2

为什么我认为这是不正确的方式以及为什么我在寻找更好的解决方案R？

• 在上面给出的分布拟合方法中,假设它set.values是从分布中完整的一组实现T

• 在我原来的问题中,我知道密度曲线第一部分的点 - 我不知道它的尾部,我想估计尾部(和整个密度函数)

Answer 1

这是一个更好的尝试，就像之前它用来optim查找约束在框中的一组值（由调用中的lower和upper向量定义optim）的最佳值一样。请注意，除了威布尔分布形状参数之外，它还缩放 x 和 y 作为优化的一部分，因此我们有 3 个参数需要优化。

不幸的是，当使用所有点时，它几乎总是在约束框的边缘找到一些东西，这向我表明，威布尔可能并不适合所有数据。问题是这两点——它们太大了。您会在第一个图中看到尝试拟合所有数据。

如果我放弃前两点而只适合其余的点，我们会得到更好的拟合。您在第二个图中看到了这一点。我认为这是一个很好的选择，无论如何它都是约束框内部的局部最小值。

library(optimx)
sample <- c(60953,7787,3056,2359,1759,1819,1189,1077,1080,985,622,648,518,
            611,1037,727,489,432,371,1125,69,595,624)
t.sample <- 0:22

s.fit <- sample[3:23]
t.fit <- t.sample[3:23]

wx <- function(param) { 
  res <- param[2]*dweibull(t.fit*param[3],shape=param[1])
  return(res)
} 
minwx <- function(param){
  v <- s.fit-wx(param)
  sqrt(sum(v*v))
}

p0 <- c(1,200,1/20)
paramopt <- optim(p0,minwx,gr=NULL,lower=c(0.1,100,0.01),upper=c(1.1,5000,1))

popt <- paramopt$par
popt
rms <- paramopt$value
tit <- sprintf("Weibull - Shape:%.3f xscale:%.1f  yscale:%.5f rms:%.1f",popt[1],popt[2],popt[3],rms)

plot(t.sample[2:23], sample[2:23], type = "p",col="darkred")
lines(t.fit, wx(popt),col="blue")
title(main=tit)