The HST PanCET Program: Hints of Na i and Evidence of a Cloudy Atmosphere for the Inflated Hot Jupiter WASP-52b

We obtained time series spectroscopy during three transits of WASP-52b with HST/STIS on UT 2016 November 1 and UT 2016 November 29 using the G430L grating (2892–5700 Å) and UT 2017 May 11 using the G750L grating (5240–10270 Å). Two additional transits were observed with Spitzer/IRAC as part of GO program 13038 (PI Stevenson) in the 3.6 μm channel on UT 2016 October 18 and in the 4.5 μm channel on UT 2018 March 22. Table 1 summarizes the transit observations and instrument settings for each visit.

2.1.1. HST/STIS

2.1.2. Spitzer/IRAC

2.2.1. HST/STIS

2.2.2. Spitzer/IRAC

We acquired 135 out-of-transit R-band images over the 2014–2015, 2015–2016, and 2016–2017 observing seasons with Tennessee State University's 14 inch Celestron Automatic Imaging Telescope. These observations, however, do not include the epochs of the Spitzer and ground-based (Chen et al. 2017; Louden et al. 2017) transit observations. We therefore used 740 out-of-transit V-band images from Ohio State University's All-sky Automated Survey for Supernovae¹⁶ program (Shappee et al. 2014; Kochanek et al. 2017) for more detailed activity monitoring coverage. The ASAS-SN observations were taken over the 2013–2014, 2014–2015, 2015–2016, and 2016–2017 observing seasons. The number of observations, date range, mean magnitude, and standard deviation of the data taken in each observing season are included in Table 2.

To quantify the level of activity in WASP-52, we used observations from the Advanced CCD Imaging Spectrometer (ACIS) on the Chandra X-ray telescope. These data were taken from the Chandra public archive (GO 15728; PI Wolk). Standard CIAO¹⁷ tasks were applied to reduce the data. The software ISIS (Houck & Denicola 2000) was used to fit the spectrum with a metallicity fixed to photospheric values ([Fe/H] = 0.03), and the interstellar matter (ISM) absorption was assumed to be N_H = 4 × 20 cm⁻³ (based on target location and distance), consistent with the fit. The resulting 1 − T fit has $\mathrm{log}\,{\text{}}T({\rm{K}})={6.54}_{-0.29}^{+1.06}$ and log EM (cm⁻³) = ${51.14}_{-0.63}^{+0.33}$. We measured an X-ray luminosity of L_x = (3.5 ± 1.0) × 10²⁸ erg s⁻¹ in the 0.12–2.48 keV band (J. Sanz-Forcada et al. 2018, in preparation). The light curve shows variability, but poor statistics hamper a detailed analysis. The corresponding log(L_x/L_bol) = −4.7 indicates that WASP-52 is a moderately active star. Further details will be included in J. Sanz-Forcada et al. (2018, in preparation).

We initially adopted a simple sinusoidal model of period P = 11.8 ± 3.3 days based on the Hébrard et al. (2013) rotational period, but we found that this approach did not accurately model the variations in amplitude and period of the AIT and ASAS-SN variability monitoring data over all observing seasons. Discrepancies between the sinusoidal model and the photometric monitoring data are likely due to different spot configurations at different observed epochs (Dang et al. 2018), suggesting that a quasi-periodic model would more accurately fit the data.

We therefore jointly modeled the AIT and ASAS-SN ground-based stellar activity data using a GP regression, a framework that has been shown to accurately disentangle stellar activity signals from planetary signals in radial velocity data (e.g., Aigrain et al. 2012; Haywood et al. 2014) and photometry (e.g., Pont et al. 2013; Aigrain et al. 2016; Angus et al. 2018). We ran a GP optimization routine (Ambikasaran et al. 2015) using the George package for Python. We used a three-component kernel to model the quasi-periodicity of the ground-based activity monitoring data, irregularities in the amplitude of the ground-based photometry, and stellar noise. See Appendix A.1 for further details regarding the functional form of our chosen kernel. We use a gradient-based optimization routine to find the best-fit hyperparameters and set the 11.8 ± 3.3 day rotation period from Hébrard et al. (2013) as a uniform prior with an uncertainty three times larger than the literature value. Figure 1 shows the ground-based variability monitoring data overplotted with the GP model for all seasons modeled jointly and for each observing season separately.

Zoom In Zoom Out Reset image size

The GP regression model with a quasi-periodic kernel reproduces the observed flux variations well for epochs during which we have stellar activity data, but has large uncertainties outside of a given season. The model accurately predicts the stellar activity behavior for each observing season when fit separately, but this model prediction has a significantly larger 1σ uncertainty when fitting the data from all seasons together. We therefore modeled the data for each observing season separately, and we used the amplitude of the photometric variation given by the GP model for each epoch to correct the transmission spectrum for the effects of stellar activity (see Section 3.4).

