# Towards precise stellar ages: combining isochrone fitting with empirical
gyrochronology

Ruth Angus11affiliation: American Museum of Natural History, Central Park
West, Manhattan, NY, USA affiliationmark:  22affiliation: Center for
Computational Astrophysics, Flatiron Institute, 162 5th Avenue, Manhattan, NY,
USA affiliationmark:  , Timothy D. Morton33affiliation: Department of
Astronomy, University of Florida, Gainesville, FL, USA affiliationmark:
22affiliation: Center for Computational Astrophysics, Flatiron Institute, 162
5th Avenue, Manhattan, NY, USA affiliationmark:  , Daniel Foreman-
Mackey22affiliation: Center for Computational Astrophysics, Flatiron
Institute, 162 5th Avenue, Manhattan, NY, USA , Jennifer van
Saders44affiliation: Institute for Astronomy, University of Hawai’i at Mānoa,
Honolulu, HI, USA , Jason Curtis55affiliation: Department of Astronomy,
Columbia University, Manhattan, NY, USA , Stephen R. Kane66affiliation:
Department of Earth and Planetary Sciences, University of California,
Riverside, CA 92521, USA , Megan Bedell22affiliation: Center for Computational
Astrophysics, Flatiron Institute, 162 5th Avenue, Manhattan, NY, USA , Rocio
Kiman77affiliation: Department of physics, CUNY Graduate Center, City
University of New York, Manhattan, NY, USA affiliationmark:  11affiliation:
American Museum of Natural History, Central Park West, Manhattan, NY, USA ,
David W. Hogg88affiliation: Center for Cosmology and Particle Physics, New
York University, Manhattan, NY, USA affiliationmark:  22affiliation: Center
for Computational Astrophysics, Flatiron Institute, 162 5th Avenue, Manhattan,
NY, USA affiliationmark:  99affiliation: Center Data Science, New York
University, Manhattan, NY, USA affiliationmark:  1010affiliation: Max-Planck
Institute for Astronomy, Königstuhl, Heidelberg, Germany , John
Brewer1111affiliation: Yale center for astronomy and astrophysics, Yale
University, New Haven, CT, USA

###### Abstract

We present a new age-dating technique that combines gyrochronology with
isochrone fitting to infer ages for FGKM main-sequence and subgiant field
stars. Gyrochronology and isochrone fitting are each capable of providing
relatively precise ages for field stars in certain areas of the Hertzsprung-
Russell diagram: gyrochronology works optimally for cool main-sequence stars,
and isochrone fitting can provide precise ages for stars near the main-
sequence turnoff. Combined, these two age-dating techniques can provide
precise and accurate ages for a broader range of stellar masses and
evolutionary stages than either method used in isolation. We demonstrate that
the position of a star on the Hertzsprung-Russell or color-magnitude diagram
can be combined with its rotation period to infer a precise age via both
isochrone fitting and gyrochronology simultaneously. We show that
incorporating rotation periods with 5% uncertainties into stellar evolution
models improves age precision for FGK stars on the main sequence, and can, on
average, provide age estimates up to three times more precise than isochrone
fitting alone. In addition, we provide a new gyrochronology relation,
calibrated to the Praesepe cluster and the Sun, that includes a variance model
to capture the rotational behavior of stars whose rotation periods do not
lengthen with the square-root of time, and parts of the Hertzsprung-Russell
diagram where gyrochronology has not been calibrated. This publication is
accompanied by an open source Python package, stardate, for inferring the ages
of main-sequence and subgiant FGKM stars from rotation periods, spectroscopic
parameters and/or apparent magnitudes and parallaxes.

##  1 Introduction

Age is the most difficult stellar property to measure, and the difficulty of
age-dating is particularly acute for low mass (GKM) stars on the main sequence
(MS). Using conventional dating methods, uncertainties on the ages of these
stars can be as large as the age of the Universe. GKM dwarfs are difficult to
age-date because most of their physical and observable properties do not
change rapidly. This is represented in the spacing of isochrones on a
Hertzsprung-Russell diagram (HRD) or color-magnitude diagram (CMD). On the MS,
isochrones are tightly spaced and, even with very precise measurements of
effective temperature and luminosity, the position of a MS star on the HRD may
be consistent with range of isochrones spanning several billion years (see
Soderblom, 2010, for a review of stellar ages). At the main-sequence turnoff
however, isochrones are spread further apart, so that sufficiently precisely
measured temperatures and luminosities can yield ages that are extremely
precise (with minimum statistical uncertainties on the order of 5-10%) (e.g.
Pont and Eyer, 2004). The classical method for measuring stellar ages is
isochrone placement, or isochrone fitting, where surface gravity changes
resulting from fusion in the core (usually observed via luminosity, $L$, and
effective temperature, $T_{\mathrm{eff}}$, or absolute magnitude and colour)
are compared with a set of models that trace stellar evolution across the HRD,
or CMD. CMD/HRD position has been thoroughly mapped with physical models, and
can be used to calculate relatively accurate (but not necessarily precise)
ages, barring some systematic variations between different models, (e.g. Yi et
al., 2001; Dotter et al., 2008; Dotter, 2016). On the MS itself, there is
little differentiation between stars of different ages in the $L$ and
$T_{\mathrm{eff}}$ plane, so ages tend to be very imprecise. The method of
inferring a star’s age from its rotation period, called ‘gyrochronology’, is
much better suited for measuring ages on the MS because MS stars spin down
relatively rapidly.

Magnetic braking in MS stars was first observed by Skumanich (1972) who,
studying young clusters and the Sun, found that the rotation periods of Solar-
type stars decay with the square-root of time. It has since been established
that the rotation period of a MS star depends, to first order, only on its age
and effective temperature or color (e.g. Barnes, 2003). The convenient
characteristic of stars that allows their ages to be inferred from their
current rotation periods and independently of their primordial ones, comes
from the steep dependence of spin-down rate on rotation period (Kawaler,
1989). Stars spinning with high angular velocity will experience a much
greater angular momentum loss rate than slowly spinning stars and for this
reason, no matter the initial rotation period, Solar type stars will have the
same rotation period after around the age of the Hyades, 500-700 million years
(Irwin and Bouvier, 2009; Gallet and Bouvier, 2015). After this time, the age
of a star can be inferred, to first order, from its dust-corrected color (e.g.
B-V or $G_{BP}-G_{RP}$) and current rotation period alone (See Epstein and
Pinsonneault, 2014, for an analysis of how initial conditions effect
gyrochronal ages).

The relation between age, rotation period and mass has been studied in detail,
and several different models have been developed to capture the rotational
evolution of Sun-like stars. Some of these models are theoretical and based on
physical processes; modeling angular momentum loss as a function of stellar
properties as well as the properties of the magnetic field and stellar wind
(e.g. Kawaler, 1988, 1989; Pinsonneault et al., 1989; van Saders and
Pinsonneault, 2013; Matt et al., 2015; van Saders et al., 2016). Other models
are empirical and capture the behavior of stars from a purely observational
standpoint, using simple functional forms that can reproduce the data (e.g.
Barnes, 2003, 2007; Mamajek and Hillenbrand, 2008; Angus et al., 2015). Both
types of model, theoretical and empirical, must be calibrated using
observations. Old calibrators are especially important because new evidence
suggests that rotational evolution goes through a transition at old age or,
more specifically, at a large Rossby number, $Ro$ (the ratio of rotation
period to the convective overturn timescale). For example, stars shown to be
old from Kepler asteroseismic data rotate more rapidly than expected given
their age (Angus et al., 2015; van Saders et al., 2016). A new physically
motivated gyrochronology model, capable of reproducing these data, was
recently introduced (van Saders et al., 2016). It relaxes magnetic braking at
a critical Rossby number of around 2, approximately the Solar value. This
model predicts that, after stellar rotation periods lengthen enough to move
stars across this $Ro$ threshold, stars conserve angular momentum and maintain
a nearly constant rotation period from then until they evolve off the MS.

The gyrochronology models that capture post $Ro$-threshold, rotational
evolution (van Saders et al., 2016) are the current state-of-the-art in
rotation dating. These models can be computed over a grid of stellar
parameters, and interpolated-over to predict the age of a star. The process of
measuring the age of a field star with these models is similar to inferring an
age using any set of isochrones, with the difference being that rotation
period is an additional observable dimension. Ages calculated using these
models are therefore likely to be much more precise than using rotation-free
isochrones since rotation period provides an additional anchor-point for the
age of a star. We present here a complementary method that combines isochrones
with an empirical gyrochronology model. The methodology is related to the van
Saders et al. (2016) model in that both use a combination of rotation periods
and other observable properties that track stellar evolution on the HRD in
concert. The main difference is that the gyrochronology model used here is an
entirely empirically calibrated one, as opposed to a physically derived one.

One major advantage of using a physically motivated gyrochronology model is
the ability to rely on physics to interpolate or extrapolate over parts of
parameter space with sparse data coverage. However, rotational spin-down is a
complex process that is not yet fully understood and currently no physical
model can accurately reproduce all the data available. For this reason, even
physically motivated gyrochronology models cannot always be used to reliably
extrapolate into unexplored parameter space. Physical models, when calibrated
to data, can provide insight into the physics of stars, however, if accurate
and precise prediction of stellar properties is desired, empirical models can
have advantages over physical ones. For example, the data may reveal complex
trends that cannot be reproduced with our current understanding of the
physical processes involved, but may be captured by more flexible data-driven
models. In addition, it is relatively straightforward to build an element of
stochasticity into empirical models, i.e. to allow for and incorporate
outliers or noisy trends. This is particularly important for stellar spin down
because rotation periods can be affected by additional confounding variables
which are not always observed (having a binary or planetary companion, for
example). A further advantage of empirical models is that inference is more
tractable: it can be extremely fast to fit them to data.

