gyrointerp package

gyrointerp.gyro_posterior module

Main drivers:: gyro_age_posterior

gyro_age_posterior_list

gyro_age_posterior_mcmc

gyrointerp.gyro_posterior.gyro_age_posterior(Prot, Teff, Prot_err=None, Teff_err=None, age_grid=np.linspace(0, 3000, 500), interp_method='pchip_m67', bounds_error='4gyrlimit', n=None, N_grid='default', age_scale='default', popn_parameters='default', verbose=False)[source]

Given a single star’s rotation period and effective temperature, and (optionally) their uncertainties, what is the posterior probability distribution for its age?

The answer returned by this code assumes that the gyrointerp.models.slow_sequence_residual model holds, which is the probability distribution described in BPH23 for the distribution of rotation periods at any given age and temperature.

If Prot_err and Teff_err are not specified, they are assumed to be 1% relative and 100 K, respectively. Spectroscopic temperature are acceptable, though the preferred effective temperature scale is implemented in gyrointerp.teff, in the given_dr2_BpmRp_AV_get_Teff_Curtis2020 function. This requires an accurate estimate for the reddening. Whatever your effective temperature scale, it should ideally be compared against that in Appendix A of Curtis+2020.

Parameters:

Prot (int or float) – Rotation period in units of days.
Prot_err (int or float) – Rotation period uncertainty in units of days.
Teff – (int or float): Effective temperature in units of degrees Kelvin. Must be between 3800 and 6200 K.
Teff_err (int or float) – Effective temperature uncertainty in units of degrees Kelvin.
age_grid (np.ndarray) – Grid over which the age posterior is evaluated; units here and throughout are fixed to be megayears (10^6 years). A fine choice is 500 points uniformly distributed between 0 and 3000 Myr: np.linspace(0, 3000, 500), assuming that the default choices of bounds_error == '4gyrlimit' and interp_method == 'pchip_m67' are being used.
interp_method (str) – How will you interpolate between the polynomial fits to the reference open clusters? “pchip_m67” is the suggested default method, which uses Piecewise Cubic Hermite Interpolating Polynomials (PCHIP) to interpolate over not only 0.08-2.6 Gyr, but also sets the gradient in Prot vs Time in the 1-2.6 Gyr interval based on the observations of M67 from Barnes+2016 and Dungee+2022. This yields an evolution of the rotation period envelope that is smooth and by design fits the cluster data from the age of alpha-Per through M67. Other available interpolation methods include “skumanich_vary_n”, “alt”, “diff”, “skumanich_fix_n_0.XX”, “1d_linear”, “1d_slinear”, “1d_quadratic”, and “1d_pchip”, some of which are described in Appendix A of BPH23. Unless you know what you are doing, “pchip_m67” is recommended.
bounds_error (str) – How will you extrapolate at <0.08 Gyr and >2.6 Gyr? Available options are “nan”, “limit” or “4gyrlimit”. If “limit”, then extrapolate by returning the fixed limiting rotation period at the oldest or youngest cluster given in models.slow_sequence.reference_model_ids. If “4gyrlimit” (the suggested default), extrapolate out to 4 Gyr by also including M67. Past 4Gyr, use the same behavior as “limit”. If one is interested in obtaining unbiased ages near the recommended 2.6 Gyr limit of this code, use “4gyrlimit”, otherwise “limit” will overestimate the probability density beyond 2.6 Gyr. Finally, if “nan”, ages above or below the minimum reference age return NaNs.
n (None, int, or float) – Power-law index analogous to the Skumanich braking index, but different in detail (read the source code to see how). This is used only if interp_method == "alt" or interp_method == "diff", neither of which is recommended for most users. Default is None.
N_grid (str or int) – The number of grid points in effective temperature and residual-period over which the integration is performed to evaluate the posterior (Equation 1 of BPH23). “default” sets it to N_grid = (Prot_grid_range)/Prot_err, where “Prot_grid_range” is set to 20 days, the range of the grid used in the integration. This default ensures convergence, because numerical tests show convergence down to ~0.7x smaller grid sizes. If an integer is passed, will use that instead. For most users, “default” is best.
age_scale (str) – “default”, “1sigmaolder”, or “1sigmayounger”. Shifts the entire age scale appropriately, based on the user’s beliefs about what ages of reference clusters are correct. The scale is as described in the systematic uncertainty tests of BPH23, and defined in gyrointerp.age_scale.
popn_parameters (str) – “default”, or (dict) containing the population-level free parameters. Keys of “a0”, “a1”, “k0”, “k1”, “y_g”, “l_hidden”, and “k_hidden” must all be specified. Wrapped by gyro_age_posterior_mcmc, for users who wish to do the population-level hyperparameter sampling described by BPH23.