The temperature difference between the stellar photosphere and the spotted region introduces a slope in the planet spectrum (e.g., Kreidberg 2017), so it is necessary to correct for this effect. Using the results of the GP regression model described in Section 3.3, we followed the prescription of Huitson et al. (2013) to correct the spectrophotometric light curves for the effect of unocculted stellar spots of a fixed size. This method involves estimating the variability amplitude for a broadband wavelength value (determined by the photometric filter of the activity monitoring data), which can then be used as an anchor for the wavelength-dependent flux correction on the HST and Spitzer data.

We first converted the ASAS-SN and AIT photometric variability monitoring data to relative flux. Excluding in-transit measurements, we estimated the nonspotted stellar flux to be F_⋆ = max(F) + kσ, where F is the variability monitoring data, σ is the dispersion of the photometric measurements, and k is a factor fixed to unity (Aigrain et al. 2012; see also Appendix A.2). The variability monitoring data were then normalized to the nonspotted stellar flux to estimate the amount of dimming. Table 3 gives the dimming values for each transit observation. We also computed the amplitude of the spot corrections at the variability monitoring wavelength Δf₀ = 1 − f_norm, where f_norm is the mean of the normalized flux array.

To derive the wavelength-dependent flux correction, we used a stellar flux model (T_eff = 5000 K, log(Z) = −1.5, log(g) = 4.5) and a spot model (T_eff = 4750 K, log(Z) = −1.5, log(g) = 4.5) from the Kurucz (1993) 1D ATLAS grid¹⁸ of stellar atmospheric models. The stellar flux model was chosen based on the effective temperature, metallicity, and surface gravity of WASP-52 given in Hébrard et al. (2013). The spot model is the same as the stellar model but 250 K cooler (Berdyugina 2005). To select this spot model, we tested different cold spots ranging from 3500 to 4750 K (corresponding to temperature differences between the spot and stellar photosphere from 1500 to 250 K). Figure 2 shows that spots at colder temperatures exhibit less flux dimming at longer wavelengths and a weaker slope in the optical. To correct for the effect of unocculted spots, we therefore use a spot model at 4750 K for which cold spots give the strongest slope to account for the maximum possible contribution from starspots in the data.

Zoom In Zoom Out Reset image size

We interpolated the stellar model to the spot model grid and computed the wavelength-dependent correction factor derived in Sing et al. (2011b):

$$\begin{eqnarray}&&f(\lambda ,T)=\left(1-\displaystyle \frac{{F}_{\lambda ,{T}_{\mathrm{spot}}}}{{F}_{\lambda ,{T}_{\mathrm{star}}}}\right)/\left(1-\displaystyle \frac{{F}_{{\lambda }_{o},{T}_{\mathrm{spot}}}}{{F}_{{\lambda }_{o},{T}_{\mathrm{star}}}}\right),\end{eqnarray} \tag{ 1 }$$

where ${F}_{\lambda ,{T}_{\mathrm{spot}}}$ is the stellar model flux at temperature T_spot and at the wavelength of the transit observations, ${F}_{\lambda ,{T}_{\mathrm{star}}}$ is the stellar model flux at the wavelength of the transit observations and at temperature T_star, ${F}_{{\lambda }_{o},{T}_{\mathrm{spot}}}$ is the stellar model flux at temperature T_spot at the reference wavelength of the activity monitoring data, and ${F}_{{\lambda }_{o},{T}_{\mathrm{star}}}$ is the stellar model flux at temperature T_star at the activity monitoring reference wavelength. The final flux dimming correction was then calculated as Δf = Δf₀ × f(λ, T).

We computed the average flux correction over each bandpass (see Table 3) and applied the correction to each spectrophotometric light curve using

$$\begin{eqnarray}&&{y}_{\mathrm{corrected}}=y+\displaystyle \frac{{\rm{\Delta }}f}{(1-{\rm{\Delta }}f)}\overline{{y}_{\mathrm{oot}}},\end{eqnarray} \tag{ 2 }$$

where y_corrected is the corrected light curve flux, y is the original (uncorrected) light curve, and $\overline{{y}_{\mathrm{oot}}}$ is the mean of the out-of-transit exposures. We then fit analytic transit light curve models (Mandel & Agol 2002) to the stellar activity-corrected light curves, as detailed in Section 4.

We produced the STIS broadband (wavelength-integrated) light curves by summing the time series over the complete wavelength range of the bandpass (2892–5700 Å for the G430L grating; 5240–10270 Å for the G750L grating). The white light curves for each of the STIS visits are shown in Figure 3. Photometric uncertainties were derived based on pure photon statistics. The raw white light curves exhibited instrumental systematics related to the orbital motion of the spacecraft (Gilliland et al. 1999; Brown 2001). In particular, the HST focus is known to experience significant variations on the spacecraft orbital timescale resulting from thermal expansion/contraction during the spacecraft's day/night orbital cycle.

