ExoSystem Object

Functions for fitting data on exoplanet systems.

class exosystem.ExoSystem(direc, init_planets='init_planets.txt', init_star='init_star.txt', init_lcs='init_lcs.txt', init_rv='init_rv.txt', init_ld='init_ld.txt', init_priors='init_priors.txt', init_ttvs='init_ttvs.txt')

A class for fitting exoplanet systems.

Initialize an ExoSystem object with a directory and init files. The init files must all be created using alfred.

Parameters:

direc (str) – The directory in which the init files, data, outputs, and plots of the ExoSystem are stored.
init_planets (str, optional) – Name of the planet parameter init file. Default is init_planets.txt.
init_star (str, optional) – Name of the star parameter init file. Default is init_star.txt.
init_lcs (str, optional) – Name of the light curve data init file. Default is init_lcs.txt.
init_rv (str, optional) – Name of the RV data init file. Default is init_rv.txt.
init_ld (str, optional) – Name of the limb darkening parameter init file. Default is init_ld.txt.
init_priors (str, optional) – Name of the priors init file. Default is init_priors.txt.
init_ttvs (str, optional) – Name of the TTV init file. Default is init_ttvs.txt.

fit(name, nburn, nrun, fit_transit, fit_rv, fit_star, nwalk=0, fit_ld=False, use_priors=False, rv_bkg_order=0, star_run=None, save_samples=False, sigma_clip=5, lc_supersample_size=600, show_plots=True, order_a=False, skip_state_check=False)

Fit the light curve, RV, and/or stellar data for this ExoSystem using MCMC.

Parameter optimization is done with scipy minimize before running MCMC. Additionally, when fitting transits to the light curves, will perform up to 10 rounds of sigma clipping to get rid of outliers in the lightcurve, stopping when the number of clipped points in a round is less than 10. Clipped points are masked out, and masks are saved in the Masks folder. During each round, masked points are shown in a plot before continuing. MCMC is performed with a joint likelihood function on all of the specified data.

When fitting stellar parameters simultaneously with planet parameters, runs a star-only fit first to get physical starting values for the stellar parameters.

Flattened and thinned MCMC chains are stored to pickle files ending in _res.p in the Output folder, as well as chains of derived parameters in a pickle file ending in _dres.p. A human readable table of median derived parameters, as well as their upper and lower uncertainties, is saved to the Results folder and printed at the end. Plots of the MCMC chains and corner plots are generated and saved to a folder in the Plots folder, along with plots of the best fit models to the data. Some run settings (e.g. fit_transit) are saved in a pickle file ending in _run_settings.p.

Parameters:

name (str) – Name of the fit. This name will be attached to all ouputs, including pickle files of data, human readable results tables, and the folder in Plots in which this run’s plots will be saved. This name is also used for loading results back in to be manipulated or plotted again.
nburn (int) – Number of burn-in steps for the MCMC. These steps are thrown out before saving the results and making plots. The burn-in allows the chains to settle into the maximum likelihood.
nrun (int) – Number of sampling steps for the MCMC. These steps are saved and used for results and making plots. They do not include the burn-in steps.
fit_transit (bool) – Whether or not to fit transits to the light curve data.
fit_rv (bool) – Whether or not to fit a model to the RV data.
fit_star (bool) – Whether or not to fit stellar parameters. Can be fit on their own, or if transit data is also being fit (with or without RV data as well). Cannot be fit with just RV data.
nwalk (int) – Number of walkers to use for the MCMC. This needs to be at least 2 times the number of free parameters. If nwalk is less than that value, or if nwalk isn’t provided, nwalk will be set to exactly 2 times the number of free parameters.
fit_ld (bool, optional) – Whether or not to fit quadratic limb darkening coefficients as free parameters. Will be overriden to False if fit_star is True, as limb darkening coefficients are generated at each step in that case. Default is False.
use_priors (bool, optional) – Whether or not to use user specified priors from the init_priors file specified during ExoSystem initialization. Default is False.
rv_bkg_order (int) – Polynomial order for the background trend fit to the RV data. Can be 0 (for a flat line), 1 (for a slope), or 2 (for a quadratic). Default is 0.
star_run (str) – Name of a previous fit run that included stellar parameters. If specified, provides starting values for the stellar parameters in this fit. If not specified, or it doesn’t exist, runs a star only fit before continuing to get physical starting parameters. Default is None.
save_samples (bool, optional) – Whether or not to save the full, unflattened, un-thinned MCMC chains to a pickle file. This can be handy if you expect to want to remove problematic walkers that wandered off from the final results. Most of the time this isn’t necessary, and these take up additional storage. Default is False.
sigma_clip (float) – Sets how strict the sigma clipping is. During each round, in each light curve, points farther than sigma_clip times the root median squared of the residual are masked out. If you see in the plots that points are being masked out that probably shouldn’t be, increase this parameter and run it again. To turn off sigma clipping set this parameter to 0. Default is 5.
lc_supersample_size (int) –
Sets the exposure time, in seconds, to supersample the transit light curves to. For each data set, the number of supersamples will be that data set’s exposure time divided by this value, rounded to the nearest integer.

