An 11 Earth-mass, Long-period Sub-Neptune Orbiting a Sun-like Star

To date, nearly four thousand exoplanets have been discovered, but over three quarters of them orbit their host star with periods of less than 50 days (NASA Exoplanet Archive;²⁵ accessed 2019 April 13). However, this is the result of observational biases rather than a feature of the underlying exoplanet population. Bias to short periods is especially strong for the transit method, the most common method of exoplanet detection. Nevertheless, Petigura et al. (2018) finds that from 1 to 24 R_⊕, the planet occurrence rate either increases or plateaus as a function of period out to many hundreds of days. Therefore, despite the estimated abundance of long-period planets (i.e., planets with periods longer than 50 days²⁶ ), our understanding of them is still very incomplete. Relative to the short-period population, there are very few long-period exoplanets (particularly in the low-mass regime) with precise and accurate densities and compositions, and even fewer with atmospheric characterization.

Thus, a larger sample of masses for long-period planets would allow us to address a number of interesting questions. For example, it would allow us to study the rocky to gaseous planet transition and how it depends on stellar flux. We could also investigate planet compositions in or near the habitable zone of Sun-like stars.

Another interesting feature to study would be the planet radius occurrence gap detected by Fulton et al. (2017) and Fulton & Petigura (2018). Owen & Wu (2017) and Van Eylen et al. (2018) have proposed that photoevaporation strips planets near their host stars down to the core, thus creating the gap. Lopez & Rice (2018) have investigated the period dependence of the gap position and Zeng et al. (2017) have analyzed the relationship between gap position and stellar type. More long-period planets, with or without planet masses, would provide new insights into the nature and cause of this radius occurrence gap.

In this paper, we characterize the long-period exoplanet Kepler-538b, the only known planet in the Kepler-538 system, first validated by Morton et al. (2016). There is a possible second transiting planet candidate with a period of 117.76 days, but its existence is very much in question; we briefly discuss this candidate in Section 5.3. We determine the properties of the host star, a G-type star slightly smaller than the Sun. We also determine properties of the exoplanet including the orbital period, mass, radius, and density by modeling transit photometry, radial velocity (RV) data, and stellar activity indices. We find that Kepler-538b is the smallest long-period planet to date with both a measured radius and RV mass.

The format of this paper is as follows. In Section 2, we detail our photometric and spectroscopic observations of the planet and its host star. We then discuss stellar parameterization in Section 3 and modeling of photometry and spectroscopy in Section 4. Our results are then presented and discussed in Section 5. Finally, we summarize and conclude our paper in Section 6.

Photometric observations of the Kepler-538 system were collected with the Kepler spacecraft (Borucki et al. 2008) across 17 quarters beginning in 2009 May and ending in 2013 May. Kepler collected both long-cadence and short-cadence observations of this system. Short-cadence observations (in quarters 3, 7–12 and 17) were collected every 58.89 s, and long-cadence observations (in all other quarters) were collected every 1765.5 s (∼29.4 minutes). In particular, we used pre-search data conditioning (PDC) light curves from these quarters downloaded from the Mikulski Archive for Space Telescopes.

Although Kepler-538 was not validated until Morton et al. (2016), it was flagged as a Kepler Object of Interest well before that. As a result, we have conducted a great deal of spectroscopic follow-up on Kepler-538 since it was identified as a candidate host star by the Kepler mission.

First, we collected two spectra with the Tillinghast Reflector Echelle Spectrograph (TRES; Fűrész 2008), an R = 44,000 spectrograph on the 1.5 m Tillinghast reflector at the Fred Lawrence Whipple Observatory (located on Mt. Hopkins, Arizona). These spectra were collected on the nights of 2010 May 28 and 2010 July 5 and had exposure times of 12 and 15 minutes, respectively.

We also downloaded RVs from 26 spectra collected with the High Resolution Echelle Spectrometer (HIRES) instrument (Vogt et al. 1994) at the Keck I telescope from 2010 July 25 to 2014 July 11. These spectra were originally collected as part of the Kepler Follow-up Observing Program. The standard California Planet Search setup was used (Howard et al. 2010) and the C2 decker was utilized to conduct sky subtraction. Exposure times averaged 1800 s.

Finally, we gathered 83 spectra with the High Accuracy Radial velocity Planet Searcher in North hemisphere (HARPS-N) instrument (Cosentino et al. 2012, 2014) on the 3.6 m Telescopio Nazionale Galileo (TNG) on La Palma. These observations were made from 2014 June 20 to 2015 November 7, all with exposure times of 30 minutes. They were collected as part of the HARPS-N Collaboration's Guaranteed Time Observations (GTO) program. Using the technique described in Malavolta et al. (2017), we confirmed that none of these spectra suffered from Moon contamination.

Stellar atmospheric parameters (effective temperature, metallicity, and surface gravity) were determined in two different ways. First, we combined the two TRES spectra and used the Stellar Parameter Classification tool (SPC; Buchhave et al. 2012). SPC compares an input spectrum against a library grid of synthetic spectra from Kurucz (1992), interpolating over the library to find the best match as well as uncertainties on the relevant stellar parameters. This method provides a measure for the rotational velocity as well.

Second, we used ARES+MOOG on the combination of our 83 HARPS-N spectra. More details about this method, based on equivalent widths (EWs), are found in Sousa (2014) and references therein. In short, ARESv2 (Sousa et al. 2015) automatically calculates the EWs of a set of neutral and ionised iron lines (Sousa et al. 2011). These are then used as input in MOOG²⁷ (Sneden 1973), assuming local thermodynamic equilibrium and using a grid of ATLAS plane-parallel model atmospheres (Kurucz 1993). Following Sousa et al. (2011), we added systematic errors in quadrature to our errors. The value for surface gravity was corrected for accuracy following Mortier et al. (2014). The results from SPC and ARES+MOOG agreed well within uncertainties.

We then estimated stellar mass, radius, and thus density with the isochrones package, a Python routine for inferring model-based stellar properties from known observations (Morton 2015). We supplied the spectroscopic effective temperature, metallicity, the Gaia DR2 parallax (Gaia Collaboration et al. 2016, 2018), and multiple photometric magnitudes (B, V, J, H, K, W1, W2, W3, and G) as input. Note that we did not use the surface gravity as an input parameter as this parameter is not well determined spectroscopically (e.g., Mortier et al. 2014). We ran isochrones four times, using the two different sets of spectroscopic parameters and two sets of isochrones, Modules for Experiments in Stellar Astrophysics (MESA) Isochrones and Stellar Tracks (MIST) and Dartmouth.²⁸

