The SPInS program

A module which contains the main program for SPInS as well as various classes which intervene when calculating the priors and likelihood function:

  • Distribution: a class which represents a probability distribution

  • Prior_list: a class with a list of priors

  • Likelihood: a class used to represent the likelihood function

  • Probability: a class which groups the priors and likelihood function together

This module relies on the emcee package to apply an MCMC algorithm which will return a representative sample of models for a given set of seismic an classic constraints.

Warning

In various places in this module, for instance in the Prior_list and Likelihood classes, various methods return what is described as a \(\chi^2\) value. Technically, these are not \(\chi^2\) values, but rather \(-\chi^2/2\), i.e. the argument of the exponential function which intervenes in the Gaussian probability distribution.

class SPInS.Distribution(_type, _values)[source]

A class which represents a probability distribution, and can yield its value for a given input parameter, or provide a random realisation.

Note

Derived from a class originally written by G. Davies.

Parameters:

  • _type (string) – type of probability function (current options include “Uniform”, “Gaussian”, “Truncated_gaussian”, “IMF1”, “IMF2”, “Above”, “Below”, “Uninformative”, or “Casa2011”)

  • _values (list of floats) – list of parameters relevant to the probability function

property error_bar

Returns an error bar based on the distribution. This does not necessarily correspond to the one-sigma value but rather to what is the most convenient value.

Returns:

the error bar

Return type:

float

property mean

Returns the mean value of the probability distribution.

Returns:

the mean value of the probability distribution

Return type:

float

property nparams

Return the number of relevant parameters for a given distribution.

Returns:

the number of relevant parameters

Return type:

int

print_me()[source]

Print type and parameters of probability distribution.

re_centre(value)[source]

Re-centre the probability distribution around the input value.

Parameters:

value (float) – new value around which to centre the distribution

re_normalise(value)[source]

Re-normalise the probability distribution so that its characteristic width corresponds to the input value.

Parameters:

value (float) – new value around for the chacteristic width

realisation(size=None)[source]

Return random values which statistically follow the probability distribution.

Parameters:

size (int or tuple of ints) – shape of random variates (NOTE: None corresponds to 1 value)

Returns:

a set of random realisations

Return type:

float

to_string()[source]

Produce nice string representation of the distribution.

Returns:

nice string representation of the distribution

Return type:

string

type

Type of probability function (“Uniform”, “Gaussian”, “Truncated_gaussian”, “IMF1”, “IMF2”, “Above”, “Below”, “Uninformative”, or “Casa2011”)

values

List of parameters relevant to probability function

class SPInS.Likelihood[source]

A class which describes the likelihood function and allows users to evaluate it.

add_constraint(a_constraint)[source]

Add a supplementary constraint to the list of constraints.

Parameters:

constraint ((string, Distribution)) – supplementary constraint

add_constraints_from_prior(prior_dict)[source]

Extract constraints from a priors dictionary with the name of the parameters as keys, and distributions as values.

Parameters:

prior_dict ({str : Distribution} dictionary) – dictionary with parameters as keys and distributions as values

constraints

List of constraints which intervene in the likelihood function.

evaluate(my_model)[source]

Calculate ln of likelihood function (i.e. a \(\chi^2\) value multiplied by -1/2) for a given model.

Parameters:

my_model (model.Model) – model for which the \(\chi^2\) value (multiplied by -1/2) is being calculated

Returns:

the \(\chi^2\) value (multiplied by -1/2), i.e. ln(likelihood)

Return type:

float

Note

  • This avoids model interpolation and can be used to gain time.

  • This includes a correction to account for the use of a non-dimensional age parameter

read_constraints(filename)[source]

Read a file with constraints.

Parameters:

filename (string) – name of file with constraints

class SPInS.Prior_list(prior_dict)[source]

A class which contains a list of priors as well as convenient methods for adding priors and for evaluating them.

init_priors(prior_dict)[source]

Initialise the priors using a python dictionary with the name of the parameters as keys, and distributions as values.

Parameters:

prior_dict ({str : Distribution} dictionary) – dictionary with parameters as keys and distributions as values

print_me()[source]

Print contents of prior list, including the names of the parameters, the functions applied to these and the distributions.

priors

