雋立然　王亞東

摘 要：面向多基因組的研究，以建模多個體關(guān)系和比較個體差異為主要研究內(nèi)容。多基因組可視化可以幫助研究者依據(jù)多個體關(guān)系，有目的地分析、比較多基因組之間的差異。多個基因組遺傳變異層面的比較，因為變異數(shù)量巨大、并且絕大部分變異并無信息性，故而很難在有限的顯示空間內(nèi)可視化。本文根據(jù)多基因組可視化的需求，分析了多基因組可視化的數(shù)據(jù)降維策略，提出了基于LDA模型及KL散度的多基因組相似度求解方法，建立了基于MDS算法的多基因組可視化降維方法，并使用千人基因組第三階段的基因組變異數(shù)據(jù)，驗證上述方法的可靠性。

關(guān)鍵詞：可視化；多基因組；LDA模型；MDS算法

中圖分類號：TP18 文獻(xiàn)標(biāo)識號：A 文章編號：2095-2163（2015）01-

Abstract： Multiple genome research focuses on modelling multiple individual relationships and comparing individual differences. Visualization of multiple genomes can help researchers to analyze and compare differences among multiple individual genomes intuitively. However， multiple genomes are difficult to be visualized in limited displaying space owing to the huge number of variants. And most of genome variants are less informative. This paper analyzed the data dimension reducing strategy for multiple genome visualization， proposes multiple genome similarity solving algorithm based on LDA model and KL-divergence， and designs multiple genome visualization dimention reducing method. The phase3 datasets of 1000 genomes project are used to verify the effectiveness and reliability of the above methods.

Keywords： Visualization； Multiple Genomes； LDA model； MDS Algorithm

0 引言

近年來隨著個體基因組測序的普及，測序的個體基因組數(shù)量大大增加，基于多個個體基因組數(shù)據(jù)的研究也日漸增多。遺傳學(xué)、人類學(xué)、社會學(xué)等許多學(xué)科都高度重視對人群的研究，早期的HapMap計劃[1]就有意識地搜集世界各地不同人群的基因組數(shù)據(jù)，作為其延續(xù)，2014年年中，千人基因組計劃公開發(fā)布了第三階段的數(shù)據(jù)，共包括屬于5個超級人群、26個人群的2 504個個體的基因組變異數(shù)據(jù)。在遺傳疾病的研究中，對多個疾病樣本與多個正常樣本的基因組進(jìn)行比照分析、對多個疾病亞型的樣本基因組進(jìn)行測試分析，均為常見的科學(xué)手段。因此，同時對多個個體的基因組進(jìn)行比較、分析即已成為生命科學(xué)和醫(yī)學(xué)研究中的重要需求。多基因組可視化能夠顯著提升多個個體基因組的比較和分析效率，也是重要的研究課題。

多基因組可視化并非多個個人基因組可視化的簡單集成，特別是當(dāng)需要可視化的個人基因組數(shù)量較多時，簡單集成的方法無法直觀地表達(dá)多個基因組之間的異同。多基因組可視化關(guān)注的是多個基因組之間的關(guān)系，也并非是基因組的一般性特征，這又不同于一般性的基因組可視化。多個基因組遺傳變異層面的比較，因為變異數(shù)量巨大、并且絕大部分變異并無信息性，故而很難在有限的顯示空間內(nèi)可視化，也即使研究者很難從大量變異數(shù)據(jù)中篩選出重要的變異。通過幫助研究者們查看多個基因組在遺傳變異層面的比較結(jié)果，并且尋找多個基因組中對研究有用的變異，則是多基因組可視化的主要目的。

本文根據(jù)多基因組可視化的需求，探討了多基因組可視化面臨的主要問題，分析了多基因組可視化的數(shù)據(jù)降維策略；提出了基于LDA模型及KL散度的多基因組相似度求解方法，其中，LDA模型由于可以給出相似基因組之間的共同潛在特征相關(guān)的變異列表和概率分布，將更加有利于對研究者所關(guān)心的變異進(jìn)行識別；并且建立了基于MDS算法的多基因組可視化降維方法；最后，本文使用千人基因組第三階段的基因組變異數(shù)據(jù)，分析和測試了上述方法的有效性。

