Allow bounds on params

.github/workflows/R-CMD-check.yaml

@@ -8,27 +8,16 @@ name: R-CMD-check
 
 jobs:
   R-CMD-check:
-    runs-on: ${{ matrix.config.os }}
-    name: ${{ matrix.config.os }} (R ${{ matrix.config.r }})
-
-    strategy:
-      fail-fast: false
-      matrix:
-        config:
-          - {os: macos-latest,   r: 'release'}
-          - {os: windows-latest, r: 'release'}
-          - {os: ubuntu-latest,  r: 'release'}
-          - {os: ubuntu-latest,  r: 'devel'}
+    runs-on: ubuntu-latest
+    env:
+      GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
 
     steps:
       - uses: actions/checkout@v4
 
-      - uses: r-lib/actions/setup-pandoc@v2
-
       - uses: r-lib/actions/setup-r@v2
         with:
-          r-version: ${{ matrix.config.r }}
-          use-public-rspm: true
+          r-version: 'release'
 
       - uses: r-lib/actions/setup-r-dependencies@v2
         with:
@@ -37,4 +26,5 @@ jobs:
 
       - uses: r-lib/actions/check-r-package@v2
         with:
+          args: 'c("--no-manual", "--as-cran")'
           upload-snapshots: true

.gitignore

@@ -24,3 +24,7 @@ tests/testthat/_snaps/
 # devtools / roxygen
 /doc/
 /Meta/
+
+# local R session config (machine-specific)
+.Rprofile
+.Renviron

DESCRIPTION

@@ -19,5 +19,4 @@ Imports:
     parallel
 Suggests:
     testthat (>= 3.0.0)
-Roxygen: list(markdown = TRUE)
-RoxygenNote: 7.1.2
+RoxygenNote: 7.3.2

R/GetAlgoParams.R

@@ -18,6 +18,8 @@
 #' @param give_up_init An integer for how many failed initialization attempts before stopping the optimization routine. 100 is the default value.
 #' @param stop_check An integer for how often to check the convergence criterion. The default is 10 iterations.
 #' @param stop_tol A convergence metric must be less than value to be labeled as converged. The default is 1e-4.
+#' @param lower A numeric scalar or n_params-dimensional vector specifying lower bounds for each parameter. Default is -Inf (no lower bound).
+#' @param upper A numeric scalar or n_params-dimensional vector specifying upper bounds for each parameter. Default is Inf (no upper bound).
 #' @param converge_crit A string denoting the convergence metric used, valid metrics are 'stdev' (standard deviation of population weight in the last stop_check iterations) and 'percent' (percent improvement in median particle weight in the last stop_check iterations). 'stdev' is the default.
 #' @return A list of control parameters for the optim_SQGDE function.
 #' @export
