1. 引言
下機(jī)的FASTQ數(shù)據(jù)通常需要進(jìn)行質(zhì)控和預(yù)處理,以保證下游分析輸?shù)臏?zhǔn)確性。fastp軟件僅掃描數(shù)據(jù)文件一次,就可以完成FASTQC + cutadapt + Trimmomatic 的功能;而且它使用C++開發(fā),利用高效算法,并且支持多線程,加快了處理速度。
2. 特點(diǎn)
對(duì)輸入和輸出數(shù)據(jù)進(jìn)行詳盡的質(zhì)量剖析,生成人性化的報(bào)告。(quality curves, base contents, KMER, Q20/Q30, GC Ratio, duplication, adapter contents...)
過濾掉低質(zhì)量,太短和太多N的"bad reads"。
從頭部(5')或尾部(3')計(jì)算滑動(dòng)窗內(nèi)的堿基質(zhì)量均值,并切除低質(zhì)量堿基。(類似Trimmomatic,但是非??欤?。
頭尾剪裁。(上下機(jī)序列可能不穩(wěn)定)
去除接頭。(可以不用輸入接頭序列,算法可以自動(dòng)識(shí)別接頭序列并進(jìn)行剪裁)
不匹配堿基對(duì)矯正。(對(duì)于雙端測(cè)序(PE)的數(shù)據(jù),軟件能夠矯正重疊區(qū)域的不匹配堿基對(duì)。如:重疊區(qū)域內(nèi),一個(gè)堿基質(zhì)量很高,另一個(gè)非常低,fastp會(huì)依據(jù)質(zhì)量值進(jìn)行校正。
修剪尾部(3')的polyG。(修剪ployX尾,可以去掉不想要的ploy。如:mRNA-Seq 中的ployA。)
處理使用了唯一分子標(biāo)識(shí)符(UMI)的數(shù)據(jù),并將UMI轉(zhuǎn)換為序列名稱。
不僅生成質(zhì)控和過濾結(jié)果的HTML報(bào)告(存儲(chǔ)在fastp.html,可用-h可指定),還生成了程序可讀性非常強(qiáng)的JSON報(bào)告(fastp.json,可用-j指定)。
為了適合并行處理,fastp將輸出進(jìn)行分拆。有兩種模式可選,a. 指定拆分的個(gè)數(shù); b. 指定拆分后每個(gè)文件的行數(shù)。
fastp支持gzip的輸入和輸出,同時(shí)支持SE和PE數(shù)據(jù)處理;支持PacBio/Nanopore的長(zhǎng)reads數(shù)據(jù);也支持交錯(cuò)輸入。
3. 安裝fastp
Bioconda 安裝
# note: the fastp version in bioconda may be not the latest
conda install -c bioconda fastp
Linux 系統(tǒng)安裝
# this binary was compiled on CentOS, and tested on CentOS/Ubuntu
wget http://opengene.org/fastp/fastp
chmod a+x ./fastp
源下載安裝
# get source (you can also use browser to download from master or releases)
git clone https://github.com/OpenGene/fastp.git
# build
cd fastp
make
# Install
sudo make install
注:如編譯時(shí)遇到(undefined reference to gzbuffer),需要更新zlib。
4. 輸入和輸出
-
fastp支持雙端(PE)數(shù)據(jù)和單端(SE)數(shù)據(jù)的輸入和輸出
對(duì)于SE數(shù)據(jù),使用指令-i 或 --in1輸入,-o 或--out1輸出。
對(duì)于PE數(shù)據(jù),使用指令-I 或 --in2輸入,-O或--out2輸出。
如果不指定輸出文件名,而且輸入名以.gz結(jié)尾,QC會(huì)執(zhí)行完并壓縮。
-
STDIN輸入
--stdin 指明使用STDIN輸入。
--interleaved_in 指定interleaved paired-end stream.
-
STDOUT輸出 - fastp將結(jié)果傳到STDOUT,以便傳遞到bwa和bowtie2。
- --stdout
對(duì)于PE數(shù)據(jù),輸出將與FASTQ交錯(cuò),因此會(huì)有 record1-R1 -> record1-R2 -> record2-R1 -> record2-R2 -> record3-R1 -> record3-R2 ...記錄。
store the unpaired reads for PE data
store the reads that fail the filters
process only part of the data
do not overwrite exiting files
split the output to multiple files for parallel processing
merge PE reads
5. 功能介紹
-
質(zhì)量過濾 - 質(zhì)量過濾功能默認(rèn)開啟
-Q 或 disable_quality_filtering關(guān)閉。
-n 或 --n_base_limit限制堿基N的數(shù)量
-q 或 --qualified_quality_phred指定合格的質(zhì)量值。(默認(rèn)為 -q 15,即質(zhì)量值大于等于Q15的即為合格)
-u 或 --unqualified_percent_limit指定不合格堿基的最高百分比。(默認(rèn)40%)
-e 或--average_qual通過質(zhì)量百分比過濾,如果質(zhì)量值低于平均值,就被丟棄。
-
長(zhǎng)度過濾 - 長(zhǎng)度過濾默認(rèn)開啟
-L 或 --disable_length_filtering關(guān)閉。
-l 或 --length_required指定所需最小長(zhǎng)度。(-l 30表示低于30個(gè)堿基的read會(huì)被丟棄)
--length_limit指定最長(zhǎng)read長(zhǎng)度。(mall RNA sequencing可用)
結(jié)合起來就可以實(shí)現(xiàn)常用的discard模式,以保證所有輸出的序列都一樣長(zhǎng)。
-
低復(fù)雜度過濾 - 默認(rèn)關(guān)閉
- -y 或 --low_complexity_filter開啟。
復(fù)雜度由的定義為reads中與下一個(gè)堿基不同的堿基的百分比。
# a 51-bp sequence, with 3 bases that is different from its next base
seq ='AAAATTTTTTTTTTTTTTTTTTTTTGGGGGGGGGGGGGGGGGGGGGGCCCC'
complexity = 3/(51-1) = 6%
- -Y 或 --complexity_threshold 指定閾值,默認(rèn)為30%。
在fastp的報(bào)告中,第一個(gè)Summary表格很清楚地顯示了過濾的統(tǒng)計(jì)信息。
- 接頭處理 - 默認(rèn)開啟,自動(dòng)識(shí)別
-
-A 或 --disable_adapter_trimming關(guān)閉。
對(duì)于SE數(shù)據(jù),可以用a 或 --adapter_sequence參數(shù)來輸入接頭。(fastp通過分析前1M的reads估計(jì)接頭序列,對(duì)于SE數(shù)據(jù),估計(jì)可能不準(zhǔn)確,而且如果接頭序列是特制的,也可能導(dǎo)致估計(jì)出錯(cuò)。)
對(duì)于PE數(shù)據(jù),自動(dòng)檢測(cè)接頭功能是關(guān)閉的,--detect_adapter_for_pe開啟。亦可以用--adapter_sequence和 --adapter_sequence_r2指定R1, R2接頭序列。
--adapter_fasta 還可以用file 的形式指定需要修剪的任何序列,序列長(zhǎng)度大>6bp。
>Illumina TruSeq Adapter Read 1
AGATCGGAAGAGCACACGTCTGAACTCCAGTCA
>Illumina TruSeq Adapter Read 2
AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT
>polyA
AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
修剪順序:--adapter_sequence | --adapter_sequence_r2 | --adapter_fasta
修剪的接頭序列可以在fastp報(bào)告中找到.
- 根據(jù)質(zhì)量修剪每條read - 有三種選擇,可以全選
-5或 --cut_front指定從5'開始丟棄低于質(zhì)量低的,同時(shí)多N堿基也會(huì)被去除;cut_front_window_size設(shè)定滑窗大??; cut_front_mean_quality指定閾值。
-3或 --cut_tail, cut_tail_window_size, cut_tail_mean_quality.
-r或 --cut_right, cut_right_window_size, cut_right_mean_quality.(滑窗從頭向尾移動(dòng),如果窗內(nèi)平均質(zhì)量低于閾值,丟棄滑窗內(nèi)和右邊所有序列。)
*all these three operations will interfere deduplication for SE data, and --cut_front or --cut_right may also interfere deduplication for PE data. *
If --cut_right is enabled, then there is no need to enable --cut_tail, since the former is more aggressive. If --cut_right is enabled together with --cut_front, --cut_front will be performed first before --cut_right to avoid dropping whole reads due to the low quality starting bases.
cut_front will interfere deduplication for both PE/SE data, and --cut_tail will interfere deduplication for SE data.
- PE數(shù)據(jù)的堿基矯正 - 默認(rèn)關(guān)閉
通過overlap檢索,矯正后堿基質(zhì)量和對(duì)照同等。
-c or --correction 開啟. 可調(diào)參數(shù)overlap_len_require (default 30), overlap_diff_limit (default 5) and overlap_diff_limit_percent (default 20%).
- 全局裁剪 - 統(tǒng)一剪裁reads的頭部和尾部
用于統(tǒng)一去掉低質(zhì)量cycle。如:Illumina最后一個(gè)cycle通常質(zhì)量是非常低,通過-t 1 or --trim_tail1=1去掉。
For read1 or SE data, the front/tail trimming settings are given with -f, --trim_front1 and -t, --trim_tail1.
For read2 of PE data, the front/tail trimming settings are given with -F, --trim_front2 and -T, --trim_tail2. (如read2無指定,read2裁剪參數(shù) = read1。)
-b或 --max_len1指定read1長(zhǎng)度, -B或 --max_len2指定read2長(zhǎng)度。(如read2無指定,read2裁剪參數(shù) = read1。)
1, UMI preprocessing (--umi)
2, global trimming at front (--trim_front)
3, global trimming at tail (--trim_tail)
4, quality pruning at 5' (--cut_front)
5, quality pruning by sliding window (--cut_right)
6, quality pruning at 3' (--cut_tail)
7, trim polyG (--trim_poly_g, enabled by default for NovaSeq/NextSeq data)
8, trim adapter by overlap analysis (enabled by default for PE data)
9, trim adapter by adapter sequence (--adapter_sequence, --adapter_sequence_r2. For PE data, this step is skipped if last step succeeded)
10, trim polyX (--trim_poly_x)
11, trim to max length (---max_len)
- polyG裁剪 - 默認(rèn)對(duì)Illumina NextSeq/NovaSeq data開啟
對(duì)于兩色發(fā)光法的Illumina設(shè)備(NextSeq /NovaSeq),在沒有光信號(hào)情況下base calling的結(jié)果會(huì)返回G,所以在序列的尾端可能會(huì)出現(xiàn)較多的polyG,需要被去除。
fastp會(huì)通過機(jī)器ID自動(dòng)化地識(shí)別NextSeq / NovaSeq的數(shù)據(jù),然后進(jìn)行polyG識(shí)別和裁剪。
-g or --trim_poly_g 向其他數(shù)據(jù)開啟;-G or --disable_trim_poly_g關(guān)閉。
<poly_g_min_len>設(shè)定最小長(zhǎng)度,默認(rèn)值10。
- polyX尾裁剪 - 默認(rèn)關(guān)閉
-x or --polyX; <poly_x_min_len>, 默認(rèn)值10。
同時(shí)開啟時(shí)的裁剪順序:polyG | polyA。
- 分子標(biāo)簽(UMI)處理 - 消除重復(fù)和糾正錯(cuò)誤
-U 或 -umi 開啟。
--umi_loc指定UMI所在的位置,可以是(index1、 index2、 read1、 read2、 per_index、 per_read )中的一種;分別表示UMI是在index位置上,還是在插入片段中。如果指定了是在插入序列中,還需要使用 --umi_len 參數(shù)來指定UMI所占的堿基長(zhǎng)度。
#初始序列
@NS500713:64:HFKJJBGXY:1:11101:1675:1101 1:N:0:TATAGCCT+GACCCCCA
AAAAAAAAGCTACTTGGAGTACCAATAATAAAGTGAGCCCACCTTCCTGGTACCCAGACATTTCAGGAGGTCGGGAAA
+
6AAAAAEEEEE/E/EA/E/AEA6EE//AEE66/AAE//EEE/E//E/AA/EEE/A/AEE/EEA//EEEEEEEE6EEAA
?
fastp -i R1.fq -o out.R1.fq -U --umi_loc=read1 --umi_len=8
?
@NS500713:64:HFKJJBGXY:1:11101:1675:1101:AAAAAAAA 1:N:0:TATAGCCT+GACCCCCA
GCTACTTGGAGTACCAATAATAAAGTGAGCCCACCTTCCTGGTACCCAGACATTTCAGGAGGTCGGGAAA
+
EEE/E/EA/E/AEA6EE//AEE66/AAE//EEE/E//E/AA/EEE/A/AEE/EEA//EEEEEEEE6EEAA
- 輸出拆分 - 兩種模式
-s or --split指定file數(shù)目。
-S or --split_by_lines指定file中的行數(shù)。
- **過表達(dá)序列(Overrepresented sequence analysis) **- 默認(rèn)關(guān)閉
-p or --overrepresentation_analysis開啟。( count 長(zhǎng)度為10bp, 20bp, 40bp, 100bp or (cycles - 2 ), 默認(rèn)值20兼顧速度與精確性)
在fastp的報(bào)告中,不僅提供overrepresented sequence的序列個(gè)數(shù)和占比,還提供了其在測(cè)序cycles中的分布情況。
- PE數(shù)據(jù)reads拼接
對(duì)于PE數(shù)據(jù),fastp通過-m/--merge選項(xiàng)實(shí)現(xiàn)拼接模式。