如何拟合看起来像高斯的数据？

机器算法验证 r 正态分布
2022-03-02 11:09:36

我对统计很陌生，所以请原谅我可能使用了错误的词汇。

我有一些数据在绘制时看起来（对我来说）像高斯。

数据是从 jpeg 图像中提取的。这是从图像中截取的一条垂直线，仅使用红色数据（来自 RGB）。

以下是完整数据（27 个数据点）：

> r
 [1] 0.003921569 0.031372549 0.023529412 0.015686275 0.003921569 0.027450980
 [7] 0.003921569 0.015686275 0.031372549 0.105882353 0.305882353 0.490196078
[13] 0.560784314 0.615686275 0.592156863 0.505882353 0.364705882 0.227450980
[19] 0.050980392 0.031372549 0.019607843 0.054901961 0.031372549 0.015686275
[25] 0.027450980 0.003921569 0.011764706

> dput(r)
c(0.00392156862745098, 0.0313725490196078, 0.0235294117647059, 
0.0156862745098039, 0.00392156862745098, 0.0274509803921569, 
0.00392156862745098, 0.0156862745098039, 0.0313725490196078, 
0.105882352941176, 0.305882352941176, 0.490196078431373, 0.56078431372549, 
0.615686274509804, 0.592156862745098, 0.505882352941176, 0.364705882352941, 
0.227450980392157, 0.0509803921568627, 0.0313725490196078, 0.0196078431372549, 
0.0549019607843137, 0.0313725490196078, 0.0156862745098039, 0.0274509803921569, 
0.00392156862745098, 0.0117647058823529)
plot(r)

在此处输入图像描述

我想找到一个尽可能接近绘图/数据的高斯。

我尝试使用 R 包 mixtools 中的 normalmixEM。

> fit = normalmixEM(r)

但这似乎在默认情况下试图适应两个高斯的混合。

我尝试使用参数 k 指定只有一个高斯：

> fit = normalmixEM(r, k = 1)
Error in normalmix.init(x = x, lambda = lambda, mu = mu, s = sigma, k = k,  : 
  arbmean and arbvar cannot both be FALSE

我怎样才能拟合数据？

2个回答

我建议使用非线性最小二乘法进行此分析。

# First present the data in a data-frame
tab <- data.frame(x=seq_along(r), r=r)
#Apply function nls
(res <- nls( r ~ k*exp(-1/2*(x-mu)^2/sigma^2), start=c(mu=15,sigma=5,k=1) , data = tab))

从输出中，我能够获得以下拟合的“高斯曲线”：

v <- summary(res)$parameters[,"Estimate"]
plot(r~x, data=tab)
plot(function(x) v[3]*exp(-1/2*(x-v[1])^2/v[2]^2),col=2,add=T,xlim=range(tab$x) )

在此处输入图像描述

合身并不惊人......不会x↦sin⁡(x)/x功能是更好的模型？

拟合高斯分布和拟合高斯密度曲线是有区别的。正在做normalmixEM的是前者。你想要的是（我猜）后者。

粗略地说，拟合分布是如果您制作数据的直方图并尝试查看它具有什么样的形状，您会做什么。相反，您正在做的只是绘制一条曲线。这条曲线恰好在中间有一个驼峰，就像你通过绘制高斯密度函数得到的那样。

为了得到你想要的，你可以使用类似的东西optim来拟合你的数据曲线。下面的代码将使用非线性最小二乘法找到给出最佳拟合高斯曲线的三个参数：m是高斯均值，s是标准差，k是任意缩放参数（因为高斯密度被约束为积分为 1，而你的数据不是）。

x <- seq_along(r)

f <- function(par)
{
    m <- par[1]
    sd <- par[2]
    k <- par[3]
    rhat <- k * exp(-0.5 * ((x - m)/sd)^2)
    sum((r - rhat)^2)
}

optim(c(15, 2, 1), f, method="BFGS", control=list(reltol=1e-9))
其它你可能感兴趣的问题
上一篇Boruta 特征选择（在 R 中）是否考虑了变量之间的相关性？ 下一篇变化点分析