All four results were consistent, so we followed Malavolta et al. (2018) and derived our final set of parameters and uncertainties from the 16th, 50th, and 84th percentile values of the combined posteriors, minimizing systematic biases from using different spectroscopic methods or isochrones. The results of this analysis are listed in Table 1.

Table 1.  Stellar Parameters of Kepler 538

Parameter	Unit	SPC	ARES+MOOG	Combined
Stellar parameters
Effective temperature Teff	K	5547 ± 50	5522 ± 72	⋯
Surface gravity log g	g cm−2	4.51 ± 0.10	4.55 ± 0.12	⋯
Metallicity [m/H]	dex	−0.03 ± 0.08	⋯	⋯
Metallicity [Fe/H]	dex	⋯	−0.15 ± 0.05	⋯
Radius R*	R⊙	${0.8707}_{-0.0060}^{+0.0063}$	${0.8727}_{-0.0062}^{+0.0063}$	${0.8717}_{-0.0061}^{+0.0064}$
Mass M*	M⊙	${0.925}_{-0.036}^{+0.034}$	0.870 ± 0.024	${0.892}_{-0.035}^{+0.051}$
Density ρ*	ρ⊙	${1.404}_{-0.068}^{+0.061}$	1.31 ± 0.052	${1.349}_{-0.0716}^{+0.089}$
Distance	pc	${156.67}_{-0.70}^{+0.71}$	${156.65}_{-0.68}^{+0.70}$	${156.66}_{-0.69}^{+0.71}$
Age	Gyr	${3.8}_{-2.0}^{+2.1}$	${6.7}_{-1.6}^{+1.8}$	${5.3}_{-3.0}^{+2.4}$
Projected rotational velocity v sin i	km s−1	1.1 ± 0.5	⋯	⋯

Download table as:  ASCIITypeset image

As a useful check, we find that our estimates of stellar effective temperature, stellar radius, and distance are all within 1σ of the Gaia DR2 revised Kepler stellar parameters (Berger et al. 2018).

3.1. Consistency with Stellar Activity and Gyrochronology

3.2. Possible Binarity of Kepler-538

Our analysis of photometric and spectroscopic data included a simultaneous fit to both data types. Therefore, we first describe the data reduction process and model components of photometry and spectroscopy separately, then discuss the combined model afterward.

4.1. Photometric Data

We cleaned and reduced the photometric Kepler data using the lightkurve Python package (Barentsen et al. 2019). Each quarter was cleaned and reduced separately. For a given quarter, observation times without a corresponding flux were removed. Then, a crude light curve model based on the exoplanet parameters reported in the NASA Exoplanet Archive³⁰ (accessed 2019 February 16) was subtracted from the light curve so that in-transit data would not be clipped or flattened out in the next steps. Next, we flattened the light curve using the lightkurve flatten function, which uses a Savitzky–Golay filter. A window length of 615 or 41 was selected (i.e., 615 or 41 consecutive data points) for short-cadence and long-cadence data respectively, which is approximately three times the ratio between the transit duration and the observation cadence. Then, we clipped outlier data points discrepant from the median flux by more than 5σ. Lastly, we added the transit model from the earlier step back to the light curve. The reduced data can be seen in Figure 1, plotted in time and also phase-folded to the period of Kepler-538b.

Zoom In Zoom Out Reset image size

We modeled the light curve with the BATMAN Python package (Kreidberg 2015), which is based on the Mandel & Agol (2002) transit model. The model included a baseline offset parameter, a white noise parameter (to allow for instrumental and systematic noise in the data), two quadratic limb-darkening parameters (using the Kipping 2013a parameterization), the transit time (i.e., reference epoch), orbital period, planet radius relative to stellar radius, transit duration, impact parameter, eccentricity, and longitude of periastron.

We assumed uniform, Jeffreys, or modified Jeffreys priors for most of the parameters in this model, which are listed in Table 2. A Jeffreys prior is less informative than a uniform prior when the prior range is large and the scale of the parameter is unknown. A modified Jeffreys prior has the following form (Gregory 2007):

$$\begin{eqnarray*}&&p(X)=\displaystyle \frac{1}{X+{X}_{0}}\displaystyle \frac{1}{\mathrm{ln}\left(\tfrac{{X}_{\max }+{X}_{0}}{{X}_{\min }+{X}_{0}}\right)},\end{eqnarray*}$$

where X_min and X_max are the minimum and maximum prior value and X₀ is the location of a knee in the prior. A modified Jeffreys prior behaves like a Jeffreys prior above the knee at X₀ and behaves likes a uniform prior below the knee; this is useful when the prior includes zero (creating an asymptote for a conventional Jeffreys prior). A Jeffreys prior is simply a modified Jeffreys prior with the knee at X₀ = 0.

Table 2.  Transit and RV Parameters of Kepler 538b

