Removed Rsolnp/csolnp (not ready for primetime). Some fixes for crashes due to RTMB and NaN parameters

Author: alexiosg

DESCRIPTION

@@ -7,7 +7,7 @@ Maintainer: Alexios Galanos <alexios@4dscape.com>
 Description: Unobserved components time series model using the linear innovations state space representation (single source of error) with choice of error distributions and option for dynamic variance. Methods for estimation using automatic differentiation, automatic model selection and ensembling, prediction, filtering, simulation and backtesting. Based on the model described in Hyndman et al (2012) <doi:10.1198/jasa.2011.tm09771>.
 Depends: R (>= 4.1.0), Rcpp (>= 0.12.9), tsmethods (>= 1.0.0)
 LinkingTo: Rcpp (>= 0.12.9), TMB, RcppEigen
-Imports: TMB (>= 1.7.20), methods, tsaux, tsdistributions, zoo, xts, copula, flextable, data.table, sandwich, nloptr, Rsolnp, RTMB, viridisLite, future, future.apply, progressr
+Imports: TMB (>= 1.7.20), methods, tsaux, tsdistributions, zoo, xts, copula, flextable, data.table, sandwich, nloptr, RTMB, viridisLite, future, future.apply, progressr
 Suggests: 
     rmarkdown,
     tstests,

NAMESPACE

@@ -63,7 +63,6 @@ importFrom(RTMB,eigen)
 importFrom(RTMB,matrix)
 importFrom(Rcpp,evalCpp)
 importFrom(Rcpp,loadModule)
-importFrom(Rsolnp,csolnp)
 importFrom(TMB,MakeADFun)
 importFrom(copula,P2p)
 importFrom(copula,normalCopula)

NEWS.md

@@ -1,8 +1,8 @@
 # tsissm 1.0.2
 
-* The SLSQP solver (from nloptr) sometimes will return NaN for the parameters. Added a stop
-clause else it may crash the process (due to RTMB eigen crashing).
-* Added option for solnp solver (csolnp from Rsolnp)
+* The SLSQP solver (from nloptr) sometimes will return NaN for the parameters during
+the optimization phase. Added a stop clause else it may crash the process (due 
+to RTMB eigen crashing).
 * Added a solver control function (issm_control)
 
 # tsissm 1.0.1

R/backtest.R

@@ -20,7 +20,8 @@
 #' @param seed an value specifying if and how the random number generator should 
 #' be initialized (\sQuote{seeded}). Either NULL or an integer that will be used in a 
 #' call to set.seed before simulating the response vectors.
-#' @param solver either \dQuote{nlortr} or \dQuote{solnp}.
+#' @param solver only \dQuote{nlortr} currently supported.
+#' @param control optional control parameters.
 #' @param trace whether to show the progress bar. The user is expected to have
 #' set up appropriate handlers for this using the \dQuote{progressr} package.
 #' @param ... not used.
@@ -57,8 +58,12 @@
 #'
 tsbacktest.tsissm.autospec <- function(object, start = floor(length(object$y)/2), end = length(object$y),
                                    h = 1, estimate_every = 1, rolling = FALSE, weights_scheme = c("AIC","BIC","U"), weights = NULL,
-                                   seed = NULL, solver = "nloptr", trace = FALSE, ...)
+                                   seed = NULL, solver = "nloptr", control = NULL, trace = FALSE, ...)
 {
+    solver <- "nloptr"
+    if (is.null(control)) {
+        control <- issm_control(solver = "nloptr", algorithm = "SLSQP", trace = 0)
+    }
     if (object$top_n > 1) {
         top_n <- object$top_n
         weights_scheme <- match.arg(weights_scheme[1], c("AIC","BIC","U"))
@@ -68,10 +73,10 @@ tsbacktest.tsissm.autospec <- function(object, start = floor(length(object$y)/2)
             weights_scheme <- "U"
         }
         out <- .backtest_ensemble(object, start = start, end = end, h = h, estimate_every = estimate_every, rolling = rolling, weights_scheme = weights_scheme, 
-                                  weights = weights, seed = seed, solver = solver, trace = trace)
+                                  weights = weights, seed = seed, solver = solver, control = control, trace = trace)
     } else {
         out <- .backtest_top(object, start = start, end = end, h = h, estimate_every = estimate_every, rolling = rolling, seed = seed, solver = solver, 
-                             trace = trace)
+                             control = control, trace = trace)
     }
     return(out)
 }
