Binary Companions of Evolved Stars in APOGEE DR14: Search Method and Catalog of ∼5000 Companions

For every source in the sample of APOGEE stars defined in Section 2, we obtain a posterior sampling in binary-system parameter space, treating it as a single-lined (SB1) spectroscopic binary system with a single companion. This sampling is performed with The Joker (Price-Whelan et al. 2017) under a relatively uninformative prior pdf, and the resulting posterior samplings are used to discover and characterize individual binary star systems and generate a catalog of companions. We now describe the assumptions and method used to generate samplings for individual systems:

• no multiple companions: All variations in radial velocity of the primary star are induced by a single companion. This is motivated by the idea that triple star systems are usually hierarchical so that the period of the inner binary is typically much shorter than the orbital period of the outer body (e.g., Tokovinin 2018). At present, we ignore the possibility of higher-order multiple systems.
• Kepler: All variations in velocity of the primary are gravitational, and we therefore ignore the possibility of coherent intrinsic variation from, e.g., stellar oscillations.
• SB1: All spectra are single-lined; that is, we assume that the secondary is significantly fainter and is thus undetected in the spectra. This assumption is motivated by the fact that we expect RGB stars to be substantially more luminous than their typical companion. However, there are known main-sequence double-lined binary stars in the APOGEE catalog (El-Badry et al. 2018), and an expected but unknown fraction of RGB–RGB binaries.
• simple noise model: Measurements are unbiased and noise estimates are correct up to an unknown excess variance. All noise contributions result in Gaussian uncertainties on the individual radial velocity measurements.

The Joker is a custom-built Monte Carlo sampler designed to produce independent posterior samples in Keplerian orbital parameters, given radial velocity measurements under the assumptions listed above. Our parameterization of the orbital elements is similar to the notation in Murray & Correia (2010): the radial velocity v at time t is given by

$$\begin{eqnarray}&&v(t;{\boldsymbol{\theta }})={v}_{0}+K\,[\cos (\omega +f(t;e,P,{M}_{0},{t}_{0}))+e\,\cos \omega ]\end{eqnarray} \tag{ 1 }$$

where ${\boldsymbol{\theta }}=(P,e,{M}_{0},\omega ,K,{v}_{0})$—period, eccentricity, mean anomaly at a reference time t₀, argument of pericenter, velocity semi-amplitude, barycentric velocity—and the true anomaly, f, is a function of time and the specified parameters (see Section 2 of Price-Whelan et al. 2017 or Equation (63) in Murray & Correia 2010). In addition to the orbital parameters listed above, The Joker can also generate samples in an "excess variance" parameter, s², that is added to the per-visit measurement variances. This parameter allows us to test whether the visit RV uncertainties are underestimated For upper RGB stars the inferred excess variance will be a combination of extra systematic uncertainty and true astrophysical surface jitter (e.g., Hekker et al. 2008).

The Joker was designed for the extremely multimodal pdfs expected when the number of radial velocity measurements of a source is small, or the data are sparse (in phase coverage) or noisy. While other Markov chain Monte Carlo (MCMC) methods have difficulty producing independent samples with such data, The Joker succeeds by brute force: after generating an initial (very large) library of prior samples from an assumed prior pdf (see below), the (typically multimodal) likelihood is evaluated at each sample and used to reject all prior samples above the likelihood surface (the rejection sampling step). In practice, given a number of requested samples for each star, the sampling proceeds iteratively: since it is easier to accept samples when the data are sparse or noisy, far more prior sample draws (and thus likelihood evaluations) must occur under very constraining data.

Here we describe the specifics of generating posterior samplings for all of the APOGEE targets in our parent sample. We execute the full procedure twice for different goals (as described in Section 5), and only here outline the key steps in the pipeline.

For all 96,231 APOGEE stars with ≥3 good visits (see Section 2), we use The Joker to generate posterior samplings for each star under the assumptions listed above (see Section 3.1). We start by generating a library of 536,870,912 prior samples generated under a prior similar to that defined in Price-Whelan et al. (2017):

• 1.  
  uniform or isotropic in angle parameters,
• 2.  
  uniform in log-period over the domain [1, 32768] day,
• 3.  
  a beta distribution over eccentricity (using parameters from Kipping 2013).