In this work we calibrated a new empirical gyrochronology relation, fit to the
Praesepe open cluster and the Sun, in Gaia $G_{BP}-G_{RP}$ color. Gaia $G$,
$G_{BP}$ and $G_{RP}$ apparent magnitudes are now the most abundant
photometric measurements, available for more than a billion stars (Gaia
Collaboration et al., 2018). In fact, an important point of context for this
work is the new availability of data relevant to gyrochronology and stellar
ages. Gaia now provides broad-band photometry and parallaxes for over a
billion stars, and Kepler, K2 and TESS are providing rotation periods for
hundreds of thousands of stars. Gyrochronology is becoming one of the most
readily available age-dating methods, so continuing to improve gyrochronology
relations and methods is important. Like any other gyrochronology relation,
this new Praesepe-based model does not perfectly reproduce all the observed
data and some simple modifications could make significant improvements, for
example, by including a mixture model to account for outliers and binaries,
and by removing the period-age, period-color separability to account for
different period-color shapes seen in clusters of different ages. We leave
these improvements for a future project and, for now, test this new Praesepe-
based model, which is built into an open source Python package for
gyrochronology called stardate. stardate provides the framework for
simultaneous gyrochronology and isochrone fitting, and, because stardate is
modular, it would be straightforward to update the gyrochronology relation in
the future.

This paper is laid out as follows. In section 2 we describe our new age-dating
model and its implementation, in section 3 we test this model on simulated
stars and cluster stars, and in section 4 we discuss the implications of these
tests and future pathways for development. Throughout this paper we use the
term ‘observables’ to refer to the following observed properties of a star:
$T_{\mathrm{eff}}$, log(g), observed bulk metallicity ($[\mathrm{Fe/H}]$),
parallax ($\bar{\pi}$), photometric colors in different passbands (${\bf
m_{x}}=[m_{J},m_{H},m_{K},m_{B},m_{V},m_{G},m_{GBP},m_{GRP}...]$, etc) and
rotation period ($P_{\mathrm{rot}}$). The term ‘parameters’ refers to the
physical properties of that star: age ($t$), equivalent evolutionary phase
(EEP), model metallicity ($\mathrm{M/H}$), distance ($D$) and V-band
extinction ($A_{V}$). These are the properties that generate the observables.

##  2 Method

###  2.1 A new empirical gyrochronology relation

We fit a broken power law to the 650 Myr Praesepe cluster in order to
calibrate a gyrochronology relation that takes advantage of new data available
from the K2 and Gaia surveys. This relation captures the detailed shape of
Praesepe’s rotation period-color relation, and is calibrated to Gaia
$G_{BP}-G_{RP}$ color. Praesepe is the ideal calibration cluster because it
has the largest number of members with precisely measured rotation periods,
over a large range of colors (spanning spectral types A0 through M6, Rebull et
al., 2017), of any open cluster. It is relatively tightly clustered on the sky
and many of its members were targeted in a single K2 campaign, from which
rotation periods have been measured via light curve frequency analysis
(Douglas et al., 2017; Rebull et al., 2017). We compiled rotation periods,
Gaia photometry and Gaia parallaxes for members of the 650 Myr Praesepe
cluster (Fossati et al., 2008) by identifying Praesepe members with measured
rotation periods from Douglas et al. (2017) in the K2-Gaia crossmatch catalog
provided at https://gaia-kepler.fun/. This catalog cross-matched the EPIC
catalog (Huber et al., 2016) with the Gaia DR2 catalog (Gaia Collaboration et
al., 2018), using a 1” search radius. The result was a sample of 757 stars
with rotation periods, parallaxes and Gaia $G$, $G_{BP}$ and $G_{RP}$-band
photometry, shown in figure 1111This analysis was performed in a Jupyter
notebook available here:
https://github.com/RuthAngus/stardate/blob/master/paper/code/Praesepe.ipynb.
Although this new model does not perfectly describe the rotation period of
every star, it provides a better representation of rotational evolution than
previous models such as the Angus et al. (2015) model (shown as the blue solid
line in figure 1).

Figure 1:  The rotation periods of Praesepe members (Douglas et al., 2016),
vs. their Gaia colors ($G_{BP}-G_{RP}$) with a broken power law model, fit to
these data. The dashed line and shaded region shows the mean and variance
model described in equation 3 at 650 Myrs (lower model) and 4.56 Gyrs (upper
model). The solid blue line shows the (Angus et al., 2015) gyrochronology
relation at 650 Myrs (lower model) and 4.56 Gyrs (upper model). Shaded regions
show the 1$\sigma$ range of the rotation period model (equation 3). The
rotation periods of stars bluer than around 0.56 dex and redder than around
2.7 dex in $G_{BP}-G_{RP}$ are modeled as a broad log-normal distribution with
a standard deviation of 0.5 dex, added to observational uncertainties. This
figure was generated in a Jupyter notebook available at
https://github.com/RuthAngus/stardate/blob/master/paper/code/Fitting_Praesepe.ipynb

In order to fit a relation to Praesepe, we removed rotational outliers bluer
than $G_{BP}-G_{RP}$ = 2.7 via sigma-clipping and fit a 5th-order polynomial
to the remaining FGK and early M stars. We found that a 5th order polynomial
provided a substantially better fit than lower-order polynomials, which were
not able to capture the sharp ‘elbow’ in the rotation period-color relation.
Additional orders provided either a worse fit, tending towards extreme values
at the boundaries, or diminishing returns in goodness-of-fit. We also fit a
straight line to the late M dwarfs ($G_{BP}-G_{RP}$ $>$ 2.7), to capture the
mass-dependent initial rotation periods of low mass stars (Somers et al.,
2017). We fit a separable straight line function to the period-age relation
using the ages of Praesepe and the Sun222The fitting process was performed in
a Jupyter notebook available at
https://github.com/RuthAngus/stardate/blob/master/paper/code/Fitting_Praesepe.ipynb..
This new Praesepe-calibrated gyrochronology relation is,

$\log_{10}(P_{\mathrm{rot}})=c_{A}\log_{10}(t)+\sum_{n=0}^{4}c_{n}[\log_{10}(G_{BP}-G_{RP})]^{n}$
(1)

for stars with $G_{BP}-G_{RP}$ $<$ 2.7 and

$\log_{10}(P_{\mathrm{rot}})=c_{A}\log_{10}(t)+\sum_{m=0}^{1}b_{m}[\log_{10}(G_{BP}-G_{RP})]^{m}$
(2)

for stars with $G_{BP}-G_{RP}$ $>$ 2.7, where $P_{\mathrm{rot}}$ is rotation
period in days and $t$ is age in years. Best-fit coefficient values are shown
in table 1.

Table 1: Coefficient values for equations 1 and 2. Coefficient | Value
---|---
$c_{A}$ | 0.65 $\pm$ 0.05
$c_{0}$ | -4.7 $\pm$ 0.5
$c_{1}$ | 0.72 $\pm$ 0.05
$c_{2}$ | -4.9 $\pm$ 0.2
$c_{3}$ | 29 $\pm$ 2
$c_{4}$ | -38 $\pm$ 4
$b_{0}$ | 0.9 $\pm$ 0.5
$b_{1}$ | -13.6 $\pm$ 0.1

We calibrated this new relation using Praesepe and the Sun alone, without
including other clusters, because different open clusters have slightly
different period-color relationships (Agüeros et al., 2018; Agüeros, 2018;
Curtis and Agüeros, 2018; Curtis et al., 2019) and, given that we used an age-
color separable relation, adding in extra clusters was unlikely to improve the
calibration. We did not include asteroseismic stars because most are slightly
evolved and would require a gyrochronology relation that depends on $\log g$.
This new relation, fit to Praesepe and the Sun, does not perfectly predict the
rotation periods of stars at all colors and ages, but it provides several
improvements over previous empirical gyrochronology relations. Firstly, it
uses new K2 rotation period measurements to model the period-color relation of
Praesepe in detail, secondly, it includes a model for the rotational behavior
of M dwarfs, and thirdly, it is calibrated to Gaia $G_{BP}-G_{RP}$ color: a
directly observable quantity and the most widely available photometric color
index.

Equation 1, describing the rotational evolution of FGK and early M stars, is
most closely analogous to previously calibrated empirical gyrochronology
relations (e.g. Barnes, 2003, 2007; Mamajek and Hillenbrand, 2008; Barnes,
2010; Angus et al., 2015). It describes stars with Sun-like magnetic dynamos
that follow a ‘Skumanich-like’ magnetic braking law, i.e. their rotation
period increases with the square root of their age. It does not describe stars
hotter than around 6250 K which have a thin convective layers and a weak
magnetic dynamo, nor does it describe fully convective stars which take a long
time to converge onto the Skumanich (1972) braking law (Krishnamurthi et al.,
1997). In addition, stars with Rossby numbers larger than around 2 do not show
Skumanich-like magnetic braking (van Saders et al., 2016, 2018), so equation 1
does not describe these stars. It also does not describe the rotation periods
of subgiants or giants, whose rotation periods are influenced by their
expanding radii, changing winds, and core-to-surface differential rotation
(e.g. van Saders and Pinsonneault, 2013; Tayar and Pinsonneault, 2018).
Finally, this relation does not describe the rotation periods of dynamically
or magnetically interacting binaries which often rotate more rapidly than
isolated stars at the same age and color (Douglas et al., 2016). In order to
include non-Skumanich type stars in our combined gyrochronal and isochronal
model, we designed a composite gyrochronology relation which describes a mean
and variance model for the rotation period distributions across the HRD.
Rotation periods are modeled differently for stars of different photometric
color, EEP, Rossby number, age and metallicity. The rotation period model for
each of these groups is described below.

  * •

The rotation periods of late F, GK, and early M dwarfs (0.56 $<$
$G_{BP}-G_{RP}$ $<$ 2.7), with Rossby numbers less than 2, are modeled as a
log-normal distribution, with mean given by equation 1, and variance given by
the squared inverse of their observational uncertainties.

  * •

The rotation periods of fully convective stars ($G_{BP}-G_{RP}$ $>$ 2.7) are
modeled with a log-normal distribution with mean given by equation 2, and an
extra standard deviation of 0.5 dex, added to observational uncertainties.
This distribution reflects the observed rotation periods of late M dwarfs in
the Praesepe cluster which span a broad range at every color. Late M dwarfs
with masses $\lesssim$ 0.3 M⊙, temperatures $\lesssim$ 3500 and
$G_{BP}-G_{RP}$$\gtrsim$ 2.7 exhibit weak magnetic braking until at least
after the age of Praesepe ($\sim$650 million years).

  * •

