Model dispersion with PRISM; an alternative to MCMC for rapid analysis of models

This is the documentation for the PRISM package, an efficient and rapid alternative to MCMC methods for optimizing and analyzing scientific models. PRISM was made by Ellert van der Velden (@1313e) as part of a Ph.D under supervision of A/Prof. Alan Duffy at Swinburne University of Technology. It is written in pure Python 2/3 and publicly available on GitHub.

Warning

This is the documentation of the v1.0.x branch of PRISM, which still supports Python 2.7. Starting with v1.1.0, PRISM no longer supports Python 2.7, and this branch will not be maintained.

The documentation of PRISM is spread out over several sections:

Introduction

Rapid technological advancements allow for both computational resources and observational/experimental instruments to become better, faster and more precise with every passing year. This leads to an ever-increasing amount of scientific data being available and more research questions being raised. As a result, scientific models that attempt to address these questions are becoming more abundant, and are pushing the available resources to the limit as these models incorporate more complex science and more closely resemble reality.

However, as the number of available models increases, they also tend to become more distinct, making it difficult to keep track of their individual qualities. A full analysis of every model would be required in order to recognize these qualities. It is common to employ Markov chain Monte Carlo (MCMC) methods and Bayesian statistics for performing this task. However, as these methods are meant to be used for making approximations of the posterior probability distribution function, there must be a more efficient way of analyzing them.

PRISM tries to tackle this problem by using the Bayes linear approach, the emulation technique and history matching to construct an approximation (‘emulator’) of any given model. The use of these techniques can be seen as special cases of Bayesian statistics, where limited model evaluations are combined with advanced regression techniques, covariances and probability calculations. PRISM is designed to easily facilitate and enhance existing MCMC methods by restricting plausible regions and exploring parameter space efficiently. However, PRISM can additionally be used as a standalone alternative to MCMC for model analysis, providing insight into the bahvior of complex scientific models. With PRISM, the time spent on evaluating a model is minimized, providing developers with an advanced model analysis for a fraction of the time required by more traditional methods.

Why use PRISM?

  • Written in pure Python 2/3, for versatility;
  • Stores results in HDF5-files, allowing for easy user-access;
  • Can be executed in serial or MPI, on any number of processes;
  • Compatible with Windows, Mac OS and Unix-based machines;
  • Accepts any type of model and comparison data;
  • Built as a plug-and-play tool: all main classes can also be used as base classes;
  • Easily linked to any model by writing a single custom ModelLink subclass;
  • Capable of reducing relevant parameter space by factors over 100,000 using only a few thousand model evaluations;
  • Can be used alone for analyzing models, or combined with MCMC for efficient model parameter estimations.

Getting started

Installation

PRISM can be easily installed by either cloning the repository and installing it manually:

$ git clone https://github.com/1313e/PRISM
$ cd PRISM
$ pip install .

or by installing it directly from PyPI with:

$ pip install prism

PRISM can now be imported as a package with import prism. For using PRISM in MPI, mpi4py >= 3.0.0 is required (not installed automatically).

The PRISM package comes with two ModelLink subclasses. These ModelLink subclasses can be used to experiment with PRISM to see how it works. Using PRISM has several examples explaining the different functionalities of the package.

Example usage

See Minimal example for a documented explanation on this example.

# Imports
from prism import Pipeline
from prism.modellink import GaussianLink

# Define model data and create ModelLink object
model_data = {3: [3.0, 0.1], 5: [5.0, 0.1], 7: [3.0, 0.1]}
modellink_obj = GaussianLink(model_data=model_data)

# Create Pipeline object
pipe = Pipeline(modellink_obj)

# Construct first iteration of the emulator
pipe.construct()

# Create projections
pipe.project()

The PRISM pipeline

The *PRISM* pipeline

The structure of the PRISM pipeline.

The overall structure of PRISM can be seen in Fig. 1 and will be discussed below. The Pipeline object plays a key-role in the PRISM framework as it governs all other objects and orchestrates their communications and method calls. It also performs the process of history matching and refocusing (as explained in ref{subsec:History matching}). It is linked to the model by a user-written ModelLink object (see ModelLink: A crash course), allowing the Pipeline object to extract all necessary model information and call the model. In order to ensure flexibility and clarity, the PRISM framework writes all of its data to one or several HDF5-files using h5py, as well as numpy.

The analysis of a provided model and the construction of the emulator systems for every output value, starts and ends with the Pipeline object. When a new emulator is requested, the Pipeline object creates a large Latin-Hypercube design (LHD) of model evaluation samples to get the construction of the first iteration of the emulator systems started. To ensure that the maximum amount of information can be obtained from evaluating these samples, a custom Latin-Hypercube sampling code was written. This produces LHDs that attempt to satisfy both the maximin criterion as well as the correlation criterion. This code is customizable through PRISM and publicly available in the e13Tools Python package.

This Latin-Hypercube design is then given to the Model Evaluator, which through the provided ModelLink object evaluates every sample. Using the resulting model outputs, the Active Parameters for every emulator system (individual data point) can now be determined. Next, depending on the user, polynomial functions will be constructed by performing an extensive Regression process for every emulator system, or this can be skipped in favor of a sole Gaussian analysis (faster, but less accurate). No matter the choice, the emulator systems now have all the required information to be constructed, which is done by calculating the Prior Expectation and Prior Covariance values for all evaluated model samples (\(\mathrm{E}(D_i)\) and \(\mathrm{Var}(D_i)\)).

Afterward, the emulator systems are fully constructed and are ready to be evaluated and analyzed. Depending on whether the user wants to prepare for the next emulator iteration or create a projection (see Projections), the Emulator Evaluator creates one or several LHDs of emulator evaluation samples, and evaluates them in all emulator systems, after which an Implausibility Check is carried out. The samples that survive the check can then either be used to construct the new iteration of emulator systems by sending them to the Model Evaluator, or they can be analyzed further by performing a Projection. The Pipeline object performs a single cycle by default (to allow for user-defined analysis algorithms), but can be easily set to continuously cycle.

In addition to the above, PRISM also features a high-level Message Passing Interface (MPI) implementation using the Python package mpi4py. All emulator systems in PRISM can be constructed independently from each other, in any order, and only require to communicate when performing the implausibility cut-off checks during history matching. Additionally, since different models and/or architectures require different amounts of computational resources, PRISM can run on any number of MPI processes (including a single one in serial to accommodate for OpenMP codes) and the same emulator can be used on a different number of MPI processes than it was constructed on (e.g., constructing an emulator using 8 MPI processes and reloading it with 6). More details on the MPI implementation and its scaling can be found in MPI implementation.

In Using PRISM and ModelLink: A crash course, the various components of PRISM are described more extensively.

MPI implementation

Given that most scientific models are either already parallelized or could benefit from parallelization, we had to make sure that PRISM allows for both MPI and OpenMP coded models to be connected. Additionally, since individual emulator systems in an emulator iteration are independent of each other, the extra CPUs required for the model should also be usable by the emulator. For that reason, PRISM features a high-level MPI implementation for using MPI-coded models, while the Python package NumPy handles the OpenMP side. A mixture of both is also possible (due to the worker_mode context manager).

Here, we discuss the MPI scaling tests that were performed on PRISM. For the tests, the same GaussianLink class was used as in Minimal example, but this time with \(32\) emulator systems (comparison data points) instead of \(3\). In PRISM, all emulator systems are spread out over the available number of MPI processes as much as possible while also trying to balance the number of calculations performed per MPI process. Since all emulator systems are stored in different HDF5-files, it is possible to reinitialize the Pipeline using the same Emulator class and ModelLink subclass on a different number of MPI processes. To make sure that the results are not influenced by the variation in evaluation rates, we constructed an emulator of the Gaussian model and used the exact same emulator in every test.

The tests were carried out using any number of MPI processes between \(1\) and \(32\), and using a single OpenMP thread each time for consistency. We generated a Latin-Hypercube design of \(3\cdot10^6\) samples and measured the average evaluation rate of the emulator using the same Latin-Hypercube design each time. To take into account any variations in the evaluation rate caused by initializations, this test was performed \(20\) times. As a result, this Latin-Hypercube design was evaluated in the emulator a total of \(640\) times, giving an absolute total of \(1.92\cdot10^9\) emulator evaluations.

MPI scaling test results.

Figure showing the MPI scaling of PRISM using the emulator of a simple Gaussian model with \(32\) emulator systems. The tests involved analyzing a Latin-Hypercube design of \(3\cdot10^6\) samples in the emulator, determining the average evaluation rate and executing this a total of \(20\) times using the same sample set every time. The emulator used for this was identical in every instance. Left axis: The average evaluation rate of the emulator vs. the number of MPI processes it is running on. Right axis: The relative speed-up factor vs. the number of MPI processes, which is defined as \(\frac{f(x)}{f(1)\cdot x}\) with \(f(x)\) the average evaluation rate and \(x\) the number of MPI processes. Dotted line: The minimum acceptable relative speed-up factor, which is always \(1/x\). Dashed line: A straight line with a slope of \({\sim}0.645\), connecting the lowest and highest evaluation rates. The tests were performed using the OzSTAR computing facility at the Swinburne University of Technology, Melbourne, Australia.

In Fig. 2, we show the results of the performed MPI scaling tests. On the left y-axis, the average evaluation rate vs. the number of MPI processes the test ran on is plotted, while the relative speed-up factor vs. the number of MPI processes is plotted on the right y-axis. The relative speed-up factor is defined as \(f(x)/(f(1)\cdot x)\) with \(f(x)\) the average evaluation rate and \(x\) the number of MPI processes. The ideal MPI scaling would correspond to a relative speed-up factor of unity for all \(x\).

In this figure, we can see the effect of the high-level MPI implementation. Because the emulator systems are spread out over the available MPI processes, the evaluation rate is mostly determined by the runtime of the MPI process with the highest number of systems assigned. Therefore, if the number of emulator systems (\(32\) in this case) cannot be divided by the number of available MPI processes, the speed gain is reduced, leading to the plateaus like the one between \(x=16\) and \(x=31\). Due to the emulator systems not being the same, their individual evaluation rates are different such that a different evaluation rate has a bigger effect on the average evaluation rate of the emulator the more MPI processes there are. This is shown by the straight dashed line drawn between \(f(1)\) and \(f(32)\), which has a slope of \({\sim}0.645\).

The relative speed-up factor shows the efficiency of every individual MPI process in a specific run, compared to using a single MPI process. This also shows the effect of the high-level MPI implementation, giving peaks when the maximum number of emulator systems per MPI process has decreased. The dotted line shows the minimum acceptable relative speed-up factor, which is always defined as \(1/x\). On this line, the average evaluation rate \(f(x)\) for any given number of MPI processes is always equal to \(f(1)\).

Using PRISM

Here, various different aspects of how the PRISM package can be used are described.

Minimal example

A minimal example on how to initialize and use the PRISM pipeline is shown here. First, one has to import the Pipeline class and a ModelLink subclass:

>>> from prism import Pipeline
>>> from prism.modellink import GaussianLink

Normally, one would import a custom-made ModelLink subclass, but for this example one of the two ModelLink subclasses that come with the PRISM package is used (see Writing a ModelLink subclass for the basic structure of writing a custom ModelLink subclass).

Next, the ModelLink should be initialized, which is the GaussianLink class in this case. In addition to user-defined arguments, every ModelLink subclass takes two optional arguments, model_parameters and model_data. The use of either one will add the provided parameters/data to the default parameters/data defined in the class. Since the GaussianLink class does not have default data defined, it is required to supply it with some data during initialization (using an array, dict or external file):

>>> # f(3) = 3.0 +- 0.1, f(5) = 5.0 +- 0.1, f(7) = 3.0 +- 0.1
>>> model_data = {3: [3.0, 0.1], 5: [5.0, 0.1], 7: [3.0, 0.1]}
>>> modellink_obj = GaussianLink(model_data=model_data)

Here, the GaussianLink class was initialized by giving it three custom data points and using its default parameters. One can check this by looking at the representation of this GaussianLink object:

>>> modellink_obj
GaussianLink(model_parameters={'A1': [1.0, 10.0, 5.0], 'B1': [0.0, 10.0, 5.0],
                               'C1': [0.0, 5.0, 2.0]},
             model_data={7: [3.0, 0.1], 5: [5.0, 0.1], 3: [3.0, 0.1]})

The Pipeline class takes several optional arguments, which are mostly paths and the type of Emulator class that must be used. It also takes one mandatory argument, which is an instance of the ModelLink subclass to use. Since it has already been initialized above, the Pipeline class can be initialized:

>>> pipe = pipeline(modellink_obj)
>>> pipe
Pipeline(GaussianLink(model_parameters={'A1': [1.0, 10.0, 5.0], 'B1': [0.0, 10.0, 5.0],
                                        'C1': [0.0, 5.0, 2.0]},
                      model_data={7: [3.0, 0.1], 5: [5.0, 0.1], 3: [3.0, 0.1]}),
         working_dir='prism_0')

Since no working directory was provided to the Pipeline class and none already existed, it automatically created one (prism_0).

PRISM is now completely ready to start emulating the model. The Pipeline allows for all steps in a full cycle (see PRISM pipeline) to be executed automatically:

>>> pipe.run()

which is equivalent to:

>>> pipe.construct(analyze=False)
>>> pipe.analyze()
>>> pipe.project()

This will construct the next iteration (first in this case) of the emulator, analyze it to check if it contains plausible regions and make projections of all active parameters. The current state of the Pipeline object can be viewed by calling the details() method (called automatically after most user-methods), which gives an overview of many properties that the Pipeline object currently has.

This is all that is required to construct an emulator of the model of choice. All user-methods, with one exception (evaluate()), solely take optional arguments and perform the operations that make the most sense given the current state of the Pipeline object if no arguments are given. These arguments allow for one to modify the performed operations, like reconstructing/reanalyzing previous iterations, projecting specific parameters, evaluating the emulator and more.

Projections

After having made an emulator of a given model, PRISM can show the user the knowledge it has about the behavior of this model by making projections of the active parameters in a specific emulator iteration. These projections are created by the project() method, which has many different properties and options. For showing them below, the same emulator as the one in Minimal example is used.

Properties

Projections (and their figures) are made by analyzing a large set of evaluations samples. For 3D projections, this set is made up of a grid of proj_res x proj_res samples for the plotted (active) parameters, where the values for the remaining parameters in every individual grid point are given by an LHD of proj_depth samples. This gives the total number of analyzed samples as proj_res x proj_res x proj_depth.

Every sample in the sample set is then analyzed in the emulator, saving whether or not this sample is plausible and what the implausibility value at the first cut-off is (the first value in impl_cut). This yields proj_depth results per grid point, which can be used to determine the fraction of samples that is plausible and the minimum implausibility value at the first cut-off in this point. Doing this for the entire grid and interpolating them, creates a map of results that is independent of the values of the non-plotted parameters. For 2D projections, it works the same way, except that only a single active parameter is plotted.

Note

When using a 2D model, the projection depth used to make a 2D projection will be proj_depth, which is to be expected. However, when using an nD model, the projection depth of a 2D projection is equal to proj_res x proj_depth. This is to make sure that for an nD model, the density of samples in a 2D projection is the same as in a 3D projection.

The project() method solely takes optional arguments. Calling it without any arguments will produce six projection figures: three 2D projections and three 3D projections. One of each type is shown below.

2D projection figure of model parameter A1.

2D projection figure of model parameter \(A_1\). The vertical dashed line shows the parameter estimate of \(A_1\), whereas the horizontal red line shows the first implausibility cut-off value.

3D projection figure of model parameters A1 and B1.

3D projection figure of model parameters \(A_1\) and \(B_1\). The dashed lines show the estimates of both parameters.

