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

We present an optical to near-infrared transmission spectrum of the inflated hot Jupiter WASP-52b using three transit observations from the Space Telescope Imaging Spectrograph (STIS) mounted on the Hubble Space Telescope, combined with Spitzer/Infrared Array Camera (IRAC) photometry at 3.6 microns and 4.5 microns. Since WASP-52 is a moderately active (log(Lx/Lbol) = -4.7) star, we correct the transit light curves for the effect of stellar activity using ground-based photometric monitoring data from the All-Sky Automated Survey for Supernovae (ASAS-SN) and Tennessee State University's Automatic Imaging Telescope (AIT). We bin the data in 38 spectrophotometric light curves from 0.29 to 4.5 microns and measure the transit depths to a median precision of 90 ppm. We compare the transmission spectrum to a grid of forward atmospheric models and find that our results are consistent with a cloudy spectrum and evidence of sodium at 2.3-sigma confidence, but no observable evidence of potassium absorption even in the narrowest spectroscopic channel. We find that the optical transmission spectrum of WASP-52b is similar to that of the well-studied inflated hot Jupiter HAT-P-1b, which has comparable surface gravity, equilibrium temperature, mass, radius, and stellar irradiation levels. At longer wavelengths, however, the best fitting models for WASP-52b and HAT-P-1b predict quite dissimilar properties, which could be confirmed with observations at wavelengths longer than ~1 micron. The identification of planets with common atmospheric properties and similar system parameters will be insightful for comparative atmospheric studies with the James Webb Space Telescope.


