Into the UV: The Atmosphere of the Hot Jupiter HAT-P-41b Revealed

We begin our exploration of HAT-P-41b's atmosphere by comparing our observed spectra with a grid of 1D forward models. The library of "generic" exoplanet transmission spectra of Goyal et al. (2019) spans a broad range of atmospheric compositions, temperature, and aerosol properties, under the assumption of equilibrium chemistry. Our fitting procedure follows a similar methodology to Wakeford et al. (2018), repeated for each WFC3 UVIS G280 data reduction. The resulting best-fitting models for the jitter detrending and marginalization reductions are shown in Figure 2 (top panel). Though the forward model grid can fit most of the observations, mismatches occur in several spectral regions, especially near 1.5 and 3.6 μm. As a result, the ${\chi }_{\nu }^{2}$ values for this equilibrium model were 2.24 and 2.60 for the jitter and marginalization cases, respectively.

Zoom In Zoom Out Reset image size

We find slight variations between the preferred forward models for each data reduction. The best-fitting model for the marginalization data reduction estimates the atmospheric temperature at the pressures probed via transmission as 1400 K, with a solar metallicity, solar ratio of carbon to oxygen (C/O, 0.56), and a significant cloud deck. In contrast, the best-fitting model for the jitter detrending reduction approach prefers a colder atmosphere (900 K), 10× solar metallicity, subsolar C/O (0.35), and a significant cloud deck. The fits are broadly consistent within the limits provided by their pseudo-probability distributions. ¹⁹ However, note that the best-fitting model in the marginalization case allows for the presence of VO, which can be present under chemical equilibrium assumptions at 1400 K at relevant atmospheric pressures, to explain the spectral structure shortward of 1 μm (Figure 2, top panel). This illustrates a degree of sensitivity to the chosen data reduction technique—we further quantify this in Section 3.3.

Our forward model comparison indicates that clouds could play a key role in shaping HAT-P-41b's transmission spectrum. However, our grid includes clouds via a simple cloud deck pressure decoupled from the atmospheric dynamics and thermochemical structure. We therefore verified the plausibility of the formation and advection of clouds in HAT-P-41b's atmosphere by running a general circulation model (GCM) using the SPARC/MITgcm (e.g., Showman et al. 2009; Kataria et al. 2016). We assume the same physical parameters of HAT-P-41b used throughout this study and an atmospheric composition of solar metallicity consistent with the preferred 1D atmospheric models. The resulting globally and spatially averaged temperature and vertical mixing profiles for the dayside, nightside, and each terminator are shown in Figure 3. The temperature profiles of the west (morning; green profiles) and east (evening; purple profiles) terminators cross multiple condensation curves, indicating that a broad range of cloud species may be present in HAT-P-41b's observable atmosphere. The ∼200 K difference between the terminators may drive differing cloud properties on each limb of the planet. With their plausibility and constraints on vertical mixing established, we turn to detailed microphysical modeling to investigate the physical nature and composition of clouds in the atmosphere of HAT-P-41b.

Zoom In Zoom Out Reset image size

We simulate cloud distributions for HAT-P-41b using the Community Aerosol and Radiation Model for Atmospheres (CARMA). CARMA is a one-dimensional bin-scheme aerosol microphysics model that computes vertical and size distributions of aerosol particles. CARMA solves a discretized continuity equation that accounts for aerosol nucleation, condensation, evaporation, and transport (see Gao et al. 2018 and references therein). The specific microphysical setup that we use for modeling condensational clouds in this work is described by Powell et al. (2019) and has been shown to reproduce trends in hot Jupiter cloudiness across a broad range of parameter space. In this work we do not vary the microphysical parameters of the condensate species and we assume that the volatile species in the atmosphere have a solar abundance. We post-process these results to calculate transmission spectra using a modified version of Exo-Transmit that includes a Mie theory prescription for the cloud opacities as described by Powell et al. (2019).

To model the observations of HAT-P-41b we must choose input planetary properties. As the atmospheric metallicity and microphysical parameters are held constant, the remaining tunable parameters are the temperature profile and the amount of atmospheric vertical mixing. For both parameters we begin by considering the globally averaged 3D GCM temperature and vertical mixing profile shown in Figure 3. We then vary both profiles by a constant factor, because both parameters are not well constrained, to derive cloud properties that give rise to simulated spectra that best match the observations.

