Computing occurrence rates
Note
This page is currently under construction!
With the Architectures of Exoplanetary Systems IV paper, SysSimPyPlots is now updated to include user-friendly functions for calculating planet occurrence rates from the existing models!
The multi-planet models enables us to compute both:
planet occurrence rates (the mean number of planets per star), and
fractions of stars with planets (the fraction of stars with at least one such planet).
These are computed from a large number of simulated catalogs from a given model by counting the number of planets within a given set of period, radius, and/or mass bounds provided by the user. Our functions described below also automatically calculates the uncertainties in the estimated quantities arising from uncertainties in the model parameters and stochastic noise during the simulation process.
A few simple examples
First, we need to load the simulated catalogs from a given model. We have a faster function for loading just the periods, radii, and masses from the physical catalogs, which are sufficient for computing simple occurrence rates:
# Load modules and catalogs as before:
from syssimpyplots.general import *
from syssimpyplots.load_sims import *
from syssimpyplots.compute_occurrence import *
load_dir = '/path/to/a/simulated/catalog/' # replace with your path!
# Load a simulated physical catalog:
sssp_per_sys_all = load_planets_periods_radii_masses_as_summary_stats_per_sys_many_fast(load_dir)
Note
The function load_planets_periods_radii_masses_as_summary_stats_per_sys_many_fast will load the first 100 simulated catalogs in the directory provided by load_dir, by default, but this can be changed by setting the optional argument runs=.
This is also much faster than loading each full catalog by repeatedly calling the function compute_summary_stats_from_cat_phys, as was done in the “Loading physical catalogs” section.
The output sssp_per_sys_all is a list of dictionaries (one for each simulated physical catalog), each containing the basic summary statistics (mainly things like P_all for all the orbital periods, radii_all for all the planet radii, and mass_all for all the planet masses).
To compute the occurrence rates, we can call the following function providing the list of simulated catalogs and simple bounds in period, radius, and/or mass. Some examples are given below:
\(R_p = 2\)-\(4 R_\oplus\) and \(P = 3\)-\(100\) days:
occurrence_all, fswp_all = compute_occurrence_rates_and_fswps_in_period_radius_mass_bounds_many_catalogs(sssp_per_sys_all, radius_bounds=(2,4), period_bounds=(3,100))
# Printed output:
No mass bounds given; assuming the full simulation range of (0.1,1000.0) Earth masses.
Occurrence rate = 0.422_{-0.030}^{+0.049}
f_swp = 0.295_{-0.019}^{+0.023}
\(M_p = 50\)-\(1000 M_\oplus\) and \(P = 3\)-\(11\) days:
occurrence_all, fswp_all = compute_occurrence_rates_and_fswps_in_period_radius_mass_bounds_many_catalogs(sssp_per_sys_all, mass_bounds=(50,1000), period_bounds=(3,11))
# Printed output:
No radius bounds given; assuming the full simulation range of (0.5,10.0) Earth radii.
Occurrence rate = 0.008_{-0.004}^{+0.006}
f_swp = 0.007_{-0.004}^{+0.005}
\(R_p = 0.8\)-\(1.2 R_\oplus\), \(M_p = 0.65\)-\(0.98 M_\oplus\), and \(P = 180\)-\(270\) days (e.g., Venus analogs):
occurrence_all, fswp_all = compute_occurrence_rates_and_fswps_in_period_radius_mass_bounds_many_catalogs(sssp_per_sys_all, radius_bounds=(0.8,1.2), mass_bounds=(0.65,0.98), period_bounds=(180,270))
# Printed output:
Occurrence rate = 0.006_{-0.002}^{+0.003}
f_swp = 0.006_{-0.002}^{+0.003}
The outputs of the function calls are an array of occurrence rates (occurrence_all) and an array of fractions of stars with planets (fswp_all), with each value in the arrays corresponding to a single simulated catalog. The function will also print the credible regions (central 68%) of each quantity.
Tip
For bounds that are not provided, the full simulation range will be assumed.