1

我使用我知道的正态分布的权重随机生成 1000 个数据点。现在我试图最小化 -log 似然函数来估计 sig^2 的值和权重。我从概念上理解了这个过程,但是当我尝试对其进行编码时,我就迷失了。

这是我的模型:

p(y|x, w, sig^2) = N(y|w0+w1x+...+wnx^n, sig^2)

我已经在谷歌上搜索了一段时间,我了解到 scipy.stats.optimize.minimize 函数对此很有用,但我无法让它正常工作。我尝试过的每个解决方案都适用于我从中获得解决方案的示例,但我无法将其推断为我的问题。

x = np.linspace(0, 1000, num=1000)
data = []
for y in x:
        data.append(np.polyval([.5, 1, 3], y))

#plot to confirm I do have a normal distribution...
data.sort()
pdf = stats.norm.pdf(data, np.mean(data), np.std(data))
plt.plot(test, pdf)
plt.show()

#This is where I am stuck.
logLik = -np.sum(stats.norm.logpdf(data, loc=??, scale=??))

我发现方程 error( w ) = .5*sum(poly(x_n, w ) - y_n)^2 与最小化权重误差有关,因此可以最大化我对权重的可能性,但我没有不明白如何编写代码...我发现 sig^2 有类似的关系,但有同样的问题。有人可以澄清如何做到这一点来帮助我的曲线拟合吗?也许去发布我可以使用的伪代码?

4

2 回答 2

3

是的,实现可能性拟合minimize很棘手,我花了很多时间。这就是我包装它的原因。如果我可以无耻地插入我自己的包symfit,你的问题可以通过这样做来解决:

from symfit import Parameter, Variable, Likelihood, exp
import numpy as np

# Define the model for an exponential distribution
beta = Parameter()
x = Variable()
model = (1 / beta) * exp(-x / beta)

# Draw 100 samples from an exponential distribution with beta=5.5
data = np.random.exponential(5.5, 100)

# Do the fitting!
fit = Likelihood(model, data)
fit_result = fit.execute()

我不得不承认我并不完全了解你的发行版,因为我不了解你的角色w,但也许以这段代码为例,你会知道如何适应它。

如果没有,请告诉我您模型的完整数学方程式,以便我进一步帮助您。

有关更多信息,请查看文档。(有关幕后发生的事情的更多技术描述,请阅读此处此处。)

于 2016-10-25T07:32:55.947 回答
1

我觉得你的设置有问题。通过最大似然,您可以获得最大化观察数据的概率的参数(给定特定模型)。您的模型似乎是:

在此处输入图像描述

其中 epsilon 是 N(0, sigma)。

所以你最大化它:

在此处输入图像描述

或等效地获取日志以获取:

在此处输入图像描述

在这种情况下, f 是对数正态概率密度函数,您可以使用 得到stats.norm.logpdf。然后,您应该使用scipy.minimize最大化一个表达式,该表达式将是 stats.norm.logpdf在每个 i 点处评估的总和,从 1 到您的样本大小。

如果我对您的理解正确,则您的代码缺少y向量加x向量!向我们展示这些向量的示例,我可以更新我的答案以包含用于估计具有该日期的 MLE 的示例代码。

于 2016-10-26T20:46:16.467 回答