The best-fitting model, shown in Figure 2 (bottom panel), is 400 K hotter than the globally averaged temperature profile from the 3D GCM and has two orders of magnitude less global mixing. The increase in temperature limits the supersaturation of the condensible species and thus the formation of clouds. This increase in temperature is consistent with temperature variations seen in the GCM, especially limb-to-limb variations (e.g., Kataria et al. 2016). Reducing the amount of vertical mixing further limits the formation of clouds as well as the vertical extent in the atmosphere of the cloud particles (see Powell et al. 2019). The level of reduction in the vertical mixing calculated from the GCM is consistent with studies of GCM tracer transport, where the derived mixing of tracers is commonly two orders of magnitude less efficient than the transport of atmospheric gases (e.g., Parmentier et al. 2013). Both of these effects give rise to a population of large aluminum (Al₂O₃) clouds at ∼10⁻¹–10⁻³ bar that dominate the cloud opacity. These clouds are optically thick at short wavelengths and optically thin at longer wavelengths (Vahidinia et al. 2014). The results of this one-dimensional globally averaged cloud modeling lead to a similar spectral fit to the first-order forward models considered in Section 3.1 (e.g., the jitter fit in the top panel of Figure 2). However, the inclusion of wavelength-dependent cloud opacities improves the fit at 3.6 μm. Nevertheless, this cloud model fails to reproduce the full shape of the H₂O feature centered at 1.4 μm. The corresponding ${\chi }_{\nu }^{2}$ values for this cloud model were 2.67 and 2.85 for the jitter and marginalization cases, respectively. As an additional sanity check, we ran simple Mie theory model comparisons (Wakeford & Sing 2015), which required large (∼0.5–10 μm) particles that result in gray opacities out to 2 μm and miss key molecular absorption features. In seeking a model capable of explaining the observations over the full wavelength range, we turn now to retrieval analyses.

Atmospheric retrievals relax many aforementioned assumptions, such as chemical equilibrium, opting instead to parameterize the atmospheric state. Bayesian sampling techniques explore millions of potential states, comparing their resultant transmission spectra with observations to derive posterior probability distributions for each model parameter. This inverse approach allows atmospheric properties, such as chemical abundances, temperature profiles, and cloud properties, to be retrieved directly from the data.

Our retrieval philosophy employs two central principles to ensure robust atmospheric inferences. First, we conduct retrievals for each data reduction separately, establishing the sensitivity of derived atmospheric properties to different reduction techniques. Second, three different retrieval codes are independently applied to each data reduction: POSEIDON (MacDonald & Madhusudhan 2017), NEMESIS (Irwin et al. 2008; Barstow et al. 2017; Krissansen-Totton et al. 2018), and ATMO (Amundsen et al. 2014; Tremblin et al. 2015; Wakeford et al. 2017). We can thereby quantify the impact on our results of differing modeling choices (e.g., molecular line lists, temperature profiles, cloud parameterizations).

We consider an extensive range of potential atmospheric components. Our investigations include opacity due to the following chemical species: H₂, He, H⁻, Na, K, Li, TiO, VO, AlO, SiO, TiH, CrH, FeH, AlH, CaH, SiH, H₂O, CH₄, CO, CO₂, NH₃, HCN, NO, H₂S, SH, PH₃, and C₂H₂. Isothermal (NEMESIS and ATMO) and non-isothermal temperature structures (POSEIDON, Madhusudhan & Seager 2009) were considered. Three cloud models were examined: (i) an opaque cloud deck with a vertically uniform haze (ATMO, Wakeford et al. 2017); (ii) a single cloud deck, with variable top and base pressures, and power-law extinction with a variable index (NEMESIS, Barstow et al. 2017); and (iii) patchy clouds and haze around the terminator (POSEIDON, MacDonald & Madhusudhan 2017). Although these cloud prescriptions do not explicitly include Mie scattering calculations, our model fitting exercise with CARMA (see Section 3.2) demonstrates that more complex microphysical prescriptions for clouds do not fully explain the observed spectrum of HAT-P-41b. The influence of stellar contamination was also considered in the retrieval process (NEMESIS, e.g., Pinhas et al. 2018). Iterative expansion of the considered molecules was performed between the three retrieval codes until a minimal basis set was identified. All three codes use nested sampling to explore the parameter space, via either MultiNest (NEMESIS and POSEIDON, Feroz & Hobson 2008; Feroz et al. 2009, 2019; Buchner et al. 2014) or dynesty (ATMO, Speagle 2020).