For the excess variance parameter, s², we use a Gaussian over the transformed parameter $y=\mathrm{ln}{s}^{2}$ with the mean and standard deviation (μ_y, σ_y) indicated where the runs are described (see Section 5). The reference time for each star is set to the first visit epoch; M₀ then becomes the mean anomaly at the first visit observation. Table 1 contains descriptions of all parameters and priors used.

Table 1.  Summary and Description of Parameters

Name	Prior	Description
P	$\mathrm{ln}P\sim { \mathcal U }(1,32768)\,\mathrm{day}$	period
e	e ∼ Beta(0.867, 3.03)	eccentricity
t0	fixed	reference time
M0	${M}_{0}\sim { \mathcal U }(0,2\pi )\,\mathrm{rad}$	mean anomaly at reference time
ω	$\omega \sim { \mathcal U }(0,2\pi )\,\mathrm{rad}$	argument of pericenter
s2	$\mathrm{ln}{s}^{2}\sim { \mathcal N }({\mu }_{y},{\sigma }_{y}^{2})$	extra variance added to each visit variance
K	${ \mathcal N }(0,{\sigma }_{v}^{2})\,\mathrm{km}\,{{\rm{s}}}^{-1}$	velocity semi-amplitude
v0	${ \mathcal N }(0,{\sigma }_{v}^{2})\,\mathrm{km}\,{{\rm{s}}}^{-1}$	system barycentric velocity

Note. Beta(a, b) is the beta distribution with shape parameters (a, b), ${ \mathcal U }(a,b)$ the uniform distribution over the domain (a, b), and ${ \mathcal N }(\mu ,{\sigma }^{2})$ is the normal distribution with mean μ and variance σ². For the systemic velocity and semi-amplitude, (v₀, K), we use a broad Gaussian prior that is formally inconsistent between The Joker and follow-up MCMC sampling (see Section 3.2.1): in The Joker we assume Gaussian priors with σ_v much larger than the measurement uncertainty, so they can be neglected (${\sigma }_{v}\approx \infty$), whereas when running MCMC we fix ${\sigma }_{v}={10}^{3}\,\mathrm{km}\,{{\rm{s}}}^{-1}$.

Download table as:  ASCIITypeset image

We request 256 posterior samples for each source. Depending on the data quality and phase coverage of the visits, The Joker will require different numbers of prior samples in order to rejection-sample down to the requested number of posterior samples: for few-epoch or noisy RV data, many prior samples will pass the rejection step, whereas for very precise or many-epoch RV data, The Joker may need to process the full library of prior samples. We therefore generate the posterior samples using an iterative procedure that adaptively predicts how many prior samples to test for each star. For sources with very constraining data, The Joker may return fewer than the requested number of samples (as few as one sample). When this occurs, we continue sampling either using standard MCMC or by increasing the size of the prior cache and continuing rejection sampling with The Joker.

3.2.1. "Needs MCMC": Following up The Joker with MCMC

If just one posterior sample is returned after exhausting the full library of prior samples, or if multiple (but fewer than 256) are returned that all lie within a single mode of the posterior pdf, the posterior pdf over orbital parameters is treated as effectively unimodal: these stars are flagged "needs MCMC." In this case, we use the location of the returned sample (if only one is returned) or a randomly chosen sample from those returned (if multiple samples are returned within one mode) to generate a small Gaussian ball of initial conditions and use standard MCMC to continue sampling until we obtain 256 samples.

In detail, we use an ensemble MCMC sampling algorithm (Goodman & Weare 2010) implemented in Python (emcee; Foreman-Mackey et al. 2013) to perform the samplings. emcee uses an ensemble of "walkers" (separate Markov chains) that naturally adapt to the geometry of the parameter space by generating proposal steps based on the locations of other walkers (i.e., the ensemble chains are not independent). We transform the standard Keplerian orbital parameters to a safer parameterization, $(\mathrm{ln}P,\sqrt{K}\,\cos ({M}_{0}+\omega )$, $\sqrt{K}\,\sin ({M}_{0}+\omega )$, $\sqrt{e}\,\cos \omega$, $\sqrt{e}\,\cos \omega$, $\mathrm{ln}{s}^{2},{v}_{0})$, for MCMC sampling. This reparameterization is safer and more efficient for sampling with emcee, which expects parameters to be components of a vector so that linear operations can be applied (see, e.g., Hogg & Foreman-Mackey 2018); the angle variables (ω, M₀) do not meet this requirement in the standard parameterization. We use the same prior pdfs as in The Joker when running MCMC (see Table 1).

