2015-04-18 5 views
1

Когда я пытаюсь вычислить расстояние Махаланобиса со следующим кодом python, я получаю некоторые записи Nan в результате. У вас есть представление о том, почему это происходит? Мой data.shape = (181, 1500)Scipy - Nan при расчете расстояния Mahalanobis

from scipy.spatial.distance import pdist, squareform 

data_log = log2(data + 1) # A log transform that I usually apply to my data 
data_centered = data_log - data_log.mean(0) # zero centering 
D = squareform(pdist(data_centered, 'mahalanobis')) 

Я также попытался:

data_standard = data_centered/data_centered.std(0, ddof=1) 
D = squareform(pdist(data_standard, 'mahalanobis')) 

Также получил пренебрежимо малых. Ввод не поврежден, и другие расстояния, такие как расстояние корреляции, могут быть рассчитаны просто отлично. По какой-то причине, когда я уменьшаю количество функций, я перестаю получать Nans. Например следующие примеры не получает Nan:

D = squareform(pdist(data_centered[:,:200], 'mahalanobis')) 
D = squareform(pdist(data_centered[:,180:480], 'mahalanobis')) 

в то время как те, другие получают Nans:

D = squareform(pdist(data_centered[:,:300], 'mahalanobis')) 
D = squareform(pdist(data_centered[:,180:600], 'mahalanobis')) 

Любые подсказки? Является ли это ожидаемым поведением, если какое-либо условие для ввода не выполняется?

+1

Вы можете поделиться куском с вашими данными? – marmeladze

+0

0,9,0,0,0,1,0,0,0,2,0,0,0,0,0,0,0,2,7,0,0,1,0,0,1 , 0,0,0,0,0,0,0,0,12,0,0,1,0,0,0 0,9,0,0,0,0,0,0,0,0 , 0,0,0,0,2,0,0,0,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,7,2 , 0,1,0,0,0 0,5,1,0,0,1,5,0,0,2,0,0,0,0,8,0,0,0,4,0 , 0,0,7,4,0,0,0,0,0,0,0,1,0,9,0,1,3,0,0,1 0,4,0,0,0 , 0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,1,2,0,0,0 , 0,0,0,5,0,0,0,0,0,0 0,13,0,0,0,0,1,0,0,0,0,2,6,0,3 , 0,0,0,6,0,1,0,0,3,2,0,0,0,0,0,0,0,0,5,4,0,0,0,0,0 1,17,0,0,5,2,0,0,0,0,0,0,7,0,0,0,0,0,1,1,0,0,0,1,1 , 0,0,0,0,1,0,0,0,5,0,0,3,2,0,0 0,8,0,0,3,0,0,0,0,0 , 0,2,10,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0 , 0,0,0,6,0 – Gioelelm

+0

Это довольно маленькое окно, которое соответствует комментариям. Как бы вы хотели, чтобы я разделил его? – Gioelelm

ответ

2

У вас меньше наблюдений, чем признаков, поэтому матрица ковариации V, вычисленная кодом scipy, является единственной. Код не проверяет это и слепо вычисляет «обратную» ковариационную матрицу. Поскольку этот численно рассчитанный обратный в основном мусор, продукт (x-y)*inv(V)*(x-y) (где x и y - это наблюдения) может оказаться отрицательным. Тогда квадратный корень этого значения приводит к nan.

Например, этот массив также приводит к nan:

In [265]: x 
Out[265]: 
array([[-1. , 0.5, 1. , 2. , 2. ], 
     [ 2. , 1. , 2.5, -1.5, 1. ], 
     [ 1.5, -0.5, 1. , 2. , 2.5]]) 

In [266]: squareform(pdist(x, 'mahalanobis')) 
Out[266]: 
array([[ 0.  ,   nan, 1.90394328], 
     [  nan, 0.  ,   nan], 
     [ 1.90394328,   nan, 0.  ]]) 

Вот расчет Махаланобиса сделано "вручную":

In [279]: V = np.cov(x.T) 

В теории, V вырождена; следующее значение эффективно 0:

In [280]: np.linalg.det(V) 
Out[280]: -2.968550671342364e-47 

Но inv не видит проблему, и возвращает обратное:

In [281]: VI = np.linalg.inv(V) 

Давайте вычислить расстояние между x[0] и x[2] и убедиться, что мы получаем то же самое не-нан значение (1,9039), возвращаемый pdist, когда мы используем VI:

In [295]: delta = x[0] - x[2] 

In [296]: np.dot(np.dot(delta, VI), delta) 
Out[296]: 3.625 

In [297]: np.sqrt(np.dot(np.dot(delta, VI), delta)) 
Out[297]: 1.9039432764659772 

Вот что происходит, когда мы пытаемся вычислить расстояние между x[0] и x[1]:

In [300]: delta = x[0] - x[1] 

In [301]: np.dot(np.dot(delta, VI), delta) 
Out[301]: -1.75 

Тогда квадратный корень из этой величины дает nan.


В SciPy 0,16 (будет выпущен в июне 2015 года), вы получите сообщение об ошибке вместо нан или мусора.Сообщение об ошибке описывает проблему:

In [4]: x = array([[-1. , 0.5, 1. , 2. , 2. ], 
    ...:  [ 2. , 1. , 2.5, -1.5, 1. ], 
    ...:  [ 1.5, -0.5, 1. , 2. , 2.5]]) 

In [5]: pdist(x, 'mahalanobis') 
--------------------------------------------------------------------------- 
ValueError        Traceback (most recent call last) 
<ipython-input-5-a3453ff6fe48> in <module>() 
----> 1 pdist(x, 'mahalanobis') 

/Users/warren/local_scipy/lib/python2.7/site-packages/scipy/spatial/distance.pyc in pdist(X, metric, p, w, V, VI) 
    1298          "singular. For observations with %d " 
    1299          "dimensions, at least %d observations " 
-> 1300          "are required." % (m, n, n + 1)) 
    1301     V = np.atleast_2d(np.cov(X.T)) 
    1302     VI = _convert_to_double(np.linalg.inv(V).T.copy()) 

ValueError: The number of observations (3) is too small; the covariance matrix is singular. For observations with 5 dimensions, at least 6 observations are required. 
Смежные вопросы