2024-06-16

李恒大神的又一杰作：bioawk让你轻松搞定生物序列分析

原创小果生信果 2023-11-30 19:00:45

收录于话题

大家好，小果今天想和大家分享一个很有用的工具，叫做bioawk。

Bioawk是生物信息学大神李恒基于awk编辑的。在此必须膜拜以下！做过基因组二代测序分析的生信人几乎都用过BWA，他就是BWA的第一作者，他还是SAMtools、MAQ、TreeSoft和TreeFam等著名生物信息软件的核心作者！

李恒大神不仅在软件开发方面有着卓越的能力，还在生物信息学理论和应用方面有着深厚的造诣。他发表了超过200篇高水平的论文，其中包括Nature、Science、Cell等顶级期刊。他的论文被引用了超过2万次，H指数达到了70。他还获得了很多荣誉和奖励，例如：中国青年科技奖、国家杰出青年科学基金、长江学者特聘教授、中国生物信息学会理事长等。

bioawk是一个基于awk的文本处理程序，专门用于处理生物信息学相关的数据格式，比如FASTA，FASTQ，GFF，SAM、BED、VCF等和带有列名的TAB分隔格式。bioawk可以让你快速地对这些数据进行过滤，转换，统计和分析，而不需要写复杂的脚本或者使用其他的软件。bioawk的语法和awk基本一致，但是增加了一些针对生物信息学数据的特性，比如可以直接读取序列的长度，质量，注释等信息，也可以方便地对序列进行反向互补，翻译等操作。bioawk还支持多种输出格式，比如CSV，TSV，JSON等，方便你将结果导入其他的程序或者平台进行进一步的分析。

接下来，小果将给大家介绍一些bioawk的基本用法和常见的应用场景。如果你对bioawk感兴趣，或者想要提高你的文本处理能力，那么请继续阅读吧。

bioawk的安装

bioawk是一个开源的软件，你可以从GitHub上克隆它的源代码，并且按照说明进行编译和安装。你也可以使用conda或者brew等包管理器来安装bioawk。

源代码下载网址：https://github.com/lh3/bioawk

安装命令：

$ git clone git://github.com/lh3/bioawk.git $ cd bioaw$ make或者connda安装：$：connda install -p ~/YOUR/conda path/ bioawk -y

bioawk的基本语法

Bioawk是一个命令行工具，它的用法是：./bioawk [-F fs] [-v var=value] [-c fmt] [-tH] [-f progfile | ‘prog’] [file …]。下面是各个选项的解释：

·-F fs：指定输入字段分隔符。默认为“空格”。

·-v var=value：在程序开始运行之前，将变量var设置为值value。

·-c fmt：指定输入格式。当使用此选项时，bioawk将根据格式或输入的第一行无缝添加命名字段的变量，具体取决于fmt。此选项还使bioawk能够读取gzip压缩的文件。

·-t：使用此选项等同于bioawk -F’t’ -v OFS=”t”。

·-H：打印输出标题行。

·-f progfile：从文件progfile中读取程序。

·’prog’：直接在命令行中指定程序。

·[file …]：指定要处理的文件。如果未指定，则从标准输入读取。

我们可以用bioawk -c -h命令列出支持的格式：

bed:     1:chrom 2:start 3:end 4:name 5:score 6:strand 7:thickstart 8:thickend 9:rgb 10:blockcount 11:blocksizes 12:blockstartssam:    1:qname 2:flag 3:rname 4:pos 5:mapq 6:cigar 7:rnext 8:pnext 9:tlen 10:seq 11:qualvcf:    1:chrom 2:pos 3:id 4:ref 5:alt 6:qual 7:filter 8:infogff:    1:seqname 2:source 3:feature 4:start 5:end 6:score 7:filter 8:strand 9:group 10:attributefastx:       1:name 2:seq 3:qual 4:comment

bioawk的基本语法和awk一样，都是`bioawk ‘pattern {action}’ file`。pattern是一个匹配条件，用于选择符合条件的记录（record）。action是一个执行动作，用于对选择的记录进行处理。file是一个输入文件，通常是一个生物信息学相关的数据文件。如果没有指定file，那么bioawk会从标准输入中读取数据。

bioawk的常见应用

下面我将给大家展示一些bioawk在生物信息学中的常见应用。这些只是一些简单的例子，你可以根据你自己的需求和想象力来发挥bioawk的更多功能。

1. 统计文件信息。

当生信小白辛辛苦苦终于拿到我们的测序结果时，会发现得到了一堆看不懂的东西/(ㄒoㄒ)/~~，如何将这一堆让人脑壳痛的数据转换成我们可以看懂的数据呢？bioawk可以实现fastx数据的统计。如果我们有一个名为test.fastq的FASTQ文件，内容如下：

