遺伝率推定とMxとMixed Linear Model - ryamadaの遺伝学・遺伝統計学メモ

遺伝率とその推定のためのACEモデル・アプリケーション"Mx"
- MxはStructural Equation Model
Mixed Linear Model
いずれにしろ、次のようなことが基本になる
- 表現型の分散共分散行列 $V$ が、 $V=\sigma_G^2 G + \sigma_{F}^2 F$ のように分解できて、この分散共分散行列に従う多次元正規分布の点(多次元ベクトル〜個体数の長さのベクトル)が観測データであると考える
- そして、そのときの遺伝率が $\frac{\sigma_G^2}{\sigma_G^2+\sigma_F^2}$ であることから
- 観測データ(フェノタイプのデータベクトル)に、想定されるGを持ち込んで、 $\frac{\sigma_G^2}{\sigma_G^2+\sigma_F^2}$ の推定値を求めることが、遺伝率の推定
- ここで、双生児解析ならGには、「一卵性双生児の相関係数は１、二卵性は0.5」を決め打ちに使うし、同胞やその他の血縁関係の場合にも、血縁関係に基づく値を用いてGを作ることもあれば、計測ジェノタイプデータがあるならそれから、遺伝要因の重なり具合を推定してGの値とすることもある。血縁関係にない人を集めてきた場合にゲノムワイドSNPデータなどから、個人間のにかより具合を計算してGを作れば、ゲノムワイドSNPデータに基づく遺伝率推定となる
どうしてもわからないといけないのは
- 多変量正規分布
- 分散共分散行列
それを使って
- 一卵性双生児同士の表現型は、遺伝率が１ならいつも一致すること
- 遺伝率が０のときは、血縁関係の程度によらず、表現型はばらばらなこと
- 遺伝率が０から１の間のときには、血縁関係が強いほど表現型が似ていること
をシミュレーションしておこう

# install.packages("mvtnorm")
library(mvtnorm)

# 遺伝的な関係の強さを表す行列を作る
n <- 6 # 個体数
G <- diag(rep(1,n)) # 自身との相関は１
G
# 血縁関係に応じて相関の値を与える
G[1,2] <- G[2,1] <- 1 # s1,s2は一卵性双生児
G[3,4] <- G[4,3] <- 0.5 # s3,s4は二卵性双生児
G[5,6] <- G[6,5] <- 0.02 # s5,s6は少し血縁関係がある
# (s1,s2),(s3,s4),(s5,s6)以外の関係は全くの無関係
G

# こんな６個体をn.iter回、集めてみることにする
n.iter <- 10^4
# 多変量正規分布に従う乱数を発生させる

# 遺伝要因以外の寄与分は、独立(非対角成分は０)とする
E <- diag(rep(1,n))

# 遺伝率 k
# 遺伝率は、GとEとの割合のこと
par(mfrow=c(2,2))
ks <- c(1,0.8,0.5,0) # ４通りの遺伝率を試す
Vs <- Xs <- list()
for(i in 1:length(ks)){
	Vs[[i]] <- ks[i]*G + (1-ks[i])*E
	Xs[[i]] <- rmvnorm(n.sample,sigma=Vs[[i]])
	plot(Xs[[i]][,5],Xs[[i]][,6],xlim=range(Xs[[i]]),ylim=range(Xs[[i]]),cex=0.1,pch=20,xlab="one",ylab="the other",main=paste("k=",ks[i]))
	points(Xs[[i]][,3],Xs[[i]][,4],cex=0.1,pch=20,col=3)
	points(Xs[[i]][,1],Xs[[i]][,2],cex=0.1,pch=20,col=2)

}
par(mfcol=c(1,1))