For each star that is flagged "needs MCMC," we run emcee with 1024 walkers for 16,384 steps, take the final walker positions, and downsample at random until we have 256 samples. We compute the Gelman–Rubin convergence statistic, ${\hat{R}}_{j}$, (Gelman & Rubin 1992) for each parameter j and include these values in the catalogs below when standard MCMC is run. We also provide a binary flag, "converged," for each sampling continued with MCMC that is set to true if

$$\begin{eqnarray}&&\mathop{\mathrm{mean}}\limits_{j}({\hat{R}}_{j})\lt 1.1.\end{eqnarray} \tag{ 2 }$$

3.2.2. "Needs More Prior Samples": Extending The Joker Sampling

We construct a comparison sample of simulated data with no RV variability to assess our false-positive rates in the selections below (see Section 5). We randomly pick 16,384 stars from the parent sample used in this work and replace the visit velocity measurements with simulated data. For each star n, we randomly sample an excess variance parameter value from the prior, s_n. For each visit k, we then sample a new velocity v_nk by drawing from a Gaussian with mean equal to the mean of the real data visit velocities, ${\bar{v}}_{n}$, and variance equal to the sum of the visit variance (uncertainty), σ_nk, and the excess variance,

$$\begin{eqnarray}&&{v}_{{nk}}\sim { \mathcal N }({\bar{v}}_{n},{\sigma }_{{nk}}^{2}+{s}_{n}^{2}).\end{eqnarray} \tag{ 3 }$$

When we save the comparison data, we only store the visit uncertainty and "forget" the fact that the simulated data are generated with excess variance (to be inferred with The Joker).

We do not expect a sharp transition in inferred orbital parameters between stars with and without companions: there is a continuum of companion properties. For example, the companions can have low mass or long periods, or be at high inclination, all of which will make the velocity semi-amplitude, K, small for a given system. Therefore, there is no simple cut on radial velocity data alone that would select a complete sample of stars with companions. It is nevertheless possible to define a cut that selects stars with high-confidence companions.

We select a sample of stars that confidently have companions using percentiles computed from the log of the posterior samples in the velocity semi-amplitude parameter, $\mathrm{ln}K$. Figure 4 shows the distribution of first percentiles of $\mathrm{ln}K$ computed for all stars with posterior samplings from The Joker (filled, blue histogram), along with the same for posterior samplings for the null comparison sample (solid, dark line). To select stars with companions, we use a cut in the first percentile of the $\mathrm{ln}K$ samples such that 1% of the comparison sample is selected, i.e., our estimated false-positive rate is ≈1%. We use a threshold of $\mathrm{ln}K=-0.2$ to meet this constraint (vertical, dashed line in Figure 4); 4898 stars pass this cut. For brevity below, we refer to this sample as the "high-K" stars, and the complementary sample as the "low-K" stars.

Zoom In Zoom Out Reset image size

Table 2 contains a list of all stars in the parent sample and the value of the first percentile over the $\mathrm{ln}K$ samples for each star. Figures 5 and 6 show eight examples of high-K stars, i.e., stars that likely have companions, with different numbers of visits, N, indicated on the panels. Table 3 contains a list of all APOGEE_IDs for stars in the high-K sample. Figures 7 and 8 shows the same for eight low-K stars, i.e., stars that have either no companions, low-mass companions, or companions at long periods or high inclination. The majority of the stars in the high-K sample have poorly constrained orbital parameters because the posterior pdfs are very multimodal. The high-K sample selected from Table 2 includes all stars in the catalogs described in the next two subsections, but the following two catalogs are mutually exclusive.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Starting from the high-K sample, we select stars that have posterior period samples that fall within a single mode. We define a period resolution for each star ${\rm{\Delta }}=4{P}_{min}^{2}/(2\pi T)$, where P_min is the minimum period sample value and T is the epoch span of the data for that star. We consider a sampling to be unimodal when

$$\begin{eqnarray}&&{P}_{\max }-{P}_{\min }\lt {\rm{\Delta }}\end{eqnarray} \tag{ 4 }$$

