structure 2.3 - 中文_第1頁
structure 2.3 - 中文_第2頁
structure 2.3 - 中文_第3頁
structure 2.3 - 中文_第4頁
structure 2.3 - 中文_第5頁
已閱讀5頁,還剩40頁未讀 繼續免費閱讀

下載本文檔

版權說明:本文檔由用戶提供并上傳,收益歸屬內容提供方,若內容存在侵權,請進行舉報或認領

文檔簡介

1、Structure 2.3使用手冊Jonathan K. PritchardaXiaoquan WenaDaniel Falushb 1 2 3a芝加哥大學人類遺傳學系b牛津大學統計學系軟件來自/structure.html 2010年2月2日1我們在Structure項目中的其他的同事有Peter Donnelly、Matthew Stephens和Melissa Hubisz。2開發這個程序的第一版時作者(JP、MS、PD)在牛津大學統計系。3關于Structure的討論和問題請發給在線的論壇上:structure-software

2、。 在郵遞問題之前請查對這個文檔并搜索以前的討論。1 引言程序Structure使用由不連鎖的標記組成的基因型數據實施基于模型的聚類方法來推斷群體結構。這種方法由普里查德(Pritchard)、斯蒂芬斯(Stephens)和唐納利(Donnelly)(2000a)在一篇文章中引入,由Falush、斯蒂芬斯(Stephens)和普里查德(Pritchard)(2003a,2007)在續篇中進行了擴展。我們的方法的應用包括證明群體結構的存在,鑒定不同的遺傳群體,把個體歸到群體,以及鑒定移居者和摻和的個體。簡言之,我們假定有K個群體(這里K可能是未知的)的一個模型,每個群體在每個位點上由一組等位基因

3、頻率來刻畫。樣本內的個體被(按照概率)分配到群體,或共同分配到兩個或更多個群體,如果它們的基因型表明它們是混和的。假定在群體內,位點處于哈迪-溫伯格平衡和連鎖平衡。不精確地講,個體被按達到這一點那樣的方法指定到群體。我們的模型不假定一個特別的突變過程,并且它可以應用于大多數通常使用的遺傳標記,包括微衛星(microsatellites)、SNP和RFLP。模型假定在亞群體內標記不處于連鎖不平衡(LD),因此我們不能處理極其靠近的標記。從2.0版開始,我們現在能夠處理弱連鎖的標記。雖然這里實現的計算方法是相當強有力的,但是為了保證明智的答案,在運行程序的過程中還是需要謹慎。例如,不可能從理論上確

4、定合適的運行長度(時間),這需要用戶自己做一些實驗。這份資料描述軟件的使用和解釋,并補充發表的文章,這些文章提供了對方法的更正式的描述和評價。1.1 概述軟件包Structure由幾個部分組成。程序的計算部分用C語言編寫。我們發布源碼和用于各種平臺(目前有蘋果機,Windows,Linux,Sun)的可執行文件。C可執行文件讀取用戶提供的一個數據文件。還有一個Java前端為用戶提供各種有幫助的工具,包括對輸出的簡單的處理。你也可以從命令行調用Structure而不是使用前端。這份資料包括關于怎樣格式化數據文件、怎樣選擇合適的模型、以及怎樣解釋結果的信息。它也有關于使用兩種界面(命令行和前端)的

5、細節以及各種用戶定義的參數的匯總。1.2 在2.3版中有哪些更新?2.3版(2009年4月發布)引入了新的模型用于改進數據集結構的推論,其中(1)數據對于通常的結構模型來說信息不夠,不足以提供準確的推論,但是(2)抽樣的地點與群體歸屬關系(population membership)相關。在這種情形下,通過明確利用抽樣地點信息,我們使結構得到改善,經常允許性能提高很多(Hubisz et al., 2009)。我們希望在下幾個月釋放更進一步的改進。 表1:實例數據文件。這里MARKERNAMES = 1, LABEL = 1, POPDATA = 1, NUMINDS = 7, NUMLOCI

6、 = 5, MISSING = -9, POPFLAG = 0, LOCDATA = 0, PHENOTYPE = 0, EXTRACOLS = 0。第2列顯示個體的地理取樣位置。我們也可以把數據存儲為每個個體一行(ONEROWPERIND = 1),在這種情況下第一行為“George 1 -9 -9 145 -9 66 64 0 0 92 94”。Loc_a Loc_b Loc_c Loc_d Loc_e 喬治1-9 14566092喬治1-9 -9 64094保拉110614268192保拉110614864094馬修2110145-9 092馬修2110148661-9 鮑勃210814

7、264194鮑勃2-9 142-9 094Anja 1112142-9 1-9 Anja 111414266194彼得1-9 145660-9 彼得1110145-9 1-9 卡斯坦2108145620-9 卡斯坦2110145641922 數據文件的格式基因型數據的格式顯示在表2中(表1顯示一個例子)。基本上,整個數據集被作為一個矩陣安排在單個文件里,其中個體的數據在行里,位點在列里。用戶能對格式做出若干選擇,大多數這些數據(除基因型外!)是可選擇的。對于一個二倍體生物,每個個體的數據可以是作為連續的2行被儲存,其中每個位點在一列,或者在一行中,其中每個位點在連續的兩列。除非你打算使用連鎖模