A projection figure is made up of two subplots. The upper subplot shows a map of minimum implausibility values that can be reached for any given value (combination) of the plotted parameter(s). The lower subplot gives a map of the fraction of samples that is plausible in a specified point on the grid (called “line-of-sight depth” due to the way it is calculated). Another way of describing this map is that it gives the probability that a parameter set with given plotted value(s) is plausible.

Both projection types have a different purpose. A 3D projection gives insight into what the dependencies (or correlations) are between the two plotted parameters, by showing where the best (top) and most (bottom) plausible samples can be found. On the other hand, a 2D projection is quite similar in meaning to a maximum likelihood optimization performed by MCMC methods, with the difference being that the projection is based on expectations rather than real model output. A combination of both subplots allows for many model properties to be derived, especially when they do not agree with each other.

Options

The project() method takes two (optional) arguments, emul_i and proj_par. The first controls which emulator iteration should be used, while the latter provides the model parameters of which projections need to be made. Since it only makes sense to make projections of active parameters, all passive parameters are filtered out of proj_par. The remaining parameters are then used to determine which projections are required (which also depends on the requested projection types). For example, if one wishes to only obtain projections of the \(A_1\) and \(B_1\) parameters (which are both active) in iteration 1, then this can be done with:

>>> pipe.project(1, ('A1', 'B1'))

This would generate the figures shown above, as well as the 2D projection figure of \(B_1\). By default, the last constructed emulator iteration and all model parameters are requested.

The remaining input arguments can only be given as keyword arguments, since they control many different aspects of the project() method. The proj_type argument controls which projection types to make. For 2D models, this is always ‘2D’ and cannot be modified. However, for nD models, this can be set to ‘2D’ (only 2D projections), ‘3D’ (only 3D projections) or ‘both’ (both 2D and 3D projections). By default, it is set to ‘both’.

The figure argument is a bool, that determines whether or not the projection figures should be created after calculating the projection data. If True, the projection figures will be created and saved, which is done by default. If False, the data that is contained within the projection figures will be calculated and returned in a dict. This allows the user to either let PRISM create the projection figures using the standard template or create the figures themselves.

The align argument controls the alignment of the subplots in every projection figure. By default, it aligns the subplots in a column (‘col’), as shown in the figures above. Aligning the subplots in a row (‘row’) would give Fig. 3 as the figure below.

2D projection figure of model parameter A1 with the 'row' alignment.

2D projection figure of model parameter \(A_1\) with the ‘row’ alignment.

The smooth argument is also a bool, that determines what to do if a grid point in the projection figure contains no plausible samples, but does contain a minimum implausibility value below the first non-wildcard cut-off. If False, which is the default, these values are kept in the figure, which may show up as artifact-like features. If True, these values are set to the first cut-off, basically removing them from the projection figure. This may however also remove interesting features. Below are two identical projections, one that is smoothed and one that is not, to showcase this difference (these projections are from the second iteration, since this effect rarely occurs in the first iteration).

Non-smoothed 3D projection figure of model parameters A1 and B1.

Non-smoothed 3D projection figure of model parameters \(A_1\) and \(B_1\).

Smoothed 3D projection figure of model parameters A1 and B1.

Smoothed 3D projection figure of model parameters \(A_1\) and \(B_1\).

In these figures, one can see that the non-smoothed projection shows many features in the upper subplot that look like artifacts. These features are however not artifacts, but caused by a sample (or samples) having its highest implausibility value being below the first implausibility cut-off, but still being implausible due to failing a later cut-off. For example, if the implausibility cut-offs are [4.0, 3.7, 3.5] and a sample has implausibility values [3.9, 3.8, 3.2], it is found implausible due to failing to meet the second cut-off. However, since the first value is still the highest implausibility value, that value is used in the projection figure. Smoothing figures usually allows for 3D projections (2D projections rarely show this) to become less crowded, but they do throw away information. It should therefore only be used when necessary.

The force argument is a bool, which controls what to do if a projection is requested for which data already exists. If False (default), it will use the previously acquired projection data to create the projection figure if it does not exist, skip if it does or return the figure data if figure is False. If True, the projection data and all associated projection figures will be deleted, and the projection will be recalculated.

The remaining six arguments are keyword argument dicts, that need to be passed to the various different plotting functions that are used for creating the projection figures. The fig_kwargs dict is passed to the figure() function when creating the projection figure instance. The impl_kwargs_2D and los_kwargs_2D dicts are passed to the plot() function when making the minimum implausibility and line-of-sight depth subplots, respectively, for the 2D projections. Similarly, the impl_kwargs_3D and los_kwargs_3D dicts are passed to the hexbin() function for 3D projections. And, finally, the line_kwargs dict is passed to the draw() function for drawing the parameter estimate lines.

Dual nature (normal/worker mode)

PRISM features a high-level MPI implementation, as described in MPI implementation: all user-methods and most major methods are to be executed by all MPI ranks at the same time, and PRISM will automatically distribute the work among the available ranks within this function/method. This allows for PRISM to be used with both serial and parallel models, by setting the MPI_call flag accordingly, while also allowing for the same code to be used in serial and parallel. However, given that the emulator of PRISM can be very useful for usage in other routines, like Hybrid sampling, an external code will call PRISM’s methods. In order to use PRISM in parallel with a parallelized model, this code would have to call PRISM with all MPI ranks simultaneously at all times, which may not always be possible (e.g., when using MCMC methods).

Therefore, PRISM has a dual execution/call nature, where it can be switched between two different modes. In the default mode, PRISM works as described before, where all MPI ranks call the same user-code. However, by using the worker_mode context manager, all code within will be executed in worker mode. When in worker mode, all worker ranks are continously listening for calls from the controller rank, made with the _make_call() method. They will continue to do so until the controller makes the pipe._make_call(None) call, which is done automatically when the controller exits worker_mode. Manually making this call should solely be done in advanced use-cases.

In worker_mode, one uses the following structure (assuming that the Pipeline instance is called pipe):

# Code to be executed in default mode

with pipe.worker_mode:
    if pipe.is_controller:
        # Code to be executed in worker mode

# More code to be executed in default mode

Note

All code that is inside the worker_mode context manager should solely be executed by the controller rank. If not, all worker ranks will execute this code after the controller ranks exits the context manager. Currently, it is not possible to make a context manager handle this automatically (the rejected PEP 377 describes this perfectly).

The _make_call() method accepts almost anything that can be called. It can also be used when not in worker_mode, in which case it works the exact same way for all MPI ranks. Its sole limitation is that all supplied arguments must be pickleable (e.g., compiled code objects are NOT pickleable due to safety reasons), both when used in worker_mode and outside of it.

The worker_mode can be used in a variety of ways, as described below. It can be used to access any attribute of the Pipeline instance:

with pipe.worker_mode:
    if pipe.is_controller:
        # Construct first emulator iteration
        pipe._make_call('construct', 1)

        # Print latest constructed emulator iteration
        print(pipe._make_call('emulator._get_emul_i', 1, 0))

        # Make a specific projection with the 'row' alignment
        pipe._make_call('project', 1, (0, 1), align='row')

which is equivalent to:

# Construct first emulator iteration
pipe.construct(1)

# Print latest constructed emulator iteration
print(pipe.emulator._get_emul_i(1, 0))

# Make a specific projection with the 'row' alignment
pipe.project(1, (0, 1), align='row')

The above two code snippets are equal to each other, and the worker_mode will most likely be used very rarely in this fashion. However, by supplying the _make_call() method with a callable function (that can be pickled), externally defined functions can be executed:

# Enable worker mode
with pipe.worker_mode:
    if pipe.is_controller:
        # Import print function that prepends MPI rank to message
        from prism._internal import rprint

        # Make call to use this function
        # Equivalent to 'rprint("Reporting in.")'
        pipe._make_call(rprint, "Reporting in.")

This is especially useful when one combines a serial code with PRISM, but wants PRISM to execute in MPI. An application example of this is Hybrid sampling.

Hybrid sampling

A common problem when using MCMC methods is that it can often take a very long time for MCMC to find its way on the posterior probability distribution function, which is often referred to as the burn-in phase. This is because, when considering a parameter set, there is usually no prior information that this parameter set is (un)likely to result into a desirable model realization. This means that such a parameter set must first be evaluated in the model before any probabilities can be calculated. However, by constructing an emulator of the model, one can use it as an additional prior for the posterior probability calculation. Therefore, although PRISM is primarily designed to make analyzing models much more efficient and accessible than normal MCMC methods, it is also very capable of enhancing them. This process is called hybrid sampling, which can be performed easily with the utils module and will be explained below.

Algorithm

Hybrid sampling allows one to use PRISM to first analyze a model’s behavior, and later use the gathered information to speed up parameter estimations (by using the emulator as an additional prior in a Bayesian analysis). Hybrid sampling works in the following way:

  1. Whenever an MCMC walker proposes a new sample, it is first passed to the emulator of the model;
  2. If the sample is not within the defined parameter space, it automatically receives a prior probability of zero (or \(-\infty\) in case of logarithmic probabilities). Else, it will be evaluated in the emulator;
  3. If the sample is labeled as implausible by the emulator, it also receives a prior probability of zero. If it is plausible, the sample is evaluated in the same way as for normal sampling;
  4. Optionally, a scaled value of the first implausibility cut-off is used as an additional (non-zero) prior probability. This can be enabled by using the hybrid input argument for the get_lnpost_fn() function.

Since the emulator that PRISM makes of a model is not defined outside of the parameter space given by par_rng, the second step is necessary to make sure the results are valid. There are several advantages of using hybrid sampling over normal sampling:

  • Acceptable samples are guaranteed to be within plausible space;
  • This in turn makes sure that the model is only evaluated for plausible samples, which heavily reduces the number of required evaluations;
  • No burn-in phase is required, as the starting positions of the MCMC walkers are chosen to be in plausible space;
  • As a consequence, varying the number of walkers tends to have a much lower negative impact on the convergence probability and speed;
  • Samples with low implausibility values can optionally be favored.

Usage

In order to help the user with combining PRISM with MCMC to use hybrid sampling, the utils module provides two functions: get_walkers() and get_lnpost_fn(). The get_walkers() function analyzes a set of proposed init_walkers and returns the positions that are valid (and the number of positions that are valid). By default, it uses the available impl_sam of the last constructed iteration, but it can also be supplied with a custom set of proposed walkers or an integer stating how many proposed walkers the function should check:

>>> # Use impl_sam if it is available
>>> n, p0 = get_walkers(pipe)

>>> # Request 2000 proposed samples
>>> n_walkers = 2000
>>> n, p0 = get_walkers(pipe, init_walkers=n_walkers)

>>> # Use custom init_walkers
>>> from numpy.random import rand
>>> init_walkers = rand(n_walkers, pipe.modellink.n_par)
>>> n, p0 = get_walkers(pipe, init_walkers=init_walkers)

Given that most sampling methods operate in unit space (all parameter values are between zero and unity), the get_walkers() function automatically assumes that all samples are defined in unit space (and therefore using rand() works). It can however be disabled by setting the unit_space input argument to False. One has to keep in mind that, because of the way the emulator works, there is no guarantee for a specific number of valid walkers to be obtained. Having the desired emulator iteration already analyzed may give an indication how many samples in total need to be proposed to be left with a specific number.

When the initial positions of the MCMC walkers have been determined, one can use them in an MCMC parameter estimation algorithm, avoiding the burn-in phase. This in itself can already be very useful, but it does not allow for hybrid sampling yet. Most MCMC methods require the definition of an lnpost function, which takes a parameter set and returns the corresponding natural logarithm of the posterior probability. In order to do hybrid sampling, this lnpost function must have the algorithm described above implemented.

The get_lnpost_fn() function factory provides exactly that. It takes a user-defined ext_lnpost function and a Pipeline object, and returns a function definition get_lnpost(par_set, *args, **kwargs). This get_lnpost function first analyzes a proposed par_set in the emulator, passes par_set (along with any additional arguments) to ext_lnpost if the sample is plausible, or returns \(-\infty\) if it is not. The return-value of the ext_lnpost function is then returned by the get_lnpost function as well. To make sure that the get_lnpost function can be used in both execution modes (see Dual nature (normal/worker mode)), all parallel calls to the Pipeline object are done with the _make_call() method.

The use of a function factory here allows for all input arguments to be validated once and then saved as local variables for the get_lnpost function. Not only does this avoid that all arguments have to be provided and validated for every individual call, but it also ensures that the same arguments are used every time, as local variables of a function cannot be modified by anything. Since users most likely use get_walkers() and get_lnpost_fn() frequently together, the get_walkers() function allows for the ext_lnpost argument to be supplied to it. This will automatically call the get_lnpost_fn() function factory using the provided ext_lnpost and the same input arguments given to get_walkers(), and return the obtained get_lnpost function in addition to the initial MCMC walkers.

Application

Using the information above, using hybrid sampling on a model of choice can be done quite easily. For performing the MCMC analysis, we will be using the emcee package in this example.

Assume that we want to first analyze and then optimize the Gaussian model given by the GaussianLink class. So, we first have to make an emulator of the model:

>>> from prism import Pipeline
>>> from prism.modellink import GaussianLink
>>> model_data = {3: [3.0, 0.1], 5: [5.0, 0.1], 7: [3.0, 0.1]}
>>> modellink_obj = GaussianLink(model_data=model_data)
>>> pipe = Pipeline(modellink_obj)
>>> pipe.construct()

Using the constructed emulator, we can perform a model parameter optimization using hybrid sampling. For this, we need to define an ext_lnpost function, for which we will use a simple Gaussian probability function:

def ext_lnpost(par_set, pipe):
    # Create parameter dict for call_model
    par_dict = dict(zip(pipe.modellink.par_name, par_set))

    # Use wrapped model to obtain model output
    mod_out = pipe.modellink.call_model(pipe.emulator.emul_i,
                                        par_dict,
                                        pipe.modellink.data_idx)

    # Get the model and data variances
    # Since the value space is linear, the data error is centered
    md_var = pipe.modellink.get_md_var(pipe.emulator.emul_i,
                                       par_dict,
                                       pipe.modellink.data_idx)
    data_var = [err[0]**2 for err in pipe.modellink.data_err]

    # Calculate the posterior probability and return it
    sigma_2 = md_var+data_var
    diff = pipe.modellink.data_val-mod_out
    return(-0.5*(np.sum(diff**2/sigma2)))

Since the Pipeline object already has the model wrapped and linked, we used that to evaluate the model. The GaussianLink class has a centered data error, therefore we can take the upper bound for every error when calculating the variance. However, for more complex models, this is probably not true.

Next, we have to obtain the starting positions for the MCMC walkers. Since we want to do hybrid sampling, we can obtain the get_lnpost function at the same time as well:

>>> from prism.utils import get_walkers
>>> n, p0, get_lnpost = get_walkers(pipe, unit_space=False,
                                    ext_lnpost=ext_lnpost, hybrid=True)

By setting hybrid to True, we use the implausibility cut-off value as an additional prior. Now we only still need the EnsembleSampler class and NumPy (for the ext_lnpost function):

>>> import numpy as np
>>> from emcee import EnsembleSampler

Now we have everything that is required to perform a hybrid sampling analysis. In most cases, MCMC methods require to be executed on only a single MPI rank, so we will use the worker_mode:

# Activate worker mode
with pipe.worker_mode:
    if pipe.is_controller:
        # Create EnsembleSampler object
        sampler = EnsembleSampler(n, pipe.modellink.n_par,
                                  get_lnpost, args=[pipe])

        # Run mcmc for 1000 iterations
        sampler.run_mcmc(p0, 1000)

        # Execute any custom operations here
        # For example, saving the chain data or plotting the results