(see Section 2 of Price-Whelan et al. 2017); 320 stars pass this cut. Unlike the multimodal stars, the posterior pdfs for these stars can be approximated using point estimates and standard deviations of their respective samplings. We report maximum a posteriori (MAP) sample values for the orbital parameters, along with other computed quantities for all 320 stars in this high-K, unimodal sample. For stars with measured primary masses, M₁, from other work (Ness et al. 2015), we compute the minimum companion mass, ${M}_{2,\min }$, and include both masses in this catalog. We additionally join with the APOGEE DR14 allStar catalog and provide all columns from this catalog for convenience. Table 4 contains descriptions and units for all columns in the high-K, unimodal sample catalog, available online.¹⁸

Table 4.  High-K, Unimodal Systems

Column name	Unit/format	Description
APOGEE_ID		identifier used by APOGEE
P	days	P, period
P_err	days
M0	rad	M0, phase at reference epoch
M0_err	rad
e		e, eccentricity
e_err
omega	rad	ω, argument of pericenter
omega_err	rad
jitter	km s−1	s, excess variance parameter
jitter_err	km s−1
K	km s−1	K, velocity semi-amplitude
K_err	km s−1
v0	km s−1	v0, systemic velocity
v0_err	km s−1
t0	barycentric MJD	t0, reference epoch
converged		binary flag indicating whether the sampling converged
Gelman–Rubin		Gelman–Rubin statistic for each MCMC parameter
M1	M⊙	primary mass estimate (Ness et al. 2015)
M1_err	M⊙
M2_min	M⊙	M2,min, minimum M2 mass
M2_min_err	M⊙
clean_flag		[0 = good, 1 = suspicious, 2 = bad], score from by-eye vetting
q_min		minimum mass ratio
q_min_err
R1	${R}_{\odot }$	radius of the primary star
R1_err	${R}_{\odot }$
a_sini	au	projected separation
a_sini_err	au
a2_sini	au	projected semimajor axis of the companion orbit
a2_sini_err	au
DR14RC	True/False	if the star is included in the APOGEE DR14 red clump catalog
TINGRC	True/False	if the star is included in red clump catalog of Ting et al. (2018)
⋮		all columns from Ness et al. (2015)
⋮		all columns from APOGEE DR14 allStar file
(320 rows)

Note. Description of the data table containing summary information for all stars in the high-K, unimodal sample: stars that likely have companions and have well-determined orbital parameters. All orbital parameter values are from the maximum a posteriori (MAP) posterior sample (either from The Joker or from emcee). All columns ending in _err are estimates of the standard deviation of the posterior samples, σ, computed using the median absolute deviation, MAD, as σ ≈ 1.5 × MAD.

(This table is available in its entirety in FITS format.)

Download table as:  ASCIITypeset image

We have also visually inspected the inferred orbits for all systems in this sample and have flagged systems with questionable or invalid fits. The value of the flag is: 0, when the orbits look reasonable; 1, when the orbits look reasonable but the value of the inferred excess variance is large; and 2, when the orbits are clearly poor fits. Many of the stars flagged as "2" look like they may be triple systems, because they tend to have variations in radial velocity over two distinct timescales that are never well-fit by a single Keplerian orbit. The stars flagged as "1" tend to be either upper RGB stars, where astrophysical surface jitter can lead to RV modulations, or SB2 systems, where absorption lines from the secondary confuse the RV pipeline and lead to strange RV signals. This flag is included in the catalog (Table 4) as the column clean_flag; to select a clean sample of companions that have been successfully vetted by eye, select only stars with clean_flag == 0.

A subset of the high-K sample (Section 5.1) have multimodal posterior distributions over orbital parameters that are limited to a few qualitatively different solutions. For these stars, one or a few more RV measurements would likely lead to uniquely determined orbital parameters, making these stars prime candidates for follow-up efforts. As examples of such cases, we include an additional catalog of systems that have effectively bimodal posterior samplings in orbital period. We identify these systems using k-means clustering (Lloyd 1982) of the posterior samples in period. In detail, for all stars in the high-K sample, we compute $\mathrm{ln}P$ for all period samples, then use k-means clustering with k = 2 as implemented in the scikit-learn package (Pedregosa et al. 2011) to identify two clusters of samples. We initialize the cluster positions at either end of the range of possible periods. For each cluster, we then ask whether the samples assigned to that cluster are unimodal (Equation (4)). If the samples in each respective cluster are unimodal, we call the sampling "bimodal"; 106 systems are identified as bimodal.

