Learn R Programming

cpgen (version 0.1)

cGBLUP: Genomic BLUP

Description

This function allows fitting a mixed model with one random effect besides the residual using clmm. The random effect $\mathbf{a}$ follows some covariance-structure $\mathbf{G}$

Usage

cGBLUP(y,G,X=NULL, scale_a = 0, df_a = -2, scale_e = 0, df_e = -2,
          niter = 10000, burnin = 5000, seed = NULL, verbose=TRUE)

Arguments

y
vector of phenotypes
G
Relationship matrix / covariance structure for random effects
X
Optional Design Matrix for fixed effects. If omitted a column-vector of ones will be assigned
scale_a
prior scale parameter for $a$
df_a
prior degrees of freedom for $a$
scale_e
prior scale parameter for $e$
df_e
prior degrees of freedom for $e$
niter
Number of iterations
burnin
Burnin
seed
Seed
verbose
Prints progress to the screen

Value

  • List of 6:
  • var_ePosterior mean of the residual variance
  • var_aPosterior mean of the random-effect variance
  • bPosterior means of the fixed effects
  • aPosterior means of the random effects
  • posterior_var_ePosterior of the residual variance
  • posterior_var_uPosterior of the random variance

Details

Kang et al. (2008): $$\mathbf{y} = \mathbf{Xb} + \mathbf{a} + \mathbf{e} \textrm{ with: } \mathbf{a} \sim MVN(\mathbf{0},\mathbf{G}\sigma^2_a)$$

By finding the decomposition: $\mathbf{G = UDU'}$ and premultiplying the model equation by $\mathbf{U'}$ we get: $$\mathbf{U'y = U'Xb + U'a + U'e}$$ with: $$Var(\mathbf{U'y}) = \mathbf{U'G'U} \sigma^2_a + \mathbf{U'U} \sigma^2_e$$ $$\mathbf{U'UDU'U}\sigma^2_a + \mathbf{I}\sigma^2_e$$ $$\mathbf{D}\sigma^2_a + \mathbf{I}\sigma^2_e$$

After diagonalization of the variance-covariance structure the transformed model is being fitted by passing $\mathbf{D}^{1/2}$ as the design matrix for the random effects to clmm. The results are subsequently backtransformed and returned by the function.

References

Kang, H. M., N. A. Zaitlen, C. M. Wade, A. Kirby, D. Heckerman, M. J. Daly, and E. Eskin. "Efficient Control of Population Structure in Model Organism Association Mapping." Genetics 178, no. 3 (February 1, 2008): 1709-23. doi:10.1534/genetics.107.080101.

See Also

clmm, cgrm, cGWAS.emmax

Examples

Run this code
# generate random data
rand_data(500,5000)

# compute a genomic relationship-matrix
G <- cgrm(M,lambda=0.01)

# run model
mod <- cGBLUP(y,G)

Run the code above in your browser using DataLab