如何利用特定步长的窗口提取短序列？

Question

下面的代码在窗口大小为4的每个序列中提取短序列。如何按步长2移动窗口并提取4个碱基对？在

示例代码

from Bio import SeqIO

with open("testA_out.fasta","w") as f:
        for seq_record in SeqIO.parse("testA.fasta", "fasta"):
            i = 0
            while ((i+4) < len(seq_record.seq)) :
              f.write(">" + str(seq_record.id) + "\n")
              f.write(str(seq_record.seq[i:i+4]) + "\n")
              i += 2

示例输入种皮在

^{pr2}$

测试输出示例

>human1
ACCC
>human1
CCGA
>human1
GATT

这个输出的问题是有一个T遗漏了，所以在这个例子中我希望也包括它。我怎样才能得到这个输出呢？用反向提取的方式，包括从头到尾提取时可能遗漏的碱基对。有人能帮我吗？在

预期输出

>human1
ACCC
>human1
CCGA
>human1
GATT
>human1
ATTT
>human1
CGAT    
>human1
CCCG

Answer 1

对于任何窗口大小和任何步长：

fasta='ACCCGATTT'
windowSize=4
step=1
i=0
while (i+windowSize)<=len(fasta):
    currentWindow=fasta[i:i+windowSize]
    print(currentWindow)
    i+=step

windowSize=4，step=2的输出：

^{pr2}$

windowSize=4，step=1的输出：

ACCC                       
CCCG
CCGA
CGAT
GATT
ATTT

最后一个与“预期输出”完全相同，排序不同。在

Answer 2

您可以将for循环与range一起使用，对range使用第三个step参数。这样，它比使用while循环要干净一些。如果数据不能除以块大小，那么最后一个块将更小。在

data = "ACCCGATTT"
step = 2
chunk = 4
for i in range(0, len(data) - step, step):
    print(data[i:i+chunk])

输出是

^{pr2}$

Answer 3

您的特定示例可以通过将步长改为1来解决。但你的问题似乎是在问，“如果序列中没有足够的字符，我如何从序列末尾以相同的窗口大小重复？”。因此，一个这样做会产生影响的例子可能是

AAAATTT

窗口大小为6，步长为2，其中AAAATT来自“前进”方向，而{}来自“反向”方向，但没有其他子序列。在

显然，向前和向后运行一次代码可以做到这一点，但它会引入重复，这通常不是一件好事。但是，您可以重构问题，以便将步骤划分为多个步骤对。在

对于一个步长为x的序列，您可以将y分为x%y和y-（x%y），然后按这些成对的步骤向前移动。（当x%y==0时跳过对的第一个成员。）

我只发布字符串处理函数，因为这些函数都不是针对基因序列的。在

^{pr2}$