在Python中使用numpy/scipy拟合高斯KDE

Question

我正在给一个叫“diff”的变量拟合一个高斯核密度估计器，这个变量是两个向量的差。具体的做法是这样的：gaussian_kde_covfact(diff, smoothing_param) —— 其中gaussian_kde_covfact的定义如下：

class gaussian_kde_covfact(stats.gaussian_kde):
    def __init__(self, dataset, covfact = 'scotts'):
        self.covfact = covfact
        scipy.stats.gaussian_kde.__init__(self, dataset)

    def _compute_covariance_(self):
        '''not used'''
        self.inv_cov = np.linalg.inv(self.covariance)
        self._norm_factor = sqrt(np.linalg.det(2*np.pi*self.covariance)) * self.n

    def covariance_factor(self):
        if self.covfact in ['sc', 'scotts']:
            return self.scotts_factor()
        if self.covfact in ['si', 'silverman']:
            return self.silverman_factor()
        elif self.covfact:
            return float(self.covfact)
        else:
            raise ValueError, \
                'covariance factor has to be scotts, silverman or a number'

    def reset_covfact(self, covfact):
        self.covfact = covfact
        self.covariance_factor()
        self._compute_covariance()

这个方法是有效的，但有一种特殊情况，就是当diff是一个全为0的向量时。我会遇到这样的错误：

 File "/srv/pkg/python/python-packages/python26/scipy/scipy-0.7.1/lib/python2.6/site-packages/scipy/stats/kde.py", line 334, in _compute_covariance
    self.inv_cov = linalg.inv(self.covariance)
  File "/srv/pkg/python/python-packages/python26/scipy/scipy-0.7.1/lib/python2.6/site-packages/scipy/linalg/basic.py", line 382, in inv
    if info>0: raise LinAlgError, "singular matrix"
numpy.linalg.linalg.LinAlgError: singular matrix

有没有什么办法可以解决这个问题呢？在这种情况下，我希望它返回一个密度，基本上是完全集中在0的差值上，而在其他地方没有任何质量。

谢谢。

Answer 1

一个质量集中在某一点的密度分布不是高斯分布，所以严格来说，你想做的事情是没有定义的（这样的分布没有有限的协方差）。

现在，在你的情况下，如果一个向量的所有值都是零，你可以特别处理这个情况，绕过整个系统。检测这种情况的一个简单方法是计算差值的最大值，并将其与一个很小的数（eps，具体来说是numpy.finfo(x.dtype).eps，针对向量x）进行比较。你也可以通过捕捉LinalgError来简单检测，但你需要小心区分协方差不明确的情况和零值的情况。