Python bitarray 反向补码

Question

我正在使用 Python 的 bitarray module 将写入二进制文件的 DNA 序列转换为其反向补码。每个核苷酸由以下格式的两位表示：

A - 00, C - 01, G - 10, T - 11.

例如，
AGCTACGG (00 10 01 11 00 01 10 10) 的反向补码将是 CCGTAGCT (01 01 10 11 00 10 01 11)。

这个序列正好占用 16 位（2 个字节），但是长度为 9 的序列将占用 18 位，并且它被填充以占用24 位（3 字节）。

目前我使用 for 循环进行转换，但这个解决方案非常慢。

def reverse_complement( my_bitarray, seq_length ):

    for i in range(0, 2 * seq_length - 1, 2):

        if my_bitarray[i] == my_bitarray[i + 1]:

            if my_bitarray[i] == 0:
                my_bitarray[i], my_bitarray[i + 1] = 1, 1

            else:
                my_bitarray[i], my_bitarray[i + 1] = 0, 0

    #padding if the bitarray is not a multiple of 8 bits in length
    if seq_length / 4 != int():
        my_bitarray.reverse()
        my_bitarray.fill()
        my_bitarray.reverse()

    return my_bitarray

a = bitarray()
a.frombytes(seq[::-1])
b = a[int(seq_start)::] # seq without padding
b.reverse()

reverse_complement(b, seq_length)

关于如何加快此过程的任何提示？

Answer 1

如果您不介意从 PyPI 安装 boltons 包，您可以执行以下操作：

from itertools import chain

from bitarray import bitarray
from boltons.iterutils import pairwise

original = bitarray('0010011100011010')
complement = ~original
reverse_complement = bitarray(chain.from_iterable(reversed(pairwise(complement))))
assert reverse_complement == bitarray('0101101100100111')

更新：

截至boltons v16.2.0，pairwise做了别的事情，所以答案应该改成使用chunked：

from boltons.iterutils import chunked
reverse_complement = bitarray(chain.from_iterable(reversed(chunked(complement, 2))))

Answer 2

您提供的代码没有给出您指定的答案。

这是给出正确答案的代码。也许它也足够快：

def reverse_complement(my_bitarray):
    # First reverse by twos
    my_bitarray = zip(my_bitarray[0::2], my_bitarray[1::2])
    my_bitarray = reversed(list(my_bitarray))
    my_bitarray = (i for t in my_bitarray for i in t)
    my_bitarray = bitarray(my_bitarray)

    # Then complement
    my_bitarray.invert()
    return my_bitarray

请注意，您不必担心填充。 bitarray.bitarray() 为您管理所有这些。