Community targets for JWST's early release science program: evaluation of WASP-63b

TL;DR

This study evaluates WASP-63b as a JWST community target, analyzing HST data to detect water vapor and assess atmospheric properties, aiding in the selection of optimal exoplanets for JWST's early science observations.

Contribution

It provides the first near-infrared spectrum of WASP-63b, demonstrating water detection and exploring atmospheric composition and cloud effects for JWST target vetting.

Findings

Detection of water absorption feature with 3.5-4.0 sigma significance

Muted water feature suggests clouds or composition effects

Possible indication of super-solar HCN abundance

Abstract

We present observations of WASP-63b by the Hubble Space Telescope (HST) as part of "A Preparatory Program to Identify the Single Best Transiting Exoplanet for JWST Early Release Science". WASP-63b is one of the community targets under consideration for the James Webb Space Telescope (JWST) Early Release Science (ERS) program. We present a spectrum derived from a single observation by HST Wide Field Camera 3 in the near infrared. We engaged groups across the transiting exoplanet community to participate in the analysis of the data and present results from each. There is general agreement amongst all results that we find an H2O absorption feature with 3.5-4.0 sigma significance. However, the feature is muted in comparison to a clear atmosphere at solar composition. Although the detection of the water feature is robust, the reasons for the muting of this feature are ambiguous due to a degeneracy between clouds and composition. The data does not yield robust detections of any molecular species other than H2O. The group was motivated to perform an additional set of retrieval exercises to investigate an apparent bump in the spectrum at ~ 1.55 um. We explore possible disequilibrium chemistry and find this feature is consistent with super-solar HCN abundance but it is questionable if the required mixing ratio of HCN is chemically and physically plausible. The ultimate goal of this study is to vet WASP-63b as a potential community target to best demonstrate the capabilities and systematics of JWST instruments for transiting exoplanet science. In the case of WASP-63b, the presence of a detectable water feature indicates that WASP-63b remains a plausible target for ERS observations.

Tables

Table 1. Table 2 Transmission spectrum of WASP-63b measured with HST WFC3 G141 grism.

BMK			KBS			HRW
Wavelength($\mu$m)	Rp/R⋆	Unc (ppm)	Wavelength($\mu$m)	Rp/R⋆	Unc (ppm)	Wavelength($\mu$m)	Rp/R⋆	Unc (ppm)
1.1425	0.07791	310	1.1425	0.07819	290	1.1352	0.07810	460
1.1775	0.07732	290	1.1775	0.07801	275	1.1677	0.07740	450
1.2125	0.07717	300	1.2125	0.07785	290	1.2002	0.07762	430
1.2475	0.07753	300	1.2475	0.07770	270	1.2327	0.07839	420
1.2825	0.07812	280	1.2825	0.07848	260	1.2651	0.07832	420
1.3175	0.07789	270	1.3175	0.07846	285	1.2977	0.07782	420
1.3525	0.07847	290	1.3525	0.07853	265	1.3301	0.07832	420
1.3875	0.07891	300	1.3875	0.07861	270	1.3626	0.07921	420
1.4225	0.07832	290	1.4225	0.07826	290	1.3951	0.07880	425
1.4575	0.07839	300	1.4575	0.07849	280	1.4275	0.07848	440
1.4925	0.07773	330	1.4925	0.07791	290	1.4600	0.07862	445
1.5275	0.07865	330	1.5275	0.07831	295	1.4925	0.07811	450
1.5625	0.07866	370	1.5625	0.07833	320	1.5250	0.07908	455
1.5975	0.07816	350	1.5975	0.07752	310	1.5575	0.07892	485
1.6325	0.07751	360	1.6325	0.07718	320	1.5899	0.07821	490
						1.6224	0.07789	505

Full text

††slugcomment: Accepted for publication in AJ June 18, 2018

1. Introduction

The James Webb Space Telescope will revolutionize transiting exoplanet atmospheric science due to a combination of its capability for continuous, long duration observations and its larger collecting area, spectral coverage, and resolution compared to existing space-based facilities. We previously outlined a plan in Stevenson et al. (2016) to fully demonstrate the capabilities of the JWST instruments during the Early Release Science (ERS) program allowing the community to plan more efficient and successful transiting exoplanet characterization programs in later cycles.

Stevenson et al. (2016) identified a set of “community targets” which meet a certain set of criteria for ecliptic latitude, period, host star brightness, well constrained orbital parameters, and predicted strength of spectroscopic features. WASP-63b was identified as one of the strongest transmission spectroscopy candidates for JWST Early Release Science. It is an inflated planet ($1.43$ RJ) with a low mass of only $0.38$ MJ (Hellier et al., 2012) resulting in a large atmospheric scale height. It orbits a bright (11.2 Vmag) star. Additionally, WASP-63b occupies an important, underexplored, region of transmission spectroscopy due to its mass. Most exoplanets studied in detail with transmission spectroscopy are either hot Jupiters of mass ($\sim$1–3 MJ)(e.g. Deming et al., 2013; Kreidberg et al., 2014a, 2015a; Line et al., 2016; Sing et al., 2016) or Super Earth-to-Neptune mass planets ($\sim$0.01–0.1 MJ) (e.g. Kreidberg et al., 2014b; Fraine et al., 2014; Knutson et al., 2014). In order to understand formation and evolution processes, it’s important to understand the composition of atmospheres over a full and continuous range of masses (e.g. Mordasini et al., 2016).

