高通量数据处理流程

数据处理流程计算机集群测序控制PC测序仪1.控制测序过程决定测序长度、填加试剂、控制温度、控制反应时间、拍照2.图像分析对测序仪拍照的图片进行图像分析，得到亮点的光强度和坐标3.basecalling由光强度得到碱基序列4.数据传输将basecalling结果（二进制文件bcl）传输到计算机集群的存储上5.数据处理将bcl文件转化为后续信息分析所使用的文本文件（fastq,qseq）6.index拆分7.数据质量分析8.数据备份9.后续信息分析......图像分析及basecalling基本原理1234CCCCGGGTTAAACycle1Cycle2Cycle3对A发出的光拍照对C发出的光拍照对G发出的光拍照对T发出的光拍照图像分析及basecalling基本原理1234TCCAATGCACGGCycle1Cycle2Cycle3由4个cluster得到4条序列：①ATA...②CCT...③GCG...④GAC...图像分析对每个图片独立的处理图像锐化对图片进行快速傅里叶变换(FFT)，在傅里叶空间乘以滤波函数后反变换识别cluster/亮点信噪比(亮度/背景值)大于阈值的亮点区域计算亮点光强度和位置坐标在亮点区域，对光强度进行二维插值，求出最大光强度，以及最大光强度对应的位置坐标将同一个tile的所有图片中的亮点坐标对齐重叠不同图片之间存在偏移/拉伸/压缩(offset)系统、稳定的：4种光折射率不同、滤波片不同、光路不同，所以造成成像的偏移/拉伸/压缩，可利用crosstalk，计算出偏移/拉伸/压缩的数值(offset参数)偶然、随机的：flowcell表面不平、自动调整焦距、机械移动不够精确、随机振动，可利用crosstalk解决CrosstalkAC光谱间有交叠，GT光谱间有交叠，所以：碱基A的图片中包含C发出的光碱基C的图片中包含A发出的光碱基T的图片中包含G发出的光不利：不能直接比较光强度大小而得到碱基有利：利用图片中共同的亮点，将所有图片对齐重叠，解决offset问题ACAC光谱光谱滤波片图像分析结果Cluster坐标Cycle1Cycle2Cycle3......(1053,1543)A1000A500A23......C800C900C26......G20G-18G500......T24T-12T300......(1923,1723)A-13A1000A33......C-12C800C24......G40G25G500......T700T20T300......(1032,1231)A500A50A1000......C900C40C800......G26G500G-20......T32T300T-20....................................Crosstalk固有发光能力不同，需要归一化BasecallingCrosstalk校正4种光强度归一化（用DNA样品计算参数）Phasing/Prephasing校正（用DNA样品计算参数）对于每个cluster：在每个cycle中，比较4种光强度，光强度最大的就是当前cycle测到的碱基，各cycle测到的碱基连起来组成这个cluster的碱基序列；计算每个碱基的质量值ACGTphasingSequencingprimerprephasingBasecalling结果：qseq文件•每一行表示一条reads(一个cluster)•每行有11列，tab分隔：机器编号、run序号、Lane号、Tile号、X坐标、Y坐标、index标志、read1/read2标志、碱基序列、质量序列、是否通过默认的质量筛选标准•Single-end(SE)测序：1个qseq文件•Pair-end(PE)测序：2个qseq文件分别存放read1和read2的数据；2个文件的同一行属于同一个cluster•每条序列(reads)长度=上机测序循环(cycle)数量；•测序cycle数量受测序试剂盒的试剂量限制，对于GA有：36SE、36+7/8SEindex、45PE、36+7+45PEindex、76PE、74+7+76PEindex、73+8+76PEindex、101PE、101+7/8+101PEindex等对于Hiseq：91PE、91+8+91PEindex、101PE、101+8+101PEindexBasecalling结果：qseq文件列意义表示法1机器编号2run序号整数3Lane号整数1到84Tile号整数5X坐标整数6Y坐标整数7index标志index序列或者“0”8read1/read2标志1表示read1；2表示read29碱基序列大写ACGT和.：ACCCAACTCATCTGAAACA10质量序列每个碱基有一个质量值，用字符表示：字符的ASCII码值-64=质量值bbbc`bb_bb_aSa`V]`\11是否通过默认的质量筛选标准1表示通过；0表示不通过，质量差fastq文件•每4行表示一条reads(一个cluster)第一行：@序列ID，包含index序列及read1或read2标志：第二行：碱基序列，大写“ACGTN”第三行：“+”，省略了序列ID第四行：质量值序列：字符的ASCII码值-64=质量值•Single-end(SE)测序：1个fastq文件•Pair-end(PE)测序：2个fastq文件分别存放read1和read2的数据；Read1的fastq文件*1.fq中第一条reads：@FC61FL8AAXX:1:17:1012:19200#GCCAAT/1CCACTGTCATGTGAACATCACAGAGACATTTCTTGA+bbbbbbbbbbabbbbbbbbbbbbbbaaaaaaaaa_\Read2的fastq文件*2.fq中第一条reads：@FC61FL8AAXX:1:17:1012:19200#GCCAAT/2AAAATTAGCCAGGCAATGGTGGTGCATGCCTTTAATCCCAGCTA+`QVVV``V``````YVYWWYPWYYTYYWUYYYVV```````WW`•质量值@FC61FL8AAXX:1:17:1012:19200#GCCAAT/1CCACTGTCATGTGAACATCACAGAGACATTTCTTGA+bbbbbbbbbbabbbbbbbbbbbbbbaaaaaaaaa_\表示方法Illumina：字符的ASCII值-64=质量值(Sanger：字符的ASCII值-33=质量值)范围GAIllumina1.3+(09年3月之后):[2,35][B,c]GAIllumina1.0(09年3月之前):[-5,40][;,h]Hiseq:[2:38][B,f]质量值与错误率理论关系：Q=-10log10(e)质量值计算方法：根据光强信号信噪比、光强度衰减、GC含量等参数，计算质量值fastq文件fastq文件•每条序列(reads)长度read1和read2分别去除了最后一个碱基，即：36SE有效长度为35101PE有效长度为100(read1)+100(read2)101+8+101PEindex有效长度为100(read1)+100(read2)Read1中所有reads长度相同，Read2中所有reads长度相同，但是Read1和Read2长度可以不相同，取决于上机测序循环(cycle)数量•质量筛选(PF)：Illumina标准流程输出的fastq文件，去除了qseq文件中没有通过默认质量筛选标准的低质量序列(reads)GA正常PF比例：DNA80~90%，RNA70~85%•每个Lane的正常产量范围：GA20~30MPFreads—Read1和Read2各有20~30MHiseq60~80MPFreads—Read1和Read2各有60~80M碱基总产量=Read1的产量+Read2的产量=reads数量＊(Read1的长度+Read2的长度)fastq文件产量（GA）样品类型上机测序类型平均产量(Gbp/Lane)正常产量(Mreads/Lane)SmallRNA36SE0.7020（2个样品混合上1个lane）表达谱36SE0.7020（4个样品混合上1个lane）Chip-seq36SE0.8825MeDIP-seq45PE2.2025转录组76PE3.8125（2个样品混合上1个lane）Denovo101PE6.0030Denovo76PE4.5030Denovo45PE2.6430外显子76PE4.5030重测序101PE6.0030重测序76PE4.5030重测序45PE2.6430Meta76PE4.5030甲基化76PE3.7525甲基化45PE2.2025fastq文件产量（HiseqvsGA）GAHiSeqFC面积mm^2/FC5101440tile面积mm^2/tile0.531255.625tile/Lane12032rawcluster(万/tile)28250270310350PF85%87%87%87%85%PFcluster(万/tile)23.8218235270298cluterK/mm^2527444480551622readsM/Lane28.669.675.286.395.2101PEGb/Lane5.713.915.017.319.091PEGb/Lane12.513.515.517.1文库质控问题1：Pair-end关系800bp及以下文库5'3'3'5'PCRprimer1反向互补(包含5'adapter反向互补)PCRprimer1(包含5'adapter)PCRprimer2(包含3'adapter反向互补)PCRprimer2反向互补(包含3'adapter)Read1测序Read2测序5'5'3'3'•与参考序列比较或者：总之，Read1，Read2与参考序列比对结果：一正(F)一反(R)，且F的位点坐标小于R的位点坐标Read1Read2参考序列正向5'5'3'3'参考序列反向互补Read2Read1参考序列正向5'5'3'3'参考序列反向互补FR参考序列正向5'5'3'3'参考序列反向互补•总之，Read1，Read2于参考序列比对结果：一正(F)一反(R)，且F的位点坐标小于R的位点坐标•文库插入片段长度FR参考序列正向5'5'3'3'参考序列反向互补Insert-sizeInsert-size5'PCRprimer1PCRprimer23'3'5'文库分子长度2100检测报告文库长度分布与参考序列比对得到insert-size分布正常insert-size分布基因组DNA外显子PCR-free文库异常insert-size分布2K及以上文库文库质控问题1：Pair-end关系•与参考序列比较或者：总之，Read1，Read2于参考序列比对结果：一正(F)一反(R)，且F的位点坐标大于R的位点坐标Read1Read2参考序列正向5'5'3'3'参考序列反向互补Read2Read1参考序列正向5'5'3'3'参考序列反向互补FR参考序列正向5'5'3'3'参考序列反向互补正常insert-size分布2K5~6KPCR-free文库10K异常insert-size分布文库问题2：adapter污染•空载：adapter与adapter直接连接，中间没有插入片段，导致read1测到3'adapter，read2测到5'adapter的反向互补reads尾部测到adapter•插入片段过短插入片段长度小于上机测序循环(cycle)数，导致read1尾部测到3'adapter，read2尾部测到5'adapter的反向互补5'3'3'5'PCRprimer1反向互补(包含5'adapter反向互补)PCRprimer1(包含5'adapter)PCRprimer2(包含3'adapter反向互补)PCRprimer2反向互补(包含3'adapter)adapt