二代测序数据质控 QC

NGS-QC

illumina测序原理

• 高通量测序（High-Throughput Sequencing）又名二代测序|下一代测序（Next Generation Sequencing，NGS），是相对于传统的桑格测序|一代测序（Sanger Sequencing）而言的。相对于Sanger测序而言，二代测序可以提供中等的读长和适中的价格，适合de novo 测序、转录组测序、宏基因组研究等。
• Solexa的测序原理是可逆终止化学反应。Solexa是一种基于边合成边测序技术(Sequencing-By-Synthesis，SBS)的新型测序方法。通过利用单分子阵列实现在小型芯片(Flow Cell)上进行桥式PCR反应。由于新的可逆阻断技术可以实现每次只合成一个碱基，并标记荧光基团，再利用相应的激光激发荧光基团，捕获激发光，从而读取碱基信息。
• 二代测序建库测序大致流程
  DNA片段经末端修复、加ployA尾、加测序接头、纯化、PCR扩增等步骤完成整个文库制备。构建好的文库通过illumina HiSeqTM PE150进行测序。文库构建完成后，先使用Qubit2.0进行初步定量，稀释文库至1ng/μl，随后使用Agilent 2100对文库的insert size进行检测，insert size符合预期后，使用Q-PCR方法对文库的有效浓度进行准确定量（文库有效浓度>2nM），以保证文库质量。

二代测序数据拆分

• 原始下机数据睡bcl文件，根据前面建库的index信息，进行数据的拆分，除非是包lane或者包run，否则二代测序公司是不会提供该文件的
• 外包测序返回的是拆分后的rawdata及质控后的cleandata，由rawdata到cleandata的数据过滤过程称为质控

二代测序数据质控

• 质控主要进行低质量，含N，含adpter的过滤
• 过滤主要考虑的参数：
  1. 数据有效数据利用率，一般要求高于95%，现在正常项目大多在99%
  2. 数据量，数据量所有样品，高于约定数据量的95%，看合同签订的是raw还是clean
  3. Q20一般要>90%（illunima官方承诺85%）
  4. Q30一般要>85%（illunima官方承诺80%）
  5. GC含量，一般波动不大，5%波动以内，群体复杂的要特殊考虑
  6. GC波动情况（WGS几乎无波动，简化基因组及panel的另行考虑）
  7. NT比对情况，要求无污染，现在公司不会直接提供，GC波动大时，可以要求测序公司提供，以排除污染。

FQ数据格式

• 高通量测序(如Illumina HiSeqTM/MiseqTM)得到的原始图像数据文件经CASAVA碱基识别(Base Calling)分析转化为原始测序序列（Sequenced Reads），我们称之为 Raw Data或Raw Reads，结果以 FASTQ (简称为fq)文件格式存储，其中包含测序序列（reads）的序列信息以及其对应的测序质量信息。
  FASTQ格式文件中每个read由四行描述，如下：

@HWI-ST1276:71:C1162ACXX:1:1101:1208:2458 1:N:0:CGATGT
NAAGAACACGTTCGGTCACCTCAGCACACTTGTGAATGTCATGGGATCCAT
+
#55???BBBBB?BA@DEEFFCFFHHFFCFFHHHHHHHFAE0ECFFD/AEHH

• 其中：
  第一行以“@”开头，随后为Illumina 测序标识别符(Sequence Identifiers)和描述文字(选择性部分)；
  第二行是碱基序列；
  第三行以“+”开头，随后为Illumina 测序标识别符(选择性部分)；
  第四行是对应碱基的测序质量，该行中每个字符对应的 ASCII 值减去 33，即为对应第二行碱基的测序质量值。

原始数据过滤

• 测序得到的原始测序序列（Sequenced Reads）或者 raw reads，里面含有带接头的、低质量的reads。为了保证信息分析质量，必须对raw reads过滤，得到clean reads，后续分析都基于 clean reads。数据处理的条件如下（非标准条件，可参考，比较松的条件,这个是诺禾的过滤条件，大家比例会有所调整，但是都是过滤的这三项）：
  • 去除带接头(adapter)的reads pair；
  • 当单端测序read中含有的N的含量超过该条read长度比例的10%时，需要去除此对paired reads；
  • 当单端测序read中含有的低质量（Q ≤ 5）碱基数超过该条read长度比例的 50% 时，需要去除此对paired reads。

