python - 我怎样才能在给定可能元素的每个可能的长度 L 数组的数组上快速运行一个函数？

Question

我有一个函数 calc_dG，对于对应于短 DNA 序列(3 到 15 个碱基左右)的任何阵列，它会给出该序列的结合能。实际上，它只是一个数组查找。 nndG 是相邻碱基对的结合能数组，因此当使用 a 时，结合能可以用 nndG[4*S[:-1]+S[1:]] 计算,g,c,t -> 0,1,2,3 数字表示序列的方式:这意味着可以在 numpy 中非常快速地一次计算许多序列的数组。

我需要找到长度为 L 的每个序列，这些序列既适合某个模板又会产生特定范围内的结合能值。

使用迭代器很容易做到这一点:只需遍历每个可能的数组输入，计算结合能，然后记录范围内的数组。然而，这在 Python 中实现时太慢了(对于长度为 15 且每个元素有 4 个可能值的数组，有 4**15 个可能的数组，等等)。我可以使用 Wea​​ve 或其他一些方法在 C 中实现它，但我更愿意找到一个简单快速的基于数组的解决方案。

例如，如果每个元素都具有相同的可能值(例如 [0,1,2,3])，则可以使用 lambda 生成具有这些值的每个可能长度 L 一维数组的数组x: 索引(重复([4],L)).reshape((L,-1)).transpose();然后我就可以做 calc_dG( result )，并使用 result[results that are in the desired range] 得到我想要的数组作为最终结果。这比使用 Python 迭代器快得多，并且可能几乎与使用 C 迭代器一样快，甚至更快。不幸的是，它不适用于任意模板，对于较长的序列，它会耗尽内存，因为它必须在计算值之前将每个可能的数组存储在内存中。

有没有办法在不求助于 C 的情况下完成所有这些工作？

Answer 1

如果我对你的问题的理解正确，那么你正在最大化一个函数 f(i_1, i_2, ..., i_n) 在集合 {0, 1, 2, 3} 中的整数。

您可以结合使用迭代和矢量化索引。

import numpy as np
import itertools

def cartesian_chunked(n, n_items=4, chunk_dim=3):
    if n > chunk_dim:
        p = n - chunk_dim
        q = chunk_dim
        outer = itertools.product(*([range(n_items)] * (n - chunk_dim)))
    else:
        p = 0
        q = n
        def outer_iter():
            yield ()
        outer = outer_iter()

    chunk = np.zeros([n_items**q, n], dtype=int)
    chunk[:,p:] = np.indices(np.repeat([n_items], q)).reshape(q, -1).T
    for seq in outer:
        chunk[:,:p] = seq
        yield chunk

def compute_energy(indices):
    base_energies = np.array([-1, 4, 8, 2.4])
    return (base_energies[indices]).sum(axis=1) 

max_energy = 0
max_config = None

# try out 4**10 ~ 1e6 combinations, in chunks of 4**8
for chunk in cartesian_chunked(n=10, n_items=4, chunk_dim=8):
    energies = compute_energy(chunk)
    j = np.argmax(energies)
    if energies[j] > max_energy:
        max_energy = energies[j]
        max_config = chunk[j].copy() # copy! the chunk is modified

print max_energy
print max_config