Figure 9 shows a few examples of systems that meet this criterion. In many of these cases, there are certain times at which a future observation would be far more informative. Visually, those times correspond to periods when the predicted RVs have large variance based on the posterior samples at hand. For example, in the bottom left panel of Figure 9, an observation at t = 60 days would rule out one of the two possible classes of orbits.

Zoom In Zoom Out Reset image size

Table 5 contains a list of all APOGEE targets with bimodal samplings identified in this way. The table contains two rows for each source, with the values of period, eccentricity, and velocity semi-amplitude from the MAP sample from each period mode. When available, this table also includes an estimate of the primary mass, M₁, from Ness et al. (2015) and an estimate of the minimum companion mass, M_2,min, for each mode.

We again visually inspect all systems in this sample and assign a flag based on the apparent quality of the inferred orbits (see Section 5.2). Again, to select a clean sample of companions from this catalog that have been successfully vetted by eye, select only stars with clean_flag == 0.

Figure 10 shows the companion and stellar properties of the unimodal and bimodal samples: the upper left panel shows inferred period and eccentricity of the binary orbit, the lower left panel shows inferred period and surface gravity of the primary star, and the lower right panel shows the stellar parameters of the primary star.

Zoom In Zoom Out Reset image size

The period–eccentricity plot (upper left panel in Figure 10) shows a clear signature of tidal circularization: such plots for main-sequence stars typically show a more distinct circularization period, below which the majority of stars have eccentricities close to zero. Here, we see that close to P ≈ 10 days the majority of companion orbits appear to be circular, but between 10 and 100 days there appears to be a gradual decrease in eccentricity rather than a sharp transition. This is likely because the primary stars in this sample have a large range in surface gravities and therefore a large range in stellar radii (see lower right panel); this result is explored in more detail in a companion work (Price-Whelan & Goodman 2018).

In the lower left panel of Figure 10, the diagonal (thick) line shows the orbital period of a hypothetical massless companion at the surface of a $1.35\,{M}_{\odot }$ primary star (the median of our sample with measured masses) with the given surface gravity. The vertical line in this panel shows the same for a $1.35\,{M}_{\odot }$ primary with $\mathrm{log}g=0$, i.e., the orbital period at the surface of such a star close to the tip of the RGB. The horizontal line in this panel shows the approximate upper bound of the red clump, above which most stars should be fully convective. Interestingly, there is a dearth of short-period systems for stars above the red clump ($\mathrm{log}g\lesssim 2.3$) in the triangle defined by the three qualitative lines, suggestive of possible companion engulfment as the primary star envelope encases the companion.

The lower right panel of Figure 10 shows the stellar parameters of all primary stars in the unimodal and bimodal samples. Relative to Figure 2, upper RGB stars appear to be under-represented in these companion catalogs, another tentative signature of engulfment or depletion of companions on shorter-period orbits.

Most of the companions we find have minimum masses ${M}_{2,\min }\lt 1\,{M}_{\odot }$. The left panel of Figure 11 shows estimates of primary and companion masses for systems with unimodal (circle markers) and bimodal (square markers, blue lines) period samplings for the subset of 69/320 unimodal and 25/106 bimodal systems with primary mass estimates, again using primary masses from Ness et al. (2015). The upper dashed line shows the curve M_2,min = M₁, and the lower dashed line shows the hydrogen-burning limit, ${M}_{2,\min }\,=0.08\,{M}_{\odot }$. For the systems with bimodal samplings, we compute the minimum companion mass for each period mode and connect these estimates with a line. Here we restrict the selection to stars with $\mathrm{log}g\gt 2$ to avoid upper RGB stars that may have surface oscillations that mimic RV modulations from companions.

Zoom In Zoom Out Reset image size

The right panel of Figure 11 shows the ratio of the primary stellar radius to the projected separation of the two bodies, ${R}_{1}/(a\,\sin i)$, and the minimum mass ratio, q_min = M_2,min/M₁, for the same systems. Again, circle (black, gray) markers indicate systems from the unimodal sample, and square (blue) markers and lines connect the two estimates for systems with bimodal period samplings. The upper dashed line in this panel indicates the Roche radius, computed using an approximate functional form (Eggleton 1983)