8、型(見下面),否則單個個體的等位基因的次序并不重要。預基因型(pre-genotype)數據列(見下面)對每個體記錄兩次。(更一般地,對于n倍體生物來說,每個個體的數據被儲存在n個連續的行中,除非ONEROWPERIND選項被使用。) 2.1 數據文件的組成部分:輸入文件的要素如下所列。如果給出,它們一定按以下順序,然而大多數是可選的并且可以被完全刪除。用戶必須指明哪些數據被給出,或者在前端里(front end),或者(當從命令行運行Structure時)在一個單獨的文件mainparams里。同時,用戶也要指定個體和位點的數目。2.2 行1. 標記名稱(可選擇;字符串) 文件的第一行可以包

9、含數據集里的每個標記的標識符的一個列表。這一行包含整數或字母的L個字符串,其中L是位點的數目。2. 隱性等位基因(僅用于有顯性的標記數據;整數)SNP或者微衛星數據一般將不包括這一行。但是如果選項RECESSIVEALLELES被設置為1,則程序要求有這一行來表明每個標記上哪個等位基因(如果有的話)是隱性的。關于更多的信息請參閱第4.1節。該選項用于象AFLP那樣的數據,以及用于多倍體的情形,其中基因型可能是含糊的。3. 標記之間的距離(可選擇;實數)文件里的下一行是一個標記之間距離的集合,供有連鎖的位點使用。這些應該是遺傳距離(例如,厘摩),或者是這種距離的一些替代,基于(例如)物理距離。如

10、果標記距離(粗略地)與重組率成正比,則距離的實際單位不是那么重要 。前端從數據估計一個合適的尺度,但是命令行版本的用戶必須在文件extraparams里設置LOG10RMIN、LOG10RMAX和LOG10RSTART。標記必須按照連鎖群中的圖譜次序排列。當連續的標記來自不同的連鎖群(例如,不同的染色體)時,這應該用數值-1注明。第一個標記也被賦值為-1。所有其他的距離都是非負的。這一行包含L個實數。4. 連鎖相信息(可選擇;僅用于二倍體數據;在范圍0, 1內的實數)。這只供連鎖模型使用。這是L個概率的一行,出現在每個個體的基因型數據之后。如果連鎖相是完全知道的,或者沒有連鎖相信息可用,則這些

11、行是不必要的。當有來自家系數據的部分連鎖相信息,或者當來自雄性的單倍體X染色體數據和二倍體常染色體數據被一起輸入時,它們可能是有用的。對于連鎖相信息有兩種可選擇的表示:(1)個體的兩行數據被假設為分別與父本的和母本的相對應。連鎖相行表明當前標記上的排序正確的概率(設置MARKOVPHASE = 0);(2)連鎖相行表明與以前的等位基因有關的一個等位基因的連鎖相是正確的概率(設置MARKOVPHASE = 1)。第一項應該填入0.5,以便把這行填寫到L項。例如下列數據輸入表示來自一個男性的信息,有5個連鎖相未知的常染色體微衛星位點,后面是3個X染色體位點,使用母本/父本相模型:102156165

12、101143105104101100148163101143 -9 -9 -90.5 0.5 0.5 0.5 0.5 1.0 1.0 1.0其中-9表示“缺失數據”,這里缺失是由第二X染色體缺乏造成的,0.5表明常染色體位點的連鎖相是未知的,1.0表明X染色體位點由母本遺傳的概率為1.0,因此其連鎖相是已知的。相同的信息可以用markovphase模型來描述。這樣的話輸入文件將讀為:102156165101143105104101 100148163101143-9-9-9 0.5 0.5 0.5 0.5 0.5 0.5 1.0 1.0這里,2 1.0 s 表明那個第1 和第2,其次和第3 個

13、X染色體位點彼此完全同相。注意站點以站點產量在這些2 模式下將不同。在第一例子中,Structure將輸出母親和父親染色體的任務可能發生的事件。在第2 個情況下,它將輸出在輸入文件里列舉的每等位基因的可能發生的事件。5. 個體/ 基因型數據(必需的)取樣的每一個個體的數據象在下面描述的那樣安排成一行或多行。2.3 個體/基因型數據個體數據的每一行包含下列要素。這些形成數據文件里的列。1. Label(標簽)(可選擇;字符串) 一串整數或者字母,用來指明樣本中的每個個體。2. PopData(可選擇;整數)一個整數,指明一個用戶定義的群體,從其中獲得個體(例如這些整數可以指明個體取樣的地理位置)

14、。在默認的模型中,這個信息不被聚類算法使用,但是能用來幫助組織輸出(例如,將來自相同的預定義群體的個體彼此緊挨著繪圖)。3. PopFlag(可選擇;0或者1)一個布爾標簽,表明使用學習樣本時是否使用PopData(見USEPOPINFO,在下面)。(注:布爾(Boolean)變量(標簽)是取值為TRUE或FALSE的變量,在這里分別用整數1(使用PopData)和0(不使用PopData)表示。) 4. LocData(可選擇;整數)一個整數,為每個個體指明一個用戶定義的取樣地點(或者其他特性,例如一個分享的表現型)。當LOCPRIOR模型被打開時,這個信息用來幫助聚類。如果你僅僅希望使用L

15、OCPRIOR模型的PopData,那么你可以省略LocData列,并設置LOCISPOP = 1(這告訴程序使用PopData來設置地點)。5. Phenotype(可選擇;整數) 一個整數,為每個個體指明一個所關心的表現型的值(表中的f(i))。(表現型信息實際上沒有用于Structure。這里用來允許與關聯作圖程序STRAT有一個平滑的接口。) 6. Extra Columns(可選擇;字符串) 用戶把被程序忽略的附加數據包括在輸入文件里可能是方便的。這些數據就在這里輸入,可以是由整數或字符組成的串。7. Genotype Data(必需的;整數) 一個給定位點上的每個等位基因應該由一個