如果我们想统计测序数据的基本信息，我们可以使用下面的命令来打印出每条序列的名称、GC比和长度：

$ bioawk -c fastx '{ print $name, gc($seq),length($seq) }' test.fastq

输出结果如下：

这里我们使用了`-c fastx`选项来告诉bioawk我们的输入文件是一个FASTA或者FASTQ格式的文件。这样bioawk就会自动将每条序列作为一个记录，并且提供一些内置的变量来访问序列的信息，比如$name, $seq, $qual等。我们还使用了length()函数来计算序列的长度，并且用print语句来打印出结果。

我们还可以对gff文件进行信息统计。比如说，如果我们想要计算一个GFF文件中每种特征（feature）的个数和总长度。我们可以使用下面的命令来实现：

bioawk -c gff '{count[$3]++; len[$3] += $5 - $4 + 1} END {for (f in count) print f, count[f], len[f]}' test.gff

输出结果：

这里我们使用了-c gff选项来告诉bioawk我们的输入文件是一个GFF格式的文件（gff3或者gtf为文件）。这样bioawk就会自动将每条注释作为一个记录，并且提供一些内置的变量来访问注释的信息，比如$1, $2, $3等。我们还使用了数组（array）来存储每种特征的个数和总长度，并且用END块来在处理完所有记录后打印出结果。
2. 过滤序列

我们得到的测序数据是存在参差的，为例得到更好的结果，根据一些条件来过滤掉一些不需要的序列是很必要的。比如说，我们想要过滤掉长度小于100bp或者质量低于20的序列。我们可以使用下面的命令来实现：

$ bioawk -c fastx 'length($seq) >= 100 && meanqual($qual) >= 35' test.fastq > filtered_test.fastq

这里我们使用了meanqual()函数来计算序列的平均质量，仅保留平均质量不低于35的序列，并且用&&符号来表示逻辑与（and）的关系。我们还使用了重定向符号（>）来将结果输出到一个新的文件中。
3. 转换格式

在生信分析的过程中，，我们很多时候需要将一种格式的文件转换成另一种格式的文件。比如说，我们想要将一个FASTQ文件转换成一个CSV文件，每行包含序列的名称，长度和GC含量。我们可以使用下面的命令来实现：

bioawk -c fastx -v OFS=',' '{print $name, length($seq), gc($seq)}' test.fastq > test.csv

结果如下：

这里我们使用了-v选项来设置输出字段分隔符（OFS）为逗号（,）。我们还使用了print语句来打印出结果。

Bioawk新增加了一个功能，可以直接从fasta/q文件中获得序列的反向互补序列，命令如下：

bioawk -c fastx '{ print ">"$name ORS revcomp($seq) }' test.fastq > rtest.fasta

结果如下：

在此命令中，'{ print “>”$name ORS revcomp($seq) }’ 是一个 awk 脚本，它用于提取每个记录的名称和序列，并计算序列的互补序列。$name 变量表示记录的名称，$seq 变量表示序列。revcomp 函数用于计算序列的互补序列。ORS 是输出记录分隔符，它的默认值是换行符。

Bioawk也可以从BAM文件中提取FASTA序列。还可以添加一些筛选条件，获得我们想要的高质量序列。比如如果想要得到FLAG为16的序列的方向互补序列，则使用revcomp。命令如下：

samtools view test.bam | bioawk -c sam '{s=$seq; if(and($flag, 16)) {s=revcomp($seq)} print ">"$qname"n"s}'

结果文件：

总之，bioawk是一个非常强大的文本处理工具，它可以帮助我们快速地处理各种生物信息学相关的数据格式。它的语法简单易学，功能强大，可以让我们轻松地完成各种复杂的数据分析任务。如果您是一名生物信息学工作者，那么bioawk一定是您不可或缺的工具之一。希望这篇文章能够帮助你提高你的文本处理能力，让你在生物信息学中更加游刃有余。如果你想要了解更多关于bioawk的信息，你可以访问它的GitHub页面，或者查看它的手册。如果你有任何问题或者建议，欢迎在评论区留言或者联系我。谢谢你的阅读，下次再见。

参考文档：https://github.com/lh3/bioawk/blob/master/README.md

往期推荐

1.搭建生信分析流水线，如工厂一样24小时运转Snakemake——进阶命令

2.比blast还优秀的序列比对工具？HMMER来了

3.对单细胞分析毫无头绪？让popsicleR领你入门

4.小果带你绘制ROC曲线评估生存预测能力

5.软件包安装、打怪快又好，1024G存储的生信服务器；还有比这更省钱的嘛！！！