For example, if a data set has an exposure time of 1800 seconds and this parameter is set to 600 seconds, each point in the data set will be modeled as the average of three points with an exposure time of 600 seconds.

Decreasing this parameter will increase the number of supersamples, potentially improving model accuracy, but at the cost of slowing down the fit.

Default is 600 seconds.
show_plots (bool, optional) – Whether or not to show plots at the end of the run. Plots are saved regardless. Default is True.
order_a (bool, optional) – Whether or not to put an order prior on the planetary semi-major axes in a transit fit with multiple transiting planets. If True, restricts the semi-major axis of each planet to be greater than those of planets with shorter orbital periods. If fitting the stellar parameters, overrides this to False, since the stellar mass is instead used to set semi-major axes. Default is False.
skip_state_check (bool, optional) – Passed to the emcee sampler. Whether or not to skip checking whether the initial parameters can fully explore the space. Only set to True if you keep getting initial state check errors after burn in. Default is False.

Return type:

None

initial_star_fit()

Runs an initial minimization on just the stellar parameters using a minimizer to get them in the ballpark before running MCMC. Starts at solar parameters. Called by fit.

Do not run this by itself.

run_sigma_clip(name, show_plots)

Not meant to be run on its own. Performs sigma clipping on the light curve data, with a threshold set by the sigma_clip parameter in fit. Initial fits are done using scipy minimize on the log likelihood function. Stops after 10 iterations or when the number of clipped points in an iteration is under 10 per data set. Saves plots to Plots/name/sigma_clip/.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save the plots to.
show_plot (bool) – Whether or not to show the plots.

initialize_chains()

Not meant to be run on its own. Initializes random positions for the MCMC chains around the initial best fits from optimization. Ensures that parameters stay within their bounds when applicable using truncated Gaussian distributions. Otherwise, just uses Gaussian distributions to add a little variance to the chains. Also builds a dict, self.parnames, along the way which keeps track of the index of each parameter in the pos array.

Returns:: The initial walker positions.
Return type:: ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]

run_sampler(pos, skip_state_check)

Not meant to be run on its own. Runs the MCMC sampler for burn-in and sample steps using emcee.

Parameters:

pos (ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]) – The initial positions of the chains generated by initialize_chains.
skip_state_check (bool) – Whether or not to skip the initial state check after the burn-in.

save_samples(name)

Saves the full sample chains out to a pickle file (not flattened or thinned). This is only automatically run if save_samples was set to True during fitting. Otherwise, can be run manually after a fit. This file takes up much more storage space than the normal thinned and flattened samples which are always saved. The pickle file is saved to Output/name_samples.p.

Parameters:: name (str) – Name of the run. Sets the name of the pickle file to name_samples.p.

continue_run(name, nrun, save_samples=False, show_plots=True, skip_state_check=False)

Continues running the MCMC sampler from where it left off, without a burn in. Remakes all results and plots, and saves to the provided name (does not have to be the same name as the previous run).

Can only be run if a run has already been performed with this ExoSystem object. Will not work on a loaded in run.

Parameters:

name (str) – Name of the fit. This name will be attached to all ouputs, including pickle files of data, human readable results tables, and the folder in Plots in which this run’s plots will be saved. This name is also used for loading results back in to be manipulated or plotted again.
nrun (int) – Number of sampling steps for the MCMC. These steps are saved and used for results and making plots. They do not include the burn-in steps.
save_samples (bool, optional) – Whether or not to save the full, unflattened, un-thinned MCMC chains to a pickle file. This can be handy if you expect to want to remove problematic walkers that wandered off from the final results. Most of the time this isn’t necessary, and these take up additional storage. Default is False.
show_plots (bool, optional) – Whether or not to show plots at the end of the run. Plots are saved regardless. Default is True.
skip_state_check (bool, optional) – Passed to the emcee sampler. Whether or not to skip checking whether the initial parameters can fully explore the space. Only set to True if you keep getting initial state check errors after burn in. Default is False.

flatten_chains(): Flattens and thins by a factor of 20 ExoSystem.samples, ExoSystem.log_likes, and ExoSystem.blobs (if it exists, stores linear regression periods and Tcs of TTV planets from each sample). Stores these as ExoSystem.flat_samples, ExoSystem.flat_log_likes, and ExoSystem.flat_blobs, respectively.

make_results(name)

Converts flat_samples and flat_log_likes into a dict that stores the chains of parameters that were directly fit, as well as the log likelihood of each sample. Saves these to the pickle file Output/name_res.p. Accessible through ExoSystem.chains.

