表1:利用隨機遊走模型從首次複製的DBP基因組關係矩陣中得到前10個SNPs和對應的原始p值
全文
Burak Karacaoren＊
土耳其安塔利亞Akdeniz大學農學院動物科學係*通訊作者:Burak Karacaören，土耳其安塔利亞Akdeniz大學農學院動物科學係，E-mail: burakkaracaoren@akdeniz.edu.tr
背景
譜係信息線性混合模型是橫斷麵基因組學研究中常用的遺傳關係檢測和校正方法。本研究的主要目的是利用隨機步行-卡爾曼濾波方法對GAW18數據集進行動態關聯映射分析。通過對動態基因和環境效應的自回歸結構對模型進行了擴展，將平穩過程納入模型。
方法
我們使用的是隨機遊走模型，如下圖所示
$${y_t} = {a_t} + {\varepsilon _t}，{\varepsilon _t} \propto N\左({O，\sigma _e^2} \右);{a_{t + 1}} = {a_t} + {\eta _t}，{\eta _t} \propto {\mkern 1mu} N\左({O,A\sigma _n^2} \右)$$在(1)中，第一個方程稱為觀測方程，第二個方程稱為狀態方程。我們假設觀察，yt，取決於不可觀測的量，αt，我們的目標是做統計推斷αt(州)
結果
用基因組關係矩陣預測舒張壓的誤差、遺傳和永久環境方差分量分別為17.3(0.0006)、10.9(0.0007)和8.2(0.0006);用譜係關係矩陣預測舒張壓的誤差、遺傳和永久環境方差分量分別為18.0(0.0007)、9.2(0.0008)和8.4(0.0006)。時間點1,2,3的Rs 11711953被發現與MAP4基因相關。
結論
可能由於時間點較少，該模型沒有檢測到所有真實的基因組信號。基因組關係矩陣給出了較好的膨脹因子。隨機遊走是一個非平穩過程，本文通過調整∆參數擴展了平穩情況下的模型。在基因組研究中，如果不考慮縱向基因和環境隨時間的影響，可能會導致無法檢測到真實信號和/或由於隨機誤差也可能導致假陽性。
協會的映射;隨機漫步;卡爾曼濾波器;吉布斯抽樣
結合譜係信息的混合模型方法通常用於橫斷麵基因組學研究中遺傳關係的檢測和校正。由於基因表達可能隨時間而改變，重複測量對檢測相關基因組信號[2]更有用。動態關聯研究通常采用功能映射方法。然而，解釋非參數函數回歸係數的結果在生物學上可能是困難的。動態關聯映射[3]采用隨機回歸係數模型。然而，隨機回歸模型在獲得軌跡開始和結束的精確估計[4]方麵有局限性。
此外，上述兩種模型(以及動態關聯映射文獻中的其他模型)都需要一整套的觀察結果才能進行預測:這可能也會產生等待數月(如果不是數年)才能獲得預測的問題。在本研究中，我們假設基因組信號隨時間的變化可以通過狀態空間形式的隨機遊走-卡爾曼濾波模型進行跟蹤，以獲得縱向殘差。由於卡爾曼濾波，我們不必等待收集整個數據集來進行模型評估，因此隻要進行測量，就可以進行估計。並且由於關聯映射中使用了縱向殘差:如果信號是真實的，生物推理也可以很容易地推導出來。
最近，我們在貝葉斯環境[5]中擴展了[1]的語法模型。本文采用[6]模型進行動態關聯映射，采用隨機步行-卡爾曼濾波方法對GAW18數據集進行動態關聯映射。通過對動態基因和環境效應的自回歸結構對模型進行了擴展，將平穩過程納入模型。
GAW 18在3個時間點對849個個體的200個重複提供模擬表型。我們使用舒張壓(DBP)表型進行關聯映射。我們在第一次複製的3個時間點使用相關的849個個體分析了來自3號染色體的65519個snp。
質量控製
我們使用了來自3號染色體的849個有65519個SNPs的係譜個體進行關聯映射。由於小等位基因頻率<1%，我們排除了7229個SNPs, Hardy Weinberg檢驗(p<0.001)排除了208個SNPs，遺漏檢驗(p>0.1)排除了2個SNPs，在分析[7]中留下58080個SNPs。我們排除了44個基因型過低的個體，隻剩下數據集中的805個個體。Kolmogrow-Smirnow檢驗用於評價反應變量的正態性。時間、性別、吸煙狀況、年齡和係譜數作為固定效應納入後續分析，其基礎是利用預測和觀察結果之間的相關性進行初步分析。
隨機遊走模型
我們使用的是隨機遊走模型，如下圖所示
$ $ {y_t} ={\α_t} + {\ varepsilon _t}, {\ varepsilon _t} \ propto N(0 \σ_e ^ 2) $ $
$ ${\α_ {t + 1}} ={\α_t} +{\埃塔_t},{\埃塔_t} \ propto N(0 \σ_n”^ 2 )............( 1) $ $
在(1)中，第一個方程稱為觀測方程，第二個方程稱為狀態方程。我們$ $ \ eqalign{&假定\ \ {\ rm {}}, {\ rm{}} \觀察,{y_t}, \, rm{}}{\ \,取決於\ {\ rm {}}, {\ rm{}},難以察覺的\ \ {\ rm{}}數量,{\α_t}, {\ rm {}} \,, \ \ cr & \ {\ rm {}}, {\ rm {}} {\ rm{}} \目標,是\ \ {\ rm {}}, {\ rm{}} \,統計\ {\ rm {}}, {\ rm{}},推理\ \ {\ rm{}}, \,α_t}{\ \({州}\右)。\,我們\ {\ rm {}}, {\ rm {}} \, {\ rm{}}常數\ \ cr & {\ rm{}}方差\ \ {\ rm {}}, {\ varepsilon _t} \,{\埃塔_t} \ \,, \ \σ_e ^ 2和\ \σ_n”^ 2 \ \,分別與\ {\ rm {}}, {\ rm{}}獨立\ cr} $ $相同和正態分布隨機變量與零的意思。我們假設基因影響和永久的環境影響。
性狀遺傳分析采用混合模型;
$ $ y = X \β+ {Z_a} + {Z_p} p + e ................( 2) $ $
其中，y為觀測向量，β為固定效應向量，a為隨機效應向量，p為隨機永久環境效應向量，X, Z一個, Zp為設計矩陣，e為隨機$$\eqalign{& residual\，{\rm{}}效應的向量。\ \,σ_a ^ 2 \ \,, \σ_p ^ 2 \ {\ rm {}} \, {\ rm{}}, \ \σ_e ^ 2 \; \, \ {\ rm {}}, {\ rm{}}遺傳、\,{\ rm{}}永久\環境\ {\ rm {}}, {\ rm {}} \, {\ rm{}}錯誤\ \ cr & {\ rm{}}差異。\，A\，是\，{\rm{}} \，{\rm{}}加性\cr} $$個體遺傳關係矩陣;是一個單位矩陣。結合個體的基因型和譜係，通過共祖矩陣的係數得到A。
下麵，基於貝葉斯原理，我們展示了KF-RW方法中使用的一般假設。(2)根據以下遞歸關係[8]給出(3)中無常數項的關節後驗比例分布;
$ $ p \左({{\θ_t} \左|{{\θ_{\離開({- t} \右)}},{最大}}\。} \) \ propto p \左({{\θ_t} \左|{{\θ_{\離開({t - 1} \右)}},{Y_ {t - 1}}} \。左({}\右)p \{\θ_ {t + 1}} \左|{{\θ_t}, {Y_ {t - 1}}} \。左({}\右)p \ {{\ rm {Y}} _t} \左|{{\θ_t}, {Y_ {t - 1}}} \。} \右),$ $
$ $ f \離開({y / b、p \σ_a ^ 2 \σ_p ^ 2 \σ_e ^ 2} \) \ propto{\離開({\σ_e ^ 2} \右)^ {- {1 \ / 2}N}} {\ mkern 1μ}\ exp \離開({-{1 \ / 2}{{\離開({y - Xb - {Z_a} - {Z_p} p} \右)}^ \ '}\離開({y - Xb - {Z_a} - {Z_p} p} \右)/ \σ_e ^ 2} \右)$ $
$ $ {\ rm {X}} \三角洲{\離開({\σ_a ^ 2} \右)^ -}^ {{1 \ / 2}{N_a}} \ exp \離開({{1 \ / 2}{{{rm \{‘}}}_1}{{\ rm{一}}^ {- 1}}{a_1} / \σ_a ^ 2}左\)\[{\刺激\ limits_ {t = 2} ^ t{{{\離開({\σ_a ^ 2} \右)}^ {- {1 \ / 2}{N_a}}}} \ exp \離開({-{1 \ / 2}{{\離開({{{rm \{一}}_ {\ rm {t}}}, {{\ rm{一}}_ {{\ rm {t}} - 1}}} \右)}^ \ '}{{\ rm{一}}^{- 1}}\離開({{{rm \{一}}_ {\ rm {t}}}, {{\ rm{一}}_ {{\ rm {t}} - 1}}} \右)/ \σ_a ^ 2} \右)}\右]$ $
$ $ {\ rm {X}} \三角洲{\離開({\σ_p ^ 2} \右)^ {- {1 \ / 2}{N_p}}} \ exp \離開({{1 \ / 2}{{{rm \ {p '}}} _1} {{rm \ p{}} _1} / \σ_e ^ 2}左\)\[{\刺激\ limits_ {t = 2} ^ t{{{\離開({\σ_p ^ 2} \右)}^ {- {1 \ / 2}{N_p}}}} \ exp \離開({-{1 \ / 2}{{\離開({{{rm \ p {}} _ {\ rm {t}}}, {{rm \ p {}} _ {{\ rm {t}} - 1}}} \右)}^ \ '}{{\ rm{我}}^{- 1}}\離開({{{rm \ p {}} _ {\ rm {t}}}, {{rm \ p {}} _ {{\ rm {t}} - 1}}} \右)/ \σ_p ^ 2} \右)}\右]$ $
$ $ {\ rm {X}}{\離開({\σ_e ^ 2} \右)^{- \離開({{{{\ν_e}} \ / 2} + 1} \右)}}\ exp \離開了 [ { - {{{\ ν_e} {S_e}} \ /{2 \σ_e ^ 2}}} \右]{\離開({\σ_a ^ 2} \右)^ -}^{\離開({{{{\ν_a}} \ / 2} + 1} \右)}\ exp \離開了 [ { - {{{\ ν_a} {S_a}} \ /{2 \σ_a ^ 2}}} \右]{\離開({\σ_p ^ 2} \右)^{- \離開({{{\νp} \ / 2} + 1} \右)}}\ exp \離開了 [ { - {{{\ ν_p} {S_p}} \ /{2 \σ_p ^ 2}}} \ ]...................( 3) $ $
(3)的最後一行是方差參數先驗假設的比例倒卡方分布密度的乘積。對於隨機漫步模型，假設∆為1。經過代數運算後，條件分布可以寫成:
$ $ {b_t} \左|{\σ_a ^ 2 \σ_p ^ 2 \σ_e ^ 2, {a_t}, {p_t}, {y_t} \ sim N \離開({{{\離開({X 'X} \右)}^ {- 1}}X ' \離開({y - {Z_a} - {Z_p} {p_t}} \右),{{\離開({X 'X} \右)}^{- 1}}\σ_e ^ 2} \右)}\ $ $
$ $ {a_t} \左|{\σ_a ^ 2 \σ_p ^ 2 \σ_e ^ 2, {b_t}, {p_t}, {y_t}} \。\ sim $ $
$ $ N \左({{{\左({{1 \ /{\σ_e ^ 2}} {{{rm \ {Z '}}} _ {{a_t}}} {{\ rm {Z}} _ {{a_t}}} +{2 \ /{\σ_a ^ 2}} {{\ rm{一}}^{- 1}}}\右)}^{- 1}}\離開({{1 \ /{\σ_e ^ 2}} {\ rm {Z '}} \離開({{{rm \ {y}} _t}, {{\ rm {X}} _ {\ rm{一}}}{{\ rm {b}} _t}, {{\ rm {Z}} _ {{{rm \ p {}} _ {\ rm {t}}}}} {{rm \ p {}} _ {\ rm {t}}}} \右)+{1 \ /{\σ_a ^ 2}} {{\ rm{一}}^{- 1}}{現代{t + 1}}} \右),{{\離開({{1 \ /{\σ_e ^ 2}} {{{rm \ {Z '}}} _ {{{rm \{一}}_ {\ rm {t}}}}} {{\ rm {Z}} _ {{{rm \{一}}_ {\ rm {t}}}}} +{2 \ /{\σrm _a ^ 2}}{{\{一}}^{- 1}}}\右)}^{- 1}}}\右)$ $
$ $ {p_t} \左|{\σ_a ^ 2 \σ_p ^ 2 \σ_e ^ 2, {a_t}, {b_t}, {y_t}} \。\ sim $ $
$ $ N \左({{{\左({{1 \ /{\σ_e ^ 2}} {{{rm \ {Z '}}} _ {{p_t}}} {{\ rm {Z}} _ {{p_t}}} +{2 \ /{\σ_p ^ 2}} {{\ rm{我}}^{- 1}}}\右)}^{- 1}}\離開({{1 \ /{\σ_e ^ 2}} {\ rm {Z '}} \離開({{{rm \ {y}} _t}, {{\ rm {X}} _ {\ rm{一}}}{{\ rm {b}} _t}, {{\ rm {Z}} _ {{{rm \ p {}} _ {\ rm {t}}}}} {{\ rm{一}}_ {\ rm {t}}}} \右)+{1 \ /{\σ_p ^ 2}} {{\ rm{我}}^{- 1}}{現代{t + 1}}} \右),{{\離開({{1 \ /{\σ_e ^ 2}} {{{rm \ {Z '}}} _ {{{rm \ p {}} _ {\ rm {t}}}}} {{\ rm {Z}} _ {{{rm \ p {}} _ {\ rm {t}}}}} +{2 \ /{\σrm _p ^ 2}}{{\{我}}^{- 1}}}\右)}^{- 1}}}\右)$ $
$ $ \σ_a ^ 2 \左|{\σ_p ^ 2 \σ_e ^ 2, {p_t}, {a_t},} \。{b_t}, {y_t} \ sim{{\離開({{Q_a} +{\ν_a} {S_a}} \右)}\ /{\氣_ {DF} ^ 2}} $ $
$ $ \σ_p ^ 2 \左|{\σ_a ^ 2 \σ_e ^ 2, {p_t}, {a_t},} \。{b_t}, {y_t} \ sim{{\離開({{Q_p} +{\ν_p} {S_p}} \右)}\ /{\氣_ {DF} ^ 2}} $ $
$ $ \σ_e ^ 2 \左|{\σ_a ^ 2 \σ_p ^ 2, {p_t}, {a_t},} \。{b_t}, {y_t} \ sim{{\離開({{Q_e} +{\ν_e} {S_e}} \右)}\ /{\氣_ {DF} ^ 2}} $ $
其中後兩行Qa, Qp和Qe分別表示誤差項和DF自由度的二次形式。我們運行了10000次迭代的模型，使用了5000次DBP的磨合周期。為了減少自相關性，我們每10次迭代采樣一次。我們嚐試了不同參數的反Wishart先驗分布，以獲得殘差。
我們使用一個混合模型，使用R軟件[10]進行全基因組關聯分析[9,7]:
$ $ {\ rm {y =}} \, {\ rm {Xb}} \, {\ rm { + }}\,{\ rm {e }}..........{\ rm {(4)}} $ $
其中y包含殘差或來自(3)的隨機效應，b表示固定效應(SNP)， X和為發生率$$矩陣，而\，e\，是\，a\，矢量\，包含\，殘差\，\，和\，假設\，正常\，分布\，與\，I\，\sigma _e^2。\，I\，是\，一個\，單位\，矩陣，\sigma _e^2$$是殘差方差。我們使用5%的錯誤發現閾值來檢測關聯映射中的基因組信號。我們還采用了橫截麵語法[1]方法，按每個時間點進行比較。利用基因組共祖矩陣[9]和譜係信息估計DBP的遺傳力分別為0.299和0.259。
在不了解底層仿真模型的情況下進行分析。但是，我們在討論結果時使用了GAW18的答案。
我們用Kolmogrow-Smirnov檢驗證實正常。然而，由於我們采用貝葉斯殘差:所有響應變量轉化為正態分布(P > 0.01)。時間、性別、吸煙狀況、年齡和係譜數作為固定效應納入後續分析，其基礎是利用預測和觀察結果之間的相關性進行初步分析。我們發現，預測和觀察結果之間的相關性最高可達0.15。用基因組關係矩陣預測誤差、遺傳和永久環境方差分量分別為17.3(0.0006)、10.9(0.0007)和8.2(0.0006);用譜係關係矩陣預測DBP誤差分量分別為18.0(0.0007)、9.2(0.0008)和8.4(0.0006)。DBP的模擬遺傳率為0.317，而基因組親緣關係估計值更接近其真實值(表1-3)。
關於基因進化和長期環境影響的假設可能很重要。與隨機遊走模型相比，一定程度的自回歸結構可能更加真實。我們根據(3)中參數空間的限製，使用∆對模型進行簡單的優化。我們考慮了偏離隨機漫步的兩種極端情況(∆=0.1和∆=0.9)。隨機遊走假設基因效應會隨著時間的推移在上下方向緩慢變化，∆=1.0，。通過調整∆，可以引入基因效應和環境效應的自回歸結構，得到平穩分布。在這裏，時間序列將圍繞平均軌跡分布。
誤差、遺傳和永久環境方差成分預測為10.26(0.0006)、1.98(0.0001)和1.98(0.0003)，∆=0.1和14.2(0.0002)，∆=0.9為DBP，預測為6.4(0.0009)和5.5(0.0006)。用∆=1.0，∆=0.1和∆=0.9預測遺傳力分別為0.299,0.139和0.246。與自回歸結構相比，隨機漫步的結果更好(DBP用0.317模擬)。我們假設:增加時間點應降低基因組膨脹因子[7]，因為隨著時間的推移，相關性和子結構信息的積累。我們在混合模型中同時使用了基於基因組和基於譜係的關係矩陣(3)。在DBP的三個時間點上，基因組關係矩陣給出的基因組膨脹因子分別為1.40、1.33和1.63，而基於譜係的關係矩陣為1.57、1.83和1.59。由於時間點較少(t=3)，我們仍然獲得了高水平的基因組膨脹因子(λ > 1)。從表1和表2可以看出，基因組關係和譜係關係矩陣在不同時間點檢測到的snp集合大多不同。
然而，小樣本量和小時間點都可能導致假陽性和假陰性。尤其對於第一個時間點，這可能是正確的:基因組關係矩陣在5%的假發現率(FDR)下檢測到154個snp(在5%的FDR下，時間點2和3分別檢測到134和56個snp)，譜係關係矩陣在5%的FDR下檢測到216個snp(在5%的FDR下，時間點2和3分別檢測到300和96個snp)。由於基因組膨脹因子較小，我們研究了因果snp的基因組關係矩陣結果。時間點1,2,3的rs11711953被發現與MAP4基因相關。
我們使用語法方法對每個時間點(以及它們的平均值)進行橫斷麵分析(表3)。然而，在多次假設修正後，我們沒有檢測到任何基因組信號。盡管rs1948722在ARHGEF3附近有來自時間點2的信號(p<0.00012)，但經過FDR的多次假設修正後，SNP變得不顯著。語法p值的幅度(表3)比隨機遊走模型的p值(表1,2)要大。這清楚地表明，隨著時間的推移，基因和環境的縱向影響需要考慮到適當的方法。否則，由於基因組信號將受到隨機誤差的汙染，這可能導致信號未被檢測到或也可能導致假陽性。隨機遊走是一個非平穩過程，本文通過調整∆參數擴展了平穩情況下的模型。然而，對於基因動力學和永久環境效應的非平穩假設是否有用，還需要進行理論和經驗的動態關聯研究。
與譜係信息相比，基因組關係矩陣給出了更好的膨脹因子和遺傳力估計。由於隨機遊走模型具有卡爾曼濾波的遞歸結構，在實際應用中可用於較長的時間序列。當縱向觀測可用時(如每日或每月)，模型可以通過卡爾曼濾波順序預測在線基因組信號。在基因組研究中，如果不考慮縱向基因和環境隨時間的影響，可能會導致無法檢測到真實信號和/或由於隨機誤差也可能導致假陽性。
遺傳分析講習班由美國國立衛生研究院(NIH)資助，由美國國立普通醫學科學研究所授予R01 GM031575贈款。這項工作得到了Akdeniz大學106號項目研究基金的支持。作者希望感謝與Luc Janss博士關於自動回歸結構的有益討論。
表2:DBP首次複製的譜係關係矩陣采用隨機遊走模型得到前10個SNPs和對應的原始p值
表3:使用第一次複製DBP的語法獲得的前10個snp和相應的原始p值
沒有相互競爭的利益。
- Aulchenko SY, de Koning DJ, Haley C(2007)使用混合模型和回歸的全基因組快速關聯:一種基於全基因組譜係的數量性狀位點關聯分析的快速簡單方法。遺傳學177:577 - 585。［Ref。］
- 吳銳，林敏(2006)功能定位——動態複雜性狀遺傳結構的定位與研究。Nat Rev Genet 7: 229- 237。［Ref。］
- Macgregor S, Knott SA, White I, Visscher PM(2005)複雜譜係中縱向定量性狀數據的數量性狀位點分析。遺傳學171:1365 - 1376。［Ref。］
- Faro EL, Cardoso VL, Albuquerque LG(2008)應用隨機回歸模型對卡拉蘇小母牛(牛牛科)測試日產奶量進行方差分量估計。Genet Mol Biol 31;［Ref。］
- Karacaören B (2014) F2小鼠數據集中使用祖先信息標記的生長相關性狀的混合映射。J Bioinform編譯生物學［Ref。］
- Karacaören B(2015)映射動態數量性狀的貝葉斯隨機遊走方法。J應用了非線性動力學。
- Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MA，等(2007)PLINK:全基因組關聯和基於群體的連鎖分析工具集。簡網81:559-575。［Ref。］
- West M, Harrison J(1997)貝葉斯預測和動態模型。施普林格。［Ref。］
- Endelman JB(2011)用R包rrBLUP進行基因組選擇的嶺回歸和其他核。植物Gen 4: 250-255。［Ref。］
- R開發核心團隊(2014)統計計算語言和環境。R統計計算基金會，維也納，奧地利。［Ref。］
在此下載臨時PDF
Aritcle類型:研究文章
引用:Karacaören B(2015)基於GAW 18數據集的卡爾曼濾波模型的動態關聯映射。Int J Mol基因et基因Ther 1 (1): http://dx.doi.org/10.16966/2471-4968.101
版權:©2015 Karacaören B.這是一篇開放獲取的文章，根據創作共用署名許可協議發布，該協議允許在任何媒體上不受限製地使用、發布和複製，前提是注明原作者和來源。
出版的曆史: