When using normal likelihood,how to define sigma

Questions v5
1

I want to ask a silly question. When using a normal distribution as likelihood function, fit a curve, the values on the curve is denoted as obs like this

with pm.Model() as background_model:
    amp = pm.Uniform('amp',lower=-100,upper=-0)
    rate = pm.Uniform('rate', lower=0,upper= 200)
    t_start = pm.Uniform('t_start',lower=0,upper=500)
    b1=pm.Uniform('b1',lower=0,upper=4)    
    y_observed=pm.Normal(
        "y_observed",
        mu=amp*np.exp(-(t_back-t_start)/rate)+b1,
        sigma=noise,
        observed=obs,
        )
    output = pm.Deterministic('output',amp*np.exp(-(t_back-t_start)/rate)+b1)
    prior = pm.sample_prior_predictive()
    posterior_b = pm.sample(draws = Samples, target_accept = 0.9,chains=4,cores=1)
    posterior_only_background = pm.sample_posterior_predictive(posterior_b)

how do we determine the selection of the parameter sigma? What factors are usually related to it?

2

The equation

amp*np.exp(-(t_back-t_start)/rate)+b1

describes the average trend which relate the predictor variables (t_back I assume in this case) to outcome variables which is the observed. However the data will not exactly match this trend, more likely it will be distributed with some noise around this data. sigma in this particular model represents how much noise you allow. Normally, if there is enough data this is also left as a parameter to be estimated with something like this

noise = pm.HalfNormal("noise", 5)

Here I wrote 5 arbitrarily, it should be chosen to represent the natural scale of data in your case but this is still much more flexible than trying to guess an exact sigma for your problem. If you standardise your data in some way, this is generally close to 1. You might get a bit more in depth idea on the choice of particular priors for sigma here:

Prior distributions for variance parameters in hierarchical models

1 Like