Also creates a dict of flat chains of derived parameters (e.g. period from log(period), or planet equilibrium temperature). Saves these to the pickle file Output/name_dres.p. Accessible through ExoSystem.derived_chains.

Combines results and derived_results into a single human-readable ascii table which is saved to Results/name.txt. The astropy table can be directly accesed through ExoSystem.results.

Writes out specific run settings (fit_transit, fit_rv, fit_star, fit_ld, use_priors, lc_supersample_size, and fixed parameters if applicable) to the pickle file Output/name_run_settings.p. These are loaded back in with ExoSystem.load_results and are used to tell plotting functions what to plot.

Parameters:: name (str) – Name of the run. Sets the names of the output files from this function.

make_plots(name, show_plots=True)

Generates plots for the run, including chains plots, corner plots, and plots of best fit models. When plotting best fit models, also runs the corresponding function to generate those best fit models.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save the plots to.
show_plots (bool) – Whether or not to show plots. Plots are saved regardless. Default is True.

property results

Astropy table storing median and 1 sigma (16th and 84th percentile) values for fit parameters and derived parameters.

Use the chains object or derived chains object for flattened chains of these parameters, as well as a description of them in the docstrings.

property chains

Dict which stores the flattened chains of fit parameters. Keys are the parameter names.

Potential fit parameters:

log(P) x: Natural log of the orbital period, in days, for planet number x.
Tc x: Time of conjunction, in BJD-2457000, for planet number x.
ror x: Planet to star radius ratio, for planet number x.
log(a/rs) x: Natural log of the orbital semi-major axis to stellar radius ratio, for planet number x.
cos(i) x: Cosine of the orbital inclination, restricted between 0 and 90 degrees, for planet number x.
log(K) x: Natural log of the RV semi-amplitude, in m/s, for planet number x.
secw x: Square root of the eccentricity times the cosine of the argument of periastron, for planet number x.
sesw x: Square root of the eccentricity times the sine of the argument of periastron, for planet number x.
e x: Orbital eccentricity, for planet number x. Only fit if w is fixed.
w x: Orbital argument of periastron, in radians, for planet number x. Only fit if e is fixed.
fp x: Planet to star flux ratio, for planet number x.
TT x y: Transit time (when fitting for TTVs), in BJD-2457000, for planet number x transit number y.
F0 x: Baseline flux value, for transit data set x.
log(rho_gp) x: Natural log of the GP period, in days, for transit data set x.
log(sigma_gp) x: Natural log of the GP standard deviation, for transit data set x.
gamma: RV systemic velocity, in m/s.
gamma_dot: 1st derivative of the RV systemic velocity, in m/s/day. Relative to time of first RV measurement.
gamma_ddot: 2nd derivative of the RV systemic velocity, in m/s/day^2. Relative to time of first RV measurement.
rv_offset x: RV offset relative to first data set, in m/s, for RV data set x.
u1 x: Linear coefficient for quadratic limb darkening, for filter x.
u2 x: Non-linear coefficient for quadratic limb darkening, for filter x.
eep: Equivalent evolutionary phase, used to interpolate across MIST models for stellar fitting.
log10(age): Log base 10 of the stellar age, in yr.
feh: Stellar metallicity, in dex.
distance: Distance to the system, in pc.
AV: Visual band extinction to the system, in mags.

property derived_chains

Dict which stores flat chains of derived parameters. Keys are the parameter names.

Potential derived parameters:

P x: Orbital period, in days, for planet number x.
Rp x: Planetary radius, in earth radii, for planet number x.
a/rs x: Orbital semi-major axis to stellar radius ratio, for planet number x.
a x: Orbital semi-major axis, in AU, for planet number x.
i x: Orbital inclination, in degrees, for planet number x (restricted between 0 and 90 degrees).
e x: Orbital eccentricity, for planet number x.
w x: Orbital argument of periastron, in radians, for planet number x.
K x: RV semi-amplitude, in m/s, for planet number x.
Mp x: Planetary mass, in earth masses, for planet number x.
Mpsini x: Planetary mass times the sine of the orbital inclination, in earth masses, for planet number x (used when transit wasn’t fit).
rhop x: Planetary mean density, in g/cm^3, for planet number x.
teq x: Planetary equilibrium temperature (zero-albedo blackbody), in K, for planet number x.
sinc x: Insolation flux reaching the planet, in units of Earth insolation flux, for planet number x.
b x: Transit impact parameter, for planet number x.
depth x: Transit depth, in ppm, for planet number x.
dur x: Transit duration, in hours, for planet number x.
rhos x: Implied stellar density from transit of planet number x, in g/cm^3.
TSM x: Transmission spectroscopy metric, for planet number x.
rstar: Stellar radius, in solar radii.
mstar: Stellar mass, in solar masses.
mstar_init: Initial stellar mass at the zero-age main sequence point, in solar masses.
rhostar: Mean stellar density from stellar fitting, in g/cm^3.
Tstar: Stellar effective temperature, in K.
loggstar: log10 of the stellar surface gravity, in cm/s^2.
Lstar: Stellar bolometric luminosity, in solar luminosities.
Mbolstar: Stellar bolometric magnitude, in mags.