16、獨特的整數來編碼(例如微衛星重復得分)。2.4 缺失的基因型數據缺失數據應該用沒在數據中的其他地方出現過的一個數字來標明(按照慣例經常使用-9)。這個數字也可以用于有單倍體和二倍體數據混合的地方(例如男性中的X和常染色體位點)。缺失數據值是與描述數據集特性的其它參數一起被設置的。2.5 格式化的錯誤。我們已經進行了相當仔細的錯誤檢查,以保證數據集的格式正確,并且程序將試圖提供一些關于存在的任何問題的性質的提示。前端要求在每行的結束回車,不允許在行內回車;Structure的命令行版本以與處理空格或制表符(Tab)同樣的方式處理回車。可能出現的一個問題是,在將數據導入Structure之前用來組

17、裝數據的編輯程序可能引入隱藏的格式化字符,經常在行的末尾,或者在文件的末尾。前端能自動除去大多數這些錯誤,但是當數據文件好像處于正確的格式時,這類問題可能對錯誤負責。如果你正在把數據導入到一個Unix系統,dos2unix功能可能對徹底清理這些錯誤有幫助。3 用戶的建模決策3.1 祖先模型個體的祖先有4個主要模型:(1) 非混合模型(個體離散地來自一個群體或者另一個群體);(2)混合模型(每個個體從K個群體中的每一個抽取他/她的基因組的一部分);(3)連鎖模型(象混合模型一樣,但是連鎖的位點更可能來自相同的群體);(4)有先驗信息的模型(允許Structure使用關于取樣地點的信息:或者幫助用

18、弱的數據進行的聚類,發現遷移者,或者預定義一些群體)。關于模型1、2 、4的詳情見Pritchard等(2000a)和Hubisz 等(2009),關于模型3的詳情見Falush等(2003a)。1. 非混合模型。每個體完全來自K個群體之一。輸出報告個體i來自群體k的后驗概率。每個群體的先驗概率是1 / K。這個模型適合于研究完全離散的群體,并且經常比混合模型在檢測微妙的結構方面更強有力。2. 混合模型。個體可能具有混合的祖先。這可以表述為個體i從群體k中的祖先那里繼承了他的/她的基因組的一部分。輸出記錄這些比例的后驗平均估計值。以祖先向量q(i)為條件,每個等位基因的起源是獨立的。我們推薦這

19、個模型作為大多數分析的起始點。這是處理真實群體的大多數復雜性的一個相當靈活的模型。混合是真實數據的一個普通特征,如果你使用非混合模型,你或許不會發現它。混合模型也能以一種自然的方式處理混合的區域(hybrid zones)。表2:數據文件的格式,為兩行的格式。大多數這些組成部分是可選的(欲了解詳細信息,參見正文)。Ml是標記l的標識符。rl表明哪個等位基因,如果有的話,在每個標記上是隱性的(僅針對顯性的基因型數據)。Di,i+1是標記i和i + 1之間的距離。ID(i)是個體i的標簽,g(i)是個體i的一個預先定義的群體索引(PopData);f(i)是一個被用來合并學習樣品的標簽(PopFl

20、ag);l(i)是個體i的取樣地點(LocData);f(i)可以儲存個體i的表現型;y1(i), ., yn(i)用于儲存額外的數據(這些數據會被程序忽略);(xli,1, xli,2)儲存個體i在位點l上的基因型。pi(l)是個體i中的標記l的連鎖相的信息。3. 連鎖模型。這實質上是將混合模型推廣,來處理“混合連鎖不平衡”,即,在最近混和的群體中的連鎖標記之間出現的相關性。Falush等(2003a)描述了該模型和更詳細的計算。基本的模型是,過去的t個世代,有一次混合事件,將K個群體混合了。如果你考慮單個染色體,它由一系列“塊(chunk)”組成,這些“塊”是從混合時的祖先那里作為離散的單

21、位遺傳來的。出現混合LD是因為連鎖的等位基因經常在相同的塊上,因此來自相同的祖先群體。塊的大小被假設為獨立的指數隨機變量,具有平均長度1/t(以摩爾根為單位)。在實踐中我們估計“重組率”r,所用的數據對應于從現在的塊切換到新的塊的比率。個體i里的每個塊以概率qk(i)獨立地來自群體k,其中qk(i)是那個個體的祖先來自群體k的比例。總起來,新模型保留了混合模型的主要要素,但是在單個塊上的全部等位基因必須來自相同的群體。新的MCMC算法結合了可能的塊大小和斷點。它對于每個體報告總的祖先,考慮連鎖,并且也能報告染色體的每一點兒的起源的可能性,如果用戶想要的話。當使用連鎖的位點來研究混合的群體時,這

22、個新模型表現得比原先的混合模型更好。它得到對祖先向量的更準確的估計,并且能從數據中抽出更多的信息。這對混合作圖應該是有用的。該模型不是用于處理非常緊密連鎖的標記之間的背景LD的。顯然,這個模型是大多數混合群體的復雜現實的大大的簡化。不過,混合的主要的效應是在連鎖的標記之間建立長遠的相關性,因此我們這里的目的是在一個相當簡單的模型中將那個特征包括進來。計算比混合模型的要慢一點,特別對于大的K和不知道連鎖相的數據。不過,它們對于數千個位點和個體以及多個群體來說還是切實可行的。如果有關于標記的相對位置的信息(通常是一張遺傳圖譜),則只能使用該模型。4. 使用先驗的群體信息。Structure的默認模