Returns:

p_ages

Numpy array containing posterior probabilities at each point of the age_grid. Values are NaN if Teff outside of [3800, 6200] K.

Return type:

np.ndarray

gyrointerp.gyro_posterior.gyro_age_posterior_list(cache_id, Prots, Teffs, Prot_errs=None, Teff_errs=None, star_ids=None, age_grid=np.linspace(0, 3000, 500), N_grid='default', bounds_error='4gyrlimit', interp_method='pchip_m67', nworkers=None)[source]

Given rotation periods and effective temperatures for many stars, run them in parallel. This is a thin wrapper to gyro_age_posterior assuming default parameters.

Parameters:

cache_id (str) – The output posteriors will be cached to ~/.gyrointerp_cache/{cache_id} (required).
Prots (np.ndarray) – 1-D array of rotation periods
Teffs (np.ndarray) – 1-D array of temperatures, same length as Prots
Prot_errs (np.ndarray) – 1-D array of rotation period uncertainties. If None, assumes 1% relative uncertainties by default.
Teff_errs (np.ndarray) – 1-D array of effective temperature uncertainties. If None, assumes 100K by default.
star_ids (np.ndarray of strings) – Arbitrary strings naming your stars; optional. For example, if you give “TIC1234567”, posteriors will be written to CSV files with a pattern matching TIC1234567_ProtXX.XXXX_TeffYYYY.Y_limitgrid_defaultparameters.csv. If None, then the identifier is omitted.
nworkers (int or None) – Number of workers to thread over. By default, will be taken to be all available CPU cores.

Returns:

List of paths to all available output posteriors at ~/.gyrointerp_cache/{cache_id}. If you re-use your cache_id, this means you will get more than you asked for!

gyrointerp.gyro_posterior.gyro_age_posterior_mcmc(Prot, Teff, Prot_err=None, Teff_err=None, age_grid=np.linspace(0, 3000, 500), verbose=False, bounds_error='4gyrlimit', interp_method='pchip_m67', N_grid='default', n=None, age_scale='default', N_pop_samples=512, N_post_samples=10000, cachedir=None)[source]

Given the rotation period and effective temperature of a single star, sample over the population-level hyperparameters a1/a0, ln k0, ln k1, and y_g to determine the posterior probability distribution of the age. These are the dotted lines in the upper panel of Fig3 in BPH23.

Parallelization is done over the hyperparameters. However, the computational cost for a given star is about 1000x that of running the best-fit hyperparameters, as implemented in gyro_age_posterior. Use of this function is therefore generally not recommended, unless you have an understood need for doing things this way.

Arguments are as in gyro_age_posterior, but with four additions:

Parameters:

cache_id (str) – The output posteriors will be cached to ~/.gyrointerp_cache/{cache_id} (required).
N_pop_samples (int) – The number of draws from the posteriors for the aforementioned parameters to average over.
N_post_samples (int) – For each of the above draws, the number of draws from the resulting age posterior to cache before concatenating them all together.
cachedir (str) – Path to directory where individual posteriors will be cached (for every star, N_post_samples files will be written here!). It is highly recommended that you specify this.

Returns:

p_ages

Numpy array containing posterior probabilities at each point of the age_grid. Values are NaN if Teff outside of [3800, 6200] K.

Return type:

np.ndarray

gyrointerp.models module

Functions to fit rotation versus effective temperature sequences, or to quickly return the results of those fits (including their interpolations!)

Contents:: reference_cluster_slow_sequence

slow_sequence

slow_sequence_residual
Helper functions:: teff_zams

g_lineardecay

gyrointerp.models.g_lineardecay(age, bounds_error='4gyrlimit', y_g=1 / 2)[source]: Function g(t) from BPH23 defining the linear rate at which the uniform component amplitude from models.slow_sequence_residual decreases with age. Unity at <=120 Myr, decreasing linearly to y_g by 300 Myr (eg. 1/2, 1/4, 1/6). Decreases linearly thereafter, and once it reaches zero, it stays at zero.