A list of priors or tight ball distributions. Each contain the following items (in this order - the numbers correspond to the indices):

  1. the name of the parameter (including functions applied to the parameter).

  2. a list of functions (not strings) which when applied to the corresponding MCMC parameter produces the parameter in item 0.

  3. the derivatives of the functions in item 1. These derivatives are necessary for correctly calculating the probability in terms of the corresponding MCMC parameter rather than the parameter given in item 0.

  4. the inverted list of functions. When applied to the parameter in item 0, it produces one of the MCMC parameters.

  5. a probability distribution.

  6. a multiplicative coefficient on ln(prior) to account for common parameters when fitting multiple stars.

realisation(size=None)[source]

Return an array with realisations for each prior. The last dimension will correspond to the different priors.

Parameters:

size (int or tuple of ints) – shape of random variates (for each prior)

Returns:

a set of realisations

Return type:

numpy float array

Note

Realisations of priors do not take into account the multiplicative coefficient of common parameters for multiple stars as a single realisation will be used for all of the stars.

class SPInS.Probability(_priors, _likelihoods)[source]

A class which combines the priors and likelihood function, and allows the the user to evalute ln of the product of these.

Parameters:

  • _priors (Prior_list) – input set of priors

  • _likelihood (Likelihood) – input likelihood function

evaluate(my_model, star_index)[source]

Evalulate the ln of the product of the priors and likelihood function, i.e. the probability, for a given model and a given star, to within an additive constant.

Parameters:

  • my_model (model.Model) – input model

  • star_index (int) – index of star for which we would like to evaluate the probability

Returns:

the ln of the probability

Return type:

float

Note

This avoids model interpolation and can be used to gain time.

is_outside(params)[source]

Test to see if the given set of parameters lies outside the grid of models. This is done by evaluating the probability and seeing if the result indicates this.

Parameters:

params (array-like) – input set of parameters

Returns:

True if the set of parameters corresponds to a point outside the grid.

Return type:

boolean

like_function(params)[source]

Evaluate the combined likelihood function for all of the stars.

Parameters:

params (array-like) – input set of parameters

Returns:

the ln of the likelihood

Return type:

float

likelihoods

The likelihood functions (one for each star).

prior_function(params)[source]

Evaluate the combined prior function for all of the stars.

Parameters:

params (array-like) – input set of parameters

Returns:

the ln of the prior

Return type:

float

priors

The set of priors (the same priors apply to each star).

SPInS.autocorr_time = []

integrated autocorrelelation time (this is useful for testing convergence)

SPInS.best_MCMC_params = None

best parameters from the MCMC run

SPInS.best_MCMC_result = -1e+300

best ln(probability) result from the MCMC run

SPInS.best_grid_params = None

best parameters after scanning the grid (once per star)

SPInS.best_grid_result = -1e+300

ln(probability) obtained for best_grid_params

SPInS.check_configuration()[source]

Test the version of the EMCEE package and set the mc2 variable accordingly.

Test the values of the variables in check_configuration to make sure they’re acceptable. If an unacceptable value is found, then this will stop SPInS and explain what variable has an erroneous value.

SPInS.dont_run_emcee(p0)[source]

Dummy subroutine which simply returns the initial distribution without doing any emcee iterations.

Parameters:

p0 (np.array) – the initial set of walkers

Returns:

the dummy sampler for the MCMC run, dummy array with walker percentiles

Return type:

dummy sampler object, np.array

SPInS.find_a_blob(params)[source]

Find a blob (i.e. supplementary output parameters) for a given set of parameters (for one model). The blob also includes the log(P) value as a first entry.

Parameters:

params (array-like) – input set of parameters

Returns:

list of supplementary output parameters

Return type:

list of floats

SPInS.find_best_model(star_index)[source]

Scan through grid of models to find “best” model for a given probability function for an individual star (i.e. the product of priors and a likelihood function).

Parameters:

star_index (int) – index of star being fitted

Returns:

the best model for the star, the number of accepted parameters (i.e. models), the number of rejected parameters (i.e. models)

Return type:

float array, int, int

SPInS.find_best_model_in_track(track_star)[source]

Scan through an evolutionary track to find “best” model for prob, the probability function (i.e. the product of priors and a likelihood function) for a given observed star.

Parameters:

track_star ((int, int)) – number of the evolutionary track and index of the observed star.

Returns:

the ln(probability) value, and the “best” model, the parameters of the accepted models, and those of the rejected models