1 LDA模型的基本理論

Latent Dirichlet Allocation（LDA）模型[3]是無監(jiān)督學(xué)習(xí)的概率主題模型，該模型假設(shè)每個文檔有不同概率的多個主題，而文檔中的詞匯則通過這些主題以一定概率具體選擇而生成。因此，通過學(xué)習(xí)語料庫中的文檔，LDA模型可以獲取語料庫中潛在的主題，并得到每個文檔的主題混合分布，以及每個主題中的主題-詞匯概率。由于是無監(jiān)督學(xué)習(xí)算法，使得LDA并不需要輸入標(biāo)注后的語料庫，并且對于每一個主題，都可以找出一個詞的集合對其進(jìn)行描述。LDA模型基于bag-of-words假設(shè)，即不考慮詞在文檔中的順序，并且也不考慮文檔之間的順序。但LDA并不假設(shè)詞匯或文檔是獨立同分布的。LDA模型可以用于文本主題識別，文本分類和文本相似度計算等問題。

LDA模型的主要任務(wù)是尋找使語料庫文檔具有較高概率的模型，并使語料庫以外的其他類似文檔也具有較高概率，以識別非語料庫文檔的主題。因此，LDA模型的基本策略是用一組隨機(jī)混合的隱含主題分布表示文檔，并使用詞匯的概率分布來描述每個主題。通過觀測到的語料庫內(nèi)文檔中詞匯的分布估計文檔-主題向量和主題-詞匯矩陣等參數(shù)，從而獲得主題的詞匯描述、文本的主題相似度等信息。

LDA模型最早通過變分貝葉斯期望最大化算法（Variational Bayesian Expectation Maximization，VBEM）估計參數(shù)[3]，也可以使用較快的吉布斯采樣（Gibbs Sampling）方法估計參數(shù)[4]。在LDA模型基礎(chǔ)上，D.Blei等和D. Ramage等隨后又提出了有監(jiān)督的LDA模型sLDA[5]和L-LDA[6]。

2基于LDA模型的多基因組相似度計算方法

基于此，一般化地考慮LDA模型，LDA模型試圖從可觀測到的離散數(shù)據(jù)單元和離散數(shù)據(jù)單元的無序集合的關(guān)聯(lián)關(guān)系中，為這些集合學(xué)習(xí)到一個有意義的隱含屬性，該屬性是集合包含其內(nèi)容數(shù)據(jù)的標(biāo)志，也是集合與集合之間進(jìn)行語義性比較的基礎(chǔ)，并且該屬性還可以使用一些離散數(shù)據(jù)單元所描述或定義。但是LDA模型并不限定該屬性和數(shù)據(jù)之間或?qū)傩院图现g有因果關(guān)系。

綜上理解可知，LDA模型可以應(yīng)用于多基因組的相似度計算和比較研究中。人類基因組之間有高達(dá)99.9%的相似性，個體基因組一般被表示為相對于一個標(biāo)準(zhǔn)的參考基因組的一組變異信息。因此基因組可以表示為許多變異的集合，對于基因組而言，這些變異只有分子位置上的順序關(guān)系，語義關(guān)聯(lián)上的順序關(guān)系可以被忽略。而根據(jù)不同的問題背景，該基因組可能具有不同的屬性，如不同人群、超級人群，或者疾病-正常、疾病的不同亞型等等。本文以屬于不同人群的多個基因組為例，應(yīng)用LDA模型計算多基因組的相似程度，本例中，個體的基因組事實上可能是多個人群的混合，如混血。圖1顯示了多基因組相似度計算問題的數(shù)據(jù)與LDA模型術(shù)語間的映射關(guān)系。

本文根據(jù)先驗知識的變異預(yù)篩選策略能夠使多基因組相似度計算問題的規(guī)模降低到LDA模型的求解算法能夠求解的范圍內(nèi)，并更好地識別有意義的人群-變異關(guān)系和基因組相似特征變異。

與一般的LDA模型解類似，人群（多基因組的子類別）的相似度以及多個個體基因組之間的相似度可以使用Kullback-Leibler Divergence，即KL散度（KL距離）[7]來刻畫，由于KL散度的不對稱性，也可以使用對稱KL散度，即KL散度的算術(shù)平均數(shù)、幾何平均數(shù)、調(diào)和平均數(shù)，或者JS散度（Jensen-Shannon divergence）及其平方根[8， 9]，本文將使用JS散度的平方根作為兩個個體基因組之間的相似性度量。