@@ -80,8 +85,12 @@ tsbacktest.tsissm.autospec <- function(object, start = floor(length(object$y)/2)
 #' @rdname tsbacktest
 #' @export
 tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_orig)/2), end = length(object$target$y_orig),
-                                   h = 1, estimate_every = 1, rolling = FALSE, seed = NULL, solver = "nloptr", trace = FALSE, ...)
+                                   h = 1, estimate_every = 1, rolling = FALSE, seed = NULL, solver = "nloptr", control = NULL, trace = FALSE, ...)
 {
+    solver <- "nloptr"
+    if (is.null(control)) {
+        control <- issm_control(solver = "nloptr", algorithm = "SLSQP", trace = 0)
+    }
     parameter <- b <- forecast_dates <- NULL
     data <- xts(object$target$y_orig, object$target$index)
 
@@ -152,7 +161,7 @@ tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_
                                sampling = object$target$sampling, sample_n = object$variance$sample_n, 
                                init_garch =  object$variance$init_variance, garch_order = object$variance$order, 
                                variance = object$variance$type, distribution = object$distribution$type)
-        mod <- try(estimate(spec, solver = solver, scores = FALSE), silent = TRUE)
+        mod <- try(estimate(spec, solver = solver, control = control, scores = FALSE), silent = TRUE)
         model_coef <- coef(mod)
         log_lik <- as.numeric(logLik(mod))
         aic <- as.numeric(AIC(mod))
@@ -242,7 +251,7 @@ tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_
 }
 
 .backtest_top <- function(object, start = floor(length(object$y)/2), end = length(object$y), h = 1,
-                          estimate_every = 1, rolling = FALSE, seed = NULL, solver = "nloptr", 
+                          estimate_every = 1, rolling = FALSE, seed = NULL, solver = "nloptr", control = NULL, 
                           trace = FALSE, ...)
 {
     parameter <- b <- forecast_dates <- NULL
@@ -305,7 +314,7 @@ tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_
                                sampling = object$sampling, sample_n = object$sample_n, 
                                init_garch =  object$init_garch, garch_order = object$garch_order, 
                                variance = object$variance, distribution = object$distribution, top_n = 1)
-        mod <- try(estimate(spec, solver = solver, trace = FALSE), silent = TRUE)
+        mod <- try(estimate(spec, solver = solver, control = control, trace = FALSE), silent = TRUE)
         model_coef <- coef(mod)
         log_lik <- as.numeric(logLik(mod))
         aic <- as.numeric(AIC(mod))
