时间戳

首页 > 解决方案 > 如何将多个fasta序列切片为fasta格式范围内长度的子字符串？

python - 如何将多个fasta序列切片为fasta格式范围内长度的子字符串？

问题描述

给定一个包含多个不同长度的fasta格式蛋白质序列的文件，如何生成长度为i的fasta格式蛋白质片段（子串）（i分别大于或等于5但不超过每个蛋白质的长度）？

例如，蛋白质序列文件：

>NP_12345.1
ACDEFGH
>XP_543211.2
KLMNOP
...

我想要输出fasta格式：

>NP_12345.1(1-5)
ACDEF
>NP_12345.1(1-6)
ACDEFG
>NP_12345.1(1-7)
ACDEFGH
>NP_12345.1(2-6)
CDEFG
>NP_12345.1(2-7)
CDEFGH
>NP_12345.1(3-7)    
DEFGH
>XP_543211.2(1-5)
KLMNO
>XP_543211.2(1-6)
KLMNOP
>XP_543211.2(2-6)
LMNOP
...

谁能帮助我？提前致谢。

注意：我可以使用

Seqkit sliding -s 1 -W 5

要生成只有 5 个氨基酸的蛋白质片段或肽的 fasta 格式，但如果是 6 个氨基酸的肽，我必须修改参数 -W 6。还有其他一站式解决方案吗？

标签： pythonbashperlawk

解决方案

使用 Perl，请尝试：

perl -e '
$/ = "\xd\xa";  # required if input.txt is created with DOS newline format
while (<>) {
    chomp($name = $_);
    chomp($seq = <>);
    $len = length($seq);
    for ($i = 0; $i <= $len - 5; $i++) {
        for ($j = $i + 5; $j <= $len; $j++) {
            printf("%s(%d-%d)\n", $name, $i+1, $j);
            print substr($seq, $i, $j-$i), "\n";
        }
    }
}' input.txt

产生：

>NP_12345.1(1-5)
ACDEF
>NP_12345.1(1-6)
ACDEFG
>NP_12345.1(1-7)
ACDEFGH
>NP_12345.1(2-6)
CDEFG
>NP_12345.1(2-7)
CDEFGH
>NP_12345.1(3-7)
DEFGH
>XP_543211.2(1-5)
KLMNO
>XP_543211.2(1-6)
KLMNOP
>XP_543211.2(2-6)
LMNOP

希望这可以帮助。

推荐阅读