And that is basically all that is required for using PRISM together with MCMC. For a normal MCMC approach, the same code can be used, except that one has to use ext_lnpost instead of get_lnpost (and, obtain the starting positions of the walkers in a different way).

General usage rules

Below is a list of general usage rules that apply to PRISM.

  • Unless specified otherwise in the documentation, all input arguments in the PRISM package that accept…

    • a bool (True/False) also accept 0/1 as a valid input;

    • None indicate a default value or operation for obtaining this input argument. In most of these cases, the default value depends on the current state of the PRISM pipeline, and therefore a small operation is required for obtaining this value;

      Example

      Providing None to pot_active_par, where it indicates that all model parameters should be potentially active.

    • the names of model parameters also accept the internal indices of these model parameters. The index is the order in which the parameter names appear in the par_name list or as they appear in the output of the details() method;

    • a sequence of integers, floats and/or strings will accept (almost) any formatting including most special characters as separators as long as they do not have any meaning (like a dot for floats or valid escape sequences for strings);

      Example

      The following sequences are equal:
      • A, 1, 2.0, B;
      • [A,1,2.,B];
      • “A 1 2.0 B”;
      • “’[“ (A / }| \n; <1{}) \,,”>2.000000 !! \t< )?%\B .
    • the path to a data file (PRISM parameters, model parameters, model data) will read in all the data from that file as a Python dict, with a colon : acting as the separator between the key and value.

  • Depending on the used emulator type, state of loaded emulator and the PRISM parameter values, it is possible that providing values for certain PRISM parameters has no influence on the outcome of the pipeline. This can be either because they have non-changeable default values or are simply not used anywhere (given the current state of the pipeline);

    Examples

    • If method != ‘gaussian’, it causes sigma to have no use in the pipeline;
    • Switching the bool value for use_mock while loading a constructed emulator has no effect, since the mock data is generated (or not) when constructing a new emulator and cannot be changed or swapped out afterward.
  • All docstrings in PRISM are written in RTF (Rich Text Format) and are therefore best viewed in an editor that supports it (like Spyder).

External data files

When using PRISM, there are three different cases where the path to an external data file can be provided. As mentioned in General usage rules, all external files are read-in as a Python dict, with the colon being the separator between the key and value. Additionally, all lines are read as strings and converted back when assigned in memory, to allow for many different mark-ups to be used. Depending on which of the three files is read-in, the keys and values have different meanings. Here, the three different files are described.

PRISM parameters file

This file contains the non-default values that must be used for the PRISM parameters. These parameters control various different functionalities of PRISM. It is provided as the prism_file argument when initializing the Pipeline class. When certain parameters are read-in depends on their type:

  • Emulator parameters: Whenever a new emulator is created;
  • Pipeline parameters: When the Pipeline class is initialized;
  • Implausibility parameters: When the analyze() method is called (saved to HDF5) or when an emulator iteration is loaded that has not been analyzed yet (not saved to HDF5);
  • Projection parameters: When the project() method is called.

The default PRISM parameters file can be found in the prism/data folder and is shown below:

n_sam_init          : 500                   # Number of initial model evaluation samples
proj_res            : 25                    # Number of projected grid points per model parameter
proj_depth          : 250                   # Number of emulator evaluation samples per projected grid point
base_eval_sam       : 800                   # Base number for growth in number of model evaluation samples
sigma               : 0.8                   # Gaussian sigma/standard deviation (only required if method == 'gaussian')
l_corr              : 0.3                   # Gaussian correlation length(s)
impl_cut            : [0.0, 4.0, 3.8, 3.5]  # List of implausibility cut-off values
criterion           : 'multi'               # Criterion for constructing LHDs
method              : 'full'                # Method used for constructing the emulator
use_regr_cov        : False                 # Use regression covariance
poly_order          : 3                     # Polynomial order for regression
n_cross_val         : 5                     # Number of cross-validations for regression
do_active_anal      : True                  # Perform active parameter analysis
freeze_active_par   : True                  # Active parameters always stay active
pot_active_par      : None                  # List of potentially active parameters
use_mock            : False                 # Use mock data

In this file, the key is the name of the parameter that needs to be changed, and the value what it needs to be changed to. PRISM itself does not require this default file, as all of the default values are hard-coded, and is therefore never read-in. An externally provided PRISM parameters file is only required to have the non-default values.

Model parameters file

This file contains the non-default model parameters to use for a model. It is provided as the model_parameters input argument when initializing the ModelLink subclass (a dict or array-like can be provided instead as well). Keep in mind that the ModelLink subclass may not have default model parameters defined.

An example of the various different ways model parameter information can be provided is given below:

# name      : lower_bndry       upper_bndry     estimate
  A         : 1                 5               3
  Bravo     : 2                 7               None
  C42       : 3                 6.74

In this file, the key is the name of the model parameter and the value is a sequence of integers or floats, specifying the lower and upper boundaries of the parameter and, optionally, its estimate. The contents of this file is equal to providing the following as model_parameters during initialization of a ModelLink subclass:

# As a dict
model_parameters = {'A': [1, 5, 3],
                    'Bravo': [2, 7, None],
                    'C42': [3, 6.74]}

# As an array_like
model_parameters = [['A', [1, 5, 3]],
                    ['Bravo', [2, 7, None]],
                    ['C42', [3, 6.74]]]

# As two array_likes zipped
model_parameters = zip(['A', 'Bravo', 'C42'],
                       [[1, 5, 3], [2, 7, None], [3, 6.74]])

Providing None as the parameter estimate or not providing it at all, implies that no parameter estimate (for the corresponding parameter) should be used in the projection figures.

Model data file

This file contains the non-default model comparison data points to use for a model. It is provided as the model_data input argument when initializing the ModelLink subclass (a dict or array-like can be provided instead as well). Keep in mind that the ModelLink subclass may not have default model comparison data defined.

An example of the various different ways model comparison data information can be provided is given below:

# data_idx  : data_val          data_err        data_spc
  1, 2      : 1                 0.05    0.05    'lin'
  3.0       : 2                 0.05            'log'
  ['A']     : 3                 0.05    0.15
  1, A, 1.0 : 4                 0.05

Here, the key is the full sequence of the data identifier of a data point, where any character that is not a letter, number, minus/plus or period acts as a separator between the elements of the data identifier. The corresponding value specifies the data value, data error(s) and data value space. Braces, parentheses, brackets and many other characters can be used as mark-up in the data identifier, to make it easier for the user to find a suitable file lay-out. A full list of all characters that can be used for this can be found in prism._internal.aux_char_list.

Similarly to the model parameters, the following is equal to the contents of this file:

# As a dict
model_data = {(1, 2): [1, 0.05, 0.05, 'lin'],
              3.0: [2, 0.05, 'log'],
              'A': [3, 0.05, 0.15],
              (1, 'A', 1.0): [4, 0.05]}

# As an array_like
model_data = [[(1, 2), [1, 0.05, 0.05, 'lin']],
              [3.0, [2, 0.05, 'log']],
              ['A', [3, 0.05, 0.15]],
              [(1, 'A', 1.0), [4, 0.05]]]

# As two array_likes zipped
model_data = zip([(1, 2), 3.0, 'A', (1, 'A', 1.0)],
                 [[1, 0.05, 0.05, 'lin'], [2, 0.05, 'log'], [3, 0.05, 0.15], [4, 0.05]])

It is necessary for the data value to be provided at all times. The data error can be given as either a single value, where it assumed that the data point has a centered \(1\sigma\)-confidence interval, or as two values, where they describe the upper and lower bounds of the \(1\sigma\)-confidence interval. The data value space can be given as a string or omitted, in which it is assumed that the value space is linear.

Note

The parameter value bounds are given as [lower bound, upper bound], whereas the data errors are given as [upper error, lower error]. The reason for this is that, individually, the order for either makes the most sense. Together however, it may cause some confusion, so extra care needs to be taken.

Descriptions

Terminology

Below is a list of the most commonly used terms/abbreviations in PRISM and their meaning. If a term is also a variable within the PRISM code, its name is given in brackets.


Active emulator system [active_emul_s]
An emulator system that has a data point assigned to it.
Active parameters [active_par/active_par_data]
The set of model parameters that are considered to have significant influence on the output of the model and contribute at least one polynomial term to one/the regression function.
Adjusted expectation [adj_exp]
The prior expectation of a parameter set, with the adjustment term taken into account. It is equal to the prior expectation if the emulator system has perfect accuracy.
Adjusted values [adj_val]
The adjusted expectation and variance values of a parameter set.
Adjusted variance [adj_var]
The prior variance of a parameter set, with the adjustment term taken into account. It is zero if the emulator system has perfect accuracy.
Adjustment term
The extra term (as determined by the BLA) that is added to the prior expectation and variance values that describes all additional correlation knowledge between model realization samples.
Analysis/Analyze
The process of evaluating a set of emulator evaluation samples in the last emulator iteration and determining which samples should be used to construct the next iteration.
BLA
Abbreviation of Bayes linear approach.
Construction/Construct
The process of calculating all necessary components to describe an iteration of the emulator.
Construction check [ccheck]
A list of keywords determining which components of which emulator systems are still required to finish the construction of a specified emulator iteration.
Controller (rank)
An MPI process that controls the flow of operations in PRISM and distributes work to all workers and itself. By default, a controller also behaves like a worker, although is not identified as such.
(Inverted) Covariance matrix [cov_mat/cov_mat_inv]
The (inverted) matrix of prior covariances between all model realization samples and itself.
Covariance vector [cov_vec]
The vector of prior covariances between all model realization samples and a given parameter set.
Data error [data_err]
The \(1\sigma\)-confidence interval of a model comparison data point, often a measured/calculated observational error.
Data (point) identifier [data_idx]
The unique identifier of a model comparison data point, often a sequence of integers, floats and strings that describe the operations required to extract it.
Data point [data_point]
A collection of all the details (value, error, space and identifier) about a specific model comparison data point that is used to constrain the model with.
Data (value) space [data_spc]
The value space (linear, logarithmic or exponential) in which a model comparison data point is defined.
Data value [data_val]
The value of a model comparison data point, often an observed/measured value.
Emulation method [method]
The specific method (Gaussian, regression or both) that needs to be used to construct an emulator.
Emulator [emul/emulator]
The collection of all emulator systems together, provided by an Emulator object.
Emulator evaluation samples [eval_sam_set]
The sample set (to be) used for evaluating the emulator.
(Emulator) Iteration [emul_i]
A single, specified step in the construction of the emulator.
Emulator system [emul_s]
The emulated version of a single model output/comparison data point in a single iteration.
Emulator type [emul_type]
The type of emulator that needs to be constructed. This is used to make sure different emulator types are not mixed together by accident.
Evaluation/Evaluate
The process of calculating the adjusted values of a parameter set in all emulator systems starting at the first iteration, determining the corresponding implausibility values and performing an implausibility check. This process is repeated in the next iteration if the check was successful and the requested iteration has not been reached.
Evaluation set
Same as sample set.
External model realization set [ext_real_set]
A set of externally calculated and provided model realization samples and their outputs.
Frozen (active) parameters
The set of model parameters that, once considered active, will always stay active if possible.
FSLR
Abbreviation of forward stepwise linear regression.
Gaussian correlation length [l_corr]
The maximum distance between two values of a specific model parameter within which the Gaussian contribution to the correlation between the values is still significant.
Gaussian sigma [sigma]
The standard deviation of the Gaussian function. It is not required if regression is used.
HDF5
Abbreviation of Hierarchical Data Format version 5.
Implausibility (cut-off) check [impl_check]
The process of determining whether or not a given set of implausibility values satisfy the implausibility cut-offs of a specific emulator iteration.
Implausibility cut-offs [impl_cut]
The maximum implausibility values an evaluated parameter set is allowed to generate, to be considered plausible in a specific emulator iteration.
(Univariate) Implausibility value [uni_impl_val]
The number of sigmas an emulator system expects the (real) model output corresponding to a given parameter set, to be away from the data point it is compared against, given its adjusted values. It takes into account all variances associated with the parameter set, which are the observational variance (given by data_err), adjusted emulator variance (adj_var) and the model discrepancy variance (md_var).
Implausibility wildcard
A maximum implausibility value, preceding the implausibility cut-offs, that is not taken into account during the implausibility cut-off check. It is denoted as \(0\) in provided implausibility cut-off lists.
LHD
Abbreviation of Latin-Hypercube design.
Master (HDF5) file [hdf5_file]
(Path to) The HDF5-file in which all important data about the currently loaded emulator is stored. A master file is usually accompanied by several emulator system (HDF5) files, which store emulator system specific data and are externally linked to the master file.
MCMC
Abbreviation of Markov chain Monte Carlo.
Mock data
The set of comparison data points that has been generated by evaluating the model for a random parameter set and perturbing the output by the model discrepancy variance.
Model

A black box that takes a parameter set, performs a sequence of operations and returns a unique collection of values corresponding to the provided parameter set.

Note

This is how PRISM ‘sees’ a model, not the used definition of one.

2D model
A model that has/takes 2 model parameters.
2+D/nD model
A model that has/takes more than 2 model parameters.
ModelLink (subclass) [modellink/modellink_obj]
The user-provided wrapper around the model that needs to be emulated, provided by a ModelLink object.
Model discrepancy variance [md_var]
A user-defined value that includes all contributions to the overall variance on a model output that is created/caused by the model itself. More information on this can be found in Model discrepancy variance (md_var).
Model evaluation samples [add_sam_set]
The sample set (to be) used for evaluating the model.
Model output(s) [mod_out/mod_set]
The model output(s) corresponding to a single (set of) model realization/evaluation sample(s).
Model realization samples
Same as model evaluation samples.
Model realizations (set) [mod_real_set]
The combination of model realization/evaluation samples and their corresponding model outputs.
MPI
Abbreviation of Message Passing Interface.
MPI rank
An MPI process that is used by any PRISM operation, either being a controller or a worker.
MSE
Abbreviation of mean squared error.
OLS
Abbreviation of ordinary least-squares.
Parameter set [par_set]
A single combination/set of model parameter values, used to evaluate the emulator/model once.
Passive parameters
The set of model parameters that are not considered active, and therefore are considered to not have a significant influence on the output of the model.
(PRISM) Pipeline [pipe/pipeline]
The main PRISM framework that orchestrates all operations, provided by a Pipeline object.
Plausible region
The region of model parameter space that still contains plausible samples.
Plausible samples [impl_sam]
A subset of a set of emulator evaluation samples that satisfied the implausibility checks.
Polynomial order [poly_order]
Up to which order polynomial terms need to be taken into account for all regression processes.
Potentially active parameters [pot_active_par]
A user-provided set of model parameters that are allowed to become active. Any model parameter that is not potentially active will never become active, even if it should.
PRISM
The acronym for Probabilistic Regression Instrument for Simulating Models. It is also a one-word description of what PRISM does (splitting up a model into individually emulated model outputs).
PRISM (parameters) file [prism_file]
(Path to) The text-file that contains non-default values for the PRISM parameters that need to be used for the currently used Pipeline instance. It is None if no such file is used.
Prior covariance [prior_cov]
The covariance value between two parameter sets as determined by an emulator system.
Prior expectation [prior_exp]
The expectation value of a parameter set as determined by an emulator system, without taking the adjustment term (from the BLA) into account. It is a measure of how much information is captured by an emulator system. It is zero if regression is not used, as no information is captured.
Prior variance [prior_var]
The variance value of a parameter set as determined by an emulator system, without taking the adjustment term (from the BLA) into account.
Projection/Project
The process of analyzing a specific set of active parameters in an iteration to determine the correlation between the parameters.
Projection figure
The visual representation of a projection.
Regression
The process of determining the important polynomial terms of the active parameters and their coefficients, by using an FSLR algorithm.
Regression covariance(s) [poly_coef_cov]
The covariances between all polynomial coefficients of the regression function. By default, they are not calculated and it is empty if regression is not used.
Residual variance [rsdl_var]
The variance that has not been captured during the regression process. It is empty if regression is not used.
Root directory [root_dir]
(Path to) The directory/folder on the current machine in which all PRISM working directories are located. It also acts as the base for all relative paths.
Sample [sam]
Same as a parameter set.
Sample set [sam_set]
A set of samples.
Worker (rank)
An MPI process that receives its calls/orders from a controller and performs the heavy-duty operations in PRISM.
Working directory [working_dir]
(Path to) The directory/folder on the current machine in which the PRISM master file, log-file and all projection figures of the currently loaded emulator are stored.
Worker mode [worker_mode]
A mode initialized by worker_mode, where all workers are continuously listening for calls made by the controller rank and execute the received messages. This allows for serial codes to be combined more easily with PRISM. See Dual nature (normal/worker mode) for more information.