@@ -394,7 +403,7 @@ tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_
 
 .backtest_ensemble <- function(object, start = floor(length(object$y)/2), end = length(object$y),
                                h = 1, estimate_every = 1, rolling = FALSE, weights_scheme = c("U","AIC","BIC"), weights = NULL,
-                               seed = NULL, solver = "nloptr", trace = FALSE, ...)
+                               seed = NULL, solver = "nloptr", control = NULL, trace = FALSE, ...)
 {
     if (weights_scheme == "U") {
         wfun <- function(x) {
@@ -470,7 +479,7 @@ tsbacktest.tsissm.spec <- function(object, start = floor(length(object$target$y_
                                sampling = object$sampling, sample_n = object$sample_n, 
                                init_garch =  object$init_garch, garch_order = object$garch_order, 
                                variance = object$variance, distribution = object$distribution, top_n = top_n)
-        mod <- try(estimate(spec, solver = solver, trace = FALSE), silent = TRUE)
+        mod <- try(estimate(spec, solver = solver, control = control, trace = FALSE), silent = TRUE)
         w <- wfun(mod)
         L <- index_table[[i]]
         M <- split(L, by = "filter_date")

R/estimation.R

@@ -3,7 +3,7 @@
 #' @description Estimates a model given a specification object using
 #' maximum likelihood.
 #' @param object an object of class \dQuote{tsissm.spec} or \dQuote{tsissm.autospec}.
-#' @param solver a choice or either \code{nloptr} or \code{solnp}.
+#' @param solver only \code{nloptr} currently supported.
 #' @param control solver control parameters (see \code{\link{issm_control}}).
 #' @param scores whether to calculate the analytic scores (Jacobian) of the
 #' likelihood. This is not available for the \dQuote{tsissm.autospec} object.
@@ -90,13 +90,11 @@ estimate.tsissm.autospec <- function(object, solver = "nloptr", control = NULL,
         if (solver == "nloptr") {
             control <- issm_control(solver = "nloptr", algorithm = "SLSQP", trace = 0) 
         } else {
-            control <- issm_control(solver = "solnp", trace = 0) 
+            stop("\nonly nloptr supported")
         }
     } else {
         if (solver == "nloptr") {
             control$print_level <- 0
-        } else {
-            control$trace <- 0
         }
     }
     distribution <- object$distribution
@@ -698,13 +696,6 @@ tmb_inputs_issm_constant <- function(spec)
                           lb = spec_list$lower, ub = spec_list$upper, opts = control, fun = fun, issmenv = issmenv) 
         }
         pars <- sol$solution
-    } else if (solver == "solnp") {
-        cfun <- make_constraint_solnp(object, fun, issmenv)
-        sol <- csolnp(pars = fun$par, fn = spec_list$llh_fun, gr = spec_list$grad_fun, ineq_fn = cfun$ineq_fn,
-                      ineq_jac = cfun$ineq_jac, ineq_lower = cfun$ineq_lower, ineq_upper = cfun$ineq_upper, 
-                      lower = spec_list$lower, upper = spec_list$upper, control = control, fun = fun, issmenv = issmenv)
-        sol$status <- sol$convergence
-        pars <- sol$pars
     }
 
     spec_list$data$xseed <- as.numeric(fun$env$data$xseed)
@@ -799,13 +790,6 @@ tmb_inputs_issm_constant <- function(spec)
         }
         spec_list$data$xseed <- as.numeric(fun$env$data$xseed)
         pars <- sol$solution
-    } else if (solver == "solnp") {
-        cfun <- make_constraint_solnp(object, fun, issmenv)
-        sol <- csolnp(pars = fun$par, fn = spec_list$llh_fun, gr = spec_list$grad_fun, ineq_fn = cfun$ineq_fn,
-                      ineq_jac = cfun$ineq_jac, ineq_lower = cfun$ineq_lower, ineq_upper = cfun$ineq_upper, 
-                      lower = spec_list$lower, upper = spec_list$upper, control = control, fun = fun, issmenv = issmenv)
-        sol$status <- sol$convergence
-        pars <- sol$pars
     }
 
     if (return_scores) {

R/solvers.R

@@ -1,28 +1,21 @@
 #' ISSM solver control parameters
 #'
-#' @param solver choice of \sQuote{nloptr} and \sQuote{solnp}
+#' @param solver only \sQuote{nloptr} at present
 #' @param algorithm (nloptr) the algorithm to use
 #' @param trace (integer) controls print level information
 #' @param xtol_rel (nloptr) relative tolerance on optimization parameters
 #' @param maxeval (nloptr) number of function evaluations to stop on
 #' @param xtol_abs (nloptr) absolute tolerances on optimization parameters
-#' @param max_iter (solnp) maximum number of major (outer) iterations
-#' @param min_iter (solnp) maximum number of minor (inner) iterations per outer iteration
-#' @param rho (solnp) onitial penalty parameter for the augmented Lagrangian
-#' @param tol (solnp) convergence tolerance
 #' @details
-#' The function provides defaults for use when using either the \code{\link{nloptr}[nloptr]}
-#' or \code{\link{csolnp}[Rsolnp]}. For the 
-#' former, additional control parameters may be appended to the list if the user so wishes 
-#' (\sQuote{nloptr} has many more options).
-#' @returns a list with the options which is then passed to the appropriate solver.
+#' The function provides defaults for use, but additional control parameters may 
+#' be appended to the list if the user so wishes (\sQuote{nloptr} has many more options).
+#' @returns a list with the options which is then passed to the solver.
 #' @rdname solver_control
 #' @export
 issm_control <- function(solver = "nloptr", algorithm = c("SLSQP","AUGLAG/MMA","AUGLAG/CCSAQ"), 
-                         trace = 0, xtol_rel = 1e-14, maxeval = 1000, xtol_abs = 1e-12,
-                         rho = 1, max_iter = 400, min_iter = 800, tol = 1e-8)
+                         trace = 0, xtol_rel = 1e-14, maxeval = 1000, xtol_abs = 1e-12)
 {
-    solver <- match.arg(solver[1], c("nloptr","solnp"))
+    solver <- match.arg(solver[1], c("nloptr"))
     if (solver == "nloptr") {
         algorithm <- match.arg(algorithm[1], c("SLSQP","AUGLAG/MMA","AUGLAG/CCSAQ"))
         maxeval <- max(1, maxeval)
@@ -37,12 +30,5 @@ issm_control <- function(solver = "nloptr", algorithm = c("SLSQP","AUGLAG/MMA","
                                      local_opts = list(algorithm = "NLOPT_LD_CCSAQ", maxeval = maxeval, xtol_rel = xtol_rel))
         )
         return(ctrl)
-    } else if (solver == "solnp") {
-        max_iter <- max(1, max_iter)
-        min_iter <- min(1, min_iter)
-        trace <- max(0, trace)
-        tol <- abs(tol)
-        ctrl <- list(trace = trace, max_iter = max_iter, min_iter = min_iter, tol = tol)
-        return(ctrl)
     }
 }