gyrointerp.models.reference_cluster_slow_sequence(Teff, model_id, poly_order=7, verbose=True)[source]

Given a set of temperatures, get the rotation periods implied by a polymial fit to the Prot-Teff slow sequence of a particular cluster between 3800 and 6200 K.

Parameters:

Teff (np.ndarray / float / list-like iterable) – Effective temperature in Kelvin. Curtis+2020 Gaia DR2 BP-RP scale, or spectroscopic effective temperatures, preferred above all other options.
model_id (str) – String identifying the desired reference cluster. Can be any of ['α Per', 'Pleiades', 'Blanco-1', 'Psc-Eri', 'NGC-3532', 'Group-X', 'Praesepe', 'NGC-6811', '120-Myr', '300-Myr', 'NGC-6819', 'Ruprecht-147', '2.6-Gyr', 'M67'], where ‘120-Myr’ will concatenate of Pleiades, Blanco-1, and Psc-Eri into one polynomial fit, and ‘300-Myr’ will concatenate NGC-3532 and Group-X.
poly_order (int) – Integer order of the polynomial fit.

Returns:

Prot_model

Numpy array containing rotation periods for each requested temperature.

Return type:

np.ndarray

gyrointerp.models.slow_sequence(Teff, age, poly_order=7, reference_model_ids=['α Per', '120-Myr', '300-Myr', 'Praesepe', 'NGC-6811', '2.6-Gyr'], reference_ages=[80, 120, 300, 670, 1000, 2600], verbose=True, bounds_error='4gyrlimit', interp_method='pchip_m67', n=None)[source]

Given an age and a set of temperatures, return the implied slow sequence rotation periods, as derived from interpolation using the reference clusters with known ages. This function is the “mean gyrochronal model”, i.e., it assumes slow sequence evolution.

Parameters:

age (int or float) – An integer or float corresponding to the age for which we want a rotation period. Units: Myr (=10^6 years).
Teff (float or iterable of floats) – Effective temperature(s) of the sample to be dated. Units: Kelvin. Must be between 3800 and 6200 K.
reference_model_ids (list) – This list can include any of ['α Per', 'Pleiades', 'Blanco-1', 'Psc-Eri', 'NGC-3532', 'Group-X', 'Praesepe', 'NGC-6811', '120-Myr', '300-Myr', '2.6-Gyr', 'NGC-6819', 'Ruprecht-147', 'M67'] The default is set as described in the manuscript, to enable gyro-age derivations between 0.08-2.6 Gyr. Note that “120-Myr” and “300-Myr” are concenations of the relevant clusters.
reference_ages (iterable of floats) – Ages (units of Myr) corresponding to reference_model_ids.
verbose (bool) – True or False to choose whether to print error messages. Default is False
interp_method (str) –
How will you interpolate between the polynomial fits to the reference open clusters? “pchip_m67” is the suggested default method, which uses Piecewise Cubic Hermite Interpolating Polynomials (PCHIP) to interpolate over not only 0.8-2.6 Gyr, but also sets the gradient in Prot vs Time in the 1-2.6 Gyr interval based on the observations of M67 from Barnes+2016 and Dungee+2022. This yields an evolution of the rotation period envelope that is smooth and by design fits the cluster data from the age of alpha-Per through M67. Other available interpolation methods include “skumanich_vary_n”, “alt”, “diff”, “skumanich_fix_n_0.XX”, “1d_linear”, “1d_slinear”, “1d_quadratic”, and “1d_pchip”, some of which are described in Appendix A of BPH23. Unless you know what you are doing, “pchip_m67” is recommended.
bounds_error (str) – How will you extrapolate at <0.08 Gyr and >2.6 Gyr? Available options are “nan”, “limit” or “4gyrlimit”. If “limit”, then extrapolate by returning the fixed limiting rotation period at the oldest or youngest cluster given in models.slow_sequence.reference_model_ids. If “4gyrlimit” (the suggested default), extrapolate out to 4 Gyr by also including M67. Past 4Gyr, use the same behavior as “limit”. If one is interested in obtaining unbiased ages near the recommended 2.6 Gyr limit of this code, use “4gyrlimit”. Finally, if “nan”, ages above or below the minimum reference age return NaNs.
n (None, int, or float) – Power-law index analogous to the Skumanich braking index, but different in detail (see the implementation). This is used only if interp_method == "alt" or interp_method == "diff", neither of which is recommended for most users. Default is None.