Zoom In Zoom Out Reset image size

As in past studies (e.g., Huitson et al. 2013; Sing et al. 2013; Nikolov et al. 2014), we accounted for and detrended these instrumental systematic effects by fitting a fourth-order polynomial to the flux dependence on HST orbital phase. We excluded the first orbit and the first exposure of each subsequent orbit in accordance with common practice, because these data have unique, complex systematics (see Figure 3) and were taken while the telescope was thermally relaxing into its new pointing position. We applied orbit-to-orbit flux corrections to the STIS data by fitting a polynomial of the spacecraft orbital phase (ϕ_t), drift of the spectra on the detector (x and y), the shift of each stellar spectrum cross-correlated with the first spectrum of the time series (ω), and time (t).

We then generated systematics models spanning all possible combinations of detrending variables (see Appendix B.1 and Table 7 for details) and performed separate fits using each systematics model included in the two-component function. For each model, we fixed P, i, and a/R_⋆ to the values given in Hébrard et al. (2013), assumed zero eccentricity, and fit for T₀, R_p/R_⋆, stellar baseline flux, and instrument systematic trends. We determined the best-fit parameters of the two-component model using a Levenberg–Marquardt least-squares fitting routine (Markwardt 2009) to derive the system parameters and calculated the Akaike Information Criterion (AIC; Akaike 1974) for each function. The results of these fits are shown in Table 4. The STIS stellar spectra and raw and detrended white light curves for each visit are shown in Figure 3.

We marginalized over the entire set of functions following the framework outlined in Gibson (2014). Marginalization over multiple systematics models assumes equally weighted priors for each model tested. Table 8 shows the results for each systematics model using the STIS white light curves. The reduced chi-squared χ_r from the fits can be used as a proxy for the photon noise level of a given data set (Nikolov et al. 2014). We selected the systematics model for use in detrending the STIS white light curves based on the lowest AIC value because the light curves show no significant correlation to the quantified systematic parameters (Nikolov et al. 2014).

We produced the spectrophotometric light curves by dividing the spectra into bins of width 17–400 Å and integrating the flux from each bandpass. The width of the bins is determined primarily by the need to achieve a given photometric precision to be sensitive to features in the planet's transmission spectrum that are comparable in amplitude to an atmospheric scale height. The smaller bin sizes were chosen to be centered at specific absorption features, such as Na i at 5893 Å and K i at 7665 Å.

We modeled the systematic errors using two methods. In the first approach, we independently fit each of the binned light curves with the same family of transit+systematics models (see Appendix B.1) as the broadband light curve. The only differences are that we fixed the midtransit time T₀ and the scaled semimajor axis a/R_⋆ to the white light curve best-fit values. We also fixed the limb-darkening coefficients (computed following the procedure outlined in Section 4.4) to the derived theoretical values. In the second approach, we performed a common mode correction to remove color-independent systematic trends from each spectral bin. Then, we fit the common-mode-corrected spectroscopic light curves by fitting for residuals with a parameterized model of six fewer free parameters (c₁ − c₄, T₀, and a/R_⋆) and marginalizing over the entire set of functions defined in Appendix B.1. The common mode trends are computed by dividing the raw flux of the white light curve in each grating by the best-fit analytic transit model. We applied the common mode correction by dividing each binned light curve by the derived common-mode flux. Removing common mode trends is known to reduce the amplitude of the observed HST breathing systematics, since these trends are similar in wavelength across the detector (Sing et al. 2013; Nikolov et al. 2014). The common-mode-corrected light curves are shown in Figure 3.

Both methods produced similar results (i.e., consistent baseline in R_p/R_⋆). Since the common mode correction produces lower dispersion in the spectrophotometric light curves and smaller R_p/R_⋆ uncertainties (Nikolov et al. 2015), we report the common-mode-corrected results for the fitted R_p/R_⋆ (before and after applying the correction for unocculted spots discussed in Section 3.4) and nonlinear limb-darkening coefficients for each spectroscopic channel in Table 5. The raw and detrended STIS spectrophotometric light curves for the G430L and G750L gratings are shown in Figures 4–6.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We modeled the 3.6 and 4.5 μm IRAC transit photometry in accordance with the methods outlined in Sing et al. (2011b, 2013) and Nikolov et al. (2014). To correct for flux variations from intrapixel sensitivity, we fit a polynomial to the stellar centroid position (Charbonneau et al. 2005, 2008; Reach et al. 2005; Knutson et al. 2008). This technique is effective on short timescales (<10 hr) and for small (<0.2 pixels) variations in the stellar centroid position (Lewis et al. 2013). We corrected for systematic effects using a model given by the linear combination