PRISM parameters

Below are descriptions of all the parameters that can be provided to PRISM in a text-file when initializing the Pipeline class (using the prism_file input argument):

n_sam_init (Default: 500)
Number of model evaluation samples that is used to construct the first iteration of the emulator. This value must be a positive integer.
proj_res (Default: 25)
Number of emulator evaluation samples that is used to generate the grid for the projection figures (it defines the resolution of the projection). This value must be a positive integer.
proj_depth (Default: 250)
Number of emulator evaluation samples that is used to generate the samples in every projected grid point (it defines the accuracy/depth of the projection). This value must be a positive integer.
base_eval_sam (Default: 800)
Base number of emulator evaluation samples that is used to analyze an iteration of the emulator. It is multiplied by the iteration number and the number of model parameters to generate the true number of emulator evaluations, in order to ensure an increase in emulator accuracy. This value must be a positive integer.
sigma (Default: 0.8)
The Gaussian sigma/standard deviation that is used when determining the Gaussian contribution to the overall emulator variance. This value is only required when method == ‘gaussian’, as the Gaussian sigma is obtained from the residual variance left after the regression optimization if regression is included. This value must be non-zero.
l_corr (Default: 0.3)
The amplitude(s) of the Gaussian correlation length. This number is multiplied by the difference between the upper and lower value boundaries of the model parameters to obtain the Gaussian correlation length for every model parameter. This value must be positive and either a scalar or a list of n_par scalars (where the values correspond to the sorted list of model parameters).
impl_cut (Default: [0.0, 4.0, 3.8, 3.5])
A list of implausibility cut-off values that specifies the maximum implausibility values a parameter set is allowed to have to be considered ‘plausible’. A zero can be used as a filler value, either taking on the preceding value or acting as a wildcard if the preceding value is a wildcard or non-existent. Zeros are appended at the end of the list if the length is less than the number of comparison data points, while extra values are ignored if the length is more. This must be a sorted list of positive values (excluding zeros).
criterion (Default: ‘multi’)
The criterion to use for determining the quality of the LHDs that are used, represented by an integer, float, string or None. This parameter is the only non-PRISM parameter. Instead, it is used in the lhd()-function of the e13Tools package. By default, ‘multi’ is used to give equal priority to maximizing minimum distances and minimizing the maximum correlation between pair-wise samples.
method (Default: ‘full’)

The method to use for constructing the emulator. ‘gaussian’ will only include Gaussian processes (no regression), which is much faster, but also less accurate. ‘regression’ will only include regression processes (no Gaussian), which is more accurate than Gaussian only, but underestimates the emulator variance by multiple orders of magnitude. ‘full’ includes both Gaussian and regression processes, which is slower than Gaussian only, but by far the most accurate both in terms of expectation and variance values.

‘gaussian’ can be used for faster exploration especially for simple models. ‘regression’ should only be used when the polynomial representation of a model is important and enough model realizations are available. ‘full’ should be used by default.

Warning

When using PRISM on a model that can be described completely by the regression function (anything that has an analytical, polynomial form up to order poly_order like a straight line or a quadratic function), use the ‘gaussian’ method unless unavoidable (in which case n_sam_init and base_eval_sam must be set to very low values).

When using the regression method on such a model, PRISM will be able to capture the behavior of the model perfectly given enough samples, in which case the residual (unexplained) variance will be approximately zero and therefore sigma as well. This can occassionally cause floating point errors when calculating emulator variances, which in turn can give unexplainable artifacts in the adjustment terms, therefore causing answers to be incorrect.

Since PRISM’s purpose is to identify the characteristics of a model and therefore it does not know anything about its workings, it is not possible to automatically detect such problems.

use_regr_cov (Default: False)
Whether or not the regression variance should be taken into account for the variance calculations. The regression variance is the variance on the regression process itself and is only significant if a low number of model realizations (n_sam_init and base_eval_sam) is used to construct the emulator systems. Including it usually only has a very small effect on the overall variance value, while it can slow down the emulator evaluation rate by as much as a factor of 3. This value is not required if method == ‘gaussian’ and is automatically set to True if method == ‘regression’. This value must be a bool.
poly_order (Default: 3)
Up to which order all polynomial terms of all model parameters should be included in the active parameters and regression processes. This value is not required if method == ‘gaussian’ and do_active_anal is False. This value must be a positive integer.
n_cross_val (Default: 5)
Number of (k-fold) cross-validations that must be used for determining the quality of the active parameters analysis and regression process fits. If this parameter is zero, cross-validations are not used. This value is not required if method == ‘gaussian’ and do_active_anal is False. This value must be a non-negative integer and not equal to 1.
do_active_anal (Default: True)
Whether or not an active parameters analysis must be carried out for every iteration of every emulator system. If False, all potentially active parameters listed in pot_active_par will be active. This value must be a bool.
freeze_active_par (Default: True)
Whether or not active parameters should be frozen in their active state. If True, parameters that have been considered active in a previous iteration of an emulator system, will automatically be active again (and skip any active parameters analysis). This value must be a bool.
pot_active_par (Default: None)
A list of parameter names that indicate which parameters are potentially active. Potentially active parameters are the only parameters that will enter the active parameters analysis (or will all be active if do_active_anal is False). Therefore, all parameters not listed will never be considered active. If all parameters should be potentially active, then a None can be given. This must either be a list of parameter names or None.
use_mock (Default: False)
Whether or not mock data must be used as comparison data when constructing a new emulator. Mock data is calculated by evaluating the model for a randomly chosen set of parameter values, and adding the model discrepancy variances as noise to the returned data values. When using mock data for an emulator, it is not possible to change the comparison data in later emulator iterations. This value must be a bool.

HDF5

Whenever PRISM constructs an emulator, it automatically stores all the calculated data for it in an HDF5-file named 'prism.hdf5' in the designated working directory. This file contains all the data that is required in order to recreate all emulator systems that have been constructed for the emulator belonging to this run. If the Pipeline class is initialized by using an HDF5-file made by PRISM, it will load in this data and return a Pipeline object in the same state as described in the file.

Below is a short overview of all the data that can be found inside a PRISM master HDF5-file. HDF5-files can be viewed freely by the user using the HDFView application made available by The HDFGroup.


The general file contains:
  • Attributes (9/10): Describe the general non-changeable properties of the emulator, which include:

    • Emulator type and method;
    • Gaussian parameters;
    • Name of used ModelLink subclass;
    • Used PRISM version;
    • Polynomial order;
    • Bools for using mock data or regression covariance;
    • Mock data parameters if mock data was used.
  • Every emulator iteration has its own data group with the iteration number as its name. This data group stores all data/information specific to that iteration.


An iteration data group ('i') contains:
  • Attributes (9): Describe the general properties and results of this iteration, including:

    • Active parameters for this emulator iteration;
    • Implausibility cut-off parameters;
    • Number of emulated data points, emulator systems, emulator evaluation samples, plausible samples and model realization samples;
    • Bool stating whether this emulator iteration used an external model realization set.
  • 'emul_n': The data group that contains all data for a specific emulator system in this iteration. The value of 'n' indicates which emulator system it is, not the data point. See below for its contents;

  • 'impl_sam': The set of emulator evaluation samples that survived the implausibility checks and will be used to construct the next iteration;

  • 'proj_hcube': The data group that contains all data for the (created) projections for this iteration, if at least one has been made. See below for its contents;

  • 'sam_set': The set of model realization samples that were used to construct this iteration. In every iteration after the first, this is the 'impl_sam' of the previous iteration;

  • 'statistics': An empty data set that stores several different types of statistics as its attributes, including:

    • Size of the MPI communicator during various construction steps;
    • Average evaluation rate/time of the emulator and model;
    • Total time cost of most construction steps (note that this value may be incorrect if a construction was interrupted);
    • Percentage of parameter space that is still plausible within the iteration.

An emulator system data group ('i/emul_n') contains:
  • Attributes (5+): List the details about the model comparison data point used in this emulator system, including:

    • Active parameters for this emulator system;
    • Data errors, identifiers, value space and value;
    • Regression score and residual variance if regression was used.
  • 'cov_mat': The pre-calculated covariance matrix of all model evaluation samples in this emulator system. This data set is never used in PRISM and stored solely for user-convenience;

  • 'cov_mat_inv': The pre-calculated inverse of 'cov_mat';

  • 'exp_dot_term': The pre-calculated second expectation adjustment dot-term (\(\mathrm{Var}\left(D\right)^{-1}\cdot\left(D-\mathrm{E}(D)\right)\)) of all model evaluation samples in this emulator system.

  • 'mod_set': The model outputs for the data point in this emulator system corresponding to the 'sam_set' used in this iteration;

  • 'poly_coef' (if regression is used): The non-zero coefficients for the polynomial terms in the regression function in this emulator system;

  • 'poly_coef_cov' (if regression and regr_cov are used): The covariances for all polynomial coefficients 'poly_coef';

  • 'poly_idx' (if regression is used): The indices of the polynomial terms with non-zero coefficients if all active parameters are converted to polynomial terms;

  • 'poly_powers' (if regression is used): The powers of the polynomial terms corresponding to 'poly_idx'. Both 'poly_idx' and 'poly_powers' are required since different methods of calculating the polynomial terms are used depending on the number of required terms and samples;

  • 'prior_exp_sam_set': The pre-calculated prior expectation values of all model evaluation samples in this emulator system. This data set is also never used in PRISM.


A projections data group ('i/proj_hcube') contains individual projection data groups ('i/proj_hcube/<name>'), which contain:
  • Attributes (4): List the general properties with which this projection was made, including:

    • Implausibility cut-off parameters (they can differ from the iteration itself);
    • Projection depth and resolution.
  • 'impl_los': The calculated line-of-sight depth for all grid points in this projection;

  • 'impl_min': The calculated minimum implausibility values for all grid points in this projection.

FAQ

How do I contribute?

Contributing to PRISM is done through pull requests in the repository. If you have an idea on what to contribute, it is recommended to open a GitHub issue about it, such that the maintainers can give their advice or help. If you want to contribute but do not really know what, then you can take a look at the ToDos that are scattered throughout the code. When making a contribution, please keep in mind that it must be compatible with all Python versions that PRISM supports (2.7/3.5+), and preferably with all operating systems as well.

How do I report a bug/problem?

By opening a GitHub issue and using the Bug report template.

What does the term mean?

A list of the most commonly used terms in PRISM can be found on the Terminology page.

What OS are supported?

PRISM should be compatible with all Windows, Mac OS and UNIX-based machines, as long as they support one of the required Python versions. Compatibility is currently tested for Linux 16.04 (Xenial)/Mac OS-X using Travis CI and Windows 32-bit/64-bit using AppVeyor.

Community guidelines

PRISM is an open-source and free-to-use software package (and it always will be), provided under the BSD-3 license (see below for the full license).

Users are highly encouraged to make contributions to the package or request new features by opening a GitHub issue. If you would like to contribute to the package, but do not know what, then there are quite a few ToDos in the code that may give you some inspiration. As with contributions, if you find a problem or issue with PRISM, please do not hesitate to open a GitHub issue about it.

And, finally, if you use PRISM as part of your workflow in a scientific publication, please consider citing the PRISM pipeline paper (and starring the repository):

Note

Put reference to PRISM paper here!

License

BSD 3-Clause License

Copyright (c) 2019, Ellert van der Velden
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

* Redistributions of source code must retain the above copyright notice, this
  list of conditions and the following disclaimer.

* Redistributions in binary form must reproduce the above copyright notice,
  this list of conditions and the following disclaimer in the documentation
  and/or other materials provided with the distribution.

* Neither the name of the copyright holder nor the names of its
  contributors may be used to endorse or promote products derived from
  this software without specific prior written permission.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

Additions

Below are some bigger ideas/improvements that may be added to PRISM if there is demand:

  • Low-level MPI implementation (probably by using D2O);

    With 6 emulator systems and 4 processes, the three different MPI levels would be:

    • No level: 6-0-0-0;
    • High-level: 2-2-1-1;
    • Low-level: 1.5-1.5-1.5-1.5.
  • Dynamic implausibility cut-offs;

  • Allow for a master projection figure to be made (kind of like a double corner plot);

  • Set PRISM parameters properties manually (using a property setter) instead of only reading it in from file. This will additionally allow for external functions to be executed that determine the values of certain properties;

  • Allow for user-provided methods in the ModelLink subclass to be executed at specific points in the emulator construction;

  • Implement multi-variate implausibilities;

  • Allow for no ModelLink object to be provided, which blocks construction but enables everything emulator-only related;

  • Allow for old PRISM master files to be provided when making a new emulator, recycling work done previously;

  • If MPI_call is False for the ModelLink subclass, use all MPI ranks to evaluate a part of sam_set in the model simultaneously. This will require a check or flag that the model can be called in multiple instances simultaneously (to accommodate for models that, for example, need to read files during evaluations). Added benefit of this is that it would become possible to add the option for the user to set a preferred number of MPI processes calling the model (in MPI), allowing PRISM to split up the available processes if more efficient;

  • GPU acceleration;

  • Adding the theory behind PRISM to the docs.

Pipeline

class prism.Pipeline(modellink_obj, root_dir=None, working_dir=None, prefix=None, prism_file=None, emul_type=None, comm=None)[source]

Defines the Pipeline class of the PRISM package.

The Pipeline class is the main user class of the PRISM package and provides a user-friendly environment that gives access to all operations within the package.

_Projection__analyze_proj_hcube(hcube)

Analyzes an emulator projection hypercube hcube.

Parameters:hcube (1D array_like of int of length {1, 2}) – Array containing the internal integer identifiers of the main model parameters that require a projection hypercube.
Returns:
  • impl_min_hcube (1D ndarray object) – List containing the lowest implausibility value that can be reached in every single grid point on the given hypercube.
  • impl_los_hcube (1D ndarray object) – List containing the fraction of the total amount of evaluated samples in every single grid point on the given hypercube, that still satisfied the implausibility cut-off criterion.
_Projection__draw_2D_proj_fig(hcube, impl_min, impl_los, proj_res)

Draws the 2D projection figure for the provided hcube, given the impl_min and impl_los values.

Parameters:
  • hcube (1D array_like of int of length 1) – Array containing the internal integer identifiers of the main model parameters that require a projection figure.
  • impl_min (1D ndarray object) – List containing the lowest implausibility value that can be reached in every single grid point on the given hypercube.
  • impl_los (1D ndarray object) – List containing the fraction of the total amount of evaluated samples in every single grid point on the given hypercube, that still satisfied the implausibility cut-off criterion.
  • proj_res (int) – Number of emulator evaluations used to generate the grid for the given hypercube.
