最大似然估计伪代码

Question

我需要写一个最大似然估计器，用来估算一些玩具数据的均值和方差。我有一个包含100个样本的向量，是用numpy.random.randn(100)生成的。这些数据应该是均值为零、方差为1的高斯分布。

我查了维基百科和其他一些资料，但因为我没有统计学背景，所以有点困惑。

有没有什么伪代码可以参考一下最大似然估计器？我对最大似然估计的直觉有了，但就是不知道从哪里开始编码。

维基百科上说要取对数似然的最大值。我理解的是：我需要用不同的参数来计算对数似然，然后选择那些给出最大概率的参数。但我不明白的是：我到底要在哪里找到这些参数？如果我随机尝试不同的均值和方差来获得高概率，那我应该什么时候停止尝试呢？

Answer 1

你需要一个数字优化的方法。虽然不确定Python里有没有现成的实现，但如果有的话，应该是在numpy或scipy这些库里。

可以看看像“尼尔德-梅德算法”或者“BFGS”这样的东西。如果实在找不到，可以用Rpy这个工具，调用R语言里的“optim()”函数。

这些函数的工作原理是搜索函数的空间，试图找出最大值在哪里。想象一下你在雾里找山顶。你可能会选择一直朝着最陡的方向走，或者可以让一些朋友带着对讲机和GPS去探路，帮你测量一下。无论哪种方法，都可能会让你误以为找到了最高点，所以通常需要从不同的地方多试几次。否则，你可能会认为南边的山顶是最高的，其实北边还有一个更高的山顶在遮挡着它。

Answer 2

我刚看到这个，虽然有点老旧，但我希望能帮助到其他人。之前的评论对机器学习优化的描述还不错，但没有人提供实现的伪代码。Python的Scipy库里有一个最小化器可以做到这一点。下面是线性回归的伪代码。

# import the packages
import numpy as np
from scipy.optimize import minimize
import scipy.stats as stats
import time

# Set up your x values
x = np.linspace(0, 100, num=100)

# Set up your observed y values with a known slope (2.4), intercept (5), and sd (4)
yObs = 5 + 2.4*x + np.random.normal(0, 4, 100)

# Define the likelihood function where params is a list of initial parameter estimates
def regressLL(params):
    # Resave the initial parameter guesses
    b0 = params[0]
    b1 = params[1]
    sd = params[2]

    # Calculate the predicted values from the initial parameter guesses
    yPred = b0 + b1*x

    # Calculate the negative log-likelihood as the negative sum of the log of a normal
    # PDF where the observed values are normally distributed around the mean (yPred)
    # with a standard deviation of sd
    logLik = -np.sum( stats.norm.logpdf(yObs, loc=yPred, scale=sd) )

    # Tell the function to return the NLL (this is what will be minimized)
    return(logLik)

# Make a list of initial parameter guesses (b0, b1, sd)    
initParams = [1, 1, 1]

# Run the minimizer
results = minimize(regressLL, initParams, method='nelder-mead')

# Print the results. They should be really close to your actual values
print results.x

这个对我来说效果很好。虽然这只是基础部分，它没有对参数估计进行分析或提供置信区间，但这是一个开始。你也可以使用机器学习的方法来找到，比如说，常微分方程（ODE）和其他模型的估计值，具体我在这里有描述。

我知道这个问题有点老，希望你们到现在已经搞明白了，但我也希望其他人能从中受益。

Answer 3

如果你要进行最大似然估计的计算，第一步就是要假设一个依赖于某些参数的分布。因为你是自己生成数据的（你甚至知道你的参数），所以你可以“告诉”程序假设是高斯分布。不过，你并不告诉程序你的参数（比如0和1），而是把它们先留作未知，之后再计算出来。

现在，你有一个样本向量（我们叫它 x，它的元素是 x[0] 到 x[100]），你需要对它进行处理。为此，你需要计算以下内容（f 表示高斯分布的概率密度函数）：

f(x[0]) * ... * f(x[100])

正如我在链接中所示，f 使用了两个参数（希腊字母 µ 和 σ）。现在，你需要计算 µ 和 σ 的值，使得 f(x[0]) * ... * f(x[100]) 的值尽可能大。

当你完成这个计算后，µ 就是你对均值的最大似然估计值，而 σ 是对标准差的最大似然估计值。

需要注意的是，我并没有明确告诉你 如何 计算 µ 和 σ 的值，因为这个过程相当数学化，我也不太懂；我只是告诉你获取这些值的技术，这种技术也可以应用于其他分布。

因为你想要最大化原始的表达式，你可以“简单地”最大化原始表达式的对数——这样可以避免处理所有这些乘法，把原始表达式转化为一些加法。

如果你真的想计算，可以做一些简化，得到以下表达式（希望我没有搞错）：

现在，你需要找到 µ 和 σ 的值，使得上面的表达式达到最大值。这个过程是一个非常复杂的任务，叫做非线性优化。

你可以尝试的一个简化方法是：固定一个参数，然后计算另一个参数。这样可以避免同时处理两个变量。