load_results(name)

Loads in the results from a previous run. Flattened chains will be accessable in ExoSystem.chains and ExoSystem.derived_chains.

Also loads the median results table

If the run includes a transit fit, loads in and applies masks to the light curves.

Attempts to load in best fit transit models and RV models (if they exist) to ExoSystem.fm and ExoSystem.rvm, respectively.

Attempts to load in SED best fit magnitudes and errors (if they exist) into ExoSystem.magfit and ExoSystem.magfiterr, respectively.

Attempts to load in unflattened, un-thinned MCMC chains if they exist. Paramater name map is saved to ExoSystem.parnames, chains are saved to ExoSystem.samples, log likelihood values are saved to ExoSystem.log_likes_full.

Also loads in some run settings from the run.

Parameters:: name (str) – Name of the previous run to load in.

load_results_tab(name)

Loads in the median results table from Results/name.txt as an astropy table accessible through ExoSystem.results.

Parameters:: name (str) – Name of the previous run to load in.

load_lcs_mod(name)

Loads in the best fit transit model of a previous run into ExoSystem.fm.

Parameters:: name (str) – Name of the previous run to load in.

load_rv_mod(name)

Loads in the best fit RV model of a previous run into ExoSystem.rvm.

Parameters:: name (str) – Name of the previous run to load in.

load_magfit(name)

Loads in the best fit stellar SED model of a previous run into ExoSystem.magfit and ExoSystem.magfiterr.

Parameters:: name (str) – Name of the previous run to load in.

load_samples(name)

Loads in unflattened, un-thinned MCMC chains. Paramater name map is saved to ExoSystem.parnames, chains are saved to ExoSystem.samples, log likelihood values are saved to ExoSystem.log_likes. If any TTVs were fit, linear regression periods and Tcs for each sample are saved to Exosystem.blobs.

Parameters:: name (str) – Name of the previous run to load in.

load_run_settings(name)

Loads in run settings from the run, specifically fit_transit, fit_rv, fit_star, fit_ld, use_priors, lc_supersample_size, and any parameters which were fixed. These are used to tell plotting functions what they need to plot.

Parameters:: name (str) – Name of the previous run to load in.

delete_run(name)

Searches for files from a run of the given name and deletes them.

Parameters:: name (str) – Name of the previous run to delete.

initialize_ttvs(name)

Loads in an Init_ttvs file if it exists, otherwise prompts the user to create one. Then, sets up initial variables for fitting ttvs. These are: - self.ttvs0: A dict with keys being planet numbers, and the values being a sorted array of all observed transit times. - self.ttvsectors: A dict with keys being the planet number followed by the transit number (1-indexed), and the values being the index of the transit data set which the transit occurs in. - self.ttvi: A dict with keys being planet numbers, and the values being an array of integers representing the number of periods occurring between the first observed transit and each transit being fit.

Parameters:: name (str) – Name of the Init_ttvs file.

plot_pl_chains(name, show_plot=True)

Plots the MCMC chains for the planetary parameters. Saves the plot as pl_chains.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_ttv_chains(name, show_plot=True)

Plots the MCMC chains for the planetary transit times of planets with TTVs. Saves the plot as tt_chains.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_det_chains(name, show_plot=True)

Plots the MCMC chains for the light curve detrending GP hyper parameters. Saves the plot as gp_chains.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_star_chains(name, show_plot=True)

Plots the MCMC chains for the stellar parameters. Saves the plot as star_chains.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_big_chains(name, show_plot=False)

Plots the MCMC chains for every fit parameter. Saves the plot as bigchains.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is False.

plot_pl_corner(name, show_plot=True)

Creates corner plots of the MCMC posteriors for each planet’s parameters. Saves the plots as p#_corner.png (where # is the planet number) to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_star_corner(name, show_plot=True)

Creates a corner plot of the MCMC posterior for the stellar parameters. Saves the plot as star_corner.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_big_corner(name, show_plot=False)

Creates a corner plot of the MCMC posterior of every parameter. Saves the plot as bigcorner.png to the folder for this run, set by the name parameter. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is False.

gen_lcs(name)

Generates model light curves from the best fit parameters. Saves these to the pickle file Output/name_lcm.p. These model light curves are accessible through ExoSystem.lc_mod. See the documentation on ExoSystem.lc_mod for details on the formatting of the object.

