THREE SUPER-EARTHS ORBITING HD 7924^{* †}

We collected 599 new RV measurements for HD 7924 over the last five years since the discovery of HD 7924b (H09) using the HIRES spectrograph on the Keck I telescope (Vogt et al. 1994). When this new data is combined with the data from H09 and APF we have over 10 yr of observational baseline (see Figure 1(a)). Our data collection and reduction techniques are described in detail in H09. On Keck/HIRES we observe the star through a cell of gaseous iodine in order to simultaneously forward model the instrumental line broadening function (PSF) and the subtle shifts of the stellar lines relative to the forest of iodine lines. The HIRES detector was upgraded in August of 2004 and the RV zero-point between the pre-and post-upgrade data may not necessarily be the same. For this reason, we allow for separate RV zero-points for the pre-upgrade, post-upgrade, and APF data sets. All RV measurements and associated Ca ii H and K ${{S}_{{\rm HK}}}$ activity indices (the ratio of the flux in the cores of the H and K lines to neighboring continuum levels) are listed in Table 2.

Zoom In Zoom Out Reset image size

Table 2.  Radial Velocities of HD 7924

${\rm BJ}{{{\rm D}}_{{\rm TDB}}}$	RV	Uncertainty	Instrumenta	${{S}_{{\rm HK}}}$
(–2440000)	(m s−1)	(m s−1)
12307.77162	5.23	1.31	k	⋯
12535.95639	3.24	1.16	k	⋯
13239.08220	−3.00	0.91	j	0.218
13338.79766	3.83	1.08	j	0.227
16505.99731	−7.97	1.46	a	⋯
16515.90636	3.88	1.63	a	⋯

Note.

^ak—pre-upgrade Keck/HIRES, j—post-upgrade Keck/HIRES, a—APF

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

The APF is a 2.4 m f/15 Cassegrain telescope built by Electro-optical Systems Technologies housed in an IceStorm-2 dome located at Lick Observatory atop Mount Hamilton, 20 miles east of San Jose, CA. The telescope operates completely unattended using a collection of Python, Tcl, bash, and csh scripts that interact with the lower level software operating on the Keck Task Library keyword system (Lupton & Conrad 1993). Every scheduled observing night, our automation software queries an online Google spreadsheet (that serves as our target database) and creates two observing plans: one for good conditions and the other for poor. Weather permitting, shortly after sunset the observatory opens automatically and determines whether to use the observing plan designed for good or poor conditions by monitoring a bright star. Once the observing plan is started, the high-level software continues to monitor the conditions and adjusts the schedule or switches observing plans if necessary. Observations of our targets (primarily for radial velocity measurements) continue throughout the night until conditions deteriorate too much to remain open or the until morning 9° twilight.

The Levy Spectrograph is a high-resolution slit-fed optical echelle spectrograph mounted at one of the two Nasmyth foci of the APF designed specifically for the detection and characterization of exoplanets (Burt et al. 2014; Radovan et al. 2014; Vogt et al. 2014b). Each spectrum covers a continuous wavelength range from 3740 to 9700 Å. We observed HD 7924 using a 10 wide decker for an approximate spectral resolution of R = 100,000. Starlight passes through a cell of gaseous iodine that serves as a simultaneous calibration source for the instrumental PSF and wavelength reference. In addition, we collected a high signal-to-noise ratio (S/N) spectrum through the 05 wide decker (R = 150,000) with the iodine cell out of the light path. This spectrum serves as a template from which we measure the relative Doppler shifts of the stellar absorption lines with respect to the iodine lines while simultaneously modeling the PSF and wavelength scale of each spectrum. APF guides on an image of the star before the slit using a slanted, uncoated glass plate that deflects 4% of the light to a guide camera. Photon-weighted times of mid-exposure are recorded using a software-based exposure meter based on guide camera images that monitors the sky-subtracted light entering the slit during an exposure (Kibrick et al. 2006).

We measure relative RVs using a Doppler pipeline descended from the iodine technique in Butler et al. (1996). For the APF, we forward-model 848 segments of each spectrum between 5000 and 6200 Å. The model consists of a stellar template spectrum, an ultra high-resolution Fourier transform spectrum of the iodine absorption of the Levy cell, a spatially variable PSF, a wavelength solution, and RV. We estimate RV uncertainties (Table 2) as the uncertainty on the mean RV from the large number of spectral segments. Well-established RV standard stars have a RV scatter of ∼2–3 m s$^{-1}$based on APF measurements. This long-term scatter represents a combination of photo-limited uncertainties, stellar jitter, and instrumental systematics. We collected a total of 109 RV measurements of HD 7924 on 80 separate nights (post outlier rejection) over a baseline of 1.3 yr with a typical per-measurement uncertainty of 2.0 m s⁻¹.

For both the APF and Keck data, Julian dates of the photon-weighted exposure mid-times were recored during the observations, then later converted to Barycentric Julian date in the dynamical time system (${\rm BJ}{{{\rm D}}_{{\rm TDB}}}$) using the tools of Eastman et al. (2010).⁸