Return type:

float, model.Model, tuple of float arrays, tuple of float arrays

SPInS.find_best_models()[source]

Find best parameters for each model treated individually then combine parameters to obtain a global set of optimal parameters. If there are common parameters, these are averaged between the different solutions.

Returns:

the age ranges for the evolutionary tracks with the best models, and the common age range

Return type:

1D numpy float array, float

SPInS.find_blobs(samples)[source]

Find blobs (i.e. supplementary output parameters) from a set of samples (i.e. for multiple models).

Parameters:

samples (list/array of array-like) – input set of samples

Returns:

set of supplementary output parameters

Return type:

np.array

SPInS.find_output_params(my_model)[source]

This finds the supplementary output params as specified by the config.output_params variable.

Parameters:

my_model (np.array) – a model for which we wish to obtain the supplementary output parameters

Returns:

a list with the supplementary parameters

Return type:

list of floats

SPInS.find_param_map()[source]

Finding mapping between input set of parameters and parameters for each star, taking into account the fact that some of the parameters may be common parameters (such as the age for a binary system/cluster)

Returns:

the mapping from the reduced to the extended set of parameters

Return type:

a 2D int array

SPInS.find_weights_hist(array, nbins)[source]

This finds appropriate weights for created a normalised histogram (thus avoiding the use of the keywords “normed” or “density”).

Parameters:

  • array (numpy float array) – the array to be normalised

  • nbins (int) – number of bins involved in the histogram

Returns:

the weight array

Return type:

numpy float array

SPInS.grid = None

grid of models

SPInS.grid_params_MCMC = ()

parameters used in the MCMC run (dimension = ndims)

SPInS.grid_params_MCMC_fundamental = ()

fundamental parameters used in the MCMC run

SPInS.grid_params_MCMC_index = ()

array with indices of the parameters used in the MCMC run

SPInS.init_walkers()[source]

Initialise the walkers used in emcee.

Returns:

array of starting parameters

Return type:

np.array

SPInS.initial_distributions = None

Prior_list type object which stores the distributions for the initial tight ball used to initialise the MCMC run (provided tight_ball`==``True`).

SPInS.interpolation_tests(filename)[source]

Carry out various interpolation tests and write results in binary format to file.

Parameters:

filename (string) – name of file in which to write test results.

Note

The contents of this file may be plotted using methods from plot_interpolation_test.py.

SPInS.labels = None

labels for plots

SPInS.like_function(params)[source]

Wrapper function which allows prob.like_function to be pickled and thus used in a parallelised context.

Parameters:

params (array-like) – input set of parameters

Returns:

the ln of the likelihood

Return type:

float

SPInS.load_binary_data(filename)[source]

Read a binary file with a grid of models.

Parameters:

filename (string) – name of file with grid in binary format

Returns:

the grid of models

Return type:

model.Model_grid

SPInS.log0 = -1e+300

a large negative value used to represent ln(0)

SPInS.mc2 = None

True if emcee’s version is 2 or prior

SPInS.my_map = None

pointer to the map function (either the parallel or sequential versions)

SPInS.naccepted_parameters = []

number of sets of parameters associated with accepted models for each star

SPInS.ncommon = 0

number of common parameters when fitting multiple stars

SPInS.ndims = 0

number of dimensions for MCMC parameters per star

SPInS.ndims_total = 0

total number of dimensions for MCMC parameters

SPInS.nrejected_parameters = []

number of sets of parameters associated with rejected models for each star

SPInS.nstars = 0

number of stars being fitted simultaneously

SPInS.output_folder = None

folder in which to write the results

SPInS.output_params_funcs = ()

tuple with functions applied to supplementary output parameters

SPInS.output_params_index = ()

array with indices of supplementary output parameters

SPInS.plot_distrib_iter(percentiles, labels, names, folder)[source]

Plot individual distribution of walkers as a function of iterations.

Parameters:

  • percentiles (np.array) – array with percentiles from emcee run

  • labels (list of strings) – labels for the different variables in parameters space

  • names (list of strings) – names of the different variables in parameters space (used for filenames)

  • folder (string) – specify name of file in which to save plots of walkers.

Warning

This method must be applied before the samples are reshaped, and information on individual walkers lost.

SPInS.plot_histograms(samples, names, fancy_names, truths=None)[source]