$$\begin{eqnarray}&&f(t)={a}_{0}\,+{a}_{1}x\,+{a}_{2}{x}^{2}\,+{a}_{3}y\,+{a}_{4}{y}^{2}\,+{a}_{5}{xy}\,+{a}_{6}t,\end{eqnarray} \tag{ 3 }$$

where f(t) is the stellar flux as a function of time, the coefficients a₀ through a₆ are free fitting parameters, x and y are the detector positions of the stellar centroid, and t is time. We generate all possible model combinations of Equation (3), which we marginalize over using the Gibson (2014) procedure as detailed in Section 4.1 and Appendix B.1. Using the WASP-52 system parameters from Hébrard et al. (2013) as priors, we jointly fit for all parameters and find that our results (see Table 4) agree within 1σ with the STIS white light curve analysis. To measure R_p/R_⋆ for the transmission spectrum, we fixed the orbital period P, normalized planet semimajor axis a/R_⋆, inclination i, and central transit time T₀ to the best-fit values from the joint fit. The limb-darkening coefficients were also fixed to their theoretical values based on 3D stellar atmosphere models (see Section 4.4). The measured planetary radius and limb-darkening coefficients are included in Table 5. Figures 7 and 8 show the raw and detrended Spitzer transit light curves.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We modeled the limb darkening of WASP-52b using the four-parameter nonlinear limb-darkening law (Claret 2000) given by

$$\begin{eqnarray}&&\displaystyle \frac{I(\mu )}{I(1)}=1\,-\,\displaystyle \sum _{n=1}^{4}{c}_{n}(1-{\mu }^{n/2}),\end{eqnarray} \tag{ 4 }$$

where I(1) is the intensity at the center of the stellar disk, c_n (n = 1–4) are the limb-darkening coefficients, and μ = cos(θ), where θ is the angle between the normal to the stellar surface and the line of sight.

To derive the stellar limb-darkening coefficients, we followed the procedure described in Sing (2010) and initially used values for the four limb-darkening coefficients based on 1D ATLAS theoretical stellar models (Kurucz 1993). We then derived the limb-darkening coefficients from 3D stellar models (Magic et al. 2015) and compared to the 1D results to eliminate the known wavelength-dependent degeneracy of limb darkening with transit depth (Sing et al. 2008a). This approach reduces the number of free parameters in the fit (typically four parameters per grating), but may cause an underestimation of errors in the derived spectrum. The derived nonlinear 3D limb-darkening coefficients for each spectrophotometric light curve are shown in Table 5.

We inspect the presence and significance of the Na i and K i features in the WASP-52b transmission spectrum using a grid of spectrophotometric channels ranging in width from 30 to 255 Å in steps of 15 Å, centered on the Na i (5893 Å) and K i (7665 Å) resonance doublets. The minimum bin size (30 Å) is defined to include both doublet lines. If a planetary atmospheric signal is present, this binning scheme should demonstrate a gradual decay in the measured transit depth for larger bin sizes. The results of this analysis are shown in Figure 11.

Zoom In Zoom Out Reset image size

For the Na i feature, we note a gradual decrease in the measured transit depth for bins 30–100 Å wide. For wider bins, the signal is largely washed out and the transit depth remains unchanged within the uncertainties. This trend is expected when observing a narrow absorption peak and is consistent with the presence of a cloud deck. By measuring the difference in transit depth between the narrowest (30 Å) bin and the flat transit depth baseline, we detect the core of the Na i doublet at 2.3σ confidence. In the case of K i, we do not see evidence of absorption even in the narrowest spectroscopic channel, suggesting that this feature is either masked by a thick cloud deck or not as abundant in the atmosphere of WASP-52b. For further discussion, see Section 6.1.

We fit the combined STIS+Spitzer transmission spectrum to the publicly available grid of forward model transmission spectra (Goyal et al. 2018) produced using the ATMO 1D radiative–convective equilibrium model (Amundsen et al. 2014; Tremblin et al. 2015, 2016; Drummond et al. 2016). The models are generated for the parameters (e.g., mass, radius, gravity) of WASP-52b.

The grid¹⁹ includes 3920 model transmission spectra of WASP-52b for five temperatures (1015, 1165, 1315, 1465, 1615 K), seven metallicities (0.005, 0.1, 1, 10, 50, 100, 200× solar), seven C/O ratios (0.15, 0.35, 0.56, 0.70, 0.75, 1.0, and 1.5), four values of the haziness parameter α_haze (1, 10, 150, and 1100), and four values of the cloudiness parameter α_cloud (0, 0.06, 0.2, and 1). The parameter α_haze is a proxy for the haze enhancement factor of small scattering aerosol particles suspended in the atmosphere, where α_haze = 1 indicates no haze and α_haze = 1100 indicates thick hazes. The cloudiness parameter α_cloud gives the strength of gray scattering due to H₂ at 350 nm, with α_cloud = 0 corresponding to no clouds and α_cloud = 1 corresponding to a thick cloud deck. See Goyal et al. (2018) and references therein for further details. The transmission spectra are computed assuming isothermal pressure−temperature (P − T) profiles and condensation without rainout (local condensation).