Parameter	Unit	This Paper	Priors
Transit parameters
Period P	day	81.73778 ± 0.00013	Unif(81.73666, 81.73896)
Time of first transit	BJD-2454833	${211.6789}_{-0.0011}^{+0.0010}$	Unif(211.6671, 211.6901)
Orbital eccentricity e	⋯	${0.041}_{-0.029}^{+0.034}$ (<0.11)a	Beta(0.867, 3.03)b,c
Longitude of periastron ω	degree	${140}_{-90}^{+140}$	Unif(0, 360)
Impact parameter b	⋯	${0.41}_{-0.21}^{+0.10}$	Unif(0, 1)
Transit duration t14	hr	${6.62}_{-0.13}^{+0.21}$	Unif(0, 24)
Radius ratio ${R}_{{\rm{p}}}/{R}_{* }$	⋯	${0.02329}_{-0.00033}^{+0.00039}$	Jeffreys(0.001, 1)
Quadratic limb-darkening parameter q1	⋯	${0.164}_{-0.042}^{+0.067}$	Unif(0, 1)
Quadratic limb-darkening parameter q2	⋯	${0.74}_{-0.22}^{+0.16}$	Unif(0, 1)
Normalized baseline offset	ppm	$-{2.1}_{-2.8}^{+2.7}$	Unif(−100, 100)
Photometric white noise amplitude	ppm	${112.2}_{-2.4}^{+2.5}$	ModJeffreys(1, 1000, 234)
RV parameters
Semi-amplitude K	m s−1	${1.69}_{-0.38}^{+0.39}$	ModJeffreys(0.01, 10, 2.1)
HIRES RV white noise amplitude	m s−1	${3.25}_{-0.48}^{+0.56}$	ModJeffreys(0, 10, 2.1)
HARPS-N RV white noise amplitude	m s−1	${2.24}_{-0.27}^{+0.29}$	ModJeffreys(0, 10, 2.1)
HARPS-N FWHM white noise amplitude	m s−1	${6.71}_{-0.46}^{+0.52}$	Jeffreys(0.01, 10)
HIRES RV offset amplitude	m s−1	$-{0.50}_{-0.87}^{+0.78}$	Unif(−5, 5)
HARPS-N RV offset amplitude	m s−1	$-{37322.07}_{-0.73}^{+0.58}$	Unif(−37330, −37315)
HARPS-N FWHM offset amplitude	m s−1	${6655.4}_{-8.6}^{+7.5}$	Unif(6600, 6700)
GP RV convective blueshift amplitude Vc	m s−1	${0.86}_{-0.54}^{+0.75}$	ModJeffreys(0, 15, 2.1)
GP RV rotation modulation amplitude Vr	m s−1	${4.0}_{-3.0}^{+5.7}$	ModJeffreys(0, 15, 2.1)
GP FWHM amplitude Fc	m s−1	${13.3}_{-4.9}^{+5.9}$	Jeffreys(0.01, 25)
GP stellar rotation period P*	day	${25.2}_{-1.2}^{+6.5}$ d	Unif(20, 40)
GP inverse harmonic complexity λp	⋯	${5.2}_{-2.5}^{+2.8}$	Unif(0.25, 10)
GP evolution timescale λe	day	${370}_{-140}^{+200}$	Jeffreys(1, 1000)
Derived parameters
Planet radius Rp	R⊕	${2.215}_{-0.034}^{+0.040}$	⋯
System scale $a/{R}_{* }$	⋯	${87.5}_{-1.6}^{+1.5}$	⋯
Planet semimajor axis a	au	${0.4669}_{-0.0090}^{+0.0087}$	⋯
Orbital inclination i	degree	${89.73}_{-0.06}^{+0.14}$	⋯
Planet mass Mp	${M}_{\oplus }$	${10.6}_{-2.4}^{+2.5}$	⋯
Planet mean density ρp	ρ⊕	0.98 ± 0.23	⋯
Planet mean density ρp	g cm−3	5.4 ± 1.3	⋯
Planet insolation flux Sp	${S}_{\oplus }$	${2.99}_{-0.27}^{+0.31}$	⋯
Planet equilibrium temperature Teq (albedo = 0.3)	K	380	⋯
Planet equilibrium temperature Teq (albedo = 0.5)	K	350	⋯

Notes.

^a95% confidence limit. ^bBeta distribution parameter values from Kipping (2013b). ^cPrior also truncated to exclude e > 0.95. ^dRotation period uncertainties are highly asymmetric because the posterior includes a large peak at 25 days and a smaller peak at 31 days.

Download table as:  ASCIITypeset image

The only parameter with a different prior was orbital eccentricity. We applied a beta prior to orbital eccentricity using the values recommended by Kipping (2013b); we also truncated the prior to exclude e > 0.95.

Additionally, we also applied a stellar density prior. This was done given the fact that stellar density can be measured in two distinct ways: from photometry for a transiting exoplanet and from a stellar spectrum combined with stellar evolutionary tracks (we used the latter method in Section 3). Specifically, stellar density can be calculated via the following equation (Seager & Mallén-Ornelas 2003; Sozzetti et al. 2007):

$$\begin{eqnarray}&&{\rho }_{* }=\displaystyle \frac{3\pi }{{{GP}}^{2}}{\left(\displaystyle \frac{a}{{R}_{* }}\right)}^{3},\end{eqnarray} \tag{ 1 }$$

where the orbital period (P) and the normalized semimajor axis (a/R_*) are exoplanet properties that can be derived from the light curve. We applied a Gaussian prior to the exoplanet-derived stellar density using the density (and corresponding uncertainties) derived from spectra and stellar evolutionary tracks.

4.2. RV Data

Our RV analysis of Kepler-538b included not only the RV values determined from our HIRES and HARPS-N spectra, but also a number of indicators of stellar activity estimated from these spectra. For HARPS-N, these included the cross-correlation function (CCF) bisector span inverse slope (hereafter BIS), the CCF FWHM, and $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$. Our data reduction was performed with the data reduction software (DRS) 3.7 HARPS-N pipeline which applied a G2 stellar type mask. For HIRES, RVs are estimated with an iodine cell rather than cross correlation, so $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$ was calculated but not BIS or FWHM.

The RV and FWHM observations (and the corresponding model fit) can be seen in Figure 2. Additionally, all RV, FWHM, BIS, and $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$ values are listed in Table 3. There is a clear long-term trend in the FWHM observations (and to a lesser extent in the BIS and $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$ observations). However, we could not determine whether these trends have a stellar or instrumental origin, nor why there are no similar trends in the RV observations. On the one hand, when we checked three standard stars observed by HARPS-N during the same period of time, only one showed a similar FWHM trend. On the other hand, a FWHM trend in HARPS-N observations was also reported by Benatti et al. (2017) due to a defocusing problem, but that issue was corrected in 2014 March, before our first HARPS-N observations began. Still, perhaps a similar but slower and smaller drift affected our observations.

Zoom In Zoom Out Reset image size

Table 3.  RV Observations and Activity Indicators, Determined from the DRS

BJD	RV	RV Error	FWHM	BIS	${\mathrm{log}}_{10}({R}_{{HK}}^{{\prime} })$	${\mathrm{log}}_{10}({R}_{{HK}}^{{\prime} })$ Error	Instrument
	(m s−1)	(m s−1)	(km s−1)	(km s−1)