23、式只使用遺傳學的信息來了解群體結構。不過,經常有可以與聚類相關的附加信息(例如,取樣的個體的物理特性或者取樣的地理位置)。目前,Structure可以用3種方式使用這種信息: LOCPRIOR模型:利用取樣位置作為先驗信息來輔助聚類用于結構信號比較弱的數據集。有一些數據集,其中有真實的群體結構(例如,取樣位置之間的顯著的FST),但是信號太弱,標準的Structure模型不能發現。對于標記很少、個體很少或者非常弱的Structure,經常是這樣的情況。在這種情形下,為了提高性能,Hubisz等(2009)發展了新模型,利用地點信息來輔助聚類。對于這樣的數據集,其中結構的信號太弱以致使用標準的S

24、tructure模型不能被發現,新模型經常能提供群體結構和個體祖先的準確的推斷。簡言之,LOCPRIOR模型的基本原理如下。通常,Structure假定個體的所有部分都大約是先驗等可能的。因為可能的部分的數目非常巨大,對于Structure來說,需要信息非常豐富的數據來斷定個體的任何特定的部分被聚類到群具有強的統計支持。相反,LOCPRIOR模型認為實際上,來自相同的取樣位置的個體經常來自相同的群體。因此,建立LOCPRIOR模型以期望取樣的位置可能關于祖先是信息豐富的。如果數據表明位置是信息豐富的,那么LOCPRIOR模型允許Structure使用這種信息。Hubisz等(2009)發展了一

25、對LOCPRIOR模型:一種用于沒有混合的情況,一種用于有混合的情況。在兩種情況中,內在的模型(以及似然函數)與標準版本相同。關鍵的差別是允許structure使用地點信息來幫助聚類(即,通過修改先驗信息來得到與位置有關的更偏愛的聚類解決方案)。LOCPRIOR模型具有合乎需要的特性:(i)當不存在結構時,它們不傾向于發現結構;(ii)當個體的祖先與取樣位置不相關時,他們能夠忽視取樣的信息;(iii)當群體結構的信號非常強大時,舊模型和新模型基本上給出相同的答案。因此,我們建議在大多數數據數量非常有限的情形下使用新模型,特別是當標準的Structure模型不提供一個Structure的清晰信號

26、時。但是,因為現在已經積累了標準的Structure模型的很多經驗,我們建議對于信息非常豐富的數據集將基本模型作為默認(Hubisz 等等,2009)。為了運行LOCPRIOR模型,用戶必須首先為每個個體指定“取樣地點”,作為一個整數編碼。即,我們假定樣品是在一組分離的位置收集的,并且我們不使用關于地點的任何空間信息。(我們認識到,在一些研究中,每個個體可能在一個不同的地點收集,因此將個體塞進一套更小的分離的地點可能不是對數據的理想的代表。) “地點”也可以代表一個表現型、生態型(ecotype)或者民族團體(ethnic group)。地點被鍵入到輸入文件中,要么在PopData列(設置LO

27、CISPOP = 1)中,要么作為一個單獨的LocData列(參閱第2.3節)。為了使用LOCPRIOR 模型,你必須首先指定或者用混合模型用非混合的模型。如果你使用的是圖形用戶界面版本,則勾選“use sampling locations as prio”(用取樣位置作為先驗信息)框。如果你使用的是命令行版本,則設置LOCPRIOR = 1。(注意,LOCPRIOR與連鎖模型不兼容。) 我們迄今的經驗是當不存在結構時,LOCPRIOR模型不偏向于檢測到假的結構。你可以把相同的診斷用于是否有真的結構,當你沒使用LOCPRIOR時。另外查看r的值可能有幫助,它確定由位置攜帶的信息的數量。r的值接

28、近1,或者<1,表明位置是信息豐富的。r的更大的值表明或者沒有群體結構,或者結構不依賴位置。USEPOPINFO模型:使用取樣位置來對移居者或者雜交種進行檢驗供信息非常豐富的數據數據集使用。在一些數據集里,用戶可能發現預確定的組(例如取樣位置)幾乎正好與結構聚類相對應,除了少數似乎被錯誤歸類的個體以外。Pritchard等(2000a)提出了正式的Bayesian檢驗,用于評價是否在這個樣品內的任何個體是他們認為的群體的移民,或者具有新近的移民祖先。注意這個模型假定被預先規定的群體通常是正確的。它采用十分強大的數據來克服先驗的錯誤分類。在使用USEPOPINFO模型之前,你也應該在沒有群

29、體信息的情況下運行程序,以保證預確定的群體與遺傳學的信息粗略一致。為了使用這模型,把USEPOPINFO設置為1,并且選擇MIGRPRIOR的一個值(在Pritchard等(2000a)中它是v)。你可以在0.001到0.1的范圍內為v選擇一個值。每個個體的預確定的群體被設置在輸入數據文件中(見PopData)。用這種方式,在輸入文件里被分配到群體k的個體在Structure算法中將被分配到群k。因此,被預先規定的群體應該是在1和MAXPOPS (K)之間的整數。如果任何個體的PopData超出這個范圍,它們的q將按正常的方式被更新(即沒有先驗的群體信息,根據將被使用的模型,如果USEPOPI

30、NFO被關上的話。)USEPOPINFO模型:預先指定一些個體的起源的群體來幫助未知起源的個體的祖先估計。使用USEPOPINFO模型的第二個方法是定義“學習樣本”(learning samples),它被預定義為來自特定的群。然后用Structure來聚類剩下的個體。注意:在前端里,這個選項使用“Update allele frequencies using only individuals with POPFLAG=1”選項被打開,位于“Advanced Tab”標簽下。學習樣品是利用數據文件里的PopFlag列實現的。預先規定的群體被用于那些個體,它們的PopFlag = 1(并且它們的

31、PopData在(1.K)中)。對于PopFlag = 0的個體,PopData值被忽略。如果數據文件里沒有PopFlag列,那么當USEPOPINFO被開啟時,PopFlag被為全部個體設置為1。具有PopFlag = 0的或者PopData不在(1.K)中的個體的祖先,根據混合或者沒有混合的模型被更新,象由用戶指定的那樣。如上所述,如果有很少的個體沒有預先規定的群體,將a設置成一個明智的值來可能是有幫助的。USEPOPINFO的應用可能在幾個方面有幫助。例如,可能有一些個體的來源是已知的,我們的目標是對未知來源的另外的個體進行歸類。例如,我們可能從一群已知品種(編號為1 . . .K)的狗

