Accessing target in transformed parameters: stanc2 vs stanc3

Developers
1

stanc2 (the C++ Stan-to-C++ compiler) allowed the following model:

functions { 
  void test_lp(real r) {
     target += normal_lpdf(r | 0, 1);
  }
}
parameters { 
  real y; 
}
transformed parameters { 
  test_lp(5.0);  
}
model {
  y ~ normal(0, 1); 
}

while

parameters { 
  real y; 
}
transformed parameters { 
  target += normal_lpdf(r | 0, 1);
}
model {
  y ~ normal(0, 1); 
}

was not allowed (with the target is not allowed in transformed parameters error). In stanc3 none of these are allowed.
To me this seems like a stanc2 bug, just wanted to confirm.

This came up in a recent stanc3 issue.

cc: @seantalts @nhuurre @rybern et al

EDIT:

The example from the issue is from https://rstudio-pubs-static.s3.amazonaws.com/435225_07b4ab5afa824342a4680c9fb2de6098.html

2 Likes
2

You might want to tag @Bob_Carpenter, he’ll know if this is a stanc2 bug.

2 Likes
3

The allowance of a distribution call to lp through a function seems like a bug. However, it brings up the case that there should be an allowed lp increment for log-abs-det-jacobian adjustments in transformed parameters. I’m thinking a special function like increment_target() such as

parameters {
  real<lower=0> y_inv;
}
transformed parameters {
  real<lower=0> y;
  y = 1 / y_inv;  // change variables
  increment_target(-2 * log(y_inv));  //  Jacobian adjustment;
}
model {
  y ~ gamma(2,4);
}

It seems unnecessary in this situation but there are more complicated situations where I’m doing an adjustment in a function in transformed parameters and then I have to figure out a way to pass it to a dummy variable. This is difficult especially when one has matrices as outputs. So now I have an array of 2 matrices where matrix is my output and the other is mostly zeros except for the jacobian adjustment. Such as

// what's done today
matrix[] f (int a){
matrix[a, a] out[2];
...
out[1] = my_output;
out[2, 1, 1] = log(abs(jacobian(y));
return out;
}
// proposal
matrix f (int a){
matrix[a, a] out;
...
out = my_output;
jac_adj = log(abs(jacobian(y));
increment_target(jac_adj);
return out;
}
3 Likes