2455402.854339	−8.78	1.32	⋯	⋯	−5.078	⋯	HIRES
2455414.971547	−1.24	1.33	⋯	⋯	−5.003	⋯	HIRES
2455486.859621	−0.83	1.50	⋯	⋯	−4.985	⋯	HIRES
2455544.719659	−5.95	2.14	⋯	⋯	−4.971	⋯	HIRES
2455760.087400	4.58	1.39	⋯	⋯	−4.982	⋯	HIRES
2455796.934228	−3.56	1.25	⋯	⋯	−4.974	⋯	HIRES
2455797.920566	−1.33	1.22	⋯	⋯	−4.972	⋯	HIRES
2455799.056115	−5.98	1.48	⋯	⋯	−5.005	⋯	HIRES
2456114.931149	−6.02	1.27	⋯	⋯	−4.944	⋯	HIRES
2456133.896429	−2.53	1.31	⋯	⋯	−4.933	⋯	HIRES
2456147.919325	0.44	1.31	⋯	⋯	−4.967	⋯	HIRES
2456163.912379	6.63	1.31	⋯	⋯	−4.922	⋯	HIRES
2456164.801818	−0.23	1.30	⋯	⋯	−4.924	⋯	HIRES
2456166.047782	1.05	1.37	⋯	⋯	−4.927	⋯	HIRES
2456166.759374	−2.55	1.44	⋯	⋯	−4.913	⋯	HIRES
2456167.990464	−3.55	1.51	⋯	⋯	−4.912	⋯	HIRES
2456451.100822	2.47	1.48	⋯	⋯	−4.953	⋯	HIRES
2456483.086583	5.47	1.94	⋯	⋯	−4.943	⋯	HIRES
2456486.833662	0.39	1.31	⋯	⋯	−4.938	⋯	HIRES
2456488.822611	7.76	1.19	⋯	⋯	−4.926	⋯	HIRES
2456494.987645	5.79	1.29	⋯	⋯	−4.926	⋯	HIRES
2456506.780605	1.27	1.21	⋯	⋯	−4.931	⋯	HIRES
2456507.968496	−0.32	1.23	⋯	⋯	−4.947	⋯	HIRES
2456532.877110	2.73	1.14	⋯	⋯	−4.940	⋯	HIRES
2456830.887055	−0.74	1.39	⋯	⋯	−4.962	⋯	HIRES
2456850.049952	2.59	1.26	⋯	⋯	−4.954	⋯	HIRES
2456828.616553	−37327.64	6.07	6.63443	−0.02220	−4.9522	0.0699	HARPS-N
2456828.651774	−37320.76	1.50	6.66839	−0.03107	−4.9819	0.0106	HARPS-N
2456829.664594	−37319.24	1.56	6.66376	−0.03206	−4.9788	0.0111	HARPS-N
2456830.665375	−37319.62	1.83	6.67197	−0.03353	−4.9882	0.0146	HARPS-N
2456831.690035	−37319.77	1.70	6.66691	−0.03495	−4.9634	0.0128	HARPS-N
2456832.615999	−37314.68	1.80	6.66854	−0.03182	−4.9975	0.0150	HARPS-N
2456833.672301	−37322.65	2.12	6.67240	−0.03957	−4.9553	0.0182	HARPS-N
2456834.581042	−37315.83	2.01	6.66876	−0.03718	−4.9764	0.0169	HARPS-N
2456834.677908	−37321.70	1.98	6.66877	−0.03961	−4.9768	0.0169	HARPS-N
2456835.587887	−37318.54	2.07	6.65839	−0.02603	−5.0090	0.0197	HARPS-N
2456845.576470	−37322.26	1.44	6.66668	−0.03079	−4.9739	0.0097	HARPS-N
2456846.662015	−37327.96	2.10	6.65639	−0.02891	−4.9885	0.0184	HARPS-N
2456847.656794	−37321.88	2.67	6.65969	−0.03820	−4.9938	0.0273	HARPS-N
2456848.652903	−37327.73	1.68	6.66082	−0.03667	−4.9900	0.0131	HARPS-N
2456849.657878	−37326.37	2.17	6.66712	−0.02913	−4.9558	0.0182	HARPS-N
2456850.660745	−37324.98	2.17	6.66481	−0.03397	−4.9575	0.0182	HARPS-N
2456851.654237	−37323.89	1.66	6.66590	−0.03297	−4.9625	0.0121	HARPS-N
2456852.655703	−37316.98	2.55	6.67237	−0.03670	−4.9994	0.0260	HARPS-N
2456853.657053	−37318.56	1.62	6.67355	−0.03209	−4.9818	0.0122	HARPS-N
2456865.684262	−37320.32	1.67	6.66629	−0.02961	−4.9818	0.0130	HARPS-N
2456866.681774	−37323.10	3.48	6.66279	−0.04100	−4.9629	0.0388	HARPS-N
2456883.639193	−37324.45	1.89	6.65864	−0.03904	−5.0004	0.0164	HARPS-N
2456884.647365	−37324.61	1.84	6.66539	−0.02953	−4.9982	0.0166	HARPS-N
2456885.644031	−37322.71	1.86	6.66503	−0.03266	−4.9753	0.0159	HARPS-N
2456886.642561	−37346.88	11.81	6.66764	−0.08202	−5.0609	0.2193	HARPS-N
2456887.651622	−37322.33	1.83	6.67033	−0.03390	−4.9883	0.0152	HARPS-N
2456888.580937	−37321.19	2.74	6.65114	−0.02926	−4.9577	0.0266	HARPS-N
2456889.585275	−37324.05	2.15	6.65998	−0.03761	−4.9716	0.0189	HARPS-N
2456903.541993	−37318.65	1.41	6.66569	−0.03568	−4.9851	0.0095	HARPS-N
2456919.514886	−37322.96	2.44	6.66569	−0.03738	−4.9735	0.0224	HARPS-N
2456922.547287	−37323.20	1.67	6.65621	−0.04098	−5.0131	0.0146	HARPS-N
2456923.501548	−37320.74	1.62	6.66346	−0.03370	−5.0026	0.0127	HARPS-N
2456924.510113	−37318.66	2.43	6.65831	−0.03765	−4.9763	0.0237	HARPS-N
2456936.514073	−37326.01	1.72	6.65525	−0.03770	−4.9864	0.0143	HARPS-N
2456939.418861	−37323.19	1.36	6.65561	−0.03499	−4.9986	0.0093	HARPS-N
2456969.402685	−37323.75	3.10	6.65044	−0.03308	−5.0094	0.0377	HARPS-N
2457106.734166	−37320.49	2.77	6.66677	−0.03337	−4.9801	0.0289	HARPS-N
2457116.717298	−37324.15	1.51	6.64918	−0.03136	−5.0282	0.0120	HARPS-N
2457118.706394	−37327.91	1.87	6.64875	−0.03660	−5.0338	0.0181	HARPS-N
2457121.726137	−37324.71	1.61	6.64485	−0.03692	−5.0174	0.0141	HARPS-N
2457153.685174	−37323.84	2.87	6.64814	−0.04127	−5.0063	0.0321	HARPS-N
2457156.714776	−37324.15	12.94	6.66453	−0.05998	−4.9161	0.1624	HARPS-N
2457159.642662	−37323.22	2.14	6.65282	−0.03562	−4.9989	0.0202	HARPS-N
2457160.638323	−37320.99	2.03	6.65589	−0.04059	−5.0130	0.0185	HARPS-N
2457161.626357	−37324.63	1.70	6.65201	−0.04164	−5.0095	0.0140	HARPS-N
2457180.658376	−37322.50	1.60	6.64953	−0.03798	−5.0160	0.0127	HARPS-N
2457181.686408	−37322.02	1.75	6.65033	−0.03871	−5.0016	0.0146	HARPS-N
2457182.670828	−37322.60	1.58	6.64797	−0.03663	−5.0014	0.0122	HARPS-N
2457183.652886	−37321.34	1.81	6.64944	−0.03550	−4.9890	0.0145	HARPS-N
2457184.643705	−37324.94	2.35	6.65649	−0.03450	−5.0265	0.0237	HARPS-N
2457185.662466	−37324.52	1.51	6.65258	−0.03701	−5.0096	0.0115	HARPS-N
2457186.662672	−37325.81	1.52	6.65510	−0.03815	−5.0237	0.0118	HARPS-N
2457188.679310	−37328.26	1.58	6.65099	−0.03864	−5.0010	0.0123	HARPS-N
2457189.672084	−37322.40	2.63	6.64860	−0.03852	−5.0649	0.0316	HARPS-N
2457190.685669	−37325.59	1.61	6.64984	−0.03343	−5.0211	0.0135	HARPS-N
2457191.685746	−37323.07	1.57	6.65434	−0.03606	−5.0243	0.0128	HARPS-N
2457192.684342	−37324.65	1.67	6.65251	−0.03217	−5.0374	0.0147	HARPS-N
2457193.684869	−37323.90	1.44	6.65760	−0.03912	−5.0097	0.0107	HARPS-N
2457195.594752	−37318.06	1.65	6.65239	−0.03455	−4.9930	0.0132	HARPS-N
2457221.626801	−37324.84	1.37	6.64834	−0.04333	−5.0237	0.0101	HARPS-N
2457222.569536	−37324.99	1.78	6.65155	−0.04076	−4.9982	0.0150	HARPS-N
2457223.579194	−37322.72	3.17	6.63767	−0.04077	−5.0049	0.0375	HARPS-N
2457225.522395	−37319.94	3.63	6.64970	−0.03534	−5.0279	0.0478	HARPS-N
2457226.582966	−37317.08	3.17	6.66790	−0.02615	−5.0750	0.0424	HARPS-N
2457226.606265	−37321.98	2.58	6.65109	−0.03404	−5.0105	0.0277	HARPS-N
2457227.627333	−37316.52	1.67	6.66113	−0.03348	−5.0084	0.0136	HARPS-N
2457228.630703	−37320.68	2.72	6.66677	−0.02712	−4.9272	0.0251	HARPS-N
2457229.584812	−37316.53	2.51	6.65780	−0.02953	−5.0054	0.0266	HARPS-N
2457230.528273	−37320.18	3.34	6.66154	−0.02550	−4.9889	0.0385	HARPS-N
2457254.631919	−37325.16	2.41	6.66627	−0.03812	−5.0093	0.0260	HARPS-N
2457256.398580	−37314.12	2.36	6.64870	−0.04351	−4.9894	0.0239	HARPS-N
2457257.413714	−37315.46	1.88	6.65391	−0.03681	−4.9651	0.0149	HARPS-N
2457267.435503	−37324.84	2.04	6.64591	−0.03502	−5.0136	0.0193	HARPS-N
2457268.492540	−37322.86	2.21	6.64810	−0.04496	−5.0215	0.0224	HARPS-N
2457269.418733	−37320.32	1.69	6.65385	−0.03696	−5.0211	0.0145	HARPS-N
2457270.407599	−37321.14	1.36	6.65014	−0.03809	−5.0301	0.0099	HARPS-N
2457271.408119	−37325.01	1.44	6.65177	−0.03459	−5.0336	0.0110	HARPS-N
2457273.426969	−37326.84	1.46	6.64860	−0.04486	−5.0183	0.0109	HARPS-N
2457301.384627	−37319.51	1.47	6.65573	−0.03632	−5.0121	0.0111	HARPS-N
2457302.383904	−37318.58	1.79	6.65019	−0.03970	−5.0308	0.0159	HARPS-N
2457321.426080	−37318.14	2.40	6.64231	−0.03293	−5.0135	0.0243	HARPS-N
2457330.417736	−37318.77	1.89	6.65257	−0.03583	−5.0053	0.0167	HARPS-N
2457334.397358	−37321.22	1.64	6.65030	−0.04224	−5.0085	0.0131	HARPS-N

