The AIMS
program
A module which contains the main program for AIMS as well as various classes which intervene when calculating the priors and likelihood function:
Distribution
: a class which represents a probability distributionPrior_list
: a class with a list of priorsMode
: a class used to represent observed modesCombination
: a class used to represent a linear frequency combinationCombination_function
: a class used to represent functions of frequency combinationsLikelihood
: a class used to represent the likelihood functionProbability
: 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 AIMS.Combination[source]
A class which contains indices and coefficients which represent a linear combination of frequencies.
- add_coeff(j, coeff)[source]
Append the given index and coefficient to the list of indices and coefficients.
- Parameters
j (int) – index of the mode
coeff (float) – coefficient used in the frequency combination
- coeff
Coefficients of the linear combination of frequencies.
- find_values(freq)[source]
Find the value associated with this frequency combination using the input frequencies.
- Parameters
freq (float array like) – list of frequencies
- index
Indices in of the linear combination of frequencies.
- value
Value of the linear combination of frequencies.
- class AIMS.Combination_function(_name, _function, _gradient_function, _offset=0.0)[source]
A class which applies functions to linear frequency combinations.
- add_combination(combination)[source]
Append a frequency combination to this combination function.
- Parameters
combination (
Combination
) – a linear frequency combination
- combinations
The frequency combinations which intervene in the function.
- find_values(freq)[source]
Find the value and gradient of this frequency combination function for the input set of frequencies.
- Parameters
freq (float array like) – list of frequencies
- function
The function which is applied to the frequency combinations.
- gradient
The gradient of the frequency combination function at the observed frequencies.
- gradient_function
The gradient of the function which is applied to the frequency combinations.
- name
Human-readable name of the frequency combination function
- offset
Additive offset for the frequency combination function.
- value
Value of the frequency combination function at the observed frequencies.
- class AIMS.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 “Gaussian”, “Truncated_gaussian”, “Uniform”, “IMF1”, “IMF2”, “Uninformative”, “Above”, “Below”)
_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
- 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
- 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”, “Uninformative”, “Above”, or “Below”)
- values
List of parameters relevant to probability function
- class AIMS.Likelihood[source]
A class which described the likelihood function and allows users to evaluate it.
- Rkkinv
Inverse of R matrix from the WhoSGlAd method (see Eqs. (D.5) and (D.10) from Farnir et al. (2019))
Note
The first index corresponds to the p vector whereas the second index corresponds to the normalised q vector.
- add_combinations(string, num_list, den_list=[], target_ell=None)[source]
This finds the indices of modes which intervene in a frequency combination or ratio, as specified by the mandatory and optional arguments. These indices, the relevant coefficients, the numerator, the denominator, and the resultant value of the combination are stored in the
combination_functions
variable.- Parameters
string (string) – name of the frequency combination
num_list (list of (int,int,float)) – list of relative mode identifications and coefficients used to define a frequency combination or the numerator of a frequency ratio. This list contains tuples of the form (delta n, delta l, coeff).
den_list (list of (int,int,float)) – list of relative mode identifications and coefficients used to define the denominator of a frequency ratio. If absent, then, it is assumed that a linear combination of frequencies is represented. The form is the same as for
num_list
.target_ell (int) – this is used to impose a specific l value on the first selected mode.
- add_constraint(constraint)[source]
Add a supplementary constraint to the list of constraints.
- Parameters
constraint ((string,
Distribution
)) – supplementary constraint
- add_dnu_constraint(string, l_targets=[0])[source]
Add the large frequency separation as a contraint. The coefficients are obtained via a least-squares approach. The approach taken here has two advantages:
Correlations between the large frequency separation and other seismic constraints will be taken into account.
The same modes will be used in the same way, both for the observations and the models.
- Parameters
string (string) – name of this seismic constraint
l_targets (list of int) – specifies for which l values the large frequency separation is to be calculated. If
None
is supplied, all modes will be used.
Note
This uses an analytical approach and is therefore the prefered method.
- add_dnu_constraint_matrix(string, l_targets=[0])[source]
Add the large frequency separation as a contraint. The coefficients are obtained via a least-squares approach. The approach taken here has two advantages:
Correlations between the large frequency separation and other seismic constraints will be taken into account.
The same modes will be used in the same way, both for the observations and the models.
- Parameters
string (string) – name of this seismic constraint
l_targets (list of int) – specifies for which l values the large frequency separation is to be calculated. If
None
is supplied, all modes will be used.
Note
This uses a matrix approach and is therefore not the prefered method.
- add_nu_min_constraint(string, target_ell=0, pos=0, min_n=False)[source]
Add the minimun frequencies/modes of a specific ell value as a seismic constraint. Typically, such constraints are used as an “anchor” when combined with constraints based on frequency ratios.
- Parameters
string (string) – name of the seismic constraint
target_ell (int) – ell value of the minimum frequency/mode
pos (int) – position of desired frequency (0 = lowest, 1 = second lowest …)
min_n (boolean) – if
False
, look for minimum observational frequency. IfTrue
, look for minimum radial order.
- add_seismic_constraint(string)[source]
Add seismic contraints based on the keyword given in
string
.- Parameters
string (string) –
keyword which specifies the type of constraint to be added. Current options include:
nu
: individual frequenciesnu0
: individual frequencies (radial modes only)nu_min0
: radial mode with minimum frequencynu_min1
: radial mode with second lowest frequencynu_min2
: radial mode with third lowest frequencyr02
: \(r_{02}\) frequency ratiosr01
: \(r_{01}\) frequency ratiosr10
: \(r_{10}\) frequency ratiosdnu
: individual large frequency separations (using all modes)dnu0
: individual large frequency separations (using radial modes only)d2nu
: second differences (using all modes)avg_dnu
: average large frequency separation (using all modes)avg_dnu0
: average large frequency separation (using radial modes only)whosglad_dnu
: WhoSGlAd version of average large frequency separation (all modes)whosglad_dnu0
: WhoSGlAd version of average large frequency separation (l=0 modes)whosglad_dnu1
: WhoSGlAd version of average large frequency separation (l=1 modes)whosglad_r01
: WhoSGlAd version of average r01 frequency ratiowhosglad_r02
: WhoSGlAd version of average r02 frequency ratiowhosglad_Delta01
: WhoSGlAd Delta01 seismic constraintwhosglad_Delta02
: WhoSGlAd Delta02 seismic constraintwhosglad_eps0
: WhoSGlAd average frequency offset (l=0 modes)whosglad_eps1
: WhoSGlAd average frequency offset (l=1 modes)whosglad_AHe
: WhoSGlAd amplitude indicator for helium contentwhosglad_Abcz
: WhoSGlAd amplitude indicator for the base of the convection zone
- add_whosglad_AHe_constraint(string, maxpower=2)[source]
Add a seismic constraint which provides a normalised He glitch amplitude using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
maxpower (int) – maximum power in Pj(n)=n^j polynomial
- add_whosglad_Abcz_constraint(string, maxpower=2)[source]
Add a seismic constraint which provides a normalised glitch amplitude for the base of the convection zone using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
maxpower (int) – maximum power in Pj(n)=n^j polynomial
- add_whosglad_Delta_constraint(string, l_target, maxpower=2)[source]
Add a ratio of large frequency separations as a contraint using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
l_target (int) – specifies for which l to calculate the average Delta constraint.
maxpower (int) – maximum power in Pj(n)=n^j polynomial
- add_whosglad_dnu_constraint(string, l_targets=[0], maxpower=2)[source]
Add an average large frequency separation as a contraint using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
l_targets (list of int) – specifies for which l values the large frequency separation is to be calculated. If
None
is supplied, all modes will be used.maxpower (int) – maximum power in Pj(n)=n^j polynomial
- add_whosglad_epsilon_constraint(string, l_target, maxpower=2)[source]
Add an average frequency offset as a contraint using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
l_target (int) – specifies for which l to calculate the average epsilon offset.
maxpower (int) – maximum power in Pj(n)=n^j polynomial
- add_whosglad_ratio_constraint(string, l_target, maxpower=2)[source]
Add an average frequency ratio as a contraint using the WhoSGlAd method (Farnir et al. 2019). This approach has the advantage of leading to (near) orthogonal seismic indicators.
- Parameters
string (string) – name of this seismic constraint
l_target (int) – specifies for which l to calculate the average frequency ratio.
maxpower (int) – maximum power in Pj(n)=n^j polynomial
- apply_constraints(my_model)[source]
Calculate a \(\chi^2\) value for the set of constraints (excluding seismic constraints based on mode frequencies).
- Parameters
my_model (
model.Model
) – model for which the \(\chi^2\) value is being calculated- Returns
the \(\chi^2\) value deduced from classic constraints
- Return type
float
- assign_n(my_model)[source]
Assign the radial orders based on proximity to theoretical frequencies from an input model.
- Parameters
my_model (
model.Model
) – input model
- classic_weight
Absolute weight to be applied to classic constraints (incl. nu_max constraint).
- clear_seismic_constraints()[source]
This clears the seismic constraints. Specifically, the list of seismic combinations, and associated covariance matrix and its inverse are reinitialised.
- coeff
2D float array with the coefficients for each frequency combination. The indices are:
The index of the term
The cumulative index of the frequency combination
- combination_functions
This contains functions of frequency combinations.
- compare_frequency_combinations(my_model, mode_map, a=[])[source]
This finds a \(\chi^2\) value based on a comparison of frequency combination functions, as defined in the
combination_function
variable.- Parameters
my_model (
model.Model
) – model for which the \(\chi^2\) value is being calculatedmode_map (list of int) – a mapping which relates observed modes to theoretical ones
a (array-like) – parameters of surface correction terms
- Returns
the \(\chi^2\) value for the seismic constraints
- Return type
float
Note
I’m assuming none of the modes are missing (i.e. that mode_map doesn’t contain the value -1)
- constraints
List of constraints which intervene in the likelihood function.
- construct_whosglad_basis(maxpower=2, glitchpower=- 5)[source]
Construct orthogonal basis for the smooth component of a pulsation spectrum using the Gram-Schmidt orthonormalisation method as described in the WhoSGlAd method (Farnir et al. 2019).
- Parameters
maxpower (int) – maximum power in Pj(n)=n^j polynomial
glitchpower (int) – lower power used in Helium glitch component
- cov
Covariance matrix which intervenes when calculating frequency combinations.
- create_combination_arrays()[source]
Create array form of frequency combination functions to be used with a fortran based routine for calculating the seismic chi^2 value.
- create_mode_arrays()[source]
Create arrays with mode parameters (n, l, freq), which can be interfaced with fortran methods more easily.
- dot_product_whosglad(v1, v2)[source]
Dot product used for constructing an orthonormal set of seismic indicators using the WhoSGlAd method (Farnir et al. 2019).
- Parameters
v1 (float np.array) – a vector of frequencies (or some polynomial of the radial order)
v2 (float np.array) – a vector of frequencies (or some polynomial of the radial order)
- Returns
the dot product between v1 and v2
- Return type
float
- evaluate(my_model)[source]
Calculate ln of likelihood function (i.e. a \(\chi^2\) value) for a given model.
- Parameters
my_model (
model.Model
) – model for which the \(\chi^2\) value is being calculated- Returns
the \(\chi^2\) value, optionally the optimal surface amplitudes (depending on the value of
AIMS_configure.surface_option
), integers to indicate whether the model was rejected due to classic or seismic contraints- Return type
float, np.array (optional), int, int
Note
This avoids model interpolation and can be used to gain time.
- find_covariance()[source]
This prepares the covariance matrix and its inverse based on the frequency combinations in
combination_functions
.Warning
This method should be called after all of the methods which add to the list of frequency combinations.
- find_l_list(l_targets, npowers=2)[source]
Find a list of l values with the following properties:
each l value only occurs once
each l value given in the parameter
l_targets
is in the result l list, except if there is less than npowers modes with this l valueif the parameter
l_targets
isNone
, look for all l values with npowers or more modes associated with them
- Parameters
l_targets (list of int) – input list of l values
npowers (int) – the number of powers in the fit
- Returns
new list of l values with the above properties
- Return type
list of int
- find_map(my_model, use_n)[source]
This finds a map which indicates the correspondance between observed modes and theoretical modes from
my_model
.- Parameters
my_model – model for which the \(\chi^2\) value is being calculated
use_n (boolean) – specify whether to use the radial order when finding the map from observed modes to theoretical modes. If
False
, the map is based on frequency proximity.
- Returns
the correspondance between observed and theoretical modes from the above model, and the number of observed modes which weren’t mapped onto theoretical modes
- Return type
list of int, int
Note
a value of -1 is used to indicate that no theoretical mode corresponds to a particular observed mode.
only zero or one observed mode is allowed to correspond to a theoretical mode
- find_vec(a_combination_function)[source]
This finds a set of coefficients which intervene when constructing the coviance matrix for frequency combination functions.
- Parameters
a_combination_function (
Combination_function
) – variable which specifies the frequency combination function.- Returns
the above set of coefficients
- Return type
np.array
- find_weights()[source]
Find absolute weights for seismic and classic constraints based on options in
AIMS_configure.py
.
- find_weights_new()[source]
Find absolute weights for seismic and classic constraints based on options in
AIMS_configure.py
.
- fvalues
Array with the observed frequencies
- get_optimal_surface_amplitudes(my_model, mode_map)[source]
Find optimal surface correction amplitude, for the surface correction specified by
surface_option
.- Parameters
my_model (
model.Model
) – the model for which we’re finding the surface correction amplitudemode_map (list of int) – a mapping which relates observed modes to theoretical ones
- Returns
optimal surface correction amplitudes
- Return type
np.array
- guess_dnu(with_n=False)[source]
Guess the large frequency separation based on the radial modes.
- Parameters
with_n (boolean) – specifies whether to use the n values already stored with each mode, when calculating the large frequency separation.
- Returns
the large frequency separation
- Return type
float
- guess_n()[source]
Guess the radial order of the observed pulsations modes.
This method uses the large frequency separation, as calculated with
guess_dnu()
, to estimate the radial orders. These orders are subsequently adjusted to avoid multiple modes with the same identification. The resultant radial orders could be off by a constant offset, but this is not too problematic when computing frequency combinations or ratios.
- indices
2D int array with the mode indices for each frequency combination. The indices are:
The index of the term
The cumulative index of the frequency combination
- invcov
Inverse of covariance matrix,
Likelihood.cov
.
- lvalues
Array with the l values of the observed modes
- ncoeff
The number of terms in each frequency combination. The index is a cumulative index that uniquely designates one of the frequency combinations from one of the frequency combination functions. It is also the same as the second index in the
coeff
andindices
variables.
- ncomb
1D int array which specifies the number of frequency combinations within each frequency combination function
- nvalues
Array with the n values of the observed modes
- qk
Normalised q vectors from the WhoSGlAd method (Farnir et al. 2019).
- read_constraints(filename, factor=1.0)[source]
Read a file with pulsation data and constraints.
- Parameters
filename (string) – name of file with pulsation data.
factor (float) – multiplicative factor for pulsation frequencies. Can be used for conversions.
- seismic_weight
Absolute weight to be applied to seismic constraints
- class AIMS.Mode(_n, _l, _freq, _dfreq)[source]
A class which describes an observed pulsation mode.
- Parameters
_n (int) – radial order of observed mode
_l (int) – harmonic degree of observed mode.
_freq (float) – pulsation frequency (in \(\mathrm{\mu Hz}\)).
_dfreq (float) – error bar on pulsation frequency (in \(\mathrm{\mu Hz}\)).
Warning
Negative values are not accepted for
_l
,_freq
, or_dfreq
.- dfreq
Error bar on pulsation frequency (in \(\mathrm{\mu Hz}\)).
- freq
Pulsation frequency (in \(\mathrm{\mu Hz}\)).
- l
Harmonic degree of observed mode.
- match(a_mode)[source]
Check to see if input mode has the same (n,l) values as the current mode.
- Parameters
a_mode (
Mode
) – input mode which is being compared with current mode.- Returns
True
if the input mode has the same (n,l) values as the current mode.- Return type
boolean
- n
Radial order of observed mode.
- class AIMS.Prior_list[source]
A class which contains a list of priors as well as convenient methods for adding priors and for evaluating them.
- add_prior(aPrior)[source]
Add a prior to the list.
- Parameters
aPrior (
Distribution
) – prior which is to be added to the list.
- priors
A list of probability distributions which correspond to priors.
- class AIMS.Probability(_priors, _likelihood)[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)[source]
Evalulate the ln of the product of the priors and likelihood function, i.e. the probability, for a given model, to within an additive constant.
- Parameters
my_model (
model.Model
) – input model- Returns
the ln of the probability, integers indicating if the model has bee rejected based on classic constraints, seismic constraints, and/or priors
- Return type
float, int, int, int
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 evaluate 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
- likelihood
The likelihood function.
- priors
The set of priors.
- AIMS.accepted_parameters = []
list of parameters associated with accepted models
- AIMS.append_osm_parameter(config_osm, name, value, step, rate, bounds)[source]
Add a parameter in xlm format in the file with the classic constraints for OSM.
- Parameters
config_osm (lxml.etree._Element) – XLM etree element to which to add the parameter
name (string) – name of the parameter
value (float) – value of the parameter
step (float) – parameter step (this intervenes when numerically calculating derivatives with respect to this parameter)
rate (float) – parameter rate (this corresponds to a tolerance on this parameter)
bounds (float tuple) – bounds on the parameter
- AIMS.append_osm_surface_effects(modes_osm, name, numax, values)[source]
Add a method with which to calculate surface effects to the OSM contraint file.
- Parameters
modes_osm (lxml.etree._Element) – XML element to which to add the surface effects method
name (string) – name of the method
numax (float) – value of numax
values (float tuple) – values which intervene in the method
- AIMS.autocorr_time = []
integrated autocorrelelation time (this is useful for testing convergence)
- AIMS.best_MCMC_model = None
best model from the MCMC run
- AIMS.best_MCMC_params = None
parameters for the model
best_MCMC_model
- AIMS.best_MCMC_result = -1e+300
ln(probability) result for the model
best_MCMC_model
- AIMS.best_age_range = 0.0
Age range on track with
best_model_model
- AIMS.best_grid_model = None
best model from a scan of the entire grid
- AIMS.best_grid_params = None
parameters for the model
best_grid_model
- AIMS.best_grid_result = -1e+300
ln(probability) result for the model
best_grid_model
- AIMS.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 AIMS and explain what variable has an erroneous value.
- AIMS.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
- AIMS.find_best_model()[source]
Scan through grid of models to find “best” model for a given probability function (i.e. the product of priors and a likelihood function).
- AIMS.find_best_model_in_track(ntrack)[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).- Parameters
ntrack (int) – number of the evolutionary track
- Returns
the ln(probability) value, the “best” model, the parameters of accepted models, the parameters of rejected models, the age range of the track, and the numbers of models rejected due to classic constraints, seismic contraints, and priors
- Return type
(float,
model.Model
, 2D float list, 2D float list, float, int, int, int)
- AIMS.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
- AIMS.grid = None
grid of models
- AIMS.grid_params_MCMC = ()
parameters used in the MCMC run (excluding surface correction parameters)
- AIMS.grid_params_MCMC_with_surf = ()
parameters used in the MCMC run (including surface correction parameters)
- AIMS.init_walkers()[source]
Initialise the walkers used in emcee.
- Returns
array of starting parameters
- Return type
np.array
- AIMS.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
.
- AIMS.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
- AIMS.log0 = -1e+300
a large negative value used to represent ln(0)
- AIMS.mc2 = None
True if emcee’s version is 2 or prior
- AIMS.my_map = None
pointer to the map function (either the parallel or sequential versions)
- AIMS.nreject_classic = 0
Number of models rejected because of classic constraints
- AIMS.nreject_prior = 0
Number of models rejected based on priors
- AIMS.nreject_seismic = 0
Number of models rejected because of seismic constraints
- AIMS.nsurf = 0
number of surface term parameters
- AIMS.output_folder = None
folder in which to write the results
- AIMS.plot_distrib_iter(percentiles, labels, 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 dimensions in parameters space
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.
- AIMS.plot_distrib_iter_old(samples, labels, folder)[source]
Plot individual distribution of walkers as a function of iterations.
- Parameters
samples (np.array) – samples from the emcee run
labels (list of strings) – labels for the different dimensions in parameters space
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.
- AIMS.plot_echelle_diagram(my_model, my_params, model_name)[source]
Produce an echelle diagram for input model.
- Parameters
my_model (
model.Model
) – model for which we’re making an echelle diagrammy_params (array-like) – parameters of the model
model_name (string) – name used to describe this model. This is also used when naming the file with the echelle diagram.
- AIMS.plot_frequencies(grid)[source]
Plot frequencies along an evolutionary track. For debugging purposes only.
- Parameters
grid (
Model_grid
) – grid of models
- AIMS.plot_frequency_diff(my_model, my_params, model_name, scaled=False)[source]
Make a diagram with frequency differences.
- Parameters
my_model (
model.Model
) – model for which we’re plotting frequencie differencesmy_params (array-like) – parameters of the model
model_name (string) – name used to describe this model. This is also used when naming the file with the plot.
scaled (boolean) – if True, scale frequency differences by frequency error bars
- AIMS.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
- AIMS.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.
- AIMS.pool = None
pool from which to carry out parallel computations
- AIMS.prob = None
Probability
type object that represents the probability function which includes the likelihood and priors
- AIMS.rejected_parameters = []
list of parameters associated with rejected models
- AIMS.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
- AIMS.statistical_model = None
model corresponding to statistical parameters
- AIMS.statistical_params = None
parameters for the model
statistical_model
- AIMS.statistical_result = -1e+300
ln(probability) result for the model
statistical_model
- AIMS.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
- AIMS.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
- AIMS.tight_ball_distributions = None
Prior_list
type object with the distributions for the initial tight ball
- AIMS.write_LEGACY_summary(filename, KIC, labels, samples)[source]
Write a one line summary of the statistical properties based on a sequence of realisations to a file. The format matches that of the LEGACY project.
The results include:
average values for each variable (statistical mean)
error bars for each variable (standard mean deviation)
- Parameters
filename (string) – name of file in which to write the statistical properties
KIC (string) – KIC number of the star
labels (list of strings) – names of relevant variables
samples (np.array) – samples for which statistical properties are calculated
- AIMS.write_SPInS_file_cgs(filename)[source]
Write list file compatible with SPInS. The quantities mass, luminosity, and radius are provided in cgs.
Note
The header in the output file might need to be edited by hand.
- AIMS.write_SPInS_file_solar(filename)[source]
Write list file compatible with SPInS. The quantities mass, luminosity, and radius are provided solar units.
Note
The header in the output file might need to be edited by hand.
- AIMS.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
- AIMS.write_combinations(filename, samples)[source]
Produce a list of linear combinations of grid models (based on interpolation) corresponding to the provided model parameters.
- Parameters
filename (string) – name of the file to which to write the model combinations
samples (np.array) – set of model parameters for which we would like to obtain the grid models and interpolation coefficients
- AIMS.write_list_file(filename)[source]
Write list file from which to generate binary grid. Various filters can be included to reduce the number of models.
Note
This code is intended for developpers, not first time users.
- AIMS.write_model(my_model, my_params, my_result, model_name, extended=False)[source]
Write text file with caracteristics of input model.
- Parameters
my_model (
model.Model
) – model for which we’re writing a text filemy_params (array-like) – parameters of the model
my_result (float) – ln(P) value obtained for the model
model_name (string) – name used to describe this model. This is also used when naming the text file.
extended – if set to True, all of the theoretical modes are saved in the text file, including those not matched to observations
- AIMS.write_osm_don(filename, my_model)[source]
Write file with choice of physical ingredients to be used by CESAM or CESTAM and OSM.
- Parameters
filename (string) – name of file which will contain the physical ingredients
my_model (
model.Model
) – model from which which is derived various physical constraints/settings
Note
Written by B. Herbert.
- AIMS.write_osm_frequencies(filename, my_model)[source]
Write file with frequencies for Optimal Stellar Model (OSM), written by R. Samadi.
- Parameters
filename (string) – name of file which will contain the frequencies
my_model (
model.Model
) – model from which are derived the radial orders
Note
Written by B. Herbert.
- AIMS.write_osm_xml(filename, my_params, my_model)[source]
Write file with classic constraints for OSM
- Parameters
filename (string) – name of file with classic constraints
my_model (
model.Model
) – model used in deriving some of the constraints
Note
Originally written by B. Herbert. Includes some modifications.
- AIMS.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
- AIMS.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
- AIMS.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
- AIMS.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