python中有效地提取稀疏数据的移动平均和阈值以上的过滤

2024-03-29 13:57:29 发布

您现在位置:Python中文网/ 问答频道 /正文

我正忙着做一些基因组分析,有点卡住了。我有一些非常稀疏的数据,需要找到移动平均值超过某个阈值的地方,将每个点标记为1或0。数据的类型是唯一的,因此我无法使用可用的程序进行分析。在

每个点代表人类基因组上的一个点(碱基对)。对于每个数据集,有200000个潜在点。数据本质上是一个约12000个索引/值对的列表,其中所有其他点都假定为零。我需要做的是在整个数据集中取一个移动平均值,并返回平均值高于阈值的区域。在

我目前正在从数据集中按顺序读取每个点,并围绕找到的每个点构建一个数组,但对于大窗口大小,这是非常缓慢的。有没有一种更有效的方法来做到这一点,也许是用scipy或熊猫?在

编辑:杰米下面的魔术代码很好用(但我还不能投赞成票)!我非常感激。在


Tags: 数据标记程序区域类型基因组列表顺序
1条回答
网友
1楼 · 发布于 2024-03-29 13:57:29

你可以用numpy来矢量化整件事。我构建了一个随机数据集,包含0到1999999之间的12000个索引,以及一个同样长的0到1之间的随机浮动列表:

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

然后,我在每个indices周围构造一个总窗口大小2*win+1的索引数组,以及一个相应的数组,其中包含该点对移动平均值的贡献:

^{pr2}$

剩下的就是计算出重复指数,并将对移动平均线的贡献相加:

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

相关问题 更多 >