## INVESTIGATING TRENDS IN ATMOSPHERIC COMPOSITIONS OF COOL GAS GIANT PLANETS USING SPITZER SECONDARY ECLIPSES

Nicole L. Wallack , Heather A. Knutson , Caroline V. Morley 1 1 2 3 , , Julianne I. Moses , Nancy H. Thomas , 4 1 Daniel P. Thorngren , Drake Deming , Jean-Michel D´ esert , Jonathan J. Fortney , and Joshua A. Kammer 5 6 7 8 9 1 Division of Geological &amp; Planetary Sciences, California Institute of Technology, Pasadena, CA 91125, USA 2 Department of Astronomy, Harvard University, Cambridge, MA 02138, USA 3 Department of Astronomy, University of Texas at Austin, Austin, TX 78712, USA 4 Space Science Institute, 4750 Walnut Street, Suite 205, Boulder, CO 80301, USA 5 Department of Physics, University of California, Santa Cruz, CA 95064, USA 6 Department of Astronomy, University of Maryland at College Park, College Park, MD 20742, USA 7 Anton Pannekoek Institute for Astronomy, University of Amsterdam, 1090 GE Amsterdam, The Netherlands 8 Department of Astronomy and Astrophysics, University of California, Santa Cruz, CA 95064, USA

9 Southwest Research Institute, San Antonio, TX 78238, USA

Draft version August 2, 2019

## ABSTRACT

We present new 3.6 and 4.5 µ m secondary eclipse measurements for five cool (T glyph[lessorsimilar] 1000 K) transiting gas giant planets: HAT-P-15b, HAT-P-17b, HAT-P-18b, HAT-P-26b, and WASP-69b. We detect eclipses in at least one bandpass for all planets except HAT-P-15b. We confirm and refine the orbital eccentricity of HAT-P-17b, which is also the only planet in our sample with a known outer companion. We compare our measured eclipse depths in these two bands, which are sensitive to the relative abundances of methane versus carbon monoxide and carbon dioxide, respectively, to predictions from 1D atmosphere models for each planet. For planets with hydrogen-dominated atmospheres and equilibrium temperatures cooler than ∼ 1000 K, this ratio should vary as a function of both atmospheric metallicity and the carbon-to-oxygen ratio. For HAT-P-26b, our observations are in good agreement with the low atmospheric metallicity inferred from transmission spectroscopy. We find that all four of the planets with detected eclipses are best matched by models with relatively efficient circulation of energy to the nightside. We see no evidence for a solar-system-like correlation between planet mass and atmospheric metallicity, but instead identify a potential (1 9 . σ ) correlation between the inferred CH /(CO+CO ) ratio and stellar metallicity. 4 2 Our ability to characterize this potential trend is limited by the relatively large uncertainties in the stellar metallicity values. Our observations provide a first look at the brightness of these planets at wavelengths accessible to the James Webb Space Telescope , which will be able to resolve individual CH 4 , CO, and CO 2 bands and provide much stronger constraints on their atmospheric compositions.

## 1. INTRODUCTION

Observations of the ever-expanding ensemble of exoplanetary systems provide unique statistical insights into the formation and evolution of planetary systems. This is perhaps best illustrated by the classic correlation between gas giant planet frequency and host star metallicity (Fischer &amp; Valenti 2005), which suggests that these planets most likely formed via core accretion (e.g., Pollack et al. 1996; Johansen &amp; Lambrechts 2017). Observations of the masses and radii of extrasolar gas giant planets also indicate that, like the giant planets in our solar system, the average densities of these planets tend to increase with decreasing mass (Miller &amp; Fortney 2011; Thorngren et al. 2016). These trends are consistent with a picture in which Jovian-mass planets were able to accrete substantially more gas from the protoplanetary disk than Neptune-mass planets, either because their cores reached the critical mass for gas accretion earlier or because they formed in a region of the disk with a higher gas surface density (e.g., at shorter orbital periods).

growing planet's atmosphere will enrich its bulk metallicity and leave a unique compositional fingerprint that will vary according to its formation location and epoch (e.g., Oberg et al. 2011; Espinoza et al. 2017). ¨ In the solar system, Jupiter has both a smaller core mass fraction and a lower atmospheric carbon-to-hydrogen ratio than Neptune (e.g., Lodders 2003). However, with only one planetary system it is difficult to determine the relative importance of the formation location in determining atmospheric metallicity. If we instead consider the broader population of exo-Jupiters and exo-Neptunes, which presumably originate from a variety of formation locations, we can ask whether exo-Neptunes consistently exhibit higher atmospheric metallicities than exo-Jupiters (and therefore whether such enhancement is largely independent of formation location) or whether both populations span a wide range of atmospheric metallicities that reflect their varied formation locations and accretion histories (Humphries &amp; Nayakshin 2018).

The atmospheric compositions of gas giant planets should in theory allow us to distinguish between these two scenarios, as the incorporation of solids into the nwallack@caltech.edu

In principle we can directly determine the mean molecular weights and corresponding metallicities of transiting planet atmospheres by measuring their wavelengthdependent transit depths or transmission spectra (e.g., Seager &amp; Sasselov 2000). However, a majority of the gas giant planets observed to date have clouds in their day-

Table 1 System Properties

|                                                                                                                 | HAT-P-15b                                                                                                                                 | HAT-P-17b                                                                                                        | HAT-P-18b                                                                                                | HAT-P-26b                                                                                                                            | WASP-69b                                                                                             |
|-----------------------------------------------------------------------------------------------------------------|-------------------------------------------------------------------------------------------------------------------------------------------|------------------------------------------------------------------------------------------------------------------|----------------------------------------------------------------------------------------------------------|--------------------------------------------------------------------------------------------------------------------------------------|------------------------------------------------------------------------------------------------------|
| T ∗ (K) Mass (M Jup ) Radius (R Jup ) T eq (K) a e b ω (deg) b Period (days) c T c (BJD-2,450,000) c References | 5568 ± 90 1.949 +0 . 08 - 0 . 078 1.072 ± 0.043 902 ± 27 0.200 +0 . 026 - 0 . 028 262.5 +2 . 4 - 2 . 9 10.863502(37) 4638.56094(48) 1,2,3 | 5246 ± 80 0.537 ± 0.017 1.010 ± 0.029 791 ± 17 0.3417 ± 0.0036 200.5 ± 1.3 10.3385230(90) 4801.17018(20) 3,4,5,6 | 4803 ± 80 0.200 ± 0.019 0.995 ± 0.052 822 ± 22 < 0.087( < 0.16) ... 5.5080291(42) 4715.02254(39) 3,7,8,9 | 5011 ± 55 0.059 ± 0.007 0.565 +0 . 072 - 0 . 032 1028 ± 21 0.14 +0 . 12 - 0 . 08 46 +33 - 71 4.2345023(15) 5304.65218(25) 6,10,11,12 | 4715 ± 50 0.250 ± 0.023 1.057 ± 0.047 961 ± 20 < 0.11( < 0.23) ... 3.8681382(17) 5748.83422(18) 3,13 |

## Notes.

a Calculated assuming planet-wide heat circulation and zero albedo. Uncertainties on the temperature calculated from the uncertainties on a/R ∗ and T . ∗ The semi-major axis and R ∗ values and their corresponding uncertainties for HAT-P-15b, HAT-P-17b, and WASP-69b are taken from reference 3 and for HAT-P-18b from reference 9. R ∗ and semi-major axis for HAT-P-26b are from Wakeford et al. (2017).

b The orbital eccentricity e and longitude of periapse ω are derived from fits to radial velocity data. c Uncertainties on the last two digits are parenthesized.

References. (1) Kov´cs et al. (2010), (2) Torres et al. (2012), (3) Bonomo et al. (2017), (4) Howard et al. a (2012), (5) Fulton et al. (2013), (6) Mortier et al. (2013), (7) Hartman et al. (2011a), (8) Kirk et al. (2017), (9) Seeliger et al. (2015), (10) Hartman et al. (2011b), (11) Knutson et al. (2014b), (12) Stevenson et al. (2016), (13) Anderson et al. (2014)

night terminator region that attenuate the amplitude of the expected absorption features (e.g., Sing et al. 2016; Barstow et al. 2017), leading to degeneracies between cloud-top pressure and atmospheric metallicity for observations with low signal-to-noise detections (e.g., Benneke &amp; Seager 2012). This problem is especially acute for exo-Neptunes, which typically have smaller planet-star radius ratios and higher surface gravities than their Jovian counterparts, both of which make it more challenging to detect atmospheric absorption during the transit. There are currently only three exo-Neptunes with published transmission spectra (GJ 436b, Knutson et al. 2014a; HAT-P-11b, Fraine et al. 2014; and HAT-P-26b, Wakeford et al. 2017), and of these three HAT-P-26b is the only one with a relatively clear atmosphere and correspondingly strong constraints on its atmospheric metallicity. Interestingly, this planet appears to have an atmospheric metallicity substantially lower than that of Neptune (Wakeford et al. 2017).

## 2017).