We computed the mean model prediction for the wavelength range of each spectroscopic channel (see Table 5) and performed a least-squares fitting of the band-averaged model to the spectrum. For the fitting procedure, we allowed the vertical offset in R_p/R_⋆ between the spectrum and model to vary while holding all other parameters fixed in order to preserve the model shape. The number of degrees of freedom for each model is n − m, where n is the number of data points and m is the number of fitted parameters. Since n = 38 and m = 1, the number of degrees of freedom for each model is constant. From the fits, we computed the χ² statistic to quantify our model selection.

Figure 9 shows the best-fit model, representative clear and hazy models, and a flat model compared to the observed transmission spectrum. The best-fitting model (χ² = 39.3) is cloudy (α_cloud = 1.0) and slightly hazy (α_haze = 10) with a 2.3σ signature of Na i absorption, a temperature of T = 1315 K, solar metallicity ([M/H] = 0.0), and slightly supersolar C/O (C/O = 0.70). The selected clear model (χ² = 49.5) has a lower temperature (T = 1015 K) and no clouds (α_cloud = 0.00) or hazes (α_haze = 1). The representative hazy model (χ² = 43.3) is similar to the clear model, but with extreme haziness (α_haze = 1100). The flat model (χ² = 52.2) represents a featureless (gray) spectrum.

The χ² contour plot for the model grid fits is shown in Figure 12. The grid provides constraints on the atmospheric and physical parameters of WASP-52b, with the 1σ confidence region favoring high cloudiness, slight haziness, solar metallicity, slightly supersolar C/O ratio (>0.56), and an equilibrium temperature of 1315 K.

Zoom In Zoom Out Reset image size

In addition to the combined STIS+Spitzer transmission spectrum we report here, there are two ground-based optical transmission spectrum measurements for WASP-52b. Most recently, Louden et al. (2017) used spectroscopy between 4000 and 8750 Å for two transit observations from the ACAM instrument mounted on the William Herschel Telescope (WHT). Chen et al. (2017) observed one transit with the Gran Telescopio Canarias's (GTC) OSIRIS instrument in the 5220–9030 Å wavelength range. We show these transmission spectra compared to our STIS+Spitzer results in Figure 9.

With the WHT/ACAM observations, Louden et al. (2017) modeled spot-crossing events via GPs, adopted a harmonic analysis of ground-based photometric monitoring, and found varying levels of activity over time with evidence of differential rotation. These results reveal a flat transmission spectrum attributed to an optically thick cloud deck. The GTC/OSIRIS observations indicate a cloudy atmosphere with a ∼3σ detection of Na i and a weaker detection of K i. Calculations of the integrated absorption depth for the Na i and K i signals suggest an inverted temperature structure for the upper atmosphere of WASP-52b (Chen et al. 2017).

Our R_p/R_⋆ baseline is consistent within 1σ with the ground-based transmission spectrum from WHT/ACAM (Louden et al. 2017). The most significant difference is that the WHT spectrum does not show any variation in R_p/R_⋆ around 5893 Å. If there is a weak signal of Na i absorption from the planet, the resolution of the spectrum, composed of equally sized spectrophotometric bins of width 250 Å, may be washing it out, as illustrated in Figure 11.

The baseline of our spectrum matches less well with the ground-based GTC/OSIRIS transmission spectrum (Chen et al. 2017). We note a constant offset in the absolute measured transit depths of our STIS+Spitzer spectrum and the GTC measurement, with a difference in the R_p/R_⋆ baseline of ∼3σ. The authors attribute their shallower transit depth measurement compared to previous studies (e.g., Hébrard et al. 2013; Kirk et al. 2016; Mancini et al. 2017) to the effects of stellar activity. We find evidence of a Na i signal that is consistent with the GTC detection within 1σ, but our spectrum shows no evidence of K i absorption, contrary to the Chen et al. (2017) result. This discrepancy could be due to the different methods used to correct for the effects of stellar activity (see Section 3 and see Chen et al. 2017).

We compared WASP-52b to the well-studied inflated hot Jupiter HAT-P-1b (see Nikolov et al. 2014) because both planets have comparable system parameters and atmospheric properties. HAT-P-1b and WASP-52b have overlapping surface gravity, equilibrium temperature, mass, radius, and stellar irradiation, and the transmission spectra of both planets are marginally flat with evidence of Na i absorption but no observable K i absorption. Table 6 compares the stellar and planetary parameters for WASP-52 (Hébrard et al. 2013; Table 4) and HAT-P-1 (Nikolov et al. 2014).

