Testing non-genetic predictors in polygenic mixed model

Questions about GenABEL (aka *ABEL) suite of packages
Forum rules
Please remember not to post any sensitive data on this public forum.
The first few posts of newly registered users will be moderated in order to filter out any spammers.

When get a solution to the problem you posted, please change the topic name (e.g. from "how to ..." to "[SOLVED] how to ..."). This will make it easier for the community to follow the posts yet to be attended.
Posts: 2
Joined: Fri Jan 31, 2014 5:04 pm

Testing non-genetic predictors in polygenic mixed model

Postby NickM » Fri Jan 31, 2014 5:15 pm


I need to fit a model that tests non-genetic predictors of interest while accounting for the relatedness between individuals. For example, in the simplified simulated example below, I have a response variable, Y, and a fixed effect, X, measured on 20 subjects:

N <- 20
Subj <-factor(1:N)
X <- rnorm(N, 0.2, 4)
Y <- rnorm(N, X*0.4, 6)
Data <- data.frame(Subj, X, Y)

I also have a hypothetical relatedness matrix (Rmat):

Rmat <- diag(N)
corrs <- runif(N*(N-1)/2, 0, 1)
Rmat[lower.tri(Rmat)] <-corrs
Rmat <- t(Rmat)
Rmat[lower.tri(Rmat)] <-corrs
rownames(Rmat)<-Subj; colnames(Rmat)<-Subj

Can I use GenABEL to fit a mixed model that tests for the association between X and Y while accounting for the relatedness/ kinship matrix? I seem to be getting close using the polygenic function, but this doesn't provide standard errors/p-values for the fixed effect(s), while if I try to run grammar it seems to want genetic/SNP data as my predictors of interest.

Any help/ pointers greatly appreciated!

Many thanks,


Nicola Pirastu
GenABEL senior expert
GenABEL senior expert
Posts: 151
Joined: Wed Feb 09, 2011 3:24 pm

Re: Testing non-genetic predictors in polygenic mixed model

Postby Nicola Pirastu » Sun Feb 09, 2014 8:43 am


you have basically three choiches:

1) use polygenic in combiantion with the likelihood ratio test, but you will not get standard errors and if you add more than one variable at a time you will not get individual variables p-values. In this case you should first fit a null model like:

Code: Select all

nullmod=polygenic(y~sex+age,data=yourdata, kinship=yourkinship)

Then you fit the full model:

Code: Select all

yourmod=polygenic(y~sex+age+ yourtest.variable,data=yourdata, kinship=yourkinship)

then you can use the likelyhood ratio test.

Code: Select all

p=pchisq(nullmod$h2an$minimum - yourmod$h2an$minimum, df=1, lower=F)

Beta estimates are in yourmod$h2an$estimate. They do not present names but they are ordered in the order you put them. The vector is 2 elements longer since heritability and residual variance are there as well. Remember that if you want this to make sense you have to use exactly the same people in the 2 models (model.frame() is great to achieve this). Also if you add more than 1 variable in the test model, remember to add degrees of freedom when you estimate the pvalue.

2) use polygenic-hglm which will give you the effects and the standard errors and everything else, but is a hierarchical mixed model which is not exactly the same.

3) use GRAMMAR+, in this case you first regress out kinship, you get the enviromental corrected residuals ($grresidualY) and then use this to run a normal regression with glm.

Hope this solves your problems.



Posts: 2
Joined: Fri Jan 31, 2014 5:04 pm

Re: Testing non-genetic predictors in polygenic mixed model

Postby NickM » Fri Feb 14, 2014 10:35 am

Hi Nicola,

Thanks very much for your response. The first two solutions you offered seem to address my problem perfectly.

For your third proposal, if I were to fit the "enviromental corrected residuals" in a regular regression, wouldn't this assume that each individual is unrelated (i.e. independent)? I can see how regressing out the kinship would account for any confounding (such as that due to population structure), but am struggling to understand how any relatedness/ non-independence of individuals would be addressed using this approach.

Best regards,


Return to “GenABEL”

Who is online

Users browsing this forum: No registered users and 1 guest