We rejected a handful of low S/N spectra (S/N $\lt \;70$ per pixel) and measurements with uncertainties greater than nine times the median absolute deviation of all measurement uncertainties relative to the median uncertainty for each instrument. This removed a total of 11 RVs out of the 906 total measurements in the combined data set (Keck plus APF).

After outlier rejection we bin together any velocities taken less than 0.5 days apart on a single telescope. Since data taken in short succession are likely affected by the same systematic errors (e.g., spectrograph defocus), these measurements are not truly independent. Binning helps to reduce the effects of time-correlated noise by preventing multiple measurements plagued by the same systematic errors from being given too much weight in the Keplerian analysis. When the data are binned together and the uncertainty for the resulting data point is added in quadrature with the stellar jitter (an additional error term that accounts for both stellar and instrumental systematic noise), the binned data point receives only as much weight as a single measurement. While this likely reduces our sensitivity slightly we accept this as a tradeoff for more well-behaved errors and smoother ${{\chi }^{2}}$ surfaces. An independent analysis using the unbinned data in Section 4.4 finds three planets having the same orbital periods, eccentricities, and masses within 1σ as those discovered by analyzing the binned data.

We identify significant periodic signals in the RVs using an iterative multi-planet detection algorithm based on the two-dimensional Keplerian Lomb–Scargle (2DKLS) periodogram (O'Toole et al. 2009). Instead of fitting sinusoidal functions to the RV time series (Lomb 1976; Scargle 1982), we create a periodogram by fitting the RV data with Keplerian orbits at many different starting points on a two-dimensional (2D) grid over orbital period and eccentricity. This technique allows for relative offsets and uncertainties between different data sets to be incorporated directly into the periodogram and enhances the sensitivity to moderate and high eccentricity planets.

We fit the Keplerian models using the Levenberg–Marquardt (LM) ${{\chi }^{2}}$-minimization routine in the RVLIN IDL package (Wright & Howard 2009). Multi-planet models are sums of single-planet models, with planet–planet gravitational interactions neglected. Such an approximation is valid since the interaction terms are expected to be $\ll$1 m s⁻¹ for non-resonant, small planets. We define a grid of search periods following the prescription of Horne & Baliunas (1986) and at each period we seed an LM fit at five evenly spaced eccentricity values between 0.05 and 0.7. Period and eccentricity are constrained to intervals that allow them to vary only half the distance to adjacent search periods and eccentricities. All other model parameters are free to vary, including the parameters of any previously identified planets. The period and eccentricity for previously identified planets are constrained to be within ±5% and $_{-10}^{+5}$%, respectively, of their initial values but all other parameters are unconstrained. This prevents slightly incorrect fits of the first detected planets from injecting periodic residuals that could mimic further planetary signals. The 2DKLS periodogram power at each point in the grid is

$$\begin{eqnarray}&&Z(P,e)=\frac{{{\chi }^{2}}-\chi _{B}^{2}}{\chi _{B}^{2}},\end{eqnarray} \tag{ 1 }$$

where ${{\chi }^{2}}$ is the sum of the squared residuals to the current $N+1$ planet fit, and $\chi _{B}^{2}$ is the sum of the squared residuals to the best N-planet fit. In the first iteration of the planet search (comparing a one-planet model to a zero-planet model), $\chi _{B}^{2}$ is the sum of the squared error-normalized residuals to the mean (the B subscript stands for baseline). The 2D periodogram is collapsed into a one-dimensional periodogram as a function of period by taking the maximal Z for each period searched (i.e., the best fit eccentricity for every period).

We start the iterative planet search by comparing a zero-planet model (flat line) to a grid of one-planet models. A strong signal at a period of 5.4 days is detected at very high significance in the first iteration. This planet was initially discovered by H09 using only 22% of the data presented in this work. When we calculate two versus one planet and three versus two planet periodograms, we find two additional highly significant signals with periods of 15.3 and 24.5 days respectively (see Figure 1). These two periodic signals are best fit by Keplerian orbital models with semi-amplitudes of 2.3 and 1.7 m s⁻¹, respectively, with no significant eccentricity. See Tables 3 and 4 for the full orbital solution.

Table 3.  Orbital Parameters

Parameter	Value	Units
Modified DE-MCMC Step Parametersa
log(Pb)	0.73223 $\pm \;2e-05$	log(days)
$\sqrt{{{e}_{b}}}{\rm cos} {{\omega }_{b}}$	0.15 $_{-0.17}^{+0.13}$	⋯
$\sqrt{{{e}_{b}}}{\rm sin} {{\omega }_{b}}$	−0.09 $_{-0.15}^{+0.18}$	⋯
log(Kb)	0.555 $_{-0.024}^{+0.023}$	m s−1
log(Pc)	1.184663 $_{-9.3e-05}^{+9.2e-05}$	log(days)
$\sqrt{{{e}_{c}}}{\rm cos} {{\omega }_{c}}$	0.20 $_{-0.24}^{+0.17}$	⋯
$\sqrt{{{e}_{c}}}{\rm sin} {{\omega }_{c}}$	0.11 $_{-0.20}^{+0.17}$	⋯
log(Kc)	0.364 $_{-0.040}^{+0.037}$	m s$^{-1}$
log(Pd)	1.3883 $_{-0.00031}^{+0.00027}$	log(days)
$\sqrt{{{e}_{d}}}{\rm cos} {{\omega }_{d}}$	0.31 $_{-0.24}^{+0.16}$	⋯
$\sqrt{{{e}_{d}}}{\rm sin} {{\omega }_{d}}$	0.02 $_{-0.37}^{+0.30}$	⋯
log(Kd)	0.219 $_{-0.057}^{+0.052}$	m s$^{-1}$
Model Parameters
Pb	5.39792 ± 0.00025	days
${{T}_{{\rm conj},b}}$	2455586.38 $_{-0.110}^{+0.086}$	${\rm BJ}{{{\rm D}}_{{\rm TDB}}}$
eb	0.058 $_{-0.040}^{+0.056}$	⋯
${{\omega }_{b}}$	332 $_{-50}^{+71}$	degrees
Kb	3.59 $_{-0.19}^{+0.20}$	m s$^{-1}$
Pc	15.299 $_{-0.0033}^{+0.0032}$	days
${{T}_{{\rm conj},c}}$	2455586.29 $_{-0.47}^{+0.40}$	${\rm BJ}{{{\rm D}}_{{\rm TDB}}}$
ec	0.098 $_{-0.069}^{+0.096}$	⋯
${{\omega }_{c}}$	27 $_{-60}^{+52}$	degrees
Kc	2.31 $_{-0.20}^{+0.21}$	m s$^{-1}$
Pd	24.451 $_{-0.017}^{+0.015}$	days
${{T}_{{\rm conj},d}}$	2455579.1 $_{-0.9}^{+1.0}$	${\rm BJ}{{{\rm D}}_{{\rm TDB}}}$
ed	0.21 $_{-0.12}^{+0.13}$	⋯
${{\omega }_{d}}$	119 $_{-97}^{+210}$	degrees
Kd	1.65 ± 0.21	m s$^{-1}$
${{\gamma }_{{\rm post}-{\rm upgrade}\,{\rm Keck}}}$	−0.19 ± 0.16	m s$^{-1}$
${{\gamma }_{{\rm pre}-{\rm upgrade}\,{\rm Keck}}}$	2.0 $_{-1.2}^{+1.1}$	m s$^{-1}$
${{\gamma }_{{\rm APF}}}$	0.28 $_{-0.47}^{+0.46}$	m s$^{-1}$
${{\sigma }_{{\rm jitt}}}$	2.41 $_{-0.10}^{+0.11}$	m s$^{-1}$

Note.

^aMCMC jump parameters that were modified from the physical parameters in order to speed convergence and avoid biasing parameters that must physically be finite and positive.

Download table as:  ASCIITypeset image

We find a fourth, significant signal at ∼2400 days, but we also find a very similar signal (approximately sinusoidal with the same period and phase) upon inspection of the time series of ${{S}_{{\rm HK}}}$ stellar activity measurements (see Figure 2). Although we fit for this periodic signal as an additional Keplerian we do not interpret this as an additional planet. Instead we interpret this as the signature of the stellar magnetic activity cycle. We searched for a fifth periodic signal and found two additional marginally significant peaks at 17.1 and 40.8 days. Due to their marginal strengths and the fact that the two periods are related by the synodic month alias ($1/17.1\approx 1/40.8+1/29.5$) we are especially cautious in the interpretation of these signals. We scrutinize these candidate signals in depth in Sections 4.4 and 4.5. We also see a peak at around 40 days in a periodogram of the ${{S}_{{\rm HK}}}$ time series after the long-period signal from the stellar magnetic activity cycle is removed that is presumably the signature of rotationally modulated star spots. We conclude that the 17.1 and 40.8 day signals are most likely caused by rotational modulation of starspots, but the 5.4, 15.3, and 24.5 day signals are caused by three planetary companions with minimum (${{M}_{p}}{\rm sin} {{i}_{p}}$) masses of 8.7 ${{M}_{\oplus }}$, 7.9 ${{M}_{\oplus }}$, and 6.4 ${{M}_{\oplus }}$.

Zoom In Zoom Out Reset image size

We determine the orbital parameters and associated uncertainties for the three planet system using the ExoPy Differential-evolution Markov Chain Monte Carlo (DE-MCMC, Ter Braak 2006) engine described in Fulton et al. (2013) and Knutson et al. (2014). We treat the total RV model as a sum of Keplerian orbits each parameterized by orbital period (P_i), time of inferior conjunction (${{T}_{{\rm conj},i}}$), eccentricity (e_i), argument of periastron of the star's orbit (${{\omega }_{i}}$), and velocity semi-amplitude (K_i) where i is an index corresponding to each planet (b–d). An RV "jitter" term (${{\sigma }_{{\rm jitt}}}$) is added in quadrature with the measurement uncertainties at each step in the MCMC chains. We also fit for independent RV zero-points for the APF, pre-upgrade Keck, and post-upgrade Keck data. The long-period RV signal presumably caused by the stellar magnetic activity cycle is treated as an additional Keplerian orbit with the same free parameters as for each of the three planets. In order to speed convergence and avoid biasing parameters that must physically be finite and positive we step in the transformed and/or combinations of parameters listed in Table 3. We add a ${{\chi }^{2}}$ penalty for large jitter values of the following form:

$$\begin{eqnarray}&&\chi _{{\rm new}}^{2}={{\chi }^{2}}+2\mathop{\sum }\limits_{n}{\rm ln} (\sqrt{2\pi (\sigma _{{\rm vel},n}^{2}+\sigma _{{\rm jitt}}^{2})}),\end{eqnarray} \tag{ 2 }$$

where ${{\sigma }_{{\rm vel},n}}$ is the velocity uncertainty for each of the n measurements (Johnson et al. 2011).

We fit for the 23 free parameters by running 46 chains in parallel continuously checking for convergence using the prescription of Eastman et al. (2013). When the number of independent draws (T_z as defined by Ford 2006) is greater than 1000 and the Gelman–Rubin statistic (Gelman et al. 2003; Holman et al. 2006) is within 1% of unity for all free step parameters we halt the fitting process and compile the results in Table 3. We also list additional derived properties of the system in Table 4 that depend on the stellar properties listed in Table 1.

Table 4.  Derived Properties

Parameter	Value	Units
${{e}_{b}}{\rm cos} {{\omega }_{b}}$	0.0073 $_{-0.0076}^{+0.0210}$	⋯
${{e}_{b}}{\rm sin} {{\omega }_{b}}$	−0.0031 $_{-0.0190}^{+0.0063}$	⋯
ab	0.05664 $_{-0.00069}^{+0.00067}$	AU
${{M}_{b}}{\rm sin} {{i}_{b}}$	8.68 $_{-0.51}^{+0.52}$	${{M}_{\oplus }}$
Sba	113.7 $_{-3.6}^{+3.7}$	${{S}_{\oplus }}$
${{T}_{{\rm eq},b}}$ b	825.9 $_{-6.5}^{+6.6}$	K
${{e}_{c}}{\rm cos} {{\omega }_{c}}$	0.017 $_{-0.018}^{+0.051}$	⋯
${{e}_{c}}{\rm sin} {{\omega }_{c}}$	0.008 $_{-0.013}^{+0.037}$	⋯
ac	0.1134 $_{-0.0014}^{+0.0013}$	AU
${{M}_{c}}{\rm sin} {{i}_{c}}$	7.86 $_{-0.71}^{+0.73}$	${{M}_{\oplus }}$
Sca	28.35 $_{-0.89}^{+0.92}$	${{S}_{\oplus }}$
${{T}_{{\rm eq},c}}$ b	583.6 $_{-4.6}^{+4.7}$	K
${{e}_{d}}{\rm cos} {{\omega }_{d}}$	0.059 $_{-0.054}^{+0.084}$	⋯
${{e}_{d}}{\rm sin} {{\omega }_{d}}$	0.001 $_{-0.074}^{+0.076}$	⋯
ad	0.1551 $_{-0.0019}^{+0.0018}$	AU
${{M}_{d}}{\rm sin} {{i}_{d}}$	6.44 $_{-0.78}^{+0.79}$	${{M}_{\oplus }}$
Sda	15.17 $_{-0.48}^{+0.49}$	${{S}_{\oplus }}$
${{T}_{{\rm eq},d}}$ b	499 ± 4	K

Notes.

^aStellar irradiance received at the planet relative to the Earth. ^bAssuming a bond albedo of 0.32; the mean total albedo of super-Earth-size planets (Demory 2014).

Download table as:  ASCIITypeset image

We attempted to empirically determine the probability that Gaussian random noise in the data could conspire to produce an apparent periodic signal with similar significance to the periodogram peaks corresponding to each of the planets. We calculated 1000 2DKLS periodograms, each time scrambling the velocities in a random order drawn from a uniform distribution. We located and measured the height of the global maxima of each periodogram and compare the distribution of these maxima to the periodogram peak heights at the periods of the three planets in the original periodograms. The power for each periodogram within the set of 1000 is a ${\Delta }{{\chi }^{2}}$ between a two-planet plus activity model to a three-planet model assuming the other two planets as "known." Figure 3 shows that the distribution of maxima from the periodograms of the scrambled data are clearly separated from the original peaks. None of the trials produce periodogram peaks anywhere near the heights of the original peaks corresponding to the three planets. Because the 2DKLS periodogram allows eccentric solutions, we explore the scrambled RVs for non-sinusoidal solutions and are therefore sensitive to a wide variety of false alarm signals. We conclude that the false alarm probabilities for all three planets are $\lt 0.001$.

Zoom In Zoom Out Reset image size

The high cadence of the APF data set allows us to explore short-period orbital solutions. Traditional observations on large telescopes such as Keck often yield only a few nights of data per year, making it difficult to determine whether a short orbital period is an alias of a longer period, or a true physical signal. Although planets with orbital periods shorter than one day are uncommon in the galaxy (around 0.83 ± 0.18% of K dwarfs; Sanchis-Ojeda et al. 2014), 11 ultra-short period planets have been found⁹ due to their high detectability. To search as carefully as possible for short period signals, we use the unbinned data sets, which consist of 797 RVs from Keck and 109 RVs from the APF. Our use of the unbinned data in this section also provides independent confirmation of the results obtained with the binned data above.

Dawson & Fabrycky (2010) outline a rigorous procedure to distinguish between physical and alias periods. Our method for finding the orbital periods and distinguishing aliases is as follows:

• 1.  
  Determine the window function of the data to understand which aliases are likely to appear.
• 2.  
  Compute the periodogram of the data, determining the power and phase at each input frequency.
• 3.  
  If there is a strong peak in the periodogram, fit an N-planet Keplerian (starting with N = 1), using the periodogram peak as the trial period.
• 4.  
  Subtract the N-planet Keplerian from the data.
• 5.  
  Compute the periodogram of the Nth planet in the model and compare it to the periodogram of the Nth planet in the data minus the model of the other planets.
• 6.  
  If a second peak in the periodogram has similar height to the tallest peak and is located at an alias period, repeat steps 3–5 using that trial period.
• 7.  
  If you explored an alias period, choose the model that minimizes ${{\chi }^{2}}$ and best reproduces the observed periodogram. Subtract this model from the RVs.
• 8.  
  Treat the residuals as the new data set and go back to step 2. Examine the residuals from the N-planet fit for additional planets, and continue until there are no more signals in the periodogram.

The window functions of the individual and combined Keck and APF RV time series are shown in Figure 4. The window function is given by

$$\begin{eqnarray}&&W(\nu )=\frac{1}{N}\mathop{\sum }\limits_{j=1}^{N}{\rm exp} (-2\pi i\nu {{t}_{j}}),\end{eqnarray} \tag{ 3 }$$

where ν is the frequency in units of days⁻¹ and t_j is the time of the jth observation. The data sets are complementary: 109 RVs from APF over the past year are well-distributed over the months and the year, and so are only susceptible to the daily aliases, whereas the 797 RVs from Keck over the last decade are distributed in a way that gives some power to daily aliases as well as longer-period aliases. The power in the combined window function illustrates that we might be susceptible to daily aliases, and weak signals in the periodogram might even be susceptible to monthly or yearly aliases.

Zoom In Zoom Out Reset image size

To take the periodogram of the time series, we use a version of fasper (Press & Rybicki 1989) written for Python. We find the same peaks in the periodogram of the data and residuals at 5.4, 15.3, and 24.5 days that we interpret as planets. The periodograms of the data and periodograms of the Keplerian models associated with these periods are shown in Figure 5. Again, we recover a fourth peak at 2570 days that we attribute to long-term stellar activity due to the correlation between the RVs and the ${{S}_{{\rm HK}}}$ values. The three-planet plus stellar activity model found with the alias-search method yields periods, eccentricities, and planet masses within 1σ of the results quoted in Table 3 for all three planets. ${{\omega }_{b}}$ is consistent within 1σ, but inconsistent for planets c and d. However, the arguments of periastron for all three planets are poorly constrained due to their nearly circular orbits. We find a fifth peak at 40.8 days, which is also prominent in the periodogram of the stellar activity and is likely the rotation period of the star. The 40.8 days signal has a strong alias at 17.1 days and so we test Keplerian models at both 40.8 and 17.1 days to discriminate which is the true signal and which is the alias (see Figure 6). We find that the periodogram of the model 40.8 days signal better matches the alias structure in the periodogram of the data, and so we prefer 40.8 days as a candidate stellar rotation period.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Although this star is relatively inactive with $R{{^{\prime} }_{{\rm HK}}}$ values in the literature between −4.89 (Isaacson & Fischer 2010) and −4.85 (Canto Martins et al. 2011), some low-level chromospheric activity is detectable with our high-precision RVs. The most obvious feature is the long period signal with a period of ≈6.6 years and an amplitude of ≈5 m s$^{-1}$. Although this signal looks promising as a long-period sub-Jupiter-mass planet candidate, upon inspection of the ${{S}_{{\rm HK}}}$ time series we notice that this activity indicator is highly correlated with this long-period signal in the velocities (see Figure 2). We interpret this signal as the signature of the stellar magnetic activity cycle of which we have observed nearly two full cycles. If we subtract this long-period signal from the ${{S}_{{\rm HK}}}$ values and make a periodogram of the residuals we find a marginally significant peak at ≈41 days. This is most likely caused by rotational modulation of starspots, because 41 days is near the expected rotation period (38 days, Isaacson & Fischer 2010) for a star of this spectral type and age.

To gain confidence in the source of these Doppler signals (activity or planets) we employed a time-dependent RV–${{S}_{{\rm HK}}}$ decorrelation that probes evolution of the RV periodogram on timescales of the stellar rotation. Simple linear decorrelation is limited by the natural phase offset between RV and activity signals. At any given time, the measured RV and activity both depend on the flux-weighted fractional coverage of the magnetically active region on the stellar disk. The measured ${{S}_{{\rm HK}}}$ peaks when the magnetically active regions are closest to the disk center, where its projected area is largest and the star would otherwise be brightest. However, the RV approaches zero at the center of the stellar disk. The activity-induced RV shift is maximized when the active region is closer to the limb, ∼1/8 phase from center (e.g., Queloz et al. 2001; Aigrain et al. 2012). In addition, the magnetically active regions continuously evolve on timescales of several stellar rotation periods weakening the RV–${{S}_{{\rm HK}}}$ correlations over long timescales.

Motivated by these limitations, our alternative RV–${{S}_{{\rm HK}}}$ decorrelation method is effective on timescales of both the long-term magnetic activity cycle and the stellar rotation period. All APF data are excluded from this analysis because activity measurements are not presently available.

Starting with the same outlier-removed time-series as described above, we bin together any radial velocity and ${{S}_{{\rm HK}}}$ measurements taken on the same night. The period of the long-term stellar magnetic activity cycle is then identified from a Lomb–Scargle (LS) periodogram of the activity time-series. That period serves as the initial guess in a subsequent single-Keplerian fit to the RV time-series using rv_mp_fit. Next, planet candidates are identified using an iterative approach in which the number of fitted planets N increases by one starting at N = 1. $N+1$ Keplerians are fit at each step, one accounting for the long-term stellar magnetic activity cycle:

• 1.  
  Determine period of $N{\rm th}$ planet candidate: identify the period, P_N, of the highest peak in the LS periodogram of the residuals to the N-planet Keplerian fit to the RV time-series.
• 2.  
  Obtain best N-planet solution: perform an $(N+1)$-planet Keplerian fit with initial guesses being the N Keplerians from the previous fit plus a Keplerian with initial period guess P_N.
• 3.  
  Repeat: steps 1–2 for $N=N+1$.

We wish to subsequently remove residual RV signatures that can be confidently attributed to magnetically active regions rotating around the stellar surface. Since the lifespan of such regions is typically a few rotation periods, each measured RV was decorrelated with ${{S}_{{\rm HK}}}$ activity measurements within ±100 days. We chose ±100 day windows because such windows were long enough to allow us to detect the 41 day rotation period in each window, and short enough to capture long-term variability in the power spectrum on time-scales of hundreds of days. However, before discussing our decorrelation technique, we investigate the degree to which our RV and ${{S}_{{\rm HK}}}$ measurements are correlated and the reliability of our planet candidate signals.

Figure 7 shows a time-series of LS periodograms corresponding to all RVs measured within ±100 days. Each panel contains more than 100 individual columns corresponding to each periodogram derived from all data within ±100 days of the times for each column. These periodograms are only displayed at times when a minimum data cadence requirement (discussed below) is met and are grayed out elsewhere. Individual periodograms have been normalized independently to facilitate inter-epoch comparison of the relative distribution of power across all periods. In each panel from top to bottom, an additional planet signal has been removed from the data. The long-term magnetic activity cycle is removed from the data in all three panels. The fitted longterm magnetic activity-induced component has also been removed in all cases. After removing planet b (top panel), the majority of epochs/columns show significant power at 15.3 days, corresponding to planet c. The few exceptions are epochs with the lowest cadence. The higher cadence epochs also have significant power near 24.5 days, corresponding to planet d. The 24.5 day signal is more apparent after planet c is also removed (middle panel). The presence of 15.3 and 24.5 day peaks across all high-cadence epochs is consistent with the coherent RV signature of a planet.

Zoom In Zoom Out Reset image size

Figure 7 also reveals significant RV power at ∼30–50 days periods which can be attributed to magnetically active regions (i.e., spots and plages) rotating around the stellar surface at the ∼40.8 days stellar rotation period. The activity-induced RV shifts dominate the periodograms after the planets b, c, and d have been removed (bottom panel). A number of epochs also show substantial power at ∼17.1 days. This is unlikely to be the signature of a fourth planet as the period corresponds to the monthly synodic alias of the 40.8 day rotation period. Moreover, a system of planets with 15.3 and 17.1 day periods is unlikely to be dynamically stable.

The direct correlation between RV and stellar activity on stellar rotation timescales is further highlighted in Figure 8. The top panel is the ${{S}_{{\rm HK}}}$-analog of Figure 7 in that each column is the running periodogram of all ${{S}_{{\rm HK}}}$ measurements within ±100 days. The suite of 30–50 day peaks mirror those in the RV periodograms of Figure 7, indicating that magnetic activity is responsible for RV variation on these timescales.

Zoom In Zoom Out Reset image size

The observed ∼20 days range in the period of these activity-induced peaks might seem surprisingly large, given that differential stellar rotation is likely of order several days (e.g., the observed difference in rotation period at the Sun's equator and 60° latitude is ∼5 days). However, some of this variation could also be statistical noise, simply an artifact of the limited sampling of ∼40 days signals over a relatively small time baseline. To test this hypothesis, we generated a synthetic RV time series of the exact same cadence as our real Keck observations and attempted to recover an injected signal of period 40.8 days. We first injected our three-planet model, adding Gaussian noise with standard deviation 2.36 m s⁻¹, the median velocity error in the post-upgrade Keck dataset (see Table 6). We then superposed the 40.8 days sinusoid of semi-amplitude 1.46 m s⁻¹. The chosen semi-amplitude is the median of semi-amplitudes of all Keplerian fits to activity-induced RVs as determined by our decorrelation algorithm (see below). We applied the same planet-detection algorithm described above to these simulated data and subtracted three Keplerians from the best four-Keplerian model, leaving only the 40.8 days signal.

Table 6.  APF vs. Keck Radial Velocity Precision

Instrument	rms	Median Uncertaintya	${{N}_{{\rm obs}}}$ b	Mean Cadence	VNRc
	(m s$^{-1}$)	(m s$^{-1}$)		(days)	m s$^{-1}$yr−1
HD 7924
Pre-upgrade Keck/HIRES	1.20	2.64	7	110.4	0.66
Post-upgrade Keck/HIRES	2.50	2.36	281	13.5	0.48
APF/Levy	2.80	2.86	80	5.9	0.36
HD 10700 (τ Ceti)
Pre-upgrade Keck/HIRES	2.87	2.84	84	17.5	0.63
Post-upgrade Keck/HIRES	2.32	2.22	190	19.4	0.53
APF/Levy	2.16	2.12	66	7.0	0.30
HD 9407
Pre-upgrade Keck/HIRES	1.89	1.94	12	162.8	1.26
Post-upgrade Keck/HIRES	1.89	1.86	202	18.3	0.52
APF/Levy	2.34	2.30	104	4.5	0.26

Notes.

^aStellar jitter has been added in quadrature with the binned velocities such that the ${{\chi }^{2}}$ of the velocities with respect to their median value is 1.0. ^bBinned in units of 0.5 days. ^cSee Equation (4).

Download table as:  ASCIITypeset image

The bottom plot in Figure 8 shows the resulting 200 days-running-window periodogram, akin to the bottom plot in Figure 7. The ∼10–15 days fluctuation in the recovered period of the injected 40.8 days signal is consistent with activity-induced features in our real ${{S}_{{\rm HK}}}$ and RV data. We conclude that the varying power of the RV signal with periods of 30–50 days observed in the real data is a limitation of the observing cadence and not intrinsic to the star.

The strong correlation observed between RV and stellar activity indicates tremendous potential for decorrelation to validate our planet candidates. We apply a decorrelation algorithm that iterates through each RV measurement $v({{t}_{i}})$ where t_i is the time of the ${{i}_{{\rm th}}}$ measurement:

• 1.  
  Define subsets S_i and V_i, consisting of all ${{S}_{{\rm HK}}}$ and RV measurements within 100 days of time t_i.
• 2.  
  Define ${{n}_{S}}$ and ${{n}_{V}}$ as the number of datapoints in S_i and V_i respectively. If n_S and ${{n}_{V}}\gt {{n}_{{\rm min} }}$ then proceed to next step, otherwise revert to step 1 for $i=i+1$.
• 3.  
  Identify the period, P_i, of the highest peak in the LS periodogram of S_i.
• 4.  
  Perform a single-Keplerian fit to V_i with fixed period P_i.
• 5.  
  Subtract Keplerian fit from $v({{t}_{i}})$.
• 6.  
  Repeat for $i=i+1$.

As step 2 indicates, we only decorrelated an RV measurement with ${{S}_{{\rm HK}}}$ activity if both the number of ${{S}_{{\rm HK}}}$ and RV measurements within ±100 days, n_S and n_V respectively, exceeded some minimum number, ${{n}_{{\rm min} }}$. This was done to avoid removal of spurious correlations. ${{n}_{{\rm min} }}$ was chosen by examining the dependence of correlation significance on n_S and n_V. More specifically, for each RV measurement, we computed the Pearson product-moment correlation coefficient, r, between RV and ${{S}_{{\rm HK}}}$ periodograms from all data within ±100 days. In cases of n_S and n_V between ∼10–20, r values were widely distributed between 0 and 1, suggesting correlations were spurious. Supporting this notion, the best-fitting Keplerians to the activity-correlated RV noise had unusually large semi-amplitudes, as high as 10 m s⁻¹, even in cases where $r\gt 0.5$. In contrast, when n_S and n_V exceeded 23, r ranged from 0.34 to 0.76, with a median of 0.51 and semi-amplitudes were a more reasonable 1–3 m s⁻¹. We therefore adopted ${{n}_{{\rm min} }}$ = 23 resulting in activity decorrelation of 115 of 281 RV measurements.

The top panel of Figure 9, shows the evolution of the periodogram of the entire RV time series as the activity decorrelation scheme is cumulatively applied to data points in chronological order. That is, the columnar periodogram at each given time, t, corresponds to that of the entire RV time series, with activity decorrelation of all RV measurements before and including time t. Planets b and c have been removed as well as the long term magnetic activity. The bottom panel dispays the RV periodograms before and after activity decorrelation, corresponding to the leftmost and rightmost periodograms in the top panel respectively. After decorrelation the periodogram has significant power at 24.5, 35.1, and 40.8 days and several other more moderate peaks. However, after decorrelation (rightmost column), the 24.5 signal dominates, in support of its planetary origin, while most other peaks, including those at 40.8 and 35.1 days, have been drastically reduced, consistent with manifestations of stellar activity. The single exception is the peak at 17.1 days. While it is unclear why this peak does not vanish, we remind the reader that is the monthly synodic alias of the 40.8 day stellar rotation period, and a planet with such an orbital period would likely be dynamically unstable with planet c (P = 15.3 day). The three-planet result of this analysis is consistent with the other analyses previously described. This method of decorrelating stellar noise from RV measurements will be useful for distinguishing planets from stellar activity in other high-cadence RV planet searches.

Zoom In Zoom Out Reset image size

In order to compare the architecture of the HD 7924 system with other compact multi-planet systems we compiled a catalog of similar multi-planet systems from Kepler. We restricted the Kepler multi-planet systems to those with exactly three currently known transiting planets in the system,¹⁰ a largest orbital period less than or equal to 30 days, and all three planets must have radii smaller than 4 ${{R}_{\oplus }}$. This left a total of 31 systems from the sample of confirmed Kepler systems (Borucki et al. 2011; Batalha et al. 2013; Marcy et al. 2014; Rowe et al. 2014). Their architectures are presented in Figure 13.

Zoom In Zoom Out Reset image size

We notice nothing unusual about the architecture of HD 7924 when compared with the Kepler systems. Most of these multi-planet systems contain three planets with masses between 5 and 10 ${{M}_{\oplus }}$ and semimajor axis between 0.05 and 0.3 AU. The systems of Kepler-194, Kepler-124, Kepler-219, Kepler-372, Kepler-310, Kepler-127, and Kepler-342 all host one inner planet and two outer planets that are closer to each other than to the innermost planet as in the HD 7924 system. The Kepler-372 system is particularly similar with one planet orbiting at 0.07 AU and a pair of outer planets orbiting at 0.14 and 0.19 AU. This demonstrates that the HD 7924 system is not abnormal in the context of other known compact multi-planet systems. While this is not a direct demonstration of dynamical stability, the fact that many other systems exist with very similar architectures gives strong empirical evidence that the HD 7924 system is a stable planetary configuration.

The HD 7924 planetary system is one of very few RV-detected systems hosting three super Earths. Only HD 40307 (Mayor et al. 2009a), and HD 20794 (Pepe et al. 2011) are known to host three planets with masses all below 10 ${{M}_{\oplus }}$. Both of these stars are unobservable from most northern hemisphere observatories. If we expand the mass limit to include systems with three or more planets with masses below 25 ${{M}_{\oplus }}$we find four additional systems that match these criteria, HD 69830 (Lovis et al. 2006), Gl 581 (Mayor et al. 2009b), HD 10180 (Lovis et al. 2011), and 61 Virginis (Vogt et al. 2010). These systems are difficult to detect given their small and complex RV signals and HD 7924 is the first discovery of such a system since 2011.

In order to check that the HD 7924 system is dynamically stable for many orbital periods we ran a numerical integration of the three planet system using the MERCURY code (Chambers 1999). We started the simulation using the median orbital elements presented in Table 3 and let it proceed 10⁵ years into the future. We assume that the system is perfectly coplanar with zero mutual inclination and we assume that the masses of the planets are equal to their minimum masses (${\rm sin} i=0$). No close passages between any of the three planets were found to occur ($\leqslant$1 Hill radius) during the entire simulation, suggesting that the system is stable in the current configuration.

Following the arguments of Fabrycky et al. (2014) we calculated the sum of separations between the planets (b–c, plus c–d) in units of their Hill radii. Fabrycky et al. (2014) found that the sum of the separations for the vast majority of Kepler multi-planet systems is larger than 18 Hill radii. The sum of the separations between the planets of the HD 7924 system is 36 Hill radii, providing further evidence that this system is stable and not abnormally compact relative to the many Kepler multi-planet systems.