tests/testthat/test_optimization_methods.R
In agisga/grpSLOPE: Group Sorted L1 Penalized Estimation
library(grpSLOPE)
context("proxGroupSortedL1()")
y <- 1:10
grp <- c(1,1,1,2,2,2,3,3,3,4)
sol <- c(0, 0, 0, 0.353261553, 0.441576942, 0.529892330,
1.974292890, 2.256334731, 2.538376572, 1.0000000001)
lam <- c(10, 9, 8, 7)
test_that("the prox is evaluated correctly when the groups are consequtive blocks", {
x <- proxGroupSortedL1(y = y, group = grp, lambda = lam)
expect_equal(x, sol)
})
test_that("the prox is evaluated correctly when the groups are not consequtive blocks", {
ord <- sample(y, 10)
x <- proxGroupSortedL1(y = y[ord], group = grp[ord], lambda = lam)
expect_equal(x, sol[ord], tolerance=1e-6)
})
#---
context("proximalGradientSolverGroupSLOPE()")
# set.seed(1)
# A <- matrix(runif(100, 0, 1), 10, 10)
A.vec <- c(0.26550866314209997653961,0.37212389963679015636444,0.57285336335189640522003,0.90820778999477624893188,0.20168193103745579719543,0.89838968496769666671753,0.94467526860535144805908,0.66079779248684644699097,0.62911404389888048171997,0.06178627046756446361542,0.20597457489930093288422,0.17655675252899527549744,0.68702284665778279304504,0.38410371821373701095581,0.76984141999855637550354,0.49769924208521842956543,0.71761850826442241668701,0.99190609483048319816589,0.38003517943434417247772,0.77744522131979465484619,0.93470523110590875148773,0.21214252128265798091888,0.65167376608587801456451,0.12555509596131742000580,0.26722066872753202915192,0.38611409254372119903564,0.01339033315889537334442,0.38238795707002282142639,0.86969084572046995162964,0.34034899668768048286438,0.48208011547103524208069,0.59956582542508840560913,0.49354130704887211322784,0.18621760141104459762573,0.82737331860698759555817,0.66846673819236457347870,0.79423986072652041912079,0.10794362588785588741302,0.72371094604022800922394,0.41127442964352667331696,0.82094629411585628986359,0.64706019381992518901825,0.78293276228941977024078,0.55303631164133548736572,0.52971958019770681858063,0.78935623168945312500000,0.02333120233379304409027,0.47723006503656506538391,0.73231373867020010948181,0.69273155648261308670044,0.47761962213553488254547,0.86120947683230042457581,0.43809710722416639328003,0.24479727703146636486053,0.07067904714494943618774,0.09946616017259657382965,0.31627170718275010585785,0.51863426319323480129242,0.66200507641769945621490,0.40683018718846142292023,0.91287592425942420959473,0.29360337276011705398560,0.45906572625972330570221,0.33239467418752610683441,0.65087046707049012184143,0.25801678071729838848114,0.47854524827562272548676,0.76631067064590752124786,0.08424691436812281608582,0.87532133003696799278259,0.33907293784432113170624,0.83944035018794238567352,0.34668348915874958038330,0.33377493079751729965210,0.47635124507360160350800,0.89219833584502339363098,0.86433947063051164150238,0.38998954347334802150726,0.77732069883495569229126,0.96061799721792340278625,0.43465948477387428283691,0.71251467871479690074921,0.39999436889775097370148,0.32535215187817811965942,0.75708714802749454975128,0.20269225514493882656097,0.71112122246995568275452,0.12169192102737724781036,0.24548851395957171916962,0.14330437942408025264740,0.23962941509671509265900,0.05893437727354466915131,0.64228825853206217288971,0.87626921269111335277557,0.77891467744484543800354,0.79730882588773965835571,0.45527445361949503421783,0.41008408204652369022369,0.81087024277076125144958,0.60493329027667641639709)
A <- matrix(A.vec, 10, 10)
grp <- c(0, 0, 1, 1, 2, 2, 2, 2, 2, 3)
wt <- c(2, 2, 2, 2, 5, 5, 5, 5, 5, 1)
x <- c(0, 0, 5, 1, 0, 0, 0, 1, 0, 3)
y <- A %*% x
lam <- 0.1 * (10:7)
sol <- c(0,0,3.856002988,2.080742942,0,0,0,0,0,3.512828045)
test_that("it works correctly when the groups are consequtive blocks", {
result <- proximalGradientSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_is(result$L.iter, "numeric")
expect_true(result$L.iter >= result$iter)
})
test_that("it works correctly when the groups are not consequtive blocks", {
ord <- sample(1:10, 10)
result <- proximalGradientSolverGroupSLOPE(y=y, A=A[ , ord], group=grp[ord], wt=wt[ord],
lambda=lam, dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol[ord]), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_is(result$L.iter, "numeric")
expect_true(result$L.iter >= result$iter)
})
test_that("solution is computed correctly when each group is a singleton", {
grp <- 1:10
lam <- 0.1*(10:1)
wt <- rep(1, 10)
sol <- c(0.2036292051,0.2036292051,3.8097062299,0.8573736336,0.8573736336,
0.2036292051,0.2036292051,0.7057561857,0.2036292051,2.3940891644)
result <- proximalGradientSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-6)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_true(is.null(result$L.iter))
})
#---
context("admmSolverGroupSLOPE()")
# set.seed(1)
# A <- matrix(runif(100, 0, 1), 10, 10)
A.vec <- c(0.26550866314209997653961,0.37212389963679015636444,0.57285336335189640522003,0.90820778999477624893188,0.20168193103745579719543,0.89838968496769666671753,0.94467526860535144805908,0.66079779248684644699097,0.62911404389888048171997,0.06178627046756446361542,0.20597457489930093288422,0.17655675252899527549744,0.68702284665778279304504,0.38410371821373701095581,0.76984141999855637550354,0.49769924208521842956543,0.71761850826442241668701,0.99190609483048319816589,0.38003517943434417247772,0.77744522131979465484619,0.93470523110590875148773,0.21214252128265798091888,0.65167376608587801456451,0.12555509596131742000580,0.26722066872753202915192,0.38611409254372119903564,0.01339033315889537334442,0.38238795707002282142639,0.86969084572046995162964,0.34034899668768048286438,0.48208011547103524208069,0.59956582542508840560913,0.49354130704887211322784,0.18621760141104459762573,0.82737331860698759555817,0.66846673819236457347870,0.79423986072652041912079,0.10794362588785588741302,0.72371094604022800922394,0.41127442964352667331696,0.82094629411585628986359,0.64706019381992518901825,0.78293276228941977024078,0.55303631164133548736572,0.52971958019770681858063,0.78935623168945312500000,0.02333120233379304409027,0.47723006503656506538391,0.73231373867020010948181,0.69273155648261308670044,0.47761962213553488254547,0.86120947683230042457581,0.43809710722416639328003,0.24479727703146636486053,0.07067904714494943618774,0.09946616017259657382965,0.31627170718275010585785,0.51863426319323480129242,0.66200507641769945621490,0.40683018718846142292023,0.91287592425942420959473,0.29360337276011705398560,0.45906572625972330570221,0.33239467418752610683441,0.65087046707049012184143,0.25801678071729838848114,0.47854524827562272548676,0.76631067064590752124786,0.08424691436812281608582,0.87532133003696799278259,0.33907293784432113170624,0.83944035018794238567352,0.34668348915874958038330,0.33377493079751729965210,0.47635124507360160350800,0.89219833584502339363098,0.86433947063051164150238,0.38998954347334802150726,0.77732069883495569229126,0.96061799721792340278625,0.43465948477387428283691,0.71251467871479690074921,0.39999436889775097370148,0.32535215187817811965942,0.75708714802749454975128,0.20269225514493882656097,0.71112122246995568275452,0.12169192102737724781036,0.24548851395957171916962,0.14330437942408025264740,0.23962941509671509265900,0.05893437727354466915131,0.64228825853206217288971,0.87626921269111335277557,0.77891467744484543800354,0.79730882588773965835571,0.45527445361949503421783,0.41008408204652369022369,0.81087024277076125144958,0.60493329027667641639709)
A <- matrix(A.vec, 10, 10)
grp <- c(0, 0, 1, 1, 2, 2, 2, 2, 2, 3)
wt <- c(2, 2, 2, 2, 5, 5, 5, 5, 5, 1)
x <- c(0, 0, 5, 1, 0, 0, 0, 1, 0, 3)
y <- A %*% x
lam <- 0.1 * (10:7)
sol <- c(0,0,3.856002988,2.080742942,0,0,0,0,0,3.512828045)
test_that("it works correctly when the groups are consequtive blocks", {
result <- admmSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam, rho=1,
absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
test_that("it works correctly when the groups are not consequtive blocks", {
ord <- sample(1:10, 10)
result <- admmSolverGroupSLOPE(y=y, A=A[ , ord], group=grp[ord], wt=wt[ord],
lambda=lam, rho=1, absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol[ord]), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
test_that("solution is computed correctly when each group is a singleton", {
grp <- 1:10
lam <- 0.1*(10:1)
wt <- rep(1, 10)
sol <- c(0.2036292051,0.2036292051,3.8097062299,0.8573736336,0.8573736336,
0.2036292051,0.2036292051,0.7057561857,0.2036292051,2.3940891644)
result <- admmSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
rho=1, absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
#--- comparison of FISTA vs ADMM for high and low-dimensional problems
generate_data <- function(n = NA, p = NA, n_grp = NA, seed = 2018) {
set.seed(seed)
A <- matrix(rnorm(n * p), n, p)
# group membership per predictor
grp <- sample(x = 1:n_grp, size = p, replace = TRUE)
# weights: sqrt of the group size as weight for each group
wt_per_grp <- sqrt(c(table(grp)))
wt_per_predictor <- wt_per_grp[as.character(grp)]
# a sparse vector of coefficients per predictor (i.e., the true solution)
grp_signif <- 1:10
ind_signif <- which(grp %in% grp_signif)
x <- rep(0, p)
x[ind_signif] <- rnorm(length(ind_signif), mean = 5, sd = 1)
# response variable
y <- A %*% x + rnorm(n)
# normalize
col_norms <- apply(A, 2, function(a) sqrt(sum((a - mean(a))^2)))
A <- scale(A, center = TRUE, scale = col_norms)
y <- scale(y, center = TRUE, scale = FALSE)
# regularizing parameters
lambda <- lambdaGroupSLOPE(method = "corrected", fdr = 0.1, n.obs = n,
group = grp, wt = wt_per_grp)
return(list("y" = y, "A" = A, "grp" = grp, "lambda" = lambda,
"wt_per_predictor" = wt_per_predictor))
}
# high-dimensional problem
test_that("ADMM agrees with FISTA in a high-dimensional problem", {
n <- 200
p <- 400
n_grp <- 100
dat <- generate_data(n = n, p = p, n_grp = n_grp)
# optimization
result_admm <- admmSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
rho = 1,
absolute.tol = 1e-6,
relative.tol=1e-6,
verbose = FALSE)
result_fista <- proximalGradientSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
dual.gap.tol = 1e-6,
infeas.tol = 1e-6,
verbose = FALSE)
expect_equal(result_admm$x, result_fista$x, tolerance=1e-6)
expect_identical(result_admm$status, 1)
expect_identical(result_fista$status, 1)
expect_is(result_admm$iter, "numeric")
expect_is(result_fista$iter, "numeric")
expect_true(result_admm$iter > 1)
expect_true(result_fista$iter > 1)
})
# low-dimensional problem
test_that("ADMM agrees with FISTA in a low-dimensional problem", {
n <- 400
p <- 100
n_grp <- 25
dat <- generate_data(n = n, p = p, n_grp = n_grp)
# optimization
result_admm <- admmSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
rho = 1,
absolute.tol = 1e-6,
relative.tol=1e-6,
verbose = FALSE)
result_fista <- proximalGradientSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
dual.gap.tol = 1e-6,
infeas.tol = 1e-6,
verbose = FALSE)
expect_equal(result_admm$x, result_fista$x, tolerance=1e-6)
expect_identical(result_admm$status, 1)
expect_identical(result_fista$status, 1)
expect_is(result_admm$iter, "numeric")
expect_is(result_fista$iter, "numeric")
expect_true(result_admm$iter > 1)
expect_true(result_fista$iter > 1)
})
agisga/grpSLOPE documentation built on June 3, 2023, 10:16 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
library(grpSLOPE)
context("proxGroupSortedL1()")
y <- 1:10
grp <- c(1,1,1,2,2,2,3,3,3,4)
sol <- c(0, 0, 0, 0.353261553, 0.441576942, 0.529892330,
1.974292890, 2.256334731, 2.538376572, 1.0000000001)
lam <- c(10, 9, 8, 7)
test_that("the prox is evaluated correctly when the groups are consequtive blocks", {
x <- proxGroupSortedL1(y = y, group = grp, lambda = lam)
expect_equal(x, sol)
})
test_that("the prox is evaluated correctly when the groups are not consequtive blocks", {
ord <- sample(y, 10)
x <- proxGroupSortedL1(y = y[ord], group = grp[ord], lambda = lam)
expect_equal(x, sol[ord], tolerance=1e-6)
})
#---
context("proximalGradientSolverGroupSLOPE()")
# set.seed(1)
# A <- matrix(runif(100, 0, 1), 10, 10)
A.vec <- c(0.26550866314209997653961,0.37212389963679015636444,0.57285336335189640522003,0.90820778999477624893188,0.20168193103745579719543,0.89838968496769666671753,0.94467526860535144805908,0.66079779248684644699097,0.62911404389888048171997,0.06178627046756446361542,0.20597457489930093288422,0.17655675252899527549744,0.68702284665778279304504,0.38410371821373701095581,0.76984141999855637550354,0.49769924208521842956543,0.71761850826442241668701,0.99190609483048319816589,0.38003517943434417247772,0.77744522131979465484619,0.93470523110590875148773,0.21214252128265798091888,0.65167376608587801456451,0.12555509596131742000580,0.26722066872753202915192,0.38611409254372119903564,0.01339033315889537334442,0.38238795707002282142639,0.86969084572046995162964,0.34034899668768048286438,0.48208011547103524208069,0.59956582542508840560913,0.49354130704887211322784,0.18621760141104459762573,0.82737331860698759555817,0.66846673819236457347870,0.79423986072652041912079,0.10794362588785588741302,0.72371094604022800922394,0.41127442964352667331696,0.82094629411585628986359,0.64706019381992518901825,0.78293276228941977024078,0.55303631164133548736572,0.52971958019770681858063,0.78935623168945312500000,0.02333120233379304409027,0.47723006503656506538391,0.73231373867020010948181,0.69273155648261308670044,0.47761962213553488254547,0.86120947683230042457581,0.43809710722416639328003,0.24479727703146636486053,0.07067904714494943618774,0.09946616017259657382965,0.31627170718275010585785,0.51863426319323480129242,0.66200507641769945621490,0.40683018718846142292023,0.91287592425942420959473,0.29360337276011705398560,0.45906572625972330570221,0.33239467418752610683441,0.65087046707049012184143,0.25801678071729838848114,0.47854524827562272548676,0.76631067064590752124786,0.08424691436812281608582,0.87532133003696799278259,0.33907293784432113170624,0.83944035018794238567352,0.34668348915874958038330,0.33377493079751729965210,0.47635124507360160350800,0.89219833584502339363098,0.86433947063051164150238,0.38998954347334802150726,0.77732069883495569229126,0.96061799721792340278625,0.43465948477387428283691,0.71251467871479690074921,0.39999436889775097370148,0.32535215187817811965942,0.75708714802749454975128,0.20269225514493882656097,0.71112122246995568275452,0.12169192102737724781036,0.24548851395957171916962,0.14330437942408025264740,0.23962941509671509265900,0.05893437727354466915131,0.64228825853206217288971,0.87626921269111335277557,0.77891467744484543800354,0.79730882588773965835571,0.45527445361949503421783,0.41008408204652369022369,0.81087024277076125144958,0.60493329027667641639709)
A <- matrix(A.vec, 10, 10)
grp <- c(0, 0, 1, 1, 2, 2, 2, 2, 2, 3)
wt <- c(2, 2, 2, 2, 5, 5, 5, 5, 5, 1)
x <- c(0, 0, 5, 1, 0, 0, 0, 1, 0, 3)
y <- A %*% x
lam <- 0.1 * (10:7)
sol <- c(0,0,3.856002988,2.080742942,0,0,0,0,0,3.512828045)
test_that("it works correctly when the groups are consequtive blocks", {
result <- proximalGradientSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_is(result$L.iter, "numeric")
expect_true(result$L.iter >= result$iter)
})
test_that("it works correctly when the groups are not consequtive blocks", {
ord <- sample(1:10, 10)
result <- proximalGradientSolverGroupSLOPE(y=y, A=A[ , ord], group=grp[ord], wt=wt[ord],
lambda=lam, dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol[ord]), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_is(result$L.iter, "numeric")
expect_true(result$L.iter >= result$iter)
})
test_that("solution is computed correctly when each group is a singleton", {
grp <- 1:10
lam <- 0.1*(10:1)
wt <- rep(1, 10)
sol <- c(0.2036292051,0.2036292051,3.8097062299,0.8573736336,0.8573736336,
0.2036292051,0.2036292051,0.7057561857,0.2036292051,2.3940891644)
result <- proximalGradientSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
dual.gap.tol=1e-6, infeas.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-6)
expect_identical(result$status, 1)
expect_is(result$L, "numeric")
expect_true(result$L > 0)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
expect_true(is.null(result$L.iter))
})
#---
context("admmSolverGroupSLOPE()")
# set.seed(1)
# A <- matrix(runif(100, 0, 1), 10, 10)
A.vec <- c(0.26550866314209997653961,0.37212389963679015636444,0.57285336335189640522003,0.90820778999477624893188,0.20168193103745579719543,0.89838968496769666671753,0.94467526860535144805908,0.66079779248684644699097,0.62911404389888048171997,0.06178627046756446361542,0.20597457489930093288422,0.17655675252899527549744,0.68702284665778279304504,0.38410371821373701095581,0.76984141999855637550354,0.49769924208521842956543,0.71761850826442241668701,0.99190609483048319816589,0.38003517943434417247772,0.77744522131979465484619,0.93470523110590875148773,0.21214252128265798091888,0.65167376608587801456451,0.12555509596131742000580,0.26722066872753202915192,0.38611409254372119903564,0.01339033315889537334442,0.38238795707002282142639,0.86969084572046995162964,0.34034899668768048286438,0.48208011547103524208069,0.59956582542508840560913,0.49354130704887211322784,0.18621760141104459762573,0.82737331860698759555817,0.66846673819236457347870,0.79423986072652041912079,0.10794362588785588741302,0.72371094604022800922394,0.41127442964352667331696,0.82094629411585628986359,0.64706019381992518901825,0.78293276228941977024078,0.55303631164133548736572,0.52971958019770681858063,0.78935623168945312500000,0.02333120233379304409027,0.47723006503656506538391,0.73231373867020010948181,0.69273155648261308670044,0.47761962213553488254547,0.86120947683230042457581,0.43809710722416639328003,0.24479727703146636486053,0.07067904714494943618774,0.09946616017259657382965,0.31627170718275010585785,0.51863426319323480129242,0.66200507641769945621490,0.40683018718846142292023,0.91287592425942420959473,0.29360337276011705398560,0.45906572625972330570221,0.33239467418752610683441,0.65087046707049012184143,0.25801678071729838848114,0.47854524827562272548676,0.76631067064590752124786,0.08424691436812281608582,0.87532133003696799278259,0.33907293784432113170624,0.83944035018794238567352,0.34668348915874958038330,0.33377493079751729965210,0.47635124507360160350800,0.89219833584502339363098,0.86433947063051164150238,0.38998954347334802150726,0.77732069883495569229126,0.96061799721792340278625,0.43465948477387428283691,0.71251467871479690074921,0.39999436889775097370148,0.32535215187817811965942,0.75708714802749454975128,0.20269225514493882656097,0.71112122246995568275452,0.12169192102737724781036,0.24548851395957171916962,0.14330437942408025264740,0.23962941509671509265900,0.05893437727354466915131,0.64228825853206217288971,0.87626921269111335277557,0.77891467744484543800354,0.79730882588773965835571,0.45527445361949503421783,0.41008408204652369022369,0.81087024277076125144958,0.60493329027667641639709)
A <- matrix(A.vec, 10, 10)
grp <- c(0, 0, 1, 1, 2, 2, 2, 2, 2, 3)
wt <- c(2, 2, 2, 2, 5, 5, 5, 5, 5, 1)
x <- c(0, 0, 5, 1, 0, 0, 0, 1, 0, 3)
y <- A %*% x
lam <- 0.1 * (10:7)
sol <- c(0,0,3.856002988,2.080742942,0,0,0,0,0,3.512828045)
test_that("it works correctly when the groups are consequtive blocks", {
result <- admmSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam, rho=1,
absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
test_that("it works correctly when the groups are not consequtive blocks", {
ord <- sample(1:10, 10)
result <- admmSolverGroupSLOPE(y=y, A=A[ , ord], group=grp[ord], wt=wt[ord],
lambda=lam, rho=1, absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol[ord]), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
test_that("solution is computed correctly when each group is a singleton", {
grp <- 1:10
lam <- 0.1*(10:1)
wt <- rep(1, 10)
sol <- c(0.2036292051,0.2036292051,3.8097062299,0.8573736336,0.8573736336,
0.2036292051,0.2036292051,0.7057561857,0.2036292051,2.3940891644)
result <- admmSolverGroupSLOPE(y=y, A=A, group=grp, wt=wt, lambda=lam,
rho=1, absolute.tol=1e-6, relative.tol=1e-6, verbose=FALSE)
expect_equal(result$x, as.matrix(sol), tolerance=1e-4)
expect_identical(result$status, 1)
expect_is(result$iter, "numeric")
expect_true(result$iter > 1)
})
#--- comparison of FISTA vs ADMM for high and low-dimensional problems
generate_data <- function(n = NA, p = NA, n_grp = NA, seed = 2018) {
set.seed(seed)
A <- matrix(rnorm(n * p), n, p)
# group membership per predictor
grp <- sample(x = 1:n_grp, size = p, replace = TRUE)
# weights: sqrt of the group size as weight for each group
wt_per_grp <- sqrt(c(table(grp)))
wt_per_predictor <- wt_per_grp[as.character(grp)]
# a sparse vector of coefficients per predictor (i.e., the true solution)
grp_signif <- 1:10
ind_signif <- which(grp %in% grp_signif)
x <- rep(0, p)
x[ind_signif] <- rnorm(length(ind_signif), mean = 5, sd = 1)
# response variable
y <- A %*% x + rnorm(n)
# normalize
col_norms <- apply(A, 2, function(a) sqrt(sum((a - mean(a))^2)))
A <- scale(A, center = TRUE, scale = col_norms)
y <- scale(y, center = TRUE, scale = FALSE)
# regularizing parameters
lambda <- lambdaGroupSLOPE(method = "corrected", fdr = 0.1, n.obs = n,
group = grp, wt = wt_per_grp)
return(list("y" = y, "A" = A, "grp" = grp, "lambda" = lambda,
"wt_per_predictor" = wt_per_predictor))
}
# high-dimensional problem
test_that("ADMM agrees with FISTA in a high-dimensional problem", {
n <- 200
p <- 400
n_grp <- 100
dat <- generate_data(n = n, p = p, n_grp = n_grp)
# optimization
result_admm <- admmSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
rho = 1,
absolute.tol = 1e-6,
relative.tol=1e-6,
verbose = FALSE)
result_fista <- proximalGradientSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
dual.gap.tol = 1e-6,
infeas.tol = 1e-6,
verbose = FALSE)
expect_equal(result_admm$x, result_fista$x, tolerance=1e-6)
expect_identical(result_admm$status, 1)
expect_identical(result_fista$status, 1)
expect_is(result_admm$iter, "numeric")
expect_is(result_fista$iter, "numeric")
expect_true(result_admm$iter > 1)
expect_true(result_fista$iter > 1)
})
# low-dimensional problem
test_that("ADMM agrees with FISTA in a low-dimensional problem", {
n <- 400
p <- 100
n_grp <- 25
dat <- generate_data(n = n, p = p, n_grp = n_grp)
# optimization
result_admm <- admmSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
rho = 1,
absolute.tol = 1e-6,
relative.tol=1e-6,
verbose = FALSE)
result_fista <- proximalGradientSolverGroupSLOPE(y = dat$y,
A = dat$A,
group = dat$grp,
wt = dat$wt_per_predictor,
lambda = dat$lambda,
dual.gap.tol = 1e-6,
infeas.tol = 1e-6,
verbose = FALSE)
expect_equal(result_admm$x, result_fista$x, tolerance=1e-6)
expect_identical(result_admm$status, 1)
expect_identical(result_fista$status, 1)
expect_is(result_admm$iter, "numeric")
expect_is(result_fista$iter, "numeric")
expect_true(result_admm$iter > 1)
expect_true(result_fista$iter > 1)
})
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.