1、研究背景
試劑污染和樣品交叉污染是實驗室中長期公認的一大問題，污染可能會致使實驗結果的不許確。先前研究中曾在TCGA的RNA-seq數據中發現H-HPV18（Hela細胞系衍生的HPV-18病毒），所以本文做者但願證明TCGA的RNA-seq樣本被H-HPV18病毒污染，並肯定其污染來源。web
2、分析流程
3、結果解讀
一、開發VirDetect：一種特異性病毒檢測軟件
爲了在RNA-seq數據中檢測病毒，做者開發了一個開源軟件VirDetect，該軟件使用STARV2.4 aligner進行RNA-seq reads的比對，將沒有與人類基因組比對上的數據與病毒基因組數據進行比對（圖1.A），做者在圖1中對VirDetect的檢測準確性進行了驗證。算法
爲了提升測序質量，避免資源浪費，做者將病毒基因組中的人類同源性區域和低複雜度區域進行了屏蔽（mask），做者取滑動窗口長度爲75個核苷酸，作滑動窗口算法，以93%的核苷酸類似性做爲同源片斷的標準，將這些同源片斷屏蔽（流程如圖1.B所示）。做者對屏蔽效果進行驗證，發如今對人類同源性片斷和低複雜度片斷同時屏蔽後，aligner不會再將reads比對到屏蔽片斷上（圖1.C）；當突變數<3時，中位敏感率（按讀取到病毒基因組的reads比例計算）爲99.6%；而當突變數>3時，中位敏感率呈線性降低，但圖1.E顯示中位PPV（陽性預測率，按讀取到正確基因組的比例計算）值爲97%，即便腫瘤突變負擔較高，利用VirDetect將病毒reads讀取到正確的基因組的機率仍然很高。數據庫
二、利用VirDetect查找TCGA中的病毒污染
首先做者分析了北卡羅來納大學對TCGA樣本的測序結果，繪製了多種病毒在TCGA28種腫瘤中的表達率熱圖（圖2）。圖中的大部分結果和預期一致，如HBV（乙肝病毒）在肝癌中廣泛存在，HPV16在頭頸部鱗狀細胞癌中廣泛存在等等。但與預期不符的是：做者利用VirDetect在非子宮頸部腫瘤中檢測到了HPV-18病毒（每種腫瘤樣本的具體檢測數值如圖3.A所示），其中比較突出的有透明腎細胞癌（KIRC）、肺鱗狀細胞癌（LUSC）等。做者認爲能檢測到這種數量級的病毒序列，最多是由試劑污染引發的，所以做者但願在後續的研究中繼續尋找污染的源頭。express
圖2.多種病毒在TCGA28種腫瘤中的表達率熱圖微信
先前有研究經過Hela細胞系中的23個特異性SNPs（單核苷酸多態性）證實了非子宮頸腫瘤的HPV-18是來源於Hela細胞的。做者對這一結果加以驗證，圖3.B結果說明除去CSEC（宮頸癌），和3個BLCA（膀胱癌）樣本，其他全部17個非子宮頸癌樣本的SNPs與Hela細胞的HPV-18基因組SNPs都徹底匹配，驗證了非子宮頸腫瘤的HPV-18來源於Hela細胞系app
除了HPV-18外，在96個TCGA樣本中還檢測到了XMV43病毒（如圖3.A下半部分所示），在HPV-18表達量最高的KIRC樣本中，XMV43表達量一樣最高（圖3.C），兩者的Spearman相關係數爲0.44（p=0.006），同時做者在圖3.D中統計了每種腫瘤樣本中同時存在兩種病毒的樣本數，發現他們在BRCA、HNSC、KIRC、LUSC樣本中均有明顯的共表達。所以，做者認爲TCGA樣本的污染物中可能同時含有HPV-18和XMV43兩種病毒的RNAide
圖3.C-D HPV-18和XMV43的共表達狀況ui
三、調查污染來源
爲了調查形成樣本污染的緣由，做者對UHRR（人類通用參考RNA，Universal Human Reference RNA）和大多數TCGA樣本同時進行測序，並監測文庫的構建（圖4.A、C）。大部分UHRR樣本的測序從2010年開始，當時的樣本中並無檢測到HPV-18或XMV43的序列。如圖所示，發現兩種病毒污染的高峯期發生在2011年4-6月左右，做者認爲有可能發生污染的樣本並不必定都是與其餘UHRR（+）樣本在同一天開始構建文庫的，所以不能經過圖4.A、C的結果直接判斷病毒的來源是交叉污染。url
爲了確保乳腺癌基因被加入到UHRR樣本中，做者將實驗室庫存的MCF-7和ME16C兩株乳腺癌細胞系-加入到了UHRR測序樣本中（命名爲UHRR+），在圖4.B中做者利用箱線圖對比了UHRR和UHRR+隊列中兩種病毒序列的計數，發現兩者的HPV-18水平都很高，但UHRR+隊列的XMV43水平明顯高於UHRR，這代表添加到UHRR+中區別於UHRR的兩個細胞株中的一株（MCF-7和ME16C）引入了XMV43病毒污染。spa
爲了判斷XMV43污染是由上述兩株細胞中的哪一株引入的，做者利用實驗室2013年準備的MCF-7和ME16C細胞株的RNA-seq數據檢測兩者的XMV43序列含量，發現MCF-7中的XMV43序列丰度（1%）高於ME16C（0.001%）；且整個XMV43基因組在MCF-7中的覆蓋率大於10X，而只有40%的XMV43基因組在ME16C中的覆蓋率大於10X（圖5.A）。在MCF-7細胞株中，XMV43的序列同源性更高，且所有基因組覆蓋率大於10X，所以TCGA樣本中的XMV43污染更可能來源於MCF-7。
pBABE-puro hTERT載體是用來轉導ME16C細胞株的，其中包含MMLV序列，它沒有與人類基因組顯著類似的序列和低複雜度區域，所以沒有被VirDetect屏蔽。做者發現ME16C轉錄本覆蓋了普羅黴素抗性基因和SV40啓動子（兩者均存在於pBABE-puro hTERT載體中）。正是該載體在ME16C細胞株轉導中的使用，解釋了在ME16C的RNA-seq數據中檢測到XMV43的緣由（圖5.B）
基於以上研究，做者給出了TCGA的RNA-seq樣本受到污染的整個流程：一、XMV43病毒感染MCF-7細胞株；二、將MCF-7和ME16C細胞株的RNA添加到包含Hela細胞的UHRR中；三、UHRR與TCGA同時測序，致使了TCGA樣本的污染。
圖5.比較XMV43序列在兩細胞株中的覆蓋率
四、補充一種污染物：狂犬病病毒
在10種腫瘤的共計19個樣本中檢測到了狂犬病病毒，而且集中發生在2012.11-2013.3這段時間內（圖6）。儘管病毒的表達量極低，可是在臨近的一段時間內被檢測到說明這一現象多是由該時間段的試劑污染引發的。
圖6.狂犬病病毒在樣本中的發現時間
小結
本篇文章中做者揭示了TCGA的RNA-seq數據庫污染來源，爲避免後續科研中出現數據偏差提供了幫助。在本項研究中，做者開發了能夠在RNA-seq數據中檢測病毒序列的開源軟件-Virdetect，利用該軟件檢測RNA-seq數據中的病毒表達，證明了TCGA中污染的HPV-18來自Hela細胞，XMV43主要來自MCF-7細胞株，並由此獲得了RNA-seq數據污染的流程：一、XMV43病毒感染MCF-7細胞株；二、將MCF-7和ME16C細胞株的RNA添加到包含Hela細胞的UHRR中；三、UHRR與TCGA同時測序，致使了TCGA樣本的污染。
點擊閱讀原文，便可得到小編本次與你們分享的文章。本期的分享就到這裏啦，下一期再見~
科研菌學術討論羣，在羣內能夠用本身的暱稱，廣告一概踢；其餘公衆號的宣傳也不發，就算是要發，提早和小編商量和確認，否則也是一概踢哈。
歡迎添加小編微信↑↑↑
請你們加個人時候就備註好「學術討論羣」以及本身的「單位+專業+姓名」
編輯：麻辣燙
校審：糯米飯
本文分享自微信公衆號 - 科研菌（keyanjun2020）。
若有侵權，請聯繫 support@oschina.cn 刪除。
本文參與「OSC源創計劃」，歡迎正在閱讀的你也加入，一塊兒分享。