Description Usage Arguments Details Value Author(s) Examples
To make a quadratic approximation to the likelihood function, the score and information are obtained from a prebuilt matrix of weighted second moments. This allows a parameter estimate to be obtained by one iteration of weighted least squares, or equivalently a score test. Typically the weights used to construct the prebuilt matrix correspond to the MLE under a chosen null model.
1  est.moments2(xtwx, leftvar, rightvars, n = NULL, vscale = NULL)

xtwx 
an object of class moments2, typically built using

leftvar 
name of the response variable (the left hand side of the formula). 
rightvars 
name(s) of the explanatory variables (the right hand side of the formula). 
n 
sample size, only needed for the normal linear model if there is not a single intercept
“ 
vscale 
parameter needed if 
Variables in rightvars
with nonidentifiable coefficients are
removed, with preference for keeping variables that occur earlier
rather than later in rightvars
.
When the vscale
attribute of xtwx
(or the vscale
function argument) is NULL
, this function assumes
that the xtwx
argument was calculated with unit weights and
therefore that a linear model fit is required with error variance
estimated from the data. For this application it is preferred to call
lm.moments2
, which is a wrapper for this function with
vscale=NULL
.
When the vscale
attribute of xtwx
(or the vscale
function argument) is set equal to 1, this function assumes that the
xtwx
argument was calculated with weights calculated such that
a GLS problem has been correctly set up to approximate a likelihood
function, and therefore that generalised linear model fit is required.
Values other than NULL
or 1 for the vscale
parameter may not be
what you think. Do not use other values unless you are absolutely sure what you understand
what are doing. See the source code for details.
A list with slots for the effect size estimates, standard errors, and a precision matrix.
Toby Johnson Toby.x.Johnson@gsk.com
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  data(mthfrex)
mthfrex < gls.approx.logistic(mthfrex, "HTN", c("SexC", "Age"))
xtwx < make.moments2(mthfr.params, c("HTNstar", "SexC", "Age"), mthfrex,
weightvar = "weight")
myglm < est.moments2(xtwx, "HTNstar", c("ONE", "rs6668659_T", "rs4846049_T",
"rs1801133_G", "SexC", "Age"), vscale=1)
myglm$z < myglm$betahat/myglm$se
cbind(beta = myglm$betahat, se = myglm$se, z = myglm$z,
pval = pnorm(abs(myglm$z))*2)
## Compare against results from glm
## Note have to use coded alleles used in original data
mycheck < glm(HTN ~ rs6668659_G+rs4846049_G+rs1801133_A+Sex+Age,
family="binomial", data = mthfrex$data)
coef(summary(mycheck))
## Note in results Sex factor coded differently than SexC
## Coefficients for covariates used in null model are different,
## because xtwx approximates around the fitted null model
## Look at pairwise correlations
cor(subset(mthfrex$data, select = c("rs6668659_G", "rs4846049_G",
"rs1801133_A")))^2
## SNP coefficients well approximated (given very high
## interSNP correlations) but signs ALL inverted by coded allele flips
## check error less than 10percent
stopifnot(all(1*myglm$z[2:4]/coef(summary(mycheck))[2:4,3] > 0.9))
stopifnot(all(1*myglm$z[2:4]/coef(summary(mycheck))[2:4,3] < 1.1))

Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.