2017-01-02 130 views
1

所以我有一組N個多元正態分佈,它們都具有相同的協方差。對於這些分佈中的每一個,我想計算得到值x的可能性。Python:計算N個多元正態分佈值的可能性

對於單個分佈情況,以及「X」值倍數,這是微不足道現在

from scipy.stats import multivariate_normal 
import numpy as np 

cov = [[1 ,0.1],[0.1 ,1]] 
mean = [0,0] 
Values = np.random.multivariate_normal([0,0],cov,samp) 
print multivariate_normal.pdf(Values, mean, cov) 

,如果我們扭轉這一點,假設我們只有一個值來檢查,但多個手段,每次都有相同的協方差。如下(當然在實際情況下,每次迭代的平均值不同)

means = [mean]*samples 
Value = Values[0,:] 

L = [] 
for iMean in means: 
    L.append(multivariate_normal.pdf(Value, iMean, cov)) 

print L 

有沒有更好的方法來做到這一點?如果存在任何差異,那麼假設協方差矩陣不相關也是允許的,儘管通常的解決方案是優選的。

回答

1

您可以先計算所有分佈的平方馬氏距離。 https://en.wikipedia.org/wiki/Mahalanobis_distance

然後你計算概率密度。

* https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.stats.multivariate_normal.html * https://en.wikipedia.org/wiki/Multivariate_normal_distribution

使用numpy的陣列可以避開緩慢蟒蛇循環。 我將此添加到您的示例中:

from scipy.stats import multivariate_normal 
import numpy as np 

cov = [[1 ,0.5],[0.5 ,1]] 
mean = [2,2] 

samples = 10 
means = [mean]*samples 

Value = (3,2.5) 

L = [] 
for iMean in means: 
    L.append(multivariate_normal.pdf(Value, iMean, cov)) 



mean_array = np.array(means) 
value_array = np.array(Value).astype(np.float) 
cov_array = np.array(cov) 
inv_cov_array = np.linalg.inv(cov_array) 
dim = cov_array.shape[0] 

diffs = value_array-mean_array 
maha_distances = np.sum(diffs.transpose()*np.dot(inv_cov_array,diffs.transpose()),axis=0)  
denominator = 1/np.sqrt((2*np.pi)**dim*np.linalg.det(cov_array)) 

l = denominator * np.exp(-0.5*maha_distances) 

res_dif = np.array(L) - l 
print res_dif 
+0

對不起,遲到接受!很好的答案! – Dammi