将多个 FASTA 文件拆分为单独的文件，保留其原始名称

问题描述 投票：0回答：5

我正在尝试使用之前发布在该论坛上的 AWK 脚本。我正在尝试将包含多个 DNA 序列的大型 FASTA 文件拆分为单独的 FASTA 文件。我需要将每个序列分离到它自己的 FASTA 文件中，每个新的 FASTA 文件的名称需要是原始大型 multifasta 文件中的 DNA 序列的名称（> 之后的所有字符）。

我尝试了我在 stackoverflow 上找到的这个脚本：

awk '/^>chr/ {OUT=substr($0,2) ".fa"}; OUT {print >OUT}' your_input

效果很好，但 DNA 序列直接在文件名之后开始——没有空格。 DNA 序列需要从一个新行开始（常规 FASTA 格式）。

如果能帮我解决这个问题，我将不胜感激。 谢谢！

awk split sequence fasta
5个回答
4
投票

你的意思是这样的吗？

awk '/^>chr/ {OUT=substr($0,2) ".fa";print " ">OUT}; OUT{print >OUT}' your_input

为每个“染色体/序列/事物”创建的新文件在哪里开始时有一个空行？

0
投票

我认为这应该有效。

awk '/^>chr/ {OUT=substr($0,2) ".fa"}; {print >> OUT; close(OUT)}' Input_File
0
投票

希望这个 perl 脚本能有所帮助。

#!/usr/bin/perl
open (INFILE, "< your_input.fa")
or die "Can't open file";
while (<INFILE>) {
$line = $_;
chomp $line;
if ($line =~ /\>/) { #if has fasta >
close OUTFILE;
$new_file = substr($line,1);
$new_file .= ".fa";
open (OUTFILE, ">$new_file")
or die "Can't open: $new_file $!";
}
print OUTFILE "$line\n";
}
close OUTFILE;
0
投票

.fa（或 .fasta）格式如下：

>ID1
序列
>ID2
序列

拆分 fasta 文件时，实际上不希望在其顶部插入换行符。因此Pramod的回答更为恰当。此外，可以更一般地定义 ID 以仅匹配 

>
字符。因此，完整的行将是：

awk '/^>/ {OUT=substr($0,2) ".fa"}; {print >> OUT; close(OUT)}' Input_File

如果您不想将所有拆分文件弄乱当前目录，您还可以输出到子目录中 (

subdir
)：

awk '/^>/ {OUT="subdir/" substr($0,2) ".fa"}; {print >> OUT; close(OUT)}' Input_File
0
投票

awk
将多序列 fasta 文件拆分为单独的序列文件

这个问题最好通过考虑每个序列（完整的标题）一个单一的

record
和改变
awk
的默认记录分隔符
RS
（通常是一个换行符）是唯一的（每个记录一个）
>
用于定义标题开始的符号。由于我们希望使用标题文本作为文件名，并且由于 fasta 标题不能包含换行符，因此将 
awk
的默认字段分隔符或 
FS
（通常为空白）重置为换行符也很方便.

这两个都是在一个

awk
BEGIN
块中完成的：

BEGIN{RS=">";FS="\n"}

由于文件以

>
开头，第一条记录将为空，因此必须忽略，以防止尝试写入从空记录中提取的文件名而导致的错误。因此，主要 
awk
操作块被过滤为仅处理以记录号 (
NR
) 2 开头的记录。这是通过在操作块之前放置一个条件来实现的，如下所示：

NR>1{ ... }

记录分隔符设置为

>
每条记录都是一个完整的序列，包括它的标题，每条记录在换行符处被分成字段（因为我们将字段分隔符设置为“ ")。因此，每条记录的字段 1 (
$1
) 包含我们希望用作文件名的文本。注意记录分隔符 (>) 不再是任何字段的一部分，因此整个第一个字段可用于构建文件名。在此示例中，“.fasta”已附加为文件扩展名：

fnme=$1 ".fasta";

接下来，打印 fasta 标头标记“>”，然后是整个记录（

$0
）到刚刚形成的文件名
fnme
，使用
awk
的
>
重定向：

print ">" $0 > fnme;

最后，关闭文件以防止 awk 超过允许打开文件数的系统限制，如果要写入许多文件（参见脚注）：

close(fnme);

整个过程

awk 命令

awk 'BEGIN{RS=">";FS="\n"} NR>1{fnme=$1".fasta"; print ">" $0 > fnme; close(fnme);}' example.fasta

在以下名为example.fasta

的模拟文件上测试
：

>DNA sequence 1
GCAAAAGAACCGCCGCCACTGGTCGTGAAAGTGGTCGATCCAGTGACATCCCAGGTGTTGTTAAATTGAT
CATGGGCAGTGGCGGTGTAGGCTTGAGTACTGGCTACAACAACACTCGCACTACCCGGAGTGATAGTAAT
GCCGGTGGCGGTACCATGTACGGTGGTGAAGT
>DNA sequence 2
TCCCAGCCAGCAGGTAGGGTCAAAACATGCAAGCCGGTGGCGATTCCGCCGACAGCATTCTCTGTAATTA
ATTGCTACCAGCGCGATTGGCGCCGCGACCAGGATCCTTTTTAACCATTTCAGAAAACCATTTGAGTCCA
TTTGAACCTCCATCTTTGTTC
>DNA sequence 3
AACAAAAGAATTAGAGATATTTAACTCCACATTATTAAACTTGTCAATAACTATTTTTAACTTACCAGAA
AATTTCAGAATCGTTGCGAAAAATCTTGGGTATATTCAACACTGCCTGTATAACGAAACACAATAGTACT
TTAGGCTAACTAAGAAAAAACTTT

results（终端命令和输出）

$ ls
'DNA sequence 1.fasta'  'DNA sequence 3.fasta'
'DNA sequence 2.fasta'   example.fasta

$ cat DNA\ sequence\ 1.fasta 
>DNA sequence 1
GCAAAAGAACCGCCGCCACTGGTCGTGAAAGTGGTCGATCCAGTGACATCCCAGGTGTTGTTAAATTGAT
CATGGGCAGTGGCGGTGTAGGCTTGAGTACTGGCTACAACAACACTCGCACTACCCGGAGTGATAGTAAT
GCCGGTGGCGGTACCATGTACGGTGGTGAAGT

$ cat DNA\ sequence\ 2.fasta 
>DNA sequence 2
TCCCAGCCAGCAGGTAGGGTCAAAACATGCAAGCCGGTGGCGATTCCGCCGACAGCATTCTCTGTAATTA
ATTGCTACCAGCGCGATTGGCGCCGCGACCAGGATCCTTTTTAACCATTTCAGAAAACCATTTGAGTCCA
TTTGAACCTCCATCTTTGTTC

$ cat DNA\ sequence\ 3.fasta 
>DNA sequence 3
AACAAAAGAATTAGAGATATTTAACTCCACATTATTAAACTTGTCAATAACTATTTTTAACTTACCAGAA
AATTTCAGAATCGTTGCGAAAAATCTTGGGTATATTCAACACTGCCTGTATAACGAAACACAATAGTACT
TTAGGCTAACTAAGAAAAAACTTT

脚注

“在同一个 awk 程序中连续写入大量文件。如果文件没有关闭，awk 最终可能会超过系统对一个进程中打开的文件数的限制。最好在每个进程关闭时关闭每个文件程序写完了。”

引自https://www.gnu.org/software/gawk/manual/html_node/Close-Files-And-Pipes.html

最新问题
© www.soinside.com 2019 - 2024. All rights reserved.