Constraint transforms and default priors

Modeling
,
1

I am trying to understand how to specify priors for variables that are subject to a constraint transform.

In particular, what I really want is to impose a U(0,\infty) prior on a variable.

My understanding is that when I specify

parameters {
     real<lower=0> C;
}

this transforms to a new variable which is effectively \ln C, which has Jacobian dC/d\ln C=C=\exp(\ln C) and (if I don’t explicitly give a prior to C in the model block) effectively imposes the unconstrained prior \ln C\sim U(-\infty,+\infty) (i.e., P(C)\propto1/C. First, is this correct?

If so, then, second, it seems that I have two choices to impose my desired prior. I could pick a suitably large value of Cmax and then use

model {
    C ~ uniform(0,Cmax);
}

In practice, this is fine, since there is [almost?] always enough information to find a good Cmax.
Erratum: I am pretty sure the following is incorrect – see below.

However, it occurs to me that I should be able to get the actual prior that I want by just “removing” the Jacobian of the constraint transform (and imposing no other prior on C):

model {
    target += -log(C); // note the minus sign
}

However, when I’ve tried this, the sampling hasn’t gone very well. Am I missing something?

Update: it now seems clear that I’ve got a sign error in the above. I really should be doing

model {
    target += log(C); // note the plus sign!
}

in order to actually impose the Jacobian transformation to Uniform.

Thanks, @Funko_Unko for inadvertently helping me figure this out. It would still be good to get more knowledgeable input to confirm whether this is the right thing to do (or otherwise).

Update 2 not sure this is right. See below.

1 Like
2

As far as I know, not specifying anything in the model block will yield your desired (improper?) prior.

For e.g. C~lognormal(0,1) actually yielding a lognormal prior the “previous” prior has to be uniform (I think).

3

I guess my understanding was that if you specify a prior on C, then and only then is the Jacobian applied in order to get the correct distribution – but otherwise not, hence giving the default wide-open prior \ln C\sim U(-\infty, +\infty), but it’s clear that there is some ambiguity here.

4

Oh, right that might actually be the case. Sorry for introducing additional confusion!

Edit:

Adding more flames to the confusion fire:

I think the your edited update works because you don’t have to “remove” the Jacobian adjustment because as you pointed out before, it only gets added when you write something like C~lognormal(0,1). Instead, you have to add the Jacobian adjustment that would have been added if you were actually able to specify C~uniform(0,inf).

5

However, I think even my above “update” is incorrect. By that reasoning, the following should be equivalent:

parameters {
  real<lower=0> C;
}

transformed parameters {
  real lnC = log(C);
}

model {
  target += lnC; // Is this already included implicitly?
}

and

parameters {
  real lnC;
}

transformed parameters {
  real C = exp(C);
}

model {
  target += lnC; // needed to convert to uniform in C>0.
}

However, this does not appear to be the case; the actual models (which include data and other parameters) agree when I remove the Jacobian addition for the model with parameter real<lower=0> C;.

So then it appears that the default prior is uniform in the real parameter, even with a lower limit — the Jacobian is added in by default no matter what (as @Funko_Unko said in their original response).

Can an expert comment (and/or point to appropriate documentation)?

What am I missing?

6

Yes. In other words, the Jacobian ensures that any prior (including improper flat prior) behaves as if it was applied to the constrained scale. Mike’s essay/case study Probability Theory (For Scientists and Engineers) provides the theoretical underpinning for this. Also, the Jacobian is always applied (there is an exception for the optimize method where users can switch all the jacobians on/off, but this does not depend on putting priors on parameters).

As a general rule, it is safe to assume Stan does no magic with your program - your program is just a function from parameters to target density, nothing more, nothing less. And Stan has very limited abilities to inspect the internal structure of your program. Due to Stan’s flexibility it is actually impossible to determine stuff like what statement is a “prior” in full generality - Stan is Turing complete, so you can’t even reliably determine if a statement gets exectued at all!

The relevant documentation is: 10 Constraint Transforms | Stan Reference Manual

Hope that clarifies a bit!

1 Like