数据质量统计概念：

• Raw Base(bp)：原始数据产量，测序序列的个数乘以测序序列的长度，以bp为单位。
• Clean Base(bp)：过滤之后的有效数据量，过滤后测序序列的个数乘以测序序列的长度，以bp为单位。
• Effective Rate(%)：过滤后获得clean data 与raw data的比值。
• Error Rate(%)：碱基错误率。
• GC Content(%)：碱基G和C的数量总和占总的碱基数量的百分比。
• adapter：接头，用于上机测序。建库时引入的接头序列与测序芯片（flow cell）上固定的接头相互识别。
• index：测序的标签，用于测定混合样本，通过每个样本添加的不同标签进行数据区分，鉴别测序样品。
• Q20,Q30：Phred 数值大于20、30的碱基占总体碱基的百分比，其中Phred=-10log10(e),e为错误率。
• raw data/raw reads：测序下机的原始数据。
• clean data/clean reads：对原始数据进行过滤后，剔除了低质量数据的剩余数据。后续分析均基于clean data。

分析软件及用法

Trimmix

• 进行低质量数据adpter污染数据等过滤

java -jar trimmomatic-0.36.jar PE -phred33 YM-2-1712017008_R1.fq.gz YM-2-1712017008_R2.fq.gz YM-2-1712017008.R1.clean.fastq.gz YM-2-1712017008.R1.unpaired.fq.gz YM-2-1712017008.R2.clean.fastq.gz YM-2-1712017008.R2.unpaired.fq.gz ILLUMINACLIP:TruSeq2-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:100

Trimmomatic 发表的文章至今已被引用了 2810 次，是一个广受欢迎的 Illumina 平台数据过滤工具。其他平台的数据例如 Iron torrent ，PGM 测序数据可以用 fastx_toolkit 、NGSQC toolkit 来过滤。

Trimmomatic 支持多线程，处理数据速度快，主要用来去除 Illumina 平台的 Fastq 序列中的接头，并根据碱基质量值对 Fastq 进行修剪。软件有两种过滤模式，分别对应 SE 和 PE 测序数据，同时支持 gzip 和 bzip2 压缩文件。

另外也支持 phred-33 和 phred-64 格式互相转化，现在之所以会出现 phred-33 和 phred-64 格式的困惑，都是 Illumina 公司的锅（damn you, Illumina!），不过现在绝大部分 Illumina 平台的产出数据也都转为使用 phred-33 格式了。

参数说明：

ILLUMINACLIP: 过滤 reads 中的 Illumina 测序接头和引物序列，并决定是否去除反向互补的 R1/R2 中的 R2。
SLIDINGWINDOW: 从 reads 的 5' 端开始，进行滑窗质量过滤，切掉碱基质量平均值低于阈值的滑窗。
MAXINFO: 一个自动调整的过滤选项，在保证 reads 长度的情况下尽量降低测序错误率，最大化 reads 的使用价值。
LEADING: 从 reads 的开头切除质量值低于阈值的碱基。
TRAILING: 从 reads 的末尾开始切除质量值低于阈值的碱基。
CROP: 从 reads 的末尾切掉部分碱基使得 reads 达到指定长度。
HEADCROP: 从 reads 的开头切掉指定数量的碱基。
MINLEN: 如果经过剪切后 reads 的长度低于阈值则丢弃这条 reads。
AVGQUAL: 如果 reads 的平均碱基质量值低于阈值则丢弃这条 reads。
TOPHRED33: 将 reads 的碱基质量值体系转为 phred-33。
TOPHRED64: 将 reads 的碱基质量值体系转为 phred-64。

FastQC

• 进行质控数据统计展示的软件

fastqc -o output dir -f fastq

解压后，summary.txt文件对结果做了简单的总结：

解压后，查看html格式的结果报告。结果分为如下几项