歡迎關注生信寶典:http://mp.weixin.qq.com/s/8wD14FXt7fLDo1BjJyT0ew Linux學習系列文章是生信寶典最開始主推的一塊,力圖從一個新額視角幫助初學者快速入門Linux系統,熟悉Linux下的文件和目錄,文件操作, 文件內容操作。而且教程摒棄了完美操作,列舉出常見錯誤和解決方式,管道、標準輸入輸出解惑Linux下多種信息輸出方式。
在文件排序和FASTA文件操作中簡述了awk和sed的使用,作為一個引子。本篇則詳細列舉關于awk常用的操作和一些偏門的操作。
awk基本參數解釋
awk擅長于對文件按行操作,每次讀取一行,然后進行相應的操作。
awk讀取單個文件時的基本語法格式是awk 'BEGIN{OFS=FS="\t"}{print $0, $1;}' filename。
讀取多個文件時的語法是awk 'BEGIN{OFS=FS="\t"}ARGIND==1{print $0, $1;}ARGIND==2{}' file1 file2。
awk后面的命令部分是用引號括起來的,可以單引號,可以雙引號,但注意不能與內部命令中用到的引號相同,否則會導致最相鄰的引號視為一組,引發解釋錯誤。
OFS: 文件輸出時的列分隔符 (output field separtor)
FS: 文件輸入時的列分隔符 (field separtor)
BEGIN: 設置初始參數,初始化變量
END: 讀完文件后做最終的處理
其它{}:循環讀取文件的每一行
$0表示一行內容;$1, $2, … $NF表示第一列,第二列到最后一列。
NF (number of fields)文件多少列;NR (number of rows) 文件讀了多少行: FNR 當前文件讀了多少行,常用于多文件操作時。
a[$1]=1: 索引操作,類似于python中的字典,在ID map,統計中有很多應用。
常見操作
ct@ehbio:~/sxbd$ cat ehbio.wig
variableStep chrom=chr2
300701 12.5
300702 12.5
300703 12.5
300704 12.5
300705 12.5
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}{$2=$2*10^6/(2.5*10^6); print $0}' wig
tep chrom=chr2
0
300701 4.4
300702 4.8
300703 4
300704 4.8
300705 4.8
ct@ehbio:~/sxbd$ cat count
ID Type
Pou5f1 Pluripotency
Nanog Pluripotency
Sox2 Neuron
Tet1 Epigenetic
Tet3 Epigenetic
Myc Oncogene
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}{if(FNR>1) a[$2]+=1;}END{print "Type\tCount"; for(i in a) print i,a[i];}' count
Type Count
Neuron
1
Epigenetic
2
Oncogene
1
Pluripotency
2
ct@ehbio:~/sxbd$ tail -n +
2 count | cut
-f 2 | sort | uniq -c | sed
-e 's/^ *//' -e 's/ */\t/'
2 Epigenetic
1 Neuron
1 Oncogene
2 Pluripotency
之前也提到過的列操作,從GTF文件中提取啟動子區域
sed
's/"/\t/g' GRCh38.gtf | awk
'BEGIN{OFS=FS="\t"}{if($3=="gene") {ensn=$10; symbol=$16; if($7=="+") {start=$4-1; up=start-1000; if(up<0) up=0; dw=start+500; print $1,up, dw, ensn, symbol, $7;} else if($7=="-") {start=$5-1; up=start+1000; dw=start-500; if(dw<0) dw=0; print $1,dw,up,ensn,symbol,$7}}}' | sort -k1,
1 -k2,
2 n >GRCh38.promoter.bed
ct@ehbio:~/sxbd$ cat numertic.matrix
ID A B C
a
1.002 1.234 1.999
b
2.333 4.232 0.889
ct@ehbio:~/sxbd$ awk
'{if(FNR==1) print $0; else {printf "%s%s",$1,FS; for (i=2; i<=NF; i++) printf "%.1f %s", $i, (i==NF?RS:FS)}}' numertic.matrix
ID A B C
a
1.0 1.2 2.0
b
2.3 4.2 0.9
判斷FASTQ文件中,輸出質量值的長度是與序列長度不一致的序列ID
zcat Test_2.fq.gz | awk
'{if(FNR%4==1) ID=$0; else if(FNR%4==2) seq_len=length($0); else if(FNR%4==0) {quality_len=length($0); if(seq_len!=quality_len) print ID; }}'
ct@ehbio:~/sxbd$ cat de_gene
ID log2fc padj
A
1 0.001
B -
1 0.001
C
1 0.001
D
2 0.0001
E -
0.51 0.051
F
0.1 0.1
G
1 0.1 ct@ehbio:~/sxbd$ awk
'$3<0.05 || NR==1' de_gene
ID log2fc padj
A
1 0.001
B -
1 0.001
C
1 0.001
D
2 0.0001 ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}{if(FNR==1) print $0; else {abs_log2fc=($2<0?$2*(-1):$2);if(abs_log2fc>=1 && $3<0.05) print $0;}}' de_gene
ID log2fc padj
A
1 0.001
B -
1 0.001
C
1 0.001
D
2 0.0001
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"; up="up"; dw="dw";}{if(FNR==1) {print $0 >up; print $0 >dw;} else if ($3<0.05) {if ($2>=1) print $0 >up; else if($2<=-1) print $0 >dw;}}' de_gene
ct@ehbio:~/sxbd$ head up dw
==> up <==
ID log2fc padj
A
1 0.001
C
1 0.001
D
2 0.0001 ==> dw <==
ID log2fc padj
B -
1 0.001
ID map,常用于轉換序列的ID、提取信息、合并信息等
ct@ehbio:~/sxbd$ cat id_map
ENSM Symbol Entrez
ENSG00000280516 TMEM42
693149
ENSG00000281886 TGM4
7047
ENSG00000280873 DGKD
8527
ENSG00000281244 ADAMTS13
11093
ENSG00000280701 RP11-
272 D20.
2
ENSG00000280674 ZDHHC3
51304
ENSG00000281623 Y_RNA
ENSG00000280479 CACFD1
11094
ENSG00000281165 SLC2A6
11182
ENSG00000281879 ABO
28
ENSG00000282873 BCL7A
605
ENSG00000280651 AC156455.
1 100506691
ct@ehbio:~/sxbd$ vim ensm
ct@ehbio:~/sxbd$ cat ensm
ENSG00000281244
ENSG00000281165
ENSG00000282873
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}ARGIND==1{if(FNR>1) ensm2entrez[$1]=$3;}ARGIND==2{print ensm2entrez[$1];}' id_map ensm
11093
11182
605 ct@ehbio:~/sxbd$ grep -w
-f ensm id_map | cut
-f 3
11093
11182
605
ct@ehbio:~/sxbd$ cat symbol
Tgm4
Dgkd
Aboct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}ARGIND==1{if(FNR>1) ensm2entrez[$2]=$3;}ARGIND==2{print ensm2entrez[toupper($1)];}' id_map symbol
7047
8527
28
awk
'BEGIN{OFS="\t";FS="\t"}{print log($0)/log(2)}' file
awk
'BEGIN{OFS="\t";FS="\t"}{print int($1+0.5);}' file
awk
'function abs(x){return ((x < 0.0) ? -x : x)}BEGIN{OFS="\t";FS="\t"}{pos[1]=$1;pos[2]=$2;pos[3]=$3;pos[4]=$4; len=asort(pos);for(i=len;i>1;i--) print abs(pos[i]-pos[i-1]);}' file
ct@ehbio:~/sxbd$ cat ens.bed
1 100 105
2 100 105
3 100 105
Mt
100 105
X
100 105
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}{if($1~/^[0-9XY]/) $1="chr"$1; else if($1~/M.*/) gsub(/M.*/, "chrM", $1); print $0}' ens.bed
chr1
100 105
chr2
100 105
chr3
100 105
chrM
100 105
chrX
100 105
ct@ehbio:~/sxbd$ cat trinity_id
Trinity_C1_g1_i1
Trinity_C1_g1_i2
Trinity_C1_g1_i3
Trinity_C2_g1_i1
Trinity_C3_g1_i1
Trinity_C3_g3_i2
ct@ehbio:~/sxbd$ awk
'BEGIN{OFS=FS="\t"}{count=split($1, geneL, "_"); gene=geneL[1]; for(i=2;i<count;i++) gene=gene"_"geneL[i]; print gene,$1;}' trinity_id
Trinity_C1_g1 Trinity_C1_g1_i1
Trinity_C1_g1 Trinity_C1_g1_i2
Trinity_C1_g1 Trinity_C1_g1_i3
Trinity_C2_g1 Trinity_C2_g1_i1
Trinity_C3_g1 Trinity_C3_g1_i1
Trinity_C3_g3 Trinity_C3_g3_i2
cat <<END >grade.awk
f ( avg >=
90 ) grade=
"A" ;
else if ( avg >=
80 ) grade =
"B" ;
else if (avg >=
70 ) grade =
"C" ;
else grade=
"D" ;
ENDawk
-f grade.awk grade
awk
'BEGIN{OFS="\t";FS="\t"}NR!=1{$4=$4"_"NR;print $0}' file
糅合操作
awk
'BEGIN{OFS=FS="\t"}{system("mv "$1".fq "$2".fq");}' input_mat
ct@ehbio:~/sxbd$
echo 1 | awk -v ehbio=
"shengxinbaodian" -v ehbio2=
"SXBD" '{print ehbio, ehbio2;}'
shengxinbaodian SXBD
學會了基本命令,生信分析中還有一大塊是使用已經安裝好的工具,針對軟件安裝中遇到的問題,推出了系列文章:環境變量和可執行屬性徹底釋義環境變量的概念;列舉出Linux下軟件安裝的各種方法,并針對Docker和Conda分別發文介紹。
生信寶典
Linux學習 - 常用和不太常用的實用awk命令
總結
以上是生活随笔 為你收集整理的Linux学习 - 常用和不太常用的实用awk命令 的全部內容,希望文章能夠幫你解決所遇到的問題。
如果覺得生活随笔 網站內容還不錯,歡迎將生活随笔 推薦給好友。