Output:

np.ndarray : Prot_slow_sequence

Numpy array containing period estimate for the slow sequence at the given effective temperatures and age, in units of days.

gyrointerp.models.slow_sequence_residual(age, y_grid=np.linspace(-14, 6, 1000), teff_grid=np.linspace(3800, 6200, 1001), poly_order=7, n=None, reference_model_ids=['α Per', '120-Myr', '300-Myr', 'Praesepe', 'NGC-6811', '2.6-Gyr'], reference_ages=[80, 120, 300, 670, 1000, 2600], popn_parameters='default', verbose=True, bounds_error='4gyrlimit', interp_method='pchip_m67')[source]

Given an effective temperature and an age, return the 2-D distribution of residuals around and underneath the slow sequence, sampled onto grids of y_grid X teff_grid, where y_grid is the residual of (rotation period data - mean gyrochronal model). This is Equation 7 of BPH23.

The two components of the residual are:

a gaussian in y_grid , with an age-varying cutoff in Teff (teff_zams), imposed as a logistic taper.

an age-varying and Teff-varying uniform distribution, multiplied by the inverse of the gaussian’s taper (but with independent scale length), and then truncated to ensure that stars rotate faster than zero days, and to ensure that we model only the fast sequence. This uniform distribution is also tapered by a logistic function at the “slow” end to yield a smoother transition to the gaussian.

Parameters:

teff_grid (np.ndarray) – As in models.slow_sequence.
y_grid (np.ndarray) – A grid over the residual of (rotation period - mean gyrochronal model).
poly_order (int) – As in models.slow_sequence.
reference_model_ids (list) – As in models.slow_sequence.
reference_ages (list) – As in models.slow_sequence.
interp_method (str) – As in models.slow_sequence.
bounds_error (str) – As in models.slow_sequence.
popn_parameters (str or dict) – “default”, or a dict containing the population-level free parameters. Keys of “a0”, “a1”, “k0”, “k1”, “y_g”, “l_hidden”, and “k_hidden” must all be specified.

Returns:

resid_y_Teff

2d array with dimension (N_y_grid x N_teff_grid), containing the probability distribution of the residual over the dimensions of y and Teff.

Return type:

np.ndarray

gyrointerp.models.teff_zams(age, bounds_error='limit')[source]: Physics-informed MIST effective temperature for the effective temperature a star has when it arrives on the ZAMS. Changes over 80 to 1000 Myr (the floor beyond 1000 Myr is well below the 3800 K cool limit of BPH23). For a test, see /tests/plot_teff_cuts.py.

gyrointerp.getters module

These functions return the rotation periods and effective temperatures for single members of benchmark open clusters. These are the functions used to acquire and clean the data behind Figure 1 of BPH23. Generally, most users will want to just download the data from the paper, rather than using these functions. An additional dependency if you do want to use them is the cdips module, best installed via a setup.py install from https://github.com/lgbouma/cdips.

Get all available cluster data:: _get_cluster_Prot_Teff_data
Get and clean individual cluster data:: get_alphaPer

get_Pleiades

get_Blanco1

get_PscEri

get_NGC3532

get_GroupX

get_Praesepe_Rampalli_2021

get_NGC6811

get_NGC6819

get_Ruprecht147

gyrointerp.getters.get_Blanco1(overwrite=0)[source]: Returns the Gillen+2020 (Table 1) rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_GroupX(overwrite=0)[source]: Returns the Messina+2022 rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_NGC3532(overwrite=0)[source]: Returns the Fritzewski+2021 rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_NGC6811(overwrite=0)[source]: Returns the Curtis+2019 (Table 1) rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_NGC6819(overwrite=0)[source]: Return NGC-6819 (Meibom+2015), as processed by Curtis+2020 (Table 5).

gyrointerp.getters.get_Pleiades(overwrite=0)[source]: Returns the Rebull+2016 (Table 2) rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_Praesepe_Rampalli_2021(overwrite=0)[source]: Returns the Rampalli+2021 Praesepe rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_PscEri(overwrite=0)[source]: Returns the Curtis+2019 (Table 2) rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_Ruprecht147(overwrite=0)[source]: Return Curtis+2020 (Table 1) rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ the usual Gaia DR3 columns.

Note