Plot a histogram based on a set of samples.

Parameters:

  • samples (np.array) – samples form the emcee run

  • names (list of strings) – names of the quantities represented by the samples. This will be used when naming the file with the histogram

  • fancy_names (list of strings) – name of the quantities represented by the samples. This will be used as the x-axis label in the histogram.

  • truths (list of floats) – reference values (typically the true values or some other important values) to be added to the histograms as a vertical line

SPInS.plot_walkers(samples, labels, filename, nw=3)[source]

Plot individual walkers.

Parameters:

  • samples (np.array) – samples from the emcee run

  • labels (list of strings) – labels for the different dimensions in parameters space

  • filename (string) – specify name of file in which to save plots of walkers.

  • nw (int) – number of walkers to be plotted

Warning

This method must be applied before the samples are reshaped, and information on individual walkers lost.

SPInS.pool = None

pool from which to carry out parallel computations

SPInS.prior_function(params)[source]

Wrapper function which allows prob.prior_function to be pickled and thus used in a parallelised context.

Parameters:

params (array-like) – input set of parameters

Returns:

the ln of the prior

Return type:

float

SPInS.prob = None

Probability type object that represents the probability function which includes the likelihood and priors

SPInS.run_emcee(p0)[source]

Run the emcee program.

Parameters:

p0 (np.array) – the initial set of walkers

Returns:

the emcee sampler for the MCMC run, array with walker percentiles as a function of iteration number

Return type:

emcee sampler object, np.array

SPInS.star_names = []

the names of the stars being fitted

SPInS.statistical_params = None

statistical parameters obtained from averaging the samples

SPInS.statistical_result = -1e+300

ln(probability) result for statistical_params

SPInS.string_to_title(string)[source]

Create fancy title from string.

Parameters:

string (string) – string from which the title is created.

Returns:

the fancy string title

Return type:

string

SPInS.swap_dimensions(array)[source]

Swaps the two first dimensions of an array. This is useful for handling the different conventions used to store the samples in emcee3 and ptemcee.

Parameters:

array (np.array) – input array

Returns:

array with swapped dimensions

Return type:

np.array

SPInS.threshold = -1e+290

threshold for “accepted” models. Needs to be greater than log0

SPInS.write_binary_data(infile, outfile)[source]

Read an ascii file with a grid of models, and write corresponding binary file.

Parameters:

  • infile (string) – input ascii file name

  • outfile (string) – output binary file name

SPInS.write_model(my_model, my_params, my_result, overall_result, model_name)[source]

Write text file with caracteristics of input model.

Parameters:

  • my_model (model.Model) – model for which we’re writing a text file

  • my_params (array-like) – parameters of the model

  • my_result (float) – ln(P) value obtained for the model

  • overall_result (float) – total ln(P) value obtained for the set of models

  • model_name (string) – name used to describe this model. This is also used when naming the text file.

SPInS.write_percentiles(filename, labels, samples)[source]

Write percentiles based on a sequence of realisations to a file. The results include:

+/- 2 sigma error bars (using the 2.5th and 97.5th percentiles) +/- 1 sigma error bars (using the 16th and 84th percentiles) the median value (i.e. the 50th percentile)

Parameters:

  • filename (string) – name of file in which to write the percentiles

  • labels (list of strings) – names of relevant variables

  • samples (np.array) – samples for which statistical properties are calculated

SPInS.write_readme(filename, elapsed_time)[source]

Write parameters relevant to this MCMC run.

Parameters:

  • filename (string) – name of file in which to write the statistical properties

  • elapsed_time (float) – time that elapsed during the MCMC run

SPInS.write_samples(filename, labels, samples)[source]

Write raw samples to a file.

Parameters:

  • filename (string) – name of file in which to write the samples

  • labels (list of strings) – names of relevant variables (used to write a header)

  • samples (array-like) – samples for which statistical properties are calculated

SPInS.write_statistics(filename, labels, samples)[source]

Write statistical properties based on a sequence of realisations to a file. The results include:

  • average values for each variable (statistical mean)

  • error bars for each variable (standard mean deviation)

  • correlation matrix between the different variables

Parameters:

  • filename (string) – name of file in which to write the statistical properties

  • labels (list of strings) – names of relevant variables

  • samples (np.array) – samples for which statistical properties are calculated