fix docs, penalty to integer, fix indexing for subsampling

R/ame.R

@@ -13,6 +13,7 @@
 #' estimates are returned.
 #' @param probs Vector defining the quantiles of interest. Default is 
 #' `c(0.025, 0.5, 0.975)`.
+#' @param ... Further arguments passed to specific methods.
 #' @rdname ame
 #' @export
 ame <- function(model, variable, values, ...) {

R/estimate_nhmm.R

@@ -49,6 +49,10 @@
 #' reason for the termination of optimization (e.g., reached the maximum number 
 #' of iterations) is only available by the final message show with 
 #' `verbose = TRUE`.
+#' @param penalize If `TRUE` (default), penalize the likelihood with a prior on 
+#' the coefficients
+#' @param penalty If `penalize = TRUE`, the standard deviation of the zero-mean 
+#' normal density used as a penalizing prior for the coefficients. Default is `5`.
 #' @param ... Additional arguments to [rstan::optimizing()]. Most importantly,
 #' argument `iter` defines the maximum number of iterations for optimization.
 #' The default is `2000`.

R/fit_mnhmm.R

@@ -58,7 +58,7 @@ fit_mnhmm <- function(model, inits, init_sd, restarts, threads, verbose, penaliz
           model_code, init = init,
           data = list(
             penalty = penalty,
-            penalize = penalize,
+            penalize = as.integer(penalize),
             N = model$n_sequences,
             T = model$sequence_lengths,
             max_T = model$length_of_sequences, 
@@ -110,7 +110,7 @@ fit_mnhmm <- function(model, inits, init_sd, restarts, threads, verbose, penaliz
       model_code, 
       data = list(
         penalty = penalty,
-        penalize = penalize,
+        penalize = as.integer(penalize),
         N = model$n_sequences,
         T = model$sequence_lengths,
         max_T = model$length_of_sequences,

R/fit_mnhmm2.R

@@ -62,8 +62,8 @@ fit_mnhmm2 <- function(model, inits, init_sd, restarts, threads, verbose,
         c(list(
           model_code, init = init,
           data = list(
-            penalize = penalize, 
             penalty = penalty, 
+            penalize = as.integer(penalize),
             N = model$n_sequences,
             T = model$sequence_lengths,
             max_T = model$length_of_sequences, 
@@ -102,8 +102,8 @@ fit_mnhmm2 <- function(model, inits, init_sd, restarts, threads, verbose,
         c(list(
           model_code, init = out0[[idx[i]]]$par,
           data = list(
-            penalize = penalize, 
             penalty = penalty, 
+            penalize = as.integer(penalize),
             N = model$n_sequences,
             T = model$sequence_lengths,
             max_T = model$length_of_sequences, 
@@ -154,8 +154,8 @@ fit_mnhmm2 <- function(model, inits, init_sd, restarts, threads, verbose,
     c(list(
       model_code, 
       data = list(
-        penalize = penalize, 
         penalty = penalty, 
+        penalize = as.integer(penalize),
         N = model$n_sequences,
         T = model$sequence_lengths,
         max_T = model$length_of_sequences,

R/fit_nhmm.R

@@ -54,7 +54,7 @@ fit_nhmm <- function(model, inits, init_sd, restarts, threads, verbose, penalize
           model_code, init = init,
           data = list(
             penalty = penalty,
-            penalize = penalize,
+            penalize = as.integer(penalize),
             N = model$n_sequences,
             T = model$sequence_lengths,
             max_T = model$length_of_sequences, 
@@ -102,8 +102,8 @@ fit_nhmm <- function(model, inits, init_sd, restarts, threads, verbose, penalize
     c(list(
       model_code, 
       data = list(
-        penalty = penalty,
-        penalize = penalize,
+        penalty = penalty, 
+        penalize = as.integer(penalize),
         N = model$n_sequences,
         T = model$sequence_lengths,
         max_T = model$length_of_sequences, 

inst/stan/include/loglik_multichannel.stan

@@ -18,7 +18,7 @@ array[,] vector X_o, array[] vector X_i) {
   vector[S] zeros_S = rep_vector(0, S);
   array[C + 1] int cumsum_M = append_array({0}, cumulative_sum(M));
 
-  vector[N] ll;
+  vector[N_sample] ll;
   vector[S] log_Pi;
   matrix[S, S] log_A;
   matrix[S, max_M + 1] log_B;
@@ -71,7 +71,7 @@ array[,] vector X_o, array[] vector X_i) {
       }
       log_alpha = log_alpha_new;
     }
-    ll[i] = log_sum_exp(log_alpha);
+    ll[ii] = log_sum_exp(log_alpha);
   }
   return sum(ll);
 }

inst/stan/include/loglik_multichannel_mixture.stan

@@ -21,7 +21,7 @@ array[,] vector X_o, array[] vector X_i, array[] vector X_d) {
   vector[S] zeros_S = rep_vector(0, S);
   array[C + 1] int cumsum_M = append_array({0}, cumulative_sum(M));
 
-  matrix[N, D] ll;
+  matrix[N_sample, D] ll;
   vector[S] log_Pi;
   matrix[S, S] log_A;
   matrix[S, max_M + 1] log_B;
@@ -76,14 +76,14 @@ array[,] vector X_o, array[] vector X_i, array[] vector X_d) {
         }
         log_alpha = log_alpha_new;
       }
-      ll[i, d] = log_sum_exp(log_alpha);
+      ll[ii, d] = log_sum_exp(log_alpha);
     }
   }
   vector[N_sample] ll_i;
   for(ii in 1:N_sample) {
     int i = ids[ii];
     vector[D] log_omega = log_softmax(theta * X_d[i]);
-    ll_i[i] = log_sum_exp(log_omega + ll[i, ]');
+    ll_i[ii] = log_sum_exp(log_omega + ll[ii, ]');
   }
   return sum(ll_i);
 }

inst/stan/include/loglik_singlechannel.stan

@@ -15,7 +15,7 @@ array[,] vector X_o, array[] vector X_i) {
   row_vector[K_o] zeros_row_K_o = rep_row_vector(0, K_o);
   vector[S] zeros_S = rep_vector(0, S);
 
-  vector[N] ll;
+  vector[N_sample] ll;
   vector[S] log_Pi;
   matrix[S, S] log_A;
   matrix[S, M + 1] log_B;
@@ -53,7 +53,7 @@ array[,] vector X_o, array[] vector X_i) {
       }
       log_alpha = log_alpha_new;
     }
-    ll[i] = log_sum_exp(log_alpha);
+    ll[ii] = log_sum_exp(log_alpha);
   }
   return sum(ll);
 }

inst/stan/include/loglik_singlechannel_mixture.stan

@@ -18,7 +18,7 @@ array[,] vector X_o, array[] vector X_i, array[] vector X_d) {
   row_vector[K_d] zeros_row_K_d = rep_row_vector(0, K_d);
   vector[S] zeros_S = rep_vector(0, S);
 
-  matrix[N, D] ll;
+  matrix[N_sample, D] ll;
   vector[S] log_Pi;
   matrix[S, S] log_A;
   matrix[S, M + 1] log_B;
@@ -58,14 +58,14 @@ array[,] vector X_o, array[] vector X_i, array[] vector X_d) {
         }
         log_alpha = log_alpha_new;
       }
-      ll[i, d] = log_sum_exp(log_alpha);
+      ll[ii, d] = log_sum_exp(log_alpha);
     }
   }
   vector[N_sample] ll_i;
   for(ii in 1:N_sample) {
     int i = ids[ii];
     vector[D] log_omega = log_softmax(theta * X_d[i]);
-    ll_i[i] = log_sum_exp(log_omega + ll[i, ]');
+    ll_i[ii] = log_sum_exp(log_omega + ll[ii, ]');
   }
   return sum(ll_i);
 }

src/stanExports_mnhmm.h

@@ -138,9 +138,9 @@ static constexpr std::array<const char*, 183> locations_array__ =
   " (in 'string', line 45, column 2 to column 57)",
   " (in 'string', line 46, column 9 to column 10)",
   " (in 'string', line 46, column 2 to column 39)",
-  " (in 'string', line 47, column 9 to column 10)",
-  " (in 'string', line 47, column 12 to column 13)",
-  " (in 'string', line 47, column 2 to column 18)",
+  " (in 'string', line 47, column 9 to column 17)",
+  " (in 'string', line 47, column 19 to column 20)",
+  " (in 'string', line 47, column 2 to column 25)",
   " (in 'string', line 48, column 9 to column 10)",
   " (in 'string', line 48, column 2 to column 19)",
   " (in 'string', line 49, column 9 to column 10)",
@@ -196,7 +196,7 @@ static constexpr std::array<const char*, 183> locations_array__ =
   " (in 'string', line 85, column 8 to column 34)",
   " (in 'string', line 73, column 24 to line 86, column 7)",
   " (in 'string', line 73, column 6 to line 86, column 7)",
-  " (in 'string', line 87, column 6 to column 40)",
+  " (in 'string', line 87, column 6 to column 41)",
   " (in 'string', line 64, column 26 to line 88, column 5)",
   " (in 'string', line 64, column 4 to line 88, column 5)",
   " (in 'string', line 58, column 16 to line 89, column 3)",
@@ -206,7 +206,7 @@ static constexpr std::array<const char*, 183> locations_array__ =
   " (in 'string', line 92, column 4 to column 20)",
   " (in 'string', line 93, column 11 to column 12)",
   " (in 'string', line 93, column 4 to column 54)",
-  " (in 'string', line 94, column 4 to column 48)",
+  " (in 'string', line 94, column 4 to column 50)",
   " (in 'string', line 91, column 24 to line 95, column 3)",
   " (in 'string', line 91, column 2 to line 95, column 3)",
   " (in 'string', line 96, column 2 to column 19)",
@@ -443,11 +443,12 @@ loglik_sc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
     stan::model::assign(zeros_S, stan::math::rep_vector(0, S),
       "assigning variable zeros_S");
     current_statement__ = 110;
-    stan::math::validate_non_negative_index("ll", "N", N);
+    stan::math::validate_non_negative_index("ll", "N_sample", N_sample);
     current_statement__ = 111;
     stan::math::validate_non_negative_index("ll", "D", D);
     Eigen::Matrix<local_scalar_t__,-1,-1> ll =
-      Eigen::Matrix<local_scalar_t__,-1,-1>::Constant(N, D, DUMMY_VAR__);
+      Eigen::Matrix<local_scalar_t__,-1,-1>::Constant(N_sample, D,
+        DUMMY_VAR__);
     current_statement__ = 113;
     stan::math::validate_non_negative_index("log_Pi", "S", S);
     Eigen::Matrix<local_scalar_t__,-1,1> log_Pi =
@@ -630,7 +631,7 @@ loglik_sc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
         }
         current_statement__ = 168;
         stan::model::assign(ll, stan::math::log_sum_exp(log_alpha),
-          "assigning variable ll", stan::model::index_uni(i),
+          "assigning variable ll", stan::model::index_uni(ii),
           stan::model::index_uni(d));
       }
     }
@@ -658,9 +659,9 @@ loglik_sc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
         stan::math::log_sum_exp(
           stan::math::add(log_omega,
             stan::math::transpose(
-              stan::model::rvalue(ll, "ll", stan::model::index_uni(i),
+              stan::model::rvalue(ll, "ll", stan::model::index_uni(ii),
                 stan::model::index_omni())))), "assigning variable ll_i",
-        stan::model::index_uni(i));
+        stan::model::index_uni(ii));
     }
     current_statement__ = 181;
     return stan::math::sum(ll_i);

src/stanExports_multichannel_mnhmm.h

@@ -145,9 +145,9 @@ static constexpr std::array<const char*, 206> locations_array__ =
   " (in 'string', line 48, column 2 to column 39)",
   " (in 'string', line 49, column 8 to column 13)",
   " (in 'string', line 49, column 2 to column 67)",
-  " (in 'string', line 50, column 9 to column 10)",
-  " (in 'string', line 50, column 12 to column 13)",
-  " (in 'string', line 50, column 2 to column 18)",
+  " (in 'string', line 50, column 9 to column 17)",
+  " (in 'string', line 50, column 19 to column 20)",
+  " (in 'string', line 50, column 2 to column 25)",
   " (in 'string', line 51, column 9 to column 10)",
   " (in 'string', line 51, column 2 to column 19)",
   " (in 'string', line 52, column 9 to column 10)",
@@ -219,7 +219,7 @@ static constexpr std::array<const char*, 206> locations_array__ =
   " (in 'string', line 103, column 8 to column 34)",
   " (in 'string', line 88, column 24 to line 104, column 7)",
   " (in 'string', line 88, column 6 to line 104, column 7)",
-  " (in 'string', line 105, column 6 to column 40)",
+  " (in 'string', line 105, column 6 to column 41)",
   " (in 'string', line 76, column 26 to line 106, column 5)",
   " (in 'string', line 76, column 4 to line 106, column 5)",
   " (in 'string', line 61, column 16 to line 107, column 3)",
@@ -229,7 +229,7 @@ static constexpr std::array<const char*, 206> locations_array__ =
   " (in 'string', line 110, column 4 to column 20)",
   " (in 'string', line 111, column 11 to column 12)",
   " (in 'string', line 111, column 4 to column 54)",
-  " (in 'string', line 112, column 4 to column 48)",
+  " (in 'string', line 112, column 4 to column 50)",
   " (in 'string', line 109, column 24 to line 113, column 3)",
   " (in 'string', line 109, column 2 to line 113, column 3)",
   " (in 'string', line 114, column 2 to column 19)",
@@ -482,11 +482,12 @@ loglik_mc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
       stan::math::append_array(std::vector<int>{0},
         stan::math::cumulative_sum(M)), "assigning variable cumsum_M");
     current_statement__ = 117;
-    stan::math::validate_non_negative_index("ll", "N", N);
+    stan::math::validate_non_negative_index("ll", "N_sample", N_sample);
     current_statement__ = 118;
     stan::math::validate_non_negative_index("ll", "D", D);
     Eigen::Matrix<local_scalar_t__,-1,-1> ll =
-      Eigen::Matrix<local_scalar_t__,-1,-1>::Constant(N, D, DUMMY_VAR__);
+      Eigen::Matrix<local_scalar_t__,-1,-1>::Constant(N_sample, D,
+        DUMMY_VAR__);
     current_statement__ = 120;
     stan::math::validate_non_negative_index("log_Pi", "S", S);
     Eigen::Matrix<local_scalar_t__,-1,1> log_Pi =
@@ -725,7 +726,7 @@ loglik_mc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
         }
         current_statement__ = 191;
         stan::model::assign(ll, stan::math::log_sum_exp(log_alpha),
-          "assigning variable ll", stan::model::index_uni(i),
+          "assigning variable ll", stan::model::index_uni(ii),
           stan::model::index_uni(d));
       }
     }
@@ -753,9 +754,9 @@ loglik_mc_mix(const std::vector<Eigen::Matrix<T0__,-1,-1>>& beta_i_raw,
         stan::math::log_sum_exp(
           stan::math::add(log_omega,
             stan::math::transpose(
-              stan::model::rvalue(ll, "ll", stan::model::index_uni(i),
+              stan::model::rvalue(ll, "ll", stan::model::index_uni(ii),
                 stan::model::index_omni())))), "assigning variable ll_i",
-        stan::model::index_uni(i));
+        stan::model::index_uni(ii));
     }
     current_statement__ = 204;
     return stan::math::sum(ll_i);