The rotation periods of F-type and hotter stars, with $G_{BP}-G_{RP}$ $<$ 0.56
dex, are modeled as a log-normal distribution with mean,
$\log_{10}(P_{\mathrm{rot}})=0.56$, and an additional standard deviation,
added to observational uncertainties, of 0.5 dex. Stars more massive than
around 1.25 M⊙, with a temperature $\gtrsim$ 6250 K and a $G_{BP}-G_{RP}$
$\lesssim$ 0.56 do not spin down appreciably over their main-sequence
lifetimes because they do not have the deep convective envelope needed to
generate strong magnetic fields. This model for the mean and variance of hot
star rotation periods is based on stars hotter than 6250 K in the McQuillan et
al. (2014) sample, which have a mean $\log_{10}(P_{\mathrm{rot}})$ of 0.56 dex
and a standard deviation of around 0.5 dex. Both hot and cool stars retain
rotation periods that are similar to their primordial distribution (see e.g.
Matt et al., 2012; Somers et al., 2017).

  * •

The rotation periods of stars with large Rossby numbers are modeled as
follows. The age at which a star’s Rossby number would exceed 2 is calculated
by inverting equation 1. If a star’s age is greater than this, its mean
rotation period is given by $P_{\mathrm{max}}=2/\tau$, where $\tau$ is the
convective turnover timescale, calculated via stellar mass using equation 11
from Wright et al. (2011).

  * •

The rotation periods of subgiants are described with a log-normal distribution
with mean given by equation 1, 2 or 0.56, depending on its color, and an
additional standard deviation of 5 dex. This is not an accurate model of
subgiant rotation periods (see, e.g. van Saders and Pinsonneault, 2013) but
the highly inflated variance makes it a weakly constraining one. Isochrone
fitting provides precise ages for subgiants, so by inflating the variance of
the gyrochronology relation at large EEP, we allow a star’s position on the
HRD/CMD to dominate the age information over its rotation period. This is
useful because we have not yet built an accurate rotation-age relation for
subgiants into our model.

  * •

We model stars younger than around 250 Myrs with a log-normal distribution
with mean function given by equation 1, 2 or 0.56, depending on it whether it
has an intermediate, red, or blue color respectively and a inflated standard
deviation of 0.5 dex. Rotation periods in young open clusters show a large
amount of scatter because they have not yet converged onto the Skumanich
(1972) spin-down sequence.

  * •

Finally, we model stars with very low and high metallicites (-0.2 $>$ [Fe/H]
$>$ 0.2) with a log-normal distribution with mean given by equation 1, 2 or
0.56, depending on its color, and a standard deviation of 0.5 dex.
Gyrochronology is not calibrated at these extreme metallicities due to a lack
of suitable metal poor and rich calibration stars. Rather than assume the same
gyrochronology model can be naïvely applied to these stars, we take a more
conservative approach and model them with a broad Gaussian distribution.

Inflating the variance of the rotation period distribtion for stars with non-
Skumanich magnetic braking behavior has two purposes: 1) in the cases of hot
stars and fully convective stars, it allows the broad distributions of
rotation periods observed in clusters and the field to be matched, and 2) it
down-weights the age-information provided by rotation periods in regions of
the HRD/CMD where rotation periods are not information-rich or the
gyrochronology model is inaccurate or uncalibrated. If the observational
uncertainties on rotation periods are 5% on average, which corresponds to
0.05/$\ln(10)$ $\sim$ 0.02 dex, adding 0.5 amounts to a 25$\times$ increase in
standard deviation, or a $\sim$ 600$\times$ increase in variance. So,
practically speaking, when the standard deviation is inflated to 0.5 dex or
more, ages are almost entirely inferred via isochrone fitting.

We used the following composite gyrochronology model to infer ages from
rotation periods,

$\log_{10}(P_{\mathrm{rot}})\sim\begin{cases}\mathcal{N}\left[(\ref{eqn:fgk_gyro}),(\sigma+\sigma_{P})^{2}\right],&Ro<2,0.56<G_{BP}-G_{RP}<2.7\\\
\mathcal{N}[(\ref{eqn:m_gyro}),(\sigma+\sigma_{P})^{2}],&Ro<2,G_{BP}-G_{RP}>2.7\\\
\mathcal{N}[\log_{10}(P_{\mathrm{max}}),(\sigma+\sigma_{P})^{2}],&Ro\geq 2\\\
\mathcal{N}[0.56,(\sigma+\sigma_{P})^{2}],&G_{BP}-G_{RP}<0.56,\end{cases}$ (3)

where $\sigma$ is the relative period uncertainty, divided by $\ln(10)$, on
individual rotation period measurements and $\sigma_{P}$ is an additional
scatter that is a function of EEP, age, metallicity and color. It takes a
maximal value of 0.5 for hot stars and fully convective stars, 5 for subgiants
and giants, and a minimal value of zero for late F, GK and early M dwarfs. The
variance model is shown in figure 2. Sigmoid functions were used to provide
smooth transitions between regions of low and high variance. Sharp changes in
variance would produce sharp changes in likelihood, which would cause the
posterior distributions over stellar parameters to be more difficult to
sample. The sigmoid functions shown in figure 2 reach half their maximum
values at $G_{BP}-G_{RP}$ = 0.56 and 2.5 for hot and cool stars respectively,
EEP = 454 for subgiants, age = 250 Myrs for young stars, $[$M/H$]$ = -0.2 for
metal poor stars and 0.2 for metal rich stars. The logistic growth rate, or
steepness, of the sigmoid functions is 100 dex-1 for both color transitions,
.2 EEP-1 for the EEP transition, 20 dex-1 for the age transition and 5 dex-1
for the both metallicity transitions. The additional standard deviation is
additive, so if a star is e.g. hot, evolved and metal poor, the additional
standard deviation of its rotation period rises to 6.

Figure 2:  The additional rotation period scatter, $\sigma_{P}$, added to the
observational period uncertainties in the model (see equation 3). The standard
deviation was increased for early F and hotter stars ($G_{BP}-G_{RP}$$<$
0.56), late M dwarfs ($G_{BP}-G_{RP}$ $>$ 2.7) and evolved stars (EEP
$\gtrsim$ 420) in order to down-weight the age-information supplied by
rotation periods and reproduce observed rotation period distributions. We also
increased the variance for stars younger than around 250 Myrs, because the
rotation periods of these stars typically have not yet converged onto a tight
gyrochronology sequence, and for very high and low metallicity stars (-0.2
$\gtrsim$ [Fe/H] $\gtrsim$ 0.2) because the gyrochronology relations have not
been calibrated at these extreme values. Down-weighting the gyrochronal
likelihood by the inverse variance ($1/\sigma^{2}$) allowed the ages of these
stars to be mostly inferred via isochrone fitting. Sigmoid functions were used
to provide smooth transitions between regions of low and high variance.

###  2.2 Simultaneously fitting gyrochronology and isochrones

The previous part of this section describes the model for the mean and
variance of rotation periods as a function of their ages (and colors, EEPs,
masses and metallicities). In what follows, we describe how this model was
combined with a stellar evolution model to infer stellar ages (and other
parameters) via gyrochronology and isochrone fitting simultaneously. Our goal
was to infer the age of a star from its observable properties by estimating
the posterior probability density function (PDF) over age,

$p(t|{\bf m_{x}},P_{\mathrm{rot}},\bar{\pi}),$ (4)

