User:Timothee Flutre/Notebook/Postdoc/2011/06/28
From OpenWetWare
(Difference between revisions)
(→Simple linear regression: add simulation with PVE and R code) 
m (→Simple linear regression) 

Line 88:  Line 88:  
set.seed(1859)  set.seed(1859)  
N < 100  N < 100  
  
mu < 5  mu < 5  
+  g < sample(x=0:2, size=N, replace=TRUE, prob=c(0.5, 0.3, 0.2)) # MAF=0.2  
beta < 0.5  beta < 0.5  
pve < 0.8  pve < 0.8  
beta.g.bar < mean(beta * g)  beta.g.bar < mean(beta * g)  
  sigma < sqrt((1/N) * sum((beta * g  beta.g.bar)^2) * (1pve) / pve) # 0.  +  sigma < sqrt((1/N) * sum((beta * g  beta.g.bar)^2) * (1pve) / pve) # 0.18 
y < mu + beta * g + rnorm(n=N, mean=0, sd=sigma)  y < mu + beta * g + rnorm(n=N, mean=0, sd=sigma)  
plot(x=0, type="n", xlim=range(g), ylim=range(y),  plot(x=0, type="n", xlim=range(g), ylim=range(y), 
Revision as of 09:10, 30 July 2012
Project name  Main project page Next entry 
Simple linear regression
In matrix notation: y = Xθ + ε with ε˜N_{N}(0,σ^{2}I_{N}) and θ^{T} = (μ,β)
Here is the ordinaryleastsquare (OLS) estimator of θ:
Let's now define 4 summary statistics, very easy to compute:
This allows to obtain the estimate of the effect size only by having the summary statistics available:
The same works for the estimate of the standard deviation of the errors:
We can also benefit from this for the standard error of the parameters:
with and V(ε) = σ^{2} This way, by also fixing β, it is easy to calculate the corresponding σ:
Here is some R code implementing this: set.seed(1859) N < 100 mu < 5 g < sample(x=0:2, size=N, replace=TRUE, prob=c(0.5, 0.3, 0.2)) # MAF=0.2 beta < 0.5 pve < 0.8 beta.g.bar < mean(beta * g) sigma < sqrt((1/N) * sum((beta * g  beta.g.bar)^2) * (1pve) / pve) # 0.18 y < mu + beta * g + rnorm(n=N, mean=0, sd=sigma) plot(x=0, type="n", xlim=range(g), ylim=range(y), xlab="genotypes (allele dose)", ylab="phenotypes", main="Simple linear regression") for(i in unique(g)) points(x=jitter(g[g == i]), y=y[g == i], col=i+1, pch=19) ols < lm(y ~ g) summary(ols) # muhat=5.01, betahat=0.46, R2=0.779 abline(a=coefficients(ols)[1], b=coefficients(ols)[2])