Parameters:: name (str) – Name of the run. Sets the name of the output pickle file to name_lcm.p.

property lc_mod

Stores the best fit model light curves.

fm is a dict, with a key for the phase folded models as well as for each light curve data set using their nicknames.

The ‘phase’ key contains a dict with ‘ttphase’, a 1d numpy array of phase-folded time points used for plotting the models. Also within ‘phase’ is a key for each filter/bandpass used (e.g. TESS or Kepler). These are their own dicts with the following:

fmphase: A 2d numpy array, where each row corresponds to a different transiting planet. The rows contain the best fit phase-folded model light curve for each planet, at the times of ttphase.
fmphase_err: A 3d numpy array of the 1 sigma uncertainties on fmphase. In the 0th axis, the 0 index is a 2d array corresponding to the lower 1 sigma error and the 1 index is a 2d array corresponding to the upper 1 sigma error. Each of these is formatted like fmphase, with each row corresponding to a different planet.
eclphase: The same as fmphase, but for the phase-folded secondary eclipses of any planets that were fit for this. Will be an empty array if no planets were fit for secondary eclipses.
eclphase_err: The same as fmphase_err, but for the phase-folded secondary eclipses of any planets that were fit for this. Will be an empty array if no planets were fit for secondary eclipses.

Each of the data set keys is a dict with the following:

fm: A 2d numpy array, where each row corresponds to a different transiting planet. The rows contain the best fit model light curve for each planet at the times of this data set.
fm_err: A 3d numpy array of the 1 sigma uncertainties on fm. In the 0th axis, the 0 index is a 2d array corresponding to the lower 1 sigma error and the 1 index is a 2d array corresponding to the upper 1 sigma error. Each of these is formatted like fm, with each row corresponding to a different planet.
gpf: A 1d numpy array of the GP model for detrending this light curve data set, using best fit hyper parameters.
gpf_err: A 2d numpy array of the 1 sigma uncertainties on gpf. In the 0th axis, the 0 index is a 1d array corresponding to the lower 1 sigma error and the 1 index is a 1d array corresponding to the upper 1 sigma error.

plot_full_lc(name, show_plot=True)

Plots full light curve fits for each data set. Saves to the folder in Plots set by name. Each data set plot is named transit_full_nickname.png, where nickname is the data set’s nickname. Shows the plot if show_plot is True.

If the data set was detrended, plots the GP prediction with best fit hyperparameters in the first panel. The flattened light curve and best fit transit models are plotted in the next panel. The residuals are plotted in the final panel.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_lc_phase(name, show_plot=True)

Plots phased light curve fits for each planet in each data set. Saves to the folder in Plots set by name. Each data set plot is named transits_nickname.png where nickname is the data set’s nickname. Shows the plot if show_plot is True.

The phased transit for each planet is plotted in a separate panel for each data set, with its residual below it.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

gen_rv(name)

Generates radial velocity models from the best fit parameters. Saves these to the pickle file Output/name_rvm.p. These RV models are accessible through ExoSystem.rv_mod. See the documentation on ExoSystem.rv_mod for details on the formatting of the object.

Parameters:: name (str) – Name of the run. Sets the name of the output pickle file to name_rvm.p.

property rv_mod

Stores the best fit radial velocity models.

rvm is a dict, with the following keys:

rvm: A 2d numpy array, where each row corresponds to a different RV planet. The rows contain the best fit RV model for each planet at the times of the RV data.
bkg: A 1d numpy array of the values of the background polynomial trend at the times of the RV data.
trplot: A 1d numpy array of more finely sampled times to make plotted models look smooth.
rvmplot: A 2d numpy array, where each row corresponds to a different RV planet. The rows contain the best fit RV model for each planet at the times of trplot, for smooth plotting.
bkgplot: A 1d numpy array of the values of the background polynomial trend at the times of trplot, for smooth plotting.
rvallplot: A 1d numpy array of the full best fit RV model at the times of trplot, for smooth plotting. It is the sum of all planet components and the background polynomial.
rvallplot_err: A 2d numpy array of the 1 sigma uncertainties on rvallplot. In the 0th axis, the 0 index is a 1d array corresponding to the lower 1 sigma error and the 1 index is a 1d array corresponding to the upper 1 sigma error.
trphase: A 1d numpy array of phase folded times.
rvmphase: A 2d numpy array, where each row corresponds to a different RV planet. The rows contain the best fit phase-folded RV model for each planet, at the times of trphase.
rvmphase_err: A 3d numpy array of the 1 sigma uncertainties on rvmphase. In the 0th axis, the 0 index is a 2d array corresponding to the lower 1 sigma error and the 1 index is a 2d array corresponding to the upper 1 sigma error. Each of these is formatted like rvmphase, with each row corresponding to a different planet.

plot_rv(name, show_plot=True)

