当前位置:首页>Linux>生信入门之 Linux Awk 用法详解

生信入门之 Linux Awk 用法详解

  • 2026-09-04 16:59:52
生信入门之 Linux Awk 用法详解

Awk 是生物信息学中最强大的文本处理工具之一

Awk 特别适合处理表格型数据(如 VCF、BED、GTF、表达矩阵等),是生信人员的必备技能。


一、Awk 基础语法

awk 'pattern {action}' file

常用内置变量:

  • $0:整行
  • $1, $2, $3...:第1、2、3列
  • NF:当前行的列数(Number of Fields)
  • NR:当前行号(Number of Record)
  • FS:输入分隔符(默认空格/Tab)
  • OFS:输出分隔符

二、生信实战案例

案例 1:提取 VCF 文件的特定信息

需求:从 VCF 文件中提取染色体、位置、参考碱基、变异碱基和质量值

# 跳过注释行(以#开头),提取前5列
awk '!/^#/ {print $1, $2, $4, $5, $6}' sample.vcf

# 进阶:只提取质量值 > 30 的变异
awk '!/^#/ && $6 > 30 {print $1, $2, $4, $5, $6}' sample.vcf

案例 2:计算 FASTA 文件的总长度和序列数

需求:统计 FASTA 文件中序列数量和总碱基数

awk '/^>/ {seqs++} !/^>/ {len += length($0)} END {print "序列数:", seqs, "总长度:", len}' genome.fa

# 输出示例:
# 序列数: 24 总长度: 3088286401

案例 3:处理基因表达矩阵

需求:计算每行(基因)在所有样本中的平均表达量

# 表达矩阵格式:GeneID Sample1 Sample2 Sample3 ...
awk 'NR>1 {sum=0; for(i=2;i<=NF;i++) sum+=$i; avg=sum/(NF-1); print $1, avg}' expression.tsv

# 输出示例:
# GeneA 15.3
# GeneB 23.7

案例 4:筛选 BED 文件中特定区域的特征

# 提取 chr1 上 1-1000000 区域内的所有基因
awk '$1=="chr1" && $2>=1 && $3<=1000000' genes.bed

# 统计每个染色体的基因数量
awk '{count[$1]++} END {for(chr in count) print chr, count[chr]}' genes.bed

案例 5:处理 GTF/GFF 注释文件

需求:提取所有蛋白质编码基因的基因名和染色体位置

# GTF格式:seqname source feature start end score strand frame attributes
awk '$3=="gene" && /protein_coding/ {
    match($0, /gene_name "([^"]+)"/, arr);
    print arr[1], $1, $4, $5, $7
}' annotation.gtf

# 输出示例:
# TP53 chr17 7661779 7687550 -
# BRCA1 chr17 43044295 43125364 +

案例 6:计算测序深度统计

需求:从 depth 文件中计算平均测序深度和覆盖度

# depth格式:chr pos depth
awk '{sum+=$3; count++; if($3>0) covered++} 
     END {
         print "平均深度:", sum/count; 
         print "覆盖度:", covered/count*100 "%"
     }' sample.depth

案例 7:合并多个样本的表达量数据

需求:将多个样本的 counts 文件按基因ID合并

# 假设有 sample1.count, sample2.count, sample3.count
# 每个文件格式:GeneID Count

awk '
FNR==NR {a[$1]=$2; next} 
{print $1, a[$1], $2}
' sample1.count sample2.count > merged.tsv

案例 8:处理 BLAST 输出结果

需求:提取 e-value < 1e-5 且 identity > 80% 的比对结果

# BLAST输出(outfmt 6):qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore
awk '$3>80 && $11<1e-5 {print $1, $2, $3, $11}' blast_result.txt

案例 9:计算 GC 含量

需求:计算 FASTA 文件中每条序列的 GC 含量

awk '/^>/ {if(seq) print header, gc/length*100; header=$0; gc=0; len=0; next} 
     {seq=$0; len+=length(seq); gc+=gsub(/[GgCc]/, "", seq)} 
     END {print header, gc/len*100}' sequences.fa

案例 10:处理多列注释文件并转换格式

需求:将 VEP 注释结果转换为简洁的表格

# 提取变异位置和主要影响
awk -F"	" 'NR==1 {print "Chr", "Pos", "Ref", "Alt", "Gene", "Effect"} 
     NR>1 {
         split($8, info, "|");
         print $1, $2, $4, $5, info[4], info[2]
     }' vep_output.vcf

三、常用技巧总结

技巧 命令示例
指定分隔符 awk -F"," '{print $1}'
多条件筛选 awk '$3>10 && $5=="PASS"'
正则匹配 awk '$2 ~ /chr[0-9]+/'
输出多列 awk '{print $1" "$3}'
去重计数 awk '{a[$1]++} END {for(i in a) print i, a[i]}'
求和/平均 awk '{sum+=$2} END {print sum, sum/NR}'

四、学习建议

  1. 从简单开始:先掌握基本的列提取和筛选
  2. 多练习:用真实的生信数据文件练习
  3. 善用 -F:根据文件格式设置正确的分隔符
  4. 结合其他工具:Awk + grep + sort + uniq 是黄金组合

💡 小贴士:Awk 是处理大规模生信数据的利器,熟练掌握后能大幅提升数据处理效率!


— 本文由 知联生信助手 整理发布 —

最新文章

随机文章