Download table as:  ASCIITypeset images: 1 2

We first analyzed our observations with a periodogram, then with a correlation plot, and then constructed a model for our spectroscopic data.

4.2.1. Periodogram Analysis

Before modeling our spectroscopic observations, we first investigated the frequency structure of our data. We made a generalized Lomb–Scargle periodogram (Scargle 1982; Zechmeister & Kürster 2009) of $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$, BIS, FWHM, RV, and the window function of the observation time series, all of which can be seen in Figure 3. $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$, BIS, and FWHM are indicators of stellar activity (Queloz et al. 2001; also see Haywood 2015, and references therein). The window function shows how the signals are modified by the time sampling of the measurements.

Zoom In Zoom Out Reset image size

The HARPS-N RV periodogram shows a clear peak at 82 days, the orbital period of Kepler-538b; the HIRES RV periodogram shows a weaker signal at the same period. None of the other periodograms show a similar feature, lending credence to the RV detection of Kepler-538b. The RV periodograms also exhibit two larger peaks near 0.03–0.04 days⁻¹, interpreted as the rotational frequency. Indeed, as our model fit discussed later in Section 4.3 and the results in Table 2 will show, both peaks fall within the 1σ confidence region of the stellar rotation period. (See Section 4.2.3 for a description of our rotation period estimation.) We also find that the long-term trends observed in the activity indices, combined with the spectral window, affect the periodograms, since a long-term trend is clearly noticeable (see Figure 2 and Table 3). We removed these trends and found the resulting periodograms show a peak at the rotational period, but nothing at the orbital period.