_Projection__draw_3D_proj_fig(hcube, impl_min, impl_los, proj_res)

Draws the 3D projection figure for the provided hcube, given the impl_min and impl_los values.

Parameters:
  • hcube (1D array_like of int of length 2) – Array containing the internal integer identifiers of the main model parameters that require a projection figure.
  • impl_min (1D ndarray object) – List containing the lowest implausibility value that can be reached in every single grid point on the given hypercube.
  • impl_los (1D ndarray object) – List containing the fraction of the total amount of evaluated samples in every single grid point on the given hypercube, that still satisfied the implausibility cut-off criterion.
  • proj_res (int) – Number of emulator evaluations used to generate the grid for the given hypercube.
_Projection__get_default_input_arguments()

Generates a dict containing default values for all input arguments.

Returns:kwargs_dict (dict) – Dict containing all default input argument values.
_Projection__get_default_parameters()

Generates a dict containing default values for all projection parameters.

Returns:par_dict (dict) – Dict containing all default projection parameter values.
_Projection__get_fig_path(hcube, emul_i=None)

Determines the absolute path of a projection figure corresponding to a provided projection hypercube hcube in emulator iteration emul_i and returns it.

Parameters:hcube (1D array_like of int of length {1, 2} or str) – Array containing the internal integer identifiers of the main model parameters that require a projection hypercube. If str, the name of hcube instead (_Projection__get_hcube_name()).
Other Parameters:
 emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.
Returns:
  • fig_path (str) – The absolute path to the requested projection figure.
  • fig_path_s (str) – The absolute path to the smoothed version.
_Projection__get_hcube_name(hcube)

Determines the name of a provided projection hypercube hcube and returns it.

Parameters:hcube (1D array_like of int of length {1, 2}) – Array containing the internal integer identifiers of the main model parameters that require a projection hypercube.
Returns:hcube_name (str) – The name of this projection hypercube.
_Projection__get_proj_data(hcube)

Returns the projection data belonging to the provided hypercube hcube.

Parameters:hcube (1D array_like of int of length {1, 2}) – Array containing the internal integer identifiers of the main model parameters that require a projection hypercube.
Returns:
  • impl_min_hcube (1D ndarray object) – List containing the lowest implausibility value that can be reached in every single grid point on the given hypercube.
  • impl_los_hcube (1D ndarray object) – List containing the fraction of the total amount of evaluated samples in every single grid point on the given hypercube, that still satisfied the implausibility cut-off criterion.
  • proj_res (int) – Number of emulator evaluations used to generate the grid for the given hypercube.
  • proj_depth (int) – Number of emulator evaluations used to generate the samples in every grid point for the given hypercube.
_Projection__get_proj_hcube(hcube)

Generates a projection hypercube hcube containing emulator evaluation samples The output of this function depends on the requested projection type.

Parameters:hcube (1D array_like of int of length {1, 2}) – Array containing the internal integer identifiers of the main model parameters that require a projection hypercube.
Returns:proj_hcube (3D ndarray object) – 3D projection hypercube of emulator evaluation samples. For 3D projections, the grid uses matrix indexing (second parameter varies the fastest).
_Projection__get_req_hcubes(proj_par)

Determines which projection hypercubes have been requested by the user. Also checks if these projection hypercubes have been calculated before, and depending on the value of force, either skips them or recreates them.

Parameters:proj_par (1D array_like of {int, str} or None) – For which model parameters to construct the projection figures. If 1D array_like, construct projection figures for all combinations of provided model parameters that are active, with a string referring to the name of the model parameter and an integer referring to the position in which the model parameter is shown in the details() method. If None, projection figures are made for all active model parameters.

Generates

hcubes : list of lists
List containing the parameter indices of the requested projection hypercubes.
create_hcubes : list of lists
List containing the parameter indices of the requested projection hypercubes that need to be created first.
_Projection__prepare_projections(emul_i, proj_par, **kwargs)

Prepares the pipeline for the creation of the requested projections.

Parameters:
  • emul_i (int or None) – Number indicating the requested emulator iteration.
  • proj_par (1D array_like of {int, str} or None) – For which model parameters to construct the projection figures. If 1D array_like, construct projection figures for all combinations of provided model parameters that are active, with a string referring to the name of the model parameter and an integer referring to the position in which the model parameter is shown in the details() method. If None, projection figures are made for all active model parameters.
  • kwargs (dict) – Dict of keyword arguments that were provided to project().
_Projection__process_input_arguments(emul_i, **kwargs)

Processes the input arguments given to the project() method.

Parameters:
  • emul_i (int or None) – Number indicating the requested emulator iteration.
  • kwargs (dict) – Dict of keyword arguments that were provided to project().

Generates

All default and provided arguments are saved to their respective properties.

_Projection__read_parameters()

Reads in the Projection parameters from the provided PRISM parameter file and saves them in the current Projection instance.

_Projection__save_data(data_dict)

Saves a given data dict {keyword: data} at the emulator iteration this class was initialized for, to the HDF5-file.

Parameters:

data_dict (dict) – Dict containing the data that needs to be saved to the HDF5-file.

Keyword Arguments:
 
  • keyword ({‘nD_proj_hcube’}) – String specifying the type of data that needs to be saved.
  • data ({int, float, str, array_like} or dict) – The actual data that needs to be saved at data keyword keyword. If dict, save every item individually.

Generates

The specified data is saved to the HDF5-file.

__init__(modellink_obj, root_dir=None, working_dir=None, prefix=None, prism_file=None, emul_type=None, comm=None)[source]

Initialize an instance of the Pipeline class.

Parameters:

modellink_obj (ModelLink object) – Instance of the ModelLink subclass that links the emulated model to this Pipeline instance.

Other Parameters:
 
  • root_dir (str or None. Default: None) – String containing the absolute path of the root directory where all working directories are stored. If None, root directory will be set to the directory this class was initialized in.
  • working_dir (str, int or None. Default: None) – String containing the name of the working directory of the emulator in root_dir. If int, a new working directory will be created in root_dir. If None, working directory is set to the last one that was created in root_dir that starts with the given prefix. If no directories are found, one will be created.
  • prefix (str or None. Default: None) – String containing a prefix that is used for naming new working directories or scan for existing ones. If None, all directories in root_dir are considered working directories and ‘prism_’ will be used as a prefix for new ones.
  • prism_file (str or None. Default: None) – String containing the absolute or relative path to the TXT-file containing the PRISM parameters that need to be changed from their default values. If a relative path is given, its path must be relative to root_dir or the current directory. If None, no changes will be made to the default parameters.
  • emul_type (Emulator subclass or None. Default: None) – The type of Emulator to use in this Pipeline instance. If None, use the default emulator instead.
  • comm (Intracomm object or None. Default: None) – The MPI intra-communicator to use in this Pipeline instance or MPI.COMM_WORLD if comm is None. If mpi4py is not installed, _dummyMPI is used instead.
_call_model(emul_i, par_set, data_idx)[source]

Obtain the output corresponding to the provided data_idx that is generated by the model for a given model parameter sample par_set. The current emulator iteration emul_i is also provided in case it is required by the ModelLink subclass.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • par_set (1D array_like) – Model parameter sample to calculate the model output for.
  • data_idx (list of tuples) – The list of data identifiers for which the model is requested to return the corresponding data values.
Returns:

mod_out (1D ndarray object) – Model output corresponding to given par_set.

_compile_code_snippets()[source]

Compiles all pre-defined built-in code snippets to code objects and saves them to code_objects. These code objects are used for performing standard operations in the _evaluate_sam_set() method.

_do_impl_check(emul_i, uni_impl_val)[source]

Performs an implausibility cut-off check on the provided implausibility values uni_impl_val at emulator iteration emul_i.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • uni_impl_val (1D array_like) – Array containing all univariate implausibility values corresponding to a certain parameter set for all data points.
Returns:

  • result (bool) – 1 if check was successful, 0 if it was not.
  • impl_cut_val (float) – Implausibility value at the first real implausibility cut-off.

_evaluate_model(emul_i, sam_set, data_idx)[source]

Evaluates the model for provided evaluation sample set sam_set at given data points data_idx.

This method automatically distributes the samples according to the various flags set in the ModelLink subclass.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • sam_set (1D or 2D array_like) – Parameter/sample set to evaluate in the model.
  • data_idx (list of tuples) – The list of data identifiers for which the model is requested to return the corresponding data values.
Returns:

mod_set (2D ndarray object of shape (n_sam, n_data)) – Array containing the data values corresponding to the requested data points.

_evaluate_sam_set(emul_i, sam_set, exec_code)[source]

Evaluates a provided set of emulator evaluation samples sam_set at a given emulator iteration emul_i. The provided tuple of code snippets exec_code are executed using Python’s exec() function at specific points during the analysis.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • sam_set (2D ndarray object) – Array containing model parameter value sets to be evaluated in all emulator systems in emulator iteration emul_i.
  • exec_code ({‘analyze’, ‘evaluate’, ‘hybrid’, ‘project’} or tuple) – Tuple of five code snippets (pre_code, eval_code, anal_code, post_code, exit_code) to be executed at specific points during the analysis. If string, use one of the built-in tuples in code_objects instead.
Other Parameters:
 
  • pre_code (str or code object) – Code snippet to be executed before the evaluation of sam_set starts.
  • eval_code (str or code object) – Code snippet to be executed after the evaluation of each sample in sam_set.
  • anal_code (str or code object) – Code snippet to be executed after the analysis of each sample in sam_set. This code snippet is only executed by the controller.
  • post_code (str or code object) – Code snippet to be executed after the evaluation of sam_set ends.
  • exit_code (str or code object) – Code snippet to be executed before returning the results of the evaluation of sam_set. This code snippet is only executed by the controller.
Returns:

results (object) – The object that is assigned to results, which is defaulted to None if no code snippet changes it. Preferably, the execution of post_code and/or exit_code modifies results. All MPI ranks return it.

Notes

If any of the code snippets is provided as a string, it will be compiled into a code object before starting the evaluation.

_get_default_parameters()[source]

Generates a dict containing default values for all pipeline parameters.

Returns:par_dict (dict) – Dict containing all default pipeline parameter values.
_get_eval_sam_set(emul_i)[source]

Generates an emulator evaluation sample set to be used for analyzing an emulator iteration. Currently uses the lhd() function.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
Returns:eval_sam_set (2D ndarray object) – Array containing the evaluation samples.
_get_ext_real_set(ext_real_set)[source]

Processes an externally provided model realization set ext_real_set, containing the used sample set and the corresponding data value set.

Parameters:ext_real_set (list, dict or None) – List of arrays containing an externally calculated set of model evaluation samples and its data values, a dict with keys ['sam_set', 'mod_set'] containing these arrays or None if no external set needs to be used.
Returns:
  • ext_sam_set (1D or 2D ndarray object) – Array containing the externally provided model evaluation samples.
  • ext_mod_set (1D or 2D ndarray object) – Array containing the model outputs of all specified externally provided model evaluation samples.
_get_impl_cut(impl_cut, temp)[source]

Generates the full list of impl_cut-offs from the incomplete, shortened impl_cut list.

Parameters:
  • impl_cut (1D list) – Incomplete, shortened impl_cut-offs list provided during class initialization.
  • temp (bool) – Whether the implausibility parameters should only be stored in memory (True) or should also be saved to HDF5 (False).

Generates

impl_cut : 1D ndarray object
Full list containing the impl_cut-offs for all data points provided to the emulator.
cut_idx : int
Index of the first impl_cut-off in impl_cut that is not a wildcard.
_get_impl_par(temp)[source]

Reads in the impl_cut list and other parameters for implausibility evaluations from the PRISM parameters file and saves them in the last emulator iteration.

Parameters:temp (bool) – Whether the implausibility parameters should only be stored in memory (True) or should also be saved to HDF5 (False).

Generates

impl_cut : 1D ndarray object
Full list containing the impl_cut-offs for all data points provided to the emulator.
cut_idx : int
Index of the first impl_cut-off in the impl_cut list that is not 0.
_get_iteration_data(emul_i, sam_set, ext_sam_set, ext_mod_set)[source]

Obtains the model realization data for given emulator iteration emul_i by evaluating the provided sam_set in the model and distributing model outputs to the correct emulator systems.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • sam_set (2D ndarray object) – Array containing the model evaluation samples.
  • ext_sam_set (1D or 2D ndarray object) – Array containing the externally provided model evaluation samples.
  • ext_mod_set (1D or 2D ndarray object) – Array containing the model outputs of all specified externally provided model evaluation samples.

Generates

sam_set : 2D ndarray object
Array containing the model evaluation samples for emulator iteration emul_i.
mod_set : 2D ndarray object
Array containing the model outputs of all specified model evaluation samples for emulator iteration emul_i.
_get_md_var(emul_i, par_set)[source]

Retrieves the model discrepancy variances, which includes all variances that are created by the model provided by the ModelLink instance. This method tries to call the get_md_var() method, and assumes a default model discrepancy variance of 1/6th the data value if it cannot be called. If the data value space is not linear, then this default value is calculated such to reflect that.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • par_set (1D ndarray object) – Model parameter value set to calculate the model discrepancy variances for.
Returns:

var_md (2D ndarray object) – Variance of the model discrepancy.

_get_mock_data()[source]

Generates mock data and loads it into the ModelLink object that was provided during class initialization. This function overwrites the ModelLink properties holding the parameter estimates, data values and data errors.

Generates

Overwrites the corresponding ModelLink class properties with the generated values.

_get_n_eval_sam(emul_i)[source]

This function calculates the total number of emulator evaluation samples at a given emulator iteration emul_i from base_eval_sam.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
Returns:n_eval_sam (int) – Total number of emulator evaluation samples.
_get_paths(root_dir, working_dir, prefix, prism_file)[source]

Obtains the path for the root directory, working directory and parameters file for PRISM.

Parameters:
  • root_dir (str or None) – String containing the absolute path of the root directory where all working directories are stored. If None, root directory will be set to the directory this class was initialized in.
  • working_dir (str, int or None) – String containing the name of the working directory of the emulator in root_dir. If int, a new working directory will be created in root_dir. If None, working directory is set to the last one that was created in root_dir that starts with the given prefix. If no directories are found, one will be created.
  • prefix (str or None) – String containing a prefix that is used for naming new working directories or scan for existing ones. If None, all directories in root_dir are considered working directories and ‘prism_’ will be used as a prefix for new ones.
  • prism_file (str or None) – String containing the absolute or relative path to the TXT-file containing the PRISM parameters that need to be changed from their default values. If a relative path is given, its path must be relative to root_dir or the current directory. If None, no changes will be made to the default parameters.

Generates

The absolute paths to the root directory, working directory, emulator master HDF5-file and PRISM parameters file.

_get_uni_impl(emul_i, emul_s_seq, par_set, adj_exp_val, adj_var_val)[source]

Calculates the univariate implausibility values at a given emulator iteration emul_i for specified expectation and variance values adj_exp_val and adj_var_val, corresponding to given par_set.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set (1D ndarray object) – Model parameter value set to calculate the univariate implausibility values for. Only used to pass to the get_md_var() method.
  • adj_exp_val, adj_var_val (1D array_like) – The adjusted expectation and variance values to calculate the univariate implausibility for.
Returns:

uni_impl_val (1D ndarray object) – Univariate implausibility value for all requested emulator systems.

_listen_for_calls()[source]

All worker ranks in the comm communicator start listening for calls from the corresponding controller rank and will attempt to execute the received message. Listening for calls continues until _worker_mode is set to False.

This method is automatically initialized and finalized when using the worker_mode.

_load_data()[source]