32、中收集數據,然后使用Structure為未知的(也許是雜交種)起源的另外的狗估計祖先。通過預先設置群體數目,我們可以保證Structure聚類對應于預先確定的品種,這使輸出更可解釋,并且能改進推論的準確性。(當然,如果兩個預先確定的品種在遺傳上是相同的,那么未知起源的狗可能被推斷為具有混合的祖先。USEPOPINFO的另一種用途是用于這樣一種情況:用戶想要只使用個體的一個子集來更新等位基因頻率。通常,Structure分析使用全部可得到的個體來更新等位基因頻率估計值。但是有一些情況,在那里你可能想對于一些個體估計祖先,沒有那些個體會影響等位基因頻率的估計。例如你可以有學習樣品的一個標準的收集,

33、然后周期性地你想要為新的一批基因型化的個體估計祖先。使用默認的選項,個體的祖先估計(稍微)取決于它們所在的批次。通過使用PFROMPOPFLAGONLY,你可以保證等位基因頻率估計值只依賴于PopFlag = 1的那些樣品。在不同的情況下,Murgia等(2006)想要確定一套無性系的狗瘤的起源。那些瘤如此緊密有關以至于使用的缺省設置時瘤形成它們自己的一類。通過使用PFROMPOPFLAGONLY,Murgia等迫使瘤與其他canid聚類分在一組。意見:我們建議首先運行Structure的基本的版本,以便證實被預先規定的標簽確實的確符合實際的遺傳學群體。其次,當使用學習樣品時,通過設置比0大的

34、MIGRPRIOR來允許一些錯誤的分類可能是明智的。3.2 等位基因頻率模型對于等位基因頻率有兩個基本的模型。一個模型假定每個群體內的等位基因頻率是獨立的,從一個分布中抽取,這個分布由參數l指定。那是用于Pritchard等(2000a)種的原先的模型。通常我們設置l = 1;這是缺省設置。Falush等(2003a)實施了一個模型,具有相關的等位基因頻率。這個模型標明不同群體中的頻率很可能是相似的(或許由于遷移或者由于共有的祖先)。更詳細的資料如下。獨立的模型對于很多數據集表現不錯。粗略地說,這最先說我們期望在不同的群體中的等位基因頻率彼此不同。相關的頻率模型說它們實際上可能十分相似。對于親

35、緣關系近的群體,這經常改進聚類,但是可能增加過高估計的K的危險(如下)。如果一個群體與其他群體分歧較大,則當那個群體被除去時,相關的模型有時可以取得更好的推論。估計l: 固定l = 1對于大多數數據是一個好主意,但是在一些情況下,例如SNP數據,其中大多數次要的等位基因是稀少的,這時候較小的數值可能工作得更好。對于這個原因,你可以讓程序為你的數據估計l。你可能想要這樣做一次,或許對于K = 1來說,然后將l固定在被估計的值上,因為在試圖同時國際太多的假設參數(l, a, F)時對于非識別性(non-identifiability)好像有一些問題。相關的等位基因頻率模型: 如同Falush等(2

36、003a)描述的那樣,相關的頻率模型使用一個(多維的)矢量,PA,它記錄假設的“祖先”群體中的等位基因頻率。假定在我們的樣品中代表的K個群體每個都已經經歷過與這些祖先頻率的獨立的漂移,速率分別用參數F1, F2, F3, ., FK表示。除歸因于有點不同的模型的差別和估計的差別外,被估計的Fk值應該數量上類似于FST值。此外,對于具有許多混合的數據要準確地估計Fk很難。PA被假設為具有Dirichlet先驗,具有與上面的群體頻率使用的相同的形式:pAl· D(l1, l2, . . . , ), (1)對每個l獨立。然后,群體k中的頻率的先驗為, (2)對每個k和l獨立。在這個模型里

37、,F與遺傳學距離FST有密切的關系。按照FST的標準的參數化方法,每個群體中的期望頻率由總的平均頻率給出,當等位基因的總頻率為p時,跨越亞群體的頻率的方差為p(1 p)FST。這里的模型幾乎一樣,除了我們對模型稍微做了推廣以外,通過允許每個群體以一個不通的速率(Fk)漂離祖先群體,如同群體具有不同的大小時可能被期望的那樣。我們也試圖估計“祖先頻率”,而不是使用平均的頻率。我們將獨立的先驗(prior)放于Fk上,與平均數為0.01、標準差為0.05的分布成正比(但是有PrFk ³ 1 = 0)。先驗分布的參數可以由用戶修改。一些實驗表明,0.01的先驗平均值對應于非常低細分的水平,對

38、于獨立頻率模型的數據經常導致好的表現。在其他的問題中(其中群體之間的差別更加明顯),好像數據通常壓倒了這個Fk的先驗。3.3 程序要運行多長時間? 程序從一個隨機的配置啟動,從那里采取一系列步驟穿過參數空間,每個步驟(只)依賴于前一個步驟的參數值。這個程序在運行期間引起不同的點上的Markov鏈的狀態之間的相關性。希望是通過運轉模擬足夠久,相關性將可以被忽視。有兩個問題要擔心:(1) burnin長度:在收集數據使啟動配置的影響減到最小之前模擬要運行多久,(2)在burnin以得到準確的參數估計之后模擬要運行多久。要選擇合適的burnin長度,看看由這個程序打印的歸納統計量的值是真的有幫助的(

39、例如(a,F,在群體之間的分歧距離Di,j,以及似然),以便了解它們是否已經收斂。通常10000100000的burnin非常足夠了。要選擇適當的運行長度,你需要在每個K上做幾次運行,也許長度不同,看看你是否得到一致的答案。通常,利用10000100000步運行你能得到參數(P和Q)的好的估計,但是Pr(X|K)的準確的估計可能需要更長時間的運行。實際上,你的運行時間的長度可能決定于你的計算機速度和耐心。如果你正處理極其大的數據集,并且被運行時間阻止,你可以試著修剪運行的長度和標記/個體的數量,至少為探索的分析。前端提供了幾個主要參數的時間序列曲線。在burnin階段結束之前你應該看看這些曲線

40、,以便了解這些曲線是否看起來達到了平衡。如果在burnin階段結束時數值仍然在增加或者減少,你需要增加burnin長度。如果在整個運行期間(即,不只是在burnin期間)a的估計值變化非常大,你可以通過增大ALPHAPROPSD來得到對Pr(X|K)的更準確的估計,這改進了在那種形勢下的混合。(見在第5節的一個有關的問題)。4 缺失數據,無效的等位基因和顯性標記當不斷改進Q和P時,程序忽略缺失的基因型數據。當在一個特別的位點有漏缺數據的可能性與個體在那里有什么等位基因無關時,這種方法是正確的。當具有漏缺數據的個體的Q的估計不那么準確時,沒有特別的原因阻止這樣的個體參加分析,除非他們根本幾乎沒有

41、數據。當數據以系統的方式遺漏時,出現一個嚴重的問題,如同用無效等位基因那樣。這些不適合假設的模型,即使沒有群體結構,也能夠導致明顯的違背哈迪-溫伯格。人們不應該期望假設的模型對這類破壞是穩健的。但是如果無效的等位基因可能是一個重要的問題的話,則顯性標記模型(下面)可以被使用。在樣本中有多名家庭成員也會破壞模型假定。這有時會導致K的過高估計,特別對于相關的頻率模型(Falush等,2003a),但是當K固定時,這對將個體分配給群體的影響很小。4.1 顯性標記、無效等位基因和多倍體基因型對一些類型的遺傳學標記(例如AFLP)來說,區分全部基因型是不可能的。其它類型的標記可能導致模棱兩可的基因型,如

42、果由于附近序列的變化導致PCR產物不能擴增,一部分等位基因為“無效”。從2.2版開始,我們實現了一個模型,處理與顯性標記相關的基因型的模糊性。總之,我們假定在任何特定的位點可能有對全部其他等位基因(例如A)為隱性的單個的等位基因,而全部其他的標記是共顯性的。因此AB和BB將作為“表現型”B出現在未加工的基因型數據中,AC和CC將被記錄為C,而BC將被記錄為BC。當有模糊性時,模型在可能的基因型上求和。全部的細節在Falush等(2007)里給出。 為了執行這些計算,必須告訴算法每個位點上的哪個等位基因(如果有的話)是隱性的。這通過設置RECESSIVEALLELES=1來進行,并且在輸入文件頂

43、上包括一行單L整數,在標記名稱和圖譜距離的(可選的)行之間,表明在數據集里的L個位點的每個上面的隱性等位基因。如果一個給定的位點上的全部標記是共顯性的,那么那個位點上的隱性值必須被調整成MISSING(缺失的)數據值。相反,如果隱性等位基因從未在純合狀態被觀察到,但是你認為它可能存在(例如可能有無效的等位基因),那么就把隱性值設置成在那個位點沒被觀察到的等位基因(而不是MISSING!) . 編碼基因型數據:如果表現型是不含糊的,那么它被在Structure輸入文件里按照它本來的樣子編碼。如果它是含糊的,那么它被作為顯性等位基因的純合體編碼。例如,表現型A 被編碼為AA,B被編碼為BB,BC被

44、編碼為BC,等等。如果標記是其他方面為二倍體的一個個體中的單倍體(例如男性中的X染色體),那么第2個等位基因被象以前一樣編碼為MISSING(缺失)。當A是隱性的時,基因型AB、AC等等在輸入文件里是不合法的。當RECESSIVEALLELES被用來處理無效的等位基因時,看起來是無效的純合體(homozygote null)的基因型應該作為隱性等位基因的純合體而不是作為缺失數據被輸入。在實踐中可能不確定是否一個失敗的基因型真的歸因于純合的無效等位基因。Structure應該對這些編碼為缺失的數據是穩健的,除非無效等位基因在一個位點上的頻率很高。在多倍體(PLOIDY>2)中形勢更復雜,因

45、為甚至對共顯性標記都可能有基因型的含糊。在雜合體中準確地識別出基因型經常是困難的。例如在三倍體中,表現型AB可能是AAB或者ABB。如果Structure在RECESSIVEALLELES=0的條件下運行,那么就假定沒有含糊。對于多倍體,當RECESSIVEALLELES=1時,Structure允許數據包含具有基因型模糊和不具有基因型模糊的位點。如果一些位點不含糊那么設置代碼NOTAMBIGUOUS為一個整數,這個整數不與數據內的的任何等位基因相匹配,并且不等于MISSING(缺失)。然后在輸入文件頂上的隱性等位基因的行里為不含糊的位點放置NOTAMBIGUOUS代碼。如果不是那樣,而是在一

46、個特定的位點上等位基因全部是共顯性的,但是有關于每個的數目(例如為在四倍體里的微衛星)含糊,那么就把隱性等位基因代碼設置為MISSING。最后,如果有隱性等位基因,并且還有關于每個等位基因的數目的含糊性,則設置隱性等位基因代碼來表明哪個等位基因是隱性的。存在拷貝數含糊性的等位基因的編碼與存在顯性標記的那些相似。因此,舉例來說,在四倍體中,觀察到3個共顯性位點B、C和D,這應該被編碼為B C D D或者等效地B B C D或者任何包括3個等位基因中的每一個的其他組合。它不應該被編碼為B C D (MISSING),因為這表明該特定的個體在所指的位點是三倍體。如果在這個位點上存在一個隱性等位基因A

47、,它也不能被編碼為B C D A。Pr(K)的估計: 當RECESSIVEALLELES被用于二倍體時,Markov鏈的每個步驟上的似然值是通過在可能的基因型上求和來計算的。為了便于編碼,當要么PLOIDY>2要么使用了連鎖模型時,我們以當前推算的(imputed)基因型為條件。這減小似然值,并且好像大大地擴大似然值的方差。有限的經驗表明在后一種情況下這導致對K的估計效果變差,你應該把K的這種估計看做是不可靠的。5 K(群體數目)的估計在描述這個程序的我們的文章里,我們指出這個問題應該被小心對待,由于兩個原因:(1)要獲得對Pr(X|K)的準確估計在計算上是困難的,我們的方法僅僅提供一個

48、專門的(ad hoc)近似;(2)K的生物學解釋可能不是簡單的。在我們的經驗里我們發現真正的困難在于第2個問題。我們的用于估計K的程序一般在具有少量離散的群體的數據集中計算效果較好。不過,很多現實世界的數據集并不準確地符合Structure模型(例如,由于通過距離或者近交而產生的隔離)。在那些情況里對于什么是K的正確值可能沒有一個自然的答案。或許由于這種原因,在真實的數據中我們的模型選擇標準的值隨著增加的K而繼續增加是不稀有的。那么集中于捕獲數據中的大多數結構的K的值通常是講得通的,這在生物學上似乎是合理的。5.1估計K的步驟1. (命令行版本)在文件extraparams里把COMPUTEP

49、ROBS和INFERALPHA設置為1。(前端版本)確保a允許改變。2. 對不同的MAXPOPS (K)值運行MCMC方案。最后它將輸出一行“Estimated Ln Prob of Data”。這是ln Pr(X|K)的估計。你應該對每個K獨立地運行幾次,以便證實不同運行得到的估計值是一致的。如果與不同的K獲得的估計值的變異性相比,一個給定的K的不同運行的變異性是顯著的,那么你可能需要使用更長的運行或者更長的burnin時期。如果lnPr(X|K)看起來是雙峰的(bimodal)或者多峰的(multimodal),則MCMC方案可能找到不同的答案。你可以對此進行驗證,通過比較在單個K上的不同

50、運行的Q。(參看Pritchard et al. (2000a)的數據集2A(Data Set 2A),也見下面有關多峰性(Multimodality)的部分,)。3. 計算K的后驗概率。例如,對于論文中的數據集2A(這里K是2),我們得到K ln Pr(X|K)1 -43562 -39833 -39824 -39835 -4006我們一開始可以假定一個關于K = 1, ., 5的均勻先驗分布。然后根據貝葉斯定理,Pr(K = 2)由下式給出: (3)如果我們將該式簡化為下面的公式,計算就會更容易 (4)5.2 輕微的違背模型可能導致過高估計K 當存在真正的群體結構時,這導致不連鎖的位點之間的

51、LD,以及違背哈迪溫伯格比例。粗略地說,這是被Structure算法使用的信號。但是模型的一些違背也能導致哈迪溫伯格或連鎖不平衡。這些包括近交和基因型鑒定錯誤(例如偶然的、未被發現的無效的等位基因)。即使在沒有群體結構的情況下,對于K >1,這些類型的因素也可能導致弱的統計信號。從2版本開始,我們提出相關的等位基因頻率模型(correlated allele frequency model)應該被用作默認,因為它在困難的問題上經常實現更好的執行,但是用戶應該意識到,在這樣的設置中可能更容易過高估計K,與獨立的頻率模型下相比(Falush et al. (2003a))。 下一節討論怎樣確

52、定推斷的結構是否是真實的。5.3 關于選擇K的非正式提示;結構是真實的嗎?有兩個非正式的提示,可能有助于選擇K。第一個是,對于比合適的值(有效零)更小的K,Pr(K)常常是非常小的,對于更大的K,則有或多或少的高原,如同上面顯示的數據集2A的例子中那樣。在這種情形中(其中K的幾個值給出log Pr(X|K)的相似的估計下),似乎這些估計中最小的常常是正確的。對于我們通過“或多或少的高原”所表示的東西,要提供一個堅固的規則有點難。對于小數據集來說,這可能意味著log Pr(X|K)的值在5-10的范圍內,但是Daniel Falush寫道“在非常大的數據集中,K = 3和K = 4之間的差別可能

53、是50,但是如果K = 3和K = 2之間的差別是5 000,那么我將肯定選擇K = 3”。想要使用更正式的標準(這種標準將這一點納入了考慮)的讀者可能對Evanno等(2005)的方法感興趣。我們認為考慮這一點的一種明智的方法是就模型選擇而言。即,我們可能不總是能知道K的真值,但是我們應該致力于捕獲數據里的主要結構的K的最小的值。第二提示是,如果真的有單獨的群體,那個,通常有許多有關a的值的信息,一旦Markov 鏈收斂,a通常將相對恒定(范圍經常為0.2或更少)。不過,如果沒有任何真正的結構,在運行過程中a通常變化很大。這一點的一個必然的結果是當沒有群體結構時,你將通常將看到分配給每個群體

54、的樣本的比例是大致勻稱的(每個群體中1/K),大多數個體將被公平地混和。如果一些個體被強烈地分配到一群體或者另一個,以及如果分配給每組的比例不對稱,那么這是你有真正的群體結構的強的跡象。假定你有兩個清楚的群體,但是你試圖決定是否這些中之一是更進一步再分(例如,Pr(X|K = 3)的值類似于P(X|K = 2),或者也許比P(X|K = 2)還大一點)。那么,你能嘗試的一件事情是只使用你懷疑可能被再分的群體內的個體來運行Structure,看看是否有一個如上所述的強信號。總之,你應該對根據小的Pr(K)的差別推斷的群體結構持懷疑態度,如果(1)對于分派沒有清楚的生物學解釋,(2)對全部群體的分

55、派大致勻稱,沒有個體被強烈地分配。5.4 通過距離數據的隔離通過距離的隔離指的是這樣的想法:個體可能跨越一些地區呈空間分布,帶有本地分散的。在這種形勢下,等位基因頻率跨越地區逐漸變化。基礎的Structure模型對來自這種情況的數據不很適合。當這發生時,推斷的K 的值,以及在每組中的相應的等位基因頻率可能相當任意。取決于取樣的計劃,大多數個體可能在多個組中具有混合的成員身份。即,算法將嘗試使用K的不同組分的加權平均數來對跨越地區的等位基因頻率建模。在這樣的形勢下,結果的解釋可能具有挑戰性。6 背景LD和其他miscellania 6.1 序列數據,緊密連鎖的SNP和單體型數據Structure

56、模型假定位點在群體內是獨立的(即,在群體內不處于LD)。序列數據或者來自非重組區域的(比如Y染色體或者mtDNA)的數據很可能違反這個假定。如果你有序列數據或來自多個獨立區域的密集的SNP數據,那么盡管數據不完全適合模型,Structure實際上可能表演得想當好。粗略地說,這將發生,倘若跨越不同的區域有足夠的獨立性,以至于區域內的LD不在數據中占優勢。當有足夠的獨立區域時,區域內的依賴性(dependence)的主要代價將是Structure在特別的個體的分派中低估不確定性。例如,Falush等(2003b)把Structure用于來自H. pylori的MLST(多位點序列)數據,以了解H.

57、 pylori的群體結構和遷移歷史。在那種情況下,在區域內有足夠的重組以至于群體結構的信號超過了背景LD。(關于MLST數據的更多情況,也見第10節)。在人類的應用中,Conrad等(2006)發現來自36個連鎖的區域的3000個SNP生產明智(但是嘈雜)的答案,在一個全世界的樣本中,基本上與基于微衛星的以前的結果一致見他們的Supplementary Methods Figure SM2。然而,如果數據被一個或者少數非重組的或在低重組的區域主導,那么,Structure可能被嚴重地誤導。例如,如果數據只由Y染色體數據組成,那么估計的結構大概將反映出關于Y染色體樹的某些事情,而非群體結構本身。使用這樣的數據的影響很可能是:(1)算法低估祖先估計中的不確定性的程度,在最壞的情況下,可能是有偏的或者不準確的;(2)K的估計不可能表演得好。如果你有Y或者mtDNA數據加上許多核標記,一個安全和有效的解決辦法是重新編碼來自每個連鎖區域的單體型,以至于單體型被描述為一個具有n等位基因的單個位點。如果有許多單體型,則可以把相關的單體型歸類到一起。注意連鎖模型不一定比(非)混合模型對于處理這些問題更好。連鎖模型不是設計來處理群體內的背景LD的,并且很可能被類似地干擾。6.2 多峰性Structure算法在參數空間中的一個隨機的地方開始,

溫馨提示

  • 1. 本站所有資源如無特殊說明,都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
  • 2. 本站的文檔不包含任何第三方提供的附件圖紙等,如果需要附件,請聯系上傳者。文件的所有權益歸上傳用戶所有。
  • 3. 本站RAR壓縮包中若帶圖紙,網頁內容里面會有圖紙預覽,若沒有圖紙預覽就沒有圖紙。
  • 4. 未經權益所有人同意不得將文件中的內容挪作商業或盈利用途。
  • 5. 人人文庫網僅提供信息存儲空間,僅對用戶上傳內容的表現方式做保護處理,對用戶上傳分享的文檔內容本身不做任何修改或編輯,并不能對任何下載內容負責。
  • 6. 下載文件中如有侵權或不適當內容,請與我們聯系,我們立即糾正。
  • 7. 本站不保證下載資源的準確性、安全性和完整性, 同時也不承擔用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。

評論

0/150

提交評論