@@ -27,6 +29,8 @@ GetAlgoParams = function(n_params,
                          n_iter = 1000,
                          init_sd = 0.01,
                          init_center = 0,
+                         lower = -Inf,
+                         upper = Inf,
                          n_cores_use = 1,
                          step_size = NULL,
                          jitter_size = 1e-6,
@@ -95,6 +99,28 @@ GetAlgoParams = function(n_params,
     stop('ERROR: init_center vector length must be 1 or n_params')
   }
 
+  # lower
+  lower = as.numeric(lower)
+  if(any(is.nan(lower))){
+    stop('ERROR: lower contains NaN')
+  } else if(!(length(lower) == 1 | length(lower) == n_params)){
+    stop('ERROR: lower vector length must be 1 or n_params')
+  }
+  if(length(lower) == 1) lower = rep(lower, n_params)
+
+  # upper
+  upper = as.numeric(upper)
+  if(any(is.nan(upper))){
+    stop('ERROR: upper contains NaN')
+  } else if(!(length(upper) == 1 | length(upper) == n_params)){
+    stop('ERROR: upper vector length must be 1 or n_params')
+  }
+  if(length(upper) == 1) upper = rep(upper, n_params)
+
+  if(any(lower >= upper)){
+    stop('ERROR: lower must be strictly less than upper for all parameters')
+  }
+
   # n_cores_use
   ### assign NULL value default
   if(is.null(n_cores_use)){
@@ -264,6 +290,8 @@ GetAlgoParams = function(n_params,
              'n_iter' = n_iter,
              'init_sd' = init_sd,
              'init_center' = init_center,
+             'lower' = lower,
+             'upper' = upper,
              'n_cores_use' = n_cores_use,
              'step_size' = step_size,
              'crossover_rate' = crossover_rate,

R/SQG_DE_bin_1_best.R

@@ -17,8 +17,10 @@ SQG_DE_bin_1_best=function(pmem_index,
                            objFun,
                            step_size = .8,
                            jitter_size = 1e-6,
-                           n_particles, 
+                           n_particles,
                            crossover_rate = 1,
+                           lower = -Inf,
+                           upper = Inf,
                            n_diff, ... ){
 
   # get statistics about particle
@@ -72,6 +74,7 @@ SQG_DE_bin_1_best=function(pmem_index,
       step_size*psi*(grad_approx) + # move in the direction against the gradient
       stats::runif(len_param_use,-jitter_size,jitter_size) # a little noise
   }
+  params_use = pmax(lower, pmin(upper, params_use))
   params_use = matrix(params_use,1,len_param_use)
 
   weight_proposal = NA

R/SQG_DE_bin_1_curr.R

@@ -18,6 +18,8 @@ SQG_DE_bin_1_curr=function(pmem_index,
                            jitter_size = 1e-6,
                            n_particles,
                            crossover_rate = 1,
+                           lower = -Inf,
+                           upper = Inf,
                            n_diff, ... ){
 
   # sample parents
@@ -71,6 +73,7 @@ SQG_DE_bin_1_curr=function(pmem_index,
       step_size*psi*(grad_approx) + # move in the direction against the gradient
       stats::runif(len_param_use, -jitter_size, jitter_size) # a little noise
   }
+  params_use = pmax(lower, pmin(upper, params_use))
   params_use = matrix(params_use,1,len_param_use)
 
   weight_proposal = NA

R/SQG_DE_bin_1_rand.R

@@ -19,6 +19,8 @@ SQG_DE_bin_1_rand=function(pmem_index,
                            jitter_size=1e-6,
                            n_particles,
                            crossover_rate=1,
+                           lower=-Inf,
+                           upper=Inf,
                            n_diff, ... ){
 
 
@@ -72,6 +74,7 @@ SQG_DE_bin_1_rand=function(pmem_index,
       step_size*psi*(grad_approx) + # move in the direction against the gradient
       stats::runif(len_param_use,-jitter_size,jitter_size) # a little noise
   }
+  params_use = pmax(lower, pmin(upper, params_use))
   params_use = matrix(params_use,1,len_param_use)
 
   weight_proposal = NA

R/optim_SQGDE.R

@@ -77,9 +77,11 @@ optim_SQGDE = function(ObjFun, control_params = GetAlgoParams(), ...){
   for(pmem_index in 1:control_params$n_particles){
     count = 0 # establish a count variable to avoid infinite run time
     while(weights[1,pmem_index]==Inf) {
-      particles[1, pmem_index, ] = stats::rnorm(control_params$n_params,
-                                                control_params$init_center,
-                                                control_params$init_sd)
+      particles[1, pmem_index, ] = pmax(control_params$lower,
+                                        pmin(control_params$upper,
+                                             stats::rnorm(control_params$n_params,
+                                                          control_params$init_center,
+                                                          control_params$init_sd)))
 
       weights[1, pmem_index] = ObjFun(particles[1, pmem_index, ], ...)
 
@@ -139,6 +141,8 @@ optim_SQGDE = function(ObjFun, control_params = GetAlgoParams(), ...){
                                 jitter_size = control_params$jitter_size,
                                 n_particles = control_params$n_particles,
                                 crossover_rate = control_params$crossover_rate,
+                                lower = control_params$lower,
+                                upper = control_params$upper,
                                 n_diff = control_params$n_diff, ...)),
                   control_params$n_particles,
                   control_params$n_params+1, byrow=TRUE)
@@ -152,6 +156,8 @@ optim_SQGDE = function(ObjFun, control_params = GetAlgoParams(), ...){
                                              jitter_size = control_params$jitter_size,
                                              n_particles = control_params$n_particles,
                                              crossover_rate = control_params$crossover_rate,
+                                             lower = control_params$lower,
+                                             upper = control_params$upper,
                                              n_diff = control_params$n_diff, ...)),
                   control_params$n_particles,
                   control_params$n_params+1, byrow=TRUE)

tests/testthat/test-bounds.R

@@ -0,0 +1,148 @@
+library(graDiEnt)
+
+# ── GetAlgoParams: bounds validation ─────────────────────────────────────────
+
+test_that("GetAlgoParams rejects lower >= upper", {
+  expect_error(GetAlgoParams(n_params = 2, lower = 1,      upper = 0),
+               "lower must be strictly less than upper")
+  expect_error(GetAlgoParams(n_params = 2, lower = c(0, 5), upper = c(1, 5)),
+               "lower must be strictly less than upper")
+})
+
+test_that("GetAlgoParams rejects wrong-length bounds", {
+  expect_error(GetAlgoParams(n_params = 2, lower = c(0, 0, 0)),
+               "lower vector length must be 1 or n_params")
+  expect_error(GetAlgoParams(n_params = 2, upper = c(1, 1, 1)),
+               "upper vector length must be 1 or n_params")
+})
+
+test_that("GetAlgoParams stores bounds as length-n_params vectors", {
+  cp <- GetAlgoParams(n_params = 3, lower = -1, upper = 2)
+  expect_length(cp$lower, 3)
+  expect_length(cp$upper, 3)
+  expect_equal(cp$lower, c(-1, -1, -1))
+  expect_equal(cp$upper, c(2,  2,  2))
+})
+
+test_that("GetAlgoParams accepts per-parameter bounds", {
+  cp <- GetAlgoParams(n_params = 3, lower = c(-1, 0, -5), upper = c(1, 3, 5))
+  expect_equal(cp$lower, c(-1,  0, -5))
+  expect_equal(cp$upper, c( 1,  3,  5))
+})
+
+# ── optim_SQGDE: solution stays within symmetric scalar bounds ────────────────
+
+test_that("solution stays within scalar bounds when optimum is inside", {
+  set.seed(42)
+  out <- suppressMessages(
+    optim_SQGDE(
+      ObjFun = function(x) sum((x - c(1, 2))^2),
+      control_params = GetAlgoParams(
+        n_params    = 2,
+        n_iter      = 400,
+        n_particles = 10,
+        n_diff      = 2,
+        init_sd     = 1,
+        init_center = 0,
+        lower       = -5,
+        upper       =  5
+      )
+    )
+  )
+
+  expect_true(all(out$solution >= -5))
+  expect_true(all(out$solution <=  5))
+  expect_lt(max(abs(out$solution - c(1, 2))), 0.15)
+})
+
+# ── optim_SQGDE: per-parameter asymmetric bounds ─────────────────────────────
+
+test_that("solution stays within per-parameter asymmetric bounds", {
+  set.seed(7)
+  true_opt <- c(0.5, 2.5, -0.5)
+  out <- suppressMessages(
+    optim_SQGDE(
+      ObjFun = function(x) sum((x - true_opt)^2),
+      control_params = GetAlgoParams(
+        n_params    = 3,
+        n_iter      = 500,
+        n_particles = 12,
+        n_diff      = 2,
+        init_sd     = 0.3,
+        init_center = c(0.5, 2, -0.5),
+        lower       = c( 0,  1, -1),
+        upper       = c( 1,  3,  0)
+      )
+    )
+  )
+
+  lower <- c(0, 1, -1)
+  upper <- c(1, 3,  0)
+  expect_true(all(out$solution >= lower))
+  expect_true(all(out$solution <= upper))
+  expect_lt(max(abs(out$solution - true_opt)), 0.15)
+})
+
+# ── optim_SQGDE: optimum outside bounds → clips to boundary ──────────────────
+
+test_that("solution clips to boundary when true optimum is outside bounds", {
+  set.seed(99)
+  out <- suppressMessages(
+    optim_SQGDE(
+      ObjFun = function(x) sum((x - c(-10, -10))^2),
+      control_params = GetAlgoParams(
+        n_params    = 2,
+        n_iter      = 400,
+        n_particles = 10,
+        n_diff      = 2,
+        init_sd     = 0.5,
+        init_center = 1,
+        lower       = 0,
+        upper       = 2
+      )
+    )
+  )
+
+  expect_true(all(out$solution >= 0))
+  expect_true(all(out$solution <= 2))
+  # constrained optimum is at lower boundary (0, 0)
+  expect_lt(max(abs(out$solution - c(0, 0))), 0.1)
+})
+
+# ── optim_SQGDE: unbounded by default (no behaviour change) ──────────────────
+
+test_that("default bounds (-Inf / Inf) do not alter unbounded behaviour", {
+  obj <- function(x) sum((x - c(-2, 3))^2)
+
+  set.seed(5)
+  out_nobound <- suppressMessages(
+    optim_SQGDE(
+      ObjFun = obj,
+      control_params = GetAlgoParams(
+        n_params    = 2,
+        n_iter      = 400,
+        n_particles = 10,
+        n_diff      = 2,
+        init_sd     = 1
+      )
+    )
+  )
+
+  set.seed(5)
+  out_explicitinf <- suppressMessages(
+    optim_SQGDE(
+      ObjFun = obj,
+      control_params = GetAlgoParams(
+        n_params    = 2,
+        n_iter      = 400,
+        n_particles = 10,
+        n_diff      = 2,
+        init_sd     = 1,
+        lower       = -Inf,
+        upper       =  Inf
+      )
+    )
+  )
+
+  expect_equal(out_nobound$solution, out_explicitinf$solution)
+})