最近在分析一個重測序的SV過程中出現(xiàn)了一些問題,在使用lumpyExpression分析時,獲得disordants Reads的時候出現(xiàn)異常,幾乎獲得的全部的mappingReads,原因探查過程如下:
提取disordants Reads的時候使用的方法是samtools,根據(jù)每一條reads的FLAG進行判斷,具體代碼如下:
samtools view -b -F 1294 sample.bam > sample.discordants.bam
-F代表過濾掉對應FLAG的reads,1294代表reads的FLAG情況,具體包含的reads如下,在這一步驟的核心是要過濾掉“read mapped in proper pair"的reads,即左右兩端mapping在一致區(qū)域的reads
這個html為一個小程序,如上面的截圖,輸入Flag的編號可以得到對應哪些類型的Reads
在BWA mem的軟件參數(shù)中包含一個 -P參數(shù),具體介紹如下:
-P In the paired-end mode, perform SW to rescue missing hits only but do not try to find hits that fit a proper pair.
說明如果輸入-P參數(shù),在mapping的過程中或跳過fit a proper pair的步驟,也就意味著結果不會給出read mapped in proper pair這樣的一個FLAG,而我們平臺現(xiàn)在的BWA mem代碼包含這樣一個參數(shù),如下:
具體參數(shù)測試結果如下:
包含 -P 參數(shù):
同樣的的reads的FLAG變成了83和163,代表含義如下,包含read mapped in proper pair的注釋
由此可見-P參數(shù)會影響每一條reads的FLAG,由于這種read mapped in proper pair在SV的分析過程中的判斷是十分重要的,所以建議刪除掉平臺BWA mem模塊的-P參數(shù)
經(jīng)過測試,BWA -p參數(shù)除了影響FLAG外還會影響readsMapping的位置,在CallSNP過程中會造成很大的偏差,測試結果如下:
下圖為同一條序列在不同mapping方法中的位置差異,上面為添加-p參數(shù)的結果,下面為不添加-p的結果,可以發(fā)現(xiàn),不添加-p組中,每一條reads的左右兩端均可以匹配在臨近區(qū)段,且均在第二號染色體,而在添加-p組中,reads幾乎全被mapping到9號染色體上,且左右兩端reads的距離差距非常大,IGV截圖(截圖中為chr9對應區(qū)域,上面為添加-p參數(shù),下面給為不添加-p參數(shù))中也可以發(fā)現(xiàn)添加-p參數(shù)計算出的SNP位點存在明顯的假陽性現(xiàn)象,再次證明-p參數(shù)需要被移除