Although clouds are problematic for transmission spectroscopy, observations of the thermal emission spectra of these same cloudy planets indicate the presence of strong molecular absorption features (e.g., HD 189733b, Crouzet et al. 2014, Todorov et al. 2014; GJ 436b, Morley et al. 2017). This is due in part to the shorter path length for thermal emission as compared to transmission spectroscopy, which minimizes the scattering opacity (Fortney 2005). We also expect that these tidally locked planets should exhibit day-night temperature gradients that might prevent clouds condensing in the cooler terminator region from extending into the hotter dayside region (e.g., Demory et al. 2013, Parmentier et al. 2016), although meridional advection of cloud particles may also affect the observed cloud properties (Lee et al. 2016, Lines et al. 2018). Secondary eclipse observations of the Neptune-mass planet GJ 436b ( &lt; 800 K) indicate that it has strong molecular features in its emission spectrum that can only be matched by a substantially metalenriched atmosphere (200 -1000 × solar; Stevenson et al. 2010, Moses et al. 2013, Lanotte et al. 2014, Morley et al.

In Kammer et al. (2015) we used broadband emission spectroscopy in the same 3.6 and 4.5 µ m bands to constrain the atmospheric compositions of five transiting gas giant planets with temperatures cooler than 1200 K and masses ranging between 0.3 and 3 M Jup . For these relatively cool hydrogen-rich atmospheres, models predict that the ratio of methane to carbon monoxide and carbon dioxide should act as a sensitive barometer of atmospheric metallicity and the carbon-to-oxygen ratio (Moses et al. 2013). Kammer et al. (2015) leveraged the fact that the 3.6 µ m Spitzer band is sensitive to CH 4 absorption while the 4.5 µ m band is sensitive to CO and CO 2 absorption. The ratio of the measured eclipse depths in these two bands can therefore be used to provide constraints on relative trends in atmospheric composition. Although these measurements hinted at a possible trend in atmospheric metallicity versus planet mass, our sensitivity was limited by the large measurement errors characteristic of these types of observations and our relatively small sample size.

In this study we utilize the Infra-Red Array Camera (IRAC) on board the Spitzer Space Telescope to obtain a combined total of 28 3.6 µ m and 4.5 µ m secondary eclipse observations for a sample of five additional transiting gas giant planets with temperatures below ∼ 1000 K and masses between 0.05 and 2.0 M Jup (see Table 1 for more information). Our targets in this study include HAT-P-15b (Kov´ acs et al. 2010), HAT-P-17b (Howard et al. 2012), HAT-P-18b (Hartman et al. 2011a), HATP-26b (Hartman et al. 2011b), and WASP-69b (Anderson et al. 2014). Of these five planets, HAT-P-26b is the only one with published constraints on its atmospheric metallicity from transmission spectroscopy (Wakeford et al. 2017), with a range of 0 8 . -26 × solar (1 σ ). Optical transmission spectroscopy of HAT-P-18b between 475 and 925 nm from the William Herschel Telescope indicates that it has a featureless spectrum consistent with Rayleigh scattering in this wavelength range (Kirk et al. 2017), but this result is still consistent with a wide range of at-

Table 2 Spitzer Observation Details

| Target    |   λ ( µ m) | UT Start Date   |   Length (h) |   t int (s) 1 |   t trim (h) 2 |   r pos 3 |   r phot 4 |   n bin 5 |   RMS 6 |
|-----------|------------|-----------------|--------------|---------------|----------------|-----------|------------|-----------|---------|
| HAT-P-15b |        3.6 | 2011 Nov 27     |         7.8  |           6   |            1.5 |       2.5 |        2.1 |         4 |    1.4  |
| HAT-P-15b |        3.6 | 2012 May 8      |         7.9  |           2   |            0.5 |       2.5 |        2.2 |        16 |    1.13 |
| HAT-P-15b |        3.6 | 2014 Apr 25     |        11.35 |           2   |            1   |       4   |        2   |         8 |    1.15 |
| HAT-P-15b |        4.5 | 2012 Apr 27     |         7.9  |           2   |            0   |       4   |        2.6 |         2 |    1.14 |
| HAT-P-15b |        4.5 | 2012 Nov 19     |         7.9  |           2   |            2   |       2.5 |        2.3 |         2 |    1.09 |
| HAT-P-15b |        4.5 | 2014 May 27     |        11.35 |           2   |            1.5 |       4   |        2.2 |         8 |    1.1  |
| HAT-P-17b |        3.6 | 2012 Jan 25     |         7.9  |           2   |            0   |       4   |        2   |         2 |    1.14 |
| HAT-P-17b |        3.6 | 2012 Aug 29     |         7.9  |           2   |            2   |       3   |        2.5 |        64 |    1.1  |
| HAT-P-17b |        3.6 | 2014 Sep 02     |         8.6  |           2   |            2   |       3.5 |        2.4 |         4 |    1.21 |
| HAT-P-17b |        4.5 | 2012 Feb 04     |         7.9  |           2   |            0.5 |       4   |        2.7 |        32 |    1.12 |
| HAT-P-17b |        4.5 | 2012 Sep 08     |         7.9  |           2   |            2   |       3   |        2.4 |       128 |    1.1  |
| HAT-P-17b |        4.5 | 2012 Sep 22     |         8.6  |           2   |            2   |       3.5 |        2.4 |         4 |    1.22 |
| HAT-P-18b |        3.6 | 2012 May 19     |         7.9  |           2   |            1.5 |       4   |        2   |         2 |    1.13 |
| HAT-P-18b |        3.6 | 2014 May 22     |         5.8  |           2   |            2   |       3   |        2   |         4 |    1.12 |
| HAT-P-18b |        4.5 | 2011 Aug 28     |        11.9  |           2   |            1.5 |       4   |        2.5 |        64 |    1.15 |
|           |        4.5 | 2014 Jun 2      |         5.8  |           2   |            0.5 |       4   |        2   |         4 |    1.09 |
| HAT-P-26b |        3.6 | 2014 Apr 11     |         7.4  |           2   |            1.5 |       2.5 |        2   |         2 |    1.13 |
| HAT-P-26b |        3.6 | 2014 Apr 24     |         7.4  |           2   |            1   |       3.5 |        2   |         2 |    1.16 |
| HAT-P-26b |        3.6 | 2014 Sep 10     |         7.4  |           2   |            1   |       4   |        2.3 |       512 |    1.22 |
| HAT-P-26b |        3.6 | 2014 Sep 27     |         7.4  |           2   |            1   |       3.5 |        2.2 |         2 |    1.21 |
| HAT-P-26b |        4.5 | 2014 Apr 15     |         7.4  |           2   |            2   |       4   |        2.8 |         8 |    1.14 |
| HAT-P-26b |        4.5 | 2014 May 06     |         7.4  |           2   |            1   |       3.5 |        2.4 |         2 |    1.11 |
| HAT-P-26b |        4.5 | 2014 Sep 15     |         7.4  |           2   |            0   |       3.5 |        2   |         2 |    1.2  |
| HAT-P-26b |        4.5 | 2014 Oct 02     |         7.4  |           2   |            2   |       4   |        2.1 |        64 |    1.19 |
| WASP-69b  |        3.6 | 2014 Jul 22     |         9.8  |           0.4 |            2.4 |       3.5 |        2   |         2 |    1.15 |
| WASP-69b  |        3.6 | 2014 Jul 29     |         9.8  |           0.4 |            2   |       3.5 |        2   |         8 |    1.13 |
| WASP-69b  |        4.5 | 2014 Aug 18     |         9.8  |           0.4 |            0   |       3.5 |        2.3 |         2 |    1.16 |
| WASP-69b  |        4.5 | 2015 Jan 08     |         9.8  |           0.4 |            2   |       3.5 |        2.3 |         2 |    1.13 |

- 1 Integration time
- 2 Initial trim duration
- 3 Radius of the aperture (in pixels) used to determine the location of the star on the array
- 4 Radius of the aperture (in pixels) used for the photometry
- 5 Bin size used for fits

6 Ratio of measured RMS to the photon noise limit

mospheric metallicities. Casasayas-Barris et al. (2017) detected sodium absorption at high spectral resolution in the transmission spectrum of WASP-69b, but did not place any constraints on its atmospheric metallicity.

In Section 2 we describe our photometric extraction and model fits. In Sections 3 and 4, we compare our results to atmosphere models and discuss the corresponding implications for our understanding of trends in atmospheric composition.

## 2. OBSERVATIONS AND DATA ANALYSIS

## 2.1. Photometry and Initial Model Fits

We obtained a minimum of two visits each in the IRAC 3.6 and 4.5 µ mbands (Fazio et al. 2004) for all planets in our sample, with additional observations for lower signalto-noise targets. A majority of these observations were observed in the 32 × 32 pixel subarray mode with an initial thirty-minute observation to allow for settling of the telescope followed by a peak-up pointing adjustment prior to the start of the science observation (Ingalls et al. 2012). The only exceptions are the 3.6 µ m 2011 November observation of HAT-P-15b, the 3.6 µ m 2012 January and 4.5 µ m 2012 February observations of HAT-P-17b, and the 4.5 µ m 2012 August observation of HAT-P-18b, which did not include this initial 30-minute observation and subsequent pointing adjustment. The 2011 observation of HAT-P-15b was also obtained in full-array mode instead of subarray mode. See Table 2 for additional details.

We utilize the standard Basic Calibrated Data (BCD)

images for our analysis and extract photometric fluxes as described in our previous studies (i.e., Knutson et al. 2008; Kammer et al. 2015; Morley et al. 2017). We first calculate the BJD UTC mid-exposure times for each image. We then estimate the sky background in each image by masking out a circular region with a radius of 15 pixels centered on the position of the star, iteratively trimming 3 σ outliers, and fitting a Gaussian function to a histogram of the remainder of the pixels. For the fullarray observation of HAT-P-15b, we determine the sky background using the median flux in an annulus with radii between 15 and 37 pixels centered on the position of the star. We determine the location of the star on the array using flux-weighted centroiding (e.g., Knutson et al. 2008; Deming et al. 2015) with a circular aperture. We consider aperture radii ranging between 2.5 and 4.0 pixels in 0.5 pixel steps and optimize our choice of aperture as described below. We then extract the total flux in a circular aperture centered on the position of the star using the aper routine in the DAOPhot package (Stetson 1987). We consider aperture sizes ranging from 2.0 to 3.0 pixels in steps of 0.1 pixels and from 3.0 to 5.0 pixels in steps of 0.5 pixels.

Some visits also display a ramp-like behavior at early times, which we mitigate by trimming up to two hours from the start of our time series. As discussed in Deming et al. (2015) and Kammer et al. (2015), we find that binning our data before fitting reduces the amount of time-correlated noise in the residuals. We determine the optimal flux-weighted centroiding aperture, photometric

Table 3 Best-fit Eclipse Parameters

| Target    | Band ( µ m)   | Depth (ppm)                       | Brightness Temperature (K)   | Time Offset (days) a                        | Center of Eclipse (Phase)                   | ecos( ω ) b                     |
|-----------|---------------|-----------------------------------|------------------------------|---------------------------------------------|---------------------------------------------|---------------------------------|
| HAT-P-15b | 3.6           | < 180 c c,d                       | < 971                        | 0 ± 0.0585 f                                | 0.4829 ± 0.0054 f                           | -0.0262 +0 . 0082 - 0 . 0084 g  |
| HAT-P-17b | 3.6 4.5       | 118 ± 39 < 149 c                  | 813 +49 - 61 < 708           | 0.0120 +0 . 0120 - 0 . 0130                 | 0.2997 ± 0.0012                             | -0.3146 ± 0.0018                |
| HAT-P-18b | 3.6 4.5       | 437 +146 - 144 e 326 +144 - 147 e | 1004 +78 - 94 783 +77 - 100  | 0.0091 +0 . 0054 - 0 . 0073                 | 0.5016 +0 . 0010 - 0 . 0013                 | 0.0026 +0 . 0015 - 0 . 0021     |
| HAT-P-26b | 3.6 4.5       | < 85 c 265 +68 - 72               | < 949 1087 +91 - 102         | 0.0050 ± 0.0037 0.0045 +0 . 0031 - 0 . 0038 | 0.5012 ± 0.0009 0.5011 +0 . 0007 - 0 . 0009 | 0.0019 ± 0.0014 0.0017 ± 0.0011 |
| WASP-69b  | 3.6 4.5       | 421 ± 29 463 ± 39                 | 1011 ± 17 863 +19 - 20       | 0.0033 ± 0.001                              | 0.5009 ± 0.0003                             | 0.001 ± 0.0004                  |

a Time offset from the predicted center of the eclipse.

Unless otherwise noted, we fit both channels with a common time of secondary eclipse. b Computed using the approximation for a low eccentricity orbit. Therefore, this only serves as a first-order approximation for the more highly eccentric orbit of HAT-P-17b. See P´l et al. (2010) for a detailed discussion of the correct treatment for higher eccentricity orbits. a

- c We report the 2 σ upper limit for the eclipse depth.
- d The solution for the eclipse depth at 4.5 µ m was multimodal, so we report the 2 σ upper limit corresponding to the deepest of the eclipse solutions to be conservative. See Section 2.2 for more details.
- e We report the secondary eclipse depths using a tight uniform prior (see Section 2.2 for more details).
- f Due to the nondetections in both bandpasses, we simply report the time offset and corresponding phase of the Gaussian prior derived from the RV constraints from Bonomo et al. (2017; see Section 2.2 for more details).

g We report the ecos( ω ) from Bonomo et al. (2017).

aperture, trim duration, and bin size for each visit by first fitting a combined instrumental and astrophysical model to each version of the photometry and then calculating the standard deviation of the residuals as a function of bin size stepping in powers of two as described in Kammer et al. (2015). We then calculate the least-squares difference between the measured standard deviation of the residuals and the predicted photon noise limit in each bin, which decreases as the square root of the number of points. We then select the photometric and centroiding apertures, trim duration, and bin size that minimizes this least-squares difference (i.e., the one that is closest to the photon noise at all measured timescales) for use in our subsequent analysis. Because we typically do not detect the eclipse in each individual visit, we fix the time of secondary eclipse to the predicted value (phase=0.5 for HAT-P-18b, HAT-P-26b, and WASP-69b, and using the best-fit radial velocity (RV) solution for HAT-P-15b and HAT-P-17b) during our initial optimization (see Section 2.2 for additional details). Our observations of WASP69b also showed a steep downward trend after the end of the eclipse, which was not well matched with our standard linear function of time. In order to avoid fitting a quadratic function of time, which has the potential to bias our measured eclipse depth, we also trimmed up to 2 hours from the end of the time series for all of the WASP69b observations. We optimized this trim duration in the same manner as for the initial trim duration.

Our model for each visit consists of an eclipse model and an instrumental noise model, which we fit simultaneously. We calculate our eclipse model using the batman package (Mandel &amp; Agol 2002, Kreidberg 2015), where we fix the planet-star radius ratio, orbital inclination, and the ratio of the orbital semi-major axis to the stellar radius ( a/R ∗ ) to the published values for each planet (see Table 1 for references) and allow the eclipse depth and time to vary as free parameters.

The flux we measure also depends on the position of the star on the array in each image. This is due to Spitzer 's well-documented intrapixel sensitivity variations (e.g., Charbonneau et al. 2005; Reach et al. 2005; Morales-

Calderon et al. 2006) combined with an undersampled stellar point spread function (the FWHM in the 3.6 and 4.5 µ m arrays is approximately two pixels for data taken during the post-cryogenic mission). We correct for this effect using the pixel-level decorrelation (PLD) method (Deming et al. 2015). This method uses a linear combination of individual pixel-level light curves as the instrumental noise model, and is therefore able to capture trends due to variations in both the star's position and the width of the stellar point spread function. As in Deming et al. (2015), we utilize a 3 x 3 grid of pixels centered on the location of the star in our model. We remove astrophysical flux variations in each 3 x 3 postage stamp by dividing by the summed flux across all nine pixels. Our final instrumental noise model therefore consists of nine linear coefficients corresponding to the nine individual pixel-level light curves, as well as a linear function of time to capture any long-term trends (i.e., 11 free parameters in total). For all fits we divide out our initial astrophysical model and use linear regression on the residuals to obtain an initial guess for the nine linear PLD coefficients in order to speed up convergence for these highly correlated parameters.

