ranefUvcovNew {lme4GS}R Documentation

Obtain BLUPs for new levels of random effects with user specified variance covariance-matrices.

Description

Obtain BLUPs for new levels of random effects with user specified variance covariance-matrices.

Usage

  ranefUvcovNew(object,Uvcov)

Arguments

object

is an object returned by lmerUvcov.

Uvcov

two level list with ids to be predicted and variance covariance matrix that contains information of these ids and the ids used to fit the model.

Details

Assume that the random effect uj ~ N(0, σ2jKj) and the matrix Kj is partitioned as follows:

uj=(uj1 uj2)'

and

Kj.png

The BLUP for \boldsymbol u_{j2} can be obtained as:

E(uj2|y1)=Kj21Kj11-1ûj1

Value

A list of matrix arrays one for each grouping factor

Author(s)

Paulino Perez-Rodriguez

References

Caamal-Pat D., P. Perez-Rodriguez, J. Crossa, C. Velasco-Cruz, S. Perez-Elizalde, M. Vazquez-Pena. 2021. lme4GS: An R-Package for Genomic Selection. Front. Genet. 12:680569. doi: 10.3389/fgene.2021.680569 doi: 10.3389/fgene.2021.680569

Examples



 
library(BGLR)
library(lme4GS)

data(wheat)
X<-wheat.X
Z<-scale(X,center=TRUE,scale=TRUE)
G<-tcrossprod(Z)/ncol(Z)
A<-wheat.A
rownames(G)<-colnames(G)<-rownames(A)
y<-wheat.Y[,1]

#Predict 10/100 of records selected at random. 
#The data were partitioned in 10 groups at random
#and we predict individuals in group 2.

fold<-2
y_trn<-y[wheat.sets!=fold]
y_tst<-y[wheat.sets==fold]

A_trn=A[wheat.sets!=fold,wheat.sets!=fold]
G_trn=G[wheat.sets!=fold,wheat.sets!=fold]

pheno_trn=data.frame(y_trn=y_trn,m_id=rownames(A_trn),a_id=rownames(G_trn))

#######################################################################################
#Marker based prediction
#######################################################################################
	
fm1<-lmerUvcov(y_trn~1+(1|m_id),data=pheno_trn,Uvcov=list(m_id=list(K=G_trn)))

plot(pheno_trn$y_trn,predict(fm1),xlab="Phenotype",ylab="Pred. Gen. Value")

#BLUP for individuals in the testing set
blup_tst<-ranefUvcovNew(fm1,Uvcov=list(m_id=list(K=G)))
blup_tst<-blup_tst$m_id[,1]

#Comparison
#Check the names
names(y_tst)<-rownames(G)[wheat.sets==fold]
blup_tst<-blup_tst[match(names(y_tst),names(blup_tst))]

yHat_tst<-fixef(fm1)[1]+blup_tst
points(y_tst,yHat_tst,col="red",pch=19)

#Correlation in testing set
cor(y_tst,yHat_tst)

#######################################################################################
#Pedigree based prediction
#######################################################################################

fm2<-lmerUvcov(y_trn~1+(1|a_id),data=pheno_trn,Uvcov=list(a_id=list(K=A_trn)))

plot(pheno_trn$y_trn,predict(fm2),xlab="Phenotype",ylab="Pred. Gen. Value")

#BLUP for individuals in the testing set
blup_tst<-ranefUvcovNew(fm2,Uvcov=list(a_id=list(K=A)))
blup_tst<-blup_tst$a_id[,1]

#Comparison
#Check the names
names(y_tst)<-rownames(A)[wheat.sets==fold]
blup_tst<-blup_tst[match(names(y_tst),names(blup_tst))]

yHat_tst<-fixef(fm2)[1]+blup_tst
points(y_tst,yHat_tst,col="red",pch=19)

#Correlation in testing set
cor(y_tst,yHat_tst)

#######################################################################################
#Markers + Pedigree based prediction

fm3<-lmerUvcov(y_trn~1+(1|m_id)+(1|a_id),data=pheno_trn,
               Uvcov=list(m_id=list(K=G_trn),a_id=list(K=A_trn)))

plot(pheno_trn$y_trn,predict(fm3),xlab="Phenotype",ylab="Pred. Gen. Value")

#BLUP for individuals in the testing set
blup_tst<-ranefUvcovNew(fm3,Uvcov=list(m_id=list(K=G),a_id=list(K=A)))

blup_tst_m<-blup_tst$m_id[,1]
blup_tst_a<-blup_tst$a_id[,1]

#Comparison
#Check the names
names(y_tst)<-rownames(A)[wheat.sets==fold]
blup_tst_m<-blup_tst_m[match(names(y_tst),names(blup_tst_m))]
blup_tst_a<-blup_tst_a[match(names(y_tst),names(blup_tst_a))]

yHat_tst<-fixef(fm3)[1] + blup_tst_m + blup_tst_a
points(y_tst,yHat_tst,col="red",pch=19)

#Correlation in testing set
cor(y_tst,yHat_tst)




[Package lme4GS version 0.1 Index]