The “flag_possible_binary” flag in this case is weaker than in the other comparison clusters, because there aren’t very many good stars.

gyrointerp.getters.get_alphaPer(overwrite=0)[source]: Returns the Boyle & Bouma 2023 rotation period dataframe with keys: “Prot”, “Teff_Curtis20”, “flag_benchmark_period” (for gyro calibration), _plus_ “flag_possible_binary”, “flag_ruwe_outlier”, “flag_camd_outlier”, “flag_rverror_outlier”, _plus_ the usual Gaia DR3 columns.

gyrointerp.getters.get_alphaPer_construct(overwrite=0)[source]: Constructor for getters.get_alphaPer

gyrointerp.age_scale module

This module contains a dictionary defining the default cluster age scale, as well as the assumed +1sigma and -1sigma uncertainties on those ages. These data are the same as Table 1 in BPH23.

gyrointerp.extinctionpriors module

This module defines a dictionary of mean extinction A_V values adopted for the reference open clusters when calibrating the gyrochronal model.

gyrointerp.helpers module

This module contains reusable helper functions. The most generally useful one will be get_summary_statistics.

gyrointerp.helpers.get_population_hyperparameter_posterior_samples()[source]

Access the posterior samples described in section 3.5 of BPH23. (These are generated by drivers.run_emcee_fit_gyro_model).

The returned numpy array is samples in the following parameters: a1/a0, y_g, logk0, logk1, log_f.

The notation follows Sections 3.3-3.5 of BPH23.

gyrointerp.helpers.get_summary_statistics(age_grid, age_post, N=int(100000.0))[source]

Given an age posterior probability density, age_post, over a grid over ages, age_grid, determine summary statistics for the posterior (its median, mean, +/-1 and 2-sigma intervals, etc). Do this by sampling N times from the posterior, with replacement, while weighting by the probability.

Parameters:

age_grid (np.ndarray) – Array-like of ages, in units of megayears. For instance, the default age_grid in gyro_posterior.gyro_age_posterior is np.linspace(0, 3000, 500).
age_post (np.ndarray) – Posterior probability distribution for ages; length should match age_grid. The posterior probabilities returned by``gyro_posterior.gyro_age_posterior`` and gyro_posterior.gyro_age_posterior_list are examples that would work. Generally, this helper function works for any grid and probability distribution.

Returns:

summary_statistics

Dictionary containing keys and values for median, mean, peak (mode), +/-1sigma, +/-2sigma, +/-3sigma, and +/-1sigmapct. The units of all values are megayears, except for +/-1sigmapct, which is the relative +/-1-sigma uncertainty normalized by the median of the posterior, and is dimensionless.

Return type:

dict

gyrointerp.helpers.given_dr2_get_dr3_dataframes(dr2_source_ids, runid_dr2, runid_dr3, overwrite=False)[source]

gyrointerp.helpers.left_merge(df0, df1, col0, col1)[source]

gyrointerp.helpers.prepend_colstr(colstr, df)[source]

gyrointerp.plotting module

gyrointerp.teff module

This module contains functions for calculating photometric effective temperatures.

Contents:: given_dr2_BpmRp_AV_get_Teff_Curtis2020

gyrointerp.teff.given_dr2_BpmRp_AV_get_Teff_Curtis2020(dr2_BpmRp, A_V)[source]

Empirical color-temperature relation from Appendix A of Curtis+2020. Visible in their Figure 11. Coefficients from their Table 4.

This relation was constructed using benchmark stars from Brewer+2016a, Boyajian+2012, and Mann+2015 (which should also be cited).

It is calibrated in the range 0.55 < (BP-RP)0 < 3.25, and has a scatter of about +/- 50 K.

It’s important that the passed BP-RP colors are from Gaia DR2. For FGK stars in a few test clusters (eg. NGC-3532), the typical offset is +0.02 mag and color dependent. For late K dwarf and M-dwarfs, it flips, and gets down to -0.1 mag at BP-RP of >2 (SpType>M1V). This translates to a >100 K systematic error if you use the wrong Gaia data release.

See https://www.cosmos.esa.int/web/gaia/edr3-passbands for a description of why exactly the Gaia passbands changed between reductions.

Parameters:

dr2_BpmRp (np.ndarray) – observed Gaia DR2 BP-RP array.
A_V (float) – mean reddening.

Returns:

array of effective temperatures.

Return type:

Teff (np.ndarray)