Is there an issue with mixing "~" and "target +="

Modeling
1

I’ve begun trying to make sense of some Stan code (which I think is supposed to implement parameter uncertainty inflation), which has these two lines:

y[i] ~ normal(mu_M[i], uy[i]);
target += -0.5*( V_T[i]^0.5 -fabs(mu_M[i]-y[i]) )^2 
              / uy[i]^2 ;

I notice that this combines both a sampling statement “~”, which drops constant terms in the log likelihood, and “target +=” which doesn’t, and I’m not sure if mixing the two is an issue or not. Is this a red flag, or should I not worry about it?

FYI, the full code is in these two files:

2

When in doubt appeal to the fact that the ~ notation is just syntactic sugar for a target += call. So your model is equivalent to

target += normal_lpdf(y[i] | mu_M[i], uy[i]);
target += -0.5*( V_T[i]^0.5 -fabs(mu_M[i]-y[i]) )^2 
              / uy[i]^2 ;

which may be easier to reason about.

Does that make sense? Not really given you’re defining a Gaussian distribution for y[i] twice.

3

Except that the Stan language manual makes clear that’s almost but not quite true. That’s why I was asking about it. I wanted to know if the subtle difference between the two would make a material difference.

(Also, given how V_T is defined, that’s not defining a Gaussian distribution for y[i] twice, but that’s a side note.)

4

They are equivalent for all of the algorithms currently supported in Stan. Any difference is in the normalizing constants which is superficial.

If you factor the quadratic then one of the three terms you get is the unnormalized version normal_lpdf(y[i] | mu_M[i], uy[i]) so with both lines you are adding that distribution twice, effectively halving the standard deviations.

1 Like
5

Yeah, I see what you mean. Not sure why they’re doing that, but that’s probably outside the purview of the Stan forums.

(FWIW, the paper related to those Stan files (and others) is at this DOI: 10.1002/aic.15781, if you’re interested. I’m still working on making sense of it and its associated codes.)

6

Hi,

If we use target += normal_lpdf(y[i] | mu_M[i], uy[i]) then y[i] has to be an observed data, right, I mean y[i] cannot be a missing data or any parameters?

Kind regards,
Tran.

7

No. Quite literally,

alpha ~ foo(beta1, ..., betaN);

is equivlaent to

target += foo_lpdf(alpha | beta1, ..., betaN) - const;

where the constant doesn’t depend on any of the alpha or beta that are constant. For example if you do this,

alpha ~ normal(0, 1);

where alpha is a parameter and 0 and 1 are constants, the result is equivalent (up to constant terms) to writing

target += -0.5 * alpha^2;
1 Like