7. Genetic covariances computation {SFSI}

Pairwise Genetic Covariance

Description

Pairwise genetic covariances for variables with the same experimental design

Usage

getGenCov(y, X = NULL, Z = NULL, K = NULL, trn = NULL,
          EVD = NULL, ID_geno, ID_trait, scale = TRUE,
          pairwise = TRUE, verbose = TRUE, ...)
          

Arguments

y

(numeric matrix) Response variable vector. It can be a matrix with each column representing a different response variable

X

(numeric matrix) Design matrix for fixed effects. When X=NULL a vector of ones is constructed only for the intercept (default)

Z

(numeric matrix) Design matrix for random effects. When Z=NULL an identity matrix is considered (default) thus G = K; otherwise G = Z K Z' is used

K

(numeric matrix) Kinship relationship matrix

trn

(integer vector) Which elements from vector y are to be used for training. When trn = NULL, non-NA entries in vector y will be used as training set

EVD

(list) Eigenvectors and eigenvalues from eigenvalue decomposition (EVD) of G corresponding to training data

ID_geno

(character/integer) Vector with either names or indices mapping entries of the vector y to rows/columns of matrix G

ID_trait

(character/integer) Vector with either names or indices mapping entries of the vector y to different traits

scale

TRUE or FALSE to scale each column of y by their corresponding standard deviations so the resulting variables will have unit variance

pairwise

TRUE or FALSE to calculate pairwise genetic covariances for all columns in y. When pairwise = FALSE (default) covariances of the first column in y with the remaining columns (2,...,ncol(y)) are calculated

verbose

TRUE or FALSE to whether show progress

...

Other arguments passed to the function 'fitBLUP'

Details

Assumes that both y₁ and y₂ follow the basic linear mixed model that relates phenotypes with genetic values of the form

y₁ = X b₁ + Z g₁ + e₁

y₂ = X b₂ + Z g₂ + e₂

where b₁ and b₂ are the specific fixed effects, g₁ and g₂ are the specific genetic effects of the genotypes, e₁ and e₂ are the vectors of specific environmental residuals, and X and Z are common design matrices for the fixed and genetic effects, respectively. Genetic effects are assumed to follow a Normal distribution as g₁ ~ N(0,σ²_u1K) and g₂ ~ N(0,σ²_u2K), and environmental terms are assumed e₁ ~ N(0,σ²_e1I) and e₂ ~ N(0,σ²_e2I).

The genetic covariance σ²_u1,u2 is estimated from the formula for the variance for the sum of two variables as

σ²_u1,u2 = 1/2(σ²_u3 - σ²_u1 - σ²_u2)

where σ²_u3 is the genetic variance of the variable y₃ = y₁ + y₂ that also follows the same model as for y₁ and y₂.

Likewise, the environmental covariance σ²_e1,e2 is estimated as

σ²_e1,e2 = 1/2(σ²_e3 - σ²_e1 - σ²_e2)

where σ²_e3 is the error variance of the variable y₃.

Solutions are found using the function 'fitBLUP' (see help(fitBLUP)) to sequentialy fit mixed models for all the variables y₁, y₂ and y₃.

Value

Returns a list object that contains the elements:

• varU: (vector) genetic variances.

• varE: (vector) error variances.

• covU: (vector) genetic covariances between response variable 1 and the rest.

• covE: (vector) environmental covariances between response variable 1 and the rest.

When pairwise = TRUE, varU and varE are matrices containing all variances (diagonal) and pairwise covariances (off diagonal)

Examples

  require(SFSI)
  data(wheatHTP)
  
  index = which(Y$trial %in% 1:10)    # Use only a subset of data
  Y = Y[index,]
  M = scale(M[index,])/sqrt(ncol(M))  # Subset and scale markers
  G = tcrossprod(M)                   # Genomic relationship matrix
  Y0 = scale(as.matrix(Y[,4:7]))      # Response variable
  
  ID_geno = as.vector(row(Y0))
  ID_trait = as.vector(col(Y0))
  y = as.vector(Y0)
  
  # Pairwise covariance for all columns in y
  fm = getGenCov(y=y, K=G, ID_geno=ID_geno, ID_trait=ID_trait)
  fm$varU           # Genetic variances
  fm$varE           # Residual variances
  fm$varU+fm$varE   # Phenotypic covariance
  var(Y0)           # Sample phenotypic covariance
  
  # Covariances between the first and the rest: y[,1] and y[,2:4]
  fm = getGenCov(y=y, K=G, ID_geno=ID_geno, ID_trait=ID_trait, pairwise=FALSE)
  fm$covU           # Genetic covariance between y[,1] and y[,2:4]
  
  # Containing some missing values:
  # The model is estimated from common training data
  YNA = Y0
  YNA[Y$trial %in% 1:2, 1] = NA
  YNA[Y$trial %in% 2:3, 3] = NA
  y = as.vector(YNA)
  
  fm = getGenCov(y=y, K=G, ID_geno=ID_geno, ID_trait=ID_trait)
  fm$varU 
  fm$varE