Table 6.  System Parameters and Best-fitting Model Parameters for WASP-52b and HAT-P-1b (Nikolov et al. 2014)

	WASP-52b	HAT-P-1b
Spectral type	K2V	G0V
Stellar mass M⋆  (M⊙)	0.87 ± 0.03	1.15 ± 0.05
Stellar radius R⋆  (R⊙)	0.79 ± 0.02	1.17 ± 0.03
Surface gravity log(g⋆) (cgs)	4.58 ± 0.014	4.36 ± 0.01
Stellar temperature Teff  (K)	5000 ± 100	5980 ± 50
Metallicity [Fe/H]	0.03 ± 0.12	0.13 ± 0.01
Stellar irradiation I  (erg cm−2 s−1)	6.5 ± 0.4 × 108	7.0 ± 0.4 × 108
$\mathrm{log}({R}_{{HK}}^{{\prime} })$	−4.4 ± 0.2	−4.98 ± 0.1
Mp (MJ)	0.434 ± 0.024	0.525 ± 0.019
Rp (RJ)	1.253 ± 0.027	1.319 ± 0.019
ρ (ρJ)	0.206 ± 0.009	0.213 ± 0.010
Best-fitting ATMO models
Teq (K)	1315	1322
Fe/H	0.0	1.0
C/O	0.70	0.15
αcloud	1.0	0.20
αhaze	10	10

Download table as:  ASCIITypeset image

HAT-P-1b has a precise transmission spectrum measurement from HST/STIS (Nikolov et al. 2014), which we use to compare the atmospheric properties of both planets. For this comparison, we reconstructed the HST/STIS transmission spectrum for WASP-52b using the same binning scheme of Nikolov et al. (2014) for HAT-P-1b and fit the light curves for these bins based on the methods outlined in Section 4. Figure 13 shows the HST/STIS transmission spectra of both planets with identical binning. Based on this comparison, the spectra of both planets are identical within the uncertainties with an average 1σ difference.

Zoom In Zoom Out Reset image size

As reported in Nikolov et al. (2014), the best-fit model for HAT-P-1b is a hazy spectrum with Na i absorption and an extra optical absorber to account for the observed absorption enhancement at wavelengths longer than ∼0.85 μm. To directly compare this interpretation to the WASP-52b results reported here, we fit the HAT-P-1b spectrum (Nikolov et al. 2014) to the open source ATMO grid of forward models generated for the parameters of HAT-P-1b (see Appendix C and Figure 16 for details). The best-fit model parameters for HAT-P-1b and WASP-52b are shown in Table 6.

Since WASP-52b and HAT-P-1b have similar optical transmission spectra (∼0.29–1 μm), we compare the atmospheric properties of these planets in the near-infrared to ascertain if these planets would be good candidates for comparative atmospheric analyses with JWST. Figure 14 shows the best-fit models for both planets from 0.29 to 5 μm. Beyond 1 μm, we find that the best-fitting models for HAT-P-1b and WASP-52b do not agree with each other as they do in the optical. For further discussion regarding potential reasons for this near-infrared discrepancy, see Section 6.2.

Zoom In Zoom Out Reset image size

As reported in Section 5.1, we find hints of Na i absorption but no evidence of K i in the transmission spectrum of WASP-52b. A similar trend has been observed for HD 189733b (Pont et al. 2013), HAT-P-1b (Nikolov et al. 2014), and WASP-17b (Sing et al. 2016). The reverse trend (i.e., the presence of a K i signal but no Na i) has been observed in several planets, including WASP-31b (Sing et al. 2015) and HAT-P-12b (Sing et al. 2016). We note, however, that the majority of current sodium and potassium detections in exoplanet atmospheres are low significance, and a nondetection of K i can only be interpreted as an upper limit to the abundance of that element in the atmospheric layers probed by our observations, considering their uncertainties. Based on the available data and their uncertainties, we estimate an abundance ratio of $\mathrm{ln}[\mathrm{Na}/{\rm{K}}]={8.32}_{-6.03}^{+6.51}$. This value is consistent with solar abundances but has large uncertainties and is not well constrained.

These results are consistent with our forward model analysis, which favors the presence of clouds and slight hazes that mute spectroscopic features in the planet's atmosphere. The best-fitting models give temperatures in the range T ∼ 1000–1300 K. Compared to the planet's equilibrium temperature (T_eq = 1315 K; Hébrard et al. 2013), these temperature estimates may be driven by the gradient of the Rayleigh scattering slope. Significant sodium condensation in the form of Na₂S is expected at lower temperatures (∼1000 K) and potassium condensation in the form of KCl at even lower temperatures (∼600 K; Marley et al. 2013). Tenuous Na₂S and KCl clouds could form at shallow pressures in WASP-52b's atmosphere, since the equilibrium temperature does not represent the full range of temperatures in a planet's atmosphere. The cloud deck in WASP-52b's atmosphere could therefore be composed of species other than sodium or potassium compounds, such as silicates.

