2021.3.17
持續(xù)更新中。婆芦。倦西。

?總綱

• 材料：Illumina雙末端測序原始數(shù)據(jù)（通常是兩個(gè)fastaq的.gz壓縮文件倍宾，文中為：sequence_1.fastq.gz和sequence_2.fastq.gz）
• 軟件：FastQC、Trimmomatic醇坝、SPAdes邑跪、QUAST（均可通過conda一鍵下載）

1. FastQC —— 質(zhì)量評估

FastQC是一款基于Java的軟件，它可以快速地對測序數(shù)據(jù)進(jìn)行質(zhì)量評估呼猪，并生成網(wǎng)頁版的報(bào)告画畅。

1.1 使用

fastqc <sequence_1.fastq.gz> <sequence_2.fastq.gz> -o <目錄名> -t 線程數(shù)

重要參數(shù)
1. -o：輸出文件目錄（需要提前新建一個(gè)目錄）
2. -t：選擇程序運(yùn)行的線程數(shù)
3. -f：指定輸入文件格式

1.2 輸出文件

每個(gè)原始數(shù)據(jù)文件會生成兩個(gè)文件，一個(gè)為.html網(wǎng)頁文件宋距，一個(gè)為.zip文件轴踱。打開.html可以查看序列的測序質(zhì)量情況，分為三個(gè)等級：合格項(xiàng)為綠色√谚赎，警告是黃色的淫僻！，不合格為紅色的×壶唤。（綠色越多越好）

fastqc質(zhì)量報(bào)告

2. Trimmomatic —— 過濾

Trimmomatic用于illumina二代測序數(shù)據(jù)的reads處理雳灵，主要用于對接頭(adapter)序列和低質(zhì)量序列進(jìn)行過濾。能夠識別fastq的.gz和.bz2文件闸盔。主要有兩種模式：雙端模式（PE）和單端模式（SE）悯辙。

注：二代測序下機(jī)數(shù)據(jù)一般是150bp的序列，這些序列理論上全是自己的目標(biāo)序列蕾殴，但是有可能由于測通而含有接頭序列笑撞。

2.1 使用(一條命令)

trimmomatic PE -phred33  <sequence_1.fastq.gz> <sequence_2.fastq.gz> 
 <目錄/sequence_1_paired.fastq.gz> <目錄/sequence_1_unpaired.fastq.gz> <目錄/sequence_2_paired.fastq.gz> <目錄/sequence_2_unpaired.fastq.gz> 
 ILLUMINACLIP:~/miniconda3/pkgs/trimmomatic-0.39-1/share/trimmomatic/adapters/TruSeq3-PE.fa:2:30:10:1:TRUE
 SLIDINGWINDOW:5:20 LEADING:20 TRAILING:20 MINLEN:75

重要參數(shù)
1. PE：指定為雙末端測序（單端用SE）
2. -phred33將 reads 的堿基質(zhì)量值體系轉(zhuǎn)為 phred-33，也可以設(shè)置為phred-64钓觉。
3. ILLUMINACLIP：過濾 reads 中的 Illumina 測序接頭和引物序列茴肥，并決定是否去除反向互補(bǔ)的 1/2 中的 2。ILLUMINACLIP參數(shù)后面的是.fa文件荡灾，用conda安裝的話一般在miniconda3/pkgs/trimmomatic-0.39-1/share/trimmomatic/adapters/TruSeq3-PE.fa下瓤狐，其中有兩種文件可選：TruSeq3-PE.fa和TruSeq2-PE.fa瞬铸，TruSeq3-PE.fa適用于Hiseq和Miseq機(jī)器測序，TruSeq2-PE.fa用于GAII機(jī)器測序础锐。（TruSeq和Nextera是DNA建庫常用的試劑盒）
4. SLIDINGWINDOW： 從 reads 的 5' 端開始嗓节，進(jìn)行滑窗質(zhì)量過濾。示例的意思是以5bp為窗口皆警，若這5bp堿基的平均質(zhì)量值低于20拦宣，則要進(jìn)行切除。
5. LEADING：從reads起始開始信姓，去除質(zhì)量低于閾值或?yàn)?N'的堿基鸵隧。
6. TRAILING：從 reads 的末尾開始切除質(zhì)量值低于閾值的堿基。
7. MINLEN：如果經(jīng)過剪切后 reads 的長度低于閾值則丟棄這條 reads意推。