根據(jù)LDA模型的基本理論，使用變分貝葉斯期望最大化（VBEM）算法，可以迭代求解人群分布向量。

VBEM算法引入變分參數(shù)γ和?，簡化了原來由于θ、z和w的條件關(guān)系而難以求解的概率圖模型。指定了簡化的可優(yōu)化下界的函數(shù)后，即需尋找使下界函數(shù)和真實聯(lián)合后驗分布的KL散度極小化的變分參數(shù)γ和?，具體公式為：

3 基于MDS算法的多基因組可視化降維算法

MDS方法的基本流程為：

（1）給定M個樣本的K維數(shù)據(jù)，計算每對樣本之間的相似度/距離，并存入M×M的矩陣Δ。

（2）把數(shù)據(jù)投射到低維（r維，r << K）空間，為樣本在r維空間中隨機(jī)初始位置，使用一個M×r 的矩陣X存放投影后每個樣本在r維空間的坐標(biāo)。

（3）根據(jù)樣本在低維空間的坐標(biāo)，計算每對樣本之間的距離，一般為歐氏距離，并存入M×M的矩陣D。

（4）測量Δ與D的差別，差別使用應(yīng)力值衡量，計算公式如：

（5）如果應(yīng)力值大于閾值，即低維空間的樣本距離關(guān)系還沒有足夠近似高維空間中的樣本相似關(guān)系，則移動矩陣X中的樣本坐標(biāo)，使高維空間相似度高的樣本之間的距離減小，以減小總體應(yīng)力值。

（6）重復(fù)（3）-（5），如果應(yīng)力值小于閾值，或多次循環(huán)閾值差別不大（收斂到局部最優(yōu)），算法停止。

本文使用LDA模型求解了多基因組之間的相似性，相似度量是KL散度。KL散度是不對稱的，但是在MDS算法中，作為輸入的高維空間樣本距離應(yīng)是對稱的。因此在實踐中，常常使用KL散度的算數(shù)平均數(shù)、幾何平均數(shù)、調(diào)和平均數(shù)，或者JS散度及其平方根作為相似性的度量。

JS散度（Jensen-Shannon Divergence）是對稱和平滑版本的KL散度，數(shù)學(xué)定義如下：

MDS算法的另一個關(guān)鍵問題是應(yīng)力函數(shù)的優(yōu)化方法。本文采用SMACOF算法[10]最小化應(yīng)力函數(shù)，應(yīng)力函數(shù)定義為：

（5）

也就是說，最小化應(yīng)力函數(shù)實際是盡量令。具體地，dij為r維空間上樣本i和j的歐氏距離。

4 實驗結(jié)果與結(jié)論分析

研究采用千人基因組第三階段數(shù)據(jù)作為本文方法的測試數(shù)據(jù)。

在遺傳過程中，子代個體將繼承兩個親代個體的部分變異，并產(chǎn)生少量（<100）新的變異。因此，對于子代個體來說，性狀和表型主要由親代個體遺傳信息的重新組合決定，這是子代與親代相似性的遺傳基礎(chǔ)，同時也將使親代的遺傳特征以變異為表現(xiàn)形式而保留在后代的基因組中。由于地理因素和社會因素，人類在漫長的進(jìn)化和發(fā)展過程中，總是在一定的人群（population）范圍內(nèi)通婚，這就使得每個人群中廣泛存在某些從祖先獲得的較穩(wěn)定遺傳特征，而不見或少見于其他人群。典型的特征如膚色、瞳孔、發(fā)色等等。在基因組層面，這些遺傳特征可以用一個變異或多個變異的組合進(jìn)行描述，而且這些變異在不同人群中，則呈現(xiàn)為高內(nèi)聚、低耦合的特點。

圖2就是所有個體的基因組相似度計算和可視化結(jié)果。從圖2（a）中可以看到，本文的多基因組可視化方法，盡管采用了無監(jiān)督的算法，但卻完美重現(xiàn)了全部5個超級人群的劃分：歐洲人群（European，EUR），東亞人群（East Asian，EAS），混血美洲人群（Ad Mixed American，AMR），南亞人群（South Asian，SAS）和非洲人群（African，AFR）。

