Merge branch 'master' into logistic

Author: jonathan-taylor

selectiveInference/R/RcppExports.R

@@ -34,9 +34,6 @@ randomizedLasso = function(X,
         noise_scale = 0.5 * sd(y) * sqrt(mean_diag)
     }
 
-    print(paste("ridge term", ridge_term))
-    print(paste("noise scale", noise_scale))
-    
     noise_type = match.arg(noise_type)
 
     if (noise_scale > 0) {
@@ -101,8 +98,8 @@ randomizedLasso = function(X,
     observed_unpen = result$soln[unpenalized]
     observed_subgrad = -n*result$gradient[inactive]
 
-    if (length(which(abs(observed_subgrad)>lam[1]))){
-      print("subgradient eq not satisfied")
+    if (sum(abs(observed_subgrad)>lam*(1.001)) > 0){
+      stop("subgradient eq not satisfied")
     }
 
     observed_opt_state = c(observed_unpen, observed_scalings, observed_subgrad)
@@ -355,15 +352,16 @@ conditional_density = function(noise_scale, lasso_soln) {
 randomizedLassoInf = function(X, 
                               y, 
                               lam, 
-                              family="gaussian",
-                              sampler="A",
+                              family=c("gaussian", "logistic"),
+                              sampler=c("norejection", "adaptMCMC"),
                               sigma=NULL, 
                               noise_scale=NULL, 
                               ridge_term=NULL, 
                               condition_subgrad=TRUE, 
                               level=0.9,
-			                        nsample=10000,
-			                        burnin=2000,
+                              sampler=c("norejection", "adaptMCMC"),
+                              nsample=10000,
+                              burnin=2000,
                               max_iter=100,        # how many iterations for each optimization problem
                               kkt_tol=1.e-4,       # tolerance for the KKT conditions
                               parameter_tol=1.e-8, # tolerance for relative convergence of parameter
@@ -376,6 +374,8 @@ randomizedLassoInf = function(X,
   n = nrow(X)
   p = ncol(X)
 
+  family = match.arg(family)
+  
   lasso_soln = randomizedLasso(X, 
                                y, 
                                lam, 
@@ -392,7 +392,7 @@ randomizedLassoInf = function(X,
 
   active_set = lasso_soln$active_set
   nactive = length(active_set)
-  print(paste("nactive", nactive))
+
   if (nactive==0){
     return (list(active_set=active_set, pvalues=c(), ci=c()))
   }
@@ -418,10 +418,12 @@ randomizedLassoInf = function(X,
 
   ndim = length(lasso_soln$observed_opt_state)
 
-  if (sampler =="R"){
+  sampler = match.arg(sampler)
+
+  if (sampler == "adaptMCMC"){
     S = sample_opt_variables(lasso_soln, jump_scale=rep(1/sqrt(n), ndim), nsample=nsample)
     opt_samples = as.matrix(S$samples[(burnin+1):nsample,,drop=FALSE])
-  } else if (sampler == "A"){
+  } else if (sampler == "norejection") {
     opt_samples = gaussian_sampler(noise_scale, 
                                  lasso_soln$observed_opt_state, 
                                  cur_opt_transform$linear_term,

selectiveInference/R/funs.randomized.R

@@ -1,3 +1,4 @@
+# A no-rejection MCMC algorithm Jelena and Amir have been working on
 
 log_concave_sampler = function(negative_log_density, 
                                grad_negative_log_density, 

selectiveInference/R/sampler.R renamed to selectiveInference/R/funs.sampler.R

@@ -19,6 +19,7 @@ randomizedLassoInf(X,
                    ridge_term=NULL, 
                    condition_subgrad=TRUE, 
                    level=0.9,
+                   sampler=c("norejection", "adaptMCMC"),
                    nsample=10000,
                    burnin=2000,
                    max_iter=100,       
@@ -71,6 +72,10 @@ Default is TRUE.
 \item{level}{
 Level for confidence intervals.
 }
+\item{sampler}{
+Which sampler to use -- default is a no-rejection sampler. Otherwise
+use MCMC from the adaptMCMC package.
+}
 \item{nsample}{
 Number of samples of optimization variables to sample.
 }