基因組Survey(二代測序數(shù)據(jù)質(zhì)控）

2021-06-26

一. 為什么要做基因組Survey?

Survey分析要做什么數(shù)據(jù)準備急膀？
（1）QC方法介紹
（2）NT方法介紹

1.為什么要進行Survey分析帽氓？
2.Survey分析數(shù)據(jù)準備
3.Survey 數(shù)據(jù)質(zhì)控軟件
4.重點總結(jié)

1.為什么要進行Survey分析（目的）？

1.1 Survey方案

通過質(zhì)控分瘾、NT比對朴肺，獲得高質(zhì)量的clean data，為后續(xù)分析奠定良好基礎聪轿；
基因組Survey基于小片段文庫的低深度測序數(shù)據(jù)（50X）左右；
通過K-mer分析猾浦，有效的評估基因組大小陆错、GC含量、雜合度以及重復序列的含量等信息金赦；
全面了解某一物種基因組特征的有效方法危号；
為后續(xù)的全基因denove測序的組裝策略的制定提供理論依據(jù)。

Survey.png

1.2 基因組復雜度預估

普通基因組的定義素邪？
答：單倍體、純合二倍體或者雜合度<0.5%猪半，且重復序列含量<50%兔朦，GC含量為35%到65%之間的二倍體。
復雜基因組的定義磨确？
答：雜合度>0.5%沽甥，重復序列含量>50%，多倍體乏奥，GC含量處于異常的范圍（GC含量<35%或者GC含量>65%的二倍體）摆舟。
二倍體復雜基因組進一步細分為
微雜合基因組（0.5%<雜合率<=0.8%)
高雜合基因組（雜合率>0.8%)
高重復基因組（重復序列比例>50%）
基因組大小：
基因組越大邓了，測序花錢越多

二恨诱、Survey分析數(shù)據(jù)準備？

2.1 Survey分析需要準備的數(shù)據(jù)骗炉？

二代.png

2.2 Survey測序數(shù)據(jù)質(zhì)量值說明：