where $t$ is age, ${\bf m_{x}}$ is a vector of apparent magnitudes in various
bandpasses, $P_{\mathrm{rot}}$ is the rotation period and $\bar{\pi}$ is
parallax. Spectroscopic properties ($T_{\mathrm{eff}}$, log(g) and
$\mathrm{[Fe/H]}$) and/or asteroseismic parameters ($\Delta\nu$ and
$\nu_{\mathrm{max}}$) may also be available for a star, in which case they
would appear to the right of the ‘$|$’ in the above equation since they are
observables. In order to calculate a posterior PDF over age, other stellar
parameters must be marginalized over. These parameters are distance ($D$),
V-band extinction ($A_{V}$), the inferred metallicity, $[M/H]$333The inferred
metallicity, [M/H] is a model parameter which is different to the observed
metallicity, [Fe/H] which would appear on the right side of the $`|^{\prime}$
in equation 4., and equivalent evolutionary phase (abbreviated to EEP or $E$).
EEP is a dimensionless number ranging from around 200 for M dwarfs up to
around 1600 for giants and is 355 for the Sun (see Dotter, 2016; Choi et al.,
2016). Stars are defined as subgiants when their EEP exceeds 454. Mass is
uniquely defined by EEP, age and metallicity. The marginalization involves
integrating over these extra parameters,

$\displaystyle p(t|{\bf m_{x}},P_{\mathrm{rot}},\bar{\pi})$
$\displaystyle\propto\int p({\bf
m_{x}},P_{\mathrm{rot}},\bar{\pi}|t,E,[M/H],D,A_{V})~{}p(t)\,p(E)\,p([M/H])\,p(D)\,p(A_{V})\,\mathrm{d}E\,\mathrm{d}[M/H]\,\mathrm{d}D\,\mathrm{d}A_{V}.$

This equation is a form of Bayes’ rule,

$\mathrm{Posterior}\propto\mathrm{Likelihood}\times\mathrm{Prior},$ (6)

where the likelihood of the data given the model is,

$p({\bf m_{x}},\bar{\pi},P_{\mathrm{rot}}|t,E,D,A_{V},[M/H]),$ (7)

and the prior PDF over parameters is,

$p(t)\,p(E)\,p(D)\,p(A_{V})\,p([M/H]).$ (8)

The priors we used are described in the appendix.

We assumed that the process of magnetic braking is independent of hydrogen
burning in the core, outside of the dependencies that are captured in the
model. This assumption allowed us to multiply two separate likelihood
functions together: one computed using an isochronal model and one computed
using a gyrochronal model. We assumed that the probability of observing the
measured observables given the model parameters was a Gaussian and that the
observables were identically and independently distributed. The isochronal
likelihood function was,

$\displaystyle\mathcal{L_{\mathrm{iso}}}=p({\bf
m_{x}},\bar{\pi}|t,E,[M/H],D,A_{V})$
$\displaystyle=\frac{1}{\sqrt{(2\pi)^{n}\det(\Sigma)}}\exp\left(-\frac{1}{2}({\bf
O_{I}}-{\bf I})^{T}\Sigma^{-1}({\bf O_{I}}-{\bf I})\right),$

where ${\bf O_{I}}$ is the vector of n observables: $\bar{\pi}$, ${\bf m_{x}}$
plus spectroscopic and/or asteroseismic observables if available, and $\Sigma$
is the covariance matrix of that set of observables. ${\bf I}$ is the vector
of model observables that correspond to a set of parameters: $t$, $E$,
$[M/H]$, $D$ and $A_{V}$, calculated using an isochrone model. We assumed
there is no covariance between these observables and so this covariance matrix
consists of individual parameter variances along the diagonal with zeros
everywhere else. The gyrochronal likelihood function was,

$\displaystyle\mathcal{L_{\mathrm{gyro}}}=p(P_{\mathrm{rot}}|t,E,[M/H],D,A_{V})$
(10)
$\displaystyle=\frac{1}{\sqrt{(2\pi)\det(\Sigma_{P})}}\exp\left(-\frac{1}{2}({\bf
P_{O}}-{\bf P_{P}})^{T}\Sigma_{P}^{-1}({\bf P_{O}}-{\bf P_{P}})\right),$

where ${\bf P_{O}}$ is a 1-D vector of observed logarithmic rotation periods,
and ${\bf P_{P}}$ is the vector of corresponding logarithmic rotation periods,
predicted by the model. $\Sigma_{P}$ was comprised of individual rotation
period measurement uncertainties, plus an additional variance that is a
function of EEP and $G_{BP}-G_{RP}$ color, added in quadrature. This variance
accounts for the stochastic nature of the rotation periods of very hot and
very cool stars and allowed us to predominantly use isochrone fitting to
measure the ages of subgiants. The full likelihood used in our model was the
product of these two likelihood functions,

$\mathcal{L}_{\mathrm{full}}=\mathcal{L}_{\mathrm{iso}}\times\mathcal{L}_{\mathrm{gyro}}.$
(11)

The inference processes proceded as follows. First, a set of parameters: age,
EEP, metallicity, distance and extinction, and observables for a single star
were passed to the isochronal likelihood function in equation (2.2). Then, a
set of observables corresponding to those parameters were generated from the
MIST model grid using isochrones.py (Morton, 2015) and compared to the
measured observables via the isochronal likelihood,
$\mathcal{L}_{\mathrm{iso}}$ (also computed using isochrones.py). The
parameters were also passed to the gyrochronology model (equation 3) where
$t$, $E$, $[M/H]$, $D$ and $A_{V}$ were used to calculate $G_{BP}-G_{RP}$
color and mass from the MIST model grid and, in turn, rotation period via the
gyrochronology model. EEP and $G_{BP}-G_{RP}$ were also used to calculate the
additional rotation period variance, added to the individual period
uncertainties. This model rotation period was compared to the measured
rotation period using gyrochronal likelhood function of equation 10. The
gyrochronal log-likelihood was added to the isochronal log-likelihood to give
the full likelihood, which was then added to the log-prior to produce a single
sample from the posterior PDF.

Ages were inferred with stardate using Markov Chain Monte Carlo. The joint
posterior PDF over age, mass, metallicity, distance and extinction was sampled
using the affine invariant ensemble sampler, emcee (Foreman-Mackey et al.,
2013) with 50 walkers. Samples were drawn from the posterior PDF until 100
independent samples were obtained. We actively estimated the autocorrelation
length, which indicates how many steps were taken per independent sample,
after every 100 steps using the autocorrelation tool built into emcee. The
MCMC concluded when either 100 times the autocorrelation length was reached
and the change in autocorrelation length over 100 samples was less than 0.01,
or the maximum of 500,000 samples was obtained. This method is trivially
parallelizable, since the inference process for each star can be performed on
a separate core. The age of a single star can be inferred in around 1 hour on
a laptop computer.

##  3 Results

In order to demonstrate the performance of our method, we conducted three sets
of tests. In the first we simulated observables from a set of stellar
parameters for a few hundred stars using the MIST stellar evolution models and
the gyrochronology model of equation 3. The ages predicted with our model were
compared to the true parameters used to generate the data. In the second we
tested our model by measuring the ages of individual stars in the NGC 6819
open cluster, and in the third we tested our model on Kepler asteroseismic
stars.

###  3.1 Test 1: simulated stars

For the first test we drew masses, ages, bulk metallicities, distances and
extinctions at random for 1000 stars from the following uniform distributions:

$\displaystyle\mathrm{EEP}\sim U(198,480)$ (12) $\displaystyle t\sim
U(0.5,14)\mathrm{~{}[Gyr]}$ (13) $\displaystyle[Fe/H]\sim U(-0.2,0.2)$ (14)
$\displaystyle D\sim U(10,1000)~{}\mathrm{[pc]}$ (15) $\displaystyle A_{V}\sim
U(0,0.1).$ (16)

$T_{\mathrm{eff}}$, log(g), [Fe/H], parallax, and apparent magnitudes $B$,
$V$, $J$, $H$, $K$, Gaia $G$, $G_{BP}$ and $G_{RP}$ were generated from these
stellar parameters using the MIST stellar evolution models. We added a small
amount of noise to the ‘observed’ stellar properties in order to reflect
optimistic observational uncertainties for isochrone-dating. We added Gaussian
noise with a standard deviation of 25 K to $T_{\mathrm{eff}}$, 0.01 dex to
[Fe/H] and log(g), and 10 mmags to $B$, $V$, $J$, $H$, and $K$ magnitudes.
These are just one choice of uncertainties that we could have adopted and are
extremely optimistic. The uncertainties on predicted ages will depend strongly
on all observational uncertainties, however, since this analysis is designed
to show the relative improvement in stellar age precision when rotation
periods are included, we chose to use best-case spectroscopic parameters. The
noise added to Gaia $G$-band photometry ranged from 0.3 mmag for stars
brighter than 13th magnitude, to 10 mmag for stars around 20th magnitude
(Evans et al., 2017; Gaia Collaboration et al., 2018). Noise added to Gaia
$G_{BP}$ and $G_{RP}$ bands ranged from 2 mmag for stars brighter than 13th
magnitude to 200 mmag for stars fainter than 17th. Unphysical combinations of
stellar parameters were discarded, resulting in a final sample size of 841
simulated stars. Figure 3 shows the position of these stars on an HRD (with
log(g) on the y-axis instead of luminosity to improve the visibility of the
MS), colored by their age. Rotation periods for these stars were generated
using the gyrochronology relation described in equation 3. We added 5%
Gaussian noise to all stellar rotation periods to represent realistic
measurement uncertainties of 5%. The median uncertainty on rotation periods
calculated from Kepler light curves, provided in the McQuillan et al. (2014)
catalog is 1%. However, the Aigrain et al. (2015) injection and recovery study
showed that true rotation period uncertainties are often slightly larger than
this, and noise distribution of rotation periods can be highly non-Gaussian
(e.g. Aigrain et al., 2015; Angus et al., 2018).

Figure 3:  The simulated star sample plotted on an HRD, colored by age (top
panel) and rotation period (bottom panel). HRD positions were calculated using
MIST isochrones via the isochrones.py Python package and rotation periods were
generated using equation 3. This figure was generated in a Jupyter notebook
available at
https://github.com/RuthAngus/stardate/blob/master/paper/code/Simulate_data.ipynb

Figure 4 shows the rotation periods of 841 stars generated from the
gyrochronology model.

Figure 4:  Data simulated from the rotation period model. Late F, GK and early
M dwarfs (stars with 0.56 $<$ $G_{BP}-G_{RP}$ $<$ 2.7 follow the Praesepe-
calibrated gyrochronology relation (dashed gray lines), with the exception of
old, slowly rotating stars with large Rossby numbers whose rotation periods
are fixed at 2$\times$ their convective overturn time. The rotation periods of
early F ($G_{BP}-G_{RP}$$<$ 0.56), late M dwarfs ($G_{BP}-G_{RP}$$<$ 2.7) and
subgiants (EEP $\gtrsim$ 420) were generated from a log-normal distribution
with standard deviation given by equation 3. The top panel shows the rotation
periods vs. $G_{BP}-G_{RP}$ colors of simulated stars, colored by their age
and the bottom panel shows the same stars colored by their equivalent
evolutionary phase (EEP). The gray lines describe the mean gyrochronology
model at ages 1, 3, 5, 7, 9, 11, and 13 (rotation periods rise with age). This
figure was generated in a Jupyter Notebook, available at
%The␣top␣panel␣shows␣the␣rotation␣periods␣vs.␣B-V␣colors␣of␣simulated␣stars,%9,␣11,␣and␣13␣(rotation␣periods␣rise␣with␣age).https://github.com/RuthAngus/stardate/blob/master/paper/code/Simulate_data.ipynb

We took two approaches to inferring the ages of these simulated stars: firstly
using isochrone fitting only, and secondly using isochrone fitting combined
with a gyrochronology model (stardate). Since the posterior PDFs of stars are
often multimodal, we found that the choice of initial positions of the emcee
walkers influenced the final outcome because walkers occasionally got stuck in
local minima. We found that the following set of initial parameters worked
well, though not perfectly: EEP = 330, $t=9.56$ Gyr, $[M/H]=-0.05$, $D=269$ pc
and $A_{V}=0.0$444These are the default initial parameters provided in
stardate.. Figure 5 shows the results of combining gyrochronology with
isochrone fitting for the simulated sample. The stars’ true ages are plotted
against their predicted ages, with ages inferred with gyrochronology and
isochrone fitting in color, and ages inferred using isochrone fitting only
plotted in light grey. The different panels show the results for different
types of stars: FGK dwarfs that are still undergoing magnetic braking, FGK
dwarfs that have ceased magnetic braking (their Rossby number is around 2), M
dwarfs, and evolved stars. The selection criteria for these groups are in the
panel headings. The FGK dwarfs with low Rossby numbers showed the largest
improvement: the median age precision for this group (defined as the standard
deviation of the posterior as a percentage of the median age) was 8% when
using a combination of isochrones and gyrochronology, and 22% using isochrone
fitting alone. This equates to an almost 3$\times$ improvement in age
precision. The age RMS of this group was 0.8 Gyr using isochrones and
gyrochronology, and 2 Gyr using isochrones only. Despite the fact that stars
with Rossby numbers of 2 have stopped spinning down and their rotation periods
no longer evolve with age, the ages of these stars can still be relatively
precisely constrained with isochrone fitting. The median age precision for
stars with large Rossby numbers was 11% with gyrochronology and isochrone
fitting, and 13% with isochrone fitting only. The precision of M dwarf ages
did improve overall when their rotation periods were included, but this
improvement was entirely driven by the early M dwarfs. The precision of this
group improved overall from 33% to 22% and RMS from 5.4 Gyr to 4.6 Gyr. The
precision of ages inferred for evolved stars changed very little when rotation
periods were included in the inference process. This is because the variance
of the rotation-age relation was inflated by a large amount in the
gyrochronology model (equation 3), making rotation periods almost entirely
uninformative for this group. The median age precision of subgiants from both
gyrochronology and isochrone fitting and isochrone fitting alone was 7%.

Figure 5:  The true vs. predicted ages of simulated stars. Ages calculated by
combining gyrochronology and isochrone fitting with stardate are shown in
color and ages calculated with isochrone fitting only are shown in gray. The
different panels show the results for stars with $G_{BP}-G_{RP}$ $<$ 2.2 (FGK
dwarfs) that are still braking magnetically ($Ro<2$), stars with
$G_{BP}-G_{RP}$ $<$ 2.2 that have stopped braking magnetically ($Ro\geq 2$),
stars with 2.2 $<$ $G_{BP}-G_{RP}$ (M dwarfs), and evolved stars (EEP $>$
420). Gyrochronology is highly effective for FGK stars and ages inferred with
both gyrochronology and isochrone fitting are more accurate and precise than
ages inferred via isochrone fitting only for this group. Neither
gyrochronology nor isochrone fitting can provide precise ages for M dwarfs, so
the ages of these stars are imprecise regardless of age-dating method. This
figure was generated in a Jupyter Notebook available at
%****␣master.tex␣Line␣950␣****https://github.com/RuthAngus/stardate/blob/master/paper/code/Results_plots.ipynb.

An important caveat associated with these results is that they strongly depend
on the uncertainties adopted for all observables used in the analysis:
$T_{\mathrm{eff}}$, [Fe/H], log(g), $G$, $G_{BP}$, $G_{RP}$, $J$, $H$, $K$,
$B$, $V$, and rotation period. Changing the uncertainties on these observables
will affect the uncertainties on inferred ages in different ways. This
demonstration is not intended to reflect the typical age uncertainties that
will result for all stars in reality, it merely exemplifies the increase in
stellar age precision that results from one specific choice of uncertainties.
Needless to say, estimating the uncertainties on observables accurately can be
as important as accurately measuring the observables themselves.

This simulation experiment was designed to show the theoretical improvement in
age precision when gyrochronology is incorporated into isochrone fitting.
However, it does not demonstrate the accuracy of this method because the test
data were simulated from the same model used to infer ages. The results of
this experiment are therefore extremely accurate by design. When applying this
method to real data, the results will only be accurate if the model is
accurate. In other words, stardate, like any age-dating method, provides
model-dependent ages. Stellar ages calculated with stardate depend on both the
accuracy of the MIST models and the accuracy of the gyrochronology model
(equation 3). In order to test the accuracy of stardate, we applied it to real
data, as described in the following section.

###  3.2 Test 2: Open clusters

In order to test our model on real stars with known ages, we selected a sample
of stars in the 2.5 Gyr NGC 6819 cluster. We compiled Kepler-based rotation
periods (Meibom et al., 2015), Gaia photometry and Gaia parallaxes for members
of the NGC 6819 cluster. Figure 6 shows the period-color relation of this
cluster and figure 7 shows the results of inferring the ages of individual
cluster members using a combination of gyrochronology and isochrone fitting
(via stardate) and isochrone fitting alone. The ages of F stars in this
cluster ($G_{BP}-G_{RP}$ $\sim$ 0.5-0.75) were relatively precisely
constrained by isochrone fitting alone because, at 2.5 Gyr, they are
approaching the MS turnoff. For these hot stars, ages inferred with
gyrochronology and isochrones were similar to ages inferred with isochrones
and similarly precise, showing that isochrones provide a lot of age
information for these stars and rotation periods do not add significantly more
information. The G and early K stars in this cluster ($G_{BP}-G_{RP}$
$\lesssim$ .65) were not precisely recovered from isochrone fitting alone –
the isochrone-only age posteriors tend towards the prior which is a uniform
distribution between 0 and 13.8 Gyrs. The median age of stars in the cluster
was 4.27 $\pm$ 0.48 Gyr when only isochrone fitting was used. In contrast,
including gyrochronology when inferring the ages of G and K stars in this
cluster significantly improved age precision. The median age of stars in this
cluster was 2.63 $\pm$ 0.16 Gyr when ages were inferred with a combination of
isochrone fitting and gyrochronology, using the newly calibrated Praesepe-
based gyrochronology model. The previously-calibrated Angus et al. (2015)
model resulted in a median stellar age of 2.66 $\pm$ 0.21 Gyr, which is still
consistent with the established cluster age of 2.5 Gyr. The median age of
stars in the cluster using uncorrected photometry was slightly underestimated
at 1.86 $\pm$ 0.22 Gyr. This suggests that, despite the fact that V-band
extinction is marginalized over during the inference process, correcting for
extinction before ages are estimated will reduce bias introduced by dust.

Figure 6:  The Kepler-based rotation periods of members of the 2.5 Gyr NGC
6819 open cluster. The raw $G_{BP}-G_{RP}$ colors are shown in red and the
dust-corrected colors are shown in black. The dashed line shows a
gyrochronology model that was fit to the Praesepe cluster and the Sun in this
work, interpolated to 2.5 Gyrs. The solid blue line shows a previously
calibrated gyrochronology model (Angus et al., 2015). This figure was
generated in a Jupyter Notebook available at:
https://github.com/RuthAngus/stardate/blob/master/paper/code/NGC6819.ipynb
Figure 7:  The inferred ages of members of the NGC 6819 open cluster as a
function of their $G_{BP}-G_{RP}$ color. Ages of stars inferred using a
combination of isochrone fitting and gyrochronology (Praesepe and Sun
calibration) with dereddened Gaia $G$, $G_{BP}$, and $G_{RP}$ photometry
(black circles) and uncorrected, raw, photometry (red squares). Even though
V-band extinction is marginalized over in the inference process, reddening can
still bias ages. Blue triangles, pointing up, show ages inferred using
isochrone fitting and gyrochronology, with the Angus et al. (2015)
gyrochronology model. Orange triangles, pointing down, show ages inferred
using isochrone fitting only. The ages of F stars (stars bluer than 0.7) were
precisely constrained by isochrones and including gyrochronology makes little
difference to their inferred ages. The age precision of G and K dwarfs (stars
redder than 0.7) was improved by including gyrochronology. The median age of
stars inferred using the gyrochronology model calibrated to Praesepe and the
Sun (black circles) was 2.65 $\pm$ 0.13 which is consistent with the
established cluster age (2.5 Gyr). This figure was generated in a Jupyter
Notebook available at:
%****␣master.tex␣Line␣1175␣****https://github.com/RuthAngus/stardate/blob/master/paper/code/NGC6819.ipynb

We found that 5% rotation period uncertainties resulted in the most accurate
ages for NGC 6819. The uncertainties on the measured rotation periods,
provided in Meibom et al. (2015) and shown in figure 6, were likely
underestimated for some stars. Underestimated rotation period uncertainties
can result in inaccurate age estimates. This raises the question: how should
uncertainties on rotation periods be estimated? The likelihood is weighted by
the inverse variance, so uncertainties on the rotation period control the
relative information provided by gyrochronology, isochrones, and the prior. If
rotation period uncertainties are either too large or too small, the resulting
age estimate will be imprecise and/or inaccurate. It is difficult to measure
uncertainties on rotation periods directly: standard techniques such as Lomb-
Scargle periodograms and autocorrelation functions do not provide them.
Ideally, rotation period uncertainties should capture both the measurement
precision, and the physical uncertainty introduced by the latitudinal movement
of star spots on the surface of a differentially rotating star. For example,
Donahue et al. (1996) demonstrated that the seasonal variation in measurements
of G and K star rotation periods is a function of period,
$\Delta\mathrm{P_{rot}}\propto\mathrm{P_{rot}}^{1.3\pm 0.1}$. This variation
is presumably caused by a latitudinal drift in the dominant active regions,
which traces the stellar cycle over several years, in combination with
latitudinal differential rotation. This suggests that latitudinal spot
drifting does not significantly affect stellar rotation periods when stars are
young, for example the scatter of rotation periods about the mean
gyrochronology model in Praesepe is only around 5%. However it is likely that
this effect will become more important at older ages. A thorough exploration
of how rotation period uncertainties, from both measurement uncertainty and
physical variation, affect stellar ages via gyrochronology is key to
understanding the power of gyrochronology as an age-dating method. For now, we
leave this exploration for a future study.

###  3.3 Test 3: Kepler asteroseismic stars

In order to test our method in the regime where both isochrone fitting and
gyrochronology become important, we recovered the ages of the 21 asteroseismic
stars analysed in van Saders et al. (2016). These 21 stars were observed in
Kepler’s short cadence mode and are a mixture of dwarfs and subgiants. Their
asteroseismic ages were calculated from the analysis of the frequencies of
individual oscillation modes (Mathur et al., 2012; Metcalfe et al., 2014;
Silva Aguirre et al., 2015; Ceillier et al., 2016), and their rotation periods
from their Kepler light curves (García et al., 2014). We crossmatched these
stars with the Gaia catalog to obtain parallaxes and apparant magnitudes in
the Gaia G, $\mathrm{G_{BP}}$ and $\mathrm{G_{RP}}$ band passes. We also added
J, H and K 2MASS magnitudes from the Kepler input catalog (Brown et al.,
2011), and used spectroscopic effective temperatures, spectroscopic
metallicities and rotation periods reported in table 1 of van Saders et al.
(2016). The Gaia photometry is extremely precise, and we found that
artificially inflating the uncertainties on Gaia apparent magnitudes by a
factor of 10 substantially improved the quality of fit, both in terms of MCMC
convergence and agreement with asteroseismic age measurements. In figure 8 we
show the ages of these 21 stars inferred using isochrone fitting and
gyrochronology, against their asteroseismic ages, calculated using the
Asteroseismic Modeling Portal (AMP) (Metcalfe et al., 2009, 2012, 2014). In
this figure, the colored symbols show ages inferred using a combination of
gyrochronology and isochrone fitting, implemented with the stardate Python
package. The black and grey triangles show ages inferred from gyrochronology
only, where the mass and color of the stars were not inferred, but fixed to be
the asteroseismic mass and the Gaia $G_{BP}-G_{RP}$ color. These gyrochronal
ages were calculated with age as the only free parameter, without
marginalizing over stellar mass, $G_{BP}-G_{RP}$ color or extinction. As a
result, these gyrochronal ages are not the same as the ages given by the
maximum of the gyrochronal likelihood function used in the combined age model,
they simply represent an approximation to the gyrochronal age. Grey triangles
are shown for stars where gyrochronology is not applicable because the stars
are either too metal poor, too metal rich, or too evolved. Black triangles are
shown for stars where gyrochronology is applicable. The white circles show
ages inferred from isochrones only. Dashed lines connect the three different
age measurements for the same stars.

Although the ages of all 21 stars shown in figure 8 were inferred with a joint
isochronal and gyrochronal model, most (all but 8) were either too evolved,
too metal poor, or too metal rich for gyrochronology to contribute any
information to the ages. These metal poor/rich or evolved stars lie in a
regime where the variance on their rotation period was artificially inflated
because the gyrochronology relations are not thoroughly understood or well
calibrated. The rotation periods of the remaining eight stars did contribute
to their inferred ages to some degree, however isochrones still dominated the
age information for some of them because these stars are relatively old and/or
relatively hot.

In general, there is relatively poor agreement between the asteroseismic ages
and the ages inferred using stardate. Much of this discrepancy is driven by
differences in the isochronal ages, which is likely attributable to
differences between the MIST stellar evolution models and those used in the
AMP analysis: a combination of the Aarhus stellar evolution code (ASTEC
Christensen-Dalsgaard, 2008a) and the adiabatic pulsation code (ADIPLS
Christensen-Dalsgaard, 2008b). We compared non-rotating, Solar-metallicity
MIST isochrones for middle-aged stars with Solar-metallicity BaSTI isochrones
(Pietrinferni et al., 2004; Hidalgo et al., 2018) and found that, for stars
between 4 and 8 Gyrs, in the same effective temperature range as the
asteroseismic stars, the age discrepancy between the two sets of models can be
as large as 1-2 billion years. The MIST isochrones lie above the BaSTI
isochrones on the HR diagram, leading to a systematic underprediction of ages.

The gyrochronal ages, where gyrochronology is applicable, do not show
excellent agreement with the asteroseismic ages either. The four hot stars to
the left in figure 8 are rotating more slowly than predicted by the Praesepe-
based gyrochronology models, and as a result their gyro-ages are older than
their asteroseismic ones. For these hot stars, two out of four have ages that
are still consistent, or close to consistent, with their asteroseismic ages.
The third star from the left is an anomalously slow rotator for its age and
mass and van Saders et al. (2016) also found this star to be surprisingly
slowly rotating. In contrast, the star with a black triangle symbol
(indicating that gyrochronology is applicable) furthest to the right is
rapidly rotating for its mass and age, even when weakened braking is taken
into account. This star is around Solar mass (1.0 $\pm$ 0.03 M⊙) and has an
asteroseismic age older than the Sun (7.28 $\pm$ 0.51 Gyr), yet rotates with a
period of only 19.8 $\pm$ 1.3 days. This is KIC 9098294, a single-lined
spectroscopic binary with an orbital period of around 20 days (Latham, private
communication). It is the only clear SB1 in the van Saders et al. (2016)
sample, although some others do have binary companions with long orbital
periods, for which tidal interactions are not expected to be strong.

Figure 8:  A comparison of stellar ages inferred using asteroseismic modeling
with ages inferred using a combination of isochrone fitting and
gyrochronology. Colored circles show ages inferred using isochrone fitting and
gyrochronology combined via the stardate software package. Black triangles
show the ages of all stars inferred via gyrochronology only and white circles
show ages inferred via isochrone-fitting only.

##  4 Conclusions

We present a statistical framework for measuring precise ages of MS stars and
subgiants by combining observables that relate, via different evolutionary
processes, to stellar age. Specifically, we combined HRD/CMD placement with
rotation periods, in a hierarchical Bayesian model, to age-date stars based on
both their hydrogen burning and magnetic braking history. The two methods of
isochrone fitting and gyrochronology were combined by taking the product of
two likelihoods: one that contains an isochronal model and the other a
gyrochronal one. We used the MIST stellar evolution models and computed
isochronal ages and likelihoods using the isochrones Python package. We fit a
new broken power law gyrochronology model to the Praesepe cluster and included
a modification recommended by van Saders et al. (2016) that accounts for
weakened magnetic braking at Rossby numbers larger than 2. The rotation
periods of hot stars, cool stars, evolved stars, young stars and metal poor
and rich stars were modeled with a broad log-normal distribution. We tested
stardate on simulated data and cluster stars and demonstrated that combining
gyrochronology with isochrone fitting improves the precision of age estimates
for FGK dwarfs by a factor of 3 over isochrone fitting alone, assuming 5%
measurement uncertainties on rotation periods. Incorporating rotation periods
into stellar evolution models also improves the precision of the equivalent
evolutionary phase (EEP) parameter and, since EEP, combined with age and
metallicity, determines the mass, radius and log(g) of a star, this means
rotation periods can improve the precision of all stellar parameters. Although
V-band extinction is marginalized over during inference, correcting photometry
for dust-extinction before analysis, or including it as a prior can improve
the accuracy of stellar ages measured with stardate. We also tested stardate
on a set of 21 Kepler asteroseismic stars (van Saders et al., 2016). We found
that discrepancies between ages measured with stardate and ages measured with
asteroseismology are likely produced by differences between the MIST and BaSTI
stellar evolution models. Asteroseismic and cluster stars provide an
opportunity for calibration but given the high-dimensionality of the
gyrochronology relations (i.e. rotation period depends on age, mass,
metallicity, surface gravity, etc), many stars with precise ages, spanning a
range of properties, are still needed to reliably calibrate them.

In cases where gyrochronology predicts inaccurate stellar ages it is either
because models are not correctly calibrated, because the rotation periods or
rotation period uncertainties are themselves inaccurate, or because of
rotational outliers. For example, stardate may predict inaccurate ages for
stars in close binaries whose interactions influence their rotation period
evolution. Rotational outliers are often seen in clusters (see e.g. Douglas et
al., 2016; Rebull et al., 2016; Douglas et al., 2017; Rebull et al., 2017) and
many of these fall above the main sequence on a CMD, indicating that they are
binaries. In addition, measured rotation periods may not always be accurate
and can, in many cases, be a harmonic of the true rotation period. For
example, a common rotation period measurement failure mode is to measure half
the true rotation period. The best way to prevent an erroneous or outlying
rotation period from resulting in an erroneous age measurement is to allow for
outlying rotation periods using a mixture model, a feature that could be built
into stardate in the future.

The optimal way to age-date stars is by combining all their available age-
related observables. This could ultimately include activity dating via flare
rates and chromospheric activity indices, kinematic dating and chemical
dating. Of all the established age-dating methods, gyrochronology and
isochrone fitting are two of the most complementary. The two methods are
optimal in different parts of the HRD: gyrochronology works well for FGK
dwarfs and isochrone fitting works well for subgiants and hot stars, so
combining the two methods results in consistently precise ages across a range
of masses, ages and evolutionary stages. In addition, using both methods at
once circumvents the need to decide which method to use a priori. It
eliminates the circular process of classifying a star based on its CMD
position (M dwarf, subgiant, etc), then deciding which age-dating method to
use, then inferring an age which itself depends on the classification that was
made. It is important to infer all stellar properties at once since they all
depend on each other. stardate is applicable to a large number of stars: FGKM
dwarfs and subgiants with a rotation period and broad-band photometry. This
already includes tens-of-thousands of Kepler and K2 stars and could include
millions more from TESS, LSST, WFIRST, PLATO, Gaia, and others in the future.

The code used in this project is available as a Python package called
stardate, DOI: https://doi.org/10.5281/zenodo.2712419. It is available for
download via Github555git clone https://github.com/RuthAngus/stardate.git or
through PyPI666pip install stardate_code. Documentation is available at
https://stardate.readthedocs.io/en/latest/. All code used to produce the
figures in this paper is available at https://github.com/RuthAngus/stardate.
This paper is based on code with the following Git hash:
f739562c1546e117e9bb217e1732c62b41be8061.

Some of the data presented in this paper were obtained from the Mikulski
Archive for Space Telescopes (MAST). STScI is operated by the Association of
Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555.
Support for MAST for non-HST data is provided by the NASA Office of Space
Science via grant NNX09AF08G and by other grants and contracts. This paper
includes data collected by the Kepler mission. Funding for the Kepler mission
is provided by the NASA Science Mission directorate.

This work made use of the gaia-kepler.fun crossmatch database created by Megan
Bedell

This work has made use of data from the European Space Agency (ESA) mission
Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing
and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the Gaia Multilateral Agreement.

This research was supported in part by the National Science Foundation under
Grant No. NSF PHY-1748958.

##  5 Appendix

### Priors

We used the default priors in the isochrones.py python package. The prior over
age was,

$p(A)=\frac{\log(10)10^{A}}{10^{10.5}-10^{8}},~{}~{}~{}8<A<10.5.$ (17)

where $A$, is $\log_{10}(\mathrm{Age~{}[yrs]})$. The prior over EEP was
uniform with an upper limit of 800. We found that adding this upper limit
reduced some multi-modality caused by the giant branch and resulted in better
performance. The prior over true bulk metallicity was based on the galactic
metallicity distribution, as inferred using data from the Sloan Digital Sky
Survey (Casagrande et al., 2011). It is the product of a Gaussian that
describes the metallicity distribution over halo stars and two Gaussians that
describe the metallicity distribution in the thin and thick disks:

$\displaystyle p(F)=$ $\displaystyle
H_{F}\frac{1}{\sqrt{2\pi\sigma_{\mathrm{halo}}^{2}}}\exp\left(-\frac{(F-\mu_{\mathrm{halo}})^{2}}{2\sigma_{\mathrm{halo}}}\right)$
$\displaystyle\times(1-H_{F})\frac{1}{\xi}\left[\frac{0.8}{0.15}\exp\left(-\frac{(F-0.016)^{2}}{2\times
0.15^{2}}\right)+\frac{0.2}{0.22}\exp\left(-\frac{(F-0.15)^{2}}{2\times
0.22^{2}}\right)\right],$

where $H_{F}=0.001$ is the halo fraction, $\mu_{\mathrm{halo}}$ and
$\sigma_{\mathrm{halo}}$ are the mean and standard deviation of a Gaussian
that describes a probability distribution over metallicity in the halo, and
take values -1.5 and 0.4 respectively. The two Gaussians inside the square
brackets describe probability distributions over metallicity in the thin and
thick disks. The values of the means and standard deviations in these
Gaussians are from Casagrande et al. (2011). $\xi$ is the integral of
everything in the square brackets from $-\infty$ to $\infty$ and takes the
value $\sim 2.507$. The prior over distance was,

$p(D)=\frac{3}{3000^{3}}D^{2},~{}~{}~{}0<D<3000,$ (19)

with D in kiloparsecs, and, finally, the prior over extinction was uniform
between zero and one,

$p(A_{V})=U(0,1).$ (20)

## References

  * Agüeros et al. [2018] M. A. Agüeros, E. C. Bowsher, J. J. Bochanski, P. A. Cargile, K. R. Covey, S. T. Douglas, A. Kraus, A. Kundert, N. M. Law, A. Ahmadi, and H. G. Arce.  A New Look at an Old Cluster: The Membership, Rotation, and Magnetic Activity of Low-mass Stars in the 1.3 Gyr Old Open Cluster NGC 752.  _ApJ_ , 862(1):33, Jul 2018.  10.3847/1538-4357/aac6ed.
  * Agüeros [2018] Marcel Agüeros.  The evolution of cool dwarf spin rates: Data, models, and tensions.  In _Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun_ , Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, page 17, Jul 2018.  10.5281/zenodo.1343628.
  * Aigrain et al. [2015] S. Aigrain, J. Llama, T. Ceillier, M. L. das Chagas, J. R. A. Davenport, R. A. García, K. L. Hay, A. F. Lanza, A. McQuillan, T. Mazeh, J. R. de Medeiros, M. B. Nielsen, and T. Reinhold.  Testing the recovery of stellar rotation signals from Kepler light curves using a blind hare-and-hounds exercise.  _MNRAS_ , 450(3):3211–3226, Jul 2015.  10.1093/mnras/stv853.
  * Angus et al. [2015] R. Angus, S. Aigrain, D. Foreman-Mackey, and A. McQuillan.  Calibrating gyrochronology using Kepler asteroseismic targets.  _MNRAS_ , 450:1787–1798, June 2015.  10.1093/mnras/stv423.
  * Angus et al. [2018] Ruth Angus, Timothy Morton, Suzanne Aigrain, Daniel Foreman-Mackey, and Vinesh Rajpaul.  Inferring probabilistic stellar rotation periods using Gaussian processes.  _MNRAS_ , 474(2):2094–2108, Feb 2018.  10.1093/mnras/stx2109.
  * Barnes [2003] S. A. Barnes.  On the Rotational Evolution of Solar- and Late-Type Stars, Its Magnetic Origins, and the Possibility of Stellar Gyrochronology.  _ApJ_ , 586:464–479, March 2003.  10.1086/367639.
  * Barnes [2007] S. A. Barnes.  Ages for Illustrative Field Stars Using Gyrochronology: Viability, Limitations, and Errors.  _ApJ_ , 669:1167–1189, November 2007.  10.1086/519295.
  * Barnes [2010] S. A. Barnes.  A Simple Nonlinear Model for the Rotation of Main-sequence Cool Stars. I. Introduction, Implications for Gyrochronology, and Color-Period Diagrams.  _ApJ_ , 722:222–234, October 2010.  10.1088/0004-637X/722/1/222.
  * Brown et al. [2011] T. M. Brown, D. W. Latham, M. E. Everett, and G. A. Esquerdo.  Kepler Input Catalog: Photometric Calibration and Stellar Classification.  _AJ_ , 142:112, October 2011.  10.1088/0004-6256/142/4/112.
  * Casagrande et al. [2011] L. Casagrande, R. Schönrich, M. Asplund, S. Cassisi, I. Ramírez, J. Meléndez, T. Bensby, and S. Feltzing.  New constraints on the chemical evolution of the solar neighbourhood and Galactic disc(s). Improved astrophysical parameters for the Geneva-Copenhagen Survey.  _A &A_, 530:A138, Jun 2011.  10.1051/0004-6361/201016276.
  * Ceillier et al. [2016] T. Ceillier, J. van Saders, R. A. García, T. S. Metcalfe, O. Creevey, S. Mathis, S. Mathur, M. H. Pinsonneault, D. Salabert, and J. Tayar.  Rotation periods and seismic ages of KOIs - comparison with stars without detected planets from Kepler observations.  _MNRAS_ , 456(1):119–125, Feb 2016.  10.1093/mnras/stv2622.
  * Choi et al. [2016] J. Choi, A. Dotter, C. Conroy, M. Cantiello, B. Paxton, and B. D. Johnson.  Mesa Isochrones and Stellar Tracks (MIST). I. Solar-scaled Models.  _ApJ_ , 823:102, June 2016.  10.3847/0004-637X/823/2/102.
  * Christensen-Dalsgaard [2008a] J. Christensen-Dalsgaard.  ASTEC–the Aarhus STellar Evolution Code.  _Ap &SS_, 316:13–24, August 2008a.  10.1007/s10509-007-9675-5.
  * Christensen-Dalsgaard [2008b] J. Christensen-Dalsgaard.  ADIPLS–the Aarhus adiabatic oscillation package.  _Ap &SS_, 316:113–120, August 2008b.  10.1007/s10509-007-9689-z.
  * Curtis and Agüeros [2018] Jason Curtis and Marcel Agüeros.  Problems with and Prospects for K dwarf Gyrochronology: Insights from the K2 Survey of Ruprecht 147.  In _Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun_ , Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, page 24, Jul 2018.  10.5281/zenodo.1403719.
  * Curtis et al. [2019] Jason Lee Curtis, Marcel A. Agüeros, Stephanie T. Douglas, and Søren Meibom.  A Temporary Epoch of Stalled Spin-Down for Low-Mass Stars: Insights from NGC 6811 with Gaia and Kepler.  _arXiv e-prints_ , art. arXiv:1905.06869, May 2019.
  * Donahue et al. [1996] Robert A. Donahue, Steven H. Saar, and Sallie L. Baliunas.  A Relationship between Mean Rotation Period in Lower Main-Sequence Stars and Its Observed Range.  _ApJ_ , 466:384, Jul 1996.  10.1086/177517.
  * Dotter [2016] A. Dotter.  MESA Isochrones and Stellar Tracks (MIST) 0: Methods for the Construction of Stellar Isochrones.  _ApJS_ , 222:8, January 2016.  10.3847/0067-0049/222/1/8.
  * Dotter et al. [2008] A. Dotter, B. Chaboyer, D. Jevremović, V. Kostov, E. Baron, and J. W. Ferguson.  The Dartmouth Stellar Evolution Database.  _ApJS_ , 178:89–101, September 2008.  10.1086/589654.
  * Douglas et al. [2016] S. T. Douglas, M. A. Agüeros, K. R. Covey, P. A. Cargile, T. Barclay, A. Cody, S. B. Howell, and T. Kopytova.  K2 Rotation Periods for Low-mass Hyads and the Implications for Gyrochronology.  _ApJ_ , 822:47, May 2016.  10.3847/0004-637X/822/1/47.
  * Douglas et al. [2017] S. T. Douglas, M. A. Agüeros, K. R. Covey, and A. Kraus.  Poking the Beehive from Space: K2 Rotation Periods for Praesepe.  _ApJ_ , 842:83, June 2017.  10.3847/1538-4357/aa6e52.
  * Epstein and Pinsonneault [2014] C. R. Epstein and M. H. Pinsonneault.  How Good a Clock is Rotation? The Stellar Rotation-Mass-Age Relationship for Old Field Stars.  _ApJ_ , 780:159, January 2014.  10.1088/0004-637X/780/2/159.
  * Evans et al. [2017] D. W. Evans, M. Riello, F. De Angeli, G. Busso, F. van Leeuwen, C. Jordi, C. Fabricius, A. G. A. Brown, J. M. Carrasco, H. Voss, M. Weiler, P. Montegriffo, C. Cacciari, P. Burgess, and P. Osborne.  Gaia Data Release 1. Validation of the photometry.  _A &A_, 600:A51, April 2017.  10.1051/0004-6361/201629241.
  * Foreman-Mackey et al. [2013] D. Foreman-Mackey, D. W. Hogg, D. Lang, and J. Goodman.  emcee: The MCMC Hammer.  _PASP_ , 125:306, March 2013.  10.1086/670067.
  * Fossati et al. [2008] L. Fossati, S. Bagnulo, J. Landstreet, G. Wade, O. Kochukhov, R. Monier, W. Weiss, and M. Gebran.  The effect of rotation on the abundances of the chemical elements of the A-type stars in the Praesepe cluster.  _A &A_, 483:891–902, June 2008.  10.1051/0004-6361:200809467.
  * Gaia Collaboration et al. [2018] Gaia Collaboration, A. G. A. Brown, A. Vallenari, T. Prusti, J. H. J. de Bruijne, C. Babusiaux, C. A. L. Bailer-Jones, M. Biermann, D. W. Evans, L. Eyer, and et al.  Gaia Data Release 2. Summary of the contents and survey properties.  _A &A_, 616:A1, August 2018.  10.1051/0004-6361/201833051.
  * Gallet and Bouvier [2015] F. Gallet and J. Bouvier.  Improved angular momentum evolution model for solar-like stars. II. Exploring the mass dependence.  _A &A_, 577:A98, May 2015.  10.1051/0004-6361/201525660.
  * García et al. [2014] R. A. García, T. Ceillier, D. Salabert, S. Mathur, J. L. van Saders, M. Pinsonneault, J. Ballot, P. G. Beck, S. Bloemen, T. L. Campante, G. R. Davies, J.-D. do Nascimento, Jr., S. Mathis, T. S. Metcalfe, M. B. Nielsen, J. C. Suárez, W. J. Chaplin, A. Jiménez, and C. Karoff.  Rotation and magnetism of Kepler pulsating solar-like stars. Towards asteroseismically calibrated age-rotation relations.  _A &A_, 572:A34, December 2014.  10.1051/0004-6361/201423888.
  * Hidalgo et al. [2018] Sebastian L. Hidalgo, Adriano Pietrinferni, Santi Cassisi, Maurizio Salaris, Alessio Mucciarelli, Alessandro Savino, Antonio Aparicio, Victor Silva Aguirre, and Kuldeep Verma.  The Updated BaSTI Stellar Evolution Models and Isochrones. I. Solar-scaled Calculations.  _ApJ_ , 856(2):125, Apr 2018.  10.3847/1538-4357/aab158.
  * Huber et al. [2016] Daniel Huber, Stephen T. Bryson, Michael R. Haas, Thomas Barclay, Geert Barentsen, Steve B. Howell, Sanjib Sharma, Dennis Stello, and Susan E. Thompson.  The K2 Ecliptic Plane Input Catalog (EPIC) and Stellar Classifications of 138,600 Targets in Campaigns 1-8.  _The Astrophysical Journal Supplement Series_ , 224:2, May 2016.  10.3847/0067-0049/224/1/2.
  * Irwin and Bouvier [2009] J. Irwin and J. Bouvier.  The rotational evolution of low-mass stars.  In E. E. Mamajek, D. R. Soderblom, and R. F. G. Wyse, editors, _The Ages of Stars_ , volume 258 of _IAU Symposium_ , pages 363–374, June 2009.  10.1017/S1743921309032025.
  * Kawaler [1988] S. D. Kawaler.  Angular momentum loss in low-mass stars.  _ApJ_ , 333:236–247, October 1988.  10.1086/166740.
  * Kawaler [1989] S. D. Kawaler.  Rotational dating of middle-aged stars.  _ApJ_ , 343:L65–L68, August 1989.  10.1086/185512.
  * Krishnamurthi et al. [1997] A. Krishnamurthi, M. H. Pinsonneault, S. Barnes, and S. Sofia.  Theoretical Models of the Angular Momentum Evolution of Solar-Type Stars.  _ApJ_ , 480:303–323, May 1997.  10.1086/303958.
  * Mamajek and Hillenbrand [2008] E. E. Mamajek and L. A. Hillenbrand.  Improved Age Estimation for Solar-Type Dwarfs Using Activity-Rotation Diagnostics.  _ApJ_ , 687:1264–1293, November 2008.  10.1086/591785.
  * Mathur et al. [2012] S. Mathur, T. S. Metcalfe, M. Woitaszek, H. Bruntt, G. A. Verner, J. Christensen-Dalsgaard, O. L. Creevey, G. Doǧan, S. Basu, and C. Karoff.  A Uniform Asteroseismic Analysis of 22 Solar-type Stars Observed by Kepler.  _ApJ_ , 749(2):152, Apr 2012.  10.1088/0004-637X/749/2/152.
  * Matt et al. [2012] S. P. Matt, K. B. MacGregor, M. H. Pinsonneault, and T. P. Greene.  Magnetic Braking Formulation for Sun-like Stars: Dependence on Dipole Field Strength and Rotation Rate.  _ApJ_ , 754:L26, August 2012.  10.1088/2041-8205/754/2/L26.
  * Matt et al. [2015] S. P. Matt, A. S. Brun, I. Baraffe, J. Bouvier, and G. Chabrier.  The Mass-dependence of Angular Momentum Evolution in Sun-like Stars.  _ApJ_ , 799:L23, January 2015.  10.1088/2041-8205/799/2/L23.
  * McQuillan et al. [2014] A. McQuillan, T. Mazeh, and S. Aigrain.  Rotation Periods of 34,030 Kepler Main-sequence Stars: The Full Autocorrelation Sample.  _ApJS_ , 211:24, April 2014.  10.1088/0067-0049/211/2/24.
  * Meibom et al. [2015] S. Meibom, S. A. Barnes, I. Platais, R. L. Gilliland, D. W. Latham, and R. D. Mathieu.  A spin-down clock for cool stars from observations of a 2.5-billion-year-old cluster.  _Nature_ , 517:589–591, January 2015.  10.1038/nature14118.
  * Metcalfe et al. [2009] T. S. Metcalfe, O. L. Creevey, and J. Christensen-Dalsgaard.  A Stellar Model-fitting Pipeline for Asteroseismic Data from the Kepler Mission.  _ApJ_ , 699:373–382, July 2009.  10.1088/0004-637X/699/1/373.
  * Metcalfe et al. [2012] T. S. Metcalfe, S. Mathur, G. Doğan, and M. Woitaszek.  First Results from the Asteroseismic Modeling Portal.  In H. Shibahashi, M. Takata, and A. E. Lynas-Gray, editors, _Progress in Solar/Stellar Physics with Helio- and Asteroseismology_ , volume 462 of _Astronomical Society of the Pacific Conference Series_ , page 213, Sep 2012.
  * Metcalfe et al. [2014] T. S. Metcalfe, O. L. Creevey, G. Doğan, S. Mathur, H. Xu, T. R. Bedding, W. J. Chaplin, J. Christensen-Dalsgaard, C. Karoff, and R. Trampedach.  Properties of 42 Solar-type Kepler Targets from the Asteroseismic Modeling Portal.  _ApJS_ , 214(2):27, Oct 2014.  10.1088/0067-0049/214/2/27.
  * Morton [2015] T. D. Morton.  isochrones: Stellar model grid package.  Astrophysics Source Code Library, March 2015.
  * Pietrinferni et al. [2004] Adriano Pietrinferni, Santi Cassisi, Maurizio Salaris, and Fiorella Castelli.  A Large Stellar Evolution Database for Population Synthesis Studies. I. Scaled Solar Models and Isochrones.  _ApJ_ , 612(1):168–190, Sep 2004.  10.1086/422498.
  * Pinsonneault et al. [1989] M. H. Pinsonneault, S. D. Kawaler, S. Sofia, and P. Demarque.  Evolutionary models of the rotating sun.  _ApJ_ , 338:424–452, March 1989.  10.1086/167210.
  * Pont and Eyer [2004] Frédéric Pont and Laurent Eyer.  Isochrone ages for field dwarfs: method and application to the age-metallicity relation.  _MNRAS_ , 351(2):487–504, Jun 2004.  10.1111/j.1365-2966.2004.07780.x.
  * Rebull et al. [2016] L. M. Rebull, J. R. Stauffer, J. Bouvier, A. M. Cody, L. A. Hillenbrand, D. R. Soderblom, J. Valenti, D. Barrado, H. Bouy, D. Ciardi, M. Pinsonneault, K. Stassun, G. Micela, S. Aigrain, F. Vrba, G. Somers, J. Christiansen, E. Gillen, and A. Collier Cameron.  Rotation in the Pleiades with K2. I. Data and First Results.  _AJ_ , 152:113, November 2016.  10.3847/0004-6256/152/5/113.
  * Rebull et al. [2017] L. M. Rebull, J. R. Stauffer, L. A. Hillenbrand, A. M. Cody, J. Bouvier, D. R. Soderblom, M. Pinsonneault, and L. Hebb.  Rotation of Late-type Stars in Praesepe with K2.  _ApJ_ , 839:92, April 2017.  10.3847/1538-4357/aa6aa4.
  * Silva Aguirre et al. [2015] V. Silva Aguirre, G. R. Davies, S. Basu, J. Christensen-Dalsgaard, O. Creevey, T. S. Metcalfe, T. R. Bedding, L. Casagrande, R. Handberg, M. N. Lund, P. E. Nissen, W. J. Chaplin, D. Huber, A. M. Serenelli, D. Stello, V. Van Eylen, T. L. Campante, Y. Elsworth, R. L. Gilliland, S. Hekker, C. Karoff, S. D. Kawaler, H. Kjeldsen, and M. S. Lundkvist.  Ages and fundamental properties of Kepler exoplanet host stars from asteroseismology.  _MNRAS_ , 452:2127–2148, September 2015.  10.1093/mnras/stv1388.
  * Skumanich [1972] A. Skumanich.  Time Scales for CA II Emission Decay, Rotational Braking, and Lithium Depletion.  _ApJ_ , 171:565, February 1972.  10.1086/151310.
  * Soderblom [2010] D. R. Soderblom.  The Ages of Stars.  _ARA &A_, 48:581–629, September 2010.  10.1146/annurev-astro-081309-130806.
  * Somers et al. [2017] Garrett Somers, John Stauffer, Luisa Rebull, Ann Marie Cody, and Marc Pinsonneault.  M Dwarf Rotation from the K2 Young Clusters to the Field. I. A Mass-Rotation Correlation at 10 Myr.  _ApJ_ , 850(2):134, Dec 2017.  10.3847/1538-4357/aa93ed.
  * Tayar and Pinsonneault [2018] Jamie Tayar and Marc H. Pinsonneault.  Testing Angular Momentum Transport and Wind Loss in Intermediate-mass Core-helium Burning Stars.  _ApJ_ , 868(2):150, Dec 2018.  10.3847/1538-4357/aae979.
  * van Saders and Pinsonneault [2013] J. L. van Saders and M. H. Pinsonneault.  Fast Star, Slow Star; Old Star, Young Star: Subgiant Rotation as a Population and Stellar Physics Diagnostic.  _ApJ_ , 776:67, October 2013.  10.1088/0004-637X/776/2/67.
  * van Saders et al. [2016] J. L. van Saders, T. Ceillier, T. S. Metcalfe, V. Silva Aguirre, M. H. Pinsonneault, R. A. García, S. Mathur, and G. R. Davies.  Weakened magnetic braking as the origin of anomalously rapid rotation in old field stars.  _Nature_ , 529:181–184, January 2016.  10.1038/nature16168.
  * van Saders et al. [2018] J. L. van Saders, M. H. Pinsonneault, and M. Barbieri.  Forward Modeling of the Kepler Stellar Rotation Period Distribution: Interpreting Periods from Mixed and Biased Stellar Populations.  _ArXiv e-prints_ , March 2018.
  * Wright et al. [2011] N. J. Wright, J. J. Drake, E. E. Mamajek, and G. W. Henry.  The Stellar-activity-Rotation Relationship and the Evolution of Stellar Dynamos.  _ApJ_ , 743:48, December 2011.  10.1088/0004-637X/743/1/48.
  * Yi et al. [2001] S. Yi, P. Demarque, Y.-C. Kim, Y.-W. Lee, C. H. Ree, T. Lejeune, and S. Barnes.  Toward Better Age Estimates for Stellar Populations: The Y2 Isochrones for Solar Mixture.  _ApJS_ , 136:417–437, October 2001.  10.1086/321795.