WASP-63b will be accessible to JWST approximately six months after the planned April 2019 start of Cycle 1 and ERS observations making it an ideal candidate should there be any delays in the JWST timetable. Here, we observe WASP-63b to evaluate its suitability as a prime candidate to test the capabilities of JWST. We can use the strength of the water absorption feature at 1.4 $\mu$m as a way to screen potential targets for the presence of obscuring aerosols and determine the amplitude of predicted spectral features (e.g. Deming et al., 2013; Kreidberg et al., 2014a; Sing et al., 2016; Stevenson, 2016). Ideally, a clear atmosphere with large amplitude spectroscopic features will be best suited for benchmarking the instruments and identifying their systematics.

2. Observations and Analysis

We observed the WASP-63 system using the HST Wide Field Camera 3 (WFC3) on September 19, 2016. The observations were taken as part of program GO-14642, (PI Stevenson). The observations were made using the G141 grism in forward/reverse spatial scan mode. Spatial scanning (McCullough & MacKenty, 2012) involves slewing the telescope in the cross dispersion direction during the exposure. In forward/reverse mode the telescope is exposing in both directions of the slew thus eliminating time to reset the target at the initial position on the detector between exposures. Each exposure, utilizing SPARS10, consists of 16 non-destructive reads with a total exposure time of approximately 103s which yielded peak per pixel counts near 32,000 electrons. The spectrum was read out using the 256 $\times$ 256 subarray with a scan rate of $\sim$ 0.08 arcsec/s (0.62 pixels/s). This corresponds to a total scan length of $\sim$ 8.76 arcsec which spreads the spectrum in the cross dispersion direction over $\sim$ 70 pixels. We observed WASP-63 for a total of 8 HST orbits, yielding a total of 86 time series integrations, to cover the entirety of the relatively long duration of transit ($\sim$ 5 hours).

We use the IMA files from the CalWF3 pipeline in our analysis. These files have been calibrated for dark current and zero read bias. We applied flat field corrections to each non-destructive read (NDR). Each NDR was background subtracted by considering a background window consisting of ($\sim$50) rows of pixels below the spectrum spanning the dispersion direction. A mean value for each column of the background window was taken to produce a one-dimensional background correction. The 1-D solution was then smoothed in the dispersion direction to correct for outliers. The column by column background value was then subtracted from each pixel of the image. We then extract the spectrum by taking the difference between successive NDRs. We apply a top hat filter to each NDR to limit any contribution from cosmic rays and/or overlapping spectra (Evans et al., 2016; Wakeford et al., 2017). The differences between each NDR are then summed to produce a final working image.

2.1. Band Integrated Light Curve

We perform the extraction of the band integrated light curve (white light curve) using a range of different aperture sizes in the cross-dispersion direction. The aperture sizes range from $\pm$10 pixels–$\pm$10 Each orbit includes a direct image of the star from which we calculate the stellar centroid using applying a 2-D Gaussian fit to a 5$\times$5 section of pixels centered on maximum pixel. The trace and wavelength solutions are calculated from the centroid of the undispersed image using the coefficients provided in Kuntschner et al. (2009). We extract a box from each working image with the number of columns corresponding to wavelength limits 1.125-1.65 $\mu$m and rows determined by the chosen aperture size. The band integrated light curve is the summation of all pixels within the box at each time step. The results of the initial extraction of the raw white light curve are shown in Figure 1 (top).

Fitting the white light curve requires accounting for HST systematics. We choose to follow the standard practice of discarding the first orbit as it presents different systematics from the remainder of the data (e.g. Deming et al., 2013; Stevenson et al., 2014). The raw light curve exhibits a ramp like increase in flux, commonly referred to as the ‘hook’, with each HST orbit consistent with previous observations (e.g. Berta et al., 2012; Deming et al., 2013; Fraine et al., 2014; Kreidberg et al., 2014a). The hook effect, shown in Figure 1 (middle), is generally steeper in the first frame of each HST orbit so we discard those data points. We then model the hook as an exponential plus linear function of the form $1-A\exp\left\{S(\theta-\theta_{0})\right\}+c_{1}\theta$ where $\theta$ is the HST orbital phase, $\theta_{0}$ is a reference angle for setting zero HST phase and A, S, and c1 are scaling factors. The hook model is combined multiplicatively with a second order polynomial in time over the entirety of the observation.