堿基的質(zhì)量都是以ASCII值表示的照宝，根據(jù)測序時采用的質(zhì)量方案的不同，計算十進制的質(zhì)量值的方法也有所區(qū)別句葵，常見的計算方法如下所示：
展示方式：Phred+33和Phred+64厕鹃，這里的33和64就是指ASCII值轉(zhuǎn)換為得分該減去的數(shù)值
（1）Phred+64:質(zhì)量字符的ASCII值 -64
（2）Phred+33: 質(zhì)量字符的ASCII值 -33
Illumina 測序堿基質(zhì)量值的范圍是[0,40]，即ASCII值表示為[B,h] 戒 [#,I]乍丈。
Illumina 測序錯誤率不測序質(zhì)量值簡明對應關系剂碴。具體地，如果測序錯誤率用 E 表示轻专，Illumina 堿基質(zhì)量值用 Q 表示忆矛，則有如下關系： Q = -10 log10(E)。

三请垛、 Survey數(shù)據(jù)質(zhì)控軟件洪碳？

3.1 質(zhì)控流程：

流程.png

各提出10000對比對到NT庫递览，如果都比對到同源物種，說明無污染瞳腌，如果比對到細菌真菌绞铃，可能數(shù)據(jù)有污染。

3.2 質(zhì)控軟件-trimmomatic

trimmomatic的下載
方法一：官網(wǎng)：http://www.usadellab.org/cms/index.php?page=trimmomatic 下載二進制文件嫂侍，上傳服務器

trim.png

方法二：服務器命令行下載儿捧、安裝與使用：

#wget下載二進制文件
wget http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.39.zip
#解壓縮Trimmomatic-0.39.zip文件
unzip Trimmomatic-0.39.zip
#Trimmomatic 質(zhì)控用法
#根據(jù)單端測序和雙端測序兩種模式，Trimmomatic軟件也有兩種質(zhì)控用法
#1. SE 模式
#SE模式下挑宠，只有一個輸入文件和一個質(zhì)控后的輸出文件菲盾，運行命令如下:
java -jar <path to trimmomatic jar> SE [-threads <threads>] [-phred33 | -phred64] [-trimlog <logFile>] <input> <output> <step 1> ...表示每一步的質(zhì)控參數(shù)
#2. PE模式
java -jar <path to trimmomatic.jar> PE [-threads <threads] [-phred33 | -phred64] [-trimlog <logFile>] <input 1> <input 2> <paired output 1> <unpaired output 1> <paired output 2> <unpaired output 2> <step 1>...
#實操
fq1=XX1_H3J2NDMXX_L1_1.clean.fq.gz
fq2=XX1_H3J2NDMXX_L1_2.clean.fq.gz
java -jar Trimmomatic-0.39/trimmomatic-0.39.jar \
PE \
-threads 16 \
-validatePairs $fq1 $fq2 \
#ILLUMINACLIP:BGI-SEQ-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15  AVGQUAL:20 MINLEN:36 \
ILLUMINACLIP:Trimmomatic-0.39/adapters/TruSeq3-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 AVGQUAL:20 MINLEN:36 \
-baseout ./pop_clean.fastq.gz
#輸出結(jié)果
pop_clean_1P.fastq.gz pop_clean_1U.fastq.gz pop_clean_2P.fastq.gz pop_clean_2U.fastq.gz

3.3 Trimmomatic的相關參數(shù)

Trimmomatic 過濾數(shù)據(jù)的步驟不命令行中過濾參數(shù)的順序有關，通常的過濾步驟如下：
0.threads：設置線程數(shù)目各淀，多線程運行
1.ILLUMINACLIP: 過濾 reads 中的 Illumina 測序接頭和引物序列
2.SLIDINGWINDOW: 從 reads 的 5' 端開始懒鉴，迚行滑窗質(zhì)量過濾，切掉堿基質(zhì)量平均值低于閾值的滑窗碎浇。
3.MAXINFO: 一個自動調(diào)整的過濾選項临谱，在保證 reads 長度的情況下盡量降低測序錯誤率，最大化 reads 的使用價值奴璃。
4.LEADING: 從 reads 的開頭切除質(zhì)量值低于閾值的堿基悉默。
5.TRAILING: 從 reads 的末尾開始切除質(zhì)量值低于閾值的堿基。
6.CROP: 從 reads 的末尾切掉部分堿基使得 reads 達到指定長度苟穆。
7.HEADCROP: 從 reads 的開頭切掉指定數(shù)量的堿基抄课。
8.MINLEN: 如果經(jīng)過剪切后 reads 的長度低于閾值則丟棄這條 reads。
9.AVGQUAL: 如果 reads 的平均堿基質(zhì)量值低于閾值則丟棄這條 reads雳旅。
10.TOPHRED33: 將 reads 的堿基質(zhì)量值體系轉(zhuǎn)為 phred-33跟磨。
11.TOPHRED64: 將 reads 的堿基質(zhì)量值體系轉(zhuǎn)為 phred-64。
最新的版本并不用選擇10或者11攒盈，它自己會識別格式吱晒，自動選擇

3.4 質(zhì)控報告生成軟件fastqc下載、安裝與使用

# fastqc下載 版本FastQC v0.11.9
wget https://www.bioinformatics.babraham.ac.uk/projects/fastqc/fastqc_v0.11.9.zip
unzip fastqc_v0.11.9.zip #解壓
cd FastQC
chmod 755 ./fastqc
#調(diào)用方式
./fastqc -o outdir -t threads fastq1 fastq2

QC報告展示

1.png

2.png

3.png

4.png

5.png

6.png

7.png

8.png

3.5 NT數(shù)據(jù)庫比對

1.NT庫
Partially non-redundant nucleotide from all traditional divisions of GenBank, EMBL, and DDBJ excluding GSS,STS, PAT, EST, HTG, and WGS.

image.png

非常全面的核酸數(shù)據(jù)庫
https://ftp.ncbi.nih.gov/blast/db/FASTA/

2.NT比對
軟件：BLAST
The Basic Local Alignment Search Tool (BLAST) is the most widely used sequence similarity tool. There are versions of BLAST that compare protein queries to protein databases, nucleotide queries to nucleotide databases, as well as versions that translate nucleotide queries or databases in all six frames and compare to protein databases or queries.

blastn \
-query reads_2.fa \
-db /local_data1/public_data/database/genome_DB/nt/20190417_all/nt \
-out reads_2.csv \
-outfmt "10 evalue length qseqid qlen qstart qend sacc slen sstart send pident nident sstrand qcovs qseq sseq sgi stitle" \
-num_threads 4 -evalue 1e-5 -max_target_seqs 1

3.NT比對結(jié)果文件統(tǒng)計