ZODIACAL EXOPLANETS IN TIME (ZEIT). III. A SHORT-PERIOD PLANET ORBITING A PRE-MAIN-SEQUENCE STAR IN THE UPPER SCORPIUS OB ASSOCIATION

Many known exoplanets orbit within 0.1 AU of their host star, where they are readily detected by the transit and Doppler methods (e.g., Howard et al. 2010; Fressin et al. 2013). Whether these planets formed near their present position (in situ), i.e., from circumstellar material interior to 0.5 AU (e.g., Chiang & Laughlin 2013; Ogihara et al. 2015), or accreted at distances $\gt 1\,{\rm{A}}{\rm{U}}$ and later migrated inwards (e.g., Schlaufman et al. 2009; Raymond & Cossou 2014) is actively debated. If the planets migrated, the physical mechanism(s) behind their migration is a further point of debate.

Mechanisms of planet migration can be loosely divided into three categories: interactions with the protoplanetary disk (e.g., Ida & Lin 2008; Lubow & Ida 2010), interactions between a stellar companion and the planet (the Lidov–Kozai mechanism, e.g., Wu & Murray 2003), and interactions among multiple planets (e.g., Ford & Rasio 2006; Chatterjee et al. 2008). Disk migration must occur before the protoplanetary disk dissipates/photoevaporates (≲10 Myr; e.g., Ward 1997). Migration involving angular momentum exchange with a third body typically operates on timescales much longer than disk migration (≃100 Myr to more than 1 Gyr), depending on the orbital and physical properties of the planet and perturber (Fabrycky & Tremaine 2007; Nagasawa et al. 2008).

The difference in timescales presents a possible method to distinguish between these migration mechanisms. Close-in super-Earths or Jupiter-size planets around stars younger than 100 Myr could not have migrated by a slow process such as planet–planet or planet–star interaction, and instead likely formed in situ or migrated quickly through interaction with the disk. A comparison of the occurrence rate of close-in planets over a range of ages (10–1000 Myr) would constrain the fraction of planets migrating on a given timescale.

High-precision photometry offers the best opportunity to detect the close-in planets needed to test migration theories (e.g., Janes 1996). Such planets are too close to their host star to be detected by direct imaging (e.g., Le Bouquin & Absil 2012; Bowler et al. 2015b). Starspot-induced jitter complicates the detection of the planetary reflex motion (Paulson et al. 2004; Mahmud et al. 2011), such that radial velocity (RV) surveys of young stars primarily uncover hot Jupiters around stars older than 100 Myr (e.g., Quinn et al. 2012a, 2014). Spot modulation can generate complicated variations in the light curve that make the detection of transiting planets more difficult. However, spots and transits create characteristically different patterns in a light curve, which can be separated with precise photometry. Indeed the only close-in planet (candidates) around stars $\lt 20\,{\rm{Myr}}$ old are from transit surveys (e.g., Mamajek et al. 2012; van Eyken et al. 2012; Kenworthy et al. 2015).

The repurposed Kepler mission (K2, Howell et al. 2014) has the photometric precision (tens of ppm) and observational baseline (70–80 days) required to detect small planets and rule out false-positive signals related to stellar youth (e.g., debris disks and spots). We are carrying out a search for planets transiting stars in 10–800 Myr young open clusters and OB associations using K2. Our survey, Zodiacal Exoplanets in Time, includes Upper Scorpius (≃11 Myr, Pecaut et al. 2012; Rizzuto et al. 2016), Taurus (0–5 Myr, Kenyon et al. 2008), Pleiades (≃125 Myr, Dahm 2015), and Praesepe and Hyades (650–800 Myr, Brandt & Huang 2015). Our goal is to better understand whether planets evolve from infancy (1–10 Myr) to maturity (≳1 Gyr), including how planets migrate and on what timescales.

Here we confirm and characterize a ${5.04}_{-0.37}^{+0.34}$ ${R}_{\oplus }$ planet (K2-33b) orbiting the pre-main-sequence (PMS) star K2-33 (2MASS J16101473-1919095), a member of the Upper Scorpius subgroup of the Scorpius–Centaurus (Sco–Cen) OB association. K2-33b was previously identified as a planet candidate by Vanderburg et al. (2016), but assigned inaccurate stellar and planetary parameters owing to the assumption of a main-sequence age and an unreddened spectral energy distribution for the host star. In Section 2 we describe our follow-up observations, including moderate- to high-resolution spectroscopy, adaptive optics (AO) imaging and non-redundant aperture masking, and transit photometry. Our analysis of the light-curve data is described in Section 3. In Section 4 we derive parameters for the host star. We use the available data to confirm a planetary origin of the transit signal, as we describe in Section 5. We conclude in Section 6 with a brief summary and discussion of the importance of this planet.

2.1. K2 Observations and Light-curve Extraction

From 2014 August 23 to 2014 November 13 (Campaign 2) K2 observed the core of Upper Scorpius, including K2-33. Owing to the loss of two reaction wheels, the Kepler satellite drifts. To correct the pointing, Kepler's thrusters fire every ∼6 hr. However, during the drift and subsequent thruster fire a stellar image usually moves over the detector. Combined with variations in the pixel sensitivity this generates changes in total measured flux from a given star as a function of centroid position.

Vanderburg & Johnson (2014) present a method for mitigating or removing this noise, but it is not optimized for highly variable stars, e.g., the young stars of Upper Scorpius. The transit can still be identified in the light curve of Vanderburg & Johnson (2014), but systematic trends are present, including a discontinuity in the middle of the observations when Kepler changed the direction of its roll, and large changes in the point-to-point scatter over the observing window due to poor treatment of the thruster-fire systematics. Following Becker et al. (2015) and Mann et al. (2016), we extracted a new light curve by simultaneously fitting for low-frequency variations from stellar activity, Kepler flat field, and the transits of K2-33b using a least-squares minimization. Stellar variability and the effect of errors in detector response were both modeled as splines as a function of time and centroid position with breakpoints every 0.2 days and 04, respectively. Unlike Vanderburg et al. (2016), we did not apply separate systematic corrections to the first and second halves of the K2 campaign, which removed the major discontinuity. The resulting light curve is relatively clear of visible systematic errors (see Figure 1).