We model the transit using the methods of Mandel & Agol (2002) implemented by the Python routine BATMAN Kreidberg (2015). The transit model assumes a circular orbit with a fixed period. Orbital parameters used for the transit model were taken from Hellier et al. (2012). We calculate non–linear limb darkening coefficients using the PHOENIX Code to fit theoretical spectra as described in detail in de Wit et al. (2016). We assume a stellar Teff of 5550$\pm$100 K and log(g) of 4.01 +0.02/-0.04 producing coefficients c1…c4 of (0.36627439, 0.63188403, -0.69135111, 0.23393244) for a limb darkening law of the form $I(\mu)=I_{0}[1-c_{1}\left(1-\mu^{1/2}\right)-c_{2}\left(1-\mu\right)]-c_{3}\left(1-\mu^{3/2}\right)$ $-c_{4}\left(1-\mu^{2}\right)$. During the fitting process we allow for the time of transit center (Tcen), planetary radius as a fraction of stellar radius (Rp/R⋆), a/R⋆, inclination (i), and a normalizing factor (the ratio of the amplitudes of the scan directions) to be free parameters and fit both scan directions simultaneously. All fits and uncertainty estimates are derived from the Python routine ‘emcee’ (Foreman-Mackey et al., 2013) utilizing 50 walkers with $10^{4}$ steps. Each walker is initialized at a point in parameter space determined by randomly shifting the value of each parameter by some fraction of twice the corresponding one sigma uncertainty. This forms an over–distributed sphere of starting points surrounding the initial best fit value. Given the number of observational data during both ingress and egress we place uninformative priors on the orbital properties. The first 20% of steps are removed for burn in. We test for convergence using Gelman Rubin statistics with a threshold for acceptance of 1.1 and find convergence to levels of 1.01–1.03 in all cases (Gelman & Rubin, 1992). We choose the best aperture by minimizing the scatter of the residuals of the white light curve fit. We find a best aperture of $\pm$ 44 pixels in the spatial direction centered on the spectral image centroid. We achieve a standard deviation of the normalized residuals (SDNR) of 96 parts per million (ppm). There is no indication that the selection of aperture size has a statistically significant effect on the transit depth fit. The standard deviation of the transit depth for the five apertures above and below the best aperture is 93 ppm with an average Rp/Rs of 0.077761. These results all fit well within our best fit value and uncertainty. The best fit white light curve is shown in Figure 1 (bottom) and the values and uncertainties for the planetary physical properties and ephemeris from the best fit white light curve are presented in Table 2.1.

2.2. Spectral Light Curves

The spectral light curves are extracted using the same aperture we found to minimize the SDNR of the white light curve. The range of wavelengths included in our aperture are divided into 15 bins of width $0.035$ $\mu$m. The spectrum from each frame is compared to the spectrum of the first frame using cross correlation in Fourier space to check for any shift in the wavelength–pixel solution. The shift in wavelength solution throughout the observation was on the order of a few tenths of a pixel. Each column was summed and weighted by the fraction of that pixel in the bin.

The systematics were removed from the spectral light curves using the divide white method (Deming et al., 2013; Kreidberg et al., 2014a; Stevenson et al., 2014). Removing the best fit transit from the white light curve leaves only the systematics. We divide each spectral light curve by the systematics. This assumes that the systematics are wavelength independent. We do note a linear, observation–long, wavelength dependence in the corrected spectral light curves. We use a first order polynomial to account for the wavelength dependent systematics combined with the transit model when fitting the spectral light curves. We fix the ephemerides to the white light curve solutions and use fixed, wavelength dependent, non-linear limb darkening coefficients derived in the same way as described in Section 2.1. The transit depth and normalization factors are left as the only free parameters. We fit both scan directions simultaneously. Spectral light curves and their fits are shown in Figure 2.

2.3. Transmission Spectrum

A transmission spectrum was derived from the transit depth fits of the spectral light curves. The change in the apparent planetary radius as a function of wavelength can be indicative of absorption features of molecular species in the planetary atmosphere. As a test for robustness, the spectrum of WASP-63b was extracted using multiple independent analysis pipelines in addition to the method described in detail in the previous subsections (Table 2). Figure 3 shows a comparison of the results from this methodology (BMK) with that of analysis performed using methods described in Stevenson et al. (2014) (KBS) and Wakeford et al. (2016, 2017) (HRW). The spectra are in good agreement with nearly all points from each of the methods in 1$\sigma$ agreement with each of the other methods. As a test to determine how small differences in the extracted spectrum would affect retrieval results we analyzed the KBS and HRW spectra with the Pyrat Bay model (Section 3.2.1). The KBS spectrum has a shallower bump at 1.5 $\mu$m, and thus, the retrieval only detected the water signature significantly (not HCN). The HRW atmospheric retrieval produced similar results to BMK, but with broader 1$\sigma$ constraints, due to the larger uncertainties of the dataset (compared to BMK or KBS).

We compared these results against a flat-line model using the Bayesian Information Criterion (BIC, Liddle 2007). For the atmospheric model of each dataset we searched for the combination of free parameters that minimized BIC.

For KBS, the BIC favors the atmospheric model (BIC=24.59, with three free parameters, $T$, $R_{\rm 0.1bar}$, and H$\_2$O), against the flat-line model (BIC=30.82, one free parameter). The posterior odds give the flat-line model a fractional probability of only 0.04 (Raftery 1995).

