使用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表示样本个数。