## 2.2. Simultaneous Fits and Choice of Prior on Eclipse Phase

Because a majority of the planets in our study have relatively shallow eclipse depths, we do not expect to detect the eclipse signal in fits to individual visits. We therefore carry out our initial fits to individual visits using a fixed eclipse time. For HAT-P-18b, HAT-P-26b, and WASP-69b the published RVs are consistent with a circular orbit (Knutson et al. 2014b; Bonomo et al. 2017) and we therefore fix the eclipse phase to 0.5. HAT-P-15b and HAT-P-17b have non-zero orbital eccentricities and we therefore fix the eclipse time to the predicted value from the literature (Bonomo et al. 2017). After optimizing our choice of aperture, bin size, and trim duration for each individual visit, we next carry out a joint fit to all of the visits for a given planet. In these fits we assume a common eclipse depth for each bandpass and a common

0.3

Figure 1. Phased light curves for each planet from the simultaneous fits with instrumental effects removed. We show data binned in ten-minute intervals (black filled circles) with error bars corresponding to the scatter in each bin, and overplot the best-fit eclipse model in each bandpass for comparison (red lines). The 2 σ upper limits for the best-fit eclipse depths of the 3.6 µ m and 4.5 µ m bandpasses for HAT-P-15b (see Section 2.2 for more details), the 4.5 µ m bandpass for HAT-P-17b, and the 3.6 µ m bandpass for HAT-P-26b are shown.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-8c674e14-9e8b-4d86-813d-838ec32604e8.jpg)

eclipse phase for all visits regardless of bandpass, and allow these three parameters to vary in our fits. In this case we place a uniform prior on the eclipse time spanning the range of times where the full eclipse would be visible in the data (i.e., we disallow eclipse times that are either partially or fully outside our observational window). We then fix the eclipse time to the best-fit value from the simultaneous fit and revisit our choice of optimal aperture, bin size, and trim duration for each individual visit. Lastly, we rerun the simultaneous fit using these newly optimized light curves.

We estimate the uncertainties on the best-fit eclipse parameters in these simultaneous fits using the affine- invariant Markov chain Monte Carlo (MCMC) ensemble sampler emcee (Foreman-Mackey et al. 2013). Our combined astrophysical and instrumental noise model has 13 free parameters, and we therefore use 60 walkers in our fits in order to ensure sufficient sampling of the model parameter space. We place uniform priors on all of our model parameters except where noted below. We also allow the eclipse depths to take on negative values in our fits (i.e., an increase in flux during the eclipse) in order to avoid biasing our estimate of the eclipse depth by requiring only positive values. We initialize the walkers in a tight cluster centered on the best-fit solution from a Levenberg-Marquardt minimization and carry out an

initial burn-in with a length of 10,000 steps. We then discard this initial burn-in and carry out a subsequent fit with 10 5 steps per chain. We report the median values from our MCMC chains and the corresponding 1 σ uncertainties.

We show the raw photometry for each visit with bestfit instrumental noise models from the simultaneous fit overplotted in Figures A.1 - A.3 in the appendix. Normalized light curves for these visits with best-fit eclipse light curves overplotted are shown in Figures A.4 - A.6. The standard deviation of the residuals as a function of bin size for each visit is shown in Figures A.7 - A.9, with the predicted photon noise limit for each bin size overplotted for comparison. We also combine all visits within the same bandpass and show these averaged light curves in Figure 1.

We modified our approach for HAT-P-26b, which is the only planet in our sample with four visits in each bandpass. In this case a simultaneous fit to all eight visits would require a prohibitively large model with a total of 91 free parameters. We instead elect to fit each bandpass separately, and find that the eclipse is detected at 4.5 µ mbut not at 3.6 µ m. We therefore repeat our fits to the 3.6 µ m data with a Gaussian prior on the eclipse time centered on the best-fit eclipse phase from the 4.5 µ m fits and with a width equal to the ± 1 σ uncertainty on this parameter.

HAT-P-15b was the only planet in our sample with no eclipse detection in either band. Previous RV observations of this planet indicate that it has an orbital eccentricity of 0 190 . ± 0 019 (Kov´cs et al. 2010), . a and we therefore centered our window on the predicted secondary eclipse phase rather than a phase of 0.5. Bonomo et al. (2017) subsequently reported an updated eccentricity constraint of 0 200 . +0 026 . -0 028 . with a corresponding uncertainty in the predicted eclipse phase of ± 0 058 days. . This means that our shortest observational window for this planet (7.8 hours) encompassing the entirety of the eclipse only spanned -1.5 σ to +1.3 σ , while the longest observation (11 hours) spanned -2.0 σ to +3.4 σ . Because we do not expect to detect the eclipse at a statistically significant level in a single visit, the shortest observation window becomes the limiting factor on the effective phase range of our search.

Within this range, we place an upper limit on the eclipse depth in each bandpass by carrying out fits with a Gaussian prior on the eclipse phase based on the RV fit from Bonomo et al. (2017). This results in a multimodal posterior for the best-fit eclipse center time, with one peak corresponding to a fit in which the center of the secondary eclipse occurred at the very beginning of the observations (i.e., the entirety of the eclipse is not within the data) and the other peak centered at the expected time of secondary eclipse (see Figure 2). The fitted 3.6 µ m secondary eclipse depth was consistent with zero in both peaks, but the 4.5 µ m secondary eclipse depth was bimodal. In order to be more conservative, we report the the 2 σ upper limits on the eclipse depths corresponding to the peak that is centered at the earlier time of secondary eclipse (i.e. the time that gives a positive upper limit on the 4.5 µ m eclipse depth).

We also test how the use of a Gaussian prior on the time of secondary eclipse as determined by RV measure-

Figure 2. Posterior probability distribution for the secondary eclipse center time and depth in both bands from a joint fit with a Gaussian prior derived from the RV constraints on the eclipse center time for HAT-P-15b. Contours indicate the 1 σ , 2 σ , and 3 σ bounds on these parameters.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-3bcbbe4a-a40c-427e-b768-9bc78b263ccb.jpg)

ments effects the measured eclipse depths for HAT-P17b. We find that using either a flat prior (the time of secondary eclipse must be centered between -0.1 days and 0.1 days where zero is the time of secondary eclipse predicted using the Bonomo et al. (2017) eccentricity constraints) or a Gaussian prior (with the mean of the distribution occurring at the time of the secondary eclipse corresponding to zero days and a standard deviation of 1 σ as determined from the Bonomo et al. (2017) eccentricity constraints) results in the same measured eclipse depths to within 1 σ . Therefore, we report the best-fit values using the less restrictive prior.