For HRW, the BIC favors again the atmospheric model (BIC=19.19, with four free parameters, $T$, $R_{\rm 0.1bar}$, H$\_2$O, and HCN) over the flat-line model (BIC=22.92). However, this time there is weaker evidence, with the flat-line model having a fractional probability of 0.13.

Therefore, in all three datasets, there is evidence of spectral features, but with different degrees of confidence. The BMK spectrum shows strong evidence, KBS spectrum shows moderate evidence, whereas HRW shows only weak evidence for water absorption. The HRW spectrum is based upon the marginalization technique described in Wakeford 2016 and represents the most conservative estimate of uncertainty. As such, it serves as a sort of lower bound on the detection threshold whereas the BMK and KBS methods assume we are justified in our choice of systematic model. The arguments for this choice of model are well extablished and summarized in Section 2.1 and so we will choose the BMK spectral extraction to continue on with retrieval analysis.

3. Results

The results produced by the methodology described in in detail in Section 2 (presented in Table 2 as BMK) were distributed to the members of the transiting exoplanet community who were involved with the preparation of Stevenson et al. (2016) and/or HST program GO 14642. Each was given an opportunity to provide an independent analysis of the results. Wide community involvement resulted in a number of contributions in the form of forward model comparisons and retrievals. Here we present the methods and findings from each interpretation of the BMK spectrum. All models in this section adopt the system parameters from Hellier et al. (2012).

3.1. Forward Models

3.1.1 Burrows

We apply the transit models from Howe & Burrows (2012), which adopt chemical equilibrium abundances for molecular species from Burrows & Sharp (1999) and opacities from Sharp & Burrows (2007). The atmospheric models consider an isothermal temperature profile and gray haze opacity with cross sections of 0.001–0.005 cm2g*-1* from $10^{-6}$ bar to 1 bar.

By exploring a range of temperatures, haze opacities, metallicities, and non-equilibrium CO/CH4 abundances, the best-fitting solutions pointed to solar-abundance atmospheres at a temperature of 1000 K, with a haze/cloud muting the water absorption feature at 1.4 $\mu$m (Figure 4, top panel). There is no indication of significant CH4. There is a slight degeneracy between the cloud thickness and temperature, but it is clear that the atmosphere is cloudy. These models do not show a significant metallicity dependence. Finally, a high CO abundance excess ($\sim$100 times solar) can help to fit the data at 1.5 $\mu$m, but it does not seem realistic.

3.1.2 Heng

As a complementary approach to the full retrieval calculations, we fit the data with a 3-parameter analytical model (Heng & Kitzmann, 2017). In that study, it was demonstrated that this isothermal, isobaric model matched full numerical calculations at the $\sim 0.1\%$ level over the WFC3 wavelength range. The model has three parameters: temperature, water abundance, and a constant cloud opacity. The constant cloud opacity assumes that the cloud particles are large over the wavelength range probed by WFC3 (i.e., micron-sized or larger radius). Water opacities are computed using the HELIOS-K opacity calculator (Grimm & Heng, 2015) and the HITEMP spectroscopic database. This procedure confronts the data with a simple model, which has a minimal number of parameters, to serve as a plausibility check.

Following the approach of Kreidberg et al. (2015b), we equate the reference transit radius to the white-light radius and set the reference pressure to 10 bar. We assume a hydrogen-dominated atmosphere and set the mean molecular weight to 2.4. The top panel of Figure 4 shows the best-fit model to the WFC3 WASP-63b data. Our general conclusion mirrors that of the retrieval calculations: water is present in the atmosphere of WASP-63b, but its presence is muted by a continuum, which in this case is attributed to a constant cloud opacity. The values of our fitting parameters span a temperature range from 500 to 1000 K, a water mixing ratio from $\sim 10^{-8}$ to $10^{-7}$, and a cloud opacity $\sim 10^{-8}$ to $10^{-7}$ cm2 g*-1*.

3.1.3 Morley

In order to determine the clouds that are predicted to form in the atmosphere of WASP-63b and their effect on the planet’s transmission spectrum, we ran self-consistent models including the effects of cloud condensation. These models solve for the temperature structure of the atmosphere in radiative-convective and chemical equilibrium and are more extensively described in McKay et al. (1989); Marley et al. (1996); Burrows et al. (1997); Marley et al. (1999, 2002); Fortney (2005); Saumon & Marley (2008); Fortney et al. (2008); Morley et al. (2015). Our opacity database for gases is described in Freedman et al. (2008, 2014), and we calculate the effect of cloud opacity using Mie theory, assuming spherical particles. We include iron and silicate clouds and vary the cloud sedimentation efficiency $f_{\rm sed}$ from 0.1 to 1, and find that these clouds do indeed form at high altitudes and damp the size of the signal for low sedimentation efficiencies (i.e. lofted clouds of small particles). Figure 4 top panel shows a representative transmission model for WASP-63b.

3.1.4 Parmentier

