高通量數(shù)據(jù)處理流程

1、數(shù)據(jù)處理流程計算機集群測序控制PC測序儀1.控制測序過程控制測序過程決定測序長度、填加試劑、控制溫度、控制反應(yīng)時間、拍照2.圖像分析圖像分析對測序儀拍照的圖片進行圖像分析，得到亮點的光強度和坐標3.basecalling由光強度得到堿基序列4.數(shù)據(jù)傳輸數(shù)據(jù)傳輸將basecalling結(jié)果（二進制文件bcl）傳輸?shù)接嬎銠C集群的存儲上5.數(shù)據(jù)處理數(shù)據(jù)處理將bcl文件轉(zhuǎn)化為后續(xù)信息分析所使用的文本文件（fastq,qseq）6.index拆分拆分7.數(shù)據(jù)質(zhì)量分析數(shù)據(jù)質(zhì)量分析8.數(shù)據(jù)備份數(shù)據(jù)備份9.后續(xù)信息分析后續(xù)信息分析.圖像分析及basecalling基本原理1234CCCCGGGTTAAACyc

2、le1 Cycle2 Cycle3 對A發(fā)出的光拍照對C發(fā)出的光拍照對G發(fā)出的光拍照對T發(fā)出的光拍照圖像分析及basecalling基本原理1234TCCAATGCACGGCycle1 Cycle2 Cycle3 由4個cluster得到4條序列： ATA. CCT. GCG. GAC.圖像分析對每個圖片獨立的處理圖像銳化對圖片進行快速傅里葉變換(FFT)，在傅里葉空間乘以濾波函數(shù)后反變換識別cluster/亮點信噪比(亮度/背景值)大于閾值的亮點區(qū)域計算亮點光強度和位置坐標在亮點區(qū)域，對光強度進行二維插值，求出最大光強度，以及最大光強度對應(yīng)的位置坐標將同一個tile的所有圖片中的亮點坐標對齊

3、重疊不同圖片之間存在偏移/拉伸/壓縮(offset)系統(tǒng)、穩(wěn)定的：4種光折射率不同、濾波片不同、光路不同，所以造成成像的偏移/拉伸/壓縮，可利用crosstalk，計算出偏移/拉伸/壓縮的數(shù)值(offset參數(shù))偶然、隨機的：flowcell表面不平、自動調(diào)整焦距、機械移動不夠精確、隨機振動，可利用crosstalk解決CrosstalkAC光譜間有交疊，GT光譜間有交疊，所以：堿基A的圖片中包含C發(fā)出的光堿基C的圖片中包含A發(fā)出的光堿基T的圖片中包含G發(fā)出的光不利：不能直接比較光強度大小而得到堿基有利：利用圖片中共同的亮點，將所有圖片對齊重疊，解決offset問題ACAC光譜光譜濾波片圖像分

4、析流程Template Generation利用AC之間的crosstalk、GT之間的crosstalk，將前 2個cycle的圖片，與第一個cycle的A的圖片對齊重疊，確定所有cluster的位置坐標(x,y)Registration and Intensity Extraction對于每一個cycle：將4張圖片(ACGT)中的所有亮點與cluster坐標(x,y)對應(yīng)，計算每個cluster的4種光強度對每個圖片獨立的處理圖像銳化、識別cluster/亮點、計算亮點光強度和位置坐標Template Generation利用AC之間的crosstalk、GT之間的crosstalk，將

5、前 2個cycle的圖片，與第一個cycle的A的圖片對齊重疊，確定所有cluster的位置坐標(x,y)ACGTAGTCCycle1 Cycle21243ACycle1 A Cycle1 C Cycle1 ACCACCycle2 A Cycle2 C Cycle1 ACACAC123GTCycle2 G Cycle2 T Cycle1 ACGTACGTCycle1 G Cycle1 TGT12431243圖像分析結(jié)果Cluster坐標坐標Cycle1Cycle2Cycle3.(1053,1543)A 1000A 500A 23.C 800C 900C 26.G 20G -18G 500.T

6、24T -12T 300.(1923,1723)A -13A 1000A 33.C -12C 800C 24.G 40G 25G 500.T 700T 20T 300.(1032,1231)A 500A 50A 1000.C 900C 40C 800.G 26G 500G -20.T 32T 300T -20.Crosstalk固有發(fā)光能固有發(fā)光能力不同，力不同，需要歸一化需要歸一化BasecallingCrosstalk 校正4種光強度歸一化（用DNA樣品計算參數(shù)）Phasing/Prephasing校正（用DNA樣品計算參數(shù)）對于每個cluster：在每個cycle中，比較4種光強度，光強