Regardless of the composition of these clouds, they are likely masking the K i feature and the wings of the Na i resonance core. We do not resolve the broad wings of the Na i line (Figure 11), which may suggest the presence of an extra absorber or scatterer in the atmosphere that is obscuring or masking the atmospheric Na i and K i absorption features (Seager & Sasselov 2000; Nikolov et al. 2014).

If a K i signal is truly lacking in the transmission spectrum of WASP-52b, an alternative explanation could be attributed to an underabundance of this element in the planetary atmosphere. If the cloud deck is composed of sodium or potassium compounds, however, the gas-phase abundances of these species would not reflect primordial abundances. Since K i is a weaker spectroscopic feature, it may be present in the planet's atmosphere but not detectable with the precision of the STIS data. Higher resolution, higher precision observations are necessary to confirm this idea.

Although the atmospheres of planets studied thus far appear diverse (Sing et al. 2016), we have not yet been able to identify any clear correlations between planetary atmospheric properties and other system parameters. Comparing the transmission spectra of planets with similar system parameters may therefore prove insightful in searching for common atmospheric characteristics. WASP-52b is a good target for such comparisons because the well-studied inflated hot Jupiter HAT-P-1b has comparable system parameters and atmospheric properties.

We compare the observed transmission spectrum of WASP-52b presented here to that of the well-studied inflated hot Jupiter HAT-P-1b. These two planets have similar system parameters (see Table 6), although WASP-52 is ∼0.9 times more active than HAT-P-1b (Nikolov et al. 2014). Additionally, the optical transmission spectra of both planets are marginally flat with evidence of Na i but no observable evidence of K i. We compare their transmission spectra with identical binning schemes in Figure 13, and we find that the spectra of both planets in the optical (∼0.29–1 μm) are identical within the uncertainties with an average 1σ difference.

Extending this comparison to near-infrared wavelengths, however, reveals that their transmission spectra differ considerably beyond 1 μm. The best-fit ATMO model for WASP-52b shows muted H₂O and CH₄ spectral features compared to the best-fitting HAT-P-1b model, which could indicate that WASP-52b has a higher aerosol layer. Near-infrared HST/WFC3 observations of HAT-P-1b reveal H₂O absorption at >5σ confidence (Wakeford et al. 2013), and the deep H₂O feature shown in the best-fit ATMO model matches these data. The best-fit atmospheric model for WASP-52b shows a weaker (but still observable) H₂O feature at 1.4 μm, and evidence of H₂O absorption at 1.4 μm for WASP-52b has recently been shown in Bruno et al. (2018). The activity levels of the host stars may also contribute to the divergence of the transmission spectra for these two planets at near-infrared wavelengths, although observations of the stellar UV fluxes are necessary to confirm this hypothesis. Comparative near-infrared observations with JWST can confirm the atmospheric similarities of these planets at shallower atmospheric layers compared to those probed by STIS.

We used a GP regression to model the ground-based stellar activity monitoring data (Pont et al. 2013; Haywood 2015; Aigrain et al. 2016; Angus et al. 2018). In our GP analysis, we used a three-component kernel of the form K = k₁ + k₂ + k₃. The k₁ term models the flexible (quasi-)periodicity of the ground-based activity monitoring data. It is a squared exponential kernel multiplied by an exponential sine squared kernel (k₁ = A² [k(r²) × k(x_i, x_j)]), where

$$\begin{eqnarray}&&k({r}^{2})={e}^{{r}^{2}/2}\end{eqnarray} \tag{ 5 }$$

represents the squared exponential kernel for periodic variations in the time series data parameterized by r. The exponential sine squared kernel is given by

$$\begin{eqnarray}&&k({{\boldsymbol{x}}}_{i},{{\boldsymbol{x}}}_{j})=\exp \left(-{\rm{\Gamma }}{\sin }^{2}\left[\displaystyle \frac{\pi }{P}| {x}_{i}-{x}_{j}| \right]\right),\end{eqnarray} \tag{ 6 }$$

where Γ represents the scale of the correlations and P is the period of the oscillations.

The second term, k₂, represents the irregularities in amplitude and period. It is a rational quadratic kernel of the form

$$\begin{eqnarray}&&{k}_{2}({r}^{2})={A}^{2}\left(1-\displaystyle \frac{{r}^{2}}{2\alpha }\right)\end{eqnarray} \tag{ 7 }$$

for amplitude A, Gamma distribution parameter α, and length scale r.