Zoom In Zoom Out Reset image size

2.2. Optical Spectrum from SuperNova Integral Field Spectrograph (SNIFS)

2.3. NIR Spectrum from ARCoIRIS

2.4. High-resolution NIR Spectrum

2.5. AO Imaging and Aperture Masking

On 2016 March 18 (UT), we obtained AO imaging using natural guide stars (Wizinowich et al. 2000) and non-redundant aperture masking (NRM) of K2-33 with the facility imager, NIRC2, on Keck II atop Mauna Kea. All observations were taken in vertical angle mode, using the smallest pixel scale (9.952 ± 0.002 mas pixel⁻¹). Imaging was taken with the K' filter and masking with the nine-hole mask. After AO loops closed we took eight images, each with 20 coadds and an integration time of 0.5 s per coadd. For NRM we took 10 interferograms, each with an integration time of 20 s and 1 coadd.

Each frame was linearized and corrected for geometric distortion using the NIRC2 solution from Yelda et al. (2010). Images were dark-subtracted and flat-fielded. We interpolated over "dead" and "hot" pixels. Dead pixels were identified from superflats taken in 2006–2013 as any pixel with a response of <30% in at least half of all superflats. Similarly, hot pixels were identified from a comparable set of superdarks as any pixel with $\geqslant 10$ counts in at least half of the superdarks. Pixels with flux levels $\gt 10\sigma$ above the median of the eight adjacent pixels were flagged as cosmic rays or transient hot pixels and replaced with the median.

To detect faint and wide ($\gtrsim 500$ mas) companions in the AO images we subtracted an azimuthal median model of the point-spread function (PSF) built from the smoothed PSF of K2-33. This added no additional noise at wide separations, but left the speckles in place, making it non-ideal for detecting close-in companions. To probe smaller inner working angles we instead constructed and subtracted the best-fitting PSF of another (single-star) target taken on the same night. We stacked all subtracted frames of K2-33, and identified companions by measuring the flux in apertures of 40 mas (radius) centered on every image pixel. We measured our detection limits from the standard deviation of the flux among all apertures in a 5 pixel annulus around the primary. We found no apertures that contained sufficient flux within the NIRC2 field of view to be considered as an astrophysical source.

The aperture masking observations use the complex triple product, or closure phase, to remove non-common path errors introduced by atmospheric conditions and variable optical aberrations. To remove systematics, the observation of K2-33 was paired with a calibration observation of USco J160933.8-190456, another member of Upper Scorpius (Preibisch et al. 2001). Binary system profiles can then be fit to the closure phases to produce separations and position angles and calculate contrast limits. The Appendix of Kraus et al. (2008) contains a full explanation of the data reduction and binary profile-fitting for aperture masking data. No sources were detected in the masking data.

The combination of the aperture masking and imaging observations excludes contributions of additional stars at separations 002–3'' to the light curve. Figure 2 displays the imaging and masking contrast limits as a function of separation. Owing to the edges of the detector the azimuthal coverage is not complete outside ≃3''.

Zoom In Zoom Out Reset image size

2.6. Transit Photometry from MEarth

3.1. Transit Identification

3.2. Transit Fitting

We simultaneously fit the K2 and MEarth transit data with a Monte Carlo Markov chain (MCMC) as described in Mann et al. (2016), which we briefly summarize here. We used the emcee Python module (Foreman-Mackey et al. 2013) to fit the model light curves produced by the batman package (Kreidberg 2015) using the algorithm of Mandel & Agol (2002). Following Kipping (2010) we oversampled and binned the model to match the 30 minutes K2 cadence. We used an unbinned model to fit the MEarth data due to the much lower integration time (60 s). We sampled the planet-to-star radius ratio (${R}_{{\rm{P}}}/{R}_{* }$), impact parameter (b), orbital period (P), epoch of the first transit mid-point (T₀), bulk stellar density (${\rho }_{* }$), and two limb darkening parameters (q1 and q2) for each of the two instruments (MEarth and K2). At this young age it is likely that any eccentricity was dampened by the primordial disk (e.g., Tanaka & Ward 2004; Cresswell et al. 2007), so we fix the eccentricity at zero. However, the eccentricity distribution of young planets has not been observationally constrained; we discuss lifting this assumption in Section 6.

We assumed a quadratic limb darkening law and used the triangular sampling method of Kipping (2013) in order to uniformly sample the physically allowed region of parameter space. We applied a prior on limb darkening derived from the atmospheric models of Husser et al. (2013), calculated using the LDTK toolkit (Parviainen & Aigrain 2015), which enabled us to account for errors in stellar parameters and finite grid spacing. Stellar parameters input into LDTK are derived in Section 4. Errors on the limb darkening coefficients are broadened by a factor of two to account for model uncertainties (estimated by comparing limb darkening parameters from different model grids). The filter profiles and CCD transmission functions were taken from Dittmann et al. (2016) for the MEarth bandpass and from the Kepler science center¹⁴ for the Kepler bandpass. This yielded quadratic limb darkening coefficients of ${\mu }_{1}=0.4\pm 0.1$ and ${\mu }_{2}=0.4\pm 0.1$ for Kepler and ${\mu }_{1}=0.26\pm 0.09$ and ${\mu }_{2}=0.4\pm 0.1$ for MEarth.

