最大似然估计（MLE）的一些公式与定理（python实践）

单参的情况

给定观察集$\mathcal{D}=(F_i,e_i)$，估计真实的$F_{\textrm{true}}$：

$$P(\mathcal{D}_i|F_{\textrm{true}})=\frac1{\sqrt{2\pi e_i^2}}\exp\begin{bmatrix}-\frac{(F_i-F_{\textrm{true}})^2}{2e_i^2}\end{bmatrix}$$

将$P(\mathcal{D}_i|F_{\textrm{true}})$这一条件概率（conditional probability）称为给定$F_{\textrm{true}}$（参数）下的样本$\mathcal{D}_i$发生的概率，根据上述等式其等于以$F_\textrm{true}$为均值，以$e_i$为标准差的正太分布。

在整个数据集$\mathcal{D}$（或者the entire set of measurements，观察集）上，构建似然函数（$\textbf{likelihood function}$）：

$$\mathcal{L}(\mathcal{D}|F_{\textrm{true}})=\prod_{i=1}^NP(\mathcal{D}_i|F_{\textrm{true}})=\prod_{i=1}^N\frac1{\sqrt{2\pi e_i^2}}\exp\begin{bmatrix}-\frac12\frac{(F_i-F_\textrm{true})^2}{2e_i^2}\end{bmatrix}$$

为简化计算，对数形式，也即对数似然（log_likelihood function）如下：

$$\log\,\mathcal{L}=-\frac12\sum_{i=1}^N\begin{bmatrix} \log(2\pi e_i^2)+\frac{(F_i-F_\textrm{true})^2}{2e_i^2} \end{bmatrix}$$

等式两边对$F_\textrm{true}$求偏导，$\frac{d\,\mathcal{L}}{d\,F_\textrm{true}}=0$：

$$F_\textrm{est}=\frac{\sum w_iF_i}{\sum w_i}=\frac{w^TF}{\|w\|_{\ell_1}}\\w_i=\frac1{e_i^2}\\\sigma_\textrm{est}=(\sum_{i=1}^Nw_i)^{-1/2}$$

在$e_i$相等时：

$$F_\textrm{est}=\frac{\sum_{i=1}^N F_i}N$$

def log_likelihood(theta, F, e):
    return -.5*np.sum(np.log(2*np.pi*e**2)+(F-theta[0])**2/(e**2))

双参的情况

$$F_\textrm{true}\sim\frac1{\sqrt{2\pi \sigma^2}}\exp\begin{bmatrix}-\frac12\frac{(F-\mu)^2}{\sigma^2}\end{bmatrix}$$

构建似然函数：

$$\mathcal{L}(\mathcal{D}|\theta)=\prod_{i=1}^N\frac1{\sqrt{2\pi (\sigma^2+e_i^2)}}\exp\begin{bmatrix}-\frac12\frac{(F_i-\mu)^2}{\sigma^2+e_i^2}\end{bmatrix}$$

同理，对数似然形式为：

$$\log\,\mathcal{L}(\mathcal{D}|\theta)=-\frac12\sum\begin{bmatrix}\log (2\pi(\sigma^2+e_i^2))+\frac{(F_i-\mu)^2}{\sigma^2+e_i^2}\end{bmatrix}$$

接下来的处理和单参的情况类似（取对数，求关于$\mu$和$\sigma$的偏导，置0，解析解（analytical solution，closed form solution））。

$$\mu_\textrm{est}=\frac{\sum w_iF_i}{\sum w_i}\\w_i=\frac1{\sigma^2+e_i^2}$$

def log_likelihood(theta, F, e):
    return -.5*np.sum(np.log(2*np.pi*(theta[1]**2+e**2))+(F-theta[0])**2/(theta[1]**2+e**2))

此时，我们可以使用scipy下的最优化函数，一般是最小化目标函数（objective function），所以：

from scipy import optimize

def neg_log_likelihood(theta, F, e):
    return -log_likelihood(theta, F, e)

theta_gauss = [900, 5]
theta_est = optimize.fmin(func=neg_log_likelihood, x0=theta_gauss, args=(F, e))