圖2（a）還準(zhǔn)確地定位了混血美洲人群的位置，即歐洲人群和非洲人群之間，但也延伸向東亞人群和非洲人群之間。這與中南美洲是歐洲殖民者、非裔奴隸和當(dāng)?shù)卦∶耖L期混血的事實十分吻合。特別地，其中延伸向非洲人群和東亞人群中間部分的趨勢，也與中南美洲原住民是冰川期從歐亞大陸沿白令海峽遷往美洲大陸的理論自洽。而南亞人群與歐洲人群同屬白色人種，兩者的距離比非洲人群和東亞人群更加接近。

由于各超級人群的區(qū)別十分明顯，為了進(jìn)一步討論超級人群內(nèi)部各人群的相似與區(qū)別，圖2（b）用不同形狀代表各超級人群，用顏色區(qū)分超級人群內(nèi)部的各個人群。超級人群下的各人群之間劃分也十分顯著。圖2（b）還顯示了在歐洲人群中，伊比利亞人群（IBS，圖b-1）最接近混血美洲人群中的波多黎各人群和哥倫比亞人群（PUR，CLM，圖b-2），這與中南美洲的主要殖民者是西班牙人和葡萄牙人也十分吻合。而南美洲太平洋沿岸國家秘魯?shù)娜巳海≒EL，圖b-3）則保留了較多的原住民血統(tǒng)。盡管混血美洲人群血統(tǒng)較為復(fù)雜，但本文的可視化方法仍然較為清晰地展示了混血美洲人群內(nèi)部各個人群的區(qū)別，以及這些人群與其他超級人群的聯(lián)系。

4 結(jié)束語

通過對圖2的分析，作為無監(jiān)督方法，本文的多基因組相似度計算和可視化降維算法完好地還原了千人基因組各超級人群的分布和內(nèi)部各人群的區(qū)分?？梢暬故镜母魅巳哼h(yuǎn)近關(guān)系符合生物學(xué)和歷史、地理上的各種證據(jù)，由此進(jìn)一步證實了本文多基因組可視化方法的能力和準(zhǔn)確性。

參考文獻(xiàn)：

[1] INTERNATIONAL H C. The International HapMap Project [J]. Nature， 2003， 426（6968）： 789-796.

[2] GENOMES P C， ABECASIS G R， AUTON A， et al. An integrated map of genetic variation from 1，092 human genomes [J]. Nature， 2012， 491（7422）： 56-65.

[3] BLEI D M， NG A Y， JORDAN M I. Latent Dirichlet allocation [J]. J Mach Learn Res， 2003， 3（4-5）： 993-1022.

[4] GRIFFITHS T L， STEYVERS M. Finding scientific topics [J]. Proc Natl Acad Sci U S A， 2004， 101 Suppl 1（5228-35.

[5] MCAULIFFE J D， BLEI D M. Supervised topic models； proceedings of the Advances in neural information processing systems， F， 2008 [C]// NIPS in Whistler， British Columbia， Canada，2007.

[6] RAMAGE D， HALL D， NALLAPATI R， et al. Labeled LDA： A supervised topic model for credit attribution in multi-labeled corpora； proceedings of the Proceedings of the 2009 Conference on Empirical Methods in Natural Language Processing： Volume 1-Volume 1， F， 2009 [C]// Association for Computational Linguistics. EMNLP in Singapore， 2009.

[7] KULLBACK S， LEIBLER R A. On information and sufficiency [J]. The Annals of Mathematical Statistics， 1951：79-86.

[8] JOHNSON D H， SINANOVIC S. Symmetrizing the kullback-leibler distance [J]. IEEE Transactions on Information Theory， 2001， 1（1）： 1-10.

[9] ENDRES D M， SCHINDELIN J E. A new metric for probability distributions [J]. Information Theory， IEEE Transactions on， 2003， 49（7）： 1858-60.

[10] De LEEUW J， MAIR P. Multidimensional scaling using majorization： SMACOF in R [J]. Department of Statistics， UCLA， 2011，