In order to understand how the three-dimensional structure of the planet might affect our interpretation of the planet’s transmission spectrum, we model WASP-63b with the three-dimensional global circulation model SPARC/MITgcm described in Showman et al. (2009). Our model solves for the three-dimensional temperature structure of the atmosphere assuming a cloud-free, solar-composition atmosphere. We then use the temperature map to predict the position of the clouds at the limb of the planet by comparing the partial pressure and the saturation pressure of the cloud gaseous constituents as described in Parmentier et al. (2016). The cloud top level and size of the cloud particles are free parameters representing vertical mixing and microphysics respectively. We then compute the transmission spectrum of the whole atmosphere by combining the transmission spectrum obtained with the temperature and cloud profile at each latitude around the limb (Parmentier et al., 2018).

Our global circulation model predicts a temperature difference of $400\rm K$ between the east and west limb at the 10 mbar level. As a consequence some cloud species are predicted to be condensed all over the limb of the planet whereas others should condense only on the morning limb and be evaporated on the other one, leading to a partially cloudy atmosphere (Line & Parmentier, 2016). We computed models assuming the presence of enstatite clouds ($\rm MgSiO_{3}$) or manganese sulfide clouds (MnS) corresponding to a fully cloudy and partially cloud case respectively. Our best fit spectrum with enstatite clouds is very similar to the Morley model of Fig.4. It has a cloud top pressure of 1 mbar and no constraints on the particle size. Our best fit model with MnS clouds is the Parmentier model of Fig. 4. It has a limb that is $\approx 60\%$ cloudy, resulting in a qualitatively different spectrum than the other, one-dimensional models shown here. For this model the cloud top pressure is $\approx 0.1$ mbar and the particle size is $\approx 1\mu$m. We conclude that the atmosphere of WASP-63b is unlikely to be clear with a solar-composition abundance of water. Both fully cloudy and partially cloudy atmospheres can exist, depending on the cloud composition. A higher signal to noise spectrum should be able to disentangle between the two scenarios.

3.2. Retrievals

Four groups provided atmospheric retrieval analyses for WASP-63b. The following subsections describe the retrieval procedure and the individual exploration from each group. Figures relative to the individual retrievals are included in the Appendix.

3.2.1 Blecic & Cubillos

To model the atmosphere and spectrum of WASP-63b we use the Python Radiative Transfer in a Bayesian framework (Pyrat Bay) package111http://pcubillos.github.io/pyratbay (Cubillos et al. 2018, in prep.; Blecic et al. 2018, in prep.). Pyrat Bay is an open-source reproducible code, based on the Bayesian Atmospheric Radiative Transfer package (Blecic, 2016; Cubillos, 2016). The code provides a line-by-line radiative-transfer and a thermochemical-equilibrium abundances (TEA, Blecic et al., 2016) module, which can be use in forward or retrieval mode, via a Differential-evolution MCMC sampler (Cubillos et al., 2017). The atmospheric model consist of a 1D set of concentric shell layers, in hydrostatic equilibrium, spanning from $10^{-8}$ to $100$ bar. For the temperature profile we consider either the three-channel Eddington approximation parameterization (TCEA, Line et al., 2013b) or an isothermal profile.

The Pyrat Bay code considers molecular opacities for H2O from HITEMP (Rothman et al., 2010), NH3 and CH4 from HITRAN (Rothman et al., 2013), and HCN from Exomol (Barber et al., 2014). We compressed these line-by-line data files using the repack package (Cubillos, 2017) to extract only the strong lines that dominate the spectrum, speeding up the radiative-transfer calculation. Additionally, Pyrat Bay considers collision induced absorption (CIA) from H2-H2 (Borysow et al., 2001; Borysow, 2002) and H2-He (Borysow et al., 1988, 1989; Borysow & Frommhold, 1989); and H2 Rayleigh scattering (Lecavelier Des Etangs et al., 2008). We also consider two cloud parameterizations, a simple gray-cloud opacity with constant cross section (cm*-2* molec*-1*) below $10^{-5}$ bars, and a thermal-stability cloud model based on the approach described in Ackerman & Marley (2001) and Benneke (2015), with additional flexibility (Blecic et al. 2018, in prep). We compute the opacity for either Fe or MgSiO3 condensates using Mie-scattering theory (Toon & Ackerman, 1981).

The retrieval parameterization includes free scale factors for the mixing ratios of H2O, NH3, CH4, and HCN (vertical-constant values) and the mean molecular mass of the atmosphere; either the isothermal temperature or the $\kappa$, $\gamma$, and $\beta$ parameters of the TCEA model (see Line et al., 2013b); the gray cloud cross section or the Mie-scattering cloud profile shape, condensate mole fraction, particle-size distribution, and gas number fraction just below the cloud deck.; and the planetary radius at 0.1 bar.