Our retrieved transmission spectra are compared to the observations of HAT-P-41b in Figure 4. Contrasting with the forward models of the previous sections, the retrievals prefer a combination of gas-phase optical opacity sources instead of clouds. Specifically, at least one of H⁻, AlO, CrH, and VO, ²⁰ in addition to Na, is required to explain the WFC3 UVIS G280 observations. By invoking chemical species with strong near-UV to visible opacities, but weak infrared opacities, the infrared observations can be well fit by H₂O alone. In particular, there is a clear preference for the bound–free opacity of the hydrogen anion, H⁻, which provides a smooth continuum across the UVIS range before falling off rapidly as the ionization threshold (∼1.64 μm) is approached. This continuum has similar spectral characteristics to a cloud deck across the UVIS range, potentially explaining the preference for clouds in our forward models (which do not include H⁻). The other inferred opacity sources are somewhat sensitive to the retrieval code ²¹ and chosen data reduction, as we demonstrate below. Nevertheless, all three retrieval codes agree on the interpretation of a strong near-UV to visible chemical opacity source, without the need for clouds in the observable atmosphere.

Zoom In Zoom Out Reset image size

We detect the presence of H₂O in HAT-P-41b's atmosphere at ≳5σ confidence (jitter: 6.2σ; marginalization: 5.1σ). The visible and near-UV observations additionally require at least one other prominent opacity source, with candidates identified as H⁻, CrH, AlO, VO, and Na. Besides these species, our initial retrievals—ranging in complexity from 12 to 37 free parameters—included many parameters left largely unconstrained by the present observations (e.g., temperature structure and cloud properties). Consequently, the best-fitting ${\chi }_{\nu }^{2}$ ranged from ∼1.7 to 3.0 (see Appendix B). We therefore constructed an eight-parameter "minimal" model, including only those parameters found necessary (via the Bayesian evidence) to explain the observations: an isothermal, clear atmosphere with H₂O, H⁻, Na, CrH, AlO, and VO. This best-fitting model attains ${\chi }_{\nu }^{2}$ = 1.50 and 1.72 for the jitter and marginalization cases, respectively.

The ${\chi }_{\nu }^{2}$ values obtained by our best-fitting retrieval models demonstrate a greatly improved quality of fit compared to the chemical equilibrium and self-consistent cloud models considered in Sections 3.1 and 3.2. However, we note that our chi-squared values still suggest some level of tension between the data and models. Under frequentist metrics one could still consider rejecting all the models presented here, but this argument assumes that both the data and models perfectly capture all noise sources and atmospheric physics. Tactics such as error-bar inflation (as done in studies such as Line et al. 2015; Zhang et al. 2018; Colón et al. 2020) and incorporating models of additional complexity (e.g., see discussion in Gibson 2014; Parviainen 2018) could be used to further reduce these chi-squared values, but could obfuscate our physical interpretation for HAT-P-41b from these data and associated intermodel comparisons. It is important to note that both H₂O and UV–optical opacity sources are independently required to explain the observations, irrespective of chi-squared tests, as established by our Bayesian model comparisons.

With reference to our best-fitting "minimal" model, we conducted a series of Bayesian model comparisons with POSEIDON to compute detection significances for each inferred UV–optical chemical species. The jitter reduction yields moderate evidence for H⁻ (2.9σ), weak evidence for CrH (2.7σ) and AlO (2.4σ), and a tentative hint of Na (2.0σ). The marginalization reduction yields moderate evidence for Na (2.9σ), weak evidence for H⁻ (2.6σ) and CrH (2.5σ), and tentative hints of AlO (1.9σ) and VO (1.7σ). The specific spectral features giving rise to these inferences are shown in Appendix C (see Figure 7). The differing significances for Na and VO highlight the sensitivity of some atmospheric inferences to specific data reductions. However, our most rigorous conclusion holds for both reductions: at least one of H⁻, CrH, AlO, and VO is required at >5σ (jitter: 5.1σ; marginalization: 5.5σ) to explain HAT-P-41b's transmission spectrum.

