Python 位数组逆补
Python bitarray reverse complement
我正在使用 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)
关于如何加快此过程的任何提示?
如果您不介意从 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))))
您提供的代码没有给出您指定的答案。
这是给出正确答案的代码。也许它也会足够快:
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()
为您管理所有这些。
我正在使用 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)
关于如何加快此过程的任何提示?
如果您不介意从 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))))
您提供的代码没有给出您指定的答案。
这是给出正确答案的代码。也许它也会足够快:
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()
为您管理所有这些。