Ответ 1

Многомерная нормаль теперь доступна на SciPy 0.14.0.dev-16fc0af:

from scipy.stats import multivariate_normal
var = multivariate_normal(mean=[0,0], cov=[[1,0],[0,1]])
var.pdf([1,0])

Ответ 2

Я только сделал один для своих целей, поэтому, хотя я бы поделился. Он построен с использованием "полномочий" numpy, по формуле невырожденного случая из http://en.wikipedia.org/wiki/Multivariate_normal_distribution и он проверяет ввод.

Вот код вместе с образцом прогона

from numpy import *
import math
# covariance matrix
sigma = matrix([[2.3, 0, 0, 0],
           [0, 1.5, 0, 0],
           [0, 0, 1.7, 0],
           [0, 0,   0, 2]
          ])
# mean vector
mu = array([2,3,8,10])

# input
x = array([2.1,3.5,8, 9.5])

def norm_pdf_multivariate(x, mu, sigma):
    size = len(x)
    if size == len(mu) and (size, size) == sigma.shape:
        det = linalg.det(sigma)
        if det == 0:
            raise NameError("The covariance matrix can't be singular")

        norm_const = 1.0/ ( math.pow((2*pi),float(size)/2) * math.pow(det,1.0/2) )
        x_mu = matrix(x - mu)
        inv = sigma.I        
        result = math.pow(math.e, -0.5 * (x_mu * inv * x_mu.T))
        return norm_const * result
    else:
        raise NameError("The dimensions of the input don't match")

print norm_pdf_multivariate(x, mu, sigma)

Ответ 3

В общем случае диагональной ковариационной матрицы многомерный PDF можно получить простым умножением одномерных значений PDF, возвращаемых экземпляром scipy.stats.norm. Если вам нужен общий случай, вам, вероятно, придется самому закодировать его (что не должно быть сложно).

Ответ 4

Если все еще нужно, моя реализация будет

import numpy as np

def pdf_multivariate_gauss(x, mu, cov):
    '''
    Caculate the multivariate normal density (pdf)

    Keyword arguments:
        x = numpy array of a "d x 1" sample vector
        mu = numpy array of a "d x 1" mean vector
        cov = "numpy array of a d x d" covariance matrix
    '''
    assert(mu.shape[0] > mu.shape[1]), 'mu must be a row vector'
    assert(x.shape[0] > x.shape[1]), 'x must be a row vector'
    assert(cov.shape[0] == cov.shape[1]), 'covariance matrix must be square'
    assert(mu.shape[0] == cov.shape[0]), 'cov_mat and mu_vec must have the same dimensions'
    assert(mu.shape[0] == x.shape[0]), 'mu and x must have the same dimensions'
    part1 = 1 / ( ((2* np.pi)**(len(mu)/2)) * (np.linalg.det(cov)**(1/2)) )
    part2 = (-1/2) * ((x-mu).T.dot(np.linalg.inv(cov))).dot((x-mu))
    return float(part1 * np.exp(part2))

def test_gauss_pdf():
    x = np.array([[0],[0]])
    mu  = np.array([[0],[0]])
    cov = np.eye(2) 

    print(pdf_multivariate_gauss(x, mu, cov))

    # prints 0.15915494309189535

if __name__ == '__main__':
    test_gauss_pdf()

В случае внесения будущих изменений код здесь, на GitHub

Ответ 5

Я знаю несколько пакетов python, которые используют его внутренне, с разной общностью и для разных целей, но я не знаю, предназначены ли они для пользователей.

statsmodels, например, имеет следующую скрытую функцию и класс, но не используется в statsmodels:

https://github.com/statsmodels/statsmodels/blob/master/statsmodels/miscmodels/try_mlecov.py#L36

https://github.com/statsmodels/statsmodels/blob/master/statsmodels/sandbox/distributions/mv_normal.py#L777

По существу, если вам нужна быстрая оценка, перепишите его для вашего случая использования.

Ответ 6

Я использую следующий код, который вычисляет значение logpdf, что предпочтительнее для больших размеров. Он также работает для матриц scipy.sparse.

import numpy as np
import math
import scipy.sparse as sp
import scipy.sparse.linalg as spln

def lognormpdf(x,mu,S):
    """ Calculate gaussian probability density of x, when x ~ N(mu,sigma) """
    nx = len(S)
    norm_coeff = nx*math.log(2*math.pi)+np.linalg.slogdet(S)[1]

    err = x-mu
    if (sp.issparse(S)):
        numerator = spln.spsolve(S, err).T.dot(err)
    else:
        numerator = np.linalg.solve(S, err).T.dot(err)

    return -0.5*(norm_coeff+numerator)

Код от pyParticleEst, если вы хотите, чтобы значение pdf вместо logpdf, просто введите math.exp() в возвращаемое значение

Ответ 7

Плотность можно вычислить довольно простым способом, используя функции numpy и формулу на этой странице: http://en.wikipedia.org/wiki/Multivariate_normal_distribution. Вы также можете использовать функцию правдоподобия (логарифмическая вероятность), которая менее вероятна для нижнего уровня для больших измерений и немного более простая для вычисления. Оба просто включают в себя возможность вычислить детерминант и инвертировать матрицу.

CDF, с другой стороны, является совершенно другим животным...

Ответ 8

Следующий код помог мне решить, когда задан вектор, какова вероятность того, что вектор находится в многомерном нормальном распределении.

import numpy as np
from scipy.stats import multivariate_normal

данные со всеми векторами

d= np.array([[1,2,1],[2,1,3],[4,5,4],[2,2,1]])

среднее значение данных в векторной форме, длина которых будет равна входному вектору (здесь его 3)

mean = sum(d,axis=0)/len(d)

OR
mean=np.average(d , axis=0)
mean.shape

найти коварианту векторов, которые будут иметь форму [входная векторная форма X входная векторная форма] здесь это 3x3

cov = 0
for e in d:
  cov += np.dot((e-mean).reshape(len(e), 1), (e-mean).reshape(1, len(e)))
cov /= len(d)
cov.shape

подготовка многомерного распределения Гаусса из среднего и ко-дисперсии

dist = multivariate_normal(mean,cov)

поиск функции распределения вероятностей.

print(dist.pdf([1,2,3]))

3.050863384798471e-05

Приведенное выше значение дает вероятность.