Our MCMC was allowed to explore $| b| \lt 1+{R}_{{\rm{P}}}/{R}_{* }$, P from 0 to 70 days, ${\rho }_{* }$ from 0 to $\infty$, ${R}_{{\rm{P}}}/{R}_{* }$ from 0 to 1, and T₀ from ±3 days from the initial value, all with uniform priors. All parameters were initialized to the values from our BLS search (Section 3.1), which are based on a Levenberg–Marquardt fit to the light curve (Markwardt 2009). MCMC chains were run using 150 walkers, each with 100,000 steps after a burn-in phase of 5000 steps.

We report the transit fit parameters in Table 1. For each parameter we report the median value with the errors as the 84.1 and 15.9 percentile values (corresponding to 1σ for Gaussian distributions). The model light curve with the best-fit parameters is shown in Figure 3 alongside the K2 and MEarth data. Some correlated errors are present in the MEarth light curve, primarily during ingress, which we attribute to imperfect correction of stellar variability and/or the planet crossing a spot. We also show posteriors and correlations for a subset of parameters in Figure 4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1.  Transit Fit Parameters

Parameter	Value
Period (days)	${5.424865}_{-0.000031}^{+0.000035}$
${R}_{{\rm{P}}}/{R}_{* }$	${0.0432}_{-0.0007}^{+0.0009}$
T0a (BJD – 2400000)	${56898.69288}_{-0.00120}^{+0.00118}$
Density (${\rho }_{\odot }$)	${0.51}_{-0.07}^{+0.04}$
Impact Parameter	${0.16}_{-0.11}^{+0.19}$
Duration (hr)	${4.08}_{-0.07}^{+0.07}$
$a/{R}_{* }$	${10.40}_{-0.50}^{+0.27}$
Inclination (deg)	${89.1}_{-1.1}^{+0.6}$
Eccentricity	0 (fixed)
ω (deg)	0 (fixed)
RPb (R${}_{\oplus }$)	${5.04}_{-0.37}^{+0.34}$

Notes.

^aBJD is the barycentric Julian date given in barycentric dynamical time (TBD) format. ^bPlanet radius is derived using our stellar radius (Section 4).

Download table as:  ASCIITypeset image

The transit posterior favors a low ($\lt 0.4$) impact parameter, although there is a tail in the distribution corresponding to higher impact parameter, lower ${\rho }_{* }$ ($\lt 0.3$), and larger planet radius. This region of parameter space is not ruled out by our independent stellar parameters (see Section 4), so we did apply an additional constraint in the MCMC to remove these solutions from our transit-fit posterior.

Membership in Upper Scorpius. The spatial position and kinematics of K2-33 are consistent with co-motion of the star with the ensemble space velocity of Upper Scorpius. We calculate a photometric distance to K2-33 of 140 ± 16 pc using optical and NIR photometry from the literature, and the 10 Myr isochrone from Chen et al. (2014). This is consistent with the Hipparcos distances to high-mass members of Upper Scorpius (145 ± 15 pc, de Zeeuw et al. 1999; Rizzuto et al. 2011). Using this photometric distance, proper motions from UCAC4 (−9.8 ± 1.7, −24.2 ± 1.8 mas yr⁻¹; Zacharias et al. 2012), and the mean RV from our IGRINS observations, we calculate the Galactic space velocity of K2-33 (U, V ,W) = (5.4 ± 0.5, −15.8 ± 2.2, −8.2 ± 1.2 km s⁻¹). This is consistent with the kinematic models of Chen et al. (2014) and velocity dispersion of ∼2–3 km s⁻¹ (Kraus & Hillenbrand 2008) for Upper Scorpius. Using the Bayesian method from Rizzuto et al. (2011, 2015) we calculate a probability of membership in Upper Scorpius of 96% for K2-33.

K2-33 also shows multiple indicators of youth. We measured the Na-8189 index, a well-calibrated gravity index (Slesnick et al. 2006), to be 0.946 ± 0.005, consistent with other late-type stars in Upper Scorpius and suggesting an age of 5–30 Myr. K2-33 shows a 24 μm excess (Luhman & Mamajek 2012) in Spitzer–MIPS observations, suggesting the presence of a disk, and hence an age for K2-33 of <40 Myr. Further, K2-33 was already identified as a member of the Upper Scorpius subgroup by the presence of a strong Li 6708 Å line (0.45 ± 0.15 Å, Preibisch et al. 2001), an unambiguous indicator of youth for late-type stars.

Spectral type and reddening. Following Kraus et al. (2015) and Ansdell et al. (2016), we simultaneously solved for spectral type and reddening (A_V) to account for correlations between these parameters. We compared our optical spectrum of K2-33 to a grid of 270 unreddened optical spectra of young stars from Herczeg & Hillenbrand (2014). For each template we computed the A_V value that gives the best agreement between the spectrum of K2-33 and that of the template using the reddening law from Cardelli et al. (1989) and masking out the Hα line and the strong O₂ tellurics. The resulting distribution of reduced ${\chi }^{2}$ (${\chi }_{\nu }^{2}$) values yielded a spectral type of M3.3 ± 0.2 with an A_V of ${0.75}_{-0.18}^{+0.21}$ (Figure 5). This error in spectral type does not account for systematic errors in the spectral typing scheme (which can vary by 0.5–1 spectral subtype between methods), so we instead we adopted a more conservative M3.3 ± 0.5. This does not affect the determination of A_V, but A_V could be affected if there are systematic errors in the spectrophotometric calibration of our optical spectrum or the templates of Herczeg & Hillenbrand (2014).

Zoom In Zoom Out Reset image size

To test the sensitivity to our choice of templates, we repeat the above process with M-dwarf spectral templates from Gaidos et al. (2014). These were taken with the same instrument as our spectrum of K2-33, but the spectra are predominantly from old ($\gt 1$ Gyr) stars. The templates of Gaidos et al. (2014) give a slightly earlier spectral type and higher reddening (M3.1, 0.83), but both are consistent with the values derived using templates from young stars. Comparison to the templates of Gaidos et al. (2014) also give significantly higher ${\chi }_{\nu }^{2}$ values than those from Herczeg & Hillenbrand (2014), likely because of gravity-dependent differences in the spectrum.