Plots the best fit radial velocity model to the data. Saves to the folder in Plots set by name. The plot is named rv.png. Shows the plot if show_plot is True.

The top panel shows the full RV time series, with the combined best fit model. It also plots the contribution of each planet to the RV model, as well as the background polynomial. The residual is plotted below it.

The lower panels are the RV models of each individual planet, phased to its period. The data have the contributions of the other planets and the background polynomial subtracted out. The residual is plotted below each of these panels as well.

In each plot, the data are colored according to the data set they come from, and labeled with the data set’s nickname.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

gen_magfit(name)

Generates SED magnitudes and errors from the best fit parameters. Saves these to the pickle file Output/name_magfit.p. These model magnitudes are accessible through ExoSystem.magfit and ExoSystem.magfiterr. Also saves out the bands, so that the order can be matched.

Parameters:: name (str) – Name of the run. Sets the name of the output pickle file to name_magfit.p.

plot_sed_fit(name, show_plot=True)

Plots the magnitudes of the best fit stellar parameters against the data. Residuals are shown below. Saves to the folder in Plots set by name. The plot is named sed.png. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

plot_ttvs(name, show_plot=True)

Plots the best fit transit times versus the difference from a linear ephemeris for planets with TTVs. Saves to the folder in Plots set by name. The plot is named ttvs.png. Shows the plot if show_plot is True.

Parameters:

name (str) – Name of the run. Sets the folder in Plots to save this plot to.
show_plot (bool, optional) – Whether or not to show the plot. Default is True.

examine_chains(parname)

Allows you to examine the chains of a specific parameter. Hovering over a chain with your mouse will highlight it and display the index of the chain in the 1 axis (0-indexed) of ExoSystem.samples.

Useful for identifying the index of an outlier chain so that you can remove it using ExoSystem.remove_chains.

Must have the full samples of the run available as Exosystem.samples. Either this is from having just run a fit with this ExoSystem object, or from loading in the results of a previous run which had save_samples set to True.

Parameters:: parname (str) – The name of the parameter whose chains will be plotted. See the parameters that were fit in this run with ExoSystem.parnames.

remove_chains(name, idxs, save_samples=False, show_plots=True)

Removes chains of the specified indices from samples, log_likes, and the lists of derived periods and tcs in a ttv fit. Then, reruns ExoSystem.calc_gelman_rubin, ExoSystem.flatten_chains, ExoSystem.make_results, and ExoSystem.make_plots.

One way to find outlier chains is by running ExoSystem.examine_chains first, then hovering over the plot to get indices of the outliers.

Must have the full samples of the run available as Exosystem.samples. Either this is from having just run a fit with this ExoSystem object, or from loading in the results of a previous run which had save_samples set to True.

Parameters:

name (str) – Name of the run. Does not have to be the same as the previous run. Sets the folder in Plots to save the plots to, as well as the names of the output pickle files.
idxs (ArrayLike) – Indices of the chains to be removed.
save_samples (bool) – Whether or not to save the full, unflattened, un-thinned MCMC chains to a pickle file. This can be handy if you expect to want to remove problematic walkers that wandered off from the final results. Most of the time this isn’t necessary, and these take up additional storage. Default is False.
show_plots (bool) – Whether or not to show plots. Plots are saved regardless. Default is True.

burn_steps(name, num_steps, save_samples=False, show_plots=True)

Useful if the burn-in wasn’t long enough. Cuts the first num_steps steps from samples, log_likes, and the lists of derived periods and tcs in a ttv fit. Then, reruns ExoSystem.calc_gelman_rubin, ExoSystem.flatten_chains, ExoSystem.make_results, and ExoSystem.make_plots.

Must have the full samples of the run available as Exosystem.samples. Either this is from having just run a fit with this ExoSystem object, or from loading in the results of a previous run which had save_samples set to True.

Parameters:

name (str) – Name of the run. Does not have to be the same as the previous run. Sets the folder in Plots to save the plots to, as well as the names of the output pickle files.
num_steps (int) – Number of steps to cut from the beginning of the samples.
save_samples (bool) – Whether or not to save the full, unflattened, un-thinned MCMC chains to a pickle file. This can be handy if you expect to want to remove problematic walkers that wandered off from the final results. Most of the time this isn’t necessary, and these take up additional storage. Default is False.
show_plots (bool) – Whether or not to show plots. Plots are saved regardless. Default is True.

calc_rv_bic()

Calculates the Bayesian information criterion (BIC) of the best fit RV model. Can be used to compare different RV models with different numbers of planets, background trends, planet eccentricities, etc.

Not implemented for transit models, as the out of transit points tend to dominate and so it is less reliable.

calc_transit_err_scale()

