Using distributions

uravu isn’t limited to using normally distributed ordinate values. Any distribution of ordinate values can be used, as uravu will perform a Gaussian kernel density estimation on the samples to determine a description for the distribution.

This is most easily shown in action, imagine we have some experimental data, that is distributed with a skew normal distribution, rather than the typical normal distribution. So the values for a particular -value take the (peculiar) shape,

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import skewnorm
from uravu.distribution import Distribution
from uravu.relationship import Relationship
from uravu.utils import straight_line
from uravu import plotting

np.random.seed(2)

y = skewnorm(10, 0, 0.1)

plt.hist(y.rvs(size=5000), bins=25)
plt.show()

Let’s build some synthetic data, collected by sampling some skew normal distribution across a series of values of .

x = np.linspace(1, 100, 10)
Y = []
for i in x:
    Y.append(Distribution(skewnorm.rvs(10, i*3.14, i*0.5, size=5000)+(1*np.random.randn())))

Note that the sample, in this case a series of random values from the distribution, are passed to the uravu.distribution.Distribution object and stored as a list of Distribution objects.

This list is passed to the Relationship class as shown below (note the ordinate_error keyword argument is no longer used as we are describing the distribution directly in the Y variable),

r = Relationship(straight_line, x, Y, bounds=((-10, 10), (-10, 10)))
r.max_likelihood('diff_evo')
r.variable_medians

array([ 3.37733187, -0.02284246])

plotting.plot_relationship(r)
plt.show()

It is then possible to use the standard sampling methods to investigate the distribution of the model parameters, here the model is a simple straight line relationship.

r.mcmc()

100%|██████████| 1000/1000 [01:04<00:00, 15.45it/s]

Above Markov chain Monte Carlo is used to sample the distribution of the gradient and intercept of the straight line, given the uncertainties in the ordinate values (from the distributions).

These distributions can be visualised with the plot_distribution from the uravu.plotting library.

plotting.plot_distribution(r.variables[0])
plt.show()

plotting.plot_distribution(r.variables[1])
plt.show()

We can also see how these distributions affect the agreement with the data by plotting the relationship.

plotting.plot_relationship(r)
plt.show()