7、度最大的就是當前cycle測到的堿基，各cycle測到的堿基連起來組成這個cluster的堿基序列；計算每個堿基的質(zhì)量值4種光強度歸一化2種光強度分布2種光強度分布ACGTphasingSequencing primerprephasingBasecalling結(jié)果：qseq文件 每一行表示一條reads(一個cluster) 每行有11列，tab分隔：機器編號、run序號、Lane號、Tile號、X坐標、Y坐標、index標志、read1/read2標志、堿基序列、質(zhì)量序列、是否通過默認的質(zhì)量篩選標準 Single-end(SE)測序：1個qseq文件 Pair-end(PE)測序：2個qs

8、eq文件分別存放read1和read2的數(shù)據(jù)；2個文件的同一行屬于同一個cluster 每條序列(reads)長度=上機測序循環(huán)(cycle)數(shù)量； 測序cycle數(shù)量受測序試劑盒的試劑量限制， 對于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結(jié)果：qseq文件列列意義意義表示法表示法1機器編號2ru

9、n序號整數(shù)3Lane號整數(shù)1到84Tile號整數(shù)5X坐標整數(shù)6Y坐標整數(shù)7index標志index序列或者“0”8read1/read2標志1表示read1；2表示read29堿基序列大寫ACGT和.：ACCCAACTCATCTGAAACA10質(zhì)量序列每個堿基有一個質(zhì)量值，用字符表示：字符的ASCII碼值-64=質(zhì)量值bbbcbb_bb_aSaV11是否通過默認的質(zhì)量篩選標準1表示通過；0表示不通過，質(zhì)量差fastq文件 每4行表示一條reads(一個cluster)第一行：序列ID，包含index序列及read1或read2標志：第二行：堿基序列，大寫“ACGTN”第三行：“+”，省略了序列

10、ID第四行：質(zhì)量值序列：字符的ASCII碼值-64=質(zhì)量值 Single-end(SE)測序：1個fastq文件 Pair-end(PE)測序：2個fastq文件分別存放read1和read2的數(shù)據(jù)；Read1的fastq文件 *1.fq中第一條reads：FC61FL8AAXX:1:17:1012:19200#GCCAAT/1CCACTGTCATGTGAACATCACAGAGACATTTCTTGA+bbbbbbbbbbabbbbbbbbbbbbbbaaaaaaaaa_Read2的fastq文件 *2.fq中第一條reads：FC61FL8AAXX:1:17:1012:19200#GCCAAT

11、/2AAAATTAGCCAGGCAATGGTGGTGCATGCCTTTAATCCCAGCTA+QVVVVYVYWWYPWYYTYYWUYYYVVWW 質(zhì)量值FC61FL8AAXX:1:17:1012:19200#GCCAAT/1CCACTGTCATGTGAACATCACAGAGACATTTCTTGA+bbbbbbbbbbabbbbbbbbbbbbbbaaaaaaaaa_表示方法 Illumina：字符的ASCII值 - 64 = 質(zhì)量值 (Sanger：字符的ASCII值 - 33 = 質(zhì)量值)范圍 GA Illumina1.3+(09年3月之后): 2,35 B,c GA Illumina1

12、.0 (09年3月之前): -5,40 ;,h Hiseq: 2:38 B,f質(zhì)量值與錯誤率理論關(guān)系： Q =-10 log10(e)質(zhì)量值計算方法：根據(jù)光強信號信噪比、光強度衰減、GC含量等參數(shù)，計算質(zhì)量值fastq文件fastq文件 每條序列(reads)長度read1和read2分別去除了最后一個堿基，即：36 SE 有效長度為 35101 PE 有效長度為 100(read1)+100(read2)101+8+101 PEindex 有效長度為 100(read1)+100(read2)Read1中所有reads長度相同，Read2中所有reads長度相同，但是Read1和Read2長

13、度可以不相同，取決于上機測序循環(huán)(cycle)數(shù)量 質(zhì)量篩選(PF)：Illumina標準流程輸出的fastq文件，去除了qseq文件中沒有通過默認質(zhì)量篩選標準的低質(zhì)量序列(reads)GA正常PF比例：DNA 8090%，RNA 7085% 每個Lane的正常產(chǎn)量范圍：GA 2030M PF reads Read1和Read2各有2030MHiseq 6080M PF reads Read1和Read2各有6080M堿基總產(chǎn)量 = Read1的產(chǎn)量 + Read2的產(chǎn)量 = reads數(shù)量(Read1的長度 + Read2的長度)fastq文件產(chǎn)量（GA）樣品類型樣品類型上機測序上機測序類型