Loads in all the important pipeline data into memory for the controller rank. If it is detected that the last emulator iteration has not been analyzed yet, the implausibility analysis parameters are read in from the PRISM parameters file and temporarily stored in memory.

Generates

All relevant pipeline data up to the last emulator iteration is loaded into memory.

_make_call(exec_fn, *args, **kwargs)[source]

Send the provided exec_fn to all worker ranks, if they are listening for calls, and tell them to execute it using the provided args and kwargs. All ranks that call this function will execute exec_fn as well.

If used within the worker_mode context manager, this function should only be called by the controller. If not, it should be called by all ranks that must execute exec_fn.

Parameters:
  • exec_fn (str, callable or None) – If string, a callable attribute of this Pipeline instance or a callable object that the workers should execute if not. If None, the workers stop listening for calls instead (disables worker mode).
  • args (tuple) – Positional arguments that need to be provided to exec_fn.
  • kwargs (dict) – Keyword arguments that need to be provided to exec_fn.
_multi_call_model(emul_i, sam_set, data_idx)[source]

Obtain the output set corresponding to the provided data_idx that is generated by the model for a given model parameter sample set sam_set. The current emulator iteration emul_i is also provided in case it is required by the ModelLink subclass.

This is a multi-version of _call_model().

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • sam_set (2D array_like) – Model parameter sample set to calculate the model output for.
  • data_idx (list of tuples) – The list of data identifiers for which the model is requested to return the corresponding data values.
Returns:

mod_set (2D ndarray object) – Model output set corresponding to given sam_set.

_read_parameters()[source]

Reads in the Pipeline parameters from the provided PRISM parameter file and saves them in the current Pipeline instance.

_save_data(data_dict)[source]

Saves a given data dict {keyword: data} at the last emulator iteration to the HDF5-file and as an data attribute to the current Pipeline instance.

Parameters:

data_dict (dict) – Dict containing the data that needs to be saved to the HDF5-file.

Keyword Arguments:
 
  • keyword ({‘impl_par’, ‘impl_sam’, ‘n_eval_sam’}) – String specifying the type of data that needs to be saved.
  • data ({int, float, str, array_like} or dict) – The actual data that needs to be saved at data keyword keyword. If dict, save every item individually.

Generates

The specified data is saved to the HDF5-file.

_save_statistics(emul_i, stat_dict)[source]

Saves a given statistics dict {keyword: [value, unit]} at emulator iteration emul_i to the HDF5-file. The provided values are always saved as strings.

Parameters:

emul_i (int) – Number indicating the requested emulator iteration.

Keyword Arguments:
 
  • keyword (str) – String containing the name/keyword of the statistic that is being saved.
  • value (int, float or str) – The value of the statistic.
  • unit (str) – The unit of the statistic.
analyze()[source]

Analyzes the emulator at the last emulator iteration for a large number of emulator evaluation samples. All samples that survive the implausibility checks, are used in the construction of the next emulator iteration.

Generates

impl_sam : 2D ndarray object
Array containing all emulator evaluation samples that survived the implausibility checks.
construct(emul_i=None, analyze=True, ext_real_set=None, force=False)[source]

Constructs the emulator at the specified emulator iteration emul_i, and performs an implausibility analysis on the emulator iteration right afterward if requested (analyze()).

Other Parameters:
 
  • emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the next iteration of the emulator will be constructed.
  • analyze (bool. Default: True) – Bool indicating whether or not to perform an analysis after the specified emulator iteration has been successfully constructed, which is required for constructing the next iteration.
  • ext_real_set (list, dict or None. Default: None) – List of arrays containing an externally calculated set of model evaluation samples and its data values, a dict with keys ['sam_set', 'mod_set'] containing these arrays or None if no external set needs to be used.
  • force (bool. Default: False) – Controls what to do if the specified emulator iteration emul_i already (partly) exists. If False, finish construction of the specified iteration or skip it if already finished. If True, reconstruct the specified iteration entirely.

Generates

A new HDF5-group with the emulator iteration as its name, in the loaded emulator master file, containing emulator data required for this emulator iteration.

Notes

Using an emulator iteration that has been (partly) constructed before, will finish construction or skip it if already finished when force is False; or it will delete that and all following iterations, and reconstruct the specified iteration when force is True. Using emul_i = 1 and force is True is equivalent to reconstructing the entire emulator.

If no implausibility analysis is requested, then the implausibility parameters are read in from the PRISM parameters file and temporarily stored in memory in order to enable the usage of the evaluate() and project() methods.

details(emul_i=None)[source]

Prints the details/properties of the currently loaded Pipeline instance at given emulator iteration emul_i. See Props for detailed descriptions of all printed properties.

Other Parameters:
 emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.

Props

Working directory
The relative path to the working directory of the emulator starting at the current working directory.
Emulator type
The type of this emulator, corresponding to the emul_type of the provided emul_type during Pipeline initialization. If no emulator type was provided during initialization, this is ‘default’.
ModelLink subclass
Name of the ModelLink subclass used to construct this emulator.
Emulation method
Indicates the combination of regression and Gaussian emulation methods that have been used for this emulator.
Mock data used?
Whether or not mock data has been used to construct this emulator. If so, the printed estimates for all model parameters are the parameter values used to create the mock data.
Emulator iteration
The iteration of the emulator this details overview is about. By default, this is the last (partly) constructed iteration.
Construction completed?
Whether or not the construction of this emulator iteration is completed. If not, the missing components for each emulator system are listed and the remaining information of this iteration is not printed.
Plausible regions?
Whether or not plausible regions have been found during the analysis of this emulator iteration. If no analysis has been done yet, “N/A” will be printed.
Projections available?
Whether or not projections have been created for this emulator iteration. If projections are available and analysis has been done, but with different implausibility cut-offs, a “desync” note is added. Also prints number of available projections versus maximum number of projections in parentheses.
# of model evaluation samples
The total number of model evaluation samples used to construct all emulator iterations up to this iteration, with the number for every individual iteration in parentheses.
# of plausible/analyzed samples
The number of emulator evaluation samples that passed the implausibility check out of the total number of analyzed samples in this emulator iteration. This is the number of model evaluation samples that was/will be used for the construction of the next emulator iteration. If no analysis has been done, the numbers show up as “-“.
% of parameter space remaining
The percentage of the total number of analyzed samples that passed the implausibility check in this emulator iteration. If no analysis has been done, the number shows up as “-“.
# of active/total parameters
The number of model parameters that was considered active during the construction of this emulator iteration, compared to the total number of model parameters defined in the used ModelLink subclass.
# of emulated data points
The number of data points that have been emulated in this emulator iteration.
# of emulator systems
The total number of emulator systems that are required in this emulator. The number of active emulator systems is equal to the number of data points.
Parameter space
Lists the name, lower and upper value boundaries and estimate (if provided) of all model parameters defined in the used ModelLink subclass. An asterisk is printed in front of the parameter name if this model parameter was considered active during the construction of this emulator iteration. A question mark is used instead if the construction of this emulator iteration is not finished.
evaluate(sam_set, emul_i=None)[source]

Evaluates the given model parameter sample set sam_set up to given emulator iteration emul_i. The output of this function depends on the number of dimensions in sam_set. The output is always provided on the controller rank.

Parameters:sam_set (1D or 2D array_like or dict) – Array containing model parameter value sets to be evaluated in the emulator up to emulator iteration emul_i.
Other Parameters:
 emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.
Returns:
  • impl_check (list of bool) – List of bool indicating whether or not the given samples passed the implausibility check at the given emulator iteration emul_i.
  • emul_i_stop (list of int) – List containing the last emulator iterations at which the given samples are still within the plausible region of the emulator.
  • adj_exp_val (2D ndarray object) – Array containing the adjusted expectation values for all given samples.
  • adj_var_val (2D ndarray object) – Array containing the adjusted variance values for all given samples.
  • uni_impl_val (2D ndarray object) – Array containing the univariate implausibility values for all given samples.

Prints (if ndim(sam_set) == 1)

impl_check : bool
Bool indicating whether or not the given sample passed the implausibility check at the given emulator iteration emul_i.
emul_i_stop : int
Last emulator iteration at which the given sample is still within the plausible region of the emulator.
adj_exp_val : 1D ndarray object
The adjusted expectation values for the given sample.
adj_var_val : 1D ndarray object
The adjusted variance values for the given sample.
sigma_val : 1D ndarray object
The corresponding sigma value for the given sample.
uni_impl_val : 1D ndarray object
The univariate implausibility values for the given sample.

Notes

If given emulator iteration emul_i has been analyzed before, the implausibility parameters of the last analysis are used. If not, then the values are used that were read in when the emulator was loaded.

project(emul_i=None, proj_par=None, **kwargs)

Analyzes the emulator iteration emul_i and constructs a series of projection figures detailing the behavior of the model parameters corresponding to the given proj_par. The input and output depend on the number of model parameters n_par.

Parameters:
  • emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.
  • proj_par (1D array_like of {int, str} or None. Default: None) – For which model parameters to construct the projection figures. If 1D array_like, construct projection figures for all combinations of provided model parameters that are active, with a string referring to the name of the model parameter and an integer referring to the position in which the model parameter is shown in the details() method. If None, projection figures are made for all active model parameters.
Keyword Arguments:
 
  • proj_type ({‘2D’, ‘3D’, ‘both’}. Default: ‘2D’ (2D), ‘both’ (nD)) – String indicating which projection type to create for all supplied active parameters. If n_par == 2, this is always ‘2D’ (and cannot be modified).
  • figure (bool. Default: True) – Whether or not to create the projection figures. If True, the figures are calculated, drawn and saved. If False, the figures are calculated and their data is returned in a dict.
  • align ({‘row’/’horizontal’, ‘col’/’column’/’vertical’}. Default: ‘col’) – If figure is True, string indicating how to position the two subplots. If ‘row’/’horizontal’, the subplots are positioned on a single row. If ‘col’/’column’/’vertical’, the subplots are positioned on a single column.
  • smooth (bool. Default: False) – Controls what to do if a grid point contains no plausible samples, but does contain a minimum implausibility value below the first non-wildcard cut-off. If False, these values are kept, which can show up as artifact-like features in the projection figure. If True, these values are set to the first cut-off, removing them from the projection figure. Doing this may also remove interesting features. This does not affect the projection data saved to HDF5. Smoothed figures have an ‘_s’ string appended to their filenames.
  • force (bool. Default: False) – Controls what to do if a projection hypercube has been calculated at the emulator iteration emul_i before. If False, it will use the previously acquired projection data to create the projection figure. If True, it will recalculate all the data required to create the projection figure. Note that this will also delete all associated projection figures.
  • fig_kwargs (dict. Default: {‘figsize’: (6.4, 4.8), ‘dpi’: 100}) – Dict of keyword arguments to be used when creating the subplots figure. It takes all arguments that can be provided to the figure() function.
  • impl_kwargs_2D (dict. Default: {}) – Dict of keyword arguments to be used for making the minimum implausibility (top/left) plot in the 2D projection figures. It takes all arguments that can be provided to the plot() function.
  • impl_kwargs_3D (dict. Default: {‘cmap’: ‘rainforest_r’}) – Dict of keyword arguments to be used for making the minimum implausibility (top/left) plot in the 3D projection figures. It takes all arguments that can be provided to the hexbin() function.
  • los_kwargs_2D (dict. Default: {}) – Dict of keyword arguments to be used for making the line-of-sight (bottom/right) plot in the 2D projection figures. It takes all arguments that can be provided to the plot() function.
  • los_kwargs_3D (dict. Default: {‘cmap’: ‘blaze’}) – Dict of keyword arguments to be used for making the line-of-sight (bottom/right) plot in the 3D projection figures. It takes all arguments that can be provided to the hexbin() function.
  • line_kwargs (dict. Default: {‘linestyle’: ‘–’, ‘color’: ‘grey’}) – Dict of keyword arguments to be used for drawing the parameter estimate lines in both plots. It takes all arguments that can be provided to the draw() function.
Returns:

fig_data (dict of dicts) – Dict containing the data for every requested projection figure, split up into the ‘impl_min’ and ‘impl_los’ dicts. For 2D projections, every dict contains a list with the x and y values. For 3D projections, it contains the x, y and z values. Note that due to the figures being interpolations, the y/z values can be below zero or the line-of-sight values being above unity.

Generates (if figure is True)

A series of projection figures detailing the behavior of the model. The lay-out and output of the projection figures depend on the type of figure:

2D projection figure: The output will feature a figure with two subplots for every active model parameter (n_par). Every figure gives details about the behavior of the corresponding model parameter, by showing the minimum implausibility value (top/left) and the line-of-sight depth (bottom/right) obtained at the specified parameter value, independent of the values of the other parameters.

3D projection figure (only if n_par > 2): The output will feature a figure with two subplots for every combination of two active model parameters that can be made (n_par*(n_par-1)/2). Every figure gives details about the behavior of the corresponding model parameters, as well as their dependency on each other. This is done by showing the minimum implausibility (top/left) and the line-of-sight depth (bottom/right) obtained at the specified parameter values, independent of the values of the remaining model parameters.

Notes

If given emulator iteration emul_i has been analyzed before, the implausibility parameters of the last analysis are used. If not, then the values are used that were read in when the emulator was loaded.

run(emul_i=None, force=False)[source]

Calls the construct() method to start the construction of the given iteration of the emulator and creates the projection figures right afterward if this construction was successful.

Other Parameters:
 
  • emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the next iteration of the emulator will be constructed.
  • force (bool. Default: False) – Controls what to do if the specified emulator iteration emul_i already (partly) exists. If False, finish construction of the specified iteration or skip it if already finished. If True, reconstruct the specified iteration entirely.
File

Custom File class that has added logging and automatically uses hdf5_file as the HDF5-file to open.

Type:File
base_eval_sam

Base number of emulator evaluations used to analyze the emulator systems. This number is scaled up by the number of model parameters and the current emulator iteration to generate the true number of emulator evaluations (n_eval_sam).

Type:int
code_objects

Collection of pre-compiled built-in code snippets that are used in the _evaluate_sam_set() method.

Type:dict of code objects
comm

The global MPI intra-communicator to use in this Pipeline instance. By default, this is MPI.COMM_WORLD.

Type:Intracomm
criterion

Value indicating which criterion to use in the lhd() function.

Type:str, float or None
cut_idx

The index of the first non-wildcard in a list of implausibility values. This is equivalent to the number of wildcards leading the cut-off values in impl_cut.

Type:int
do_active_anal

Whether or not to do an active parameters analysis during the construction of the emulator systems.

Type:bool
do_logging

Whether or not to save all logging messages. If False, all logging messages of level INFO and below are ignored. It also enables/disables the progress bar for making projections.

Type:bool
emulator

The Emulator instance created during Pipeline initialization.

Type:Emulator
freeze_active_par

Whether or not previously active parameters always stay active if possible.

Type:bool
hdf5_file

Absolute path to the loaded master HDF5-file.

Type:str
impl_cut

The non-wildcard univariate implausibility cut-off values for an emulator iteration.

Type:list of int
impl_sam

The model evaluation samples that will be added to the next emulator iteration.

Type:ndarray
is_controller

Whether or not this MPI process is a controller rank. If no MPI is used, this is always True.

Type:bool
is_worker

Whether or not this MPI process is a worker rank. If no MPI is used, this is always False.

Type:bool

The ModelLink instance provided during Pipeline initialization.

Type:ModelLink
n_eval_sam

The number of evaluation samples used to analyze an emulator iteration of the emulator systems. The number of plausible evaluation samples is stored in n_impl_sam. It is zero if the specified iteration has not been analyzed yet.

Type:int
n_impl_sam

Number of model evaluation samples that passed the implausibility checks during the analysis of an emulator iteration. It is zero if the specified iteration has not been analyzed yet or has no plausible samples.

