首页 > 解决方案 > 在计算多元高斯的pdf时如何利用矢量化?

问题描述

我花了几个小时在谷歌上搜索这个问题,似乎找不到任何信息。

我尝试将多元高斯 pdf 编码为:

def multivariate_normal(X, M, S):

    # X has shape (D, N) where D is the number of dimensions and N the number of observations
    # M is the mean vector with shape (D, 1)
    # S is the covariance matrix with shape (D, D)

    D = S.shape[0]
    S_inv = np.linalg.inv(S)
    logdet = np.log(np.linalg.det(S))
    log2pi = np.log(2*np.pi)
    devs = X - M

    a = np.array([- D/2 * log2pi - (1/2) * logdet - dev.T @ S_inv @ dev for dev in devs.T])

return np.exp(a)

我只通过for循环成功地计算了pdf,迭代了N次。如果我不这样做,我最终会得到一个无用的 (N, N) 矩阵。我在这里找到了另一个帖子,但该帖子已经过时并且在 matlab 中。

有没有办法利用 numpy 的矢量化?

这是我在stackoverflow上的第一篇文章,如果有什么问题请告诉我!d

标签: pythonnumpyvectorizationgaussian

解决方案


我以类似的方式遇到了这个问题,这就是我解决它的方法:

变量:

  • X = numpy.ndarray[numpy.ndarray[float]] - mxn
  • MU = numpy.ndarray[numpy.ndarray[float]] - kxn
  • SIGMA = numpy.ndarray[numpy.ndarray[numpy.ndarray[float]]] - kxnxn
  • k = 整数

其中 X 是我的特征向量,MU 是我的均值,SIGMA 是我的协方差矩阵。

为了矢量化,我根据点积的定义重写了点积:

sigma_det = np.linalg.det(sigma)
sigma_inv = np.linalg.inv(sigma)
const = 1/((2*np.pi)**(n/2)*sigma_det**(1/2))
p = const*np.exp((-1/2)*np.sum((X-mu).dot(sigma_inv)*(X-mu),axis=1))

推荐阅读