14、類型平均產(chǎn)量平均產(chǎn)量(Gbp/Lane)正常產(chǎn)量正常產(chǎn)量(Mreads/Lane)SmallRNA36SE0.70 20 （2個樣品混合上1個lane）表達譜36SE0.70 20 （4個樣品混合上1個lane）Chip-seq36SE0.88 25 MeDIP-seq45PE2.20 25 轉(zhuǎn)錄組76PE3.81 25 （2個樣品混合上1個lane）De novo101PE6.00 30 De novo76PE4.50 30 De novo45PE2.64 30 外顯子76PE4.5030重測序101PE6.00 30 重測序76PE4.50 30 重測序45PE2.64 30 Meta76

15、PE4.50 30 甲基化76PE3.75 25 甲基化45PE2.20 25 fastq文件產(chǎn)量（Hiseq vs GA）GAHiSeqFC面積面積 mm2/FC5101440tile面積面積mm2 /tile0.531255.625tile /Lane12032raw cluster(萬萬/tile)28 250 270 310 350 PF85%87%87%87%85%PF cluster(萬萬/tile)23.8 218 235 270 298 cluter K/mm2527 444 480 551 622 reads M/Lane28.6 69.6 75.2 86.3 95.2 10

16、1PE Gb/Lane5.7 13.9 15.0 17.3 19.0 91PE Gb/Lane12.5 13.5 15.5 17.1 文庫質(zhì)控問題1：Pair-end關(guān)系800bp及以下文庫5335PCR primer1反向互補(包含5adapter反向互補)PCR primer1(包含5adapter)PCR primer2(包含3adapter反向互補)PCR primer2反向互補(包含3adapter)Read1測序測序Read2測序測序5533 與參考序列比較或者：總之，Read1，Read2與參考序列比對結(jié)果：一正(F)一反(R)，且F的位點坐標小于R的位點坐標Read1Read2

17、參考序列正向5533參考序列反向互補Read2Read1參考序列正向5533參考序列反向互補FR參考序列正向5533參考序列反向互補 總之，Read1，Read2于參考序列比對結(jié)果：一正(F)一反(R)，且F的位點坐標小于R的位點坐標 文庫插入片段長度FR參考序列正向5533參考序列反向互補Insert-sizeInsert-size5PCR primer1PCR primer2335文庫分子長度2100檢測報告文庫長度分布與參考序列比對得到insert-size分布正常insert-size分布基因組DNA外顯子PCR-free文庫異常insert-size分布2K及以上文庫文庫質(zhì)控問題1：

18、Pair-end關(guān)系 與參考序列比較或者：總之，Read1，Read2于參考序列比對結(jié)果：一正(F)一反(R)，且F的位點坐標大于R的位點坐標Read1Read2參考序列正向5533參考序列反向互補Read2Read1參考序列正向5533參考序列反向互補FR參考序列正向5533參考序列反向互補正常insert-size分布2K56KPCR-free文庫10K異常insert-size分布文庫問題2：adapter污染 空載：adapter與adapter直接連接，中間沒有插入片段，導(dǎo)致read1測到3adapter，read2測到5adapter的反向互補reads尾部測到adapter 插入

19、片段過短插入片段長度小于上機測序循環(huán)(cycle)數(shù)，導(dǎo)致read1尾部測到3adapter，read2尾部測到5adapter的反向互補5335PCR primer1反向互補(包含5adapter反向互補)PCR primer1(包含5adapter)PCR primer2(包含3adapter反向互補)PCR primer2反向互補(包含3adapter)adapter空載較多導(dǎo)致堿基含量波動客戶PCR引物污染導(dǎo)致堿基含量波動文庫質(zhì)控問題3：文庫隨機性 GC含量偏差： 實驗技術(shù)（打斷、PCR、測序）本身特點，導(dǎo)致高GC和低GC區(qū)域測序覆蓋度偏低，甚至某些區(qū)域覆蓋不到； PCR-free建庫技術(shù)可減少PCR帶來的隨機性問題 duplication PCR擴增出很多一模一樣的母版分子，測序結(jié)果中很多條reads是一樣的； 基因組自身重復(fù)序列含量高導(dǎo)致duplication偏高； 數(shù)據(jù)量越大，duplication比例越高文庫質(zhì)控問題4：其它物種、樣品污染測序質(zhì)控問題 raw Clust