Effective temperature. We compared the unreddened spectrum to a grid of BT-SETTL CIFIST models¹⁵ (Allard et al. 2011), masking our regions where models poorly reproduce observed spectra and accounting for small errors in the flux and wavelength solution as detailed in Mann et al. (2013) and Gaidos et al. (2014). This method has been shown to reproduce ${T}_{\mathrm{eff}}$ values for main-sequence M dwarfs derived from interferometry (Boyajian et al. 2012), but is poorly tested on PMS stars. However, we accurately reproduced the geometric ${T}_{\mathrm{eff}}$ derived for the low-mass eclipsing binary (EB) USco CTIO5 (Kraus et al. 2015), suggesting that our method yields reasonable ${T}_{\mathrm{eff}}$ values even at young ages. To account for errors in reddening we repeated this process over the range of reddening values derived above. The effect of reddening is small, because the model comparison is driven primarily by the depth of the molecular bands instead of the overall spectral shape. We found a best-fit ${T}_{\mathrm{eff}}$ of 3540 ± 70 K.

Bolometric flux. We compiled optical BVgri photometry from the ninth data release of the AAVSO All-Sky Photometric Survey (APASS, Henden et al. 2012), NIR JHK_S photometry from The Two Micron All Sky Survey (2MASS, Skrutskie et al. 2006), griz photometry from the Sloan Digital Sky Survey (SDSS, Ahn et al. 2012), and $W1W2W3$ infrared photometry from the Wide-field Infrared Survey Explorer (WISE, Wright et al. 2010). We then scaled the (reddened) NIR and optical spectrum to the archival photometry following the procedure from Mann et al. (2015). The flux-calibrated spectrum is plotted in Figure 6. We then unreddened the calibrated spectrum. To calculate ${F}_{\mathrm{bol}}$ we integrated the spectrum over all wavelengths. As with ${T}_{\mathrm{eff}}$, this process was repeated over the range of A_V values, which effectively tripled our error on ${F}_{\mathrm{bol}}$. We found a best fit for ${F}_{\mathrm{bol}}$ of $2.25(\pm 0.26)\times {10}^{-10}$ erg cm⁻² s⁻¹.

Zoom In Zoom Out Reset image size

Stellar radius, mass, and age. We combined the distance to Upper Scorpius (145 ± 15 pc, de Zeeuw et al. 1999) with our measured ${F}_{\mathrm{bol}}$ and ${T}_{\mathrm{eff}}$ to calculate the stellar radius from the Stefan–Boltzmann relation. This gives a radius of $1.02\pm 0.13\,{R}_{\odot }$. We combined this with our transit-fit density from Section 3.2 to get a mass estimate of ${0.55}_{-0.14}^{+0.13}\,{M}_{\odot }$. However, we can significantly improve on these parameters, and also estimate the age of the system using a grid of models of PMS stellar evolution.

We compared our observables to two grids of PMS models computed with an updated version of the Dartmouth stellar evolution code (Dotter et al. 2008; Feiden & Chaboyer 2012). A number of improvements to the code that allow more accurate computation of PMS stars are summarized in Feiden (2016). Germane to K2-33, one of the grids includes effects of magnetic inhibition of convection (Feiden & Chaboyer 2012, 2013), which relieves the observed age discrepancy between early- and late-type stars in Upper Scorpius (Pecaut et al. 2012; Rizzuto et al. 2015), yielding a consistent median age of 9–10 Myr for spectral types A through M (Feiden 2016).

To infer the mass, radius, and age of K2-33, an MCMC method implemented with emcee (Foreman-Mackey et al. 2013) was used to sample the parameter space covered by our two model grids. The nonmagnetic, standard model grid covered a mass range of 0.1–0.9 ${M}_{\odot }$ with a resolution of 0.02 ${M}_{\odot }$ and a metallicity range of −0.5 to $+0.5$ dex at a resolution of 0.1 dex. The magnetic model grid covered a larger mass range of 0.1–1.7 ${M}_{\odot }$ with a resolution of 0.02 ${M}_{\odot }$, but the metallicity was restricted to [M/H] = 0.0 dex. We explored mass, metallicity (for nonmagnetic models), age, and distance to find an optimal fit to the observables, ${T}_{\mathrm{eff}}$ and ${F}_{\mathrm{bol}}$, using the likelihood function given in Mann et al. (2015).

We applied a Gassian prior on distance (145 ± 15 pc; de Zeeuw et al. 1999), a prior on ${\rho }_{* }$ drawn from our transit-fit posterior (Section 3.2), and uniform priors on mass and age. For the nonmagnetic models we applied a Gaussian prior on metallicity (0.0 ± 0.1 dex; Bubar et al. 2011; Mamajek et al. 2013). To test the robustness of the results we also ran chains with a uniform prior on ${\rho }_{* }$ or distance (restricted to 10–1000 pc) for each of the model grids.

The MCMC simulation was set up with 300 walkers at random initial starting positions and was allowed to run over 1000 iterations following a burn-in phase of 250 iterations. Convergence was diagnosed through a combination of visually monitoring trace plots for all 300 chains, monitoring the median acceptance fraction among all chains (between 25% and 50%), and monitoring the autocorrelation time for individual chains. A summary of all MCMC results is given in Table 2. Quoted values represent the median of the posterior and quoted uncertainties are the 68% Bayesian credible intervals.

Table 2.  Stellar Fit Parameters

Models	Priors		Parameters
	Distance (pc)	${\rho }_{* }$ (${\rho }_{\odot }$)	R* (${R}_{\odot }$)	M* (${M}_{\odot }$)	L* (${L}_{\odot }$)	Age (Myr)	[Fe/H]	Distance (pc)