4.2.2. Correlation Analysis

We also examined correlations between the RV observations and the other stellar activity indices. As can be seen in Figure 4, there is a slightly stronger correlation between RV and FWHM than between RV and BIS or $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$. However, there may also be useful information in the correlations between RV and BIS or $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$. In order to test this, we cross-checked results that included BIS and $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$ in the modeling against those that did not and found consistent results. For this reason, and for the sake of simplicity, in this paper we only report our analysis of RVs in conjunction with FWHM observations.

Zoom In Zoom Out Reset image size

4.2.3. General RV Modeling Approach

In order to model our RV and FWHM observations, we followed the method described in Rajpaul et al. (2015, hereafter R15), which establishes a method to characterize stellar activity that uses simultaneous regression of distinct data types (with potentially distinct time series). Here we briefly discuss Gaussian process (GP) regression and the novel approach to GPs used by R15.

In brief, a GP is a stochastic process that captures the covariance between observations and allows for the modeling of correlated noise (Rasmussen & Williams 2006). A GP is specified by a covariance matrix in which the diagonal elements are the individual observation variances and each off-diagonal element describes the covariance between two observations. The values of the off-diagonal elements are determined by a kernel function, which describes the nature of the correlated noise. GPs provide a great deal of flexibility that has made them an effective tool to account for stellar activity (Haywood et al. 2014). R15 recommended characterizing stellar activity with a quasi-periodic (QP) kernel, which balances physical motivation with simplicity. The QP kernel uses four parameters (commonly called hyperparameters) and defines the covariance matrix as follows:

$$\begin{eqnarray}&&\begin{array}{l}{K}_{\mathrm{QP}}({t}_{i},{t}_{j})\\ \ =\ {h}^{2}\exp \left(-\displaystyle \frac{{\sin }^{2}\left(\pi ({t}_{i}-{t}_{j})/{P}_{* }\right)}{2{\lambda }_{p}^{2}}-\displaystyle \frac{{\left({t}_{i}-{t}_{j}\right)}^{2}}{2{\lambda }_{e}^{2}}\right),\end{array}\end{eqnarray} \tag{ 2 }$$

where t_i and t_j are observations made at any two times, h is the amplitude hyperparameter (though not a true amplitude, as it incorporates some multiplicative constants), P_* is the period of the variability (i.e., the rotation period in the case of stellar activity), λ_p is the inverse harmonic complexity (a smoothness factor that acts as a proxy for the number of turning points and inflection points per rotation period), and λ_e is an exponential decay factor (scaling with, though not exactly equal to, the decay timescale of the spots on the star).

One of the key insights of R15 is the way in which they related multiple GPs to one another. GP regression can be used on multiple data sets by constructing a covariance matrix that describes the covariances between two observations of any type. In our case, this means any possible pairing of RV–RV, RV-FWHM, or FWHM–FWHM data points. The following equations (based on Equations (13) and 14 from R15) relate RV and FWHM:

$$\begin{eqnarray}&&{\rm{\Delta }}\mathrm{RV}={V}_{c}G(t)+{V}_{r}\dot{G}(t)\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&\mathrm{FWHM}={F}_{c}G(t).\end{eqnarray} \tag{ 4 }$$

Here, G(t) is an underlying GP directly quantifying stellar activity, and V_c, V_r, and F_c are amplitude parameters corresponding to the RV convective blueshift suppression effect, RV rotation modulation, and FWHM signal amplitude (note that this means there are three amplitude parameters instead of the single h parameter expressed in Equation (2)). Because RVs and FWHMs respond differently to the underlying stellar activity, this approach allows for more rigorous characterization of the stellar activity than methods using only RV observations, which improves the separation of the stellar and planetary signals.

We followed R15 and simultaneously modeled the HIRES RV data as well as the HARPS-N RV and FWHM data. This included a separate offset parameter and noise parameter (added in quadrature to the uncertainties) for both RV data sets and the FWHM data set (for a total of three offset parameters and three white noise parameters). Finally, the RV reflex motion due to the planet was characterized by a simple five-parameter orbital model: reference epoch, orbital period, reflex motion semi-amplitude, eccentricity, and longitude of periastron.

Because we conducted a joint fit to both photometry and spectroscopy, all orbital parameters except for reflex motion semi-amplitude are simultaneously used in our photometric model. In other words, reference epoch, orbital period, eccentricity, and longitude of periastron are used in both the photometric and spectroscopic components of our full model.

For all of the parameters used in the spectroscopic portion of the model, we assumed uniform, Jeffreys, or modified Jeffreys priors. The specific types and bounds of the priors are all listed in Table 2.

4.3. Parameter Estimation

Overall, our full model included a photometric baseline offset parameter, a photometric white noise parameter, two quadratic limb-darkening parameters, the impact parameter, the transit duration, the planet radius relative to the stellar radius, the reference epoch, the orbital period, eccentricity, longitude of periastron, the reflex motion semi-amplitude, three spectroscopic offset parameters and three spectroscopic white noise parameters (for HIRES RV, HARPS-N RV, and HARPS-N FWHM), and six GP hyperparameters (two corresponding to the two RV semi-amplitudes and one corresponding to the FWHM semi-amplitude in Equations (3) and (4), as well as the stellar rotation period, a smoothness factor, and an exponential decay factor). This yielded a total of 24 parameters, all of which are also listed in Table 2. (Note: because we only modeled the detrended and flattened photometry, our estimation of the stellar rotation period was derived solely from our spectroscopic data.)

We estimated model parameters using MultiNest (Feroz et al. 2009, 2013), a Bayesian inference tool for parameter space exploration, especially well suited for multimodal distributions. We used the following MultiNest settings for our parameter estimation: constant efficiency mode, importance nested sampling mode, multimodal mode, sampling efficiency = 0.01, 1000 live points, and evidence tolerance = 0.1.