The retrieved abundances for each inferred chemical species ²² are shown in Figure 5. All three codes provide precise H₂O abundances. Across both data reductions and all three retrieval codes, the retrieved H₂O abundance spans $\mathrm{log}({{\rm{H}}}_{2}{\rm{O}})$ ≈ −3.4 to −1.6 (with a mean precision of 0.5 dex). The H₂O abundances from ATMO are ≳0.4 dex higher than those retrieved by NEMESIS and POSEIDON. Nevertheless, the retrieved abundances are consistent to 1σ—see Table 1. These results illustrate that the data reduction technique chosen at visible wavelengths can alter retrieved H₂O abundances by ≲0.5 dex. We additionally report a precise H⁻ abundance, consistent across all three retrieval codes and both data reductions, with a mean value of $\mathrm{log}({{\rm{H}}}^{-})=-8.65\pm 0.62$. Abundance constraints for the other inferred species are relatively weak, and can vary significantly between different retrieval codes and data reductions (e.g., the VO posteriors in Figure 5). However, the derived H₂O and H⁻ abundances are robust to both modeling choices and data reduction techniques. We note that these abundances are with respect to the original "full" retrievals (see Table 1), ensuring that model uncertainties (overlapping absorption features, cloud-chemistry degeneracies, etc.) are automatically encoded in the quoted values. However, biases may still arise from neglected model complexity. In particular, retrieved H₂O and H⁻ abundances from 1D retrieval methods can be biased by up to 1 dex if compositional gradients arise between the morning and evening terminators (MacDonald et al. 2020). Such a bias would cause a slight overestimate in our retrieved H₂O abundances, and an underestimate in our H⁻ abundances (see MacDonald et al. 2020, Figure 3).

Zoom In Zoom Out Reset image size

The precise H₂O abundances we obtain can be converted into estimates of the atmospheric metallicity. By "metallicity" we refer to the atmospheric O/H ratio relative to that of its star ([Fe/H] = 0.21, Stassun et al. 2017). The retrieved molecular abundances are mapped into O/H ratios as in MacDonald & Madhusudhan (2019). We note that our metallicities should only be considered accurate to a factor of two, ²³ given that the current data are insensitive to CO. The metallicities derived by POSEIDON and NEMESIS are consistent with the stellar metallicity of HAT-P-41 for the marginalization reduction: ${1.73}_{-1.24}^{+6.26}$ × stellar and ${1.22}_{-0.92}^{+3.80}$ × stellar, respectively. The ATMO retrievals find ∼4× higher metallicities (due to the aforementioned higher H₂O abundances): ${4.35}_{-2.99}^{+10.43}$ × stellar for the marginalization reduction. Comparatively, the jitter reduction favors slightly superstellar metallicities: ${3.65}_{-2.25}^{+5.59}$ × stellar (POSEIDON), ${1.79}_{-1.33}^{+4.66}$ × stellar (NEMESIS), and ${9.57}_{-5.28}^{+8.43}$ × stellar (ATMO). Nevertheless, all derived metallicities are consistent with the stellar value to 2σ. Overall, the metallicities derived by the different retrieval codes are consistent with each other. We conclude that the metallicity of HAT-P-41b's atmosphere is consistent with being stellar, or slightly superstellar, in line with expected mass–metallicity trends.

