在 Linux 中从 Fasta 文件中提取列

在 Linux 中从 Fasta 文件中提取列

我有一个 fasta 文件,如下所示

>ENST00000632684.1 cdna chromosome:GRCh38:7:142786213:142786224:1 gene:ENSG00000282431.1 gene_biotype:TR_D_gene transcript_biotype:TR_D_gene gene_symbol:TRBD1 description:T cell receptor beta diversity 1 [Source:HGNC Symbol;Acc:HGNC:12158]
GGGACAGGGGGC
>ENST00000434970.2 cdna chromosome:GRCh38:14:22439007:22439015:1 gene:ENSG00000237235.2 gene_biotype:TR_D_gene transcript_biotype:TR_D_gene gene_symbol:TRDD2 description:T cell receptor delta diversity 2 [Source:HGNC Symbol;Acc:HGNC:12255]
CCTTCCTAC

我想提取gene_symbol 和描述。但不幸的是,描述之间有空格,我无法提取完整的描述。

我试过这个

cat Homo_sapiens.GRCh38.cdna.all.fa | grep ">" | cut -f 7,8 -d" "  > Human_Annotations

但这给了我这样的输出,其中描述已被破坏。

gene_symbol:TRBD1 description:T
gene_symbol:TRDD2 description:T

我想要这样的输出

TRBD1 T cell receptor beta diversity 1
TRDD2 T cell receptor delta diversity 2

答案1

使用awk

awk -F ':' '/^>/ { sub(" .*",    "", $10)
                   sub(" \\[.*", "", $11)
                   print $10, $11 }' file.fa

您想要提取的数据是[每个标题行的第 10 个字段中的第一个单词以及直到第 11 个字段中的所有内容(如果这些字段是:分隔的)。

该代码会删除第 10 个字段中第一个空格以及[第 11 个字段中之后的所有内容(包括[前面的空格)。

然后打印修改后的第 10 和 11 字段。

给出问题中数据的输出:

TRBD1 T cell receptor beta diversity 1
TRDD2 T cell receptor delta diversity 2

答案2

尝试这样的事情:

cat ... | sed -n '/^>/ { s/.*description: *//; s/\[.*//; p; }'

(未经测试,因为我在移动设备上。)

还有更优雅的方式;例如 Awk 循环是最灵活的。

相关内容