INTRODUCTION
Corresponding author: Munazza Alam munazza.alam@cfa.harvard.edu * National Science Foundation Graduate Research Fellow Transiting exoplanets offer unprecedented opportunities for the detection and detailed characterization of planets beyond the Solar System (Struve 1952). From transit observations, we can make inferences about the formation and evolutionary histories of these planets, their bulk compositions, and their atmospheres (Winn 2010). Atmospheric studies can be performed using transmission spectroscopy to constrain atmospheric structure and chemical composition (e.g., Charbonneau et al. 2002;Vidal-Madjar et al. 2003;, secondary eclipses to measure temperature and thermal structure (e.g., Deming et al. 2005;Charbonneau et al. 2008;Sing & López-Morales 2009;Evans et al. 2017), and orbital phase curves to probe atmospheric circulation (e.g., Knutson et al. 2007;Stevenson et al. 2017). The subject of this paper is transmission spectroscopy (Seager & Sasselov 2000;Brown 2001). During transit, light from the host star passes through the atmosphere of the planet. At wavelengths where absorption by atoms, molecules, and aerosols takes place, the planet blocks slightly more stellar flux, resulting in variations in the apparent radius of the planet as a function of wavelength. These variations in planetary radius reveal the composition of the planetary atmosphere.
The gaseous atmospheres of short-period hot Jupiters are most accessible to such observations because of their large scale heights and short orbital periods. Narrow peaks of H I (1215 Å) in the UV, Na I (5893 Å) and K I (7665 Å) in the optical, H 2 O (1.4 µm; 1.9 µm), CO (2.3 µm; 4.7 µm), CO 2 (2 µm; 15 µm), and CH 4 (2.2 µm; 7.5 µm) in the near-infrared, and scattering by molecular hydrogen (H 2 ) are expected to be prominent features in clear hot Jupiter atmospheres (Vidal-Madjar et al. 2003;Seager & Sasselov 2000;Sudarsky et al. 2003;Burrows et al. 2010;Fortney et al. 2010). Charbonneau et al. (2002) detected the first exoplanet atmosphere for the hot Jupiter HD 209458b, which was later confirmed to have absorption from Na I, H 2 Rayleigh scattering, and possible TiO/VO absorption Snellen et al. 2008;Lecavelier Des Etangs et al. 2008;Désert et al. 2008;). The first near-UV to IR (0.3-8.0 µm) transmission spectrum of the hot Jupiter HD 189733b  revealed clouds/hazes consistent with Rayleigh scattering by small condensate particles in addition to narrow peaks of Na I and K I, H 2 O absorption, and an escaping H atmosphere (Grillmair et al. 2008;Lecavelier Des Etangs et al. 2010;Bourrier et al. 2013).
To date, a diversity of hot Jupiters with a continuum of clear to cloudy atmospheres ) has been detected, with no apparent correlation between the observed spectra and other system parameters. Space-based atmospheric studies have yielded detections of Na I and K I (e.g., Charbonneau et al. 2002;Nikolov et al. 2014), elucidated the presence of thick atmospheric cloud decks (e.g., , and have provided water abundance constraints (e.g., Kreidberg et al. 2014;Wakeford et al. 2018). From groundbased transmission spectral surveys of hot Jupiters, we have detected Na I (e.g., Nikolov et al. 2016Nikolov et al. , 2018Wyttenbach et al. 2015, K I (e.g., , cloudy/hazy atmospheres (e.g., Jordán et al. 2013;Mallonn et al. 2016;Huitson et al. 2017), and Rayleigh scattering slopes (e.g., Gibson et al. 2017).
Here we present results for WASP-52b (Hébrard et al. 2013) from the Hubble Space Telescope (HST) Panchromatic Comparative Exoplanetology Treasury (PanCET) program (GO 14767;PIs Sing & López-Morales). The scientific goals of PanCET are to provide a uniform, statistically compelling ultraviolet (UV) through infrared (IR) study of clouds/hazes and chemical composition in exoplanet atmospheres, and assemble a UVOIR legacy sample of exoplanet transmission spectra that will be well-suited for follow-up with the James Webb Space Telescope (JWST).
WASP-52b is a 0.46 M Jup and 1.27 R Jup inflated hot Jupiter (T eq = 1300 K) orbiting a moderately active (log(R HK ) = −4.4 ± 0.2, Hébrard et al. 2013; log(L x /L bol ) = −4.7, Section 3.2) K2V star with a period of 1.75 days (Hébrard et al. 2013). This planet, at a spectroscopic parallax distance of 175.7±1.3 pc (Gaia Collaboration et al. 2018), is a favorable target for atmospheric studies via transmission spectroscopy due to its large scale height (H = 700 km) and deep transit (δ = 0.028). Based on its surface gravity (log(g) = 2.87 dex) and equilibrium temperature (T eq = 1315 K), WASP-52b is predicted to have a predominantly cloudy atmosphere with muted spectral features (Stevenson 2016). However, two recent ground-based atmospheric analyses of this target claim discrepant conclusions. Louden et al. (2017) cite an optically thick cloud deck to explain an observed flat transmission spectrum, but note that their results are inconsistent with deeper transit depths at longer wavelengths (Kirk et al. 2016). Conversely, Chen et al. (2017) report a cloudy atmosphere with a noticeable Na I detection and a weaker detection of K I absorption.
In this work, we measure WASP-52b's transmission spectrum over the ∼0.29−4.5 µm wavelength range by combining HST/STIS and Spitzer/IRAC observations. An outline of the paper is as follows. In Section 2, we describe the observations and data reduction techniques. In Section 3, we detail the stellar activity correction. The light curve fits and measurement of the transmission spectrum are described in Section 4. We compare our results to previously published measurements and to a grid of forward atmospheric models in Section 5, and present an interpretation of the transmission spectrum in Section 6. We summarize the paper in Section 7.

Observations
We obtained time series spectroscopy during three transits of WASP-52b with HST/STIS on UT 2016 November 01 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) The G430L and G750L STIS data have resolving power R∼500 and each consist of 37 stellar spectra taken over four, consecutive 96-minute orbits. The visits were scheduled to include the transit event in the third orbit, providing an outof-transit baseline time series before and after the transit as well as good coverage between second and third contact. We used a 128 pixel wide subarray and exposure times of 253 seconds to reduce the readout times between exposures. To minimize slit losses, the data were taken with the 52 x 2 arcsec 2 slit.

Spitzer/IRAC
Two transits of WASP-52b were observed on UT 2016 October 18 and UT 2018 March 22 with the Spitzer space telescope (Werner et al. 2004) using the 3.6 µm and 4.5 µm IRAC channels (Fazio et al. 2004). Although these observations cover the complete phase curve of the planet, for this work we only use a 6-hour portion of the phase curve centered on the transit event with enough out-of-transit baseline flux to allow for accurate analysis. Each IRAC exposure had an effective integration time of ∼2 seconds, resulting in ∼30,000 images for the portion of the phase curve corresponding to the transit.

HST/STIS
We reduced (bias-, dark-and flat-corrected) the raw 2D G430L and G750L spectra using the CALSTIS 1 pipeline (version 3.4) and the relevant calibration frames. Following the procedure detailed in Nikolov et al. (2014), we used median-combined difference images to identify and correct 1 http://www.stsci.edu/hst/stis/software/ analyzing/calibration/pipe_soft_hist/intro.html for cosmic ray events and bad pixels flagged by CALSTIS. We extracted 1D spectra from the calibrated .flt science files using IRAF's APALL task. To identify the most appropriate aperture, we extracted light curves of aperture widths ranging from 6 to 18 pixels with a step size of 1. We defined the best aperture for each grating according to the lowest photometric dispersion in the out-of-transit baseline flux. Based on this criterion, we used an aperture size of 13 pixels in our analysis. For each exposure, we computed the mid-exposure time in MJD.
Aperture extractions with no background subtraction minimize the out-of-transit standard deviation of the white light curves . Although the STIS 2D spectra are known to show negligible background sky contribution Huitson et al. 2012Huitson et al. , 2013Nikolov et al. 2015;Gibson et al. 2017), we assessed the potential bias of the sky background on the light curves by obtaining time series spectroscopy both with and without background subtraction. Comparing the light curves, we find that both data sets are fully consistent. To obtain a wavelength solution, we used the x1d files from CALSTIS to re-sample all of the extracted spectra and cross-correlate them to a common rest frame. The cross-correlation measures the shift, and the spectra are re-sampled to align them and remove sub-pixel drifts in the dispersion direction. These drifts can be associated with the different locations of the spacecraft on its orbit around the Earth (e.g., Huitson et al. 2013).

Spitzer/IRAC
We analyzed the IRAC photometry following the methodology of Nikolov et al. (2015) and Sing et al. ( , 2016. We started our analysis using the Basic Calibrated Data (.bcd) files and converted the images from flux in mega-Jansky per steradian (MJy sr −1 ) to photon counts (i.e., electrons) by multiplying each image by the gain and exposure time and then dividing by the flux conversion factor. Following the reduction procedure outlined in Knutson et al. (2012) and Todorov et al. (2013), we filtered for outliers (hot or lower pixels in the data) by following each pixel in time. We scanned the images in two passes: first by removing outliers ≥8σ away from the median value of each frame compared to the 10 surrounding images and then by removing outliers above the 4σ level following the same procedure. The total fraction of corrected pixels was 0.05%. We estimated and subtracted the sky background for each image using an iterative 3σ clipping procedure in which we excluded all pixels associated with the stellar point spread function (PSF), background stars, or hot pixels. In the last iteration, we created a histogram from the remaining pixels and determined the sky background based on a Gaussian fit to the distribution of remaining pixels. To locate the center of the PSF for each image, we used the flux-weighted centroiding method with a 5-pixel radius circular region centered on the approximate position of the star. The variation of the x and y positions of the PSF on the detectors were measured to be 0.19 and 0.24 pixels, respectively.
We extracted photometric points following fixed and timevariable photometry. In the fixed approach, we used circular apertures ranging in radius from 4 to 8 pixels in increments of 0.5. For the time-variable photometry, the aperture size was scaled by the value of the noise pixel parameter (the normalized effective background area of the IRAC point response function), which depends on the full-width-at-half-maximum (FWHM) of the stellar PSF squared (Mighell 2005;Knutson et al. 2012;Lewis et al. 2013;Nikolov et al. 2015). We identified the best results from both photometric methods by comparing the light curve residual dispersion, as well as the white and red noise components measured with the wavelet technique detailed in Carter, & Winn (2009). The time-variable method resulted in the lowest white and random red noise correlated with data points co-added in time.
Since WASP-52 is a moderately active star (log(R HK ) = −4.4±0.2, Hébrard et al. 2013), we used ground-based activity monitoring data to track stellar activity levels during the epochs of our transits. Before fitting the transit light curves, we corrected the baseline flux levels for the effect of stellar variability using a quasi-periodic Gaussian process regression model and corrected for the effect of unocculted stellar spots following the prescription of Huitson et al. (2013). Table 2 summarizes the photometric monitoring campaigns, and Table 3 includes the average flux correction for our transit observations.

Stellar Variability Monitoring
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 (AIT). These observations, however, do not include the epochs of the Spitzer and groundbased (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 2 (ASAS-SN) program (Shappee et al. 2014;Kochanek et al. 2017 Table 2.

Estimating Activity Levels
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 3 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 ISM absorption assumed to be N H = 4 × 20 cm −3 (based on target location and distance), consistent with the fit. The resulting 1-T fit has log T (K)= 6.54 +1.06 −0.29 and log EM (cm −3 ) = 51.14 +0.33 −0.63 . We measured an X-ray luminosity of L x = (3.5 ± 1.0) ×10 28 erg s −1 in the 0.12−2.48 keV band (Sanz-Forcada et al., in prep.). 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 Sanz-Forcada et al. (in prep.)

Modeling the Variability Monitoring Data
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 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 Gaussian process (GP) regression, a framework which 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;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 Gaussian process model for all seasons modeled jointly and for each observing season separately. 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 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).

Correcting for Unocculted Spots
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 intransit measurements, we estimated the non-spotted 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 was then normalized to the non-spotted 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 0 = 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 e f f = 5000 K, log(Z) = −1.5, log(g) = 4.5) and a spot model (T e f f = 4750 K, log(Z) = −1.5, log(g) = 4.5) from the Kurucz (1993) 1D ATLAS grid 4 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−4750 K (corresponding to temperature differences between the spot and stellar photosphere from 1500 K 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 star spots in the data.
We interpolated the stellar model to the spot model grid and computed the wavelength-dependent correction factor derived in : where F λ,Tspot is the stellar model flux at temperature T spot and at the wavelength of the transit observations, F λ,Tstar is the stellar model flux at the wavelength of the transit observations and at temperature T star , F λo,Tspot is the stellar model flux at temperature T spot at the reference wavelength of the activity monitoring data, and F λo,Tstar is the stellar model flux at temperature T star at the activity monitoring reference wavelength. The final flux dimming correction was then calcu- We computed the average flux correction over each bandpass (see Table 3) and applied the correction to each spectrophotometric light curve using: where y corrected is the corrected light curve flux, y is the original (uncorrected) light curve, and y 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.

LIGHT CURVE FITS
We extracted the broadband transmission spectrum and fit the spectrophotometric light curves following the methods described in Sing et al. ( , 2013 and Nikolov et al. (2014), and they are briefly summarized here. To simultaneously fit for the transit and systematic effects, we modeled each of the STIS and Spitzer transit light curves with a two-component function consisting of a transit model multiplied by a systematics model. We adopted the complete analytic transit models of Mandel & Agol (2002), which are parametrized by the mid-transit time T 0 , orbital period P, inclination i, normalized planet semi-major axis a/R , and planet-to-star radius ratio R p /R .

STIS White Light Curves
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 time-scale resulting from thermal expansion/contraction during the spacecraft's day/night orbital cycle.
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, since 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 B1 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 0 , 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 leastsquares 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 B2 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, since the light curves show no significant correlation to the quantified systematic parameters (Nikolov et al. 2014).

STIS Spectrophotometric Light Curves
We produced the spectrophotometric light curves by dividing the spectra into 17 − 400 Å width bins 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 Orbital eccentricity e 0.0 (fixed) 0.0 (fixed) 0.0 (fixed) 0.0 (fixed) Scaled semi-major axis a/R 7.38 ± 0.10 7.14 ± 0.12 7.23 ± 0.12 7.22 ± 0.07 Radius ratio Rp/R 0.1646 ± 0.0020 0.1608 ± 0.0018 0.1741 ± 0.0063 0.1639 ± 0.0005 a The values reported here are the weighted mean of fitted system parameters from the Spitzer observations. 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 mid-transit time T 0 and the scaled semi-major 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 colorindependent systematic trends from each spectral bin. Then, we fit the common mode corrected spectroscopic light curves by fitting for residuals with a parametrized model of six fewer free parameters (c 1 −c 4 , T 0 , 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 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 non-linear 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, 5, and 6.

IRAC Light Curves
We modeled the 3.6 µm and 4.5 µm IRAC transit photometry in accordance with the methods outlined in Sing et al. ( , 2013 and Nikolov et al. (2014). To correct for flux variations from intrapixel sensitivity, we fit a polynomial to the stellar centroid position (Reach et al. 2005;Charbonneau et al. 2005Charbonneau et al. , 2008Knutson et al. 2008). This technique is effective on short timescales (<10 hours) 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: where f (t) is the stellar flux as a function of time, the coefficients a 0 through a 6 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 0 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.

Limb Darkening Models
We modeled the limb darkening of WASP-52b using the four parameter non-linear limb darkening law (Claret 2000) given by: 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 ). 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.

RESULTS
The broadband STIS+Spitzer transmission spectrum for WASP-52b corrected for stellar activity and compared to theoretical forward atmospheric models (Goyal et al. 2018) and past transmission spectrum measurements (Chen et al. 2017;Louden et al. 2017) is shown in Figure 9. The transmission spectrum shows evidence of Na I absorption (5893 Å) at 2.3σ confidence and no observable detection of K I absorption.
To visualize the amplitude of the spot corrections, we also compare the raw transmission spectrum (before applying the stellar activity correction) and the spot corrected spectrum in Figure 10.

Constraints on Na I & K I
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−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.
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.

Fits to Forward Atmospheric Models
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. 2015Tremblin et al. , 2016Drummond et al. 2016). The models are generated for the parameters (e.g., mass, radius, gravity, etc.) of WASP-52b.
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 χ 2 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 (χ 2 =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 super-solar C/O (C/O=0.70). The selected clear model (χ 2 =49.5) has a lower temperature (T = 1015 K) and no clouds (α cloud =0.00) or hazes (α haze =1). The representative hazy model (χ 2 =43.3) is similar to the clear model, but with extreme haziness 5 https://bd-server.astro.ex.ac.uk/exoplanets/ WASP-52/  (α haze =1100). The flat model (χ 2 =52.2) represents a featureless (gray) spectrum. The χ 2 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 super-solar C/O ratio (>0.56) and an equilibrium temperature of 1315 K.

Comparison with Previous Results
In addition to the combined STIS+Spitzer transmission spectrum we report here, there are two ground-based optical transmission spectrum measurements for WASP-52b.  With the WHT/ACAM observations, Louden et al. (2017) modeled spot-crossing events via Gaussian processes, 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, comprised 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 c.f. Chen et al. 2017).

Comparison to HAT-P-1b
We compared WASP-52b to the well-studied inflated hot Jupiter HAT-P-1b (c.f., Nikolov et al. 2014), since 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).
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.
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 C.1 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−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 nearinfrared discrepancy, see Section 6.2. 6. DISCUSSION

Interpreting Alkali Detections in the Presence of Clouds
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 , HAT-P-1b (Nikolov et al. 2014), and WASP-17b ). 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  and HAT-P-12b . We note, however, that the majority of current sodium and potassium detections in exoplanet atmospheres are low significance and a non-detection 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 ln[Na/K] = 8.32 +6.51 −6.03 . 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 ranging from 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 2 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 2 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 comprised 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 comprised 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.

Contextualizing WASP-52b
Although the atmospheres of planets studied thus far appear diverse ), 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, since 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 that 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 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 2 O and CH 4 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 2 O absorption at >5σ confidence (Wakeford et al. 2013), and the deep H 2 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 2 O feature at 1.4 µm, and evidence of H 2 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 nearinfrared 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.

SUMMARY
Our key results are summarized as follows: • We present an optical to near-infrared transmission spectrum of WASP-52b (measured to a median precision of 90 ppm) from ∼0.29 − 5.0 µm using transit observations from HST/STIS and Spitzer/IRAC. We correct for the effects of stellar activity and fit the observed transmission spectrum to a grid of forward atmospheric models.
• Based on these fits (Figure 9), we find that our transmission spectrum measurement best matches a moderately cloudy atmospheric model with an equilibrium temperature of 1315 K, a thick cloud deck (α cloud = 1.00), a slight Rayleigh scattering slope in the blue (α haze = 10), and hints of a 2.3σ Na I signal at 5893 Å. Within the precision of our observations, we do not detect K I absorption.
• We compare the observed transmission spectra of HAT-P-1b and WASP-52b, two planetary systems with similar stellar and planetary parameters (Table 6). By constructing optical HST/STIS transmission spectra with similar binning schemes (Figure 13), we find that the spectra of these two planets are identical within the uncertainties at optical wavelengths but differ in the near-infrared ( Figure 14) based on our best fit models.
• The difference in the transmission spectra of WASP-52b and HAT-P-1b from ∼1.0 − 5.0 µm may be caused by the presence of an extra optical absorber in the atmosphere of HAT-P-1b (Nikolov et al. 2014) or uncertainties in the best fitting models, which are isothermal and therefore cannot accurately capture cloud formation.
• Comparative atmospheric observations with JWST for WASP-52b and HAT-P-1b will be key to understanding planets with similar system parameters and overlapping atmospheric properties.
In a forthcoming paper, we aim to combine the STIS+Spitzer transmission spectrum presented here with existing nearinfrared HST/WFC3 observations (Bruno et al. 2018). Using the full optical to near-infrared transmission spectrum, we will retrieve the planet's atmospheric properties (Bruno, Alam, et al., in prep.) to better constrain the atmospheric structure and chemical composition of this inflated hot Jupiter as well as precisely estimate the Na I and K I abundances in the planet's atmosphere. Such an analysis will indicate if comparative planetology and comparative atmospheric studies of WASP-52b with future JWST observations will prove insightful.
The authors thank the anonymous referee for helpful comments that greatly improved this manuscript. This paper makes use of observations from the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by the Association of Universities          Scale Heights (km) Figure 10. Comparison of the raw (gray diamonds) and spot corrected (black circles) transmission spectra for WASP-52b. The best fitting (solid blue line) and gray atmosphere (dashed black line) models from Figure 9 are shown for reference.   Figure 12. χ 2 map for WASP-52b for the no rainout ATMO model grid. The cloud, haze, and metallicity axes are log-scaled. The contours show all combinations of the grid parameters, and the colors indicate the confidence intervals corresponding to the colorbar to the right. The white regions correspond to parameter space on the grid that are not feasible given current observations at 4σ confidence and can therefore be easily ruled out.  Figure 13. Top panel: Comparison of the observed STIS transmission spectra of WASP-52b (black circles) and HAT-P-1b (purple squares; Nikolov et al. 2014) for identical wavelength bins. The HAT-P-1b spectrum is offset by an arbitrary constant such that the first spectroscopic bin is anchored to the first Rp/R value of the WASP-52b spectrum. Bottom panel: Difference between the observed WASP-52b spectrum and the offset HAT-P-1b spectrum. The spectra of both planets are identical within the uncertainties (except for one channel at ∼0.7 µm), with an average 1σ difference. The horizontal black dashed line indicates the baseline offset between the two spectra.

R p /R + Constant
HAT-P-1b WASP-52b Figure 14. Comparison of the best fit models for HAT-P-1b (purple) and WASP-52b (black) from ∼0.29 to 5 µm. The spectra are identical within the uncertainties at optical wavelengths but differ in the near-infrared.

APPENDIX
A. STELLAR ACTIVITY CORRECTION

A.1. Kernel for Gaussian Process Regression Model
We used a Gaussian process (GP) regression to model the ground-based stellar activity monitoring data Haywood 2015;Angus et al. 2018). In our GP analysis, we used a three component kernel of the form K = k 1 + k 2 + k 3 . The k 1 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 1 = A 2 [k(r 2 ) × k(x i , x j )]), where k(r 2 ) = e r 2 /2 (A1) represents the squared exponential kernel for periodic variations in the time series data parametrized by r. The exponential sine squared kernel is given by where Γ represents the scale of the correlations and P is the period of the oscillations. The second term, k 2 , represents the irregularities in amplitude and period. It is a rational quadratic kernel of the form: k 2 (r 2 ) = A 2 1 − r 2 2α (A3) for amplitude A, Gamma distribution parameter α, and length scale r. The k 3 term incorporates stellar noise in the GP model and is a squared exponential kernel (of the form shown in Equation A1) added to a white kernel of the form k 3 (x i , x j ) = cδ i j for constant c and diagonal value δ i j .

