Description Usage Arguments Value Author(s) Examples
Uses a GenABEL object and phenotype data as input. The model is fitted using the hglm
function in the hglm package.
1 2 | preFitModel(fixed = y ~ 1, random = ~1 | id, id.name = "id", genabel.data,
phenotype.data, corStruc = NULL, GRM = NULL, Neighbor.Matrix = NULL)
|
fixed |
A formula including the response and fixed effects |
random |
A formula for the random effects |
id.name |
The column name of the IDs in phen.data |
genabel.data |
An GenABEL object including marker information. This object has one observation per individual. |
phenotype.data |
A data frame including the repeated observations and IDs. |
corStruc |
A list specifying the correlation structure for each random effect. The options are: |
GRM |
A genetic relationship matrix. If not specified whilst the |
Neighbor.Matrix |
A neighborhood matrix having non-zero value for an element (i,j) where the observations i and j come from neighboring locations. The diagonal elements should be zero. |
Returns a list including the fitted hglm object fitted.hglm
, the variance-covariance matrix V
and the ratios between estimated variance components for the random effects divided by the residual variance, ratio
.
Lars Ronnegard
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 |
####### FIRST EXAMPLE USING GRM #############
data(Phen.Data) #Phenotype data with repeated observations
data(gen.data) #GenABEL object including IDs and marker genotypes
GWAS1 <- rGLS(y ~ age + sex, genabel.data = gen.data, phenotype.data = Phen.Data)
plot(GWAS1, main="")
summary(GWAS1)
#Summary for variance component estimation without SNP effects
summary(GWAS1@call$hglm)
#The same results can be computed using the preFitModel as follows
fixed = y ~ age + sex
Mod1 <- preFitModel(fixed, random=~1|id, genabel.data = gen.data,
phenotype.data = Phen.Data, corStruc=list( id=list("GRM","Ind") ))
GWAS1b <- rGLS(fixed, genabel.data = gen.data,
phenotype.data = Phen.Data, V = Mod1$V)
plot(GWAS1b, main="Results using the preFitModel function")
####### SECOND EXAMPLE USING CAR #############
# Add a fake nest variable to the data just to run the example
#In this example there are 6 nests and 60 observations per nest
Phen.Data$nest <- rep(1:6, each=60)
#A model including polygenic effects, permanent environmental effects,
#and nest effect as random
Mod2 <- preFitModel(fixed, random=~1|id + 1|nest, genabel.data = gen.data,
phenotype.data = Phen.Data, corStruc=list( id=list("GRM","Ind"), nest=list("Ind")) )
GWAS2 <- rGLS(fixed, genabel.data = gen.data, phenotype.data = Phen.Data, V = Mod2$V)
plot(GWAS2)
#Similar to previous plot because the nest effect variance component is almost 0.
###################
#Construct a fake nighbourhood matrix
D = matrix(0,6,6)
D[1,2] = D[2,1] = 1
D[5,6] = D[6,5] = 1
D[2,4] = D[4,2] = 1
D[3,5] = D[5,3] = 1
D[1,6] = D[6,1] = 1
D[3,4] = D[4,3] = 1
#The matrix shows which pair of nests that can be considered as neighbours
image(Matrix(D), main="Neighbourhood matrix")
Mod3 <- preFitModel(fixed, random=~1|id + 1|nest, genabel.data = gen.data,
phenotype.data = Phen.Data, corStruc=list( id=list("GRM","Ind"),
nest=list("CAR")), Neighbor.Matrix=D )
GWAS2b <- rGLS(fixed, genabel.data = gen.data,
phenotype.data = Phen.Data, V = Mod3$V)
plot(GWAS2b)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.