Our observations of HAT-P-18b also proved to be particularly challenging. While Hartman et al. (2011a) reported a non-zero eccentricity for HAT-P-18b with ∼ 2 σ significance, we subsequently acquired additional RV measurements and refit these data in Knutson et al. (2014b), where we found that the orbit was consistent with zero eccentricity (e= 0 11 . +0 15 . -0 08 . ). As a result, we centered the 2011 4.5 µ m and 2012 3.6 µ m observations on the predicted eclipse time for the eccentric orbit from Hartman et al. (2011a), and then centered the subsequent 2014 3.6 and 4.5 µ m observations on an orbital phase of 0.5 (i.e., a circular orbit). Although the 2011 4.5 µ m observation is not centered on a phase of 0.5, it does contain the entirety of the eclipse detected in our simultaneous fits (see Figure A.5). However, the 2012 3.6 µ mobservation only spans the first half of the eclipse.

When we tried to fit the two channels of HAT-P-18b separately as we did for HAT-P-26b, we found marginal 2.7 and 1.2 σ detections in the 3.6 µ m and 4.5 µ m bandpasses, respectively, when tight flat priors were used. However, switching to a joint fit of both bands with a flat prior allowing for the secondary eclipse to occur any time during the full range of orbital phases spanned by our observations resulted in a multimodal posterior on the best-fit eclipse center time (see Figure 3) with marginal

Figure 3. Posterior probability distribution for the secondary eclipse center time and depth in both bands from a joint fit with a relatively broad uniform prior of -0.07 to 0.07 days on the eclipse center time for HAT-P-18b. Contours indicate the 1 σ , 2 σ , and 3 σ bounds on these parameters, while the red lines show the prior constraints used in our final version of the fits for this planet.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-139f29fe-0fd2-4656-957b-77969fca4153.jpg)

eclipse depths for both channels ( &lt; σ 1 ).

The highest peak in this distribution is centered near zero, corresponding to an eclipse center time consistent with an orbital phase of 0.5. If we instead repeat our fits using a tighter uniform prior of -0 018 days to 0.027 . days, therefore excluding the other two weaker peaks, it increases the significance of the detection in both bands to 3.0 and 2.2 σ , respectively. We adopt this version of the fits as our final solution, as there is currently no evidence for a non-zero orbital eccentricity in the RV data for this planet and this peak had the highest posterior probability in our original fit.

## 3. RESULTS

We report the best-fit eclipse depths and times and their corresponding uncertainties in Table 3. We detect the eclipse in both channels with high significance for WASP-69b and with somewhat lower significance for HAT-P-18b. We detect the eclipse at 4.5 µ m but not at 3.6 µ m for HAT-P-26b, detect the eclipse at 3.6 µ m but not at 4.5 µ m for HAT-P-17b, and do not detect the eclipse in either channel for HAT-P-15b.

In order to interpret our results, we first convert the measured eclipse depth in each bandpass to a brightness temperature (e.g. Schwartz &amp; Cowan 2015). We then check for differences in brightness temperature between bands, which are indicative of changes in the shape of the planet's emission spectrum due to molecular features (see Section 4.2). We find that WASP-69b has molecular features (i.e., nonblackbody emission spectra) detected with a significance greater than 3 σ , while HAT-P-17b, HAT-P-18b, and HAT-P-26b differ from the blackbody model by less than 3 σ .

These same brightness temperatures can also be used to estimate the efficiency of heat recirculation between the planet's day- and nightsides. We find that the band- averaged brightness temperatures for all four planets with detected eclipses are consistent with their respective equilibrium temperatures (calculated assuming efficient day-night circulation and zero albedos), suggesting that they have either efficient day-night circulation, nonzero albedos, or a combination of the two (Kammer et al. 2015; Schwartz &amp; Cowan 2015, 2017).

We next use our best-fit eclipse phases to place tighter constraints on the values of ecos( ω ) for each planet. HAT-P-18b, HAT-P-26b, and WASP-69b all have time offsets that are consistent with a circular orbit (within ∼ 3 σ of the eclipse time predicted from a circular orbit). HAT-P-17b was previously known to be eccentric, and our new observations confirm and refine the published eccentricity and longitude of periastron from Bonomo et al. (2017).

## 4. DISCUSSION

## 4.1. Comparison to 1D Atmosphere Models

We compare our best-fit eclipse depths to predictions from 1D atmosphere models. Briefly, these models calculate the temperature structure of the atmosphere assuming both chemical and radiative-convective equilibrium. These models are described in more detail in Fortney et al. (2005, 2008) and Morley et al. (2013, 2017). Cross sections for molecular and atomic species are described in detail in Freedman et al. (2008, 2014), with updates to several species that are described in Marley et al.(2019, in preparation). Moderate resolution spectra are calculated using the thermal emission code described in the appendix of Morley et al. (2015). Stellar spectra are calculated using PHOENIX model atmospheres for the stellar properties given in Table 1. Surface gravities are calculated using the planet masses and radii in Table 1. To calculate the incident flux on the planet, we assume that the planet-star distance is the semimajor axis for both circular and eccentric planetary orbits. We assume that heat is either efficiently redistributed or inefficiently redistributed to the nightside. We calculate model spectra for a range of metallicities from solar to 100 × solar metallicity and a range of C/O ratios from C/O=0.15 to 1.5.

We also develop 1D thermo/photochemical kinetics and transport models for these planets to investigate the possible effects of disequilibrium chemistry (i.e., transport-induced quenching and photochemistry) on the atmospheric composition. These models use the Caltech/JPL KINETICS code (Allen et al. 1981) to solve the continuity equations for 92 neutral H-, C-, O-, and Nbearing species that interact via ∼ 1600 forward-reverse chemical reaction pairs. The reaction list is derived from Moses et al. (2013), and further details of the exoplanet disequilibrium-chemistry modeling can be found in Moses et al. (2011, 2013, 2016). Vertical transport in the models occurs through molecular and eddy diffusion, with the vertical profile of the eddy diffusion coefficient ( K zz ) assumed to be similar to that derived for HD 189733b from general circulation models (Ag´ndez u et al. 2014). Specifically, K zz = 1 × 10 7 / P ( (bar)) 0 65 . cm 2 s -1 in the radiative region in the upper troposphere and middle atmosphere (restricted to never exceeding 10 10 cm 2 s -1 in the upper atmosphere), and a constant-withaltitude value of 10 10 cm 2 s -1 in the convective region

at pressures P greater than 100 bar.

The vertical grid in the disequilibrium model consists of 198 levels separated uniformly in log pressure. The thermal structure is taken from the radiative-convective equilibrium models described above. Zero flux boundary conditions are assumed at the top and bottom boundaries, and chemical-equilibrium abundances are assumed for the initial conditions. The protosolar abundances from Table 10 of Lodders et al. (2009) are assumed to be representative of solar composition, but the models assume that 20.7% of the oxygen is removed at depth as a result of the formation of silicates and other refractory condensates (see Visscher et al. 2010). The solar spectrum at solar minimum is adopted for the stellar ultraviolet spectrum for planets with G- and K-type host stars and the composite M-type stellar spectrum described in Moses et al. (2013) is adopted for planets with M-type host stars. All fluxes are scaled to the appropriate planetstar distance.

## 4.1.1. Chemistry of Cool Hydrogen-rich Atmospheres

We first use these models to examine the effect that varying atmospheric metallicity and the carbonto-oxygen ratio has on the measured 3.6 and 4.5 µ m broadband fluxes from these planets, in order to determine the degree of degeneracy between these two parameters. This topic has been previously explored with an extensive grid of generic equilibrium chemistry models in Molaverdikhani et al. (2018), and also with both forward modeling and atmospheric retrievals for the specific case of GJ 436b (Moses et al. 2013; Morley et al. 2017). These studies indicate that GJ 436b's very low ratio of CH 4 to CO and CO 2 can only be reproduced by relatively high ( &gt; 200 × solar) metallicity atmospheres even when the C/O ratio is allowed to vary as a free parameter. In this study we focus on results using forward models, as our two broadband data points are not sufficient for a full retrieval.

We run a grid of equilibrium models with C/O ratios of [0.15, 0.6, 0.9, 1.5] and metallicities of [1, 3, 10, 30, 100] for the hottest (HAT-P-26b; T eq = 1030 ± 20 K), the coolest (HAT-P-17b; T eq = 790 ± 20 K, comparable to GJ 436b), and intermediate temperature (WASP-69b; T eq = 960 ± 20 K) planets in our sample and determine the expected 3.6 to 4.5 µ m spectral slopes from each of these models. As shown in Figures 4-6 these models indicate that the smallest 3.6 to 4.5 µ m spectral slopes (corresponding to planets that are relatively bright at 3.6 µ m and dim at 4.5 µ m) can only be achieved by models with both relatively high ( &gt; 50 × solar) metallicities and C/O ratios less than ∼ 1 5. . This is because at solar C/O the mixing ratios of H O, CH , CO, and CO 2 4 2 all increase as the atmospheric metallicity increases from 1 × to 100 × solar, leading to greater atmospheric opacity, higher temperatures, and a higher overall emission flux in the continuum regions. However, the CO and CO 2 mixing ratios increase much more rapidly with increasing metallicity than the CH 4 mixing ratio, leading to a significantly lower CH /(CO + CO ) ratio and greater 4 2 absorption in the 4.5 µ m band (where CO and CO 2 absorb) than in the 3.6 µ m band (where CH 4 absorbs).

This overall picture remains true for all but the most extreme C/O ratios. Decreasing the atmospheric C/O

ratio at 1 × solar metallicity from 0.9 to 0.15 on a relatively warm planet like HAT-P-26b has only a small effect on the dominant oxygen and carbon species H O 2 and CO, but has a much greater effect on minor constituents CH 4 and CO , with CH 2 4 being present at high C/O ratios and CO 2 at low C/O ratios. The higher C/O ratio model therefore exhibits more absorption in the 3.6 micron bandpass due to the increased presence of CH , 4 but the relative insignificance of both CH 4 and CO 2 under these conditions means that varying the C/O ratio over this range has much less of an effect on the spectral slope than changing the atmospheric metallicity.

