使用GCTA的GREML評估SNP遺傳力

• 2019 年 12 月 19 日
• 筆記

評估SNP遺傳力有兩種方法LDSC和GREML, 本文介紹下GREML評估遺傳力的方法。在GCTA軟件中，其核心就是如下所示的線性混合模型

其中y表示表型，b表示固定效應協變量的係數，u表示隨機效應自變量的係數, 這裡的隨機效應指的就是所有SNP位點對錶型的效應，e表示隨機誤差。其中u和e服從如下所示的正態分佈

表型方差可以用如下公式表示

在GCTA中，W的計算公式如下

在之前的文章中，我們介紹了GCTA對於樣本遺傳相似度的定義，公式如下

可以推導出如下關係

定義所有SNP位點的方差如下

則表型變異V的公式變換如下

定義SNP遺傳力的計算公式如下

此時只需要通過REML算法對線性混合模型的參數進行估計，然後求解上述方差，即可計算出了SNP遺傳力了。軟件實際操作的代碼如下

首先構建GRM矩陣，下游分析都是需要輸入GRM矩陣的，在實際分析中，還可以根據遺傳相似度進行過濾，去除親緣關係較近的樣本，上述代碼為了方便演示，沒有做過濾。

在進行REML分析時，有兩種選擇，第一種是默認用法，沒有協變量的校正，第二種添加了協變量的校正，注意這裡的協變量都是屬於固定效應協變量，協變量可以是離散型，也可以是連續型。

上述代碼中使用了兩個協變量，第一個為covar參數指定的離散型協變量，第二個為qcovar參數指定的連續型協變量，樣本的PCA分析的前10主成分作為協變量，即使用群體分層這個因素作為協變量。

輸出結果後綴為hsq, 內容示意如下

Variance表示方差，SE表示標準誤，其中V(G)/Vp即為SNP遺傳力，下面幾行為LRT檢驗，即似然比檢驗的結果，其中LRT = 2(logL – logL0), df表示自由度，Pval表示檢驗的p值，n表示樣本個數。