We explored several cases, obtaining qualitatively good fits in all gray-cloud, complex-cloud, and clear-atmosphere cases. As expected for transmission spectroscopy, the retrieval returned largely unconstrained parameters for the TCEA temperature model, suggesting that the data does not justify more complex models than an isothermal profile. In all cases the MCMC favors lower temperatures ($T<1000$ K) than equilibrium temperature (1500 K) at the pressures probed by the observations. We constrain the water abundance, ranging from solar to $\sim$0.1 times solar values. The observed water absorption feature is muted relative to a clear atmosphere with solar abundances. This is caused by a sub-solar water abundance, an absorbing cloud opacity, or a high mean molecular mass, which reflects in a strong correlation between the water abundance and the cloud cross-section. When we compare retrievals with the gray and complex cloud model, we find similar best-fitting spectra between the two cases. The complex-cloud retrieval does not constrain any of the cloud parameters when we set all four species abundances free. In the case when we set the water abundance as the only abundance free parameter, we find a somewhat better condensate-fraction constraint. Since the cloud opacity dominates only a limited region of the observed spectrum ($\sim$1.2–1.3 $\mu$m), we conclude that there is no need for a more complex cloud model for this study. The Reproducible Research Compendium for the Pyrat Bay models is available at https://github.com/pcubillos/KilpatrickEtal2018_WASP63b.

3.2.2 Line

We use the CHIMERA transmission model (Line et al., 2013a; Swain et al., 2014; Kreidberg et al., 2014a, 2015b; Greene et al., 2016; Line et al., 2016). For transit geometry, the code solves the radiative-transfer equation for parallel rays across the terminator of the planet (Brown, 2001; Tinetti et al., 2012). The code integrates atmospheric opacities from either correlated-K or sampled “line-by-line” absorption cross sections. To explore the parameter space, the transmission model is coupled with the PyMultiNest (Buchner et al., 2014) multimodal nested-sampling algorithm (Feroz & Hobson, 2008). The atmosphere is in hydrostatic equilibrium with height-dependent gravity, temperature, and molecular weight. The temperature profile comes from either the radiative-equilibrium model from Guillot (2010) or an isothermal profile, spanning from $10^{-7}$ to 30 bar. The atmosphere uses ’thermochemically-consistent’ molecular abundances (as defined in Kreidberg et al., 2015b), computed using the NASA CEA2 model for given elemental abundances (Lodders, 2009).

The CHIMERA code considers molecular opacities for H2O, CH4, CO, CO2, NH3, Na, K, TiO, VO, C2H2, HCN, and FeH (Freedman et al., 2014); CIA from H2-H2 and H2-He Richard et al. (2012); a Rayleigh power-law haze (Lecavelier Des Etangs et al., 2008); and either an opaque gray patchy cloud model (Line & Parmentier, 2016) or a more complex, Mie-cloud model (Lee et al., 2013).

The retrieval parameterization includes the metallicity $\rm[M/H]$, the carbon-to-oxygen ratio $\log(\rm C/O)$, and the carbon- and nitrogen-species quench pressures (Kreidberg et al., 2015b; Morley et al., 2017) to set the elemental abundances; either the isothermal temperature or the $\kappa_{\rm IR}$, $\gamma_{v}$, $T_{\rm irr}$ parameters of (Guillot, 2010); the top-pressure boundary and a ’patchy terminator’ parameter for the gray patchy cloud model, or the $Q_{0}$ and $r$ (see Lee et al., 2013) and profile shape (mixing ratio, cloud base pressure, pressure fall off index, Line et al. 2017, in prep.) for the Mie-cloud model; and a radius scale factor to set the reference altitude at 10 bar.

The ’chemically-consistent’ retrieval detects the water spectral feature at $3.6\sigma$ confidence. This is consistent with a hard upper limit on C/O near 1. The water band is muted relative to solar composition. The retrieval posterior shows two ’composition’ modes: low metallicity (${\rm[M/H]}\lesssim 1.3$ (20$\times$)) degenerate with a cloud and high metallicity (peak near $\sim$300$\times$solar). There is a turn-over degeneracy in cloud top vs. [M/H] (due to the effect on the mean molecular weight) resulting in the bi-modal marginalized metallicity distribution. Clouds can be present, but are not required to fit the spectra as given by the Bayes factor (0.45) and result in a much lower value for the low metallicity mode (¡0.1$\times$solar), while the high metallicity mode remains. The highest of the sampled metallicities (greater than $\sim$50 times solar) are possibly implausible given mass and radius of WASP-63b.

Further tests found negligible variations when imposing a temperature prior or no patchy-cloud factor. A comparison between correlated-K and line-by-line sampling opacities produced nearly identical results. Likewise, the more complex Mie-cloud model produced qualitatively similar main conclusions (with unconstrained cloud particle sizes, vertical extent, or cloud composition).

3.2.3 MacDonald & Madhusudhan

We use the nested-sampling retrieval algorithm POSEIDON (MacDonald & Madhusudhan, 2017) to analyze the WFC3 observations of WASP-63b. The code pre-computes line-by-line molecular cross sections following the methodology of Hedges & Madhusudhan (2016) and & Madhusudhan (2017). To compute detection significances we conduct nested Bayesian model comparisons. For simplicity, we model the atmosphere assuming an isothermal temperature-pressure profile, with 100 layers uniformly spaced in log-pressure from $10^{-6}$ to $100$ bar, in hydrostatic equilibrium.

The POSEIDON code considers molecular opacities for H2O from HITEMP (Rothman et al., 2010) and CH4, NH3, and HCN from Exomol (Tennyson et al., 2016); CIA from H2-H2 and H2-He (Richard et al., 2012); Rayleigh scattering (Lecavelier Des Etangs et al., 2008); and a uniform-in-altitude gray opacity cloud model.