Our full results from this analysis are presented in Table 2 and discussed in Section 5. Further, the best-fit transit model is plotted against the photometric data in Figure 1 and the phase-folded, stellar-activity-removed RV observations and model are presented in Figure 5. We find Kepler-538b to have a mass of ${M}_{p}={10.6}_{-2.4}^{+2.5}$ M_⊕, a radius of ${R}_{p}={2.215}_{-0.034}^{+0.040}$ R_⊕, a mean density of ρ_p = 0.98 ± 0.23 ρ_⊕, and negligible eccentricity (consistent with zero, <0.11 at 95% confidence). Notably, thanks to the Gaia parallax, our uncertainty on the planetary radius is less than 2%. For context, the average uncertainty, 0.037 R_⊕, is only 236 km, approximately the distance between Portland and Seattle.³¹ Finally, we also note that our estimates of transit parameters are all within 1σ of those reported in the original Kepler-538b validation paper (Morton et al. 2016).

Zoom In Zoom Out Reset image size

4.4. Model Tests

The results of our stellar characterization and light curve, RV, and FWHM modeling can be found in Tables 1 and 2.

After conducting our model fits and running the requisite follow-up tests, we found the mass of Kepler-538b to be ${M}_{p}\,={10.6}_{-2.4}^{+2.5}$ M_⊕. Combining this with the planetary radius of R_p = ${2.215}_{-0.034}^{+0.040}$ R_⊕ resulted in a planetary density of ${\rho }_{p}=0.98\,\pm 0.23{\rho }_{\oplus }$, or 5.4 ± 1.3 g cm⁻³.

Owing to its long orbital period, and its location on the mass–radius diagram, Kepler-538b likely consists of a significant fraction of ices (dominated by water ice), in addition to rocks/metals, and a small amount of gas (Zeng et al. 2018). Its host star is slightly less massive than our own Sun. Because the luminosity of a main-sequence star is a strong function of its mass (typically to the power of 3 or 4), the luminosity of the host star Kepler-538 is somewhat less than the Sun. Therefore, the snowline in the disk when this system was formed was closer in, increasing the likelihood for Kepler-538b to accrete ices during its formation.

The estimated bulk density of Kepler-538b is comparable to that of the Earth. However, this high mean density is partly due to its high mass resulting in more compression of materials under self-gravity. Its uncompressed density, as revealed by the mass–radius curves (Zeng & Sasselov 2013; Zeng et al. 2016) in Figure 6, is consistent with a composition somewhat less dense than pure-rocky and/or Earth-like rocky (1:2 iron/rock mixture). One ready explanation is that Kepler-538b is an icy core, which for some reason had not accreted as much gas as our own Uranus or Neptune (both are estimated to have a few up to ten perfect mass of gas).

Zoom In Zoom Out Reset image size

The eccentricity of Kepler-538b is small (less than 0.11 with 95% confidence). However, the planet may still have arisen from a dynamical origin, that is, inward planet migration due to planet–planet gravitational interactions (Raymond et al. 2009). Some planet formation theories have suggested the formation of multiple icy cores in relatively adjacent space near the snowline around a host star, increasing the likelihood of dynamical interactions among them and resulting in inward scatterings for some of them. If Kepler-538b were scattered inward, then its orbital eccentricity could have been higher initially, and then damped to its current value through interactions with the disk when the disk was still around. Alternatively, inward migration through planet-disk interactions may be a more likely scenario, since a disk would always keep the planet orbital eccentricity low (Chambers 2018; Morbidelli 2018) and would probably be required to damp any eccentricity from scattering.

In summary, Kepler-538b is only the tip of a huge iceberg, likely representing a class of planets common in our Galaxy, but which are not found in our own solar system. The absence of planets in between the size of the Earth and Neptune (about four Earth radii) is linked to the formation/presence of a gas giant—Jupiter (Izidoro et al. 2015; Barbato et al. 2018), and vice versa.

To date, very few exoplanets have been found on long-period orbits that also have any kind of mass measurements. In fact, according to the NASA Exoplanet Archive³² (accessed 2019 July 31), there are only 10 transiting exoplanets (excluding Kepler-538b) with an RV mass measurement and an orbital period greater than 50 days. If we look at other common methods of mass measurement (specifically transit timing variations and dynamical mass measurements of circumbinary planets), that number only increases to 37.

Further, most of those planets are quite large, more similar to Jupiter or Saturn in mass and radius than Neptune or Earth. Figure 7 demonstrates where Kepler-538b fits into this sparse region of parameter space. Kepler-538b is one of the very few small, low-mass planets well characterized to date.

Zoom In Zoom Out Reset image size

As the sample of small, long-period planets with precisely determined masses and densities grows, we will be able to address a number of fundamental questions. For example, what effect does stellar incident flux have on the size and composition of exoplanets? Since most known exoplanets have periods shorter than that of Mercury, it is difficult to analyze exoplanet composition and size for incident fluxes comparable to or less than that of Earth. Similarly, is there a relationship between the location or depth of the planet radius occurrence gap detected by Fulton et al. (2017) and a planet's mass or composition? Further characterization of this gap at longer periods would help confirm (or refute) the photoevaporation explanation of the gap and therefore provide insights about exoplanet formation.

5.1. Detection of Kepler-538b with Other Methods

5.2. Potential for Atmospheric Characterization

One interesting question to ask about Kepler-538b is whether or not it may be amenable to atmospheric characterization via transmission spectroscopy. The James Web Space Telescope (JWST; Gardner et al. 2006; Deming et al. 2009; Kalirai 2018) will devote a significant portion of its mission to the characterization of exoplanet atmospheres. The spectra shown in Figure 8 for the atmosphere of Kepler-538b were generated by the JWST Exoplanet Targeting (JET) code (C. D. Fortenbach & C. D. Dressing 2019, in preparation) assuming five observed transits. This code first takes the observed planet and system parameters (R_p, period, insolation flux, R_*, T_eff, and J-band magnitude) and then derives other key parameters (semimajor axis, T_eq, planet surface gravity, planet mass, and transit duration). In this case we used the planet mass already determined in this paper. We also assumed an optimistic low-metallicity (five times solar) planetary atmosphere with no clouds. JET then used Exo-Transmit (Kempton et al. 2017) to generate model transmission spectra and used Pandexo (Batalha et al. 2017) to generate simulated instrument spectra. We focused on the Near InfraRed Imager and Slitless Spectrograph (NIRISS) SOSS-Or1 and NIRSpec G395M instruments/modes since they are, according to Batalha & Line (2017), best suited for exoplanet transmission spectroscopy. Finally, the JET code performed a statistical analysis for multiple transits and determined if the simulated instrument spectra fit the model well enough to confirm a detection. Given current estimates of the precision (noise floor) of these JWST instruments (as well as visual inspection of the simulated spectra after five transits in Figure 8), it would likely be very difficult to detect the Kepler-538b atmosphere even with a large number of transit observations with JWST.