For cooler planets where CH 4 is expected to be the dominant carbon-bearing constituent at solar C/O ratios, extreme changes in C/O ratio can lead to significant changes in atmospheric chemistry. As the C/O ratio is reduced, these cooler planets will eventually experience a transition where the carbon shifts from CH 4 to CO 2 dominated carbon chemistries, with a corresponding major shift in spectral slope. Warmer CO-dominated planets can also transition to a different chemical regime when the C/O ratio becomes large enough that CH 4 becomes a major reservoir of carbon (see Figure 7 and Figure 8). The degree of the degeneracy between the atmospheric metallicity and the C/O ratio will be, at least in part, governed by the equilibrium temperature of the planet. In planets that are cool enough to have methane, metallicity will be the primary driver for the spectral slope (as shown by the fact that the C/O ratio does not greatly impact the spectral slope of HAT-P-17b in Figure 4). For warmer planets without significant methane, the spectral slope will be more strongly influenced by variations in the C/O ratio as shown in Figure 6 and to a lesser degree Figure 5.

As shown in Figures 7 and 8, as well as Figures A.10 - A.13, the inclusion of disequilibrium chemical processes does not appreciably change this picture (e.g., Moses et al. 2013). The disequilibrium models have slightly less CH 4 and more CO and CO 2 in their upper atmospheres as compared to the equilibrium models, but this is a relatively minor shift compared to the change in chemistry as we vary the atmospheric metallicity and C/O ratio. We therefore conclude that under the thermal conditions relevant to the planets in our sample, variations in atmospheric metallicity are the dominant factor shaping the 3.6-4.5 µ m spectral slopes, unless the atmosphere is highly enriched (C/O &gt; 1 -1 5) or depleted (C/O . &lt; 0 1) . in carbon compared to oxygen. Although a more careful consideration of the 3D coupled chemistry and dynamics may alter this picture (e.g., Cooper &amp; Showman 2005; Bordwell et al. 2018; Drummond et al. 2018; Mendon¸ ca et al. 2018; Steinrueck et al. 2018), these effects are expected to be relatively minor, as the composition in the infrared photosphere tends to be homogenized to those of the warmest dayside regions.

## 4.1.2. Model Comparison for Individual Planets

We next compare our measured secondary eclipse depths for each individual planet to a grid of four models, including atmospheric metallicities of either 1 × or 100 × solar and either full recirculation (i.e., complete redistribution of heat to the planet's nightside) or dayside-only recirculation (i.e., redistribution of heat limited to the dayside hemisphere alone). This set of four models rep-

## HAT-P-17b (T eq =790 K) at 130 mbar with Equilibrium Chemistry

Figure 4. Heat map (left) showing the spectral slope expected for equilibrium chemical models with varying C/O ratios and metallicities for HAT-P-17b. Each black point represents a forward model, where we interpolate between models to generate the heat map. We indicate the region of parameter space consistent with this planet's measured spectral slope at the 2 σ level or better with black diagonal lines. The pie charts (right) show the abundances of H 2 O, CH , CO , and CO at a pressure representative of those probed in our observations for 4 2 select models.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-b667faff-0150-4ced-9b95-2489824c3776.jpg)

## WASP-69b (T eq =960 K) at 130 mbar with Equilibrium Chemistry

Figure 5. Heat map (left) showing the spectral slope expected for equilibrium chemical models with varying C/O ratios and metallicities for WASP-69b. Each black point represents a forward model, where we interpolate between models to generate the heat map. We show the best-fit value of this planet's measured spectral slope in red and show the region of parameter space consistent at the 1 σ level or better with black diagonal lines. The pie charts (right) show the abundances of H 2 O, CH , CO , and CO at a pressure representative of those 4 2 probed in our observations for select models.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-452da324-b6bd-4fb6-909f-cc39de3605a0.jpg)

HAT-P-26b (T

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-d8f7f845-caea-4abb-809d-3624084b69e4.jpg)

eq

=1030 K) at 130 mbar with Equilibrium Chemistry

Figure 6. Heat map (left) showing the spectral slope expected for equilibrium chemical models with varying C/O ratios and metallicities for HAT-P-26b. Each black point represents a forward model, where we interpolate between models to generate the heat map. We indicate the region of parameter space consistent with this planet's measured spectral slope at the 2 σ level or better with black diagonal lines. The pie charts (right) show the relative abundances of H 2 O, CH , CO , and CO at a pressure representative of those probed in our observations 4 2 for select models.

HAT-P-17b at 130 mbar (Solar C/O)

Figure 7. Comparison of chemical abundances of H O, CH , 2 4 CO , and CO from the equilibrium chemical models (left) and dise2 quilibrium chemical models (right) for 1 × solar metallicity models (top) and 100 × solar metallicity models (bottom) all with solar C/O ratios (C/O=0.6) at 130 mbar for our coolest planet, HATP-17b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-4a5fe0da-6660-4b3e-9c70-6eb22957c1f4.jpg)

HAT-P-26b at 130 mbar (Solar C/O)

Figure 8. Comparison of the relative abundances of H 2 O, CH , 4 CO , and CO from the equilibrium chemical models (left) and dise2 quilibrium chemical models (right) for 1 × solar metallicity models (top) and 100 × solar metallicity models (bottom) all with solar C/O ratios (C/O=0.6) at 130 mbar for our hottest planet, HATP-26b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-2cdb5c43-d10a-4318-9316-a9f9de6c7a6c.jpg)

resents a reasonably compact sampling of the possible parameter space, with redistribution efficiency allowing us to make each model globally hotter or cooler while varying atmospheric metallicity serves as a simplified proxy for changes in atmospheric chemistry that can affect the spectral slope. The resulting models are shown in Figure 9. Although we include the 3.6 µ m depth for HATP-15b for completeness, we do not detect an eclipse in

Planet-Star Flux Ratio

Wavelength (μm)

Figure 9. Planet-star flux ratios as a function of wavelength from 1D atmosphere models for all five planets. We show 1 × solar and 100 × solar metallicity models with either efficient (purple and blue, respectively) or inefficient (red and orange, respectively) redistribution of energy to the planet's nightside. Our measured eclipse depths for each planet are shown as black circles, and we plot the corresponding band-integrated flux values from the models as filled squares. For HAT-P-15b (3.6 µ m only; see Section 2.2 for more details), HAT-P-17b (4.5 µ m only), and HAT-P-26b (3.6 µ m only) we show the 2 σ upper limit on the eclipse depth. We also overplot the IRAC 3.6 and 4.5 µ m response functions in gray in the bottom panel for comparison.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-d78c45d8-f449-4878-b176-389bc1b68af5.jpg)

either bandpass, and therefore refrain from any further discussion of the implications of the eclipse depths for our understanding of this planet's atmosphere.

As in Kammer et al. (2015), model-data comparisons

for all four planets with detected eclipses strongly prefer models with efficient circulation between the dayand nightside hemispheres. We find that the HAT-P26b eclipses are best matched by the 1 × solar metallicity model, in good agreement with the constraints from transmission spectroscopy presented in Wakeford et al. (2017). The WASP-69b eclipses are well matched by the 100 × solar metallicity model. Our constraints for HAT-P-17b and HAT-P-18b are somewhat weaker, but still appear to modestly favor the 100 × solar metallicity model over the 1 × solar model.

## 4.2. Model-independent Trends in Atmospheric Composition

We next consider the same model-independent metric used in Kammer et al. (2015) to search for trends in atmospheric composition. This metric is defined as the ratio of the measured brightness temperatures in the 4.5-3.6 µ m bandpasses, and should be effectively independent of the planet-star radius ratio and equilibrium planet temperature. We expect variations in this ratio to instead reflect the relative strength of CH 4 (3.6 µ m) versus CO and CO 2 (4.5 µ m) absorption features in the atmospheres of these planets. We use this ratio to search for empirical correlations with other parameters of interest including the planet's mass, bulk metallicity, and host star metallicity.

We note that the presence or absence of a temperature inversion might also alter a planet's relative brightness in these two bands in a way that mimics the shift from a CH -dominated to CO- and CO -dominated carbon 4 2 chemistry. Observations of the broader sample of hot Jupiters suggests that temperature inversions are only found in the atmospheres of the most highly irradiated planets and are most likely caused by the presence of gas phase TiO (e.g., Evans et al. 2017; Nugroho et al. 2017; Sedaghati et al. 2017; Sheppard et al. 2017). The planets in this study are too cool for TiO and VO to remain in the gas phase, and there is currently no evidence for temperature inversions in the atmospheres of planets at these temperatures.

Two of the Neptune-mass planets in our sample, GJ 436b and HAT-P-26b, have previously published constraints on their atmospheric metallicities and C/O ratios. Morley et al. (2017) re-examined all of the available secondary eclipse data for GJ 436b, which has a low 4.53.6 µ m brightness temperature ratio. They concluded that these observations were consistent with absorption features from water, carbon monoxide, and carbon dioxide corresponding to an atmospheric metallicity greater than 200 × solar and a C/O ratio consistent with solar. Similarly, transmission spectroscopy for HAT-P-26b from Wakeford et al. (2017) constrains its atmospheric metallicity to 0 8 . -26 × solar (1 σ ); this is consistent with our measurement of a relatively high 4.5-3.6 µ m brightness temperature ratio. We therefore conclude that for these two planets, differences in atmospheric metallicity are likely the primary factor driving the observed difference in 3 6 . -4 5 . µ m spectral slope, in good agreement with our predictions based on the model grids in Section 4.1.1.

We use the measured spectral slopes for our sample of planets to search for trends in spectral shape across our sample of short-period gas giant planets. If these planets have broadly solar C/O ratios and follow the same trend of increasing atmospheric metallicity with decreasing planet mass that we see for the solar system gas giants, we would expect to see a rising trend in the measured brightness temperature ratios with increasing planet mass. In Figure 10c we plot this brightness temperature ratio as a function of planet mass for the four planets with measured eclipses described herein as well as other planets with temperatures less than 1100 K and published secondary eclipse detections (see Table 4). As before, we see that planets with eccentric orbits (GJ 436b, HAT-P-17b, and WASP-8b) appear to have brightness temperature ratios that are systematically lower than those of planets on circular orbits with the same mass. Other than this apparent clustering of eccentric planets, there appears to be no obvious correlation between atmospheric composition and planetary mass for planets in this temperature regime.