The retrieval parameterization includes free scale factors for the mixing ratio of H2O, CH4, NH3, and HCN; the isothermal temperature; the gray-cloud opacity; and the reference pressure at the transit radius.

The model comparison test marginally prefer the gray-opacity case ($\chi^{2}_{\rm red}=1.21$) over a cloud-free case ($\chi^{2}_{\rm red}=1.46$), with a Bayes factor of 2.2. Adopting the gray-opacity model, we detect H2O at $4.0\sigma$ (Bayes factor = 601), HCN at $3.1\sigma$ (Bayes factor = 27.6), and ’nitrogen chemistry’ (combination of HCN and NH3) at $3.3\sigma$ (Bayes factor = 53.7). We do not detect CH4.

3.2.4 Waldmann

We retrieved the HST /WFC3 spectrum of WASP-63b using the Tau-REx atmospheric retrieval framework (Waldmann et al., 2015b, a; Waldmann, 2016). Based on the Tau code transmission forward models by Hollis et al. (2013), Tau-REx employs Nested Sampling (Feroz & Hobson, 2008) to solve the full Bayesian argument. Tau-REx can use high-resolution absorption cross-section or correlated-k tables as opacity inputs. Here we used the latter but find both to yield equivalent results for the wavelength ranges and sensitivities of the data at hand. We include pressure-dependent line broadening where such information is available, taking into account the J quantum number dependence on pressure broadening coefficients. In this study, we assume an isothermal atmospheric temperature-pressure profile, spanning from to $10^{-9}$ to 10 bar.

The Tau-REx code considers molecular opacities for H2O, CH4, NH3, HCN, TiO, and VO from Exomol (Tennyson & Yurchenko, 2012),CO and CO2 from HITEMP (Rothman et al., 2010), and O2, O3, NO2, NO, HCOOH, C2H6, and C2H2 from HITRAN (Rothman et al., 2009, 2013); CIA from H2-H2 (Borysow et al., 2001; Borysow, 2002) and H2-He (Abel et al., 2012); Rayleigh scattering (Lecavelier Des Etangs et al., 2008); and clouds using a hybrid model of gray-cloud opacities and a phenomenological Mie scattering (Lee et al., 2013).

We run two types of retrievals, a ’free’ retrieval with abundances of H2O, CH4, CO, CO2, NH3, as well as a chemical-equilibrium retrieval using an implementation of the ACE model by Agúndez et al. (2012), with the C/O ratio and atmospheric metallicity as free parameters. The rest of the retrieval parameters are the isothermal temperature; the top pressure of the gray cloud model and the particle size, composition and mixing ratio of the Mie-cloud model; and the planet reference radius at 10 bar.

We detect water with a $3.5\sigma$ significance (Bayes factor = 103) compared with a family of pure-cloud or featureless atmosphere models. We obtain $\log({\rm H}_{2}{\rm O})=-4.84^{+1.04}_{-1.53}$. We do not find any evidence of an extended Rayleigh curve due to hazes but found a gray-cloud model to be sufficient. We constrain the cloud-top pressure to $\log(p)=3.08^{+1.48}_{-0.93}$ Pa. The chemically consistent retrieval yielded two results. The first result yields a high metallicity atmosphere at 370 times solar. The second result yields a low metallicity atmosphere at 0.24 times solar. Both solutions feature comparable log-evidences and result in upper bounds of C/O at 0.49. A ratio of ${\rm C/O}<0.7$ is expected as only water is retrieved in this data set and is therefore consistent with the ’free’ retrieval approach above. The atmospheric metallicity is poorly constrained due to the presence of clouds which has the effect of muting the water feature and biasing the chemical-consistent model to either compensate with unrealistically high mean molecular weight atmospheres or unrealistically low trace gas abundances.

4. Discussion

The individual atmospheric analyses of WASP-63b agree that there is a robust water detection (3.5–4.0 $\sigma$), but with a muted absorption feature when compared to a clear solar-composition atmospheric model. It is unclear if the reason for the muting of the feature is the result of sub-solar water abundance, absorbing cloud opacity, or a high mean molecular mass. Thermochemically-consistent retrievals show a multimodal solution due to degeneracies between cloud opacity and composition (Figures: 7 and 9). Retrievals with simple temperature (isothermal) and cloud (gray opacity) models both produced fits consistent with retrievals with more complex models and thus the data does not warrant the incorporation of more complex models nor does it allow further constraints on cloud properties.

The high transit-depth values between 1.5 and 1.6 $\mu$m motivate the inclusion of HCN and the exploration of disequilibrium chemistry. Each retrieval team performed an additional retrieval exercise with a common set of assumptions to further explore the the inclusion of HCN as a means to fit the ‘bump’ in the spectrum between 1.5 and 1.6 $\mu$m. We implemented an isothermal temperature model, a gray-opacity cloud model, a free pressure–radius reference point, and opacities from H2-Rayleigh, H2-H2 and H2-He CIA, and H2O, CH4, HCN, and NH3. We adopted molecular abundances either from thermochemically-consistent calculations or from free abundances (constant vertical profiles), with the exception of the HCN abundance, which is always a free fitting variable (constant vertical profile).

