inst/md/MULTITRAIT.md

In gdlc/BGLR-R: Bayesian Generalized Linear Regression

Multi-trait prediction using SVD and univariate regressions on eigenvectors

Deriving Eigenvectors of the phenotype matrix

 library(BGLR)
 data(wheat)

 # A phenotype matrix with 4 traits
  Y=wheat.Y 

 # Molecular markers
  X=scale(wheat.X) 

 # Singular value decompositions Y=UDV'
  SVD=svd(Y)
  U=SVD$u
  D=diag(SVD$d)
  V=SVD$v

Obtaining regression coefficients

 B=matrix(nrow=ncol(X),ncol=ncol(Y))

 ETA=list(list(X=X,model='BayesB'))
 for(i in 1:ncol(Y)){
    fm=BGLR(y=U[,i],ETA=ETA,verbose=F) #use more iterations!
    B[,i]=fm$ETA[[1]]$b
 }

 # Rotating coefficients to put them in marker space
  BETA=B%*%D%*%t(SVD$v)

 # Prediction
  YHat=X%*%BETA

  # correlation in training
  for(i in 1:ncol(Y)){ print(cor(Y[,i],YHat[,i])) }

Training-testing comparing single and multi-trait

Training/Testing partition

  set.seed(195021)
  N=nrow(X)
  tst=sample(1:N,size=300)
  Y.TRN=Y[-tst,]
  Y.TST=Y[tst,]

  X.TST=X[tst,]

  ETA=list(list(X=X[-tst,],model='BayesB'))

Prediction using single-trait models

  YHAT.ST=matrix(nrow=nrow(Y.TST),ncol=4)

  for(i in 1:4){
    fm=BGLR(y=Y.TRN[,i],ETA=ETA,verbose=F) #use more iterations!
    YHAT.ST[,i]=X.TST%*%fm$ETA[[1]]$b
  } 

Prediction using eigenvectors

  SVD=svd(Y.TRN)
  U=SVD$u
  D=diag(SVD$d)
  V=SVD$v

  B=matrix(nrow=ncol(X), ncol=ncol(U))

  for(i in 1:4){
     fm=BGLR(y=U[,i],ETA=ETA,verbose=F) #use more iterations!
     B[,i]=fm$ETA[[1]]$b
  } 
  BETA=B%*%D%*%t(V)
  YHAT.EV=X.TST%*%BETA

Comparison

  for(i in 1:4){
     message( round(cor(Y.TST[,i],YHAT.ST[,i]),3),"   ",round(cor(Y.TST[,i],YHAT.EV[,i]),3))
  }

Return to examples

gdlc/BGLR-R documentation built on April 23, 2024, 11:01 p.m.