We also check to see if the atmospheric compositions of these planets are correlated with their bulk metallicities. We take published bulk metallicity values from Thorngren &amp; Fortney (2019), which used a 1D planetary model with an inert rock-ice core and a convective H/He-rockice envelope. We calculate new bulk metallicity values for GJ 3470b, which was not included in this study, using the same method (see Table 4). We find evidence for a correlation between atmospheric metallicity and bulk metallicity (Figure 10d), where planets with higher bulk metallicities have on average slightly higher 4.5-3.6 µ m brightness temperature ratios. We evaluate the significance of the proposed trend in brightness temperature ratio versus bulk metallicity using a Monte Carlo simulation where we create a series of simulated data sets by sampling from the posterior probability distributions for each point. For the brightness temperature ratios, we assume that the reported uncertainties in the 3.6 and 4.5 µ m brightness temperatures are reasonably well approximated by Gaussian distributions. We then draw 10 6 samples from these two distributions for each planet.

For planets with a nondetection in one band, we assume that the brightness temperature in the band with the non-detection follows a Gaussian distribution with a standard deviation of 100 K centered at the assumed median of the distribution as determined from the brightness temperature corresponding to the eclipse depth of the reported 2 σ upper limit (i.e., the center of the distribution is the brightness temperature corresponding to the 2 σ upper limit - 200 K); this effectively excludes solutions in which the brightness temperature for the nondetection is unphysically low (i.e., 98% of all samples are limited to temperatures within 400 K of the 2 σ upper limit). Lastly, we assume that the bulk metallicities are also well approximated by Gaussian distributions and sample from those as well.

If we exclude the three Neptune-mass planets (GJ 436b, GJ 3470b, and HAT-P-26b), we do find relatively weak (with positively sloped lines preferred over flat or negatively sloped lines at the 1.4 σ level) evidence for a linear trend with bulk metallicity. If it can be substantiated with additional measurements, this correlation would be somewhat surprising as it would suggest that gas giant planets with high bulk metallicities may have lower atmospheric metallicities. This would be the opposite of the observed trend for the solar system gas giants, in which planets with higher bulk metallicities

Table 4 Brightness Temperature Ratios

| Planet    | Mass (M Jup )         | Radius (R Jup )          | T eq (K)   | T Bright Ratio        | Bulk Z planet   | [Fe/H] ∗     | Ref        |
|-----------|-----------------------|--------------------------|------------|-----------------------|-----------------|--------------|------------|
| GJ 436b   | 0.07 ± 0.01           | 0.37 ± 0.02              | 669 ± 22   | < 0.72                | 0.83 ± 0.02     | 0.05 ± 0.141 | 1,2,3,4    |
| GJ 3470b  | 0.0432 ± 0.0051       | 0.346 ± 0.029            | 604 ± 98   | < 1.00                | 0.74 ± 0.37     | 0.27 ± 0.11  | 5,6,7      |
| HAT-P-12b | 0.21 ± 0.01           | 0.96 +0 . 03 - 0 . 02    | 963 ± 16   | 1.06 +0 . 08 - 0 . 10 | 0.32 ± 0.03     | -0.26 ± 0.06 | 1,8,9,10   |
| HAT-P-17b | 0.537 ± 0.017         | 1.010 ± 0.029            | 791 ± 17   | < 0.87                | 0.13 ± 0.04     | 0.05 ± 0.03  | 1,10,11    |
| HAT-P-18b | 0.200 ± 0.019         | 0.995 ± 0.052            | 822 ± 22   | 0.78 +0 . 08 - 0 . 09 | 0.24 ± 0.05     | 0.10 ± 0.08  | 1,10,11    |
| HAT-P-19b | 0.292 ± 0.018         | 1.132 ± 0.072            | 1010 ± 42  | 0.84 ± 0 . 06         | 0.22 ± 0.05     | 0.29 ± 0.06  | 1,10,12,13 |
| HAT-P-20b | 7.25 ± 0.19           | 0.87 ± 0.03              | 970 ± 23   | 1.00 +0 . 03 - 0 . 04 | 0.28 ± 0.05     | 0.12 ± 0.15  | 1,10,14,15 |
| HAT-P-26b | 0.059 ± 0.007         | 0.565 +0 . 072 - 0 . 032 | 1028 ± 21  | > 1.15                | 0.66 ± 0.03     | 0.01 ± 0.04  | 1,10,11    |
| WASP-8b   | 2.24 +0 . 08 - 0 . 09 | 1.04 +0 . 01 - 0 . 05    | 948 ± 22   | 0.73 ± 0 . 06         | 0.11 ± 0.04     | 0.29 ± 0.03  | 1,10,16    |
| WASP-10b  | 3.14 ± 0.27           | 1.039 +0 . 043 - 0 . 049 | 972 ± 31   | 0.94 ± 0.03           | 0.12 ± 0.02     | 0.04 ± 0.05  | 1,2,10,13  |
| WASP-67b  | 0.406 ± 0.035         | 1.091 ± 0.046            | 1003 ± 20  | > 0.97                | 0.2 ± 0.06      | 0.18 ± 0.06  | 1,10,13,17 |
| WASP-69b  | 0.250 ± 0.023         | 1.057 ± 0.047            | 961 ± 20   | 0.85 ± 0.02           | 0.21 ± 0.04     | 0.30 ± 0.06  | 1,10,11    |
| WASP-80b  | 0.54 ± 0.04           | 1.00 ± 0.03              | 825 ± 19   | 0.99 ± 0.10           | 0.17 ± 0.04     | 0.13 ± 0.11  | 1,7,9,18   |

References. (1) Thorngren &amp; Fortney (2019), (2) T eq from Southworth (2011), (3) Morley et al. (2017), (4) [Fe/H] from Rojas-Ayala et al. (2012), (5) Biddle et al. (2014), (6) Benneke et al. (2019), (7) [Fe/H] from Terrien et al. (2015), (8) Hartman et al. (2009), (9) Wong et al. (2019, in preperation), (10) [Fe/H] from Santos et al. (2013), Sousa et al. (2018), (11) This work, (12) Hartman et al. (2011a), (13) Kammer et al. (2015), (14) Bakos et al. (2011), (15) Deming et al. (2015), (16) Cubillos et al. (2013), (17) Mancini et al. (2014), (18) Triaud et al. (2015)

1050

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-a5a9a383-a875-4efa-9bb3-d53802ef1b24.jpg)

1050

Figure 10. Measured 4.5-3.6 µ m brightness temperature ratio as a function of (a) equilibrium temperature, (b) stellar metallicity, (c) planetary mass, and (d) bulk metallicity for all planets with published Spitzer eclipse depths and equilibrium temperatures less than 1100 K. We show four planets from this study (HAT-P-17b, HAT-P-18b, HAT-P-26b, and WASP-69b) as well as nine previously published planets (see Table 4 for the full list of planets and corresponding references). For planets with no eclipse detected in the 3.6 (4.5) µ m bands we plot 2 σ lower (upper) limits, respectively.

also have higher atmospheric metallicities.

Intriguingly, we do see tentative evidence for a trend in atmospheric composition with stellar metallicity (Figure 10b). This trend is not entirely surprising given that we would expect metal-rich stars to have correspondingly metal-rich disks. However, the metallicities of the host stars in our sample only vary between -0 26 and +0 29, . . corresponding to a relatively small 3.5 × change in bulk disk metallicity. In contrast, we would need to vary the atmospheric metallicities of these planets by more than two orders of magnitude to reproduce the observed change in 4.5-3.6 µ m brightness temperature ratios. If real, this trend suggests that gas giant planets in metalrich disks incorporate disproportionately more solids into their atmospheres than planets in metal-poor disks. If we assume that these planets formed in the inner disk, this might be explained by scenarios in which gas drag causes approximately centimeter-sized particles to migrate inward, increasing the concentration of solids in this region (e.g., Johansen et al. 2014).

We evaluate the significance of this proposed trend in the same manner as before (i.e., by drawing 10 6 samples from the distributions for the brightness temperatures as described above) and assuming that the reported stellar metallicities are well approximated by Gaussian distributions. We then fit the resulting 10 6 simulated data sets with a linear function and plot the corresponding posterior probability distributions for the slope and y intercept of this line (Figure 11). We find conflicting literature values for the metallicity of the M-dwarf stars in our sample (GJ 436, GJ 3470, and WASP-80), many of which also have relatively large uncertainties. We therefore opted to exclude these M-dwarf hosts when evaluating the significance of the proposed trend (Rojas-Ayala et al. 2012; Terrien et al. 2012; Lindgren &amp; Heiter 2017). We find that negative slopes are preferred, but our data are still consistent with a positive sloped or flat line at the 1 9 . σ level.

We also consider a scenario in which the points shown in Figure 11 do not follow a smooth trend, but rather comprise two distinct groups of planets- those with high atmospheric metallicities and those with low atmospheric metallicities. We follow the same approach as Schlaufman (2018) and try a hierarchical clustering algorithm, a k-means clustering algorithm, and a Gaussian-model clustering algorithm. All three algorithms are available as part of the scikit-learn package (Pedregosa et al. 2011). We find that these three methods were inconclusive as to the membership of the planets within the two groups, indicating that these measurements are consistent with a single population.

## 5. CONCLUSIONS

We present new secondary eclipse depth measurements in the 3.6 and 4.5 µ m Spitzer bands for HAT-P-17b, HAT-P-18b, HAT-P-26b, and WASP-69b and place upper limits on the eclipse depths of HAT-P-15b. Our measured times of secondary eclipse for HAT-P-18b, HATP-26b, and WASP-69b are consistent with circular orbits, and we confirm the non-zero eccentricity of HATP-17b. We compare our measured eclipse depths with 1D radiative-convective models for each planet. For HAT-P26b, which was the only planet with a well-constrained atmospheric metallicity from transmission spectroscopy,

Figure 11. Measured 4.5 µ m to 3.6 µ m brightness temperature ratio as a function of stellar metallicity. We create 10 6 simulated data sets by sampling from posterior probability distributions for the reported brightness temperatures and stellar metallicities, and show a random subset of the resulting distribution of linear fits as gray lines. We indicate the median linear solution in red, and plot the corresponding distributions in y intercept and slope below.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-b18c5c5c-c40f-4bc5-82ef-dacf6619c1dc.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-9d9c73ca-421e-426e-b2c6-0ac4332fb4ae.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-693667a3-5976-4c24-b1ab-4902a89f65e8.jpg)

our data are in good agreement with the low atmospheric metallicity reported in Wakeford et al. (2017). We find no evidence for a correlation between atmospheric composition and planetary mass. However, we do find a suggestive 1 9 . σ trend in atmospheric composition as a function of stellar metallicity. While the existence of a correlation between planetary atmospheric composition and stellar metallicity would not be surprising, the strength of the observed trend implies that short-period gas giant planets orbiting metal-rich stars may have atmospheric metallicities that are significantly higher than the bulk metallicity of the disk. However, our ability to fully understand this possible trend is limited by the availability of precise stellar metallicity measurements.

