稀疏数据的有效差异直方图

Question

我想计算一个数组 A 中所有元素与另一个数组 B 中所有元素之间差异的直方图。

所以我想要以下数据的直方图：

Delta1 = A1-B1
Delta2 = A1-B2
Delta3 = A1-B3
...
DeltaN = A2-B1
DeltaN+1 = A2-B2
DeltaN+2 = A2-B3 
...

此计算的重点是表明这些数据具有相关性，即使并非每个数据点在另一个数组中都有 "partner"，并且相关性在实践中相当嘈杂。

问题是这些文件实际上非常大，有好几 GB，并且向量的所有条目都是 64 位整数，差异非常大。为了能够使用相关函数和傅立叶变换来计算这些数据，将这些数据转换为二进制数组对我来说似乎是不可行的。

这里有一个小例子，可以让您更好地了解我正在看的内容。这种在 for 循环中使用 numpy 的 searchsorted 的实现相当慢。

import numpy as np
import matplotlib.pyplot as plt

timetagsA = [668656283,974986989,1294941174,1364697327,\
1478796061,1525549542,1715828978,2080480431,2175456303,2921498771,3671218524,\
4186901001,4444689281,5087334517,5467644990,5836391057,6249837363,6368090967,8344821453,\
8933832044,9731229532]


timetagsB = [13455,1294941188,1715828990,2921498781,5087334530,5087334733,6368090978,9731229545,9731229800,9731249954]

max_delta_t = 500
nbins = 10000 
histo=np.zeros((nbins,2), dtype = float) 
histo[:,0]=np.arange(0,nbins)   

for i in range(0,int(len(timetagsA))):
    delta_t = 0
    j = np.searchsorted(timetagsB,timetagsA[i]) 
    while (np.round(delta_t) < max_delta_t and j<len(timetagsB)):
        delta_t = timetagsB[j] - timetagsA[i] 

        if(delta_t<max_delta_t):
            histo[int(delta_t),1]+=1 

        j = j+1

plt.plot(histo[0:50,1])
plt.show()

如果有人可以帮助我找到一种更快的计算方法，那就太好了。提前致谢！

Answer 1

编辑

下面的解决方案假设您的数据太大以至于您不能将 np.substract 与 np.outer 一起使用，然后切片您想要保留的值：

arr_diff = np.subtract.outer(arrB, arrA)
print (arr_diff[(0<arr_diff ) &(arr_diff <max_delta_t)])
# array([ 14,  12,  10,  13, 216,  11,  13, 268], dtype=int64)

使用您的示例数据它可以工作，但不适用于太大的数据集

原始解决方案

让我们首先假设您的 max_delta_t 小于 timetagsB 中两个连续值之间的差值，以便于计算（然后我们可以尝试对其进行概括）。

#create the array instead of list
arrA = np.array(timetagsA)
arrB = np.array(timetagsB)
max_delta_t = np.diff(arrB).min() - 1 #here it's 202 just for the explanation

您可以以向量化的方式使用np.searchsorted：

# create the array of search 
arr_search = np.searchsorted(arrB, arrA) # the position of each element of arrA in arrB
print (arr_search)
# array([1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 4, 4, 4, 4, 6, 6, 6, 6, 7, 7, 7],dtype=int64)

将arrB与arr_search

切片，可以计算出arrB的元素对应arrA的每个元素的差值

# calculate the difference
arr_diff = arrB[arr_search] - arrA
print (arr_diff[arr_diff<max_delta_t]) # finc the one smaller than max_delta_t
# array([14, 12, 10, 13, 11, 13], dtype=int64)

所以你要找的是np.bincount

arr_bins = np.bincount(arr_diff[arr_diff<max_delta_t])
#to make it look like histo but not especially necessary
histo = np.array([range(len(arr_bins)),arr_bins]).T

现在的问题是，当 max_delta_t 大于两个连续值时，arrA 和 arrB 之间存在一些无法通过此方法获得的差异值在 arrB。这是一种方法，根据您的数据值，可能不是最有效的方法。对于 max_delta_t

的任何值

#need an array with the number of elements in arrB for each element of arrA 
# within a max_delta_t range
arr_diff_search = np.searchsorted(arrB, arrA + max_delta_t)- np.searchsorted(arrB, arrA)
#do a loop to calculate all the values you are interested in
list_arr = []
for i in range(arr_diff_search.max()+1):
    arr_diff = arrB[(arr_search+i)%len(arrB)][(arr_diff_search>=i)] - arrA[(arr_diff_search>=i)]
    list_arr.append(arr_diff[(0<arr_diff)&(arr_diff<max_delta_t)])

现在您可以 np.concatenate list_arr 并使用 np.bincount 例如：

arr_bins = np.bincount(np.concatenate(list_arr))
histo = np.array([range(len(arr_bins)),arr_bins]).T

稀疏数据的有效差异直方图

Efficient Histogram of Differences for sparse Data

python

numpy

histogram

sparse-matrix