2021.3.18
持續(xù)更新中辈毯。陨舱。。
目的：比較下目前常用的四款拼接軟件：Velvet共耍、SPAdes虑灰、SOAPdenovo、ABySS
1. FastQC —— 質(zhì)量評估
????FastQC是一款基于Java的軟件痹兜，它可以快速地對測序數(shù)據(jù)進行質(zhì)量評估穆咐，并生成網(wǎng)頁版的報告。
1.1 使用
fastqc <sequence_1.fastq.gz> <sequence_2.fastq.gz> -o <目錄名> -t 線程數(shù)
重要參數(shù)
1. -o：輸出文件目錄
2. -t：選擇程序運行的線程數(shù)
1.2 輸出文件
????每個原始數(shù)據(jù)文件會生成兩個文件字旭，一個為html網(wǎng)頁文件对湃，一個為.zip*文。打開.html可以查看序列的測序質(zhì)量情況遗淳，分為三個等級：合格項為綠色√拍柒，警告是黃色的！屈暗，不合格為紅色的×拆讯。（綠色越多越好）
2. Cutadapt —— 過濾
????Cutadapt是一款用python寫的去接頭軟件，可以去除任意指定的接頭和序列养叛，同時也可以用于質(zhì)量過濾种呐。（實際上原始數(shù)據(jù)過濾會過濾兩種序列：去接頭引物、低質(zhì)量序列）
公司返回的
rawdata不含有引物序列弃甥，可能含有接頭序列爽室，需要進行過濾后進行分析
2.1 下載
conda install cutadapt -y
注意cutadapt(v3.3)的安裝環(huán)境python版本不能高于3.9
2. 2 使用
cutadapt -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCA -A AGATCGGAAGAGCACACGTCTGAACTCCAGTCA -o sequence_filtered_1.fastq.gz -p sequence_filtered_2.fastq.gz
-q 20,20 --max-n 0.1 -m 75 sequence_1.fastq.gz sequence_2.fastq.gz
重要參數(shù)
1. -a：paired-end測序文件中正向序列文件接頭序列（文件1）
2. -A：paired-end測序文件中反向序列文件接頭序列（文件2）
3. -o：文件1去接頭后的結(jié)果文件
4. -p：文件2去接頭后的結(jié)果文件
5. -q：序列兩端堿基質(zhì)量低于某一數(shù)值時被切除，用,隔開淆攻，例如：-q 20,20阔墩。
6. -m：reads1和reads2中切除接頭后的序列長度最小值掉缺，低于這個數(shù)值則去除，例如：-m 75戈擒。
7. --max-n：N堿基占比一定比例時被去除眶明，例如--max-n 0.1，表示N堿基占read比例到10%時會被去除筐高。
8：最后是兩個原始序列文件搜囱。
可選參數(shù)
1. -- pair-filter=(any|both)：any表示read1 和 read2任何一個檢測到接頭均舍棄；both表示 read1 和 read2 全部檢測到接頭才舍棄read1 和 read2柑土。（默認(rèn)any）
2. -O ：默認(rèn)為3蜀肘，即至少三個堿基配才認(rèn)為是adapter序列。
3. -e：最大錯配比例稽屏，默認(rèn)是0.1扮宠。（解釋：cutadapt在一條read中檢測到20bp的接頭序列，那么允許該20bp的接頭序列有2個堿基的錯配）
3. Kmergenie預(yù)測最佳kmer值
????給定一組輸入狐榔，KmerGenie首先計算多個K值的k-mer豐度直方圖坛增。然后，對每個K值薄腻，它預(yù)測數(shù)據(jù)集中不同基因組k-mer的數(shù)量收捣，并返回這個數(shù)字最大化的k-mer長度。
3.1 安裝
conda create -n kmer python=2.7 -y
conda install kmergenie
3.2 使用
kmergenie <fq.list> -o <result> -s 10 -t 10
重要參數(shù)：
1. fq.list：包含需要查詢的文件庵楷，一行一個文件名（文件所在的絕對路徑）
2. -o：結(jié)果輸出的前綴名
3. -l：系統(tǒng)考慮的最小k值（默認(rèn)：15）
4. -k：系統(tǒng)考慮的最大k值（默認(rèn)：121）
5. -s：從最小k值到最大k值罢艾，每次增加的值（默認(rèn)：10）
6. -t：線程數(shù)
3.3 輸出文件
????輸出文件中主要看*report.html網(wǎng)頁文件，會推算出最合適的kmer值與基因組大小尽纽。在這里咐蚯，我的基因組得出的kmer=111。
話說弄贿，kmer值是越大越好還是越小越好呢春锋？
4. 拼接 —— contigs
4.1 ?Velvet?
????Velvet，基于kmer的序列拼接工具挎春，用于基因組的de novo組裝看疙，支持多種格式，局部拼接效果好直奋，不善于連接成scaffold能庆。
4.1.1 安裝
conda install velvet -y
如果是編譯安裝可以選擇先修改以下配置文件參數(shù)后進行安裝：
1. CATEGORIES：根據(jù)原始數(shù)據(jù)相應(yīng)增減該值的小（一般設(shè)置為10）脚线。
2. MAXKMERLENGTH：最大的kmer長度（一般設(shè)置為127）搁胆。
4.1.2 使用
分為兩步執(zhí)行：
步驟一：利用velveth對數(shù)據(jù)構(gòu)建一個hash表
velveth <output> 111 -shortPaired -fastq -separate <sequence_filtered_1.fastq.gz> <sequence_filtered_2.fastq.gz>
重要參數(shù)：
1. output：輸出文件目錄
2. 111：即hash_lenghth，用來設(shè)置k-mer的大小。也可以是31,97,2的形式渠旁，指分別拼接kmer從31到97攀例，依次增加2的序列（默認(rèn)：31）
3. -shortPaired：reads的類型。
4. -fastq：輸入文件的格式（默認(rèn)是fasta）顾腊。
5. -separate：分開兩個文件讀入粤铭。
步驟二：velvetg進行序列拼接
velvetg <output> -exp_cov auto -cov_cutoff auto
重要參數(shù)：
1. <output>：velveth生成的結(jié)果目錄
2. -exp_cov：期望的kmer覆蓋度，設(shè)置成auto用于標(biāo)準(zhǔn)的基因組測序杂靶。該參數(shù)設(shè)置成auto后梆惯，-cov_cutoff也許設(shè)置成auto。
4.2 ?SPAdes?
????SPAdes是常用的序列拼接軟件之一吗垮，支持illumina垛吗、PacBio、Nanopore烁登、Sanger怯屉、Ion Torrent等測序數(shù)據(jù)的拼接，同樣適合用于混合組裝來改善拼接效果（尤其適用于Ion Torrent數(shù)據(jù)的拼接）
4.2.1 安裝
conda install spades -y
4.2.2 使用
spades.py --pe1-1 <sequence_filtered_1.fastq.gz> --pe1-2 <sequence_filtered_2.fastq.gz> -o <spades_out>
-k 111 -t 20 --careful
重要參數(shù)：
1. --pe<#>-1：指定輸入文庫饵沧，其中#表示第幾個文庫锨络。例如，第一個pairend文庫就可以寫成--pe1-1捷泞，之后接reads1文件足删。
2. --pe<#>-2：和上面類似，如果是大片段的matepair文庫锁右，就使用--mp1-1等。
3. -t：線程數(shù)（默認(rèn)：16）讶泰。
4. --careful：減少錯誤和插入確實咏瑟，添加此項會消耗更多的時間。
5. -o：輸出目錄
可選參數(shù)：
1. -k：k-mer值列表痪署，數(shù)字必須是小于128的奇數(shù)码泞。注：若小片段文庫(150bp×2)數(shù)據(jù)量足夠(50×+)，則推薦使用-k 21,33,55,77（默認(rèn)：auto）
2. --pacbio：指定pacbio數(shù)據(jù)輸入狼犯。
3. --nanopore：指定nanopore數(shù)據(jù)輸入余寥。
4. -m：用于內(nèi)存限制（默認(rèn)：250）
5. --phred-offset：堿基質(zhì)量體系，在數(shù)據(jù)糾錯中會用到悯森，現(xiàn)在illumina數(shù)據(jù)一般采用phred 33質(zhì)量體系（默認(rèn)：auto-detect）宋舷。
4.3 ?SOAPdenovo?
????SOAPdenovo是華大基因開發(fā)的SOAP軟件包的一部分，主要用于短序列reads拼接瓢姻，尤其是illumina測序數(shù)據(jù)祝蝠。從小的細菌基因組到大的動植物基因組，人基因組都適用。特點是使用起來比較簡單绎狭，但是卻擁有很好的拼接效果细溅，尤其在基因組構(gòu)建Scaffold方面，效果很好儡嘶。
4.3.1 安裝
conda install soapdenovo2 -y
4.3.2 使用
SOAPdenovo-127mer all -s config_file -K 111 -o soap -d 1 -p 20
soapdenovo里面有兩個執(zhí)行程序：SOAPdenovo-63mer和SOAPdenovo-127mer喇聊，根據(jù)kmer值的大小進行選擇。
配置文件config_file的設(shè)置
#soapdenovo configure file
max_rd_len=150 #使用reads的最大長度蹦狂，一般設(shè)置為reads的長度誓篱。
[LIB] #文庫信息以此開頭，如有多個文庫鸥咖，每個文庫都需要有一個單獨的[LIB]標(biāo)識符
avg_ins=439 #文庫平均插入長度燕鸽，即測序文庫大小，一般取插入文庫大小的均值啼辣。
reverse_seq=0 #序列是否需要反轉(zhuǎn)啊研。如果是大片段文庫，需要設(shè)置為1鸥拧。
asm_flags=3 #組裝的標(biāo)志党远。1用于構(gòu)建contig，2用于構(gòu)建scaffold富弦，3同時用于構(gòu)建contig和scaffold沟娱，4只用于補洞，不用于拼接腕柜。（一般小片段文庫設(shè)置為3济似，大片段設(shè)置為2）
rank=1 #文庫使用的優(yōu)先級。
pair_num_cutoff=3 #可選參數(shù)盏缤，確定了reads連接成contig或scaffold的閾值砰蠢。
map_len=32 #可選參數(shù)，規(guī)定了在map過程中 reads和contig的比對長度必須達到該值（短片斷默認(rèn)：32）
q1=./reads.1.fq.gz #用于連接的數(shù)據(jù)唉铜，支持fasta和fastq格式文件台舱，也支持壓縮格式。（結(jié)尾不能有空格）
q2=./reads.2.fq.gz
重要參數(shù)：
1. -s：接配置文件
2. -o：輸出文件的文件名前綴
3. -K：輸入kmer的大小
4. -p：線程數(shù)（默認(rèn)：8）
5. -d：去除頻數(shù)不大于該值的kmer（默認(rèn)：0）
6. -D： 去除頻數(shù)不大于該值的由k-mer連接的邊（默認(rèn)：1）
7. -F：利用read對scaffold中的gap進行填補（默認(rèn)：不執(zhí)行）
4.4 ?ABySS?
????ABySS用于從頭拼接雙端測序短序列或較大基因組序列竞惋。
4.4.1 安裝
conda install abyss -y
4.4.2 使用
- 拼接一個paired-end文庫
abyss-pe k=111 name=abyss in='sequence_filtered_1.fastq.gz sequence_filtered_2.fastq.gz'
重要參數(shù)：
1. k：kmer值
2. name：輸出文件名
3. in：兩個paired-end文庫
- 拼接多個paired-end文庫
abyss-pe k=111 name=ecoli lib='pea peb' pea='pea_1.fa pea_2.fa' peb='peb_1.fa peb_2.fa'
5. 四種拼接軟件的比較
5.1 QUAST的安裝和使用
conda install quast -y
quast -o <result> <velvet_contig.fa> <abyss_contig.fa> <soapdenovo_contig.fa> <spades_contig.fa> -t 15
重要參數(shù)：
1. -o：輸出文件目錄
2. -t：選擇程序運行的線程數(shù)
3. -f：指定輸入文件格式
5.2 比較結(jié)果（contigs）
綜合contigs數(shù)量、總長度灰嫉、N50 來看拆宛，velvet拼接contig效果最好
5.2 比較結(jié)果（scaffolds）
????除了Velvet軟件之外，其他三款軟件在拼接contigs的同時也可以拼接成scaffolds胰挑。同樣綜合scaffolds數(shù)量蔓罚、總長度、N50 來看瞻颂，從paired-end reads數(shù)據(jù)一次性拼接成scaffolds效果最好的是ABySS豺谈。但是在使用abyss的過程中沒有找到調(diào)整線程數(shù)的參數(shù)贡这，用時較長茬末。
6. 總結(jié)
????Velvet拼接出的159個contigs已經(jīng)很接近公司反饋的最終結(jié)果了，后面再用報告中提到的兩款軟件試試盖矫。