Beginning in 2021, the broad infrared wavelength coverage and higher spectral resolution of the James Webb Space Telescope ( JWST ) will provide invaluable new insights into the atmospheric compositions of this population of planets (see Figure 12 for the expected JWST precision and coverage). In this study we relied on broad photometric bandpasses that span multiple absorption features including water, methane, carbon monoxide, and carbon dioxide. This necessarily results in degeneracies in the interpretation of these data, including correlations between the abundances of the various molecular species and also with the planet's dayside pressuretemperature profile. JWST will be able to resolve individual molecular bands, therefore avoiding these degeneracies and allowing for robust abundance constraints for these atmospheres. It is worth noting that these cooler planets are also ideal targets for studies examining trends in bulk metallicity as a function of planet mass (e.g., Miller &amp; Fortney 2011, Thorngren et al. 2016), as hotter planets have inflated radii that make it difficult to accurately determine their bulk metallicities (e.g., Laughlin et al. 2011, Thorngren &amp; Fortney 2018). We therefore

Figure 12. Simulated JWST observations of a single secondary eclipse observation of WASP-69b using the NIRSpec G395 grism (black points). We use our best-fit 100 × solar metallicity model for this simulation and overplot models with 1 × , 3 × , 10 × and 30 × solar metallicities for comparison (all with solar C/O ratios). The IRAC 3.6 and 4.5 µ m response functions are shown in gray and our measured eclipse depths in the two Spitzer bands are shown as filled black stars.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-e7799bcf-3395-4f09-ac95-79de64d5a094.jpg)

expect that future studies of this cool gas giant population with JWST will be able to determine for the first time whether or not the postulated correlation between planet mass, core mass fraction, and atmospheric metallicity is in fact a universal property of all gas giant planets. This in turn will tell us whether or not the interior and atmospheric compositions of gas giant planets are the inevitable outcome of the core accretion process or instead primarily reflect the diverse formation locations and disk properties of these planets.

## 6. ACKNOWLEDGEMENTS

This work is based on observations made with the Spitzer Space Telescope , which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Support for this work was provided by NASA through an award issued by JPL/Caltech. J.M.D acknowledges that the research leading to these results has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement no. 679633; ExoAtmos). J.M. acknowledges support from NASA grant NNX16AC64G.

## REFERENCES

- Ag´ undez, M., Parmentier, V., Venot, O., Hersant, F., &amp; Selsis, F. 2014, A&amp;A, 564, A73
- Allen, M., Yung, Y. L., &amp; Waters, J. W. 1981, J. Geophys. Res., 86, 3617
- Anderson, D. R., Collier Cameron, A., Delrez, L., et al. 2014, MNRAS, 445, 1114
- Bakos, G. A., Hartman, J., Torres, G., et al. 2011, ApJ, 742, 116 ´ Barstow, J. K., Aigrain, S., Irwin, P. G. J., &amp; Sing, D. K. 2017, ApJ, 834, 50
- Benneke, B., &amp; Seager, S. 2012, ApJ, 753
- Benneke, B., Knutson, H. A., Lothringer, J., et al. 2019, Nature Astronomy, 377
- Biddle, L. I., Pearson, K. A., Crossfield, I. J., et al. 2014, MNRAS, 443, 1810
- Bonomo, A. S., Desidera, S., Benatti, S., et al. 2017, A&amp;A, 602, A107
- Bordwell, B., Brown, B. P., &amp; Oishi, J. S. 2018, ApJ, 854, 8 Casasayas-Barris, N., Palle, E., Nowak, G., et al. 2017, A&amp;A, 608,

A135

- Charbonneau, D. B., Allen, L. E., Megeath, S. T., et al. 2005, ApJ, 626, 523
- Cooper, C. S., &amp; Showman, A. P. 2005, ApJ, L45
- Crouzet, N., McCullough, P. R., Deming, D., &amp; Madhusudhan, N. 2014, ApJ, 795, 166
- Cubillos, P., Harrington, J., Madhusudhan, N., et al. 2013, ApJ, 768, 42
- Deming, D., Knutson, H., Kammer, J., et al. 2015, ApJ, 805, 132 Demory, B.-O., de Wit, J., Lewis, N., et al. 2013, ApJL, 776, L25 Drummond, B., Mayne, N. J., Baraffe, I., et al. 2018, A&amp;A, 512, A105
- Espinoza, N., Fortney, J., Miguel, Y., Thorngren, D., &amp;

Murray-Clay, R. 2017, ApJL, 838, L9

- Evans, T. M., Sing, D. K., Kataria, T., et al. 2017, Nature, 548, 58 Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10 Fischer, D. A., &amp; Valenti, J. 2005, ApJ, 622, 1102

Foreman-Mackey, D., Hogg, D. W., Lang, D., &amp; Goodman, J.

2013, PASP, 125, 306

Fortney, J. J. 2005, MNRAS, 364, 649

- Fortney, J. J., Lodders, K., Marley, M. S., &amp; Freedman, R. S. 2008, ApJ, 678, 1419
- Fraine, J., Deming, D., Benneke, B., et al. 2014, Nature, 513, 526 Freedman, R. S., Lustig-Yaeger, J., Fortney, J. J., et al. 2014, ApJS, 214
- Freedman, R. S., Marley, M. S., &amp; Lodders, K. 2008, ApJS, 174, 504
- Fulton, B. J., Howard, A. W., Winn, J. N., et al. 2013, ApJ, 772 Hartman, J. D., Bakos, G. A., Torres, G., et al. 2009, ApJ, 706, ´ 785
- Hartman, J. D., Bakos, G. A., Sato, B., et al. 2011a, ApJ, 726 ´
- ´
- Hartman, J. D., Bakos, G. A., Kipping, D. M., et al. 2011b, ApJ, 728
- Howard, A. W., Bakos, G. A., Hartman, J., et al. 2012, ApJ, 749 ´
- Humphries, R. J., &amp; Nayakshin, S. 2018, MNRAS, 477, 593
- Ingalls, J. G., Krick, J. E., Carey, S. J., et al. 2012, Proceedings of SPIE, 84421Y
- Johansen, A., Blum, J., Tanaka, H., et al. 2014, Protostars and Planets VI, 547
- Johansen, A., &amp; Lambrechts, M. 2017, AREPS, 45, 359
- Kammer, J. A., Knutson, H. A., Line, M. R., et al. 2015, ApJ, 810, 118
- Kirk, J., Wheatley, P. J., Louden, T., et al. 2017, MNRAS, 3916, 3907
- Knutson, H. A., Benneke, B., Deming, D., &amp; Homeier, D. 2014a, Nature, 505, 66
- Knutson, H. A., Charbonneau, D., Allen, L. E., Burrows, A., &amp; Megeath, S. T. 2008, ApJ, 673, 526
- Knutson, H. A., Fulton, B. J., Montet, B. T., et al. 2014b, ApJ, 785
- Kov´ acs, G., Bakos, G. A., Hartman, J. D., et al. 2010, ApJ, 724, ´ 866
- Kreidberg, L. 2015, PASP, 1161
- Lanotte, A. A., Gillon, M., Demory, B.-O., et al. 2014, A&amp;A, 572, A73
- Laughlin, G., Crismani, M., &amp; Adams, F. C. 2011, ApJL, 729, 3 Lee, G., Dobbs-Dixon, I., Helling, C., Bognar, K., &amp; Woitke, P.

2016, A&amp;A, 594, A48