We begin with this common set of assumptions and then compare retrieval results from teams with differing retrieval frameworks. Notable differences include: the statistical sampler from Blecic/Cubillos (MCMC) differ from the rest (Nested sampling), the molecular opacity handling from Waldmann and Line (correlated-K) differ from the others (cross-section sampling), and the chemistry from Line (thermochemically-consistent) differs from the rest (free constant vertical profiles). Figure 4 shows the retrieved spectrum from the run using the common assumptions. All four retrievals produced consistent spectral fits, seen in the intersecting 68% confidence regions around the best-fitting models. In terms of the atmospheric characterization, these retrievals confirm the previously found water detection.

The Bayesian hypothesis testing favors the fit with HCN, improving the fit at 1.5–1.6 $\mu$m. However, the detection significance is low and inconsistent, $3.1\sigma$ (Blecic/Cubillos), $2.1\sigma$ (Line), $3.1\sigma$ (MacDonald/Madhusudhan), and $1.9\sigma$ (Waldmann). Therefore, for the currently available data, the inclusion of HCN is not statistically justified within this model parameterization. Furthermore, to reproduce the observed values requires the HCN mole fraction be $\gtrsim 10^{-5}$; much higher than thermochemical-equilibrium values at the pressures probed by the WFC3 observations (Figure 5). To produce such high HCN abundances, one would need to invoke disequilibrium-chemistry processes by either quenching or photochemistry. Quenching can occur when higher temperatures at deep layers, below the levels probed by this observation, enhance the HCN abundance without needing the high metallicities from Figure 5. If vertical mixing dominates the mid-altitude abundances of the WASP-63b atmosphere (expected at the retrieved temperatures of $\sim$1000 K), HCN could be effectively quenched, maintaining the high abundances from the deep layers throughout the probed region. Similar deviations from equilibrium chemistry have been modeled for other Jupiter-like exoplanets (Moses et al., 2011; Venot et al., 2012). None of the retrievals constrain any of the other molecular abundances that could provide additional evidence for quenching (e.g. CO, CH4, NH3). However, photochemistry could play a role in removing these other molecules from the atmosphere while enhancing the mole fraction of HCN at pressures less than a millibar. Moses et al. (2011, 2013) show that ammonia and methane can be photochemically converted to HCN at the pressure levels probed by near–IR transmission spectroscopy thus driving the retrieved abundances much higher than equilibrium values. Future observations with extended wavelength ranges and higher sensitivity, such as JWST , can help to definitively confirm or rule out the detection of HCN, and other atmospheric species, thus constraining the presence of disequilibrium chemistry.

5. Conclusions

We present the observations of one transit of the hot Jupiter WASP-63b. Observations were conducted in the near-infrared using HST WFC3 G141. This study was done as a preliminary evaluation of the suitability of WASP-63b as one of the community targets for JWST ERS science. We have detected a muted water absorption feature at $\sim$ 1.4 $\mu$m confirming WASP-63b as a potential target for ERS science. The potential presence of an absorption feature at 1.55$\mu$m is not evidence enough to make strong conclusions about the presence of other molecules in the atmosphere, however, further observations by JWST would be able to identify additional spectral features that would allow us to further constrain the atmospheric composition. The observational window for observing WASP-63b with JWST is from September 23 – April 5. JWST is currently scheduled to launch in October 2018 and ERS observations would commence in April 2019. Assuming the mission remains on schedule, WASP-63b would not be observable until several months after the ERS program window. However, if there are any delays to launch or the start of ERS observations, WASP-63b would be a prime candidate for study with multiple instruments and modes.

6. Acknowledgments

This work is based on observations made with the NASA/ESA Hubble Space Telescope, obtained from the Data Archive at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. These observations are associated with program GO-14642. BMK acknowledges funding by HST-GO-14642.047 provided by the Space Telescope Science Institute. H.R. Wakeford acknowledges support from the NASA Postdoctoral Program, administered by USRA through a contract with NASA. Jasmina Blecic is supported by NASA trough the NASA ROSES-2016/Exoplanets Research Program, grant NNX17AC03G.

The authors acknowledge the contributions and support from members of the transiting exoplanet community who have contributed to and/or supported Program GO-14642 including: E. Agol, D. Angerhausen, T. Barman, J. Barstow, N. M. Batalha, S. Birkman, D. Charbonneau, N. Cowan, N. Crouzet, S. Curry, J. M. Desert, D. Dragomir, J. Fortney, A. Garcia Munoz, N. Gibson, J. Gizis, T. Greene, J. Harrington, T. Kataria, E. Kempton, H. Knutson, L. Kreidberg, M. Lopez-Morales, M. Rocchetto, E. Schlawin, E. Shkolnik, A. Shporer, D. Sing, K. Todorov, and J. de Wit.

Appendix A Blecic & Cubillos

Appendix B Line

Appendix C MacDonald & Madhusudhan

Appendix D Waldmann