$$\begin{eqnarray}&&\displaystyle \frac{R}{a}=\displaystyle \frac{0.49\,{q}^{-2/3}}{0.6\,{q}^{-2/3}+\mathrm{ln}(1+{q}^{-1/3})}\end{eqnarray} \tag{ 5 }$$

where $q={M}_{2}/{M}_{1}$. All of these systems are consistent with being detached binaries, with the exception of one notable system (discussed below).

In Figure 11, a few interesting systems are immediately obvious: from the left panel, two systems have bimodal period samplings in which both period modes put the minimum companion mass below the hydrogen-burning limit (brown dwarf candidates), and two systems have M_2,min > M₁ (neutron star or black hole candidates). In the right panel, one system appears right on the limit of being an interacting binary (q_min ≈ 0.5), but has no significant UV flux (GALEX DR5; Bianchi et al. 2011). These systems are prime candidates for spectroscopic follow-up.

Paradoxically, there appear to exist short-period (≲100 day) companions at low surface gravities ($\mathrm{log}g\lesssim 1$), obvious in the lower left panel of Figure 10: these companions would orbit within the surface of their primary stars. Each of these systems appears fine from the APOGEE data quality flags, but could nonetheless have incorrect stellar parameters. As a test, we would ideally be able to compare the spectroscopic stellar parameters to an independent determination from, e.g., asteroseismology. Unfortunately, none of these short-period, low-$\mathrm{log}g$ stars in the unimodal or bimodal samples appears in the APOKASC (Serenelli et al. 2017) catalog, which provides asteroseismic stellar parameters for a few hundred APOGEE dwarf and giant stars. We therefore instead cross-match all high-K stars with $\mathrm{log}g\lt 1$ and 99% of their period samples below 100 days to the APOKASC catalog: one star appears in both samples, 2M19024490+4419523. In APOGEE DR14, this star has $\mathrm{log}g=0.57$, but asteroseismic $\mathrm{log}g\gt 4$ (e.g., HUBER_LOGG), consistent with being an M dwarf. Many of these systems could therefore be low-mass dwarf stars with incorrect stellar parameters in APOGEE DR14.

Another possibility is that The Joker interprets semi-coherent surface oscillations (from asteroseismic modes) with sparse time coverage as orbital RV variations (see, e.g., Hekker et al. 2008). However, for RGB stars with $\mathrm{log}g\lt 1$, these modes would likely have frequencies of ν_max ≈0.5–5 μHz (Garcia & Stello 2018), corresponding to periods of P ≈ 20–2 days and amplitudes of ≈20–200 m s⁻¹ (Huber et al. 2011; Huber 2017). Except for stars with the lowest $\mathrm{log}g$, these amplitudes would not pass our cut on K (Section 5.1).

These systems would need precise, long-term photometric follow-up to obtain asteroseismic parameters to confirm or explain their existence.

It is important to keep in mind that the companion catalogs and results presented in this article depend on the assumptions laid out in Section 3.1. We have only considered radial velocity modulations from a single massive companion—no multiple companions. This is a fundamental limitation of our search and methodology. However, stars with multiple companions will likely be included in our companion catalog anyway, but with two-body orbital solutions that either do not fit the data well or pick out the shortest-period orbit.

Related to the above, we assume that all RV variations are gravitational—Kepler. We see tentative evidence for surface oscillations at the upper giant branch (see Figure 10, lower left), which we treat at present as excess variance or intrinsic jitter. Further observations are needed to determine whether these "systems" have incorrect stellar parameters, or indeed show surface oscillations related to asteroseismic modes.

We assume that one star in each two-body system overwhelmingly dominates the luminosity and therefore the spectrum of each system—SB1. As seen in Figure 11, there are some companions with minimum masses comparable to or consistent with being larger than the masses of the observed stars. These systems are either (a) nearly edge-on MS–RGB binaries (i.e., consistent with our assumption), (b) SB2 systems where the APOGEE pipeline failed to flag the source as having broad lines or a bad fit, or (c) RGB–stellar remnant systems with a black hole or neutron star companion. A subset of these systems that look like they fall into category (c) are being followed-up to obtain further RV measurements to test whether any companions are black holes. However, The Joker can and should be extended to support sampling for SB2 systems with just one additional parameter (typically the mass ratio, or the ratio of velocity semi-amplitudes).