- Lindgren, S., &amp; Heiter, U. 2017, A&amp;A, 604
- Lines, S., Mayne, N. J., Boutle, I. A., et al. 2018, A&amp;A, 97 Lodders, K. 2003, ApJ, 1220
- Lodders, K., Palme, H., &amp; Gail, H.-P. 2009, Landolt B¨rnstein, o 712
- Mancini, L., Southworth, J., Ciceri, S., et al. 2014, A&amp;A, 568, A127
- Mandel, K., &amp; Agol, E. 2002, ApJ, 580, L171
- Mendon¸ ca, J. M., Tsai, S.-M., Matej, M., Grimm, S. L., &amp; Heng, K. 2018, ApJ, 869, 107
- Miller, N., &amp; Fortney, J. J. 2011, ApJL, 736
- Molaverdikhani, K., Henning, T., &amp; Molli`re, P. 2018, ApJ, 873, e 32
- Morales-Calderon, M., Stauffer, J. R., Kirkpatrick, J. D., et al. 2006, ApJ, 43, 544
- Morley, C. V., Fortney, J. J., Kempton, E. M., et al. 2013, ApJ, 775
- Morley, C. V., Fortney, J. J., Marley, M. S., et al. 2015, ApJ, 815, 110

- Morley, C. V., Knutson, H., Line, M., et al. 2017, AJ, 153, 1 Mortier, A., Santos, N. C., Sousa, S. G., et al. 2013, A&amp;A, 557, A70
- Moses, J. I., Visscher, C., Fortney, J. J., et al. 2011, ApJ, 737, 15
- Moses, J. I., Line, M. R., Visscher, C., et al. 2013, ApJ, 777
- Moses, J. I., Marley, M. S., Zahnle, K., et al. 2016, ApJ, 829, 66
- Nugroho, S. K., Kawahara, H., Masuda, K., et al. 2017, AJ, 154, 221
- ¨ berg, K. I., Murray-Clay, R., &amp; Bergin, E. A. 2011, ApJL, 743 O P´ al, A., Bakos, G. A., Torres, G., et al. 2010, MNRAS, 401, 2665 ´ Parmentier, V., Fortney, J. J., Showman, A. P., Morley, C. V., &amp; Marley, M. S. 2016, ApJ, 828, 1
- Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, J. Mach. Learn. Res., 12, 2825
- Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62
- Reach, W. T., Megeath, S. T., Cohen, M., et al. 2005, PASP, 117, 978
- Rojas-Ayala, B., Covey, K. R., Muirhead, P. S., &amp; Lloyd, J. P. 2012, ApJ, 748
- Santos, N. C., Sousa, S. G., Mortier, A., et al. 2013, A&amp;A, 556, A150
- Schlaufman, K. C. 2018, ApJ, 853, 37
- Schwartz, J. C., &amp; Cowan, N. B. 2015, MNRAS, 449, 4192
- -. 2017, PASP, 129, 1
- Seager, S., &amp; Sasselov, D. D. 2000, ApJ, 10, 916
- Sedaghati, E., Boffin, H. M., MacDonald, R. J., et al. 2017, Nature, 549, 238
- Sheppard, K., Mandell, A. M., Tamburo, P., et al. 2017, ApJL, 850, L32
- Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59
- Sousa, S. G., Adibekyan, V., Delgado-Mena, E., et al. 2018, A&amp;A, 58
- Southworth, J. 2011, MNRAS, 417, 2166
- Steinrueck, M. E., Parmentier, V., Showman, A. P., et al. 2018, arXiv:1808.02011v2
- Stetson, P. B. 1987, PASP, 99, 191
- Stevenson, K. B., Bean, J. L., Seifahrt, A., et al. 2016, ApJ, 817, 141
- Stevenson, K. B., Harrington, J., Nymeyer, S., et al. 2010, Nature, 464, 1161
- Terrien, R. C., Mahadevan, S., Bender, C. F., et al. 2012, ApJL, 747
- Terrien, R. C., Mahadevan, S., Deshpande, R., &amp; Bender, C. F. 2015, ApJS, 220, 16
- Thorngren, D. P., &amp; Fortney, J. J. 2018, AJ, 155, 214

-. 2019, ApJL, 874, L31

- Thorngren, D. P., Fortney, J. J., Murray-Clay, R. A., &amp; Lopez, E. D. 2016, ApJ, 831, 1
- Todorov, K. O., Deming, D., Burrows, A., &amp; Grillmair, C. J. 2014, ApJ, 796
- Torres, G., Fischer, D. A., Sozzetti, A., et al. 2012, ApJ, 757 Triaud, A. H., Gillon, M., Ehrenreich, D., et al. 2015, MNRAS, 450, 2279
- Seeliger, M., Kitze, M., Errmann, R., et al. 2015, MNRAS, 451, 4060
- Visscher, C., Lodders, K., &amp; Fegley, Jr., B. 2010, ApJ, 716, 1060 Wakeford, H. R., Sing, D. K., Kataria, T., et al. 2017, Science, 356, 628

## APPENDIX

We show the individual normalized raw light curves for each visit of HAT-P-15b, HAT-P-17b, HAT-P-18b, HAT-P26b, and WASP-69b in Figures A.1 - A.3. For each of the visits, we overplot the best-fit instrumental noise model derived from the joint fits. We then show these same individual light curves with the instrumental noise models removed in Figures A.4 - A.6 and overplot the best-fit eclipse model derived from the joint fits. We show the standard deviation of the residuals after removing both the eclipse and noise models from each visit for each planet in Figures A.7 - A.9. The predicted photon noise for each visit is overplotted for reference.

We show the mixing-ratio profiles for relevant species of interest for a 1 × and 100 × solar metallicity model for HAT-P-17b (Figures A.10 and A.11 respectively) and a 1 × and 100 × solar metallicity model for HAT-P-26b (Figures A.12 and A.13 respectively). We show both equilibrium and disequilibrium chemical models derived using the same framework as presented in Moses et al. (2011, 2013, 2016).

## FIGURES

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-3a99fd14-18cd-4118-bcc9-ffcd3d5fe012.jpg)

Time from Predicted Eclipse Center [days]

Figure A.1. Raw Spitzer photometry for each visit of HAT-P-26b. The normalized flux binned in ten-minute intervals is shown in black and the thirty-second binned flux is shown in gray. Overplotted is the best-fit instrumental model in red. Observations are shown in chronological order across each row.

Figure A.2. Raw Spitzer photometry for each visit of HAT-P-18b and WASP-69b. The normalized flux binned in ten-minute intervals is shown in black and the thirty-second binned flux is shown in gray. Overplotted is the best-fit instrumental model in red. Observations are shown in chronological order across each row.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-972bbc69-048d-4304-b1b7-136f3be7084a.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-c7ded8d7-3792-4120-b599-9dbcd0c6ce7e.jpg)

Time from Predicted Eclipse Center [days]

Figure A.3. Raw Spitzer photometry for each visit of HAT-P-15b and HAT-P-17b. The normalized flux binned in ten-minute intervals is shown in black and the thirty-second binned flux is shown in gray. Overplotted is the best-fit instrumental model in red. Observations are shown in chronological order across each row.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-c07e8ffe-3e95-4194-802e-64b13197aa64.jpg)

Time from Predicted Eclipse Center [days]

Figure A.4. Individual light curves for each visit of HAT-P-26b. The ten-minute binned normalized flux in shown in black with error bars showing the standard error of the flux in each bin. The instrumental best-fit parameters are unique to each visit and have been divided out. The red lines are the light curves with the best-fit parameters from the joint fits. Observations are shown in chronological order across each row. The 2 σ upper limits for the best-fit eclipse depths of the 3.6 µ m visits are shown.

Figure A.5. Individual light curves for each visit of HAT-P-18b and WASP-69b. The ten-minute binned normalized flux in shown in black with error bars showing the standard error of the flux in each bin. The instrumental best-fit parameters are unique to each visit and have been divided out. The red lines are the light curves with the best-fit parameters from the joint fits. Observations are shown in chronological order across each row.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-7a8fb7f4-bb52-4b61-9b1e-11eec0dd51d6.jpg)

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-3bbb777e-ebba-49b9-9184-486b0905acb8.jpg)

Time from Predicted Eclipse Center [days]

Figure A.6. Individual light curves for each visit of HAT-P-15b and HAT-P-17b. The ten-minute binned normalized flux in shown in black with error bars showing the standard error of the flux in each bin. The instrumental best-fit parameters are unique to each visit and have been divided out. The red lines are the light curves with the best-fit parameters from the joint fits. Observations are shown in chronological order across each row. The 2 σ upper limits for the best-fit eclipse depths of all visits of HAT-P-15b (see Section 2.2 for more details) and the 4.5 µ m visits of HAT-P-17b are shown.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-1f503eb4-bfef-4d8f-93c4-734d47810874.jpg)

Bins

Figure A.7. Standard deviation of the residuals of each visit of HAT-P-26b after removing the best-fit instrumental and astrophysical models as a function of bin size. The solid lines show the photon noise limits as a function of bin size scaled by 1/ √ N . Observations are shown in chronological order across each row.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-e5502935-1f8d-43ce-988b-fef6aa10e75f.jpg)

Bins

Figure A.8. Standard deviation of the residuals of each visit of HAT-P-18b and WASP-69b after removing the best-fit instrumental and astrophysical models as a function of bin size. The solid lines show the photon noise limits as a function of bin size scaled by 1/ √ N . Observations are shown in chronological order across each row.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-bd71c7f2-816c-44ea-b8e7-5723cbcfeed3.jpg)

Bins

Figure A.9. Standard deviation of the residuals of each visit of HAT-P-15b and HAT-P-17b after removing the best-fit instrumental and astrophysical models as a function of bin size. The solid lines show the photon noise limits as a function of bin size scaled by 1/ √ N . Observations are shown in chronological order across each row.

HAT-P-17b   C/O=1×solar   [Fe/H]=1×solar

Pressure (bars)

Figure A.10. Mixing ratio profiles derived using the same framework as presented in Moses et al. (2011, 2013, 2016) for relevant species of interest for a 1 × solar metallicity equilibrium (solid lines) and disequilibrium (dashed lines) model with solar C/O ratios for HAT-P-17b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-52be3276-0a21-451d-abe0-cbbf10eeb16f.jpg)

HAT-P-17b   C/O=1×solar   [Fe/H]=100×solar

Pressure (bars)

Figure A.11. Mixing ratio profiles derived using the same framework as presented in Moses et al. (2011, 2013, 2016) for relevant species of interest for a 100 × solar metallicity equilibrium (solid lines) and disequilibrium (dashed lines) model with solar C/O ratios for HAT-P-17b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-addd383a-15b4-4a4b-a1d6-f05835783ab3.jpg)

HAT-P-26b   C/O=1×solar   [Fe/H]=1×solar

Pressure (bars)

Figure A.12. Mixing ratio profiles derived using the same framework as presented in Moses et al. (2011, 2013, 2016) for relevant species of interest for a 1 × solar metallicity equilibrium (solid lines) and disequilibrium (dashed lines) model with solar C/O ratios for HAT-P-26b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-5ca4c740-4299-48de-9ec6-997ea0b7f142.jpg)

HAT-P-26b   C/O=1×solar   [Fe/H]=100×solar

Pressure (bars)

Figure A.13. Mixing ratio profiles derived using the same framework as presented in Moses et al. (2011, 2013, 2016) for relevant species of interest for a 100 × solar metallicity equilibrium (solid lines) and disequilibrium (dashed lines) model with solar C/O ratios for HAT-P-26b.

![image](https://ac.marcodsn.me/arxiv-markdown-images/1908.00014v1-figure-3cf9fe05-d5ef-4bf0-90a1-65bfd09a8825.jpg)
Nicole L. Wallack Heather A. Knutson Caroline V. Morley Julianne I. Moses Nancy H. Thomas Daniel P. Thorngren Drake Deming Jean-Michel Désert Jonathan J. Fortney Joshua A. Kammer
Investigating Trends in Atmospheric Compositions of Cool Gas Giant Planets Using Spitzer Secondary Eclipses
We present new 3.6 and 4.5 micron secondary eclipse measurements for five
cool (less than approximately 1000 K) transiting gas giant planets: HAT-P-15b,
HAT-P-17b, HAT-P-18b, HAT-P-26b, and WASP-69b. We detect eclipses in at least
one bandpass for all planets except HAT-P-15b. We confirm and refine the
orbital eccentricity of HAT-P-17b, which is also the only planet in our sample
with a known outer companion. We compare our measured eclipse depths in these
two bands, which are sensitive to the relative abundances of methane versus
carbon monoxide and carbon dioxide, respectively, to predictions from 1D
atmosphere models for each planet. For planets with hydrogen-dominated
atmospheres and equilibrium temperatures cooler than approximately 1000 K, this
ratio should vary as a function of both atmospheric metallicity and the
carbon-to-oxygen ratio. For HAT-P-26b, our observations are in good agreement
with the low atmospheric metallicity inferred from transmission spectroscopy.
We find that all four of the planets with detected eclipses are best matched by
models with relatively efficient circulation of energy to the nightside. We see
no evidence for a solar-system-like correlation between planet mass and
atmospheric metallicity, but instead identify a potential (1.9 sigma)
correlation between the inferred methane/(carbon monoxide + carbon dioxide)
ratio and stellar metallicity. Our ability to characterize this potential trend
is limited by the relatively large uncertainties in the stellar metallicity
values. Our observations provide a first look at the brightness of these
planets at wavelengths accessible to the James Webb Space Telescope, which will
be able to resolve individual methane, carbon monoxide, and carbon dioxide
bands and provide much stronger constraints on their atmospheric compositions.