Zoom In Zoom Out Reset image size

Perhaps other next-generation observatories such as the Thirty Meter Telescope (Sanders 2013), the Extremely Large Telescope (Udry et al. 2014), the Giant Magellan Telescope (Johns et al. 2012), or the Large UV/Optical/IR Surveyor (The LUVOIR Team 2018) will be able to make such a project feasible.

5.3. Possibility of a Second Planet in the System

In this paper, we analyze the Kepler-538 system in order to determine the properties of Kepler-538b, the single, known exoplanet in the system. Kepler-538 is a 0.924 M_⊙, G-type star with a visual magnitude of V = 11.27. We model the Kepler light curve and determine the orbital period of Kepler-538b to be P = 81.74 days and the planetary radius to be ${R}_{p}={2.215}_{-0.034}^{+0.040}$ R_⊕ (for reference, 0.037 = 236 km, approximately the distance between Portland and Seattle³⁶ ). These results are in agreement with previous transit fits. We also determine the planetary mass by accounting for stellar activity via a GP regression that uses information from the FWHM and RV observations simultaneously. Our model fit yields a mass estimate for Kepler-538b of ${M}_{p}={10.6}_{-2.4}^{+2.5}$ M_⊕. Combined, these results show the planet to have a density of ${\rho }_{p}=0.98\pm 0.23{\rho }_{\oplus }=5.4\pm 1.3$ g cm⁻³. This suggests a composition and atmosphere somewhere between that of Earth and Neptune, with a significant fraction of ices (dominated by water ice), in addition to rocks/metals, and a small amount of gas (Zeng et al. 2018).

To date, there have been very few precise and accurate mass measurements of long-period exoplanets. Beyond 50 days, Kepler-538b is only the 11th transiting exoplanet with an RV mass measurement (NASA Exoplanet Archive;³⁷ accessed 2019 May 4). Additional, well-constrained mass measurements of long-period planets will improve our understanding of the long-period exoplanet population. Beyond that, they will also help to answer questions about the short-period planet population, such as the nature of the planetary radius occurrence gap (Fulton et al. 2017) and the effect of stellar flux on exoplanet compositions and atmospheres.

With new, next-generation spectrographs such as the Habitable zone Planet Finder (HPF; Mahadevan et al. 2010, 2014), Keck Planet Finder (KPF; Gibson et al. 2016, 2018), EXtreme PREcision Spectrometer (EXPRES; Jurgenson et al. 2016), Echelle SPectrograph for Rocky Exoplanets and Stable Spectroscopic Observations (ESPRESSO; Mégevand et al. 2010), and NASA-NSF Exoplanet Observational Research (NN-EXPLORE) Exoplanet Investigations with Doppler Spectroscopy (NEID; Schwab et al. 2016) coming online now or in the near future, our ability to characterize long-period exoplanets will only improve. Better data will require more advanced analysis methods to extract as much information as possible. The methods used in this paper, such as GP regression, injection tests, and simultaneous modeling of RV observations and stellar activity indices, are valuable tools that strengthen the analysis of spectroscopic data, improve exoplanet characterization, and therefore better our understanding of the exoplanet population as a whole.

A.W.M. is supported by the NSF Graduate Research Fellowship grant No. DGE 1752814.

V.M.R. thanks the Royal Astronomical Society and Emmanuel College, Cambridge, for financial support.

C.D.D. acknowledges support from the NASA K2 Guest Observer program through grant 80NSSC19K0099.

This work was performed in part under contract with the California Institute of Technology (Caltech)/Jet Propulsion Laboratory (JPL) funded by NASA through the Sagan Fellowship Program executed by the NASA Exoplanet Science Institute (R.D.H.)

We acknowledge the support by INAF/Frontiera through the "Progetti Premiali" funding scheme of the Italian Ministry of Education, University, and Research.

This publication was made possible through the support of a grant from the John Templeton Foundation. The opinions expressed in this publication are those of the authors and do not necessarily reflect the views of the John Templeton Foundation.

This work was supported in part by a grant from the Carlsberg Foundation.

This paper includes data collected by the Kepler mission. Funding for the Kepler mission is provided by the NASA Science Mission directorate.

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–Hubble Space Telescope data is provided by the NASA Office of Space Science via grant NNX13AC07G and by other grants and contracts.

Some of the data presented herein were obtained at the W. M. Keck Observatory (which is operated as a scientific partnership among Caltech, UC, and NASA). The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Maunakea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.

We would like to thank the HIRES observers who carried out the RV observations with Keck.

We are grateful to Howard Isaacson for valuable proofreading and suggestions as well as providing access to HIRES $\mathrm{log}\,{R}_{{HK}}^{{\prime} }$ observations.

Based on observations made with the Italian Telescopio Nazionale Galileo (TNG) operated by the Fundación Galileo Galilei (FGG) of the Istituto Nazionale di Astrofisica (INAF) at the Observatorio del Roque de los Muchachos (La Palma, Canary Islands, Spain).

The HARPS-N project has been funded by the Prodex Program of the Swiss Space Office (SSO), the Harvard University Origins of Life Initiative (HUOLI), the Scottish Universities Physics Alliance (SUPA), the University of Geneva, the Smithsonian Astrophysical Observatory (SAO), and the Italian National Astrophysical Institute (INAF), the University of St Andrews, Queens University Belfast, and the University of Edinburgh.

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.

Facilities: Kepler - The Kepler Mission, FLWO:1.5 m (TRES) - , TNG: (HARPS-N) - , Keck:I (HIRES) - KECK I Telescope, Gaia - , NASA Exoplanet Archive - , ADS - , MAST - .

Software: MultiNest (Feroz et al. 2009, 2013), Exo-Transmit (Kempton et al. 2017), Pandexo (Batalha et al. 2017), JET (C. D. Fortenbach & C. D. Dressing 2019, in preparation), isochrones (Morton 2015), lightkurve (Barentsen et al. 2019), BATMAN (Kreidberg 2015).