The k₃ term incorporates stellar noise in the GP model and is a squared exponential kernel (of the form shown in Equation (5)) added to a white kernel of the form k₃(x_i, x_j) = cδ_ij for constant c and diagonal value δ_ij.

As described in Section 3, we account for the effects of stellar activity and unocculted starspots on the transmission spectrum by using ground-based photometry and by deriving the wavelength-dependent flux correction Δf(λ, T) (Sing et al. 2011b; see also Section 3.4, Equation (1)). This correction depends on the parameter k, which provides an assumption for the nonspotted flux on the stellar surface (Aigrain et al. 2012). We fixed k to unity based on Aigrain et al. (2012), who found this value appropriate to use for active stars. Physically, k = 1 corresponds to a spot contribution that the viewer never sees (or always sees) that is about the same as the contribution of spots that come into and out of view.

The value of this parameter is based on several assumptions; however, it is, in actuality, very ill constrained, and we warn the reader of the difficulty in selecting a value of k. To explore the potential effect of the choice of k on the final transmission spectrum, we tested different values of this parameter in the range k = 0–1 in steps of 0.2 (Figure 15). We find that larger values of k correspond (linearly) to a higher derived flux correction Δf(λ, T). When applying the activity correction to the binned light curves, we find that the light curves are rescaled to a higher R_p/R_⋆ baseline if a larger correction is applied. Thus, the parameter k is not wavelength-dependent and therefore does not affect the shape of the slope of the transmission spectrum. Rather, varying the value of k simply shifts the R_p/R_⋆ baseline of the transmission spectrum vertically.

Zoom In Zoom Out Reset image size

As described in Section 4.1, we detrended instrument systematic effects in the light curves by fitting a fourth-order polynomial to the flux dependence on HST orbital phase. Each model represents a unique linear combination of the detrending variables: orbital phase (ϕ_t), drift of the spectra on the detector (x and y), the shift (ω) of each stellar spectrum cross-correlated with the first spectrum of the time series, and time t. The variables x and y are the trace slope and an offset in the cross dispersion direction, respectively. The parameter ω is measured by cross-correlating a reference spectrum with the remaining spectra and is measured prior to resampling the spectra.

For both the G430L and G750L STIS observations, we used the 25 systematics models listed in Table 7 to detrend the light curves. After performing separate fits for each model, we marginalized over the entire set of systematics models assuming equally weighted priors to select which systematics model to use. Table 8 summarizes the selection of systematics models based on the STIS white light curves. We selected the model for detrending based on the lowest Aikake Information Criterion value.

Table 7.  White Light Curve Systematics Models

Model
G430L models
1	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t$
2	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{x}^{2}$
3	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{y}^{2}$
4	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+{x}^{2}+y$
5	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega$
6	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x$
7	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y$
8	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{3}+x$
9	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y$
10	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +x+y$
11	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}$
12	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+{\omega }^{3}$
13	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}+y$
14	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +x+{x}^{2}+{x}^{3}$
15	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y+{y}^{2}$
16	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y+{y}^{2}+{y}^{3}$
17	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+{\omega }^{3}+x$
18	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}$
19	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}+{x}^{3}$
20	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y+{y}^{2}$
21	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y+{y}^{2}+{y}^{3}$
22	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x$
23	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x+y$
24	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x+{x}^{2}+y+{y}^{2}$
25	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+{\omega }^{2}+{\omega }^{3}+x+{x}^{2}+{x}^{3}+y+{y}^{2}+{y}^{3}$
G750L models
1	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t$
2	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{x}^{2}$
3	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{y}^{2}$
4	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+{x}^{2}+y$
5	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega$
6	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x$
7	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y$
8	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +x$
9	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y$
10	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +x+y$
11	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}$
12	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega$
13	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}+y$
14	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +x+{x}^{2}+{x}^{3}$
15	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y+{y}^{2}$
16	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +y+{y}^{2}+{y}^{3}$
17	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x$
18	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}$
19	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+x+{x}^{2}+{x}^{3}$
20	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y+{y}^{2}$
21	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+y+{y}^{2}+{y}^{3}$
22	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x$
23	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x+y$
24	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x+{x}^{2}+y+{y}^{2}$
25	${\phi }_{t}+{\phi }_{t}^{2}+{\phi }_{t}^{3}+{\phi }_{t}^{4}+t+\omega +{\omega }^{2}+x+{x}^{3}+{y}^{2}+{y}^{3}$

Download table as:  ASCIITypeset image

The broadband STIS+Spitzer transmission spectrum for WASP-52b reported in Table 5 gives the weighted mean of the two STIS G430L observations (visits 52 and 53). In Tables 9 and 10, we report the spectrophotometric light curves for each G430L visit. To produce the raw transmission spectrum from the spot-corrected spectrum, the reader may simply reverse the correction given in Table 3 on each transit observation separately.