LongQC | 長讀長 Pacbio/Nanopore 測序數據質控神器
LongQC: A Quality Control Tool for Third Generation Sequencing Long Read Data
由于PacBio RS II 、PacBio sequel 和 Nanopore測序質量欠佳，如果能對長序列進行質控，則也能對測序原始數據具有初步把握。
推進LongQC的原因：1. 運行簡單 2. 結果展示美觀。
安裝
安裝依賴包
conda install h5pyconda install -c bioconda pysamconda install -c bioconda edlibconda install -c bioconda python-edlib
安裝longQC & minimap
cd /path/to/download/git clone https://github.com/yfukasawa/LongQC.gitcd LongQC/minimap2-coverage && makecd ..LongQCPath=`pwd`
由于部分python版本問題，可能還需要額外安裝sklearn的python包。
pip install sklearn
參數詳解
python $LongQCPath/longQC.py sampleqc -h
#usage: LongQC sampleqc [-h] -o OUT -x preset [-t] [-n NSAMPLE] [-s SUF] [-c TRIM] [--adapter_5 ADP5] [--adapter_3 ADP3] [-f] [-p NCPU] [-d] [-m MEM] [-i INDS] [-b] input #input 輸入待質控數據的格式，如fasta、fastq或pbbam#-o OUT 輸出文件的路徑#-x preset 輸入待質控數據的測序類型，以便自動提供數據的接頭和olvp參數，測序類型有:pb-rs2, pb-sequel, pb-hifi, ont-ligation, ont-rapid, ont-1dsq#-t 輸入轉錄本、RNA或CDA序列文件（如果有）#-n NSAMPLE 樣本序列的數量，默認為5000#-s SUF 輸出文件的前綴名字#-c TRIM 輸入去掉的接頭的保留路徑，如果沒有輸入，接頭不會保留。#--adapter_5 ADP5 輸入5'接頭#--adapter_3 ADP3 輸入3'接頭#-f 關閉一些sensitive的設置，能運行更快但是準确性下降#-p NCPU 輸入線程數，默認4個線程，線程數要求>=4#-d, --db 輸入minimap2比對的db庫#-m MEM, --mem MEM 亞數據集的内存限制，指定為千兆字節（>0和<2)，默認是0.5千兆字節#-i INDS, --index INDS 給出minimap2 (-I)的索引大小，單位為bp。在小内存的服務器上運行時這個參數要減少。默認是4g。#-b, --short 對于非常短和高錯誤的reads(500bp)，開啟高靈敏度設置。
運行
采用hifi測序原始數據試一試
python $LongQCPath/longQC.py sampleqc -x pb-hifi -o out_dir input_reads.bam
結果展示
analysis/ #minimap比對結果figs/ #分析圖片logs/ #分析日志longqc_sdust.txt #QV的統計情況QC_vals_longQC_sampleqc.jsonweb_summary.html #可視化結果
結果解讀
質控的結果概要
1. Sample name：測序樣品的名字，可通過參數-s設置2. Yield：堿基的總數量3. Number of reads：reads數量4. Q7 bases：堿基質量值（QV）大于等于7的堿基的比例 （由于是Pacbio 測序數據，所以并沒有提供QV）5. Longest read：最長reads的長度6. Estimated non-sense read fraction：評估出來的non-sense reads所占比例，應該少于30%。比例高表明: 1. 測序數據有問題 2. 覆蓋度不高
- 接頭信息
1. Number of trimmed reads：在reads末端中發現的可能是接頭（75%以上的一緻性）的序列長度；如果出乎意料地低或者沒有檢測到接頭，可能是在接頭連接步驟上有問題。2. Max seq identity：接頭序列與序列之間的最大的一緻性值。如果接頭仍然存在于數據集中，這個值應該相當高(90%)。3. Average trimmed length：接頭序列在reads中的平均末端位置
Reads長度
（Nanopore 左，PacBio 右）
這個圖展示是reads長度分布，三代測序數據會顯示單獨的高峰。如果是轉錄組測序數據，會有嚴格的大小選擇的數據（如PacBio數據）會向右移動顯示更高的峰。
1. Mean read length：樣品的reads預測的長度
2. N50：N50值
Per-reads的質量評估
（圖1 nanopore ，圖2 PacBio）
上圖将顯示了per reads的QV分布，y軸為平均QV值，x軸是reads的長度集。
理想情況下，短reads和長reads的QV分布應該相似，中值(綠色橫線)應該高于7。
圖1是開發者用DASCRUBBERR軟件計算出QV後，再進行質控的結果。
值得注意的是，圖2 由PacBio測序的序列是無法計算出QV。主要是因為PacBio的測序文件質量值是由"!"作占位符，并未有實質的QV意義。
Per read 覆蓋度
上圖為Per read 覆蓋度分布圖。如果數據沒有問題，應該是單峰（圖1）。如果是多峰值，可以考慮是否基因組或者材料的雜合情況比較高（宏基因組數據除外）。
LongQC使用GMM(用于基因組)或混合Gaussian 和lognorm 分布(用于轉錄組)自動檢測這樣的峰值，以便從背景水平中區别真實的峰值。
上圖為覆蓋度與reads長度的關系圖，用來檢查覆蓋範圍是否有意外波動。在基因組測序數據中，覆蓋範圍的波動應該在一定範圍内(默認為3 σ)。如果超出這個範圍，表明數據存在污染，文庫質量低，PacBio超載等問題。
圖1表明存在微弱波動（可忽略），圖2是正常的分布圖。
non-sense reads
(Nanopore 左，PacBio 右)
正常情況下，Normal reads 比non-sense reads的顯示更高的值。
理想情況下，正常reads的中值（橙色線）應該高于7。
有3種情況：
注意：
這個部分顯示的中值覆蓋率可能比使用參考基因組進行質控的結果要小。由于将reads比對到未更正的容易出錯的整個測序數據集上時，比對是不太敏感的，因此覆蓋率會受到這種不太敏感的結果的影響。由于上述影響，估計的基因組/轉錄組大小往往大于實際大小。一般來說，更好的數據提供更好的尺寸估計。
GC分布
上圖使用了兩種不同的方法計算GC含量，其中藍色是對end to end全長的reads進行計算，而紅色是對分塊的亞數據集進行計算。
藍色顯示更清晰的分布，因為使用更長的序列使結果受到更小的偏差。由于測序平台不同也會導緻相同的數據的reads長度不同，從而導緻全長reads不同，而紅色的分塊亞數據集對長度已經标準化，所以對測序或樣本差異更不敏感，更可靠。
Flanking region分析
可以用來檢查特定序列的連接，如接頭。如果沒有人工序列（如接頭），在兩個圖中峰值應該顯示在0的位置。如果觀察到類似于接頭的序列，則用紅色虛線标出平均長度。第一張圖表示有接頭，第二張圖表示沒有接頭
序列完整性
這裡通過DUST算法顯示序列的完整性分數
背景知識
三代測序的質控問題
二代測序采用邊合成邊測序的方法，以Illumina為例，進行橋式擴增後，在Flow cell的固相表面上獲得上百萬條成簇分布的雙鍊待測片段。然後在Flow cell中加入熒光标記的dNTP和酶，由引物起始開始合成子鍊。由于dNTP存在 3'-OH 以疊氮基團RTG修飾，這會阻礙子鍊延伸，使得每個循環隻能測得一個堿基。測完一個堿基後， Flow cell 通入液體洗掉多餘的dNTP和酶，使用顯微鏡的激光掃描特征熒光信号，根據圖像的顔色即可識别該堿基，再用化學試劑将疊氮基團與熒光基團洗去。重複此過程即可識别全部堿基。在進行堿基識别的同時，軟件還根據預先設定的規則對每一個堿基進行質量評估，為每一個堿基賦予一個質量值。不同的二代測序平台，如果其測序原理不同，則定義質量值QV所取的參數也不同。一個熒光信号被測序軟件判讀為某一種堿基，存在一定的概率(Pe)發生判斷錯誤。QV就是用來衡量Pe高低的一個定量指标，其定義為：
QV=-10 x lg Pe
根據這一定義可知，
如果堿基識别的準确率是99%，則錯誤率是1%，QV是20；
如果堿基識别的準确率是99.9%，則錯誤率是0.1%，QV是30；
如果堿基識别的準确率是99.99%，則錯誤率是0.01%，QV是40
FASTQ格式的每第四行表示這條序列的質量值，用ACSII碼表示。所以二代測序的質控軟件可以根據質量值QV來評估測序數據的質量。
而三代測序由于其随機分布的堿基錯誤率，單堿基的準确性不能直接用來衡量數據的質量，所以Pacbio Sequel測序平台的下機數據不會提供單堿基質量值QV。此外，三代測序雖然讀長長，但是錯誤率高（~15%），而且跟二代測序所用的技術原理不同。二代測序的質控軟件如FastQC（依賴QV評估測序數據）不适合三代測序的質控。當三代測序數據不提供QV，且沒有參考基因組時，要怎麼進行質控呢？
LongQC的介紹
LongQC質控軟件可以解決以上的問題。LongQC不依賴QV和參考基因組進行質控，它依靠覆蓋度（Coverage)和non-sense reads對測序數據進行評估。具體的原理和步驟如下：
LongQC中的一些解釋
non-sense reads
non-sense reads是指沒有比對到整個測序數據集的reads，這些reads在整個數據中的分布是任意且分散的，開發者認為這些reads是測序技術或建庫中導緻的錯誤reads，或者是污染的reads。通過計算這些non-sense reads在整個測序數據集所占的比例，可以評估測序數據的質量。
堿基質量值QV
由于Sequel測序數據沒有提供堿基質量信息，所以在longQC質控結果中關于QV的評估都是0。而Nanopore測序數據有Q-scores，Q-Scores是用來說明一個給定的堿基比它周圍的堿基更好還是更差的依據。開發者通過DASCRUBBER軟件對數據生成了QV，再進行質控後可以看到QV的解析結果。
最新文章
-
如何解讀世界名表級别劃分
11分鐘前 -
魏新說法話說高空抛物和墜物
12分鐘前 -
鄧中甲方劑學講稿—第四章清熱劑—清氣分熱—白虎湯
12分鐘前 -
案山朝山
15分鐘前 -
讀賈松禅《霍去病》第十五回——淺談火絨的見識和火絨的胸懷
15分鐘前 -
美國空軍一級司令部現任指揮官一覽
16分鐘前 -
推理小說五大流派，你都細品了麼
16分鐘前 -
想不到！這座2500年的江南古城，曾經比紐約還富，如今卻淪為三線
18分鐘前
熱門文章
-
威士忌高階關于 OB 與 IB，一次給你說透！
6小時前 -
每日一誦傷寒論第241條
1周前 -
老張老李侃門球之140篇
2天前 -
這個穴位可以治療多種胃痛腹痛，還可以減肥
6天前 -
診餘雜記（師傳經驗）
4天前 -
美麗中國-2870：中國最大的内陸河，塔裡木河
2天前 -
為什麼五點鐘要起床答案讓人吃驚！（現在知道還不晚）
1周前 -
治療坐骨神經痛藥酒５方
1周前
有話要說...