No modela	...	...	${1.02}_{-0.13}^{+0.13}$	${0.55}_{-0.15}^{+0.13}$	${0.15}_{-0.03}^{+0.03}$	...	...	...
	145 ± 15	uniform	${1.01}_{-0.17}^{+0.18}$	${0.56}_{-0.09}^{+0.09}$	${0.14}_{-0.04}^{+0.05}$	${10.80}_{-4.71}^{+8.93}$	0 (fixed)	${143.2}_{-21.6}^{+21.0}$
Magnetic	uniform	${0.51}_{-0.07}^{+0.04}$ b	${1.06}_{-0.08}^{+0.09}$	${0.58}_{-0.10}^{+0.10}$	${0.16}_{-0.04}^{+0.05}$	${9.31}_{-1.48}^{+1.11}$	0 (fixed)	${152.8}_{-23.5}^{+25.2}$
	145 ± 15	${0.51}_{-0.07}^{+0.04}$ b	${1.05}_{-0.07}^{+0.07}$	${0.56}_{-0.09}^{+0.09}$	${0.15}_{-0.03}^{+0.03}$	${9.34}_{-1.26}^{+1.07}$	0 (fixed)	${148.0}_{-15.9}^{+15.6}$
	145 ± 15	uniform	${0.97}_{-0.16}^{+0.18}$	${0.42}_{-0.08}^{+0.10}$	${0.14}_{-0.04}^{+0.05}$	${6.06}_{-2.83}^{+6.18}$	${0.01}_{-0.14}^{+0.13}$	${140.1}_{-20.7}^{+20.1}$
Nonmagnetic	uniform	${0.51}_{-0.07}^{+0.04}$ b	${0.94}_{-0.12}^{+0.08}$	${0.42}_{-0.09}^{+0.10}$	${0.12}_{-0.04}^{+0.04}$	${6.53}_{-1.00}^{+1.56}$	${0.00}_{-0.14}^{+0.13}$	${135.7}_{-24.3}^{+25.2}$
	145 ± 15	${0.51}_{-0.07}^{+0.04}$ b	${0.95}_{-0.09}^{+0.07}$	${0.43}_{-0.07}^{+0.09}$	${0.13}_{-0.03}^{+0.03}$	${6.57}_{-0.96}^{+1.44}$	${0.01}_{-0.14}^{+0.13}$	${139.5}_{-19.2}^{+18.0}$

Notes.

^aRadius from Stefan–Boltzmann relation, mass from radius and transit-fit density. ^bDensity prior taken from transit-fit posterior (Section 3.2).

Download table as:  ASCIITypeset image

All chains produced consistent stellar radii, luminosities, and masses. This is in large part a consequence of the observables and constraints; ${T}_{\mathrm{eff}}$, ${F}_{\mathrm{bol}}$, and distance uniquely determine the stellar radius, luminosity, and, when combined with the transit-fit density, the mass. Thus when the stronger distance and density priors were used the resulting radius/luminosity/mass is relatively independent of the model grid used. It is encouraging that even the use of a uniform distance prior yields a distance consistent with the previously established value for Upper Scorpius. Further, all fits to the magnetic models give an age consistent with the value estimated for high-mass members of Upper Scorpius (Pecaut et al. 2012).

Using a uniform prior on ${\rho }_{* }$ results in a negligible change in mass and radius, although with much larger errors. This suggests that the accuracy of our stellar parameters (but not the precision) is insensitive to the assumption of zero orbital eccentricity for K2-33b. We show a comparison of the model-based and transit-fit densities in Figure 7. Because the constraints on ${\rho }_{* }$ are relatively weak from the model comparison alone, we cannot make definitive statements about the orbital eccentricity of the system from the data alone, and instead rely on the physical argument that a planet migrating via interaction with the disk should have $\simeq 0$ orbital eccentricity.

Zoom In Zoom Out Reset image size

The nonmagnetic models favor a lower mass and radius, and younger age. Magnetic models more accurately reproduce the known distance and age of Upper Scorpius. We therefore adopted parameters from the chain utilizing magnetic models with the distance prior from de Zeeuw et al. (1999) and the density prior from our transit fit, which we use for the rest of the analysis.

Rotation period. We computed the Lomb–Scargle periodogram of the K2 light curve prior to removing the stellar variability. A strong signal is apparent at 6.29 days, which we attribute to spot coverage and the rotation period. We estimate an error of 0.17 days on the rotation period from the width of the peak in the power spectrum of the periodogram. The same period was found from the MEarth long-term monitoring (6.27 days).

Rotation velocity. We determined $v\sin {i}_{* }$ using our high-resolution IGRINS spectrum. We first determined the instrumental profile/resolution by fitting the telluric spectrum derived from the A0V standard (see Section 2) with a series of Gaussian profiles and assuming that telluric lines have negligible intrinsic width compared to the instrument resolution. Orders with fewer than four strong ($\gt 3 \%$ depth) telluric lines were ignored. We assumed that the resolution varies linearly within an order and smoothly in between orders. Our derived instrumental broadening was 0.3–0.5 Å (FWHM), consistent with the previously measured resolution of the spectrograph.

We then compared our IGRINS spectrum of K2-33 to the best-fit BT-SETTL model derived from the moderate-resolution spectra above. We broadened the model first using the instrumental profile derived from our telluric fit, then as a function of $v\sin {i}_{* }$ using the IDL code lsf_rotate (Gray 1992; Hubeny & Lanz 2011). We included seven nuisance parameters to handle normalization of the observed spectrum, small errors in wavelength calibration, and missing or overly deep/shallow lines in the atmospheric model. Orders with S/N $\lt \,20$ were removed. Each order was fit separately, each time adjusting $v\sin {i}_{* }$ and the nuisance parameters to minimize the difference between the model and IGRINS spectrum. We adopted the mean and standard error of the $v\sin {i}_{* }$ determinations across all orders as the final value and error: 8.2 ± 1.8 km s⁻¹.

