ldstatsHD (version 1.0.1)

wfrl: joint estimation of multiple regression coefficient matrices

Description

wfrl estimates jointly two regression coefficient matrices from multivariate normal distributed datasets using an ADMM based algorithm.

Usage

wfrl(D1, D2, lambda1, lambda2, automLambdas = TRUE, paired = TRUE, 
	 sigmaEstimate = "CRmad", maxiter=30, tol=1e-05, nsubset = 10000,  
	 rho = 1, rho.increment = 1, notOnlyLambda2 = TRUE)

Arguments

D1

list with the response variables. Two matrices in the list corresponding to the response variables of the two populations.

D2

list with the explanatory variables. Two matrices in the list corresponding to the explanatory variables of the two populations.

lambda1

tuning parameter for sparsity in the regression coefficients.

lambda2

tuning parameter for similarity between the regression coefficients in the two populations.

automLambdas

if TRUE the lambda's are estimated automatically with lambda1 and lambda2 being expected false positive rate levels.

paired

if TRUE, observations in D1 and D2 are assumed to be matched (\(n_1\) must be equal to \(n_2\)).

sigmaEstimate

robust method used to estimate the variance of estimated partial correlations: name that uniquely identifies "mad", "IQR" or "CRmad" (default). This measure is used to automatically select the tuning parameter (when automLambdas = TRUE).

maxiter

maximum number of iterations for the ADMM algorithm.

tol

convergence tolerance.

nsubset

maximum number of estimated partial correlation coefficients (chosen randomly) used to select lambda1 and lambda2 automatically.

rho

regularization parameter used to compute matrix inverse by eigen value decomposition (default of 1).

rho.increment

default of 1.

notOnlyLambda2

if FALSE only lambda2 is found automatically.

Value

An object of class wfrl containing the following components:

regCoef

regression coefficients.

path

non-zero structure for the regression coefficients.

diff_value

convergence control.

iters

number of iterations used.

Details

wfrl uses a weighted-fused least squares lasso maximum likelihood estimator by solving: $$ [\hat{\beta}_H,\hat{\beta}_T] = \arg\min\limits_{\beta_H,\beta_T} \left[ \frac{1}{2n}||Y-\beta_HX||^2_2 + \frac{1}{2n}||Q-\beta_TW||^2_2 +P_{\lambda_1,\lambda_2,V}(\beta)\right] $$ with $$ P_{\lambda_1,\lambda_2,V}(\beta) = \lambda_1||\beta_H||_1 + \lambda_1||\beta_T||_1 + \lambda_2||V \circ (\beta_T-\beta_H)||_1. $$

where \(\lambda_1\) is the sparsity tuning parameter, \(\lambda_2\) is the similarity tuning parameter, and \(V = [v_{ij}]\) is a \(p\times p\) matrix to weight \(\lambda_2\) for each coefficient of the differential precision matrix. If datasets are independent (paired = "FALSE"), then it is assumed that \(v_{ij} = 1\) for all pairs \((i,j)\). Otherwise (paired = "TRUE"), weights are estimated in order to account for the dependence structure between datasets in the differential network estimation. An ADMM-type recursive algorithm is used to solve the optimization problem.

See details in wfgl for transforming the selection problem of the tuning parameters \(\lambda_1\) and \(\lambda_2\).

References

Danaher, P., P. Wang, and D. Witten (2014). The joint graphical lasso for inverse covariance estimation across multiple classes. Journal of the Royal Statistical Society: Series B (Statistical Methodology) (2006), 1-20.

Boyd, S. (2010). Distributed Optimization and Statistical Learning via the Alternating Direction Method of Multipliers. Foundations and Trends in Machine Learning 3(1), 1-122.

See Also

plot.wfrl for graphical representation. wfgl for joint partial correlation estimation.

Examples

# NOT RUN {
	
# example to use of wfrl
N	<- 200
EX2 <- pcorSimulatorJoint(nobs = N, nclusters = 3, nnodesxcluster = c(30, 30,30), 
                          pattern = "pow", diffType = "cluster", dataDepend = "diag", 
                          low.strength = 0.5, sup.strength = 0.9, pdiff = 0.5, nhubs = 5, 
                          degree.hubs = 20,  nOtherEdges = 30, alpha = 2.3, plus = 0, 
                          prob = 0.05, perturb.clust = 0.2, mu = 0, diagCCtype = "dicot", 
                          diagNZ.strength = 0.6, mixProb = 0.5, probSign = 0.7,  
                          exactZeroTh = 0.05)
					 
P           <- EX2$P
q           <- 50 
BETA1       <- array(0, dim = c(P, q))
diag(BETA1) <- rep(0.35,q)
BETA2       <- BETA1
diag(BETA2)[c(1:floor(q/2))] <- 0
sigma2      <- 1.3
Q           <- scale(EX2$D1)
W           <- scale(EX2$D2)
X      	    <- Q%*%BETA1 + mvrnorm(N,rep(0,q),diag(rep(sigma2,q)))
Y      	    <- W%*%BETA2 + mvrnorm(N,rep(0,q),diag(rep(sigma2,q)))
D1     	    <- list(scale(X), scale(Y))
D2     	    <- list(scale(Q), scale(W))
## not run
#wfrl1       <- wfrl(D1, D2, lambda1 = 0.05, lambda2 = 0.05, automLambdas = TRUE, paired = FALSE, 
#                   sigmaEstimate = "CRmad", maxiter = 30, tol = 1e-05, nsubset = 10000, rho = 1, 
#                   rho.increment = 1, notOnlyLambda2 = TRUE)
#print(wfrl1)

# }