bioawk是由生信大牛李恒用C开发，不仅涵盖了awk的功能，在此基础上添加了很多额外的生物大数据文件的支持，如：BED、 GFF、 SAM、VCF、FASTA/Q等等。

安装：conda install bioawk

下面简单介绍其使用方法：

基础应用：


# 查看帮助
bioawk [command] help
# 如：bioawk -c help

# 输出：
bed:
	1:chrom 2:start 3:end 4:name 5:score 6:strand 7:thickstart 8:thickend 9:rgb 10:blockcount 11:blocksizes 12:blockstarts 
sam:
	1:qname 2:flag 3:rname 4:pos 5:mapq 6:cigar 7:rnext 8:pnext 9:tlen 10:seq 11:qual 
vcf:
	1:chrom 2:pos 3:id 4:ref 5:alt 6:qual 7:filter 8:info 
gff:
	1:seqname 2:source 3:feature 4:start 5:end 6:score 7:filter 8:strand 9:group 10:attribute 
fastx:
	1:name 2:seq 3:qual 4:comment

1）提取未mapping上的Reads（不含文件头）：bioawk -c sam 'and($flag,4)' aln.sam.gz

2）提取mapping上的Reads：bioawk -Hc sam '!and($flag,4)'

3）得到反向互补序列：bioawk -c fastx '{print ">"$name;print revcomp($seq)}' seq.fa.gz

4）从SAM创建FASTA序列


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

5）输出样本foo以及bar的基因型：grep -v ^## in.vcf | bioawk -tc hdr '{print $foo,$bar}'

更多应用：

1）统计Reads：bioawk -t -c fastx 'END {print NR}' test.fastq

2）统计GC含量：bioawk -c fastx '{print $name, gc($seq)}' test.fa

3）根据ID提取Fasta文件中相应的序列信息：


$ cat test.fa
>blah_C1
ACTGTCTGTC
ACTGTGTTGTG
ATGTTGTGTGTG
>blah_C2
ACTTTATATATT
>blah_C3
ACTTATATATATATA
>blah_C4
ACTTATATATATATA
>blah_C5
ACTTTATATATT

$ cat IDs.txt
blah_C4
blah_C5

命令：bioawk -cfastx 'BEGIN{while((getline k <"IDs.txt")>0)i[k]=1}{if(i[$name])print ">"$name"\n"$seq}' test.fa

结果如下：


>blah_C4
ACTTATATATATATA
>blah_C5
ACTTTATATATT

最后补充一句，awk的功能和语法biowak也是支持的。

参考资料：

1.https://github.com/lh3/bioawk

2.https://www.biostars.org/p/127141/

阅读: 770