Sky-projected stellar inclination. The combination of $v\sin {i}_{* }$, rotation period, and stellar radius enabled a calculation of the (sky-projected) rotational inclination (i_*) of K2-33. We first calculated the equatorial velocity ${V}_{{\rm{eq}}}=\tfrac{2\pi {R}_{* }}{{P}_{{\rm{rot}}}}$, where ${P}_{{\rm{rot}}}$ is the stellar rotation period measured above, which yielded a velocity of 8.6 ± 0.7 km s⁻¹. We assumed that effects of differential rotation are encapsulated in our ${P}_{{\rm{rot}}}$ measurement error, although this depends on where on the star the spots are located. We converted $v\sin {i}_{* }$ and ${V}_{{\rm{eq}}}$ to a posterior in $\cos ({i}_{* })$, which handles regions of the posterior where $v\sin {i}_{* }$ $\gt {V}_{{\rm{eq}}}$(see Morton & Winn 2014 for more details). The resulting posterior gives a lower limit on stellar inclination of ${i}_{* }\gt 63$° at 68.3% (1σ), suggesting that the planetary orbit is not highly misaligned with the stellar rotation.

A summary of all derived stellar parameters and errors is given in Table 3.

Table 3.  Parameters of EPIC 205117205

Parameter	Value	Source
Identifiers
EPIC 205117205		EPIC
USco J161014.7-191909		Preibisch et al. (2002)
2MASS J16101473-1919095		2MASS
Astrometry
α R.A. (hh:mm:ss J2000)	16:10:14.73	EPIC
δ decl. (dd:mm:ss J2000)	−19:19:09.38	EPIC
${\mu }_{\alpha }$ (mas yr−1)	−9.8 ± 1.7	UCAC4
$\mu \delta$ (mas yr−1)	−24.2 ± 1.8	UCAC4
Photometry
B (mag)	17.353 ± 0.138	APASS
g' (mag)	16.386 ± 0.076	APASS
r' (mag)	14.860 ± 0.072	APASS
i' (mag)	13.360 ± 0.125	APASS
z' (mag)	12.613 ± 0.004	SDSS
J (mag)	11.095 ± 0.023	2MASS
H (mag)	10.332 ± 0.021	2MASS
Ks (mag)	10.026 ± 0.019	2MASS
W1 (mag)	9.890 ± 0.023	ALLWISE
W2 (mag)	9.762 ± 0.021	ALLWISE
W3 (mag)	9.610 ± 0.049	ALLWISE
Kinematics and Distance
Barycentric RV (km s−1)	−6.70 ± 0.15	This paper
U (km s−1)	−5.4 ± 0.5	This paper
V (km s−1)	−15.8 ± 2.2	This paper
W (km s−1)	−8.2 ± 1.2	This paper
Distance (pc)	145 ± 15	de Zeeuw et al. (1999)
Physical Properties
AV	${0.75}_{-0.18}^{+0.21}$	This paper
Spectral type	M3.3 ± 0.5	This paper
Rotation period (days)	6.29 ± 0.17	This paper
${T}_{\mathrm{eff}}$ (K)	3540 ± 70	This paper
${F}_{\mathrm{bol}}$ (10−10 erg cm−2 s−1)	2.25 ± 0.26	This paper
M* (${M}_{\odot }$)	${0.56}_{-0.09}^{+0.09}$	This papera
R* (${R}_{\odot }$)	${1.05}_{-0.07}^{+0.07}$	This papera
L* (${L}_{\odot }$)	${0.15}_{-0.03}^{+0.03}$	This papera
Age (Myr)	${9.3}_{-1.3}^{+1.1}$	This papera
$v\sin {i}_{* }$ (km s−1)	8.2 ± 1.8	This paper
i* (deg)	$\gt 63$	This paper

Note.

^aAdopted parameters from our model comparison; see Table 2 for more information.

Download table as:  ASCIITypeset image

5.1. Background EB

5.2. Companion EB

5.3. Eclipsing Binary

To confirm that the transiting body is non-stellar, we use the RVs from our IGRINS spectra to put an upper limit on the mass of K2-33b. We fit the RVs by assuming a circular orbit, locking the period and argument of periapsis from the transit fit (Section 3.2), and assuming the mass derived for the host star from our model interpolation (Section 4). The RVs rule out companion masses above 3.7 Jupiter masses (M_J) at $3\sigma$ (Figure 8). The constraints are tighter ($\lt 2.8{M}_{{\rm{J}}}$) if we remove the IGRINS epoch with high telluric contamination (see Section 2). If we loosen our assumptions about the orbital eccentricity then the maximum mass increases to 5.4M_J. In all cases the RVs exclude any brown dwarf or stellar companion ($\gt 13{M}_{{\rm{J}}}$) with an orbital period matching the transit signal.

Zoom In Zoom Out Reset image size

5.4. Stellar Variability

Spots and plages on the photosphere combined with stellar rotation create 1%–3% variations in the light curve of K2-33. The amplitude of this variation is roughly an order of magnitude larger than the transit depth (≃0.26%, Figure 3). Fortunately, spots create a characteristic shape in the light curve (smoothly varying) and duration (about half the rotation period) that differs from a transit (trapezoidal shape and a duration of hours). This makes them easy to differentiate in most stars. However, improper removal of the more complicated spot patterns on young stars can sometimes generate transit-like signals over short timescales (days or weeks). Our BLS search identified many such systems, one of which we show in Figure 9.

Zoom In Zoom Out Reset image size

The combination of our MEarth and K2 light curves demonstrates that the transit signal cannot be caused by stellar variability. While spot patterns can be stable on a multi-year baseline in old M dwarfs (Newton et al. 2016), spot evolution can be seen even over the 80 day observing window for Upper Scorpius members (e.g., Figures 1, 9). If the transit signal were an artifact of signals of stellar activity, then the transit depth and shape would change or disappear over the K2 observing window and between the K2 and MEarth observations. To test this we fit each transit individually as in Section 3.2, but locked the period to the previously derived value. We found that all transits yield consistent depths and durations, including MEarth data taken 1.5 yr after the K2 observations. Further, stellar signals are generally wavelength-dependent, and the separate fits to the K2 (${\lambda }_{{\rm{mean}}}\simeq 6400$ Å) and MEarth (${\lambda }_{{\rm{mean}}}\simeq 8250$ Å) light curves give consistent parameters.

