3
$\begingroup$

I am currently doing an analysis that involves fitting a model to a 1D graph.

Following the example on the emcee documentation, I started with Maximum likelihood estimation and am now looking at using MCMC and Bayesian statistics to constrain nuisance parameters in my function.

This is all fine.

The part that I am on before being able to use MCMC is that I need to have a posterior function, which is given by the product of my likelihood function and a prior. The only part I am missing is the prior.

Based on doing background reading into this, I am starting just with the very simple example of using a flat, uniform prior, like they have in the emcee documentation (fitting a model to data tutorial is what I am referencing).

Now, I understand that for a flat uniform prior, every value within your defined range has equal probability of being chosen assigned to it and that for such a function, that probability is given by p(x) = 1/(b-a), where a, b are the min and max values of your chosen range and that, outside your range the probability is zero.

So, in the emcee example they are trying to constrain the slope of a line which they define by parameter, m. They state the range is, [a,b] = [-5.0,0.5] and therefore the probability of each parameter within that range is 1/(0.5--5.0) = 1/5.5. This all makes sense.

However, in their code block demonstrating this, which I will paste below from their page, where they then take the log of the prior and therefore, log of the probability within the function, they return the values 0, for probabilities of values within the range of the uniform function and -infinity, for values outside the range of the function.

def log_prior(theta):
m, b, log_f = theta
if -5.0 < m < 0.5 and 0.0 < b < 10.0 and -10.0 < log_f < 1.0:
    return 0.0
return -np.inf

Where 0.0 is log of the probability within the function range and - infinity is log of the probability outside the function range.

For the -infinity value, this makes sense and based on some reading, I understand that while log(0) is undefined, it is common to assign a value of -infinity to it. However, a value of zero for values within the range would imply that they took the log of 1, for values within the range.

What I don't understand is why they would take the log of 1 and not the log of 1/5.5, which is the probability assigned to each possible value in the function range. Taking the log of 1 would be the log of the total probability within the function but I can't seem to find any explanations of their reasoning anywhere.

Therefore, I was wondering the correct way to do this and whether, you should always take the log of 1 for a uniform flat prior? and if so, what is the point of determining the probability of values within the function (i.e. the 1/5.5 value in this example), if you don't use it in your prior?

I have attached the link to the emcee docs page I am using below.

Following on from this, while they have given the example for the probability of the m parameter they are trying to constrain from the fit, they are actually trying to constrain three nuisance parameters (total: m , b and f as defined in the docs). If the best approach is to return a probability of 1/(b-a), how would you do when combining multiple parameter constraints? is it just the same as a joint probability distribution when you take the product of the probabilities (or sum if you're doing log)? this is what makes sense to me but I cannot find examples of anyone explicitly giving the value for multiple constraints so I ask incase there is another approach.

Many thanks,

https://emcee.readthedocs.io/en/stable/tutorials/line/

An additional question, if anyone should know is, when doing log prior and if your parameters are outside the function range, therefore assigned a value of - infinity, when going to add this to your likelihood function to get the posterior (as we are in log space), wouldn't you always get a value of -infinity? because your likelihood value is small in comparison ? Maybe this is a naive question based on my currently limited understanding but any guidance would be appreciated.

$\endgroup$
7
  • $\begingroup$ The posterior is calculated up to a constant; see the documentation. The log(5.5) is incorporated into the constant. $\endgroup$
    – dipetkov
    Feb 22 at 9:54
  • $\begingroup$ This has been explained many times on CV but one more explanation doesn't hurt. Some related threads: 1, 2, 3. $\endgroup$
    – dipetkov
    Feb 22 at 9:58
  • 1
    $\begingroup$ Thank you for your response - I did read the part of the "calculated up to a constant" but I (mis)interpreted this as referring to the removal of the Bayesian evidence term as opposed to the 1/5.5. vs 1 factor. I appreciate you tagging relevant threads $\endgroup$
    – shram
    Feb 22 at 10:14
  • $\begingroup$ That's true -- that constant includes the evidence but it can include other components as convenient for computational purposes (as long as those are constant in the parameters). $\endgroup$
    – dipetkov
    Feb 22 at 10:16
  • $\begingroup$ @dipetkov that all makes sense, thank you for explaining and clarifying for me! $\endgroup$
    – shram
    Feb 22 at 10:29

1 Answer 1

1
$\begingroup$

There is a more formal answer to your question, but I won't give it here. In broader terms, though:

What I don't understand is why they would take the log of 1 and not the log of 1/5.5, which is the probability assigned to each possible value in the function range. Taking the log of 1 would be the log of the total probability within the function but I can't seem to find any explanations of their reasoning anywhere.

It doesn't (usually) matter.

When you're using MCMC to estimate the posterior distribution, it's not actually necessary to have a function that outputs the value of the posterior density, $P(\theta | x)$. All you need is to output a value that's proportional to this, usually written as $\propto P(\theta | x)$. This is good, because to actually calculate $P(\theta | x)$ you need to also calculate the normalising constant, and this is expensive to do, but because it's a constant, it doesn't affect the proportional value.

For this reason, it doesn't actually matter if you shift the prior density up or down by a constant amount (e.g. $1$ vs $\frac{1}{5.5}$), since everything remains in proportion.

An additional question, if anyone should know is, when doing log prior and if your parameters are outside the function range, therefore assigned a value of - infinity, when going to add this to your likelihood function to get the posterior (as we are in log space), wouldn't you always get a value of -infinity? because your likelihood value is small in comparison ? Maybe this is a naive question based on my currently limited understanding but any guidance would be appreciated.

Correct. If the prior says a value of $\theta$ has a probability of 0 (it's impossible), it should also have a probability of 0 (also impossible) in the posterior, since no amount of data can make the impossible possible.

Edit: It's also worth looking at what people actually mean when they say "up to a constant" when talking about the posterior, e.g. here.

$\endgroup$
1
  • $\begingroup$ Thank you so much for your response, that all makes a lot of sense and I appreciate you taking the time to answer! $\endgroup$
    – shram
    Feb 22 at 10:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.