Low-level Solver Interface¶
This is the low-level interface to the C++ implementation of the celerite algorithm. These methods do most of the heavy lifting but most users shouldn’t need to call these directly. This interface was built using pybind11.
-
celerite.solver.
get_kernel_value
(arg0: numpy.ndarray[float64[m, 1]], arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]], arg6: numpy.ndarray[float64]) → object¶ Get the value of the kernel for given parameters and lags
- Parameters
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
tau (array[n]) – The time lags where the kernel should be evaluated.
- Returns
The kernel evaluated at
tau
.- Return type
array[n]
-
celerite.solver.
get_psd_value
(arg0: numpy.ndarray[float64[m, 1]], arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]], arg6: numpy.ndarray[float64]) → object¶ Get the PSD of the kernel for given parameters and angular frequencies
- Parameters
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
omega (array[n]) – The frequencies where the PSD should be evaluated.
- Returns
The PSD evaluated at
omega
.- Return type
array[n]
-
celerite.solver.
check_coefficients
(arg0: numpy.ndarray[float64[m, 1]], arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]]) → bool¶ Apply Sturm’s theorem to check if parameters yield a positive PSD
- Parameters
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
- Returns
True
if the PSD is everywhere positive.- Return type
-
class
celerite.solver.
CholeskySolver
¶ A thin wrapper around the C++ CholeskySolver class
-
compute
(self: celerite.solver.CholeskySolver, arg0: float, arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]], arg6: numpy.ndarray[float64[m, 1]], arg7: numpy.ndarray[float64[m, 1]], arg8: numpy.ndarray[float64[m, n]], arg9: numpy.ndarray[float64[m, n]], arg10: numpy.ndarray[float64[m, 1]], arg11: numpy.ndarray[float64[m, 1]]) → None¶ Compute the Cholesky factorization of the matrix
- Parameters
jitter (float) – The jitter of the kernel.
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
x (array[n]) – The _sorted_ array of input coordinates.
diag (array[n]) – An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,
diag
should be the measurement _variance_ (i.e. sigma^2).
-
computed
(self: celerite.solver.CholeskySolver) → bool¶ A flag that indicates if
compute
has been executed- Returns
True
ifsolver.CholeskySolver.compute()
was previously executed successfully.- Return type
-
dot
(self: celerite.solver.CholeskySolver, arg0: float, arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]], arg6: numpy.ndarray[float64[m, 1]], arg7: numpy.ndarray[float64[m, 1]], arg8: numpy.ndarray[float64[m, n]], arg9: numpy.ndarray[float64[m, n]], arg10: numpy.ndarray[float64[m, 1]], arg11: numpy.ndarray[float64[m, n]]) → numpy.ndarray[float64[m, n]]¶ Compute the dot product of a
celerite
matrix and another arbitrary matrixThis method computes
A.b
whereA
is defined by the parameters andb
is an arbitrary matrix of the correct shape.- Parameters
jitter (float) – The jitter of the kernel.
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
x (array[n]) – The _sorted_ array of input coordinates.
b (array[n] or array[n, neq]) – The matrix
b
described above.
- Returns
The dot product
A.b
as described above.- Return type
array[n] or array[n, neq]
- Raises
ValueError – For mismatched dimensions.
-
dot_L
(self: celerite.solver.CholeskySolver, arg0: numpy.ndarray[float64[m, n]]) → numpy.ndarray[float64[m, n]]¶ Compute the dot product of the square root of a
celerite
matrixThis method computes
L.z
whereA = L.L^T
is the matrix defined incompute
.- Parameters
z (array[n] or array[n, neq]) – The matrix
z
described above.- Returns
The dot product
L.b
as described above.- Return type
array[n] or array[n, neq]
- Raises
ValueError – For mismatched dimensions.
-
dot_solve
(self: celerite.solver.CholeskySolver, arg0: numpy.ndarray[float64[m, n]]) → float¶ Solve the system
b^T . A^-1 . b
A previous call to
solver.Solver.compute()
defines a matrixA
and this method solvesb^T . A^-1 . b
for a vectorb
.- Parameters
b (array[n]) – The right hand side of the linear system.
- Returns
The solution of
b^T . A^-1 . b
.- Return type
- Raises
ValueError – For mismatched dimensions.
-
grad_log_likelihood
(self: celerite.solver.CholeskySolver, arg0: float, arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]], arg3: numpy.ndarray[float64[m, 1]], arg4: numpy.ndarray[float64[m, 1]], arg5: numpy.ndarray[float64[m, 1]], arg6: numpy.ndarray[float64[m, 1]], arg7: numpy.ndarray[float64[m, 1]], arg8: numpy.ndarray[float64[m, n]], arg9: numpy.ndarray[float64[m, n]], arg10: numpy.ndarray[float64[m, 1]], arg11: numpy.ndarray[float64[m, 1]], arg12: numpy.ndarray[float64[m, 1]]) → None¶ Compute the gradient of the log likelihood of the model using autodiff
The returned gradient is with respect to the jitter and the coefficients.
- Parameters
jitter (float) – The jitter of the kernel.
a_real (array[j_real]) – The coefficients of the real terms.
c_real (array[j_real]) – The exponents of the real terms.
a_comp (array[j_complex]) – The real part of the coefficients of the complex terms.
b_comp (array[j_complex]) – The imaginary part of the coefficients of the complex terms.
c_comp (array[j_complex]) – The real part of the exponents of the complex terms.
d_comp (array[j_complex]) – The imaginary part of the exponents of the complex terms.
x (array[n]) – The _sorted_ array of input coordinates.
y (array[n]) – The observations at
x
.diag (array[n]) – An array that should be added to the diagonal of the matrix. This often corresponds to measurement uncertainties and in that case,
diag
should be the measurement _variance_ (i.e. sigma^2).
-
log_determinant
(self: celerite.solver.CholeskySolver) → float¶ Get the log-determinant of the matrix defined by
compute
- Returns
The log-determinant of the matrix defined by
solver.CholeskySolver.compute()
.- Return type
-
predict
(self: celerite.solver.CholeskySolver, arg0: numpy.ndarray[float64[m, 1]], arg1: numpy.ndarray[float64[m, 1]]) → numpy.ndarray[float64[m, 1]]¶
-
solve
(self: celerite.solver.CholeskySolver, arg0: numpy.ndarray[float64[m, n]]) → numpy.ndarray[float64[m, n]]¶ Solve a linear system for the matrix defined in
compute
A previous call to
solver.CholeskySolver.compute()
defines a matrixA
and this method solves forx
in the matrix equationA.x = b
.- Parameters
b (array[n] or array[n, nrhs]) – The right hand side of the linear system.
- Returns
The solution of the linear system.
- Return type
array[n] or array[n, nrhs]
- Raises
ValueError – For mismatched dimensions.
-
-
class
celerite.solver.
CARMASolver
¶ A thin wrapper around the C++ CARMASolver class
This solver is parameterized following carma_pack: https://github.com/brandonckelly/carma_pack
- Parameters
log_sigma (float) – The log of the variance of the process.
arparams (array[p]) – The parameters of the autoregressive component.
maparams (array[q]) – The parameters of the moving average component.
-
get_celerite_coeffs
(self: celerite.solver.CARMASolver) → Tuple[numpy.ndarray[float64[m, 1]], numpy.ndarray[float64[m, 1]], numpy.ndarray[float64[m, 1]], numpy.ndarray[float64[m, 1]], numpy.ndarray[float64[m, 1]], numpy.ndarray[float64[m, 1]]]¶ Compute the coefficients of the celerite model for the given CARMA model
-
log_likelihood
(self: celerite.solver.CARMASolver, arg0: numpy.ndarray[float64[m, 1]], arg1: numpy.ndarray[float64[m, 1]], arg2: numpy.ndarray[float64[m, 1]]) → float¶ Compute the log likelihood using a Kalman filter
- Parameters
t (array[n]) – The input coordinates of the observations.
y (array[n]) – The observations.
yerr (array[n]) – The measurement uncertainties of the observations.