Calculates the standard deviation of the residual of the best fit transit model and compares this to the median error on the light curve. Does this light curve by light curve. Prints out the ratio of this standard deviation to the median errors. Prints out a new error scaling value to put into init_lcs for each light curve, accounting for any error scaling already in effect.

It can be useful to rescale errors to if you believe they are under or overestimated, and a good place to start is looking at residual scatter. This is an alternative to using a jitter term in other fitting packages. Rescaling the errors can be very important when detrending light curves using a GP, as underestimated errors will cause the GP to try to detrend to each individual point.

calc_rv_err_scale()

Calculates the standard deviation of the residual of the best fit RV model and compares this to the median error on the RV data. Does this data set by data set. Prints out the ratio of this standard deviation to the median errors. Prints out a new error scaling value to put into init_rv for each data set, accounting for any error scaling already in effect.

It can be useful to rescale errors to if you believe they are under or overestimated, and a good place to start is looking at residual scatter. This is an alternative to using a jitter term in other fitting packages.

calc_gelman_rubin(): Calculates the Gelman-Rubin Statistic for each parameter in the MCMC fit. Generally, if the Gelman-Rubin Statistic is below 1.1 for all parameters, then the MCMC has converged. If it looks converged, but the statistic is above 1.1 for some parameters, try running for more steps.

rv_lomb_scargle(min_freq=None, max_freq=None, freq=None, plot_periods=False, use_residual=True)

Creates a Lomb-Scargle periodogram of the RV data. Also plots the 1% and 5% false alarm levels in red and orange, respectively. Generally, a signal is significant if it crosses above the 1% false alarm level. Can be useful for looking for additional periodic signals (like planets!) in the RV data.

Parameters:

min_freq (float) – If it is not None, sets the minimum frequency (in 1/days) to search in the periodogram. Default is None.
max_freq (float) – If it is not None, sets the maximum frequency (in 1/days) to search in the periodgram. Default is None.
freq (ArrayLike) – Manually sets the frequency grid to use in the periodogram. Units are 1/days. Default is None.
plot_periods (bool, optional) – Whether or not to plot the periodgram relative to period (True) or frequency (False). Default is False.
use_residual (bool, optional) – Whether or not to make this periodogram off of just the RV data (False), or the residual of the RV data with the best fit model (True). If True, uses whatever model is currently in ExoSystem.rv_mod. Default is True.

Returns:

The frequency grid, the periodogram power at each frequency, and the Lomb-Scargle object as a tuple. Allows you to calculate the maxmimum power frequency, or make other plots.

Return type:

tuple

exosystem.calc_m_from_k(p, k, e, inc, mstar)

Calculates the planet mass from the orbital period, RV semi-amplitude, eccentricity, inclination, and stellar mass. Does not make the assumption that the planet mass is much smaller than the stellar mass, so it is still valid even in brown dwarf-M dwarf scenarios.

Parameters:

p (ArrayLike) – Orbital period in years.
k (ArrayLike) – RV semi-amplitude in m/s.
e (ArrayLike) – Orbital eccentricity.
inc (ArrayLike) – Orbital inclination in radians.
mstar (ArrayLike) – Stellar mass in solar masses.

Returns:

The planet mass in Earth masses.

Return type:

ArrayLike

exosystem.lightCurve(par, t, ld, expt, ss, eclipse=False, fp=0, rstar=0)

Generates a model transit light curve using batman.

Parameters:

par (ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]) – Array of transit parameters. Parameters must be specifically in this order. - Orbital period in same units as t. - Time of conjunction in same units as t. - Planet to star radius ratio. - Orbital semi-major axis to stellar radius ratio. - Orbital inclination in degrees. - Orbital eccentricity. Should not be over 0.9 or batman has issues. - Argument of periastron in degrees.
t (ArrayLike) – Times at which to calculate the transit model.
ld (list) – A list only containing the linear and non-linear quadratic limb darkening coefficients, in that order.
expt (float) – The exposure time of each point. Used for supersampling.
ss (int) – Number of supersamples per exposure. For each point, this many points are generated across the whole exposure time, then averaged.

Returns:

The model transit light curve at each of the input times, with a baseline of 0 flux.

Return type:

ArrayLike

exosystem.calc_t_sec(p, tc, e, w, a, r)

Calculates the time of secondary eclipse center, including light travel delay.

Parameters:

p (ArrayLike) – Orbital period in days.
tc (ArrayLike) – Time of conjunction in days.
e (ArrayLike) – Orbital eccentricity.
w (ArrayLike) – Argument of periastron in radians.
a (ArrayLike) – Semi-major axis to stellar radius ratio.
r (ArrayLike) – Stellar radius in solar radii.

Returns:

The time of secondary eclipse center.

Return type:

ArrayLike

exosystem.rvModel(par, t)

Generates model RVs using _rvModel. This is compiled with cython to run faster for large amounts of data.

Parameters:

