Greg.em: Fitting precision regression models

In berryni/mDINGO: DINGO: Differential Network Analysis in Genomics

Description

This function fits the covariance regression model by Hoff and Niu (2012) using EM algorithm with the restriction of diagonal matrix for the noise variance

Usage

Greg.em(formula, data = NULL, R = 1, tol = 1e-10, itmax = 1000, verbose = F)

Arguments

formula	an object of class "formula" used in model.frame function
data	a data frame used in model.frame function
R	rank of the model
tol	a stopping criterion
itmax	maximum number of iteration
verbose	If true, estimation results for each iteration are printed

Value

A	MLE of the baseline covariance matrix
B	MLE of the regression coefficients

Author(s)

Min Jin Ha <[e-mail]>

References

Hoff, P. D. and Niu, X. (2012) A covariance regression model. Statistica Sinica, 22, 729-753.

Examples

library(glasso)
data(gbm)
x = gbm[,1]
Y = as.matrix(gbm[,-1])
p = ncol(Y)
# Estimating inverse covariance matrix using GLasso #
S = cov(Y)
w.upper = which(upper.tri(S))

rhoarray = exp(seq(log(0.001),log(1),length=100))
BIC = rep(0,length(rhoarray))
for (rh in 1:length(rhoarray)) {
    fit.gl1 = glasso(S,rho=rhoarray[rh])
    BIC[rh] = extendedBIC(gamma=0,omegahat=fit.gl1$wi,S=S,n=nrow(Y))
}
rho = rhoarray[which.min(BIC)]
fit.gl2 = glasso(S,rho=rho)
Omega = fit.gl2$wi

# Fitting (Covariance Regression on transformed data)
diag.Omega = diag(Omega)
P = -Omega/diag.Omega
diag(P) = 0

tY = Y 
mdat = apply(tY,2,mean)
sdat = apply(tY,2,sd)
std.tY = t((t(tY) - mdat)/sdat)
smat = diag(sdat)

## rank 1 covariance regression
fit.g = Greg.em(std.tY~x,R=1) 

berryni/mDINGO documentation built on May 24, 2019, 3:04 a.m.