可選參數(shù)
1. AVGQUAL：如果 reads 的平均堿基質(zhì)量值低于閾值則丟棄這條 reads豆瘫。
2. MAXINFO：一個(gè)自動調(diào)整的過濾選項(xiàng)，在保證 reads 長度的情況下盡量降低測序錯(cuò)誤率菊值，最大化 reads 的使用價(jià)值外驱。

2.2 輸出文件

一共會輸出四個(gè)文件，其中雙端序列都保留的序列在sequence_1_paired.fastq.gz和sequence_2_paired.fastq.gz中可以用于下一步的序列拼接腻窒。只保留一條的序列在sequence_1_unpaired.fastq.gz和sequence_2_unpaired.fastq.gz中昵宇。

過濾之后，還需要用fastqc對過濾后的數(shù)據(jù)進(jìn)行質(zhì)控定页，如果不符合要求趟薄，可重新更改參數(shù)進(jìn)行過濾绽诚。

3. SPAdes —— 短序列拼接

SPAdes是常用的序列拼接軟件之一典徊，支持illumina恩够、PacBio卒落、Nanopore、Sanger蜂桶、Ion Torrent等測序數(shù)據(jù)的拼接儡毕，同樣適合用于混合組裝來改善拼接效果。

3.1 使用

spades.py --pe1-1 <sequence_1_paired.fastq.gz> --pe1-2 <sequence_2_paired.fastq.gz> 
 -t 10  -m 100 --careful --phred-offset 33 -o <spades_out>

重要參數(shù)：
1. --pe<#>-1：指定輸入文庫扑媚，其中#表示第幾個(gè)文庫腰湾。例如，第一pairend文庫就可以寫成--pe1-1疆股，之后接reads1文件
2. --pe<#>-2：和上面類似费坊，如果是大片段的matepair文庫，就使用--mp1-1等旬痹。
3. -t：線程數(shù)附井，默認(rèn)是16個(gè)讨越。
4. -m：用于內(nèi)存限制。
5. --careful：減少錯(cuò)誤和插入序列永毅，添加此項(xiàng)會消耗更多的時(shí)間把跨。
6. --phred-offset：堿基質(zhì)量體系，在數(shù)據(jù)糾錯(cuò)中會用到沼死，現(xiàn)在illumina數(shù)據(jù)一般采用phred 33着逐。
7. -o：輸出目錄

可選參數(shù)：
1. -k：k-mer值列表，數(shù)字必須是小于128的奇數(shù)意蛀。注：若小片段文庫(150bp×2)數(shù)據(jù)量足夠(50×+)滨嘱，則推薦使用-k 21,33,55,77。
2. --pacbio：指定pacbio數(shù)據(jù)輸入浸间。
3. --nanopore：指定nanopore數(shù)據(jù)輸入太雨。

3.2 輸出文件

著重關(guān)注文件*scaffolds.fasta和*contigs.fasta文件即可。

4. QUAST —— 評價(jià)組裝結(jié)果

4.1 安裝

conda install quast

在用conda安裝了FastQC魁蒜、Trimmomatic囊扳、SPAdes后，環(huán)境中的python版本是3.9兜看，安裝不上QUAST了锥咸，在一個(gè)低版本的python環(huán)境下安裝了QUAST。

4.2 使用

quast.py -o <quast_out>  -R <reference_genome.fa> -t 25 <*scaffolds.fasta> <*scaffolds.fasta>

重要參數(shù)：
1. -o：輸出目錄
2. -R：參考基因組序列
3. -t：線程數(shù)

5. 總結(jié)

自己提取基因組送去公司測了以此框架圖细移，回饋的原始數(shù)據(jù)深度達(dá)到了200×+搏予，但是自己通過這個(gè)軟件流程下來發(fā)現(xiàn)最后拼接得到的scaffold的數(shù)量有2000+！；≡雪侥！但是公司最后回饋的草圖只有100條不到的scaffold數(shù)量。說明僅靠這些步驟拼接出來的質(zhì)量并不是很高＞铩Ｋ儆А！