DeepVariant: 用卷積神經網路進行DNA序列變異位點檢測

16年12月Google旗下的子公司Verily發了一篇文章描述了一個針對全基因組測序變異位點（SNP和small indel）檢測的新演算法，這個演算法不同於一般基於統計方法的軟體，而是利用了卷積神經網路識別變異位點。應該是目前第一個把神經網路運用到生物序列變異位點檢測中的，而且效果也十分理想，超過了現有演算法，並有很靈活，能夠運用到不同生物，測序平台中。

卷積神經網路被廣泛用於圖像識別任務中，一般會給一個訓練集，比如很多包含A圖片，和很多包含B的圖片，並告訴機器哪張是A，哪張是B，通過神經網路訓練後它可以對新的圖片進行分類。卷積粗略理解的話就是說在一張圖片中，可以分割出很多小塊，這些小塊上可能會有許多共同特徵，另一種理解方式是圖片中的物體對放大縮小不敏感，應用於局部的特徵也可以應用於全局。

DeepVariant的做法就是利用了現有的圖像識別卷積神經網路，具體用的神經網路架構inception v2。它利用一種轉換方法，把測序數據與reference的比對結果轉變成RGB圖片（只取存在變異位點的位置），然後用現有的圖像識別網路來分析其中的變異位點，大致轉化規則：

紅色通道：代表ATCG四種鹼基

綠色通道：序列質量 read quality

藍色通道：正向還是反向鏈 strand

在一張圖片中的第一行的像素代表reference序列，餘下行代表各read。對每一個候選位點都生成一張這樣的圖片。

---------------------------1/19/2018 更新----------------------------

最近開源的版本中轉換方法有更新（感謝ChenX同學提醒），大致為：

ATCG四種鹼基（base）

鹼基質量（base quality）

序列比對質量（mapping quality）

正反向鏈（strand）

是否支持alternative allele

是否支持reference allele

CIGAR序列中的信息（match，mismatch， insertion，deletion等）

-----------------------------------------------------------------------

然後這些候選位點被分為三類：純合子並與參考序列相同，雜合子，純合子並與參考序列不同，利用卷積神經網路進行訓練和分類。

文中使用了Genome in a Bottle reference data 作為驗證 （這個數據集中的樣本經過13種平台測序，所以裡面的變異位點高度可靠），取染色體1-18的數據作為訓練集，20-22作為測試集。並與另一個比較成熟的pipeline GATK進行對比，DeepVariant勝出。

之後文章又探索了一下DeepVariant的拓展性：

訓練好的模型直接應用在另一個reference上，效果很好。

訓練好的模型直接用在小鼠數據上，效果也很好，（也有對小鼠重新訓練的版本，差不多）

訓練好的模型對不同平台數據重新訓練後使用，效果也很好。

看起來是個很有希望的演算法，可惜目前還沒有完全開放給公眾使用，理論上Inception模型可以找到，根據文中的方法轉換圖片的話也可以完成，不過DeepVariant訓練了80個小時，很多訓練細節可能不太好掌握。

整體上這個想法還是非常有意思的，把一個不是圖像識別的問題轉化了一下，不過似乎目前這個轉化方法其實能用上的信息還比較有限，只考慮了read quality和strand，各通道利用率也比較低，比如藍色strand只有兩個可能的值，其實還有不少其它信息也能用。另外可能這個演算法的另一優勢是能夠考慮比較複雜的候選位點周邊序列的情況，圖片轉換中實際上附帶了這些信息。在其它演算法中有時也會考慮到附近序列的情況，但多為rule based，能考慮的情況就比較有限，而CNN在處理這些信息時應該比較有優勢。這個思路在生物問題中的應用應該還有很多潛力可挖，不知下一次會被應用在哪裡呢？

DeepVariant文章鏈接：http://biorxiv.org/content/biorxiv/early/2016/12/14/092890.full.pdf

----------------------------1/19/2018 更新：使用----------------------------

目前想用deep variant的話(google/deepvariant):

1.直接使用google cloud

這個應該是最簡單的了，google cloud本身有genomic部分，包含比如GATK的序列預處理，可以直接跑現成的流程。需要下載工具gsutil進行google cloud的操作，但是上傳大量數據gsutil速度有限，目前我不太清楚礎是否有更高效的傳輸方式。

費用方面新用戶有300刀信用點，最便宜一個30x全基因組測序需3-4小時，8-9刀，全外顯子組～70分鐘，0.7刀（https://cloud.google.com/genomics/deepvariant）。速度上如果僅用cpu，比GATK還慢一倍，不過GPU一加速就好說（https://blog.dnanexus.com/2017-12-05-evaluating-deepvariant-googles-machine-learning-variant-caller/）。

2. 通過Docker image

這個也比較方便，Docker相當於把軟體和這個軟體需要的環境打包好直接給你，這樣就不需要自己安裝配置各種dependency，太美好了。安裝好docker，拉一個deepvariant image直接跑，

a)下載docker，準備輸入文件，可以是自己的，也可以下載deepvariant提供的測試數據

b)拉一個docker image：

docker pull gcr.io/deepvariant-docker/deepvariant

c)使用image，-it是交互模式，-v是mount輸入文件和模型的路徑，input就是比對好的bam文件和參考序列，模型是圖片分類模型，可從google could上下載

docker run -it -v path/to/input:/dv2/input -v path/to/models:/dv2/models

d)運行(以下例子里的文件是deepvariant的測試數據）

./opt/deepvariant/bin/make_examples

--mode calling

--ref /dv2/input/ucsc.hg19.chr20.unittest.fasta.gz

--reads /dv2/input/NA12878_S1.chr20.10_10p1mb.bam

--examples output.examples.tfrecord

--regions "chr20:10,000,000-10,010,000"

./opt/deepvariant/bin/call_variants

--outfile call_variants_output.tfrecord

--examples output.examples.tfrecord

--checkpoint /dv2/models/model.ckpt

./opt/deepvariant/bin/postprocess_variants

--ref /dv2/input/ucsc.hg19.chr20.unittest.fasta.gz

--infile call_variants_output.tfrecord

--outfile output.vcf

3.conda package，有人準備了一個conda版可以嘗試（Improving pre-built DeepVariant binaries for conda packages · Issue #29 · google/deepvariant）。裝這個的話運行可以用提供的wrapper（dv_make_examples.py，dv_call_variants.py， dv_postprocess_variants.py），選項還不是很多。我自己在server上試了下還沒搞定，各種dependency問題，淚...

4.源代碼編譯安裝：大神們自便。。。

-----------------------------------------------------------------------