부모-자식 데이터를 통한 유전율 추정
본 포스팅은 부모-자식 데이터를 통해 heritability를 추정하는 것을 포스팅한다. 아래 데이터는 parent와 child의 height를 나타내는 수치이다. 옥수수의 키라고 생각하자. 알고 싶은 것은 child의 키가 유전적 요소로 설명되는 정도이다. 우리는 이 수치를 heritability(h^2)라는 값으로 표현한다. 데이터는 이곳에서 다운로드 할 수 있다.
h2=V(G)V(Y)
library(data.table)
## Warning: package 'data.table' was built under R version 3.4.3
load("C:\\workspace\\r\\genetics\\data\\mod7_data.R")
head(Galton)
## parent child
## 1 71.40788 61.7
## 2 68.57945 61.7
## 3 64.33681 61.7
## 4 62.92260 61.7
## 5 62.21549 61.7
## 6 67.16524 62.2
데이터에 약간의 random error를 주어 아래의 그래프를 통해 표현해보자.
plot(Galton$parent+rnorm(928,0,.1),Galton$child+rnorm(928,0,.1),xlab="Father's Height",ylab="Child's Height")
abline(v=mean(Galton[,1]),col=2,lty=3,lwd=3)
abline(h=mean(Galton[,2]),col=2,lty=3,lwd=3)
abline(a=0,b=1,col=2,lty=3,lwd=3)
summary(lm(Galton$child~Galton$parent))
##
## Call:
## lm(formula = Galton$child ~ Galton$parent)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.8050 -1.3661 0.0487 1.6339 5.9264
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 36.87187 1.98827 18.55 <2e-16 ***
## Galton$parent 0.45700 0.02909 15.71 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.239 on 926 degrees of freedom
## Multiple R-squared: 0.2105, Adjusted R-squared: 0.2096
## F-statistic: 246.8 on 1 and 926 DF, p-value: < 2.2e-16
회귀계수가 0.457이다. 이 때, 유전율은 0.457*2 = 약 92%로 추정된다. 왜 그런지 살펴보기 위해 우선 beta의 추정치를 아래와 같이 회귀식의 양변에 Cov(, X)를 취해 구해보자.
Y=α+βX Cov(Y,X)=βCov(X,X) ˆβ=Cov(Y,X)Var(X)
Phenotype Y에 대해 아래와 같은 일반적인 Phenotype에 대한 모델을 세우자. 이 때, A,D,C,E 각각의 컴포넌트들에 대한 것은 나중에 살펴보고 i, j 두 개체의 Covariance만 살펴본다. A : additive effect D : dominance effect C : shared environment effect E : residuals Y=A+D+C+E r은 relatedness, k2는 genome상의 IBD2의 비율, a는 shared environment를 나타내는 수치이다. 만약 shared environment이면 a=1, 아니면 a=0이다.
이 Covariance 모델은 어디까지나 가정이다. 이 가정 하에서 각각의 Variance를 추정하고 이를 통해 heritability를 추정하게 된다.
E[Cov(Yi,Yj)]=rVar(A)+k2Var(D)+aVar(C)
이제 현재 데이터에서 위 식을 적용하자. 이 때, r=0.5, k2=0, a=0을 가정하자. (부모-자식 간의 IBD2=0이므로 k2=0이다. ) 그러면 랜덤으로 i번째 부모-자식 pair를 골랐을 때, Covariance 모델은 아래와 같다. (P : Parent, C : Child의 약자)
이러면 Var(D)가 사라지기 때문에 우리는 쉽게 Var(A)를 데이터를 통해 추정할 수 있다.
E[Cov(Pi,Ci)]=0.5Var(A) Var(A)=2E[Cov(Pi,Ci)]
print(var(Galton$parent))
## [1] 6.389121
print(var(Galton$parent, Galton$child)) # Sample Covariance
## [1] 2.919806
따라서 위 데이터에는 D(Dominance effect가 없으므로 아래와 같이 유전율을 추정할 수 있다.) h2=V(G)V(Y)=V(A)V(P)
따라서 2.91*2 = 5.82 = Var(A) 이다. Var(P) = 6.39 (Phenotype의 분산을 V(P)로 추정) 이므로 이 모델하에서 추정한 유전율은 아래와 같다.
V(A)V(P)=5.82/6.39=0.91
print(paste0("Heritability : ", 2*var(Galton$parent, Galton$child) / var(Galton$parent)))
## [1] "Heritability : 0.913992906289406"
따라서 부모-자식 데이터에 이러한 모델을 통해 유전율을 추정할 때 유전율은 회귀계수의 2배를 해주면 된다.
참고
'Domains > Genetics' 카테고리의 다른 글
genetic score 를 통해 heritability 를 구하는 방법 (1) | 2023.08.30 |
---|---|
Post-GWAS 와 Functional Genomics 의 등장 (0) | 2019.08.23 |
Solar 설치 (0) | 2018.01.29 |
유전자 발현 조절 : 유전자는 어떻게 세포의 기능을 결정하는가? (2) | 2018.01.24 |
진핵생물에서의 유전자 발현 조절 과정 (0) | 2018.01.07 |