Type:int
n_sam_init

Number of evaluation samples used to construct the first iteration of the emulator systems.

Type:int
pot_active_par

The potentially active parameters. Only parameters from this list can become active during the active parameters analysis. If do_active_anal is False, all parameters in this list will be active.

Type:list of str
prism_file

Absolute path to the PRISM parameters file or None if no file was provided.

Type:str
proj_depth

Number of emulator evaluations used to generate the samples in every grid point for the last created projection figures. Note that when making 2D projections of nD models, the used depth was this number multiplied by proj_res.

Type:int
proj_res

Number of emulator evaluations used to generate the grid for the last created projection figures.

Type:int
rank

The rank of this MPI process in comm. If no MPI is used, this is always 0.

Type:int
root_dir

Absolute path to the root directory.

Type:str
size

The number of MPI processes in comm. If no MPI is used, this is always 1.

Type:int
worker_mode

Special context manager within which all code is executed in worker mode. In worker mode, all worker ranks are continuously listening for calls from the controller rank made with _make_call().

Note that all code within the context manager is executed by all ranks, with the worker ranks executing it after the controller rank exited. It is therefore advised to use an if-statement inside to make sure only the controller rank executes the code.

Using this context manager allows for easier use of PRISM in combination with serial/OpenMP codes (like MCMC methods).

Type:_GeneratorContextManager
working_dir

Absolute path to the working directory.

Type:str

Emulator

class prism.Emulator(pipeline_obj, modellink_obj)[source]

Defines the Emulator class of the PRISM package.

Description

The Emulator class is the backbone of the PRISM package, holding all tools necessary to construct, load, save and evaluate the emulator of a model. It performs many checks to see if the provided ModelLink object is compatible with the current emulator, advises the user on alternatives when certain operations are requested, automatically takes care of distributing emulator systems over MPI ranks and more.

Even though the purpose of the Emulator class is to hold only information about the emulator and therefore does not require any details about the provided ModelLink object, it will keep track of changes made to it. This is to allow the user to modify the properties of the ModelLink subclass without causing any desynchronization problems by accident.

The Emulator class requires to be linked to an instance of the Pipeline class and will automatically attempt to do so when initialized. By default, this class should only be initialized from within a Pipeline object.

__init__(pipeline_obj, modellink_obj)[source]

Initialize an instance of the Emulator class.

Parameters:
_assign_data_idx(emul_i)[source]

Determines the emulator system each data point in the provided emulator iteration emul_i should be assigned to, in order to make sure that recurring data points have the same emulator system index as in the previous emulator iteration. If multiple options are possible, data points are assigned such to spread them as much as possible.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
Returns:
  • data_to_emul_s (list of int) – The index of the emulator system that each data point should be assigned to.
  • n_emul_s (int) – The total number of active and passive emulator systems there will be in the provided emulator iteration.

Examples

If the number of data points is less than the previous iteration:

>>> emul_i = 2
>>> self._data_idx[emul_i-1]
['A', 'B', 'C', 'D', 'E']
>>> self._modellink._data_idx
['B', 'F', 'G', 'E']
>>> self._assign_data_idx(emul_i)
([1, 3, 2, 4], 5)

If the number of data points is more than the previous iteration:

>>> emul_i = 2
>>> self._data_idx[emul_i-1]
['A', 'B', 'C', 'D', 'E']
>>> self._modellink._data_idx
['B', 'F', 'G', 'E', 'A', 'C']
>>> self._assign_data_idx(emul_i)
([1, 5, 3, 4, 0, 2], 6)

If there is no previous iteration:

>>> emul_i = 1
>>> self._data_idx[emul_i-1]
[]
>>> self._modellink._data_idx
['B', 'F', 'G', 'E', 'A', 'C']
>>> self._assign_data_idx(emul_i)
([5, 4, 3, 2, 1, 0], 6)
_assign_emul_s(emul_i)[source]

Determines which emulator systems (files) should be assigned to which MPI rank in order to balance the number of active emulator systems on every rank for every iteration up to the provided emulator iteration emul_i. If multiple choices can achieve this, the emulator systems are automatically spread out such that the total number of active emulator systems on a single rank is also balanced as much as possible.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
Returns:emul_s (list of int) – A list containing the emulator systems that have been assigned to the corresponding MPI rank by the controller.

Notes

Currently, this function only uses high-level MPI. Additional speed can be obtained by also implementing low-level MPI, which will potentially be done in the future.

_cleanup_emul_files(emul_i)[source]

Opens all emulator HDF5-files and removes the provided emulator iteration emul_i and subsequent iterations from them. Also removes any related projection figures that have default names. If emul_i == 1, all emulator HDF5-files are removed instead.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
_construct_iteration(emul_i)[source]

Constructs the emulator iteration corresponding to the provided emul_i, by performing the given emulation method and pre-calculating the prior expectation and variance values of the used model evaluation samples.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.

Generates

All data sets that are required to evaluate the emulator at the constructed iteration.

_create_new_emulator()[source]

Creates a new master HDF5-file that holds all the information of a new emulator and writes all important emulator details to it. Afterwards, resets all loaded emulator data and prepares the HDF5-file and emulator for the construction of the first emulator iteration.

Generates

A new master HDF5-file ‘prism.hdf5’ contained in the working directory specified in the Pipeline instance, holding all information required to construct the first iteration of the emulator.

_do_regression(emul_i, emul_s_seq)[source]

Performs a forward stepwise linear regression for all requested emulator systems emul_s_seq in the provided emulator iteration emul_i. Calculates what the expectation values of all polynomial coefficients are. The polynomial order that is used in the regression depends on poly_order.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.

Generates (for every emulator system)

rsdl_var : float
Residual variance of the regression function.
regr_score : float
Fit-score of the regression function.
poly_coef : 1D ndarray object
Array containing the expectation values of the non-zero polynomial coefficients.
poly_powers : 2D ndarray object
Array containing the powers of the non-zero polynomial terms in the regression function.
poly_idx : 1D ndarray object
Array containing the indices of the non-zero polynomial terms in the regression function.
poly_coef_cov : 1D ndarray object (if use_regr_cov is True)
Array containing the covariance values of the non-zero polynomial coefficients.
_evaluate(emul_i, emul_s_seq, par_set)[source]

Evaluates the emulator systems emul_s_seq at iteration emul_i for given par_set.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set (1D ndarray object) – Model parameter value set to evaluate the emulator at.
Returns:

  • adj_exp_val (1D ndarray object) – Adjusted emulator expectation value for all requested emulator systems on this MPI rank.
  • adj_var_val (1D ndarray object) – Adjusted emulator variance value for all requested emulator systems on this MPI rank.

_get_active_par(emul_i, emul_s_seq)[source]

Determines the active parameters to be used for every emulator system listed in emul_s_seq in the provided emulator iteration emul_i. Uses backwards stepwise elimination to determine the set of active parameters. The polynomial order that is used in the stepwise elimination depends on poly_order.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.

Generates (for every emulator system)

active_par_data : 1D ndarray object
Array containing the indices of all the parameters that are active in the emulator iteration emul_i.
_get_adj_exp(emul_i, emul_s_seq, par_set, cov_vec)[source]

Calculates the adjusted emulator expectation values for requested emulator systems emul_s_seq at a given emulator iteration emul_i for specified parameter set par_set and corresponding covariance vector cov_vec.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set (1D ndarray object) – Model parameter value set to calculate the adjusted emulator expectation for.
  • cov_vec (2D ndarray object) – Covariance vector corresponding to par_set.
Returns:

adj_exp_val (1D ndarray object) – Adjusted emulator expectation value for all requested emulator systems on this MPI rank.

_get_adj_var(emul_i, emul_s_seq, par_set, cov_vec)[source]

Calculates the adjusted emulator variance values for requested emulator systems emul_s_seq at a given emulator iteration emul_i for specified parameter set par_set and corresponding covariance vector cov_vec.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set (1D ndarray object) – Model parameter value set to calculate the adjusted emulator variance for.
  • cov_vec (2D ndarray object) – Covariance vector corresponding to par_set.
Returns:

adj_var_val (1D ndarray object) – Adjusted emulator variance value for all requested emulator systems on this MPI rank.

_get_cov(emul_i, emul_s_seq, par_set1, par_set2)[source]

Calculates the full emulator covariances for requested emulator systems emul_s_seq at emulator iteration emul_i for given parameter sets par_set1 and par_set2. The contributions to these covariances depend on method.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set1, par_set2 (1D ndarray object or None) – If par_set1 and par_set2 are both not None, calculate covariances for par_set1 with par_set2. If par_set1 is not None and par_set2 is None, calculate covariances for par_set1 with sam_set (covariance vector). If par_set1 and par_set2 are both None, calculate covariances for sam_set (covariance matrix). When not None, par_set is the model parameter value set to calculate the covariances for.
Returns:

cov (1D, 2D or 3D ndarray object) – Depending on the arguments provided, a covariance value, vector or matrix for requested emulator systems.

_get_cov_matrix(emul_i, emul_s_seq)[source]

Calculates the (inverse) matrix of covariances between known model evaluation samples for requested emulator systems emul_s_seq at emulator iteration emul_i.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.

Generates

cov_mat : 3D ndarray object
Matrix containing the covariances between all known model evaluation samples for requested emulator systems.
cov_mat_inv : 3D ndarray object
Inverse of covariance matrix for requested emulator systems.
_get_default_parameters()[source]

Generates a dict containing default values for all emulator parameters.

Returns:par_dict (dict) – Dict containing all default emulator parameter values.
_get_emul_i(emul_i, cur_iter)[source]

Checks if the provided emulator iteration emul_i can be requested or replaces it if None was provided. This method requires all MPI ranks to call it simultaneously.

Parameters:
  • emul_i (int or None) – Number indicating the requested emulator iteration.
  • cur_iter (bool) – Bool determining whether the current (True) or the next (False) emulator iteration is requested.
Returns:

emul_i (int) – The requested emulator iteration that passed the check.

_get_exp_dot_term(emul_i, emul_s_seq)[source]

Pre-calculates the second expectation adjustment dot-term for requested emulator systems emul_s_seq at a given emulator iteration emul_i for all model evaluation samples and saves it for later use.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.

Generates

exp_dot_term : 2D ndarray object
2D array containing the pre-calculated values for the second adjustment dot-term of the adjusted expectation for requested emulator systems.
_get_inv_matrix(matrix)[source]

Calculates the inverse of a given matrix. Right now only uses the inv() function.

Parameters:matrix (2D array_like) – Matrix to be inverted.
Returns:matrix_inv (2D ndarray object) – Inverse of the given matrix.
_get_prior_exp(emul_i, emul_s_seq, par_set)[source]

Calculates the prior expectation value for requested emulator systems emul_s_seq at a given emulator iteration emul_i for specified parameter set par_set. This expectation depends on method.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set (1D ndarray object or None) – If None, calculate the prior expectation values of sam_set. If not None, calculate the prior expectation value for the given model parameter value set.
Returns:

prior_exp (1D or 2D ndarray object) – Prior expectation values for either sam_set or par_set for requested emulator systems.

_get_regr_cov(emul_i, emul_s_seq, par_set1, par_set2)[source]

Calculates the covariances of the regression function for requested emulator systems emul_s_seq at emulator iteration emul_i for given parameter sets par_set1 and par_set2.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • emul_s_seq (list of int) – List of numbers indicating the requested emulator systems.
  • par_set1, par_set2 (1D ndarray object or None) – If par_set1 and par_set2 are both not None, calculate regression covariances for par_set1 with par_set2. If par_set1 is not None and par_set2 is None, calculate regression covariances for par_set1 with sam_set (covariance vector). If par_set1 and par_set2 are both None, calculate regression covariances for sam_set (covariance matrix). When not None, par_set is the model parameter value set to calculate the regression covariances for.
Returns:

regr_cov (1D, 2D or 3D ndarray object) – Depending on the arguments provided, a regression covariance value, vector or matrix for requested emulator systems.

_load_data(emul_i)[source]

Loads in all the important emulator data up to emulator iteration emul_i into memory.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.

Generates

All relevant emulator data up to emulator iteration emul_i is loaded into memory.

_load_emulator(modellink_obj)[source]

Checks if the provided working directory contains a constructed emulator and loads in the emulator data accordingly.

Parameters:modellink_obj (ModelLink object) – Instance of the ModelLink class that links the emulated model to this Pipeline object.
_prepare_new_iteration(emul_i)[source]

Prepares the emulator for the construction of a new iteration emul_i. Checks if this iteration can be prepared or if it has been prepared before, and acts accordingly.

Parameters:emul_i (int) – Number indicating the requested emulator iteration.
Returns:reload (bool) – Bool indicating whether or not the controller rank of the Pipeline instance needs to reload its data.

Generates

A new group in the master HDF5-file with the emulator iteration as its name, containing subgroups corresponding to all emulator systems that will be used in this iteration.

Notes

Preparing an iteration that has been prepared before, causes that and all subsequent iterations of the emulator to be deleted. A check is carried out to see if it was necessary to reprepare the requested iteration and a warning is given if this check fails.

_read_data_idx(emul_s_group)[source]

Reads in and combines the parts of the data point identifier that is assigned to the provided emul_s_group.

Parameters:emul_s_group (Group object) – The HDF5-group from which the data point identifier needs to be read in.
Returns:data_idx (tuple of {int, float, str}) – The combined data point identifier.
_read_parameters()[source]

Reads in the Emulator parameters from the provided PRISM parameter file and saves them in the current Emulator instance.

_retrieve_parameters()[source]

Reads in the emulator parameters from the provided working directory and saves them in the current Emulator instance.

_save_data(emul_i, lemul_s, data_dict)[source]

Saves a given data dict {keyword: data} at the given emulator iteration emul_i and local emulator system lemul_s to the HDF5-file and as an data attribute to the current Emulator instance.

Parameters:
  • emul_i (int) – Number indicating the requested emulator iteration.
  • lemul_s (int or None) – Number indicating the requested local emulator system. If None, use the master emulator file instead.
  • data_dict (dict) – Dict containing the data that needs to be saved to the HDF5-file.
Keyword Arguments:
 
  • keyword ({‘active_par’, ‘active_par_data’, ‘cov_mat’, ‘exp_dot_term’, ‘mod_real_set’, ‘regression’}) – String specifying the type of data that needs to be saved.
  • data ({int, float, str, array_like} or dict) – The actual data that needs to be saved at data keyword keyword. If dict, save every item individually.

Generates

The specified data is saved to the HDF5-file.

_set_mock_data()[source]

Loads previously used mock data into the ModelLink object, overwriting the parameter estimates, data values, data errors, data spaces and data identifiers with their mock equivalents.

Generates

Overwrites the corresponding ModelLink class properties with the previously used values (taken from the first emulator iteration).

Sets the ModelLink object that will be used for constructing this emulator. If a constructed emulator is present, checks if provided modellink_obj argument matches the ModelLink subclass used to construct it.

Parameters:
  • modellink_obj (ModelLink object) – Instance of the ModelLink class that links the emulated model to this Pipeline object. The provided ModelLink object must match the one used to construct the loaded emulator.
  • modellink_loaded (str or None) – If str, the name of the ModelLink subclass that was used to construct the loaded emulator. If None, no emulator is loaded.
active_emul_s

The indices of the emulator systems on this MPI rank that are active.

Type:list of int
active_par

The model parameter names that are considered active. Only available on the controller rank.

Type:list of str
active_par_data

The model parameter names that are considered active for every emulator system on this MPI rank.

Type:list of str
ccheck

The emulator system components that are still required to complete the construction of an emulator iteration on this MPI rank. The controller rank additionally lists the required components that are emulator iteration specific (‘mod_real_set’ and ‘active_par’).

