基因组Survey(二代测序数据质控）

2021-06-26

一. 为什么要做基因组Survey?

Survey分析要做什么数据准备？
（1）QC方法介绍
（2）NT方法介绍

1.为什么要进行Survey分析？
2.Survey分析数据准备
3.Survey 数据质控软件
4.重点总结

1.为什么要进行Survey分析（目的）？

1.1 Survey方案

通过质控、NT比对，获得高质量的clean data，为后续分析奠定良好基础；
基因组Survey基于小片段文库的低深度测序数据（50X）左右；
通过K-mer分析，有效的评估基因组大小、GC含量、杂合度以及重复序列的含量等信息；
全面了解某一物种基因组特征的有效方法；
为后续的全基因denove测序的组装策略的制定提供理论依据。

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分析数据准备？

2.1 Survey分析需要准备的数据？

二代.png

2.2 Survey测序数据质量值说明：

碱基的质量都是以ASCII值表示的，根据测序时采用的质量方案的不同，计算十进制的质量值的方法也有所区别，常见的计算方法如下所示：
展示方式：Phred+33和Phred+64，这里的33和64就是指ASCII值转换为得分该减去的数值
（1）Phred+64:质量字符的ASCII值 -64
（2）Phred+33: 质量字符的ASCII值 -33
Illumina 测序碱基质量值的范围是[0,40]，即ASCII值表示为[B,h] 戒 [#,I]。
Illumina 测序错误率不测序质量值简明对应关系。具体地，如果测序错误率用 E 表示，Illumina 碱基质量值用 Q 表示，则有如下关系： Q = -10 log10(E)。

三、 Survey数据质控软件？

3.1 质控流程：

流程.png

各提出10000对比对到NT库，如果都比对到同源物种，说明无污染，如果比对到细菌真菌，可能数据有污染。

3.2 质控软件-trimmomatic

trimmomatic的下载
方法一：官网：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 质控用法
#根据单端测序和双端测序两种模式，Trimmomatic软件也有两种质控用法
#1. SE 模式
#SE模式下，只有一个输入文件和一个质控后的输出文件，运行命令如下:
java -jar <path to trimmomatic jar> SE [-threads <threads>] [-phred33 | -phred64] [-trimlog <logFile>] <input> <output> <step 1> ...表示每一步的质控参数
#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
#输出结果
pop_clean_1P.fastq.gz pop_clean_1U.fastq.gz pop_clean_2P.fastq.gz pop_clean_2U.fastq.gz

3.3 Trimmomatic的相关参数

Trimmomatic 过滤数据的步骤不命令行中过滤参数的顺序有关，通常的过滤步骤如下：
0.threads：设置线程数目，多线程运行
1.ILLUMINACLIP: 过滤 reads 中的 Illumina 测序接头和引物序列
2.SLIDINGWINDOW: 从 reads 的 5' 端开始，迚行滑窗质量过滤，切掉碱基质量平均值低于阈值的滑窗。
3.MAXINFO: 一个自动调整的过滤选项，在保证 reads 长度的情况下尽量降低测序错误率，最大化 reads 的使用价值。
4.LEADING: 从 reads 的开头切除质量值低于阈值的碱基。
5.TRAILING: 从 reads 的末尾开始切除质量值低于阈值的碱基。
6.CROP: 从 reads 的末尾切掉部分碱基使得 reads 达到指定长度。
7.HEADCROP: 从 reads 的开头切掉指定数量的碱基。
8.MINLEN: 如果经过剪切后 reads 的长度低于阈值则丢弃这条 reads。
9.AVGQUAL: 如果 reads 的平均碱基质量值低于阈值则丢弃这条 reads。
10.TOPHRED33: 将 reads 的碱基质量值体系转为 phred-33。
11.TOPHRED64: 将 reads 的碱基质量值体系转为 phred-64。
最新的版本并不用选择10或者11，它自己会识别格式，自动选择

3.4 质控报告生成软件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
#调用方式
./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数据库比对

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

非常全面的核酸数据库
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比对结果文件统计