Mafft 只用 Python 创建一个文件

Question

所以我正在做一个项目来对齐序列 ID 及其代码。我得到了一个条形码文件，其中包含 DNA 序列的标签，即 TTAGG。有几个标签（ATTAC、ACCAT 等）然后从序列文件中删除并放置在序列 ID 中。 示例：

sequence file --> SEQ 01  TTAGGAACCCAAA
barcode file --> TTAGG

我想要的输出文件将删除条形码并使用它来创建一个新的 fasta 格式文件。 例子： testfile.TTAGG 打开时应该有

>SEQ01
AACCCAAA

这些文件有好几个。我想获取我创建的每个文件并通过mafft 运行它们，但是当我运行我的脚本时，它只专注于mafft 的一个文件。我上面提到的文件正常，但是当mafft 运行时，它只运行最后创建的文件。

这是我的脚本：

#!/usr/bin/python

import sys
import os

fname = sys.argv[1]
barcodefname = sys.argv[2]

barcodefile = open(barcodefname, "r")
for barcode in barcodefile:
        barcode = barcode.strip()
        outfname = "%s.%s" % (fname, barcode)
        outf = open(outfname, "w+")
        handle = open(fname, "r")
        mafftname = outfname + ".mafft"
        for line in handle:
                newline = line.split()
                seq = newline[0]
                brc = newline[1]
                potential_barcode = brc[:len(barcode)]
                if potential_barcode == barcode:
                        outseq = brc[len(barcode):]
                        barcodeseq = ">%s\n%s\n" % (seq,outseq)
                        outf.write(barcodeseq)
        handle.close()
        outf.close()
cmd = "mafft %s > %s" % (outfname, mafftname)
os.system(cmd)
barcodefile.close()

我希望这已经足够清楚了！请帮忙！我尝试更改缩进，在关闭文件时进行调整。大多数情况下，它根本不会生成 .mafft 文件，有时它会生成但不添加任何内容，但大多数情况下它只适用于最后创建的文件。

示例： 代码的开头创建了诸如 -

之类的文件

testfile.ATTAC
testfile.AGGAC
testfile.TTAGG

然后当它运行mafft 时，它只会创建 testfile.TTAGG.mafft（输入正确）

我尝试关闭 outf 文件，然后再次打开它，它告诉我我在强制它。 我已将outf 文件更改为只写，不会更改任何内容。

Answer 1

mafft 之所以只对齐最后一个文件文件，是因为它的执行是在循环之外的。

就您的代码而言，您在循环的每次迭代中创建一个输入文件名变量 (outfname)，但该变量总是在下一次迭代中被覆盖。因此，当您的代码最终到达 mafft 执行命令时，outfname 变量将包含循环的最后一个文件名。

要纠正这个问题，只需在循环中插入 mafft 执行命令：

#!/usr/bin/python

import sys
import os

fname = sys.argv[1]
barcodefname = sys.argv[2]

barcodefile = open(barcodefname, "r")
for barcode in barcodefile:
        barcode = barcode.strip()
        outfname = "%s.%s" % (fname, barcode)
        outf = open(outfname, "w+")
        handle = open(fname, "r")
        mafftname = outfname + ".mafft"
        for line in handle:
                newline = line.split()
                seq = newline[0]
                brc = newline[1]
                potential_barcode = brc[:len(barcode)]
                if potential_barcode == barcode:
                        outseq = brc[len(barcode):]
                        barcodeseq = ">%s\n%s\n" % (seq,outseq)
                        outf.write(barcodeseq)
        handle.close()
        outf.close()
        cmd = "mafft %s > %s" % (outfname, mafftname)
        os.system(cmd)

barcodefile.close()