Type:list of str
cov_mat_inv

The inverses of the covariance matrices for every emulator system on this MPI rank.

Type:list of ndarray
emul_i

The last emulator iteration that is fully constructed for all emulator systems on this MPI rank.

Type:int
emul_load

Whether or not a previously constructed emulator is currently loaded.

Type:bool
emul_s

The indices of the emulator systems that are assigned to this MPI rank.

Type:list of int
emul_s_to_core

List of the indices of the emulator systems that are assigned to every MPI rank. Only available on the controller rank.

Type:list of lists
emul_type

The type of emulator that is currently loaded.

Type:str
exp_dot_term

The second expectation adjustment dot-term values of all model evaluation samples for every emulator system on this MPI rank.

Type:list of ndarray
l_corr

The Gaussian correlation lengths for all model parameters, which is defined as the maximum distance between two values of a specific model parameter within which the Gaussian contribution to the correlation between the values is still significant.

Type:ndarray
method

The emulation method to use for constructing the emulator. Possible are ‘gaussian’, ‘regression’ and ‘full’.

Type:str
mod_set

The model outputs corresponding to the samples in sam_set for every emulator system on this MPI rank.

Type:list of ndarray
n_cross_val

Number of (k-fold) cross-validations that are used for determining the quality of the regression process. It is set to zero if cross-validations are not used. If method == ‘gaussian’ and do_active_anal is False, this number is not required.

Type:int
n_emul_s

Number of emulator systems assigned to this MPI rank.

Type:int
n_emul_s_tot

Total number of emulator systems assigned to all MPI ranks combined. Only available on the controller rank.

Type:int
n_sam

Number of model evaluation samples that have been/will be used to construct an emulator iteration.

Type:int
poly_coef

The non-zero coefficients for the polynomial terms in the regression function for every emulator system on this MPI rank. Empty if method == ‘gaussian’.

Type:list of ndarray
poly_coef_cov

The covariances for all coefficients in poly_coef for every emulator system on this MPI rank. Empty if method == ‘gaussian’ or use_regr_cov is False.

Type:list of ndarray
poly_idx

The indices for all polynomial terms with non-zero coefficients in the regression function for every emulator system on this MPI rank. Empty if method == ‘gaussian’.

Type:list of ndarray
poly_order

Polynomial order that is considered for the regression process. If method == ‘gaussian’ and do_active_anal is False, this number is not required.

Type:int
poly_powers

The powers for all polynomial terms with non-zero coefficients in the regression function for every emulator system on this MPI rank. Empty if method == ‘gaussian’.

Type:list of ndarray
rsdl_var

The residual variances for every emulator system on this MPI rank. Obtained from regression process and replaces the Gaussian sigma. Empty if method == ‘gaussian’.

Type:list of float
sam_set

The model evaluation samples that have been/will be used to construct the specified emulator iteration.

Type:ndarray
sigma

Value of the Gaussian sigma. If method != ‘gaussian’, this value is not required, since it is obtained from the regression process instead.

Type:float
use_mock

Whether or not mock data has been used for the construction of this emulator instead of actual data. If True, changes made to the data in the provided ModelLink object are ignored.

Type:bool
use_regr_cov

Whether or not to take into account the regression covariance when calculating the covariance of the emulator, in addition to the Gaussian covariance. If method == ‘gaussian’, this bool is not required. If method == ‘regression’, this bool is always set to True.

Type:bool

Utilities

Provides a collection of functions useful for using/mixing PRISM with other applications.

prism.utils.get_lnpost_fn(ext_lnpost, pipeline_obj, emul_i=None, unit_space=True, hybrid=True)[source]

Returns a function definition get_lnpost(par_set, *args, **kwargs).

This get_lnpost function can be used to calculate the natural logarithm of the posterior probability, which analyzes a given par_set first in the provided pipeline_obj at iteration emul_i and passes it to the ext_lnpost function if it is plausible.

This function needs to be called by all MPI ranks.

Parameters:
  • ext_lnpost (function) – Function definition that needs to be called if the provided par_set is plausible in iteration emul_i of pipeline_obj. The used call signature is ext_lnpost(par_set, *args, **kwargs). All MPI ranks will call this function unless called within the worker_mode context manager.
  • pipeline_obj (Pipeline object) – The instance of the Pipeline class that needs to be used for determining the validity of the proposed sampling step.
Other Parameters:
 
  • emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.
  • unit_space (bool. Default: True) – Bool determining whether or not the provided sample will be given in unit space.
  • hybrid (bool. Default: True) – Bool determining whether or not the get_lnpost function should use the implausibility values of a given par_set as an additional prior.
Returns:

get_lnpost (function) – Definition of the function get_lnpost(par_set, *args, **kwargs).

See also

get_walkers(): Analyzes proposed init_walkers and returns valid p0_walkers. worker_mode: Special context manager within which all code is executed in worker mode.

Warning

Calling this function factory will disable all regular logging in pipeline_obj (do_logging set to False), in order to avoid having the same message being logged every time get_lnpost is called.

prism.utils.get_walkers(pipeline_obj, emul_i=None, init_walkers=None, unit_space=True, ext_lnpost=None, **kwargs)[source]

Analyzes proposed init_walkers and returns valid p0_walkers.

Analyzes sample set init_walkers in the provided pipeline_obj at iteration emul_i and returns all samples that are plausible to be used as MCMC walkers. The provided samples and returned walkers should be/are given in unit space if unit_space is True.

If init_walkers is None, returns impl_sam instead if it is available.

This function needs to be called by all MPI ranks.

Parameters:

pipeline_obj (Pipeline object) – The instance of the Pipeline class that needs to be used for determining the validity of the proposed walkers.

Other Parameters:
 
  • emul_i (int or None. Default: None) – If int, number indicating the requested emulator iteration. If None, the last iteration of the emulator will be used.
  • init_walkers (2D array_like, int or None. Default: None) – Sample set of proposed initial MCMC walker positions. All plausible samples in init_walkers will be returned. If int, generate an LHD of provided size and return all plausible samples. If None, return impl_sam corresponding to iteration emul_i instead.
  • unit_space (bool. Default: True) – Bool determining whether or not the provided samples and returned walkers are given in unit space.
  • ext_lnpost (function or None. Default: None) – If function, call get_lnpost_fn() function factory using ext_lnpost and the same values for pipeline_obj, emul_i and unit_space, and return the resulting function definition get_lnpost. Any additionally provided kwargs are also passed to it.
Returns:

  • n_walkers (int) – Number of returned MCMC walkers.
  • p0_walkers (2D ndarray object) – Array containing starting positions of valid MCMC walkers.
  • get_lnpost (function (if ext_lnpost is a function)) – The function returned by get_lnpost_fn() function factory using ext_lnpost, pipeline_obj, emul_i, unit_space and kwargs as the input values.

See also

get_lnpost_fn(): Returns a function definition get_lnpost(par_set, *args, **kwargs). worker_mode: Special context manager within which all code is executed in worker mode.

Notes

If init_walkers is None and emulator iteration emul_i has not been analyzed yet, a RequestError will be raised.

Internal

Contains a collection of support classes/functions/lists for the PRISM package.

exception prism._internal.RequestError[source]

Generic exception raised for invalid action requests in the PRISM pipeline.

General purpose exception class, raised whenever a requested action cannot be executed due to it not being allowed or possible in the current state of the Pipeline instance.

exception prism._internal.RequestWarning[source]

Generic warning raised for (future) action requests in the PRISM pipeline that may not be useful.

General purpose warning class, raised whenever a requested action may not produce appropriate or expected results due to the current state of the Pipeline instance. It is also raised if an obtained result can lead to such an action in the future.

class prism._internal.CFilter(MPI_rank)[source]

Custom Filter class that only allows the controller rank to log messages to the logfile. Calls from worker ranks are ignored.

class prism._internal.CLogger(*args, **kwargs)[source]

Custom Logger class that uses the CFilter.

class prism._internal.PRISM_Comm(comm=None)[source]

Custom Intracomm class that automatically makes use of the ndarray buffers when using communications. Is functionally the same as the provided comm for everything else.

Other Parameters:
 comm (Intracomm object or None. Default: None) – The MPI intra-communicator to use in this PRISM_Comm instance. If None, use MPI.COMM_WORLD instead.
bcast(obj, root)[source]

Special broadcast method that automatically uses the appropriate method (bcast() or Bcast()) depending on the type of the provided obj.

Parameters:
  • obj (ndarray or object) – The object to broadcast to all MPI ranks. If ndarray, use Bcast(). If not, use bcast() instead.
  • root (int) – The MPI rank that broadcasts obj.
Returns:

obj (object) – The broadcasted obj.

gather(obj, root)[source]

Special gather method that automatically uses the appropriate method (gather() or Gatherv()) depending on the type of the provided obj.

Parameters:
  • obj (ndarray or object) – The object to gather from all MPI ranks. If ndarray, use Gatherv(). If not, use gather() instead.
  • root (int) – The MPI rank that gathers obj.
Returns:

obj (list or None) – If MPI rank is root, returns a list of the gathered objects. Else, returns None.

Notes

If some but not all MPI ranks use a NumPy array, this method will hang indefinitely. When gathering NumPy arrays, all arrays must have the same number of dimensions and the same shape, except for one axis.

class prism._internal.RFilter(MPI_rank)[source]

Custom Filter class that prepends the rank of the MPI process that calls it to the logging message.

class prism._internal.RLogger(*args, **kwargs)[source]

Custom Logger class that uses the RFilter if the size of the intra-communicator is more than 1.

prism._internal.check_compatibility(emul_version)[source]

Checks if the provided emul_version is compatible with the current version of PRISM. Raises a RequestError if False and indicates which version of PRISM still supports the provided emul_version.

prism._internal.check_instance(instance, cls)[source]

Checks if provided instance has been initialized from a proper cls (sub)class. Raises a TypeError if instance is not an instance of cls.

Parameters:
  • instance (object) – Class instance that needs to be checked.
  • cls (class) – The class which instance needs to be properly initialized from.
Returns:

result (bool) – Bool indicating whether or not the provided instance was initialized from a proper cls (sub)class.

prism._internal.check_vals(values, name, *args)[source]

Checks if all values in provided input argument values with name meet all criteria given in args. If no criteria are given, it is checked if all values are finite. Returns values (0 or 1 in case of bool) if True and raises a ValueError or TypeError if False.

Parameters:
  • values (array_like of {int, float, str, bool}) – The values to be checked against all given criteria in args. It must be possible to convert values to a ndarray object.
  • name (str) – The name of the input argument, which is used in the error message if a criterion is not met.
  • args (tuple of {‘bool’, ‘float’, ‘int’, ‘neg’, ‘nneg’, ‘normal’, ‘npos’, ‘nzero’, ‘pos’, ‘str’}) – Sequence of strings determining the criteria that values must meet. If args is empty, it is checked if values are finite.
Returns:

return_values (array_like of {int, float, str}) – If args contained ‘bool’, returns 0s or 1s. Else, returns values.

Notes

If values contains integers, but args contains ‘float’, return_values will be cast as float.

prism._internal.convert_str_seq(seq)[source]

Converts a provided sequence to a string, removes all auxiliary characters from it, splits it up into individual elements and converts all elements back to integers, floats and/or strings.

Parameters:seq (str or array_like) – The sequence that needs to be converted to individual elements. If array_like, seq is first converted to a string.
Returns:new_seq (list) – A list with all individual elements converted to integers, floats and/or strings.
prism._internal.delist(list_obj)[source]

Returns a copy of list_obj with all empty lists and tuples removed.

Parameters:list_obj (list) – A list object that requires its empty list/tuple elements to be removed.
Returns:delisted_copy (list) – Copy of list_obj with all empty lists/tuples removed.
prism._internal.docstring_append(addendum, join='')[source]

Custom decorator that allows a given string addendum to be appended to the docstring of the target function, separated by a given string join.

prism._internal.docstring_copy(source)[source]

Custom decorator that allows the docstring of a function source to be copied to the target function.

prism._internal.docstring_substitute(*args, **kwargs)[source]

Custom decorator that allows either given positional arguments args or keyword arguments kwargs to be substituted into the docstring of the target function.

prism._internal.get_PRISM_File(prism_hdf5_file)[source]

Returns a class definition PRISM_File(mode, emul_s=None, **kwargs).

This class definition is a specialized version of the File class with the filename automatically set to prism_hdf5_file and added logging to the constructor and destructor methods.

Parameters:prism_hdf5_file (str) – Absolute path to the master HDF5-file that is used in a Pipeline instance.
Returns:PRISM_File (class) – Definition of the class PRISM_File(mode, emul_s=None, **kwargs).
prism._internal.get_info()[source]

Prints a string that gives an overview of all information relevant to the PRISM package distribution.

prism._internal.getCLogger(name=None)[source]

Create a CLogger instance with name and return it.

prism._internal.getRLogger(name=None)[source]

Create a RLogger instance with name and return it.

prism._internal.import_cmaps(cmap_dir)[source]

Reads in custom colormaps from a provided directory cmap_dir, transforms them into LinearSegmentedColormap objects and registers them in the cm module. Both the imported colormap and its reversed version will be registered.

This function is called automatically when PRISM is imported.

Parameters:cmap_dir (str) – Relative or absolute path to the directory that contains custom colormap files. A colormap file can be a NumPy binary file (‘.npy’ or ‘.npz’) or any text file.

Notes

All colormap files in cmap_dir must have names starting with ‘cm_’. The resulting colormaps will have the name of their file without the prefix and extension.

prism._internal.move_logger(working_dir, filename)[source]

Moves the logging file filename from the current working directory to the given working_dir, and then restarts it again.

Parameters:
  • working_dir (str) – String containing the directory the log-file needs to be moved to.
  • filename (str) – String containing the name of the log-file that needs to be moved.
prism._internal.np_array(obj, *args, **kwargs)[source]

Returns np.array(obj, *args, copy=False, **kwargs).

prism._internal.raise_error(err_msg, err_type=<class 'Exception'>, logger=None)[source]

Raises a given error err_msg of type err_type and logs the error using the provided logger.

Parameters:

err_msg (str) – The message included in the error.

Other Parameters:
 
  • err_type (Exception subclass. Default: Exception) – The type of error that needs to be raised.
  • logger (Logger object or None. Default: None) – The logger to which the error message must be written. If None, the RootLogger logger is used instead.
prism._internal.raise_warning(warn_msg, warn_type=<class 'UserWarning'>, logger=None, stacklevel=1)[source]

Raises a given warning warn_msg of type warn_type and logs the warning using the provided logger.

Parameters:

warn_msg (str) – The message included in the warning.

Other Parameters:
 
  • warn_type (Warning subclass. Default: UserWarning) – The type of warning that needs to be raised.
  • logger (Logger object or None. Default: None) – The logger to which the warning message must be written. If None, the RootLogger logger is used instead.
  • stacklevel (int. Default: 1) – The stack level of the warning message at the location of this function call. The actual used stack level is increased by one.
prism._internal.rprint(*args, **kwargs)[source]

Custom print() function that prepends the rank of the MPI process that calls it to the message if the size of the intra-communicator is more than 1. Takes the same input arguments as the normal print() function.

prism._internal.start_logger(filename=None, mode='w')[source]

Opens a logging file called filename in the current working directory, opened with mode and starts the logger.

Other Parameters:
 
  • filename (str or None. Default: None) – String containing the name of the log-file that is opened. If None, a new log-file will be created.
  • mode ({‘r’, ‘r+’, ‘w’, ‘w-‘/’x’, ‘a’}. Default: ‘w’) – String indicating how the log-file needs to be opened.