Как я могу быстро запустить функцию над массивом каждого возможного массива длины L с заданными возможными элементами?

Question

У меня есть функция calc_dG, которая для любого массива, соответствующего короткой последовательности ДНК (от 3 до 15 оснований или около того), дает мне энергию связи этой последовательности. На самом деле, это просто поиск по массиву. nndG - это массив энергий связи для соседних пар оснований, и, таким образом, энергии связи могут быть вычислены с помощью nndG[4*S[:-1]+S[1:]] при использовании числового способа обозначения последовательностей a, g, c, t -> 0,1,2,3: это означает, что массивы многих последовательностей могут быть вычислены сразу очень быстро.

Мне нужно найти для длины L каждую последовательность, которая соответствует некоторому шаблону и приводит к значению энергии связи в определенном диапазоне.

Это очень легко сделать с помощью итераторов: просто переберите все возможные входные данные массива, вычислите энергию привязки, а затем запишите массивы, которые находятся в диапазоне. Это, однако, слишком медленно при реализации в Python (для длины 15 с 4 возможными значениями для каждого элемента существует 4 ** 15 возможных массивов и т. Д. И т. Д.). Я мог бы использовать Weave или какой-то другой метод реализации этого в C, но я бы предпочел найти простое и быстрое решение на основе массива.

Например, если каждый элемент имеет одинаковые возможные значения (например, [0,1,2,3]), то создание массива L 1D массива любой возможной длины с этими значениями можно выполнить с помощью lambda x: indices(repeat([4],L)).reshape((L,-1)).transpose(); тогда я могу просто сделать calc_dG( result ) и использовать результат [результаты, которые находятся в желаемом диапазоне], чтобы получить массивы, которые я хочу в качестве конечного результата. Это намного быстрее, чем использование итераторов Python, и, вероятно, почти так же быстро, если не быстрее, чем использование итераторов Си. К сожалению, он не работает для произвольных шаблонов, а для более длинных последовательностей не хватает памяти, так как он должен хранить каждый возможный массив в памяти перед вычислением значений.

Есть ли способ сделать все это, не прибегая к 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