- 一部の座位はわかっていてジェノタイプが得られる
- 一部の座位はわかっいないとする
n.marker <- 10
N <- 10^5
g <- matrix(0,N,n.marker)
f <- 0.1 + runif(n.marker) * 0.4
for(i in 1:n.marker){
g[,i] <- sample(0:2,N,replace=TRUE,prob=c(f[i]^2,2*f[i]*(1-f[i]),(1-f[i])^2))
}
n.marker.unknown <- 100
g.unknown <- matrix(0,N,n.marker.unknown)
f <- 0.1 + runif(n.marker) * 0.4
for(i in 1:n.marker){
g.unknown[,i] <- sample(0:2,N,replace=TRUE,prob=c(f[i]^2,2*f[i]*(1-f[i]),(1-f[i])^2))
}
r <- rnorm(n.marker)
n.unknown.true <- 10
r.unknown <- rep(0,n.marker.unknown)
r.unknown[1:n.unknown.true] <- rnorm(n.unknown.true)
R <- rep(0,N)
for(i in 1:n.marker){
R[which(g[,i]==1)] <- R[which(g[,i]==1)] + r[i]
R[which(g[,i]==2)] <- R[which(g[,i]==2)] + r[i]*2
}
R.unknown <- rep(0,N)
for(i in 1:n.marker){
R.unknown[which(g.unknown[,i]==1)] <- R.unknown[which(g.unknown[,i]==1)] + r.unknown[i]
R.unknown[which(g.unknown[,i]==2)] <- R.unknown[which(g.unknown[,i]==2)] + r.unknown[i]*2
}
R.total <- R + R.unknown
h <- 0.7
eV <- (1-h)/h * var(R.total)
R. <- R.total + rnorm(N,0,sqrt(eV))
pheno <- R.
lm.out.whole <- lm(pheno~cbind(g,g.unknown))
pheno.whole <- lm.out.whole$fitted
hwhole <- var(pheno.whole)/var(pheno)
hwhole
lm.out.partial <- lm(pheno~g)
pheno.partial <- lm.out.partial$fitted
hpartial <- var(pheno.partial)/var(pheno)
hpartial
- もちろん、全部の座位を使えば遺伝率がきちんと算出され、見つかっていなければ、低く計算される
> hwhole
[1] 0.6980088
>
> lm.out.partial <- lm(pheno~g)
> pheno.partial <- lm.out.partial$fitted
>
> hpartial <- var(pheno.partial)/var(pheno)
> hpartial
[1] 0.5391937