Uniform priors in easyCore#

Bayesian methods are becoming more and more popular in the analysis of model-dependent data, in particular in neutron scattering. Thinking about data probabilistically, helps towards an understanding of Bayesian modelling. Bayesian modelling is about how to include prior knowledge about the system or parameters under study in the analysis.

The same data as previously will be used again.

from easyCore import np

np.random.seed(123)

a_true = -0.9594
b_true = 7.294
c_true = 3.102

N = 25
x = np.linspace(0, 10, N)
yerr = 1 + 1 * np.random.rand(N)
y = a_true * x ** 2 + b_true * x + c_true
y += np.abs(y) * 0.1 * np.random.randn(N)

The most straightforward way to include prior knowledge in some analysis is through the use of bounded parameters. Bounded parameters cannot have values less than some lower bound or greater than some upper bound, as the probability of the parameters having these values is zero. For example, if the parameter b from the quadratic model has bounds of 0 and 10, then there is an equal probability that the value of b can be anything in between 0 and 10, and a probability of 0 outside those bounds, i.e., it has a uniform prior probability distribution.

../_images/uniform.png

Fig. 6 A uniform distribution (blue line), from 0 to 10.#

Information about bounds can be included in easyCore parameters as min and max values.

from easyCore.Objects.Variable import Parameter

a = Parameter(name='a', value=a_true, fixed=False, min=-5.0, max=0.5)
b = Parameter(name='b', value=b_true, fixed=False, min=0, max=10)
c = Parameter(name='c', value=c_true, fixed=False, min=-20, max=50)

We can then perform the analysis in the same fashion as previously, however, this time the bounds will be respected.

from easyCore.Objects.ObjectClasses import BaseObj
from easyCore.Fitting.Fitting import Fitter

def math_model(x, *args, **kwargs):
    return a.raw_value * x ** 2 + b.raw_value * x + c.raw_value

quad = BaseObj(name='quad', a=a, b=b, c=c)
f = Fitter(quad, math_model)

res = f.fit(x=x, y=y, weights=yerr)

a, b, c

(<Parameter 'a': -0.92+/-0.04, bounds=[-5.0:0.5]>,
 <Parameter 'b': 7.1+/-0.4, bounds=[0:10]>,
 <Parameter 'c': 2.1+/-0.8, bounds=[-20:50]>)

Currently easyCore only supports uniform prior probability distributions, however, there are plans to extend this in the future.