PyMC (v5+) interface¶
This celerite2.pymc
submodule provides access to the celerite2 models
within the PyTensor framework. Of special
interest, this adds support for probabilistic model building using PyMC v5 or later.
Note: PyMC v4 was a short-lived version of PyMC with the Aesara backend. Celerite2 now only supports to PyMC 5, but past releases of celerite2 might work with Aesara.
The Getting started tutorial demonstrates the use of this interface, while this
page provides the details for the celerite2.pymc.GaussianProcess
class
which provides all this functionality. This page does not include documentation
for the term models defined in PyTensor, but you can refer to the
Model building section of the Python interface documentation. All of those
models are implemented in PyTensor and you can access them using something like
the following:
import pytensor.tensor as pt
from celerite2.pymc import GaussianProcess, terms
term = terms.SHOTerm(S0=at.scalar(), w0=at.scalar(), Q=at.scalar())
gp = GaussianProcess(term)
The celerite2.pymc.GaussianProcess
class is detailed below:
- class celerite2.pymc.GaussianProcess(kernel, t=None, *, mean=0.0, **kwargs)[source]¶
- apply_inverse(y, *, inplace=False)¶
Apply the inverse of the covariance matrix to a vector or matrix
Solve
K.x = y
forx
whereK
is the covariance matrix of the GP.Note
The mean function is not applied in this method.
- Parameters:
y (shape[N] or shape[N, M]) – The vector or matrix
y
described above.inplace (bool, optional) – If
True
,y
will be overwritten with the resultx
.
- Raises:
RuntimeError – If
GaussianProcess.compute()
is not called first.ValueError – When the inputs are not valid (shape, number, etc.).
- compute(t, *, yerr=None, diag=None, check_sorted=True, quiet=False)¶
Compute the Cholesky factorization of the GP covariance matrix
- Parameters:
t (shape[N]) – The independent coordinates of the observations. This must be sorted in increasing order.
yerr (shape[N], optional) – If provided, the diagonal standard deviation of the observation model.
diag (shape[N], optional) – If provided, the diagonal variance of the observation model.
check_sorted (bool, optional) – If
True
, a check is performed to make sure thatt
is correctly sorted. AValueError
will be thrown when this check fails.quiet (bool, optional) – If
True
, when the matrix cannot be factorized (because of numerics or otherwise) the solver’sLinAlgError
will be silenced and the determiniant will be set to zero. Otherwise, the exception will be propagated.
- Raises:
ValueError – When the inputs are not valid (shape, number, etc.).
LinAlgError – When the matrix is not numerically positive definite.
- conditional(name, y, t=None, include_mean=True, kernel=None, **kwargs)[source]¶
Add a variable representing the conditional density to a PyMC model
Note
The performance of this method will generally be poor since the sampler will numerically sample this parameter. Depending on your use case, you might be better served by tracking the results of
GaussianProcess.predict()
usingDeterministic
variables and computing the predictions as a postprocessing step.- Parameters:
name (str) – The name of the random variable
y (shape[N]) – The observations at coordinates
x
fromGausianProcess.compute()
.t (shape[M], optional) – The independent coordinates where the prediction should be made. If this is omitted the coordinates will be assumed to be
x
fromGaussianProcess.compute()
and an efficient method will be used to compute the mean prediction.include_mean (bool, optional) – Include the mean function in the prediction.
kernel (optional) – If provided, compute the conditional distribution using a different kernel. This is generally used to separate the contributions from different model components.
- Returns:
A
pm.MvNormal
distribution representing the conditional density.
- dot_tril(y, *, inplace=False)¶
Dot the Cholesky factor of the GP system into a vector or matrix
Compute
x = L.y
whereK = L.L^T
andK
is the covariance matrix of the GP.Note
The mean function is not applied in this method.
- Parameters:
y (shape[N] or shape[N, M]) – The vector or matrix
y
described above.inplace (bool, optional) – If
True
,y
will be overwritten with the resultx
.
- Raises:
RuntimeError – If
GaussianProcess.compute()
is not called first.ValueError – When the inputs are not valid (shape, number, etc.).
- log_likelihood(y, *, inplace=False)¶
Compute the marginalized likelihood of the GP model
The factorized matrix from the previous call to
GaussianProcess.compute()
is used so that method must be called first.- Parameters:
y (shape[N]) – The observations at coordinates
t
as defined byGaussianProcess.compute()
.inplace (bool, optional) – If
True
,y
will be overwritten in the process of the calculation. This will reduce the memory footprint, but should be used with care since this will overwrite the data.
- Raises:
RuntimeError – If
GaussianProcess.compute()
is not called first.ValueError – When the inputs are not valid (shape, number, etc.).
- marginal(name, **kwargs)[source]¶
Add the marginal likelihood to a PyMC model
- Parameters:
name (str) – The name of the random variable.
observed (optional) – The observed data
- Returns:
A
celerite2.pymc.CeleriteNormal
distribution representing the marginal likelihood.
- predict(y, t=None, *, return_cov=False, return_var=False, include_mean=True, kernel=None)¶
Compute the conditional distribution
The factorized matrix from the previous call to
GaussianProcess.compute()
is used so that method must be called first.- Parameters:
y (shape[N]) – The observations at coordinates
t
as defined byGaussianProcess.compute()
.t (shape[M], optional) – The independent coordinates where the prediction should be evaluated. If not provided, this will be evaluated at the observations
t
fromGaussianProcess.compute()
.return_var (bool, optional) – Return the variance of the conditional distribution.
return_cov (bool, optional) – Return the full covariance matrix of the conditional distribution.
include_mean (bool, optional) – Include the mean function in the prediction.
kernel (optional) – If provided, compute the conditional distribution using a different kernel. This is generally used to separate the contributions from different model components. Note that the computational cost and scaling will be worse when using this parameter.
- Raises:
RuntimeError – If
GaussianProcess.compute()
is not called first.ValueError – When the inputs are not valid (shape, number, etc.).
PyMC (v5+) support¶
This implementation comes with a custom PyMC Distribution
that represents a
multivariate normal with a celerite covariance matrix. This is used by the
celerite2.pymc.GaussianProcess.marginal()
method documented above which
adds a marginal likelihood node to a PyMC model.
- class celerite2.pymc.distribution.CeleriteNormal(name: str, *args, rng=None, dims: str | Sequence[str | None] | None = None, initval=None, observed=None, total_size=None, transform=UNSET, default_transform=UNSET, **kwargs)[source]¶
A multivariate normal distribution with a celerite covariance matrix