R/tsissm-package.R

@@ -12,7 +12,6 @@
 #' @importFrom zoo index as.zoo zoo coredata na.locf na.fill is.zoo
 #' @importFrom sandwich estfun bwNeweyWest vcovHAC vcovOPG bread
 #' @importFrom nloptr nloptr
-#' @importFrom Rsolnp csolnp
 #' @importFrom xts xts as.xts is.xts
 #' @importFrom flextable flextable as_flextable set_caption italic fontsize separate_header add_footer_row add_footer_lines append_chunks as_chunk as_equation as_paragraph compose colformat_double set_header_labels padding bold align autofit hline width
 #' @importFrom future.apply future_lapply

R/tsprofile.R

@@ -16,8 +16,9 @@
 #' @param nsim the number of paths to generate.
 #' @param seed an object specifying if and how the random number generator
 #' should be initialized. See the simulate documentation for more details.
-#' @param solver either \dQuote{nlortr} or \dQuote{solnp}.
-#' @param trace whether to show the progress bar and additionally output verbose messages. 
+#' @param solver only \dQuote{nlortr} currently available.
+#' @param trace whether to show the progress bar and additionally output verbose messages.
+#' @param control optional control parameters.
 #' The user is expected to have set up appropriate handlers for this using 
 #' the \code{\link[progressr]{handlers}} function from the \dQuote{progressr} package.
 #' @param ... not currently used.
@@ -32,12 +33,16 @@
 #' @rdname tsprofile
 #' @export
 #'
-tsprofile.tsissm.estimate <- function(object, h = 1, nsim = 100, seed = NULL, solver = "nloptr", trace = FALSE, ...)
+tsprofile.tsissm.estimate <- function(object, h = 1, nsim = 100, seed = NULL, solver = "nloptr", control = NULL, trace = FALSE, ...)
 {
+    solver <- "nloptr"
+    if (is.null(control)) {
+        control <- issm_control(solver = "nloptr", algorithm = "SLSQP", trace = 0)
+    }
     sim <- simulate(object, seed = seed, nsim = nsim, h = length(object$spec$target$y_orig) + h)
     sim <- .check_positivity(sim)
     if (is.null(sim)) return(sim)
-    profile <- profile_fun(sim$distribution, object, h = h, solver = solver, control = NULL, trace = trace)
+    profile <- profile_fun(sim$distribution, object, h = h, solver = solver, control = control, trace = trace)
     class(profile) <- "tsissm.profile"
     return(profile)
 }

README.md

@@ -2,7 +2,7 @@
 # tsissm <img src="man/figures/logo.png" align="right" height="139" alt="" />
 
 [![R-CMD-check](https://github.com/tsmodels/tsissm/actions/workflows/rcmdcheck.yaml/badge.svg)](https://github.com/tsmodels/tsissm/actions/workflows/rcmdcheck.yaml)
-[![Last-changedate](https://img.shields.io/badge/last%20change-2025--07--05-yellowgreen.svg)](/commits/master)
+[![Last-changedate](https://img.shields.io/badge/last%20change-2025--07--12-yellowgreen.svg)](/commits/master)
 [![packageversion](https://img.shields.io/badge/Package%20version-1.0.2-orange.svg?style=flat-square)](commits/master)
 [![CRAN_Status_Badge](https://www.r-pkg.org/badges/version/tsissm)](https://cran.r-project.org/package=tsissm)