getdist.mcsamples¶
-
getdist.mcsamples.
loadMCSamples
(file_root: str, ini: Union[None, str, getdist.inifile.IniFile] = None, jobItem=None, no_cache=False, settings: Optional[Mapping[str, Any], None] = None) → getdist.mcsamples.MCSamples[source]¶ Loads a set of samples from a file or files.
Sample files are plain text (file_root.txt) or a set of files (file_root_1.txt, file_root_2.txt, etc.).
Auxiliary files file_root.paramnames gives the parameter names and (optionally) file_root.ranges gives hard prior parameter ranges.
For a description of the various analysis settings and default values see analysis_defaults.ini.
Parameters: - file_root – The root name of the files to read (no extension)
- ini – The name of a .ini file with analysis settings to use
- jobItem – an optional grid jobItem instance for a CosmoMC grid output
- no_cache – Indicates whether or not we should cache loaded samples in a pickle
- settings – dictionary of analysis settings to override defaults
Returns: The
MCSamples
instance
-
class
getdist.mcsamples.
MCSamples
(root: Optional[str, None] = None, jobItem=None, ini=None, settings: Optional[Mapping[str, Any], None] = None, ranges=None, samples: Union[numpy.ndarray, Iterable[numpy.ndarray], None] = None, weights: Union[numpy.ndarray, Iterable[numpy.ndarray], None] = None, loglikes: Union[numpy.ndarray, Iterable[numpy.ndarray], None] = None, temperature: Optional[float, None] = None, **kwargs)[source]¶ The main high-level class for a collection of parameter samples.
Derives from
chains.Chains
, adding high-level functions including Kernel Density estimates, parameter ranges and custom settings.For a description of the various analysis settings and default values see analysis_defaults.ini.
Parameters: - root – A root file name when loading from file
- jobItem – optional jobItem for parameter grid item. Should have jobItem.chainRoot and jobItem.batchPath
- ini – a .ini file to use for custom analysis settings
- settings – a dictionary of custom analysis settings
- ranges – a dictionary giving any additional hard prior bounds for parameters, e.g. {‘x’:[0, 1], ‘y’:[None,2]}
- samples – if not loading from file, array of parameter values for each sample, passed
to
setSamples()
, or list of arrays if more than one chain - weights – array of weights for samples, or list of arrays if more than one chain
- loglikes – array of -log(Likelihood) for samples, or list of arrays if more than one chain
- temperatute – temperature of the sample. If not specified will be read from the root.properties.ini file if it exists and otherwise default to 1.
- kwargs –
- keyword arguments passed to inherited classes, e.g. to manually make a samples object from
- sample arrays in memory:
- paramNamesFile: optional name of .paramnames file with parameter names
- names: list of names for the parameters, or list of arrays if more than one chain
- labels: list of latex labels for the parameters
- renames: dictionary of parameter aliases
- ignore_rows:
- if int >=1: The number of rows to skip at the file in the beginning of the file
- if float <1: The fraction of rows to skip at the beginning of the file
- label: a latex label for the samples
- name_tag: a name tag for this instance
- sampler: string describing the type of samples; if “nested” or “uncorrelated” the effective number of samples is calculated using uncorrelated approximation. If not specified will be read from the root.properties.ini file if it exists and otherwise default to “mcmc”.
-
PCA
(params, param_map=None, normparam=None, writeDataToFile=False, filename=None, conditional_params=(), n_best_only=None)[source]¶ Perform principal component analysis (PCA). In other words, get eigenvectors and eigenvalues for normalized variables with optional (log modulus) mapping to find power law fits.
Parameters: - params – List of names of the parameters to use
- param_map – A transformation to apply to parameter values; A list or string containing either N (no transformation) or L (for log transform) for each parameter. By default, uses log if no parameter values cross zero
- normparam – optional name of parameter to normalize result (i.e. this parameter will have unit power)
- writeDataToFile – True to write the output to file.
- filename – The filename to write, by default root_name.PCA.
- conditional_params – optional list of parameters to treat as fixed, i.e. for PCA conditional on fixed values of these parameters
- n_best_only – return just the short summary constraint for the tightest n_best_only constraints
Returns: a string description of the output of the PCA
-
addDerived
(paramVec, name, label='', comment='', range=None)[source]¶ Adds a new derived parameter
Parameters: - paramVec – The vector of parameter values to add. For example a combination of parameter arrays from MCSamples.getParams()
- name – The name for the new parameter
- label – optional latex label for the parameter
- comment – optional comment describing the parameter
- range – if specified, a tuple of min, max values for the new parameter hard prior bounds (either can be None for one-side bound)
Returns: The added parameter’s
ParamInfo
object
-
changeSamples
(samples)¶ Sets the samples without changing weights and loglikes.
Parameters: samples – The samples to set
-
confidence
(paramVec, limfrac, upper=False, start=0, end=None, weights=None)¶ Calculate sample confidence limits, not using kernel densities just counting samples in the tails
Parameters: - paramVec – array of parameter values or int index of parameter to use
- limfrac – fraction of samples in the tail, e.g. 0.05 for a 95% one-tail limit, or 0.025 for a 95% two-tail limit
- upper – True to get upper limit, False for lower limit
- start – Start index for the vector to use
- end – The end index, use None to go all the way to the end of the vector.
- weights – numpy array of weights for each sample, by default self.weights
Returns: confidence limit (parameter value when limfac of samples are further in the tail)
-
cool
(cool=None)[source]¶ Cools the samples, i.e. multiplies log likelihoods by cool factor and re-weights accordingly :param cool: cool factor, optional if the sample has a temperature specified.
-
copy
(label=None, settings=None)[source]¶ Create a copy of this sample object
Parameters: - label – optional lable for the new copy
- settings – optional modified settings for the new copy
Returns: copyied
MCSamples
instance
-
corr
(pars=None)¶ Get the correlation matrix
Parameters: pars – If specified, list of parameter vectors or int indices to use Returns: The correlation matrix.
-
cov
(pars=None, where=None)¶ Get parameter covariance
Parameters: - pars – if specified, a list of parameter vectors or int indices to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: The covariance matrix
-
deleteFixedParams
()¶ Delete parameters that are fixed (the same value in all samples)
-
deleteZeros
()¶ Removes samples with zero weight
-
filter
(where)¶ Filter the stored samples to keep only samples matching filter
Parameters: where – list of sample indices to keep, or boolean array filter (e.g. x>5 to keep only samples where x>5)
-
get1DDensity
(name, **kwargs)[source]¶ Returns a
Density1D
instance for parameter with given name. Result is cached.Parameters: - name – name of the parameter
- kwargs – arguments for
get1DDensityGridData()
Returns: A
Density1D
instance for parameter with given name
-
get1DDensityGridData
(j, paramConfid=None, meanlikes=False, **kwargs)[source]¶ Low-level function to get a
Density1D
instance for the marginalized 1D density of a parameter. Result is not cached.Parameters: - j – a name or index of the parameter
- paramConfid – optional cached
ParamConfidenceData
instance - meanlikes – include mean likelihoods
- kwargs –
optional settings to override instance settings of the same name (see analysis_settings):
- smooth_scale_1D
- boundary_correction_order
- mult_bias_correction_order
- fine_bins
- num_bins
Returns: A
Density1D
instance
-
get2DDensity
(x, y, normalized=False, **kwargs)[source]¶ Returns a
Density2D
instance with marginalized 2D density.Parameters: - x – index or name of x parameter
- y – index or name of y parameter
- normalized – if False, is normalized so the maximum is 1, if True, density is normalized
- kwargs – keyword arguments for the
get2DDensityGridData()
function
Returns: Density2D
instance
-
get2DDensityGridData
(j, j2, num_plot_contours=None, get_density=False, meanlikes=False, **kwargs)[source]¶ Low-level function to get 2D plot marginalized density and optional additional plot data.
Parameters: - j – name or index of the x parameter
- j2 – name or index of the y parameter.
- num_plot_contours – number of contours to calculate and return in density.contours
- get_density – only get the 2D marginalized density, don’t calculate confidence level members
- meanlikes – calculate mean likelihoods as well as marginalized density (returned as array in density.likes)
- kwargs –
optional settings to override instance settings of the same name (see analysis_settings):
- fine_bins_2D
- boundary_correction_order
- mult_bias_correction_order
- smooth_scale_2D
Returns: a
Density2D
instance
-
getAutoBandwidth1D
(bins, par, param, mult_bias_correction_order=None, kernel_order=1, N_eff=None)[source]¶ Get optimized kernel density bandwidth (in units of the range of the bins) Based on optimal Improved Sheather-Jones bandwidth for basic Parzen kernel, then scaled if higher-order method being used. For details see the notes at arXiv:1910.13970.
Parameters: - bins – numpy array of binned weights for the samples
- par – A
ParamInfo
instance for the parameter to analyse - param – index of the parameter to use
- mult_bias_correction_order – order of multiplicative bias correction (0 is basic Parzen kernel); by default taken from instance settings.
- kernel_order – order of the kernel (0 is Parzen, 1 does linear boundary correction, 2 is a higher-order kernel)
- N_eff – effective number of samples. If not specified estimated using weights, autocorrelations, and fiducial bandwidth
Returns: kernel density bandwidth (in units the range of the bins)
-
getAutoBandwidth2D
(bins, parx, pary, paramx, paramy, corr, rangex, rangey, base_fine_bins_2D, mult_bias_correction_order=None, min_corr=0.2, N_eff=None, use_2D_Neff=False)[source]¶ Get optimized kernel density bandwidth matrix in parameter units, using Improved Sheather Jones method in sheared parameters. The shearing is determined using the covariance, so you know the distribution is multi-modal, potentially giving ‘fake’ correlation, turn off shearing by setting min_corr=1. For details see the notes arXiv:1910.13970.
Parameters: - bins – 2D numpy array of binned weights
- parx – A
ParamInfo
instance for the x parameter - pary – A
ParamInfo
instance for the y parameter - paramx – index of the x parameter
- paramy – index of the y parameter
- corr – correlation of the samples
- rangex – scale in the x parameter
- rangey – scale in the y parameter
- base_fine_bins_2D – number of bins to use for re-binning in rotated parameter space
- mult_bias_correction_order – multiplicative bias correction order (0 is Parzen kernel); by default taken from instance settings
- min_corr – minimum correlation value at which to bother de-correlating the parameters
- N_eff – effective number of samples. If not specified, uses rough estimate that accounts for weights and strongly-correlated nearby samples (see notes)
- use_2D_Neff – if N_eff not specified, whether to use 2D estimate of effective number, or approximate from the 1D results (default from use_effective_samples_2D setting)
Returns: kernel density bandwidth matrix in parameter units
-
getAutocorrelation
(paramVec, maxOff=None, weight_units=True, normalized=True)¶ Gets auto-correlation of an array of parameter values (e.g. for correlated samples from MCMC)
By default, uses weight units (i.e. standard units for separate samples from original chain). If samples are made from multiple chains, neglects edge effects.
Parameters: - paramVec – an array of parameter values, or the int index of the parameter in stored samples to use
- maxOff – maximum autocorrelation distance to return
- weight_units – False to get result in sample point (row) units; weight_units=False gives standard definition for raw chains
- normalized – Set to False to get covariance (note even if normalized, corr[0]<>1 in general unless weights are unity).
Returns: zero-based array giving auto-correlations
-
getBestFit
(max_posterior=True)[source]¶ - Returns a
BestFit
object with best-fit point stored in .minimum or .bestfit fileParameters: max_posterior – whether to get maximum posterior (from .minimum file) or maximum likelihood (from .bestfit file) Returns:
-
getBounds
()[source]¶ Returns the bounds in the form of a
ParamBounds
instance, for example for determining plot rangesBounds are not the same as self.ranges, as if samples are not near the range boundary, the bound is set to None
Returns: a ParamBounds
instance
-
getCombinedSamplesWithSamples
(samps2, sample_weights=(1, 1))[source]¶ Make a new
MCSamples
instance by appending samples from samps2 for parameters which are in common. By defaultm they are weighted so that the probability mass of each set of samples is the same, independent of tha actual sample sizes. The sample_weights parameter can be adjusted to change the relative weighting.Parameters: - samps2 –
MCSamples
instance to merge - sample_weights – relative weights for combining the samples. Set to None to just directly append samples.
Returns: a new
MCSamples
instance with the combined samples- samps2 –
-
getConvergeTests
(test_confidence=0.95, writeDataToFile=False, what=('MeanVar', 'GelmanRubin', 'SplitTest', 'RafteryLewis', 'CorrLengths'), filename=None, feedback=False)[source]¶ Do convergence tests.
Parameters: - test_confidence – confidence limit to test for convergence (two-tail, only applies to some tests)
- writeDataToFile – True to write output to a file
- what –
The tests to run. Should be a list of any of the following:
- ’MeanVar’: Gelman-Rubin sqrt(var(chain mean)/mean(chain var)) test in individual parameters (multiple chains only)
- ’GelmanRubin’: Gelman-Rubin test for the worst orthogonalized parameter (multiple chains only)
- ’SplitTest’: Crude test for variation in confidence limits when samples are split up into subsets
- ’RafteryLewis’: Raftery-Lewis test (integer weight samples only)
- ’CorrLengths’: Sample correlation lengths
- filename – The filename to write to, default is file_root.converge
- feedback – If set to True, Prints the output as well as returning it.
Returns: text giving the output of the tests
Gets a list of most correlated variable pair names.
Parameters: - num_plots – The number of plots
- nparam – maximum number of pairs to get
Returns: list of [x,y] pair names
-
getCorrelationLength
(j, weight_units=True, min_corr=0.05, corr=None)¶ Gets the auto-correlation length for parameter j
Parameters: - j – The index of the parameter to use
- weight_units – False to get result in sample point (row) units; weight_units=False gives standard definition for raw chains
- min_corr – specifies a minimum value of the autocorrelation to use, e.g. where sampling noise is typically as large as the calculation
- corr – The auto-correlation array to use, calculated internally by default
using
getAutocorrelation()
Returns: the auto-correlation length
-
getCorrelationMatrix
()¶ Get the correlation matrix of all parameters
Returns: The correlation matrix
-
getCov
(nparam=None, pars=None)¶ Get covariance matrix of the parameters. By default, uses all parameters, or can limit to max number or list.
Parameters: - nparam – if specified, only use the first nparam parameters
- pars – if specified, a list of parameter indices (0,1,2..) to include
Returns: covariance matrix.
-
getCovMat
()[source]¶ Gets the CovMat instance containing covariance matrix for all the non-derived parameters (for example useful for subsequent MCMC runs to orthogonalize the parameters)
Returns: A CovMat
object holding the covariance
-
getEffectiveSamples
(j=0, min_corr=0.05)¶ Gets effective number of samples N_eff so that the error on mean of parameter j is sigma_j/N_eff
Parameters: - j – The index of the param to use.
- min_corr – the minimum value of the auto-correlation to use when estimating the correlation length
-
getEffectiveSamplesGaussianKDE
(paramVec, h=0.2, scale=None, maxoff=None, min_corr=0.05)¶ Roughly estimate an effective sample number for use in the leading term for the MISE (mean integrated squared error) of a Gaussian-kernel KDE (Kernel Density Estimate). This is used for optimizing the kernel bandwidth, and though approximate should be better than entirely ignoring sample correlations, or only counting distinct samples.
Uses fiducial assumed kernel scale h; result does depend on this (typically by factors O(2))
For bias-corrected KDE only need very rough estimate to use in rule of thumb for bandwidth.
In the limit h-> 0 (but still >0) answer should be correct (then just includes MCMC rejection duplicates). In reality correct result for practical h should depend on shape of the correlation function.
If self.sampler is ‘nested’ or ‘uncorrelated’ return result for uncorrelated samples.
Parameters: - paramVec – parameter array, or int index of parameter to use
- h – fiducial assumed kernel scale.
- scale – a scale parameter to determine fiducial kernel width, by default the parameter standard deviation
- maxoff – maximum value of auto-correlation length to use
- min_corr – ignore correlations smaller than this auto-correlation
Returns: A very rough effective sample number for leading term for the MISE of a Gaussian KDE.
-
getEffectiveSamplesGaussianKDE_2d
(i, j, h=0.3, maxoff=None, min_corr=0.05)¶ Roughly estimate an effective sample number for use in the leading term for the 2D MISE. If self.sampler is ‘nested’ or ‘uncorrelated’ return result for uncorrelated samples.
Parameters: - i – parameter array, or int index of first parameter to use
- j – parameter array, or int index of second parameter to use
- h – fiducial assumed kernel scale.
- maxoff – maximum value of auto-correlation length to use
- min_corr – ignore correlations smaller than this auto-correlation
Returns: A very rough effective sample number for leading term for the MISE of a Gaussian KDE.
-
getFractionIndices
(weights, n)[source]¶ Calculates the indices of weights that split the weights into sets of equal 1/n fraction of the total weight
Parameters: - weights – array of weights
- n – number of groups to split into
Returns: array of indices of the boundary rows in the weights array
-
getGelmanRubin
(nparam=None, chainlist=None)¶ Assess the convergence using the maximum var(mean)/mean(var) of orthogonalized parameters c.f. Brooks and Gelman 1997.
Parameters: - nparam – The number of parameters, by default uses all
- chainlist – list of
WeightedSamples
, the samples to use. Defaults to all the separate chains in this instance.
Returns: The worst var(mean)/mean(var) for orthogonalized parameters. Should be <<1 for good convergence.
-
getGelmanRubinEigenvalues
(nparam=None, chainlist=None)¶ Assess convergence using var(mean)/mean(var) in the orthogonalized parameters c.f. Brooks and Gelman 1997.
Parameters: - nparam – The number of parameters (starting at first), by default uses all of them
- chainlist – list of
WeightedSamples
, the samples to use. Defaults to all the separate chains in this instance.
Returns: array of var(mean)/mean(var) for orthogonalized parameters
-
getInlineLatex
(param, limit=1, err_sig_figs=None)[source]¶ Get snippet like: A=x\pm y. Will adjust appropriately for one and two tail limits.
Parameters: - param – The name of the parameter
- limit – which limit to get, 1 is the first (default 68%), 2 is the second (limits array specified by self.contours)
- err_sig_figs – significant figures in the error
Returns: The tex snippet.
-
getLabel
()¶ Return the latex label for the samples
Returns: the label
-
getLatex
(params=None, limit=1, err_sig_figs=None)[source]¶ Get tex snippet for constraints on a list of parameters
Parameters: - params – list of parameter names, or a single parameter name
- limit – which limit to get, 1 is the first (default 68%), 2 is the second (limits array specified by self.contours)
- err_sig_figs – significant figures in the error
Returns: labels, texs: a list of parameter labels, and a list of tex snippets, or for a single parameter, the latex snippet.
-
getLikeStats
()[source]¶ Get best fit sample and n-D confidence limits, and various likelihood based statistics
Returns: a LikeStats
instance storing N-D limits for parameter i in result.names[i].ND_limit_top, result.names[i].ND_limit_bot, and best-fit sample value in result.names[i].bestfit_sample
-
getLower
(name)[source]¶ Return the lower limit of the parameter with the given name.
Parameters: name – parameter name Returns: The lower limit if name exists, None otherwise.
-
getMargeStats
(include_bestfit=False)[source]¶ Returns a
MargeStats
object with marginalized 1D parameter constraintsParameters: include_bestfit – if True, set best fit values by loading from root_name.minimum file (assuming it exists) Returns: A MargeStats
instance
-
getMeans
(pars=None)¶ Gets the parameter means, from saved array if previously calculated.
Parameters: pars – optional list of parameter indices to return means for Returns: numpy array of parameter means
-
getName
()¶ Returns the name tag of these samples.
Returns: The name tag
-
getNumSampleSummaryText
()[source]¶ Returns a summary text describing numbers of parameters and samples, and various measures of the effective numbers of samples.
Returns: The summary text as a string.
-
getParamBestFitDict
(best_sample=False, want_derived=True, want_fixed=True, max_posterior=True)[source]¶ Gets an ordered dictionary of parameter values for the best fit point, assuming calculated results from mimimization runs in .minimum (max posterior) .bestfit (max likelihood) files exists.
Can also get the best-fit (max posterior) sample, which typically has a likelihood that differs significantly from the true best fit in high dimensions.
Parameters: - best_sample – load from global minimum files (False, default) or using maximum posterior sample (True)
- want_derived – include derived parameters
- want_fixed – also include values of any fixed parameters
- max_posterior – whether to get maximum posterior (from .minimum file) or maximum likelihood (from .bestfit file)
Returns: ordered dictionary of parameter values
-
getParamNames
()¶ Get
ParamNames
object with names for the parametersReturns: ParamNames
object giving parameter names and labels
-
getParamSampleDict
(ix, want_derived=True, want_fixed=True)[source]¶ Gets a dictionary of parameter values for sample number ix
Parameters: - ix – index of the sample to return (zero based)
- want_derived – include derived parameters
- want_fixed – also include values of any fixed parameters
Returns: ordered dictionary of parameter values
-
getParams
()¶ Creates a
ParSamples
object, with variables giving vectors for all the parameters, for example samples.getParams().name1 would be the vector of samples with name ‘name1’Returns: A ParSamples
object containing all the parameter vectors, with attributes given by the parameter names
-
getRawNDDensity
(xs, normalized=False, **kwargs)[source]¶ Returns a
DensityND
instance with marginalized ND density.Parameters: - xs – indices or names of x_i parameters
- kwargs – keyword arguments for the
getNDDensityGridData()
function - normalized – if False, is normalized so the maximum is 1, if True, density is normalized
Returns: DensityND
instance
-
getRawNDDensityGridData
(js, writeDataToFile=False, num_plot_contours=None, get_density=False, meanlikes=False, maxlikes=False, **kwargs)[source]¶ Low-level function to get unsmooth ND plot marginalized density and optional additional plot data.
Parameters: - js – vector of names or indices of the x_i parameters
- writeDataToFile – save outputs to file
- num_plot_contours – number of contours to calculate and return in density.contours
- get_density – only get the ND marginalized density, no additional plot data, no contours.
- meanlikes – calculate mean likelihoods as well as marginalized density (returned as array in density.likes)
- maxlikes – calculate the profile likelihoods in addition to the others (returned as array in density.maxlikes)
- kwargs – optional settings to override instance settings of the same name (see analysis_settings):
Returns: a
DensityND
instance
-
getRenames
()¶ Gets dictionary of renames known to each parameter.
-
getSeparateChains
() → List[getdist.chains.WeightedSamples]¶ Gets a list of samples for separate chains. If the chains have already been combined, uses the stored sample offsets to reconstruct the array (generally no array copying)
Returns: The list of WeightedSamples
for each chain.
-
getSignalToNoise
(params, noise=None, R=None, eigs_only=False)¶ Returns w, M, where w is the eigenvalues of the signal to noise (small y better constrained)
Parameters: - params – list of parameters indices to use
- noise – noise matrix
- R – rotation matrix, defaults to inverse of Cholesky root of the noise matrix
- eigs_only – only return eigenvalues
Returns: w, M, where w is the eigenvalues of the signal to noise (small y better constrained)
-
getTable
(columns=1, include_bestfit=False, **kwargs)[source]¶ Creates and returns a
ResultTable
instance. See alsogetInlineLatex()
.Parameters: - columns – number of columns in the table
- include_bestfit – True to include the bestfit parameter values (assuming set)
- kwargs – arguments for
ResultTable
constructor.
Returns: A
ResultTable
instance
-
getUpper
(name)[source]¶ Return the upper limit of the parameter with the given name.
Parameters: name – parameter name Returns: The upper limit if name exists, None otherwise.
-
getVars
()¶ Get the parameter variances
Returns: A numpy array of variances.
-
get_norm
(where=None)¶ gets the normalization, the sum of the sample weights: sum_i w_i
Parameters: where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5). Returns: normalization
-
initParamConfidenceData
(paramVec, start=0, end=None, weights=None)¶ Initialize cache of data for calculating confidence intervals
Parameters: - paramVec – array of parameter values or int index of parameter to use
- start – The sample start index to use
- end – The sample end index to use, use None to go all the way to the end of the vector
- weights – A numpy array of weights for each sample, defaults to self.weights
Returns: ParamConfidenceData
instance
-
initParameters
(ini)[source]¶ Initializes settings. Gets parameters from
IniFile
.Parameters: ini – The IniFile
to be used
-
loadChains
(root, files_or_samples: Sequence, weights=None, loglikes=None, ignore_lines=None)¶ Loads chains from files.
Parameters: - root – Root name
- files_or_samples – list of file names or list of arrays of samples, or single array of samples
- weights – if loading from arrays of samples, corresponding list of arrays of weights
- loglikes – if loading from arrays of samples, corresponding list of arrays of -log(likelihood)
- ignore_lines – Amount of lines at the start of the file to ignore, None not to ignore any
Returns: True if loaded successfully, False if none loaded
-
makeSingle
()¶ Combines separate chains into one samples array, so self.samples has all the samples and this instance can then be used as a general
WeightedSamples
instance.Returns: self
-
makeSingleSamples
(filename='', single_thin=None, random_state=None)[source]¶ Make file of unit weight samples by choosing samples with probability proportional to their weight.
If you just want the indices of the samples use
random_single_samples_indices()
instead.Parameters: - filename – The filename to write to, leave empty if no output file is needed
- single_thin – factor to thin by; if not set generates as many samples as it can up to self.max_scatter_points
- random_state – random seed or Generator
Returns: numpy array of selected weight-1 samples if no filename
-
mean
(paramVec, where=None)¶ Get the mean of the given parameter vector.
Parameters: - paramVec – array of parameter values or int index of parameter to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: parameter mean
-
mean_diff
(paramVec, where=None)¶ Calculates an array of differences between a parameter vector and the mean parameter value
Parameters: - paramVec – array of parameter values or int index of parameter to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: array of p_i - mean(p_i)
-
mean_diffs
(pars: Union[None, int, Sequence] = None, where=None) → Sequence¶ Calculates a list of parameter vectors giving distances from parameter means
Parameters: - pars – if specified, list of parameter vectors or int parameter indices to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: list of arrays p_i-mean(p-i) for each parameter
-
parLabel
(i)[source]¶ Gets the latex label of the parameter
Parameters: i – The index or name of a parameter. Returns: The parameter’s label.
-
parName
(i, starDerived=False)[source]¶ Gets the name of i’th parameter
Parameters: - i – The index of the parameter
- starDerived – add a star at the end of the name if the parameter is derived
Returns: The name of the parameter (string)
-
random_single_samples_indices
(random_state=None, thin: Optional[float, None] = None, max_samples: Optional[int, None] = None)¶ Returns an array of sample indices that give a list of weight-one samples, by randomly selecting samples depending on the sample weights
Parameters: - random_state – random seed or Generator
- thin – additional thinning factor (>1 to get fewer samples)
- max_samples – optional parameter to thin to get a specified mean maximum number of samples
Returns: array of sample indices
-
readChains
(files_or_samples, weights=None, loglikes=None)[source]¶ Loads samples from a list of files or array(s), removing burn in, deleting fixed parameters, and combining into one self.samples array
Parameters: - files_or_samples – The list of file names to read, samples or list of samples
- weights – array of weights if setting from arrays
- loglikes – array of -log(likelihood) if setting from arrays
Returns: self.
-
removeBurn
(remove=0.3)¶ removes burn in from the start of the samples
Parameters: remove – fraction of samples to remove, or if int >1, the number of sample rows to remove
-
removeBurnFraction
(ignore_frac)¶ Remove a fraction of the samples as burn in
Parameters: ignore_frac – fraction of sample points to remove from the start of the samples, or each chain if not combined
-
reweightAddingLogLikes
(logLikes)¶ Importance sample the samples, by adding logLike (array of -log(likelihood values)) to the currently stored likelihoods, and re-weighting accordingly, e.g. for adding a new data constraint
Parameters: logLikes – array of -log(likelihood) for each sample to adjust
-
saveAsText
(root, chain_index=None, make_dirs=False)¶ Saves the samples as text files, including parameter names as .paramnames file.
Parameters: - root – The root name to use
- chain_index – Optional index to be used for the filename, zero based, e.g. for saving one of multiple chains
- make_dirs – True if this should (recursively) create the directory if it doesn’t exist
-
savePickle
(filename)¶ Save the current object to a file in pickle format
Parameters: filename – The file to write to
-
saveTextMetadata
(root, properties=None)[source]¶ Saves metadata about the sames to text files with given file root
Parameters: - root – root file name
- properties – optional dictiory of values to save in root.properties.ini
-
setColData
(coldata, are_chains=True)¶ Set the samples given an array loaded from file
Parameters: - coldata – The array with columns of [weights, -log(Likelihoods)] and sample parameter values
- are_chains – True if coldata starts with two columns giving weight and -log(Likelihood)
-
setDiffs
()¶ saves self.diffs array of parameter differences from the y, e.g. to later calculate variances etc.
Returns: array of differences
-
setMeans
()¶ Calculates and saves the means of the samples
Returns: numpy array of parameter means
-
setMinWeightRatio
(min_weight_ratio=1e-30)¶ Removes samples with weight less than min_weight_ratio times the maximum weight
Parameters: min_weight_ratio – minimum ratio to max to exclude
-
setParamNames
(names=None)¶ Sets the names of the params.
Parameters: names – Either a ParamNames
object, the name of a .paramnames file to load, a list of name strings, otherwise use default names (param1, param2…).
-
setParams
(obj)¶ Adds array variables obj.name1, obj.name2 etc., where obj.name1 is the vector of samples with name ‘name1’
if a parameter name is of the form aa.bb.cc, it makes subobjects so that you can reference obj.aa.bb.cc. If aa.bb and aa are both parameter names, then aa becomes obj.aa.value.
Parameters: obj – The object instance to add the parameter vectors variables Returns: The obj after alterations.
-
setRanges
(ranges)[source]¶ Sets the ranges parameters, e.g. hard priors on positivity etc. If a min or max value is None, then it is assumed to be unbounded.
Parameters: ranges – A list or a tuple of [min,max] values for each parameter, or a dictionary giving [min,max] values for specific parameter names
-
setSamples
(samples, weights=None, loglikes=None, min_weight_ratio=None)¶ Sets the samples from numpy arrays
Parameters: - samples – The sample values, n_samples x n_parameters numpy array, or can be a list of parameter vectors
- weights – Array of weights for each sample. Defaults to 1 for all samples if unspecified.
- loglikes – Array of -log(Likelihood) values for each sample
- min_weight_ratio – remove samples with weight less than min_weight_ratio of the maximum
-
std
(paramVec, where=None)¶ Get the standard deviation of the given parameter vector.
Parameters: - paramVec – array of parameter values or int index of parameter to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: parameter standard deviation.
-
thin
(factor: int)¶ Thin the samples by the given factor, giving set of samples with unit weight
Parameters: factor – The factor to thin by
-
thin_indices
(factor, weights=None)¶ Indices to make single weight 1 samples. Assumes integer weights.
Parameters: - factor – The factor to thin by, should be int.
- weights – The weights to thin, None if this should use the weights stored in the object.
Returns: array of indices of samples to keep
-
static
thin_indices_and_weights
(factor, weights)¶ Returns indices and new weights for use when thinning samples.
Parameters: - factor – thin factor
- weights – initial weight (counts) per sample point
Returns: (unique index, counts) tuple of sample index values to keep and new weights
-
twoTailLimits
(paramVec, confidence)¶ Calculates two-tail equal-area confidence limit by counting samples in the tails
Parameters: - paramVec – array of parameter values or int index of parameter to use
- confidence – confidence limit to calculate, e.g. 0.95 for 95% confidence
Returns: min, max values for the confidence interval
-
updateBaseStatistics
()[source]¶ Updates basic computed statistics (y, covariance etc.), e.g. after a change in samples or weights
Returns: self
-
updateRenames
(renames)¶ Updates the renames known to each parameter with the given dictionary of renames.
-
updateSettings
(settings: Optional[Mapping[str, Any], None] = None, ini: Union[None, str, getdist.inifile.IniFile] = None, doUpdate=True)[source]¶ Updates settings from a .ini file or dictionary
Parameters: - settings – A dict containing settings to set, taking preference over any values in ini
- ini – The name of .ini file to get settings from, or an
IniFile
instance; by default uses current settings - doUpdate – True if we should update internal computed values, False otherwise (e.g. if you want to make other changes first)
-
var
(paramVec, where=None)¶ Get the variance of the given parameter vector.
Parameters: - paramVec – array of parameter values or int index of parameter to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: parameter variance
-
weighted_sum
(paramVec, where=None)¶ Calculates the weighted sum of a parameter vector, sum_i w_i p_i
Parameters: - paramVec – array of parameter values or int index of parameter to use
- where – if specified, a filter for the samples to use (where x>=5 would mean only process samples with x>=5).
Returns: weighted sum
-
weighted_thin
(factor: int)¶ Thin the samples by the given factor, giving (in general) non-unit integer weights. This function also preserves separate chains.
Parameters: factor – The (integer) factor to thin by
-
writeCorrelationMatrix
(filename=None)[source]¶ Write the correlation matrix to a file
Parameters: filename – The file to write to, If none writes to file_root.corr
-
exception
getdist.mcsamples.
MCSamplesError
[source] An Exception that is raised when there is an error inside the MCSamples class.
-
exception
getdist.mcsamples.
ParamError
[source] An Exception that indicates a bad parameter.
-
exception
getdist.mcsamples.
SettingError
[source] An Exception that indicates bad settings.