#  Elapsed Time: 44.18 seconds (Warm-up)

#  Elapsed Time: 44.213 seconds (Warm-up)

#  Elapsed Time: 44.481 seconds (Warm-up)

#  Elapsed Time: 44.713 seconds (Warm-up)

*   x: array of continuous predictor values
*   zgp: vector of independent standard normal variables
*   a vector to be added to the linear predictor
vector gp(data array[] vector x, real sdgp, vector lscale, vector zgp) {
cov = gp_exp_quad_cov(x, sdgp, lscale[1]);
cov = gp_exp_quad_cov(x[ : , 1], sdgp, lscale[1]);
cov = cov .* gp_exp_quad_cov(x[ : , d], 1, lscale[d]);
// deal with numerical non-positive-definiteness
/* Spectral density function of a Gaussian process
* with squared exponential covariance kernel
*   x: array of numeric values of dimension NB x D
*   lscale: vector of length-scale parameters
*   numeric values of the function evaluated at 'x'
vector spd_cov_exp_quad(data array[] vector x, real sdgp, vector lscale) {
real constant = square(sdgp) * (sqrt(2 * pi()) * lscale[1]) ^ D;
real neg_half_lscale2 = -0.5 * square(lscale[1]);
out[m] = constant * exp(neg_half_lscale2 * dot_self(x[m]));
real constant = square(sdgp) * sqrt(2 * pi()) ^ D * prod(lscale);
vector[Dls] neg_half_lscale2 = -0.5 * square(lscale);
out[m] = constant * exp(dot_product(neg_half_lscale2, square(x[m])));
/* how many elements are in a range of integers?
*   start: start of the range (inclusive)
int size_range(array[] int x, int start, int end) {
/* which elements are in a range of integers?
array[] int which_range(array[] int x, int start, int end) {
array[size_range(x, start, end)] int out;
*   start: start of the range of values in x (inclusive)
array[] int start_at_one(array[] int x, int start) {
int<lower=1> N; // total number of observations
int<lower=1> Kgp_1; // number of sub-GPs (equal to 1 unless 'by' was used)
// indices of latent GP groups per observation
array[Nsubgp_1] vector[Dgp_1] Xgp_1; // covariates of the GP
int<lower=1> N_1; // number of grouping levels
int<lower=1> M_1; // number of coefficients per level
array[N] int<lower=1> J_1; // grouping indicator per observation
int<lower=1> N_2; // number of grouping levels
int<lower=1> M_2; // number of coefficients per level
array[N] int<lower=1> J_2; // grouping indicator per observation
int prior_only; // should the likelihood be ignored?
real Intercept; // temporary intercept for centered predictors
vector<lower=0>[Kgp_1] sdgp_1; // GP standard deviation parameters
array[Kgp_1] vector<lower=0>[1] lscale_1; // GP length-scale parameters
vector[Nsubgp_1] zgp_1; // latent variables of the GP
real<lower=0> sigma; // dispersion parameter
vector<lower=0>[M_1] sd_1; // group-level standard deviations
array[M_1] vector[N_1] z_1; // standardized group-level effects
vector<lower=0>[M_2] sd_2; // group-level standard deviations
array[M_2] vector[N_2] z_2; // standardized group-level effects
vector[N_1] r_1_1; // actual group-level effects
vector[N_2] r_2_1; // actual group-level effects
real lprior = 0; // prior contributions to the log posterior
lprior += student_t_lpdf(Intercept | 3, 4.6, 2.5);
lprior += student_t_lpdf(sdgp_1 | 3, 0, 2.5)
lprior += student_t_lpdf(sigma | 3, 0, 2.5)
lprior += student_t_lpdf(sd_1 | 3, 0, 2.5)
lprior += student_t_lpdf(sd_2 | 3, 0, 2.5)
vector[Nsubgp_1] gp_pred_1 = gp(Xgp_1, sdgp_1[1], lscale_1[1], zgp_1);
// add more terms to the linear predictor
target += lognormal_lpdf(Y | mu, sigma);