Finally, we have assumed that the visit velocity error distribution is Gaussian, and that the visit velocity uncertainties could be underestimated—simple noise model. We can test this assumption using posterior samples from The Joker by computing the residuals away from our best-fit two-body orbital solutions. We find that the normalized residuals appear to be very close to Gaussian over a large dynamic range of normalized residual values, indicating that this assumption may be sufficient. Figure 12 shows the distribution of normalized residuals, R_nk, for each visit: the k visit velocities for each n star, v_nk, minus the predicted radial velocity from the best-fitting sample returned by The Joker, ${\hat{v}}_{{nk}}^{* }$, normalized by the excess-variance-included uncertainty, ${\sigma }_{{nk}}^{* }=\sqrt{{\sigma }_{{nk}}^{2}+{\hat{s}}_{n}^{2}}$, where ${\hat{s}}_{n}$ is computed from the excess variance parameter of the best-fitting sample:

$$\begin{eqnarray}&&{R}_{{nk}}=\displaystyle \frac{{v}_{{nk}}-{\hat{v}}_{{nk}}^{* }}{\sqrt{{\sigma }_{{nk}}^{2}+{\hat{s}}_{n}^{2}}}.\end{eqnarray} \tag{ 6 }$$

Zoom In Zoom Out Reset image size

There are at least two other recent catalogs of stellar systems and companions based on APOGEE data.

One of these studies focused on decomposing spectra of MS stars into mixtures of stellar spectra (El-Badry et al. 2018). Conceptually, this method works because (a) for two stars of unequal mass, unexpected absorption lines will appear superimposed on the brighter star's spectrum, and (b) for stars of close to equal mass, the line depths and ratios will not be well-matched by a single stellar model. Using this technique, they identified thousands of candidate MS binaries and trinaries, but did not consider giant stars ($\mathrm{log}g\gt 4$). This sample is therefore complementary to and non-overlapping with the catalog presented in this work.

The other recent catalog searched for stellar and substellar companions to all stars in APOGEE DR12 (including the RGB) that passed a series of quality cuts and had $\geqslant 8$ visits (Troup et al. 2016). For each star in the sample, orbits were fit to the visit RVs using a multi-step orbit-fitting procedure: it starts by identifying significant periods and a few harmonics of those periods, then fits a Keplerian orbit at each of these harmonics using least-squares fitting (De Lee et al. 2013) with a modified χ² statistic that penalizes fits in which the phase coverage of the data is poor. This procedure is not guaranteed to provide a unique orbit solution.

Of the 382 companions released as a part of the search of Troup et al. (2016), only 188 of the host stars passed the stellar parameter and quality cuts used to define the parent sample in this work (see Section 2). We have looked at all of the overlapping stars to compare the previously derived companion orbital properties to the posterior samplings derived with The Joker. We find that the comparisons fall into three categories: (1) the parameters reported in Troup et al. (2016) agree with the posterior samplings, and the period distribution appears unimodal, (2) the parameters reported in Troup et al. (2016) identify one possible mode of a likely multimodal posterior pdf over orbital parameters, and (3) the radial velocity data used in Troup et al. (2016) changed significantly between APOGEE DR12 and DR14, so no meaningful comparisons can be made; roughly 1/3 of the comparison sample falls into each class. Figure 13 shows a few representative cases in which the orbital parameters of Troup et al. (2016) (orbit shown as an orange line in the left panels, parameters shown as orange + markers in the right panels) are consistent with the orbit samples from The Joker. Figure 14 shows a few representative cases in which we find that the posterior pdf over orbital parameters is multimodal, and the orbit of Troup et al. (2016) identifies one of these modes. For completeness, Figure 15 shows two instances in which the orbital parameters from Troup et al. (2016) no longer make sense, likely because the data changed between data releases.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The main motivation for this work was to produce posterior samplings for all APOGEE DR14 stars to be used in a population inference. To use these samplings for population or hierarchical inference, the orbits of each individual system do not have to be unimodal and thus all 96,231 samplings can be used in conjunction without well-determined orbital parameters for the majority of the stars. These samplings will be useful for constraining (through hierarchical inference) the period and eccentricity distributions of evolved stars, and the occurrence rates of companions to evolved stars to test for signatures of engulfment.