5.5. Debris Disk?

Many stars in Upper Scorpius exhibit excess emission at infrared and millimeter wavelengths indicative of cooler, circumstellar material, i.e., dusty primordial or debris disks in different stages of their evolution (Luhman & Mamajek 2012). K2-33 may have excess emission at 24 μm (measured by the MIPS instrument on Spitzer) but shows no significant excess in any of the WISE bands or the Spitzer 8 and 16 μm bands (Carpenter et al. 2009; Luhman & Mamajek 2012). We confirmed this by comparing the Spitzer measurements from Luhman & Mamajek (2012) with our estimate of the photospheric flux (see Section 4), extrapolated to 30 μm using a PHOENIX BT-SETTL model (Figure 10). This excess is consistent with the WISE upper limit at 22 μm (W4 channel). This excess corresponds to that of debris disks, not an evolved disk, according to the classification of Luhman & Mamajek (2012).

Zoom In Zoom Out Reset image size

The lack of a detectable excess at wavelengths $\lt 12$ μm suggests that any disk, if it exists, must lack significant material warmer than 300 K, i.e., a central hole, or a drop in the emissivity of the grains close to the star. Using the value for L_* estimated in Section 4, we estimated that the hole extends to $\gt 0.35\,{\rm{A}}{\rm{U}}$. If the grains are small and not blackbody emitters, the emissivity will be higher at shorter wavelengths, thus strengthening this constraint. The central hole of any debris disk is significantly larger than the Keplerian orbit corresponding to the transit signal (0.051 ± 0.004 AU). This suggests that the planet and disk are physically separated and unrelated phenomena.

Disks around stars can create "dips" in the light curve as vertical structures in the disks periodically occult the host star (e.g., Cody et al. 2014). This behavior has been observed in some stars in Upper Scorpius (e.g., Ansdell et al. 2016). However, such dips are usually much deeper and are quasi-periodic or periodic, with the depth changing from dip to dip. The shapes of the dips typically do not resemble a transit; they are irregular and/or have leading or lagging tails (Ansdell et al. 2016). In contrast, the signals in the light curve of K2-33 are strictly periodic and transit-shaped, exhibiting no changes over the ∼1.5 yr interval between the K2 and MEarth observations. Finally all such "dipper" stars exhibit significant excess emission at 12 and 24 μm consistent with full or evolved disks, which is not the case for K2-33.

Young stellar associations like Upper Scorpius offer a unique view of the properties and behavior of young stars. Planets around these young stars are similarly critical probes into how planets change from formation to maturity (e.g., Aigrain et al. 2007; Quinn et al. 2012b; Bowler et al. 2015a; David et al. 2016a). In this paper, we have detailed our follow-up, characterization, and confirmation of K2-33b, a ${5.04}_{-0.37}^{+0.34}$ ${R}_{\oplus }$ planet orbiting at a period of 5.425 days around a star in the ∼11 Myr (Pecaut et al. 2012) Upper Scorpius OB association.

In combination, K2-33's proper motions, RV, and lithium and mid-IR excess in its spectrum unambiguously indicate that K2-33 is a young ($\lt 20$ Myr) PMS star associated with the Upper Scorpius association. We used moderate-resolution spectra to revise the determinations of reddening, ${T}_{\mathrm{eff}}$, and ${F}_{\mathrm{bol}}$. We fit the K2 and MEarth transit photometry, which also yielded a precise stellar density. By interpolating these constraints onto a grid of PMS models that imitate magnetic effects on the star's internal structure, we derived a precise radius (6%–7%) and mass (16%) for K2-33.

The data argue strongly that the transit signal is planetary in origin. AO and radial velocities rule out a background or bound EB as the source of the transit signal. Our MEarth transit photometry combined with the K2 photometry rules out stellar variability or a disk mimicking a transit. Further, K2-33b exhibits none of the unusual light-curve behavior of PTFO 8-8695b, the candidate whose V-shaped transit signature exhibits time-variable depth and width (van Eyken et al. 2012; Ciardi et al. 2015; Yu et al. 2015).

Few young ($\lt 20$ Myr) stars have their masses/radii determined to better than 10%, with the exception of young eclipsing binaries (e.g., Kraus et al. 2015; David et al. 2016b). The precision is in part due to our careful measurement of ${F}_{\mathrm{bol}}$ and ${T}_{\mathrm{eff}}$, and the additional constraint from the fit to the transit light curve, which provides a stellar density accurate to ≲10%. This highlights the power of transiting exoplanets to probe stellar astrophysics. While our current method relies on a stellar model, when Gaia parallaxes become available (Perryman et al. 2001; de Bruijne 2012) we can instead use transiting planets to test these models.

Our stellar density assumes that the planet has zero orbital eccentricity. This is expected for a young Neptune-mass planet where recent interactions with the nascent disk have dampened eccentricity (e.g., Tanaka & Ward 2004; Cresswell et al. 2007). However, there is a dearth of young planetary systems with known eccentricities that would be necessary to confirm this observationally. To test the sensitivity of our results to the assumption of e = 0 we reran our stellar isochrone fits, but instead used a uniform density on ${\rho }_{* }$. This gave a best-fit stellar density and radius consistent with our earlier determination, but with errors larger by a factor of two (see Table 2, Figure 7). Similarly, we derived a consistent stellar radius when using the Stefan–Boltzmann relation, which is independent of both the transit-fit density and stellar models. We conclude that the accuracy of our stellar (and therefore planetary) radius is insensitive to this assumption.

