GxM: Maximum Likelihood Estimation of Gene-by-Measured Environment...

In GxM: Maximum Likelihood Estimation for Gene-by-Measured Environment Interaction Models

Description

This function is used to produce maximum likelihood estimation for model fittings. With identification of the model fitted to the given data set with certain constrains of parameters, both numerical method and closed form based method are available to obtain the maximum likelihood estimate of parameters. Hessian matrix of log-likelihood at the estimate is provided to obtain the standard errors of the estimate.

Usage

GxM(dataset, rGvalue=0.995, modelname, 
    zeroset=character(), closedform=FALSE, K=8, coreNumber=1, 
    manualinitial=NULL, priority=1, gradientlevel=2)

Arguments

dataset	a data frame of input. This data frame should contain 6 columns, including the id of twin/sibling pairs, the rG value for the pairs, as well as the observed M1, M2, P1 and P2 values. The M1, M2, P1, P2 data would NOT be automatically standardized.
rGvalue	a correlation value of genetic influences to a pair of twins/siblings. A value close to one is used to replace exact one to avoid calculations with singular matrices. This value cannot be the same for all pairs for the same reason.
modelname	a character string input specifying which model should be fitted. Possible choices include Chol,CholGxM,CholNonLin,NLMainGxM,CorrGxM and CorrNonLin. Chol P = μ_P + a_C*A_M + c_C*C_M + e_C*E_M + a_U*A_U + c_U*C_U + e_U*E_U. CholGxM P = μ_P + (a_C+α_C*M)*A_M + (c_C+κ_C*M)*C_M + (e_C+ε_C*M)* E_M + (a_U+α_U*M)*A_U + (c_U+κ_U*M)*C_U + (e_U+ε_U*M)*E_U. CholNonLin P = μ_P + a_C*A_M + c_C*C_M + e_C*E_M + a_U*A_U + c_U*C_U + e_U*E_U + γ_1*A^2_M + γ_2*C^2_M + γ_3*E^2_M + δ_1*A_M*A_U + δ_2*C_M*C_U + δ_3*E_M*E_U. NLMainGxM P = μ_P + β_1*M + β_2*M^2 + (a_U+α_U*M)*A_U + (c_U+κ_U*M)*C_U + (e_U+ε_U*M)*E_U. CorrGxM P = μ_P + (a_P+α_P*M)*A_P + (c_P+κ_P*M)*C_P + (e_P+ε_P*M)*E_P. r_A = corr(A_P,A_M); r_C = corr(C_P,C_M); r_E = corr(E_P,E_M) CorrNonLin P = μ_P + a_P*A_P + c_P*C_P + e_P*E_P + λ_1*A_M*A_P + λ_2*C_M*C_P + λ_3*E_M*E_P. r_A = corr(A_P,A_M); r_C = corr(C_P,C_M); r_E = corr(E_P,E_M)
zeroset	a set of parameters that would be set zeros in the estimation. Character string for parameter names would be used for identification.
closedform	logical type indicating whether closed form method or numerical method should be applied to likelihood calculation.
K	an integer indicating the number of points in each dimension in numerical integration.
coreNumber	a positive integer representing whether and how parallel computation should be used. The fault value is set to be 1, and no parallel computation is necessary. A number greater than 1 indicates how many processors would be used in parallel computation.
manualinitial	a list of numerical values used for initializing the parameters excluding members in the zeroset.
priority	a positive number specifying the priority of manual initialization. It takes values in 1 and 2. Priority=2 implies that the manual initialization would override possible update built in the optimization.
gradientlevel	a positive value which serves as a threshold. A warning information would be produced if the maximum absolute value of gradient components exceeds this level.

Value

An S4 class of output with 10 slots would be returned.

loglikelihood	maximum log-likelihood value from fitting indicated model to the input data set.
BIC	Bayesian information criterion (BIC) value from the model fitting.
par	parameter estimates.
hess	Hessian matrix of log-likelihood at the maximum point with respect to parameters.
gradient	gradient components of log-likelihood at the maximum point with respect to parameters.
modelname	inherited from input.
zeroset	inherited from input.
closedform	inherited from input.
K	inherited from input.
coreNumber	inherited from input.

References

Rathouz, P. J., Van Hulle, C. A., Rodgers, J. L., Waldman, I. D., Lahey, B. B. (2008) Specification, testing, and interpretation of gene-by-measured-environment interaction models in the presence of gene-environment correlation. Behavior Genetics, 38(3), 301-315.

See Also

checkGxM, summaryGxM

Examples

id <- 1:10;
rG <- c(rep(1,5),rep(0.5,5));
M1 <- rnorm(10);
M2 <- rnorm(10);
P1 <- rnorm(10);
P2 <- rnorm(10);
Data <- data.frame(id,rG,M1,M2,P1,P2);
output <- GxM(dataset = Data, modelname='Chol', closedform=TRUE);
output;
output@BIC;
output@hess;

GxM documentation built on May 2, 2019, 3:38 p.m.