par (ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]) – Array of RV parameters. Parameters must be specifically in this order. - Orbital period in same units as t. - Time of conjunction in same units as t. - RV semi-amplitude. Returns model RV points in same units as this. - Orbital eccentricity. Should not be over 0.9 or batman has issues. - Argument of periastron in radians.
t (ArrayLike) – Times at which to calculate the RV model.

Returns:

The model RVs at each of the input times, centered at 0.

Return type:

ArrayLike

exosystem.get_transit_params(par, i, ar=None)

Generates the specific list of transit parameters needed for the lightCurve function from the parameter dict used in fitting.

Parameters:

par (dict) – Dictionary of fit parameters.
i (int) – Planet number.
ar (ArrayLike) – Direct input for the semi-major axis to stellar radius ratio if it is not being fit for (during stellar fitting). Default is None, in which case the value is pulled from par.

Returns:

Array of transit parameters that can be input directly into lightCurve.

Return type:

ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]

exosystem.get_ttv_params(par, i, ttvi, ar=None)

Generates the specific list of transit parameters needed for the lightCurve function from the parameter dict used in fitting, in the case when the planet is being fit for TTVs. Calculates the period and time of conjunction using linear regression on the transit times.

Parameters:

par (dict) – Dictionary of fit parameters.
i (int) – Planet number.
ttvi (ndarray[tuple[Any, ...], dtype[int64]]) – Array of integers representing the number of periods that have occurred between the first observed transit and each observed transit that is being fit.
ar (ArrayLike) – Direct input for the semi-major axis to stellar radius ratio if it is not being fit for (during stellar fitting). Default is None, in which case the value is pulled from par.

Returns:

Array of transit parameters that can be input directly into lightCurve.

Return type:

ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]

exosystem.get_rv_params(par, i, p=None, tc=None)

Generates the specific list of RV parameters needed for the rvModel function from the parameter dict used in fitting.

Parameters:

par (dict) – Dictionary of fit parameters.
i (int) – Planet number.
p (ArrayLike) – Direct input for the orbital period if it is not being fit for (if the planet is being fit for ttvs). Default is None, in which case the value is pulled from par.
tc (ArrayLike) – Direct input for the time of conjunction if it is not being fit for (if the planet is being fit for ttvs). Default is None, in which case the value is pulled from par.

Returns:

Array of RV parameters that can be input directly into rvModel.

Return type:

ndarray[tuple[Any, ...], dtype[TypeVar(_ScalarT, bound= generic)]]

exosystem.set_gp_params(rho, sigma, t, ferr, gp)

Updates a GP object with new parameters, then computes the factorization of the covariance matrix on the data.

Parameters:

rho (float) – The period of the GP in days.
sigma (float) – The standard deviation of the GP.
t (ArrayLike) – The times of the data points.
ferr (ArrayLike) – The flux uncertainties of the data points.
gp (GaussianProcess) – The GP to update and then return.

Returns:

Returns the input GP, with updated parameters and a recomputed factorization.

Return type:

GaussianProcess

exosystem.lnNorm(data, model, err)

Computes the natural log of a normal distribution pdf. Used for the log likelihood.

Parameters:

data (ArrayLike) – The observed values.
model (ArrayLike) – The predicted values.
err (ArrayLike) – The uncertainties in the observed values.

Returns:

The log probability of each input point.

Return type:

ArrayLike

exosystem.log_like(par_in, exs)

The log likelihood function for fitting.

Checks certain parameters to make sure they are in bounds. Calculates likelihoods from any priors. Calculates log likelihood and for stellar fitting using the isochrones MIST star model. Generates transit and RV models, then calculates log likelihood relative to data using lnNorm.

Parameters:

par_in (dict) – Dictionary of fit parameters.
exs (ExoSystem) – The exosystem object currently being used to fit. Stores all of the options, data sets, etc.

Returns:

A tuple of the log likelihood, an array of periods of any planets that are fit for ttvs, and an array of conjunction times of any planets that are fit for ttvs. The former is the primary return for running MCMC with emcee. The latter two are saved as “blobs” in the emcee sampler object, and are used to track the period and time of conjunction at each step rather than recalculating the linear regression after the fact, since these are not directly fit parameters. If any parameter is out of bounds or any portion of the log likelihood is nan, returns negative infinity and empty lists.

Return type:

tuple[float, ArrayLike, ArrayLike]

exosystem.log_like_staronly(par_in, exs)

The log likelihood function for fitting stellar parameters only. Used for the initial parameter estimation before a main fit.

Parameters:

par_in (dict) – Dictionary of fit parameters.
exs (ExoSystem) – The exosystem object currently being used to fit. Stores all of the options, data sets, etc.

Returns:

The log likelihood. Returns negative infinity if any parameter is out of bounds or any of the log likelihood functions return nan.

Return type:

float