有效地採取稀疏數據的移動平均線和蟒蛇

Question

我得到我的腳溼一些基因組分析篩選高於閾值和我有點卡住了。我有一些非常稀疏的數據，需要找到移動平均值超過某個閾值的位置，並將每個點標記爲1或0.數據是唯一類型，因此我無法使用可用的程序進行分析。

每個點表示對人類基因組的一個點（鹼基對）。對於每個數據集，都有200,000,000個潛在點。數據本質上是一個約12000個索引/值對的列表，其中所有其他點都假定爲零。我需要做的是在整個數據集中取一個移動平均值，並返回平均值超過閾值的區域。

我目前正在讀的數據集順序每一點和周圍的建築，我覺得每一個點的數組，但這是大窗口大小很慢。有沒有更高效的方法來做到這一點，也許有scipy或pandas？

編輯：下傑米的魔碼的偉大工程（但我無法給予好評還）！我非常感激。

Answer 1

你可以用numpy向量化整個事物。我已經建立的該隨機數據集（近似）12000個索引0和199999999，和隨機浮點數的0和1之間在同樣長的列表之間：

indices = np.unique(np.random.randint(2e8,size=(12000,))) 
values = np.random.rand(len(indices)) 

然後我構建總窗口大小2*win+1的索引數組圍繞每個indices，以及多少有助於通過該點的移動平均的對應陣列的：

win = 10 

avg_idx = np.arange(-win, win+1) + indices[:, None] 
avg_val = np.tile(values[:, None]/(2*win+1), (1, 2*win+1)) 

所有剩下是搞清楚重複指數和增加的貢獻的移動平均值一起：

unique_idx, _ = np.unique(avg_idx, return_inverse=True) 
mov_avg = np.bincount(_, weights=avg_val.ravel()) 

您現在可以得到指數在其中，例如列表移動平均超過0.5時，如：

unique_idx[mov_avg > 0.5] 

至於性能，第一次打開上述代碼到一個函數：

def sparse_mov_avg(idx, val, win): 
    avg_idx = np.arange(-win, win+1) + idx[:, None] 
    avg_val = np.tile(val[:, None]/(2*win+1), (1, 2*win+1)) 
    unique_idx, _ = np.unique(avg_idx, return_inverse=True) 
    mov_avg = np.bincount(_, weights=avg_val.ravel()) 
    return unique_idx, mov_avg 

這裏有一些定時幾個窗口大小，對所描述的測試數據在開始處：

In [2]: %timeit sparse_mov_avg(indices, values, 10) 
10 loops, best of 3: 33.7 ms per loop 

In [3]: %timeit sparse_mov_avg(indices, values, 100) 
1 loops, best of 3: 378 ms per loop 

In [4]: %timeit sparse_mov_avg(indices, values, 1000) 
1 loops, best of 3: 4.33 s per loop