☰

for 循环在单个命令中将多个文本从父文件 grep 到多个文件

tag-icon • tag-icon shell-script text-processing grep bioinformatics
for 循环在单个命令中将多个文本从父文件 grep 到多个文件

我有 29 个 fasta 文件（.fa 作为扩展名）根据其基因命名和存储序列。

（例如：核糖体蛋白 L1、核糖体蛋白 L6P/L9E，...）

这29个fasta文件中共有722个物种。每个序列的第一行都标有其基因和物种名称，第二行则填充其序列。

1 个物种将有超过 1 个基因序列。

我想将根据基因排序的 29 个 fasta 文件中的 722 个物种转移到单独的 722 个文件中（根据物种而不是基因对它们进行排序）。

父文件中的物种名称用方括号括起来[ ]。

如何使用for循环提取722个文件并根据其序列名称命名文件？

示例来自Ribosomal Protein L1.fa：

>gi|103486926|ref|YP_616487.1| 50S ribosomal protein L1 [Sphingopyxis alaskensis RB2256]
MAKLTKKQKALEGKVDAQKLHGVDEAIKLVRELATAKFDETLEIAMNLGVDPRHADQMVRGVVTLPAGTGKDVKVAVFAR

示例来自Ribosomal Protein L6PL9E.fa：

>gi|410479108|ref|YP_006766745.1| ribosomal protein L6P/L9E [Leptospirillum ferriphilum ML-04]
MGFTHTVEFTLPSLIKASIEKQTIITLSSPDKELLGQFAADVRSIRPPEPYKGKGIKYSGEKILRKEGKTGKK

对于第一个例子，

种名：Sphingopyxis alaskensis RB2256

基因序列：MAKLTKKQKALEGKVDAQKLHGVDEAIKLVRELATAKFDETLEIAMNLGVDPRHADQMVRGVVTLPAGTGKDVKVAVFA

我想将文件命名为Sphingopyxis alaskensis RB2256.fa并将具有该物种名称的所有序列插入到该文件中。

我正在使用 bash shell 来执行此操作。我可以用来grep完成事情：

grep -A+1 "Sphingopyxis alaskensis RB2256" *.fa >> Sphingopyxis alaskensis RB2256.fa

但我需要执行 722 次才能根据物种对序列进行排序。

for循环中的grep可以用来简化工作吗？或者有其他方法可以做到这一点？

答案1

Fasta 格式不要求所有序列都在一行上。事实上，这种情况并不常见，因为大多数生物序列都很长。因此，grep在任何情况下，如果 ID 的序列超过一行，您都会失败。此外，您的grep命令将创建一个名为 的文件Sphingopyxis，而不是一个名为Sphingopyxis alaskensis RB2256.fa.

无论如何，您可以执行类似的操作，将每个序列放入物种后的文件名中：

awk -F'[][]' '/>/{n=$2}; {print >> n".fa"}' *.fa

但是，我强烈建议您不要在文件名中使用空格，因为这只会让您的生活变得更加困难。更安全的方法是：

awk -F'[][]' '/>/{n=$2; gsub(/ /,"_",n)}; {print >> n".fa"}' *.fa

将gsub物种名称中的所有空格替换为_，从而生成以下文件：

Leptospirillum_ferriphilum_ML-04.fa  Sphingopyxis_alaskensis_RB2256.fa

请注意，上述两种方法都可以处理多行序列。

相关内容