Our retrieved abundances for visible to near-UV absorbers require an atmosphere in chemical disequilibrium. Here we consider disequilibrium mechanisms that might enhance the atmospheric abundance of AlO, CrH, and H⁻ in HAT-P-41b's atmosphere, the species common to both the jitter and marginalization data reduction analyses (see discussion in Section 3.4). At the retrieved terminator temperature ²⁴ (∼1000 ± 200 K), the equilibrium abundances of CrH and AlO for an atmosphere of solar metallicity are respectively ∼7 and 5 orders of magnitude below those we infer (Woitke et al. 2018), which represents a significant discrepancy between the retrieved abundances and those predicted by equilibrium chemistry as well as with the atmospheric metallicity we infer from the retrieved H₂O abundance. Even at depth in HAT-P-41b's atmosphere (P ∼ 1 bar, T ∼ 2000 K) the expected abundances for CrH and AlO are several orders of magnitude less than the retrieved abundances (Woitke et al. 2018), which makes it difficult for a mechanism such as vertical quenching to account for their disequilibrium abundances. Given this difficulty in identifying a plausible mechanism for a significant enhancement of CrH or AlO in HAT-P-41b's atmosphere we turn our attention to H⁻, which had the highest statistical significance (2.6–2.9σ) of the UV–optical absorbing species we considered in our retrievals. Under equilibrium thermochemistry, H⁻ is not expected to exist in appreciable quantities (≥10⁻⁹ mole fraction) for T ≲ 2500 K (Kitzmann et al. 2018). We verified the need for H⁻ opacity by conducting an additional ATMO retrieval with chemical equilibrium enforced. This retrieval raised the atmospheric temperature to ≳2500 K, far above T_eq and predictions from the GCM (Figure 3), to create H⁻ opacity with an abundance consistent with that found by our free retrievals (see the supplementary material).

We explored the potential for photochemistry to enhance the abundance of H⁻ in a hot Jupiter with T ≲ 2500 K such as HAT-P-41b. In the pressure regions probed by transits (∼1 mbar), Lavvas et al. (2014) show that electron densities will be roughly an order of magnitude higher than predictions from thermal ionization only. It is also expected that photochemistry will enhance the neutral H abundance by many orders of magnitude over than those predicted by equilibrium chemistry (e.g., Liang et al. 2003; Moses et al. 2011). This expected increase in electron densities and neutral H abundance in the atmosphere of a hot Jupiter such as HAT-P-41b increases the chance for production of H⁻ via radiative electron attachment (${\rm{H}}+{{\rm{e}}}^{-}\to {{\rm{H}}}^{-}\ +\ \gamma$), collisional electron attachment (${\rm{H}}+{{\rm{e}}}^{-}+{\rm{M}}\to {{\rm{H}}}^{-}+{\rm{M}}$), and dissociative electron attachment (${{\rm{H}}}_{2}+{{\rm{e}}}^{-}\to {\rm{H}}\ +{{\rm{H}}}^{-}$). For conditions relevant to HAT-P-41b, the last process likely dominates, due to the prevalence of H₂ and the moderately large rate coefficient of the order of 10⁻¹³ cm³ s⁻¹ (Janev et al. 2003). The destruction of H⁻ can occur through collisional detachment (e.g., ${{\rm{H}}}^{-}+{\rm{H}}\to {{\rm{H}}}_{2}+{{\rm{e}}}^{-}$ or ${{\rm{H}}}^{-}+{{\rm{H}}}_{2}{\rm{O}}\to {\mathrm{OH}}^{-}+{{\rm{H}}}_{2}$) with rate coefficients for such reactions of the order of a few × 10⁻⁹ cm³ s⁻¹ (Bruhns et al. 2010; Martinez et al. 2010). From the atmospheric density of 0.1 mbar from the global-average GCM results described above and the mixing ratios of e⁻, H, and H₂ from the model of Lavvas et al. (2014), we can estimate number densities for e⁻, H, and H₂ near the ∼0.1 mbar level of HAT-P-41b's atmosphere to be ∼10⁹, 2 × 10¹³, and 6 × 10¹⁴ cm⁻³ respectively. Assuming that in steady state the H⁻ production and loss rates balance each other, and considering production solely by H₂ dissociative electron attachment and destruction solely by collisional detachment with atomic H, we estimate a number density for H⁻ of the order of 10⁶ cm⁻³. This value corresponds to an H⁻ mixing ratio of 2 × 10⁻⁹, consistent with our retrieved value for the abundance of H⁻ and roughly six orders of magnitude larger than expectations from equilibrium chemistry for a planet such as HAT-P-41b. This order-of-magnitude estimate highlights that significant enhancement of H⁻ due to photochemical processes is likely present in HAT-P-41b and many other exoplanet atmospheres, thus shaping their UV–optical spectra. The production of H⁻ will be particularly enhanced for hot planets that receive a high extreme-UV flux from their host stars and have Na in the gas phase, which increases electron production (Lavvas et al. 2014). Such conditions are expected for HAT-P-41b (this work; Hartman et al. 2012; Linsky et al. 2014) and other hot Jupiters orbiting F-stars.