python - Python中多元正态分布的最大似然估计

Question

我正在尝试使用 MLE 拟合多元正态分布的参数。

import numpy as np
from scipy.stats import norm
from scipy.optimize import minimize
from scipy.stats import multivariate_normal as mnorm

def estimation(obs,fun,init,method='Nelder-Mead'):
    mle = lambda param: -np.sum(fun(*[obs,param])) ## negate since we will minimize
    result = minimize(mle,init,method=method)
    return result.x

拟合单变量正态分布很好：

obs = np.random.normal(1,4,50000)
ini = [0,1]
print(estimation(obs,lambda ob,p:norm.logpdf(ob,p[0],p[1]),ini))

但是遇到了一些多变量问题（将数组分配给变量时出错）：

obs_m = np.random.multivariate_normal([0,0],[[1,0],[0,100]],50000)
ini_m = [[0,0],[[1,0],[0,100]]]
print(estimation(obs_m,lambda ob,p:mnorm.logpdf(ob,p[0],p[1],ini_m))

似乎优化算法不适用于任意数组/矩阵。我必须将平均数组和协方差矩阵解压缩到一个平面数组中以最小化。

ini_m = [0,0,1,0,0,100]
print(estimation(obs_m,lambda ob,p:mnorm.logpdf(ob,[p[0],p[1]],[[p[2],p[3]],[p[4],p[5]]]),ini_m))

显然，当维度增加或对于一些更复杂的分布而没有封闭形式的解决方案时，这将很快失控。在这里做什么最好？谢谢你。