API documentation¶
Orbits¶
-
class
exoplanet.orbits.
KeplerianOrbit
(period=None, a=None, t0=0.0, incl=None, b=None, ecc=None, omega=None, m_planet=0.0, m_star=None, r_star=None, rho_star=None, m_planet_units=None, rho_star_units=None, model=None, contact_points_kwargs=None, **kwargs)¶ A system of bodies on Keplerian orbits around a common central
Given the input parameters, the values of all other parameters will be computed so a
KeplerianOrbit
instance will always have attributes for each argument. Note that the units of the computed attributes will all be in the standard units of this class (R_sun
,M_sun
, anddays
) except forrho_star
which will be ing / cm^3
.There are only specific combinations of input parameters that can be used:
- First, either
period
ora
must be given. If values are given for both parameters, then neitherm_star
orrho_star
can be defined because the stellar density implied by each planet will be computed inrho_star
. - Only one of
incl
andb
can be given. - If a value is given for
ecc
thenomega
must also be given. - If no stellar parameters are given, the central body is assumed to be
the sun. If only
rho_star
is defined, the radius of the central is assumed to be1 * R_sun
. Otherwise, at most two ofm_star
,r_star
, andrho_star
can be defined.
Parameters: - period – The orbital periods of the bodies in days.
- a – The semimajor axes of the orbits in
R_sun
. - t0 – The time of a reference transit for each orbits in days.
- incl – The inclinations of the orbits in radians.
- b – The impact parameters of the orbits.
- ecc – The eccentricities of the orbits. Must be
0 <= ecc < 1
. - omega – The arguments of periastron for the orbits in radians.
- m_planet – The masses of the planets in units of
m_planet_units
. - m_star – The mass of the star in
M_sun
. - r_star – The radius of the star in
R_sun
. - rho_star – The density of the star in units of
rho_star_units
. - m_planet_units – An
astropy.units
compatible unit object giving the units of the planet masses. If not given, the default isM_sun
. - rho_star_units – An
astropy.units
compatible unit object giving the units of the stellar density. If not given, the default isg / cm^3
.
-
get_planet_position
(t)¶ The planets’ positions in the barycentric frame
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_planet_velocity
(t)¶ Get the planets’ velocity vector
Parameters: t – The times where the velocity should be evaluated. Returns: The components of the velocity vector at t
in units ofM_sun/day
.
-
get_radial_velocity
(t, K=None, output_units=None)¶ Get the radial velocity of the star
Parameters: - t – The times where the radial velocity should be evaluated.
- K (Optional) – The semi-amplitudes of the orbits. If provided, the
m_planet
andincl
parameters will be ignored and this amplitude will be used instead. - output_units (Optional) – An AstroPy velocity unit. If not given,
the output will be evaluated in
m/s
. This is ignored if a value is given forK
.
Returns: The reflex radial velocity evaluated at
t
in units ofoutput_units
. For multiple planets, this will have one row for each planet.
-
get_relative_position
(t)¶ The planets’ positions relative to the star
Note
This treats each planet independently and does not take the other planets into account when computing the position of the star. This is fine as long as the planet masses are small.
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_star_position
(t)¶ The star’s position in the barycentric frame
Note
If there are multiple planets in the system, this will return one column per planet with each planet’s contribution to the motion. The star’s full position can be computed by summing over the last axis.
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_star_velocity
(t)¶ Get the star’s velocity vector
Note
For a system with multiple planets, this will return one column per planet with the contributions from each planet. The total velocity can be found by summing along the last axis.
Parameters: t – The times where the velocity should be evaluated. Returns: The components of the velocity vector at t
in units ofM_sun/day
.
- First, either
-
class
exoplanet.orbits.
TTVOrbit
(*args, **kwargs)¶ A generalization of a Keplerian orbit with transit timing variations
Only one of the arguments
ttvs
ortransit_times
can be given and the other will be computed from the one that was provided.Parameters: - ttvs – A list (with on entry for each planet) of “O-C” vectors for each transit of each planet in units of days. “O-C” means the difference between the observed transit time and the transit time expected for a regular periodic orbit.
- transit_times – A list (with on entry for each planet) of transit times
for each transit of each planet in units of days. These times will
be used to compute the implied (least squares)
period
andt0
so these parameters cannot also be given.
-
get_planet_position
(t)¶ The planets’ positions in the barycentric frame
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_planet_velocity
(t)¶ Get the planets’ velocity vector
Parameters: t – The times where the velocity should be evaluated. Returns: The components of the velocity vector at t
in units ofM_sun/day
.
-
get_radial_velocity
(t, K=None, output_units=None)¶ Get the radial velocity of the star
Parameters: - t – The times where the radial velocity should be evaluated.
- K (Optional) – The semi-amplitudes of the orbits. If provided, the
m_planet
andincl
parameters will be ignored and this amplitude will be used instead. - output_units (Optional) – An AstroPy velocity unit. If not given,
the output will be evaluated in
m/s
. This is ignored if a value is given forK
.
Returns: The reflex radial velocity evaluated at
t
in units ofoutput_units
. For multiple planets, this will have one row for each planet.
-
get_relative_position
(t)¶ The planets’ positions relative to the star
Note
This treats each planet independently and does not take the other planets into account when computing the position of the star. This is fine as long as the planet masses are small.
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_star_position
(t)¶ The star’s position in the barycentric frame
Note
If there are multiple planets in the system, this will return one column per planet with each planet’s contribution to the motion. The star’s full position can be computed by summing over the last axis.
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
-
get_star_velocity
(t)¶ Get the star’s velocity vector
Note
For a system with multiple planets, this will return one column per planet with the contributions from each planet. The total velocity can be found by summing along the last axis.
Parameters: t – The times where the velocity should be evaluated. Returns: The components of the velocity vector at t
in units ofM_sun/day
.
-
class
exoplanet.orbits.
SimpleTransitOrbit
(period=None, t0=0.0, b=0.0, duration=None, r_star=1.0)¶ An orbit representing a set of planets transiting a common central
This orbit is parameterized by the observables of a transiting system, period, phase, duration, and impact parameter.
Parameters: - period – The orbital period of the planets in days.
- t0 – The midpoint time of a reference transit for each planet in days.
- b – The impact parameters of the orbits.
- duration – The durations of the transits in days.
- r_star – The radius of the star in
R_sun
.
-
get_relative_position
(t)¶ The planets’ positions relative to the star
Parameters: t – The times where the position should be evaluated. Returns: The components of the position vector at t
in units ofR_sun
.
Light curve models¶
-
class
exoplanet.
StarryLightCurve
(u, r_star=1.0, model=None)¶ A limb darkened light curve computed using starry
Parameters: u (vector) – A vector of limb darkening coefficients. -
get_light_curve
(orbit=None, r=None, t=None, texp=None, oversample=7, order=0, use_in_transit=True)¶ Get the light curve for an orbit at a set of times
Parameters: - orbit – An object with a
get_relative_position
method that takes a tensor of times and returns a list of Cartesian coordinates of a set of bodies relative to the central source. This method should return three tensors (one for each coordinate dimension) and each tensor should have the shapeappend(t.shape, r.shape)
orappend(t.shape, oversample, r.shape)
whentexp
is given. The first two coordinate dimensions are treated as being in the plane of the sky and the third coordinate is the line of sight with positive values pointing away from the observer. For an example, take a look atorbits.KeplerianOrbit
. - r (tensor) – The radius of the transiting body in the same units as
r_star
. This should have a shape that is consistent with the coordinates returned byorbit
. In general, this means that it should probably be a scalar or a vector with one entry for each body inorbit
. - t (tensor) – The times where the light curve should be evaluated.
- texp (Optional[tensor]) – The exposure time of each observation.
This can be a scalar or a tensor with the same shape as
t
. Iftexp
is provided,t
is assumed to indicate the timestamp at the middle of an exposure of lengthtexp
. - oversample (Optional[int]) – The number of function evaluations to use when numerically integrating the exposure time.
- order (Optional[int]) – The order of the numerical integration
scheme. This must be one of the following:
0
for a centered Riemann sum (equivalent to the “resampling” procedure suggested by Kipping 2010),1
for the trapezoid rule, or2
for Simpson’s rule. - use_in_transit (Optional[bool]) – If
True
, the model will only be evaluated for the data points expected to be in transit as computed using thein_transit
method onorbit
.
- orbit – An object with a
-
Scalable Gaussian processes¶
-
class
exoplanet.gp.
GP
(kernel, x, diag, J=-1, model=None)¶
-
class
exoplanet.gp.terms.
Term
(**kwargs)¶ The abstract base “term” that is the superclass of all other terms
Subclasses should overload the
terms.Term.get_real_coefficients()
andterms.Term.get_complex_coefficients()
methods.
-
class
exoplanet.gp.terms.
RealTerm
(**kwargs)¶ The simplest celerite term
This term has the form
\[k(\tau) = a_j\,e^{-c_j\,\tau}\]with the parameters
a
andc
.Strictly speaking, for a sum of terms, the parameter
a
could be allowed to go negative but since it is somewhat subtle to ensure positive definiteness, we recommend keeping both parameters strictly positive. Advanced users can build a custom term that has negative coefficients but care should be taken to ensure positivity.Parameters: - a or log_a (tensor) – The amplitude of the term.
- c or log_c (tensor) – The exponent of the term.
-
class
exoplanet.gp.terms.
ComplexTerm
(**kwargs)¶ A general celerite term
This term has the form
\[k(\tau) = \frac{1}{2}\,\left[(a_j + b_j)\,e^{-(c_j+d_j)\,\tau} + (a_j - b_j)\,e^{-(c_j-d_j)\,\tau}\right]\]with the parameters
a
,b
,c
, andd
.This term will only correspond to a positive definite kernel (on its own) if \(a_j\,c_j \ge b_j\,d_j\).
Parameters: - a or log_a (tensor) – The real part of amplitude.
- b or log_b (tensor) – The imaginary part of amplitude.
- c or log_c (tensor) – The real part of the exponent.
- d or log_d (tensor) – The imaginary part of exponent.
-
class
exoplanet.gp.terms.
SHOTerm
(*args, **kwargs)¶ A term representing a stochastically-driven, damped harmonic oscillator
The PSD of this term is
\[S(\omega) = \sqrt{\frac{2}{\pi}} \frac{S_0\,\omega_0^4} {(\omega^2-{\omega_0}^2)^2 + {\omega_0}^2\,\omega^2/Q^2}\]with the parameters
S0
,Q
, andw0
.Parameters: - S0 or log_S0 (tensor) – The parameter \(S_0\).
- Q or log_Q (tensor) – The parameter \(Q\).
- w0 or log_w0 (tensor) – The parameter \(\omega_0\).
-
class
exoplanet.gp.terms.
Matern32Term
(**kwargs)¶ A term that approximates a Matern-3/2 function
This term is defined as
\[k(\tau) = \sigma^2\,\left[ \left(1+1/\epsilon\right)\,e^{-(1-\epsilon)\sqrt{3}\,\tau/\rho} \left(1-1/\epsilon\right)\,e^{-(1+\epsilon)\sqrt{3}\,\tau/\rho} \right]\]with the parameters
sigma
andrho
. The parametereps
controls the quality of the approximation since, in the limit \(\epsilon \to 0\) this becomes the Matern-3/2 function\[\lim_{\epsilon \to 0} k(\tau) = \sigma^2\,\left(1+ \frac{\sqrt{3}\,\tau}{\rho}\right)\, \exp\left(-\frac{\sqrt{3}\,\tau}{\rho}\right)\]Parameters: - sigma or log_sigma (tensor) – The parameter \(\sigma\).
- rho or log_rho (tensor) – The parameter \(\rho\).
- eps (Optional[float]) – The value of the parameter \(\epsilon\). (default: 0.01)
-
class
exoplanet.gp.terms.
RotationTerm
(**kwargs)¶ A mixture of two SHO terms that can be used to model stellar rotation
This term has two modes in Fourier space: one at
period
and one at0.5 * period
. This can be a good descriptive model for a wide range of stochastic variability in stellar time series from rotation to pulsations.Parameters: - amp or log_amp (tensor) – The amplitude of the variability.
- period or log_period (tensor) – The primary period of variability.
- Q0 or log_Q0 (tensor) – The quality factor (or really the quality factor minus one half) for the secondary oscillation.
- deltaQ or log_deltaQ (tensor) – The difference between the quality factors
of the first and the second modes. This parameterization (if
deltaQ > 0
) ensures that the primary mode alway has higher quality. - mix – The fractional amplitude of the secondary mode compared to the
primary. This should probably always be
0 < mix < 1
.
Estimators¶
-
exoplanet.
estimate_semi_amplitude
(periods, x, y, yerr=None, t0s=None)¶ Estimate the RV semi-amplitudes for planets in an RV series
Parameters: - periods – The periods of the planets. Assumed to be in
days
if not an AstroPy Quantity. - x – The observation times. Assumed to be in
days
if not an AstroPy Quantity. - y – The radial velocities. Assumed to be in
m/s
if not an AstroPy Quantity. - yerr (Optional) – The uncertainty on
y
. - t0s (Optional) – The time of a reference transit for each planet, if known.
Returns: An estimate of the semi-amplitude of each planet in units of
m/s
.- periods – The periods of the planets. Assumed to be in
-
exoplanet.
estimate_minimum_mass
(periods, x, y, yerr=None, t0s=None, m_star=1)¶ Estimate the minimum mass(es) for planets in an RV series
Parameters: - periods – The periods of the planets. Assumed to be in
days
if not an AstroPy Quantity. - x – The observation times. Assumed to be in
days
if not an AstroPy Quantity. - y – The radial velocities. Assumed to be in
m/s
if not an AstroPy Quantity. - yerr (Optional) – The uncertainty on
y
. - t0s (Optional) – The time of a reference transit for each planet, if known.
- m_star (Optional) – The mass of the star. Assumed to be in
M_sun
if not an AstroPy Quantity.
Returns: An estimate of the minimum mass of each planet as an AstroPy Quantity with units of
M_jupiter
.- periods – The periods of the planets. Assumed to be in
-
exoplanet.
lomb_scargle_estimator
(x, y, yerr=None, min_period=None, max_period=None, filter_period=None, max_peaks=2, **kwargs)¶ Estimate period of a time series using the periodogram
Parameters: - x (ndarray[N]) – The times of the observations
- y (ndarray[N]) – The observations at times
x
- yerr (Optional[ndarray[N]]) – The uncertainties on
y
- min_period (Optional[float]) – The minimum period to consider
- max_period (Optional[float]) – The maximum period to consider
- filter_period (Optional[float]) – If given, use a high-pass filter to down-weight period longer than this
- max_peaks (Optional[int]) – The maximum number of peaks to return (default: 2)
Returns: A dictionary with the computed
periodogram
and the parameters for up tomax_peaks
peaks in the periodogram.
-
exoplanet.
autocorr_estimator
(x, y, yerr=None, min_period=None, max_period=None, oversample=2.0, smooth=2.0, max_peaks=10)¶ Estimate the period of a time series using the autocorrelation function
Note
The signal is interpolated onto a uniform grid in time so that the autocorrelation function can be computed.
Parameters: - x (ndarray[N]) – The times of the observations
- y (ndarray[N]) – The observations at times
x
- yerr (Optional[ndarray[N]]) – The uncertainties on
y
- min_period (Optional[float]) – The minimum period to consider
- max_period (Optional[float]) – The maximum period to consider
- oversample (Optional[float]) – When interpolating, oversample the times by this factor (default: 2.0)
- smooth (Optional[float]) – Smooth the autocorrelation function by this factor times the minimum period (default: 2.0)
- max_peaks (Optional[int]) – The maximum number of peaks to identify in the autocorrelation function (default: 10)
Returns: A dictionary with the computed autocorrelation function and the estimated period. For compatibility with the
lomb_scargle_estimator()
, the period is returned as a list with the keypeaks
.
Distributions¶
-
class
exoplanet.distributions.
UnitVector
(*args, **kwargs)¶ A vector where the sum of squares is fixed to unity
For a multidimensional shape, the normalization is performed along the last dimension.
-
class
exoplanet.distributions.
Angle
(*args, **kwargs)¶ An angle constrained to be in the range -pi to pi
The actual sampling is performed in the two dimensional vector space
(sin(theta), cos(theta))
so that the sampler doesn’t see a discontinuity at pi.
-
class
exoplanet.distributions.
QuadLimbDark
(*args, **kwargs)¶ An uninformative prior for quadratic limb darkening parameters
This is an implementation of the Kipping (2013) reparameterization of the two-parameter limb darkening model to allow for efficient and uninformative sampling.
-
class
exoplanet.distributions.
RadiusImpact
(*args, **kwargs)¶ The Espinoza (2018) distribution over radius and impact parameter
This is an implementation of Espinoza (2018) The first axis of the shape of the parameter should be exactly 2. The radius ratio will be in the zeroth entry in the first dimension and the impact parameter will be in the first.
Parameters: - min_radius – The minimum allowed radius.
- max_radius – The maximum allowed radius.
-
exoplanet.distributions.
get_joint_radius_impact
(name='', N_planets=None, min_radius=0, max_radius=1, testval_r=None, testval_b=None, **kwargs)¶ Get the joint distribution over radius and impact parameter
This uses the Espinoza (2018) parameterization of the distribution (see
distributions.RadiusImpact
for more details).Parameters: - name (Optional[str]) – A prefix that is added to all distribution names
used in this parameterization. For example, if
name
isparam_
, vars will be added to the PyMC3 model with namesparam_rb
(for the joint distribution),param_b
, andparam_r
. - N_planets (Optional[int]) – The number of planets. If not provided, it
will be inferred from the
testval_*
parameters or assumed to be 1. - min_radius (Optional[float]) – The minimum allowed radius.
- max_radius (Optional[float]) – The maximum allowed radius.
- testval_r (Optional[float or array]) – An initial guess for the radius
parameter. This should be a
float
or an array withN_planets
entries. - testval_b (Optional[float or array]) – An initial guess for the impact
parameter. This should be a
float
or an array withN_planets
entries.
Returns: Two
pymc3.Deterministic
variables for the planet radius and impact parameter.- name (Optional[str]) – A prefix that is added to all distribution names
used in this parameterization. For example, if
Utilities¶
-
exoplanet.
eval_in_model
(var, point=None, return_func=False, model=None, **kwargs)¶ Evaluate a Theano tensor or PyMC3 variable in a PyMC3 model
This method builds a Theano function for evaluating a node in the graph given the required parameters. This will also cache the compiled Theano function in the current
pymc3.Model
to reduce the overhead of calling this function many times.Parameters: - var – The variable or tensor to evaluate.
- point (Optional) – A
dict
of input parameter values. This can bemodel.test_point
(default), the result ofpymc3.find_MAP
, a point in apymc3.MultiTrace
or any other representation of the input parameters. - return_func (Optional[bool]) – If
False
(default), return the evaluated variable. IfTrue
, return the result, the Theano function and the list of arguments for that function.
Returns: Depending on
return_func
, either the value ofvar
atpoint
, or this value, the Theano function, and the input arguments.
-
exoplanet.
get_samples_from_trace
(trace, size=1)¶ Generate random samples from a PyMC3 MultiTrace
Parameters: - trace – The
MultiTrace
. - size – The number of samples to generate.
- trace – The
-
class
exoplanet.
PyMC3Sampler
(start=75, finish=50, window=25, dense=True)¶ A sampling wrapper for PyMC3 with support for dense mass matrices
This schedule is based on the method used by as described in Section 34.2 of the Stan Manual.
Parameters: - start (int) – The number of steps to run as an initial burn-in to find the typical set.
- window (int) – The length of the first mass matrix tuning phase. Subsequent tuning windows will be powers of two times this size.
- finish (int) – The number of tuning steps to run to learn the step size after tuning the mass matrix.
- dense (bool) – Fit for the off-diagonal elements of the mass matrix.
-
extend_tune
(steps, start=None, step_kwargs=None, trace=None, step=None, **kwargs)¶ Extend the tuning phase by a given number of steps
After running the sampling, the mass matrix is re-estimated based on this run.
Parameters: steps (int) – The number of steps to run.
-
get_step_for_trace
(trace=None, model=None, regular_window=0, regular_variance=0.001, **kwargs)¶ Get a PyMC3 NUTS step tuned for a given burn-in trace
Parameters: - trace – The
MultiTrace
output from a previous run ofpymc3.sample
. - regular_window – The weight (in units of number of steps) to use when regularizing the mass matrix estimate.
- regular_variance – The amplitude of the regularization for the mass
matrix. This will be added to the diagonal of the covariance
matrix with weight given by
regular_window
.
- trace – The
-
sample
(trace=None, step=None, start=None, step_kwargs=None, **kwargs)¶ Run the production sampling using the tuned mass matrix
This is a light wrapper around
pymc3.sample
and any arguments used there (for exampledraws
) can be used as input to this method too.
-
tune
(tune=1000, start=None, step_kwargs=None, **kwargs)¶ Run the full tuning run for the mass matrix
This will run
start
steps of warmup followed by chains with exponentially increasing chains to tune the mass matrix.Parameters: tune (int) – The total number of steps to run.
-
warmup
(start=None, step_kwargs=None, **kwargs)¶ Run an initial warmup phase to find the typical set