Merge pull request #71 from andrjohns/array-syntax

Update deprecated syntax for future rstan compatibility

DESCRIPTION

@@ -18,7 +18,7 @@ Imports:
     Rcpp (>= 0.12.0),
     rlang,
     RSpectra,
-    rstan (>= 2.18.1),
+    rstan (>= 2.26.0),
     rstantools (>= 2.0.0),
     stats
 Suggests: covr, devtools, knitr, pkgdown, raster, rmarkdown, roxygen2, testthat
@@ -42,8 +42,8 @@ LinkingTo:
     Rcpp (>= 0.12.0),
     RcppArmadillo (>= 0.9.300.2.0),
     RcppEigen (>= 0.3.3.3.0),
-    rstan (>= 2.18.1),
-    StanHeaders (>= 2.18.0),
+    rstan (>= 2.26.0),
+    StanHeaders (>= 2.26.0),
     RcppParallel (>= 5.0.2)
 SystemRequirements: GNU make
 Collate: 

inst/stan/colext.stan

@@ -3,7 +3,7 @@ functions{
 #include /include/functions_priors.stan
 
 //can shortcut here I think
-vector get_pY(int[] y, vector logit_p, int nd){
+vector get_pY(array[] int y, vector logit_p, int nd){
   vector[2] out;
   out[1] = nd;
   out[2] = exp(bernoulli_logit_lpmf(y | logit_p));
@@ -29,8 +29,8 @@ matrix get_phi(matrix phi_raw, int Tstart, int Tnext){
 }
 
 //Ts = indices of primary periods when site was sampled (eg not all NA)
-real lp_colext(int[] y, int[] Tsamp, int[] J, row_vector psi, matrix phi_raw,
-               vector logit_p, int[] nd){
+real lp_colext(array[] int y, array[] int Tsamp, array[] int J, row_vector psi, matrix phi_raw,
+               vector logit_p, array[] int nd){
 
   int T = size(Tsamp);
   matrix[2,2] phi_prod = diag_matrix(rep_vector(1, 2));
@@ -56,9 +56,9 @@ real lp_colext(int[] y, int[] Tsamp, int[] J, row_vector psi, matrix phi_raw,
 }
 
 //needs fixed
-vector get_loglik_colext(int[] y, int M, int[] Tsamp, int[,] J, int[,] si,
+vector get_loglik_colext(array[] int y, int M, array[] int Tsamp, array[,] int J, array[,] int si,
                          matrix psi_raw, matrix phi_raw, vector logit_p,
-                         int[,] nd){
+                         array[,] int nd){
   vector[M] out;
   for (i in 1:M){
     out[i] = lp_colext(y[si[i,1]:si[i,2]], Tsamp[si[i,3]:si[i,4]], J[i,],
@@ -82,32 +82,32 @@ int n_fixed_col;
 int n_fixed_ext;
 int n_group_vars_col;
 int n_group_vars_ext;
-int n_random_col[has_random_col ? n_group_vars_col : 1];
-int n_random_ext[has_random_ext ? n_group_vars_ext: 1];
+array[has_random_col ? n_group_vars_col : 1] int n_random_col;
+array[has_random_ext ? n_group_vars_ext: 1] int n_random_ext;
 matrix[M*(T-1), n_fixed_col] X_col;
 matrix[M*(T-1), n_fixed_ext] X_ext;
 vector[M*(T-1)] offset_col;
 vector[M*(T-1)] offset_ext;
 
-int Zdim_col[5];
+array[5] int Zdim_col;
 vector[Zdim_col[3]] Zw_col;
-int Zv_col[Zdim_col[4]];
-int Zu_col[Zdim_col[5]];
+array[Zdim_col[4]] int Zv_col;
+array[Zdim_col[5]] int Zu_col;
 
-int Zdim_ext[5];
+array[5] int Zdim_ext;
 vector[Zdim_ext[3]] Zw_ext;
-int Zv_ext[Zdim_ext[4]];
-int Zu_ext[Zdim_ext[5]];
+array[Zdim_ext[4]] int Zv_ext;
+array[Zdim_ext[5]] int Zu_ext;
 
-int prior_dist_col[3];
-int prior_dist_ext[3];
+array[3] int prior_dist_col;
+array[3] int prior_dist_ext;
 matrix[3, (n_fixed_col+1)] prior_pars_col;
 matrix[3, (n_fixed_ext+1)] prior_pars_ext;
 }
 
 transformed data{
 
-int no_detects[M, T];
+array[M, T] int no_detects;
 int include_scale;
 int include_shape;
 

inst/stan/include/data.stan

@@ -3,19 +3,19 @@ int model_code;
 int M;
 int T;
 int Tsamp_size;
-int Tsamp[Tsamp_size];
+array[Tsamp_size] int Tsamp;
 int R;
-int J[M,T];
-int y[R];
-int si[M, 6];
+array[M,T] int J;
+array[R] int y;
+array[M, 6] int si;
 int K;
-int Kmin[M,T];
+array[M,T] int Kmin;
 int y_dist;
 int z_dist;
 int n_aux1;
 int n_aux2;
 int n_aux3;
-int aux1[n_aux1]; //Used for various auxiliary data
+array[n_aux1] int aux1; //Used for various auxiliary data
 vector[n_aux2] aux2;
 vector[n_aux3] aux3;
 
@@ -27,28 +27,28 @@ int n_fixed_state;
 int n_fixed_det;
 int n_group_vars_state;
 int n_group_vars_det;
-int n_random_state[has_random_state ? n_group_vars_state : 1];
-int n_random_det[has_random_det ? n_group_vars_det: 1];
+array[has_random_state ? n_group_vars_state : 1] int n_random_state;
+array[has_random_det ? n_group_vars_det: 1] int n_random_det;
 matrix[n_obs_state, n_fixed_state] X_state;
 matrix[n_obs_det, n_fixed_det] X_det;
 vector[n_obs_state] offset_state;
 vector[n_obs_det] offset_det;
 
-int Zdim_state[5];
+array[5] int Zdim_state;
 vector[Zdim_state[3]] Zw_state;
-int Zv_state[Zdim_state[4]];
-int Zu_state[Zdim_state[5]];
+array[Zdim_state[4]] int Zv_state;
+array[Zdim_state[5]] int Zu_state;
 
-int Zdim_det[5];
+array[5] int Zdim_det;
 vector[Zdim_det[3]] Zw_det;
-int Zv_det[Zdim_det[4]];
-int Zu_det[Zdim_det[5]];
+array[Zdim_det[4]] int Zv_det;
+array[Zdim_det[5]] int Zu_det;
 
 // Stuff for custom priors
-int prior_dist_state[3];
-int prior_dist_det[3];
-int prior_dist_shape[3];
-int prior_dist_scale[3];
+array[3] int prior_dist_state;
+array[3] int prior_dist_det;
+array[3] int prior_dist_shape;
+array[3] int prior_dist_scale;
 
 matrix[3, (n_fixed_state+1)] prior_pars_state;
 matrix[3, (n_fixed_det+1)] prior_pars_det;

inst/stan/include/functions_distsamp.stan

@@ -1,4 +1,4 @@
-real lp_distsamp(int[] y, vector db, real log_lambda, real par1, real par2,
+real lp_distsamp(array[] int y, vector db, real log_lambda, real par1, real par2,
                  int point, int keyfun, vector conv_const){
 
   real lam = exp(log_lambda);
@@ -13,7 +13,7 @@ real lp_distsamp(int[] y, vector db, real log_lambda, real par1, real par2,
   return loglik;
 }
 
-vector get_loglik_distsamp(int[] y, int M, vector db, int[,] si,
+vector get_loglik_distsamp(array[] int y, int M, vector db, array[,] int si,
                            vector log_lambda, vector trans_par1, int z_dist,
                            real trans_par2, int point, int keyfun, vector conv_const){
 

inst/stan/include/functions_keyfuns.stan

@@ -39,23 +39,23 @@ real int_negexp(real log_rate, real a, real b, int point){
   return out;
 }
 
-//real p_hazard_line(real x, real xc, real[] theta, real[] x_r, int[] x_i){
+//real p_hazard_line(real x, real xc, array[] real theta, array[] real x_r, array[] int x_i){
 //  return(1 - exp(-1 * pow((x/theta[1]), (-1*theta[2]))));
 //}
 
-real p_hazard_line(real x, real[] theta){
+real p_hazard_line(real x, array[] real theta){
   return(1 - exp(-1 * pow((x/theta[1]), (-1*theta[2]))));
 }
 
-//real p_hazard_point(real x, real xc, real[] theta, real[] x_r, int[] x_i){
+//real p_hazard_point(real x, real xc, array[] real theta, array[] real x_r, array[] int x_i){
 //  return((1 - exp(-1 * pow(x/theta[1], -1*theta[2]))) * x);
 //}
 
-real p_hazard_point(real x, real[] theta){
+real p_hazard_point(real x, array[] real theta){
   return((1 - exp(-1 * pow(x/theta[1], -1*theta[2]))) * x);
 }
 
-real trap_rule_line(real[] theta, real a, real b){
+real trap_rule_line(array[] real theta, real a, real b){
   int n = 100;
   real h = (b - a) / n;
 
@@ -66,7 +66,7 @@ real trap_rule_line(real[] theta, real a, real b){
   return h/2 * (p_hazard_line(a, theta) + 2*int_sum + p_hazard_line(b, theta));
 }
 
-real trap_rule_point(real[] theta, real a, real b){
+real trap_rule_point(array[] real theta, real a, real b){
   int n = 100;
   real h = (b - a) / n;
 
@@ -83,7 +83,7 @@ real int_hazard(real log_shape, real log_scale, real a, real b, int point){
   real out;
   real shape = exp(log_shape);
   real scale = exp(log_scale);
-  real theta[2];
+  array[2] real theta;
   theta[1] = shape;
   theta[2] = scale;
 

inst/stan/include/functions_multinomPois.stan

@@ -32,7 +32,7 @@ vector pi_fun(int pi_type, vector p, int J){
   return out;
 }
 
-real lp_multinomPois(int[] y, real log_lambda, vector logit_p, int pi_type){
+real lp_multinomPois(array[] int y, real log_lambda, vector logit_p, int pi_type){
 
   real loglik = 0.0;
   real lam = exp(log_lambda);
@@ -52,11 +52,11 @@ real lp_multinomPois(int[] y, real log_lambda, vector logit_p, int pi_type){
   return loglik;
 }
 
-vector get_loglik_multinomPois(int[] y, int M, int[,] si, vector log_lambda,
+vector get_loglik_multinomPois(array[] int y, int M, array[,] int si, vector log_lambda,
                                vector logit_p, int pi_type){
 
   vector[M] out;
-  int J = num_elements(logit_p) / M; // use %/% in future version of stan
+  int J = num_elements(logit_p) %/% M; // use %/% in future version of stan
   int pstart = 1;
   int pend;
   for (i in 1:M){

inst/stan/include/functions_occu.stan

@@ -1,4 +1,4 @@
-real lp_occu(int[] y, real logit_psi, vector logit_p, int Kmin){
+real lp_occu(array[] int y, real logit_psi, vector logit_p, int Kmin){
   real out;
   out = log_inv_logit(logit_psi) + bernoulli_logit_lpmf(y | logit_p);
   if(Kmin == 1){
@@ -7,8 +7,8 @@ real lp_occu(int[] y, real logit_psi, vector logit_p, int Kmin){
   return log_sum_exp(out, log1m_inv_logit(logit_psi));
 }
 
-vector get_loglik_occu(int[] y, int M, int[,] J, int[,] si, vector logit_psi,
-                  vector logit_p, int[] Kmin){
+vector get_loglik_occu(array[] int y, int M, array[,] int J, array[,] int si, vector logit_psi,
+                  vector logit_p, array[] int Kmin){
   vector[M] out;
   for (i in 1:M){
     out[i] = lp_occu(y[si[i,1]:si[i,2]], logit_psi[i],

inst/stan/include/functions_occuRN.stan

@@ -1,4 +1,4 @@
-real lp_rn(int[] y, real log_lambda, vector logit_r, int J, int K, int Kmin){
+real lp_rn(array[] int y, real log_lambda, vector logit_r, int J, int K, int Kmin){
 
   int numN = K - Kmin + 1;
   vector[J] q = 1 - inv_logit(logit_r);
@@ -15,8 +15,8 @@ real lp_rn(int[] y, real log_lambda, vector logit_r, int J, int K, int Kmin){
   return log_sum_exp(lp);
 }
 
-vector get_loglik_rn(int[] y, int M, int[,] J, int[,] si, vector log_lambda,
-                     vector logit_p, int K, int[] Kmin){
+vector get_loglik_rn(array[] int y, int M, array[,] int J, array[,] int si, vector log_lambda,
+                     vector logit_p, int K, array[] int Kmin){
   vector[M] out;
   for (i in 1:M){
     out[i] = lp_rn(y[si[i,1]:si[i,2]], log_lambda[i], logit_p[si[i,1]:si[i,2]],

inst/stan/include/functions_occuTTD.stan

@@ -1,4 +1,4 @@
-vector ttd_prob_exp(vector y, vector log_lam, int[] delta){
+vector ttd_prob_exp(vector y, vector log_lam, array[] int delta){
   int J = num_elements(y);
   vector[J] e_lamt;
   real lam;
@@ -9,7 +9,7 @@ vector ttd_prob_exp(vector y, vector log_lam, int[] delta){
   return e_lamt;
 }
 
-vector ttd_prob_weib(vector y, vector log_lam, int[] delta, real log_k){
+vector ttd_prob_weib(vector y, vector log_lam, array[] int delta, real log_k){
   int J = num_elements(y);
   vector[J] e_lamt;
   real k = exp(log_k);
@@ -23,7 +23,7 @@ vector ttd_prob_weib(vector y, vector log_lam, int[] delta, real log_k){
 }
 
 real lp_occuTTD(vector y, real logit_psi, vector log_lam,
-                real log_k, int[] delta, int ydist){
+                real log_k, array[] int delta, int ydist){
 
   int J = num_elements(y);
   vector[J] e_lamt;
@@ -40,8 +40,8 @@ real lp_occuTTD(vector y, real logit_psi, vector log_lam,
   return log(lik);
 }
 
-vector get_loglik_occuTTD(vector y, int M, int[,] si, vector logit_psi,
-                          vector log_lam, real log_k, int[] delta, int ydist){
+vector get_loglik_occuTTD(vector y, int M, array[,] int si, vector logit_psi,
+                          vector log_lam, real log_k, array[] int delta, int ydist){
   vector[M] out;
   for (i in 1:M){
     out[i] = lp_occuTTD(y[si[i,1]:si[i,2]], logit_psi[i], log_lam[si[i,1]:si[i,2]],

inst/stan/include/functions_pcount.stan

@@ -1,4 +1,4 @@
-real lp_pcount_pois(int[] y, real log_lambda, vector logit_p, int K, int Kmin){
+real lp_pcount_pois(array[] int y, real log_lambda, vector logit_p, int K, int Kmin){
 
   real fac = 1;
   real ff = exp(log_lambda) * prod(1 - inv_logit(logit_p));
@@ -18,8 +18,8 @@ real lp_pcount_pois(int[] y, real log_lambda, vector logit_p, int K, int Kmin){
           log(fac);
 }
 
-vector get_loglik_pcount(int[] y, int M, int[,] J, int[,] si, vector log_lambda,
-                         vector logit_p, int z_dist, real beta_scale, int K, int[] Kmin){
+vector get_loglik_pcount(array[] int y, int M, array[,] int J, array[,] int si, vector log_lambda,
+                         vector logit_p, int z_dist, real beta_scale, int K, array[] int Kmin){
   vector[M] out;
   for (i in 1:M){
     out[i] = lp_pcount_pois(y[si[i,1]:si[i,2]], log_lambda[i],

inst/stan/include/functions_priors.stan

@@ -18,7 +18,7 @@ real lp_single_prior(vector x, int dist, row_vector pars1,
 }
 
 
-real lp_priors(vector beta, int[] dist, matrix pars){
+real lp_priors(vector beta, array[] int dist, matrix pars){
 
   int idx;
   real out = 0.0;
@@ -42,7 +42,7 @@ real lp_priors(vector beta, int[] dist, matrix pars){
 }
 
 real lp_random_prior(int has_random, int n_group_vars, vector b,
-                     int[] n_random, vector sigma, int dist, matrix pars){
+                     array[] int n_random, vector sigma, int dist, matrix pars){
   int idx = 1;
   real out = 0;
   int par_idx = cols(pars);

inst/stan/spatial.stan

@@ -16,7 +16,7 @@ real theta_lpdf(vector theta, real tau, matrix Qalpha, int n_eigen){
   return 0.5*(n_eigen*log(tau) - tau*quad_form(Qalpha, theta));
 }
 
-//real lp_occu_probit(int[] y, real raw_psi, vector raw_p, int Kmin){
+//real lp_occu_probit(array[] int y, real raw_psi, vector raw_p, int Kmin){
 //  real out;
 //  real psi = Phi(raw_psi);
 //  int J = num_elements(raw_p);
@@ -28,8 +28,8 @@ real theta_lpdf(vector theta, real tau, matrix Qalpha, int n_eigen){
 //  return log(out);
 //}
 
-//vector get_loglik_occu_probit(int[] y, int M, int[,] J, int[,] si, vector raw_psi,
-//                  vector raw_p, int[] Kmin){
+//vector get_loglik_occu_probit(array[] int y, int M, array[,] int J, array[,] int si, vector raw_psi,
+//                  vector raw_p, array[] int Kmin){
 //  vector[M] out;
 //  for (i in 1:M){
 //    out[i] = lp_occu_probit(y[si[i,1]:si[i,2]], raw_psi[i],