A.2. Deriving the Wavelength-Dependent Flux Correction
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 ) ; see also Section 3.4, Equation 1). This correction depends on the parameter k, which provides an assumption for the non-spotted flux on the stellar surface (Aigrain et al. 2012). We fixed k to unity based on Aigrain et al. (2012), which 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 ranging from 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 re-scaled 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.

B.1. STIS White Light Curve Systematics Models
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 re-sampling the spectra.
For both the G430L and G750L STIS observations, we used the 25 systematics models listed in Table B1 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 B2 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 (AIC) value. . Wavelength-dependent flux correction for a fixed spot temperature (4750 K) and different values of the parameter k, a proxy for the non-spotted fraction of the stellar surface. For the activity correction described in Section 3.4, we assume k=1 (i.e., the spot contribution that the viewer never sees (or always sees) is about the same as the contribution of spots that come into and out of view).

B.2. STIS Spectrophotometric Light Curves
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 B3 and B4 below, 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.
C. HAT-P-1b FORWARD MODEL FITS Figure C.1 shows the HST/STIS transmission spectrum of the well-studied inflated hot Jupiter HAT-P-1b from Nikolov et al. (2014) compared to the best fitting theoretical forward model from the ATMO grid generated for the parameters of HAT-P-1b 6 in addition to representative clear and cloudy models. We performed these fits using the procedure described in Section 5.2. We find that the best fit model has T = 1322 K with a strong Na I signal and is slightly cloudy (α cloud =0.20) and hazy (α haze =10).