Assuming a Neptune-like density, K2-33b will have an RV amplitude of ∼20 m s⁻¹. While this is well within detection limits of current RV instruments, it is smaller than the expected spot-induced RV jitter (100–200 m s⁻¹). Moving to the NIR can reduce this noise source but not eliminate it (Mahmud et al. 2011; Crockett et al. 2012). Because the planet and stellar rotation periods are known from the light curve, it may be possible to fit each separately with sufficient data and baseline. K2-33 will be an excellent target for upcoming NIR RV spectrographs (e.g., Quirrenbach et al. 2010; Artigau et al. 2014; Kotani et al. 2014).

As with K2-25b, K2-33b is considerably larger than close-in planets found around stars of similar mass by Kepler. Most planets around M dwarfs found by Kepler are 1–2.5 ${R}_{\oplus }$(Morton & Swift 2014; Dressing & Charbonneau 2015; Mulders et al. 2015; Gaidos et al. 2016b), while K2-33b is roughly twice this size at ${5.04}_{-0.37}^{+0.34}$ ${R}_{\oplus }$. K2-33b is less of a radius outlier than K2-25b, which orbits a $0.3\,{M}_{\odot }$ star. K2-33 has a mass of $0.55\,{M}_{\odot }$ and larger planets are more common around higher-mass hosts. Further, unlike with the nearby, bright, main-sequence, and photometrically well-behaved stars in the Hyades, it is not clear whether our survey is sensitive to the more common 1–2 ${R}_{\oplus }$ planets around host stars of similar mass in Upper Scorpius. However, K2-33b fits into an emerging picture in which young planets are larger than their older counterparts. Mann et al. (2016) suggested these large radii could be due to the initial heat of formation as well as inflation and escape of the atmosphere under the influence of the young, active host star (e.g., Rogers et al. 2011; Ehrenreich et al. 2015).

The upper limit on K2-33b's age provided by its ≃11 Myr stellar host suggests that it either migrated inwards via disk migration or formed in situ, because planet–star and planet–planet interactions work on much longer timescales (Fabrycky & Tremaine 2007; Nagasawa et al. 2008), and the conditions for Kozai–Lidov evolution begin only after the disk dissipates (e.g., Martin et al. 2016). This discovery makes it unlikely that such long-term dynamical interactions are responsible for all close-in planets. However, it is difficult to draw conclusions about the dominant migration or formation mechanism for close-in planets given the sample size and incomplete understanding of our transit-search pipeline's completeness.

Selection effects may be important here. It is possible that K2-33 has an atypical history of formation and migration that also made it the easiest (and hence first) planet of this age to be identified and confirmed. A full search of all young clusters and stellar associations surveyed by the K2 mission, with proper treatment of detection completeness, is underway. This, along with improved statistics provided by the TESS and PLATO missions, will provide an estimate of the planet occurrence rate as a function of time. Trends (or a lack of trends) in this occurrence rate could set constraints on planetary migration timescales.

During the final stages of the analysis for this paper we were informed by another team that an independent analysis of this system was about to be submitted (David et al. 2016c).

We thank the anonymous referee for their thoughtful comments on the paper. The authors thank Ann Marie Cody for her help with the K2 light curves when initially identifying the transit. We also thank Charlie for his help with the data analysis and proofreading the paper.

A.W.M. was supported through Hubble Fellowship grant 51364 awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS 5-26555. This research was supported by NASA grant NNX11AC33G to E.G. K.R.C. acknowledges support provided by the NSF through grant AST-1449476. A.V. is supported by the NSF Graduate Research Fellowship, grant No. DGE 1144152.

The MEarth Team gratefully acknowledges funding from the David and Lucille Packard Fellowship for Science and Engineering (awarded to D.C.). This material is based in part upon work supported by the National Science Foundation under grants AST-0807690, AST-1109468, and AST-1004488 (Alan T. Waterman Award). 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.

The ARCoIRIS observations were conducted at Cerro Tololo Inter-American Observatory, National Optical Astronomy Observatory, which is operated by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation. This work used the Immersion Grating Infrared Spectrograph (IGRINS) that was developed under a collaboration between the University of Texas at Austin and the Korea Astronomy and Space Science Institute (KASI) with the financial support of the US National Science Foundation under grant ASTR1229522, of the University of Texas at Austin, and of the Korean GMT Project of KASI. The IGRINS pipeline package PLP was developed by Dr. Jae-Joon Lee at Korea Astronomy and Space Science Institute and Professor Soojong Pak's team at Kyung Hee University. SNIFS on the UH 2.2 m telescope is part of the Nearby Supernova Factory project, a scientific collaboration among the Centre de Recherche Astronomique de Lyon, Institut de Physique Nucléaire de Lyon, Laboratoire de Physique Nucléaire et des Hautes Energies, Lawrence Berkeley National Laboratory, Yale University, University of Bonn, Max Planck Institute for Astrophysics, Tsinghua Center for Astrophysics, and the Centre de Physique des Particules de Marseille. Some of the data presented in this paper were obtained from the Mikulski Archive for Space Telescopes (MAST). STScI is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS5-26555. Support for MAST for non-HST data is provided by the NASA Office of Space Science via grant NNX09AF08G and by other grants and contracts. This research was made possible through the use of the AAVSO Photometric All-Sky Survey (APASS), funded by the Robert Martin Ayers Sciences Fund.

The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain.

Facilities: Blanco (ARCoIRIS) - Cerro Tololo Inter-American Observatory's 4 meter Blanco Telescope, UH:2.2 m (SNIFS) - , Keck:II (NIRC2) - , Smith (IGRINS) - McDonald Observatory's 2.7m Harlan J. Smith Telescope, Kepler - The Kepler Mission.

Software: Spextool (Cushing et al. 2004), xtellcor (Vacca et al. 2003), IGRINS pipeline (Lee 2015), ARCoIRIS Spextool, batman (Kreidberg 2015), emcee (Foreman-Mackey et al. 2013).