HSDiC_ADMM: Homogeneity Detection Incorporating Prior Constraint Information by ADMM

Description

simultaneous homogeneity detection and variable selection incorporating prior constraint by ADMM algorithm. The problem turn to solving quadratic programming problems of the form min(-d^T b + 1/2 b^T D b) with the constraints A^T b >= b_0. The penalty is the pairwise fusion with p(p-1)/2 number of penalties.

Usage

HSDiC_ADMM(X, Y, A.eq, A.ge, A.lbs, A.ubs, b.eq, b.ge, b.lbs, b.ubs,
  penalty = c("MCP", "SCAD", "adlasso", "lasso"), lambda2,
  admmScale1 = 1/nrow(X), admmScale2 = 1, admmAbsTol = 1e-04,
  admmRelTol = 1e-04, nADMM = 2000, admmVaryScale = FALSE)

Arguments

n-by-p design matrix.

n-by-1 response matrix.

A.eq

equality constraint matrix.

A.ge

inequality constraint matrix.

A.lbs

low-bounds matrix on variables, see examples.

A.ubs

upper-bounds matrix on variables, see examples.

b.eq

equality constraint vector.

b.ge

inequality constraint vector.

b.lbs

low-bounds on variables, see details.

b.ubs

upper-bounds on variables, see details.

penalty

The penalty to be applied to the model. Either "lasso" (the default), "SCAD", or "MCP".

lambda2

penalty tuning parameter for thresholding function.

admmScale1

first ADMM scale parameter, 1/nrow(X) is default.

admmScale2

second ADMM scale parameter, 1 is default.

admmAbsTol

absolute tolerance for ADMM, 1e-04 is default.

admmRelTol

relative tolerance for ADMM, 1e-04 is default.

nADMM

maximum number of iterations for ADMM, 2000 is default.

admmVaryScale

dynamically chance the ADMM scale parameter, FALSE is default

Value

betahat

solution vector.

stats.ADMM_inters

number of iterations.

References

'Pairwise Fusion Approach Incorporating Prior Constraint Information' by Yaguang Li

Examples

Run this code

# NOT RUN {
## data generation
set.seed(111)
n=100
p=50
r <- 1 #0.5, 0.8, 1

beta <- r*c(sample(rep(1:2, each = 10)), rep(0,10), -sample(rep(1:2, each = 10)) )
X <- matrix(rnorm(n*p),nrow = n)
sigma = 1
Y <- X %*% beta + sigma * rnorm(n, 0, 1)


# equalities
A.eq <- rbind(rep(1,p))
b.eq <- c(0)

# inequalities
A.ge <- diag( c(rep(1,30), rep(-1,20)) )
b.ge <- rep(0,p)

# lower-bounds
A.lbs <- diag(1, p)
b.lbs <- rep(-2, p)

# upper-bounds on variables
A.ubs <- diag(-1, p)
b.ubs <- rep(-2, p)

ptm <- proc.time()
fit <- HSDiC_ADMM(X, Y, A.eq, A.ge, A.lbs, A.ubs, b.eq, b.ge, b.lbs, b.ubs,
                 penalty = "adlasso", lambda2 = 0.8, admmScale2 = 1)
proc.time() - ptm

## table(round(fit$beta,1))

plot(beta, type="p", pch = 20, cex = 1)
points(fit$beta, col = 3)


# }