seqtkはfastqをfastaに変換したり、相補鎖に変換できるツール。ランダムサンプリング機能ももち、de novo transcriptome解析でアセンブルに有利なリードデプスに間引くツールとして用いられることもある(ペーパー)。動作が非常に高速のため使いやすい。似たツールにseqkitがある。導入して機能をチェックしてみる。
インストール
mac os10.12でテストした。
brewで導入できる。
brew install seqtk
#Anaconda環境なら
conda install -c bioconda seqtk
ラン
seq - common transformation of FASTA/Q
fastqをfastaに変換。(.gzも扱える。)
seqtk seq -a input.fq > output.fa
ILLUMINA 1.3+フォーマットをfastaに変換し、同時にquality20以下を小文字にする。
seqtk seq -aQ64 -q20 input.fq > output.fa
相補鎖(reverse complement)に変換。
seqtk seq -r input.fq > output.fq
様々なオプションがある。
- -q mask bases with quality lower than INT [0]
- -X mask bases with quality higher than INT [255]
- -Q quality shift: ASCII-INT gives base quality [33]
- -M mask regions in BED or name list FILE [null]
- -L drop sequences with length shorter than INT [0]
- -c mask complement region (effective with -M)
- -r reverse complement -A force FASTA output (discard quality)
- -C drop comments at the header lines
- -N drop sequences containing ambiguous bases
- -1 output the 2n-1 reads only
- -2 output the 2n reads only
- -V shift quality by '(-Q) - 33'
- -U convert all bases to uppercases
- -S strip of white spaces in sequences
mergepe - interleave two PE FASTA/Q files
ペアリードを1つのファイルにまとめインターレースにする。
seqtk mergepe R1.fq R2.fq > paired.fq
インターレースファイルからペアじゃないリードを除去するならdropseを使う。
trimfq - trim FASTQ using the Phred algorithm
低クオリティな塩基の除去。Phread quality scoreに基づく。
eqtk trimfq input.fa > output.fa
- -q error rate threshold (disabled by -b/-e) [0.05]
- -l maximally trim down to INT bp (disabled by -b/-e) [30]
デフォルトだと5%のエラー率を超えるとトリムされる(-q 0.05)。
指定範囲の塩基のトリミング。
eqtk trimfq -b 5 -e 10 input.fa > output.fa
- -b trim INT bp from left (non-zero to disable -q/-l) [0]
- -e trim INT bp from right (non-zero to disable -q/-l) [0]
- -L retain at most INT bp from the 5'-end (non-zero to disable -q/-l) [0]
subseq - extract subsequences from FASTA/Q
fastqをランダムサンプリングする。-sで乱数発生を制御している。
seqtk sample -s100 read1.fq 10000 > sub1.fq
seqtk sample -s100 read2.fq 10000 > sub2.fq
- -s RNG seed [11]
- -2 2-pass mode: twice as slow but with much reduced memory
上記では10000リードずつ取り出される。ただし-sの値を変えるとペアじゃないリードが10000リード取り出されるので注意。ペアリードは同じパラメータで行う。
fqchk - fastq QC (base/quality summary)
クオリティの分析。
seqtk fqchk input.fq > output.txt
$ head -5 output.txt
min_len: 101; max_len: 101; avg_len: 101.00; 38 distinct quality values
POS #bases %A %C %G %T %N avgQ errQ %low %high
ALL 252500 26.2 24.3 23.7 25.8 0.0 35.9 20.4 2.8 97.2
1 2500 13.0 44.2 25.4 17.4 0.0 32.6 30.7 0.4 99.6
2 2500 21.0 22.1 21.7 35.2 0.0 32.9 31.7 0.4 99.6
3 2500 25.4 25.2 22.1 27.2 0.0 33.0 31.8 0.2 99.8
Rに持ち込んでグラフ化すれば見やすい。
他にもいくつかの機能がある。
引用
GitHub - lh3/seqtk: Toolkit for processing sequences in FASTA/Q formats
Using seqtk to trim and process reads at an insanely high speed — ANGUS 2.0 documentation
FASTA/Q sequence processing toolkit -- seqtk