1.6 用GCTA计算亲缘关系矩阵和遗传力
**第一步,**首先利用命令--make-grm
计算一个基因相关性矩阵GRM。
$ gcta --bfile 1kg_EU_BMI --autosome --maf 0.01 --make-grm --out 1kg_gcta
--bfile
,类似与plink,指定用于分析的基因型文件 的名字。
--autosome
,在计算GRM时用于筛选常染色体。
--maf
,用于筛选最小等位基因频率小小于0.01的SNPs。
--make-grm
,计算GRM
--out
,指定输出文件的名字
输出:主要输出四个文件,
1kg_gcta.grm.id,包含个体IDs;1kg_gcta.grm.bin,二进制文件;1kg_gcta.grm.N.bin,另一个二进制文件,通过GCT去描述GRM;1kg_gcta.log,log文件。
**第二步,**从样本中去除相关个体。去除相关性大于0.025,其中的一个个体被保留下来。创建了一个新的文件名1kg_rm025,去除了59个个体,仅留下320个个体。
$ gcta --grm 1kg_gcta --grm-cutoff 0.025 --make-grm --out 1kg_rm025
第三步是评估SNP遗传力,需提供表型文件。BMI_pheno.txt
$ gcta --grm 1kg_rm025 --pheno BMI_pheno.txt --reml --out 1kg_BMI_h2
生成1kg_BMI_h2.hsq文件。
$ more 1kg_BMI_h2.hsq
Source Variance SE
V(G) 3.363073 8.073311
V(e) 5.260860 8.059914
Vp 8.623934 0.683136
V(G)/Vp 0.389970 0.934748
logL -506.006
logL0 -506.067
LRT 0.122
df 1
Pval 3.6369e-01
n 320
GCTA描述BMI变异主要有2部分,V(G),加性遗传学解释,V(e),环境造成的。
遗传力评估( h S N P 2 h^2_{SNP} hSNP2)是V(G)占表型总方差V§的比例,结果是0.39。因此我们评估近40%的这种表型(BMI)的差异可归因于遗传。
一个个体的遗传力评估精确度是标准误(SE),在本例中SE达到0.93,这是非常高的,是点估计值的两倍多。这意味着估计是不精确的。GCTA执行一种名为log -似然比检验的统计检验,该检验表明遗传力评估的估计值与零假设有何不同。在本例中**,给定0.36的值**,我们不能拒绝零假设,没有证据表明BMI有遗传成分。实质上,BMI中20%-30%是由遗传导致,造成上述结果主要是因为样本量太少。
参考:
An Introduction to Statistical Genetic Data Analysis.