An Exploration of Model Degeneracies with a Unified Phase Curve Retrieval Analysis: The Light and Dark Sides of WASP-43 b

Real observations of WASP-43 b were used to test and illustrate the phase curve model presented in this paper. They were obtained from Stevenson et al. (2014, 2017) with no modifications. These consist of 15 already reduced spectra of the phase curve of WASP-43 b: 0.0625, 0.125, 0.1875, 0.25, 0.3125, 0.375, 0.4375, 0.5, 0.5625, 0.625, 0.6875, 0.75, 0.8125, 0.875, and 0.9375. The data were originally obtained by the Hubble Space Telescope in November 2013 (Program GO-13467) during three full orbital phases, three transits, and two eclipses with the WFC3 G141 grism (wavelength coverage from 1.1 to 1.6 μm). These consisted of spatially scanned images corresponding to 13–14 HST orbits for the phases and 4 HST orbits for the transit and eclipses. Complementary phase curve observations were obtained by the Spitzer Space Telescope, using the 3.6 μm (two visits, including one discarded) and 4.5 μm (one visit) photometric channels (Programs 10169 and 11001, PI: Kevin Stevenson). In addition to this, we also include the transmission spectrum from Kreidberg et al. (2014; in our model, this corresponds to phase 0.0) as this provides additional constraints on the day–night limb and allows us to extract the planetary radius with greater accuracy.

In emission, a large wavelength coverage is usually required to ensure that a sufficient pressure range of the atmosphere is probed, thus constraining the temperature structure (see Appendix D). Normally, this is done by combining the observations from HST and Spitzer to obtain spectra spanning 1.1–4.5 μm. The combination of instruments, however, can bring additional difficulties as nothing guarantees the compatibility of the observations. Sources of discrepancies can come from different instrument systematics, the use of different orbital elements, different reduction pipelines, and stellar or planet temporal variations (Yip et al. 2020, 2021; Changeat et al. 2020). As a consequence, independent studies of HST WFC3 data often lead to similar spectral shapes but different absolute depths (Changeat et al. 2020). Studies of the WASP-43 b Spitzer data (Stevenson et al. 2017; Mendonça et al. 2018a; Morello et al. 2019; Bell et al. 2021; May & Stevenson 2020) have also proved that independent reduction pipelines can obtain different results for the same Spitzer data set. In this paper, we use the Spitzer data from Stevenson et al. (2017) as is and do not investigate the potential implications of these effects. We however present a complementary retrieval in the Discussion section with the HST-only data, which highlights how the Spitzer points affect our solution. To date, this WASP-43 b data set is one of the most complete. Our integrated phase curve framework accumulates the information at all phases to extract extremely precise constraints on the atmospheric properties of this planet.

For this study, we assumed that the planet was mainly composed of hydrogen and helium with a ratio He/H₂ = 0.17. We considered collision-induced absorption of the H₂–H₂ (Abel et al. 2011; Fletcher et al. 2018) and H₂–He (Abel et al. 2012) pairs and included opacities induced by Rayleigh scattering (Cox 2015) and clouds. Because clouds are most constrained by the transmission spectrum, we model them with a fully opaque cloud layer above a given pressure, restricted to the day- and night-side regions. For the chemistry, we use the molecular line lists from the Exomol project (Tennyson et al. 2016, 2020; Chubb et al. 2021), HITEMP (Rothman & Gordon 2014), and HITRAN (Gordon et al. 2016). While many molecules are considered in the chemical equilibrium scheme, we only include molecular opacities for H₂O (Barton et al. 2017; Polyansky et al. 2018), CH₄ (Hill et al. 2013; Yurchenko & Tennyson 2014), CO (Li et al. 2015), CO₂ (Rothman et al. 2010), C₂H₂ (Wilzewski et al. 2016), C₂H₄ (Mant et al. 2018), NH₃ (Yurchenko et al. 2011), and HCN (Harris et al. 2006; Barber et al. 2013).

Section 4 aims to validate our method by presenting a retrieval on mock-up data where the true solution is known in advance. This approach allows us to check that our model is correctly implemented and will ensure that the results presented on real data, later on, are not artifacts of our method. This example is performed using Ariel simulated spectra, which also gives us the opportunity to investigate the Ariel capabilities to perform phase curve observations. The planet is simulated using the stellar and planet parameters from the literature (Bonomo et al. 2017). We first create the forward model at high resolution using our phase curve model with three distinct regions: hot spot, day side, and night side. Each region is composed of 100 layers spaced in log pressures. We assume homogeneous temperature profiles and a chemical compositions for the three regions inspired by previous works (Stevenson et al. 2014, 2017; Changeat & Al-Refaie 2020). For simplicity, the chemistry is modeled constant with altitude, and only H₂O, CH₄ and CO are included. In this example, the H₂O abundance is set at 6 × 10⁻³ for the hot spot and 1 × 10⁻⁴ for the day and night sides. For CH₄, the atmosphere contains 1 × 10⁻⁷, 1 × 10⁻⁵, and 1 × 10⁻⁴ for, respectively, the hot spot, the day side, and the night side. Finally, the CO abundance is set to 1 × 10⁻³ on the hot spot and 1 × 10⁻⁴ for the other regions. For this test, the hot-spot parameters were fixed to a hot-spot offset of 1202, consistent with Stevenson et al. (2014), and a hot-spot size of 40°, matching findings from Kataria et al. (2015) and Irwin et al. (2020). This leads to 15 spectra from phase 0.0625 to phase 0.9375. We then bin down those spectra to Ariel observations assuming Tier 2 resolution and the observation of four complete planet revolutions. The noise was simulated using the Ariel Radiometric Model (ArielRad) from L. Mugnai et al. (2019, in preparation). We do not use random noise instances of our simulated observations for the retrievals, as we aim to study the retrieval biases arising from our retrieval technique (Feng et al. 2016; Changeat et al. 2019, 2020; Mai & Line 2019). Finally, we perform our phase curve retrieval test, leaving free the planetary radius, the temperature profiles, the constant chemical abundances, and the hot-spot parameters. The parameter space is explored uniformly with large bounds (see Table 1). The results of this exercise are described in Section 4.

Section 5 presents the results for the retrievals performed on the real data. Here, we consider the reduced observations described in Section 3.1. We also use the star parameters from Bonomo et al. (2017) and fix the planet mass (Changeat et al. 2020) to its radial velocity measurement (Bonomo et al. 2017). Because we always include the transmission spectrum in our retrievals, which is very sensitive to the planetary radius, we let this parameter as free. We note that complementary tests without including the transmission spectrum did not affect the results presented here. We performed an analysis of each spectrum individually (see Appendix C) and carried on with unified analyses using three classes of models of increasing complexity:

- Two-faces free: We considered a simplified geometry with no hot spot. We remove the influence from the hot-spot region by coupling all parameters (temperature, chemistry, clouds) between the hot-spot and the day-side regions. This allows us to check the consistency of our model and compare it with previous results from the literature. This setup uses a similar geometry to Feng et al. (2020). We considered free constant with altitude abundances for the molecules H₂O, CH₄, CO, CO₂, and NH₃. In the main scenario, these were coupled between all three regions. The free temperature profiles and cloud properties were left independent between the day and night regions. For the clouds, we use a simplistic fully opaque gray cloud. Two complementary runs were also used to explore the effect of the Guillot (2010) T–p profile and decoupled chemistry (different chemistry between the night and day sides).

(1) Two-faces equilibrium: For the second scenario, we used the same geometry but increased the complexity of our model by considering equilibrium chemistry for the day- and night-side regions. We used the scheme ACE from Agúndez et al. (2014). It calculates thermochemical equilibrium abundances for H-, He-, C-, O-, and N-bearing species, which is expected for hot Jupiters between 1000 and 2000 K. Only a fraction of the calculated molecules possess cross sections, so we limit the actively absorbing species to the ones described in Section 3.2. The chemistry was decoupled between the day and night regions (different chemical profiles) but the free parameters (metallicity and C/O ratio) are shared between all regions.

(2) Full: Finally, we performed simulations on the full model by introducing the hot-spot region described in Section 2. The definition of this region requires the additional hot-spot offset (Δ) and hot-spot size (α) parameters. In the first place, we attempted to retrieve these parameters but encountered inconsistencies with previous values from the literature (see posterior distribution in Appendix E). Because these degrees of freedom were already included in the reduction steps leading to the spectra (Stevenson et al. 2014, 2017), we decided to fix the hot-spot shift to the value in Stevenson et al. (2014): −122. For the hot-spot size, we tested the values 30°, 40°, and 50°, consistent with standard predictions from recent models for this planet (Kataria et al. 2015; Irwin et al. 2020) but only present the model with 40° in the result section. The complementary retrievals are presented in the discussion section and explore the impact of hot-spot size, hot-spot clouds, and the addition of the Spitzer data. We caution the fact that these apparent inconsistencies might be linked to incompatibilities of the spectra following the reduction process (this was not highlighted by our individual analysis in Appendix C), issues with our model assumptions (different planet geometry), or the lack of information content in the considered spectra.

In all scenarios, except when stated otherwise, the parameterization of the temperature profiles was done using the N-point profile from TauREx3 (Al-Refaie et al. 2019). This heuristic profile interpolates linearly between N freely moving temperature–pressure points. The profile is then smoothed over 10 layers to avoid inflection points. For the hot-spot and day regions, we retrieve seven temperature points at fixed pressures (10⁶, 10⁵, 10⁴, 10³, 10², 1, and 0.01 Pa). Because the information content decreases on the night side due to the lower emission, we found that retrieving five points (10⁶, 10⁵,10³, 10, and 0.1 Pa) was more suitable. These choices are investigated in Appendix D, where various free temperature structures are explored. This level of complexity was chosen to maximize the Bayesian evidence in the retrievals of the WASP-43 b data from HST+Spitzer with the new phase curve model. We highlight, however, that other profiles lead to equivalent Bayesian evidence in the HST-only case.

The exploration of the parameter space is performed with uniform priors using the nested sampling algorithm MultiNest (Feroz et al. 2009; Buchner et al. 2014) with 500 live points and a log likelihood tolerance of 0.5. These choices ensure an optimal free sampling of the parameter space. All of the free parameters considered in this retrieval and their uniform priors are described in Table 1. The results of our investigation on the real data from WASP-43 b observations are presented in Section 5.

For an easier comparison of our results, we provide the model integrated values in the Spitzer bands and the χ² values corresponding to each spectrum. We note that χ² values on individual spectra might be difficult to interpret for model comparison in a Bayesian context and should be viewed with caution (Edwards et al. 1963). In particular, the data set and parameters considered in this work are not expected to be independent, while the tested models are also nonlinear (Andrae et al. 2010; Gelman 2013). The stated χ² value is calculated for the best-fit model, which might not be representative of the larger pool of solutions found via our Bayesian inference technique. Within a Bayesian framework, the Bayesian evidence, log(E), is used for model comparison and accounts for both goodness of fit and model complexity (Rougier & Priebe 2020). When comparing two models, a log difference of 3 indicates a strong preference for the model with the highest evidence (Kass & Raftery 1995; Jeffreys 1998). When testing models of increasing complexity, their Bayesian evidence should also increase as long as the additional flexibility is justified. For models of adequate complexity, the Bayesian evidence should plateau around a maximum value. Finally, if the model overfits, the Bayesian evidence should decrease.

In this run (log(E) = 2238.1), we consider only two regions: day and night sides. The geometry and the chemistry model (free chemistry) are similar to what has recently been presented in Feng et al. (2020), which allows us an easier comparison with their findings. A major difference in our model is the parameterization of the temperature structure. Here, we choose a free temperature profile, and we retrieve the temperature structure from the day and night side independently. Indeed, alternative physics-based descriptions, such as the parameterization from Guillot (2010), assume a one-dimensional atmosphere experiencing irradiation from two sources (upward and downward fluxes). While this describes well a planet day side with poor atmospheric redistribution, this geometry does not accurately represent the night side of a tidally locked planet, which are by definition not receiving direct stellar flux. For all models in this paper, except when stated otherwise, we therefore opted for a free parameterization. The best-fit spectra, the geometry, and the posteriors for our two-faces free chemistry model are described in Appendix G (blue runs). These highlight the difficulty of the model to fit the Spitzer points around the secondary eclipse. This is most likely due to the symmetry enforced in this type of geometry. The temperature profile for the day side mainly decreases with altitude, which is consistent with previous studies (Stevenson et al. 2014, 2017; Kataria et al. 2015; Irwin et al. 2020; Feng et al. 2020), but we note the preference (large 1σ uncertainties) for a thermal inversion above 10³ Pa. We ran complementary retrievals, parameterizing the temperature profile with the prescription from Guillot (2010), and found a noninverted temperature profile and chemistry (orange runs in Appendix G), similar to the findings in Feng et al. (2020). This run, however, led to a particularly low Bayesian evidence of log(E) = 2045.2, potentially highlighting the lack of flexibility in the Guillot profile. We discuss the presence of thermal inversions for this planet later, in the section describing the results of the Full model scenario. On the night side, the temperature is poorly constrained due to the low flux received at these phases. We note that the model prefers high-altitude clouds. In terms of chemical species, this retrieval finds very constrained abundances for water and ammonia (see posteriors). These match the abundances found by the joint scenario in Feng et al. (2020). As opposed to their results, the free T–p profile run does not recover constraints on CO₂, which can be explained by our use of the more flexible temperature structure. Indeed, when the Guillot profile was used, we recovered a similar abundance for CO₂ (see orange posteriors in Appendix G). This retrieval exploration using the two-faces free chemistry model confirms the findings presented in Stevenson et al. (2017), Irwin et al. (2020), and Feng et al. (2020), and also highlights different solutions when using more flexible temperature structures (also see Appendix D). We also ran a complementary retrieval where the chemistry is decoupled between the day and night sides. This retrieval led to two solutions corresponding to the green (log(E) = 2242.3) and gray (log(E) = 2240.4) runs in Appendix G. Because the gray run corresponds to an unphysically high mean molecular weight solution, we focus on the green solution. Due to the large day/night temperature contrast, differences in the chemistry can be expected. This retrieval presents a very different chemistry from the previous models with a much higher water content on the day side (around 10⁻²) and a much more complex temperature structure. The chemistry recovered here is similar to the ones presented in the next sections using decoupled equilibrium chemistry models (two-faces equilibrium and full models). On the night side, we do not recover any molecules.

One major assumption, taken in the previous scenario, is the constant chemistry with altitude and longitude. Previous studies showed that atmospheric chemistry is not constant and that assuming so could lead to observable biases (Venot et al. 2012, 2020; Agúndez et al. 2014; Drummond et al. 2018; Stock et al. 2018; Woitke et al. 2018; Changeat et al. 2020; Drummond et al. 2020). In order to provide a more realistic description of chemical properties, one could either assume more complex parameterizations (Parmentier & Crossfield 2018; Changeat et al. 2019) or use self-consistent chemical models (Agúndez et al. 2014; Stock et al. 2018; Woitke et al. 2018). While both approaches are viable with high-quality data, the low wavelength coverage of HST favors the more constrained models. We therefore assumed equilibrium chemistry for all three regions using the scheme from Agúndez et al. (2014). Because we assume that atomic elements are evenly spread across the planet atmosphere, this leaves us with only two free parameters: metallicity and C/O ratio. The impacts of these assumptions are explored more in the Discussion section. We present the best-fit spectra, chemical profiles, temperature structure, and posteriors in Appendix H. As compared with the previous run, the use of equilibrium chemistry allows us to better fit the Spitzer photometric points (log(E) = 2245.8). In this scenario, the retrieved chemistry is well constrained (see posterior distribution). The atmosphere is consistent with a slightly supersolar metallicity (${0.76}_{-0.12}^{+0.1}$) and a low C/O ratio (${0.25}_{-0.06}^{+0.07}$). In particular, the use of the chemical scheme allows large variations in the abundances of carbon-bearing species between the day and night side. There are linked to the large retrieved temperature differences between the two regions (from 2500 K down to 500 K) and are matching the findings in our free chemistry run with decoupled chemistry (green run in Appendix G). As compared to the two-faces free run, an interesting difference is the retrieved abundance for water. Here, we find that the abundance for water remains constant between the two regions (as expected from thermochemical equilibrium) but that it might be two orders of magnitude higher than the values found in the free chemistry case (log(H₂O) = −4.4). In fact, we find that a very similar solution (with high water abundance) can also be recovered from a free chemistry retrieval when the abundances between the day and night-side regions are decoupled, clearly indicating the need to account for chemical variations between those two regions. In terms of temperature structure, we find a day-side temperature profile that mainly decreases with altitude and the presence of a thermosphere, similar to the free chemistry run. The night side is much colder and includes high-altitude clouds at pressures as high as 10 Pa.

The previous models have difficulties fitting the spectra for the phases near eclipse. This is due to the symmetry imposed when considering only two faces. To increase the complexity, we split the geometry into three distinct regions: hot spot, day side, and night side. For the run presented here, the model has a fixed hot-spot size of 40°. We also model the planet without including clouds on the hot spot (the cloud-top pressure was fixed at 10⁶ Pa). This choice is justified by theoretical predictions, which suggest that the hot day side of irradiated exoplanets might be cleared due to the strong stellar irradiation (Lee et al. 2016; Parmentier et al. 2016; Lines et al. 2018; Helling 2019). The best-fit spectra for this scenario and the corresponding geometries are shown in Figure 2. For this retrieval, we obtain a log evidence of 2277.4, which indicates the relevance of adding this additional region to the model.

Zoom In Zoom Out Reset image size

We also investigated the same run on the HST spectra only. Both HST+Spitzer and HST-only runs are presented in Appendix I, which provides a zoomed version of the same run around the HST wavelength range. The HST+Spitzer run is colored while the HST-only run is gray. Both fits (with and without the Spitzer photometric points) do not show major differences in the HST wavelength range and are well fitted to the observed data. Outside the HST wavelength coverage, however, we see that the inclusion of the Spitzer points leads to large differences in the retrievals.

The temperature structure associated with the HST+Spitzer retrieval is presented in Figure 3, which as expected, clearly displays a hotter day side and a more strongly irradiated hot spot.

Zoom In Zoom Out Reset image size

As in the free temperature two-faces runs presented in previous sections, the temperature structure is consistent with the presence of thermal inversions for this planet. This contrasts with previous findings from Blecic et al. (2014), Stevenson et al. (2014), and Kreidberg et al. (2014). The retrieved temperature structure indicates the presence of a stratosphere (above 10⁵ Pa) and the presence of a thermosphere (above 10 Pa) on the planet's day side and hot spot. We note that the current data do not allow us to probe the highest altitude (thermosphere) with great accuracy, leading to large 1σ uncertainties on the retrieved profile (around ±600 K). This is also shown in the contribution functions for each region, which are plotted in Appendix J. In the full model and for most of the temperature structures explored in Appendix D, the hot spot and day side present a thermal inversion at high altitudes, below 1 Pa. This drives the thermal dissociation of all molecules at those altitudes, making the atmosphere transparent. We discuss further the interplay between this high-altitude thermal inversion and our chemical model assumptions in the discussion section. For this planet, thermal inversions are predicted (Showman et al. 2009; Kataria et al. 2015) in the presence of optical absorbers such as atomic species (K, Na) or metal hydrides and oxides (e.g., AlO, FeH, SiO, TiO, TiH, or VO). The presence of aluminum oxide (AlO) in the transmission spectrum of WASP-43 b has recently been highlighted in Chubb et al. (2020). This detection, made at the terminator region, is not expected from equilibrium chemistry models at the temperatures and pressures considered (Chubb et al. 2020; Helling et al. 2020). Chubb et al. (2020) highlighted that the presence of this molecule could come from disequilibrium or dynamical processes. The higher temperature on the day side allows for this molecule and other metal oxides/hydrides in a gaseous form (Woitke et al. 2018), which could then be transported to the terminator regions (Caldas et al. 2019; Pluriel et al. 2020b). Earlier ground-based observations of the transit (Chen et al. 2014) were also consistent with optical absorption, which were interpreted as pure Rayleigh scattering or absorption from K/Na or TiO/VO. Their observations of the eclipse suggested poor day–night contrast, as also suggested in our results, and potential high-altitude emission on the planet day side.

From a theoretical perspective, recent studies have also shown that thermal inversions could occur naturally in the upper atmosphere of hot Jupiters (Lothringer et al. 2018; Lavvas & Arfaux 2021). For example, the day side of irradiated hot Jupiters might be significantly heated by photochemical hazes and/or dissociation of the main molecules and the addition of continuum opacity from negative hydrogen (H^–). In addition, local thermal inversions are also predicted to occur if sulfur-bearing species are present (Lavvas & Arfaux 2021). While the investigations in Lothringer et al. (2018) were performed on an F-type star with a larger UV flux than WASP-43 (K-type), their fiducial hot-Jupiter simulations showed that planets at distances from their host star similar to WASP-43 b (0.1 au) and with equilibrium temperatures around 1500 K could display inversions above 100 Pa, regardless of their TiO/VO content. The results presented here might provide strong observational evidence in favor of these theoretical predictions.

In Figure 4, we also show the atmospheric structure derived from our model. The planet presents an inflated day side (and hot spot) with a large-scale height difference between the day and the night regions. This conclusion also holds for the models with only two faces, strongly suggesting that using one-dimensional transmission models to analyze transit spectra might lead to large biases for these types of planets as their current formulations lack the flexibility to properly represent such three-dimensional effects. A recent study (Skaf et al. 2020) already found observational evidence of these effects in Hubble transmission spectra and, as next-generation telescopes will be more and more precise, this indicates the importance of considering complementary phase curve studies. For irradiated targets, phase curves will bring additional constraints that will allow us to break these three-dimensional degeneracies. On the night side, we retrieve a large temperature decrease with altitude, which can be inferred from the very low emission of the spectra at these phases. Similar results were already highlighted in Stevenson et al. (2014, 2017), Irwin et al. (2020), and Feng et al. (2020).

Zoom In Zoom Out Reset image size

The abundance profiles predicted by our equilibrium chemistry scheme for this planet are shown in Figure 5. The retrieved metallicity (${1.81}_{-0.17}^{+0.19}$) and C/O ratio (${0.68}_{-0.12}^{+0.11}$) for this atmosphere, shown in the posterior distribution in Appendix K, are particularly well constrained. These correspond to a slightly carbon-enriched atmosphere with a supersolar metallicity, which contrasts with the previous two-faces equilibrium retrieval. For WASP-43 b, global climate models (GCMs) have shown that supersolar metallicities provide a better fit to the Stevenson et al. (2017) data (Kataria et al. 2015). Their work, however, explored lower metallicities than considered in our study. They also highlighted that the presence of night-side clouds provides a good explanation of the observed data. Overall, we find that water remains present in all regions of the atmosphere and that its abundance is consistent with our previous two-faces equilibrium run. Carbon is converted from mainly CO on the day side to CH₄ on the night side, with an overall higher abundance of these species compared with the previous runs. We note that NH₃, which will be a very strong absorber in the spectra from the next-generation telescopes, such as JWST (Greene et al. 2016), Twinkle (Edwards et al. 2019b), and Ariel (Tinetti et al. 2018), might reach detectable abundances on the night side of the planet.

Zoom In Zoom Out Reset image size

In terms of cloud properties, our full scenario did not recover evidence of opaque absorbing layers on the night side. Although clouds have been suggested as an explanation for the very low temperature of the night side, we highlight the fact that all spectra, even at low orbital phases, are showing some prominent water absorption features at 1.4 μm, which explains why clear atmosphere solutions are favored when enough flexibility is given to the model. While we did not detect direct evidence of night-side clouds, these are not ruled out, and because the signal on the night-side phases is lower, the retrieval does not inform us about clouds at pressures higher than 10⁴ Pa.

Overall, when compared to our previous runs without a hot spot and/or the results from Stevenson et al. (2017), Irwin et al. (2020), Feng et al. (2020), we obtain a very different picture. When a hot spot is added, our retrieval finds that thermal inversions (stratosphere and thermosphere) with a rich chemistry provides the best fit to the combined phase curve spectra.

In the three scenarios dealing with real data, we recovered particularly well-constrained solutions. By considering all spectra together, the information content contained in our analysis is increased as compared to traditional single-spectrum retrievals. This can be seen by comparing the results to our phase curve retrievals on individual spectra (see Appendices C and D) and the results from this section with our unified exploration. This allows us to extract more complex information as compared to previous studies. In all three cases, we recover a decent fit of the spectra. This is showcased in Appendix F, where we compare the best-fit spectra of our three scenarios. We note the increase in complexity between the different scenarios allows us to better reflect the observed data, as shown by the increase in log(E), especially in the Spitzer region. For example, the inclusion of the hot spot with a shift of −122 in the Full model improved the capability to fit the slightly higher flux observed from phases 0.0625–0.4375 (See also Figure 2). This is expected as the raw phase curve data from WASP-43 b exhibit day–night contrast and asymmetric emission, features that require a minimum of three regions to be properly described. Interesting different interpretations of the data can be made from the solutions recovered by our three scenarios.

More precisely, our two-faces free chemistry retrieval indicates similar results to previous studies for this planet (Stevenson et al. 2017; Feng et al. 2020), but it also highlights the need for increased complexity, matching the information content in the spectra. We find a hot day side with a large temperature decrease and a cool night side. The retrieved chemistry is also fairly similar, with a relatively low abundance of water. When we introduce variations with altitude and longitude in the chemistry (two-faces model with equilibrium chemistry), the retrieval prefers models with higher water abundances (two orders of magnitude larger). This behavior drives a supersolar metallicity and C/O ratio. In our most complex scenario, when a hot spot is added, the temperature profile displays a stratosphere with a thermal inversion around 10³ Pa on the day side and hot spot. This is associated with a different chemistry with a supersolar metallicity and a slightly higher C/O ratio. The solution also appears stable to changes in the model parameters such as clouds and hot-spot size. In all three scenarios, the best-fit T–p profiles on the day side are consistent with a thermosphere, a high-altitude thermal inversion. In essence, these tests demonstrate the potential of phase curve data in describing three-dimensional effects and the feasibility of extracting this information using retrieval frameworks. This also clearly shows the model dependence of data analysis techniques and the need to use models of adapted complexity.

When performing the retrievals with free hot-spot size and offsets, we encounter solutions that do not match the values recovered in Stevenson et al. (2014, 2017), motivating the need to fix those two values (see Appendix E). While there might be observational reasons to fix the hot-spot offset (we use the peak emission in the phase from Stevenson et al. 2014), the hot-spot size cannot be constrained easily from prior analyses. For our baseline run, it was set to 40° following suggestions from three-dimensional studies (Kataria et al. 2015; Irwin et al. 2020; Helling et al. 2020) but other values might, in fact, also explain the observed data. The impact of the hot-spot size can be investigated by running complementary retrievals with various hot-spot sizes. We considered the cases with 30° and 50° fixed hot-spot sizes. For the chemistry and the clouds (see Appendix K), comparisons of the three simulations indicate very similar results. The temperature structure also remains similar in most of the atmosphere, but there are a few differences in the top part of the hot spot and the day side where the retrieved temperature increases with a smaller hot spot. As only little differences appear, it explains the difficulties in recovering reliable information on the geometry of this planet with current data. Using similar techniques to the one presented in this paper, we however believe that data from the next-generation telescopes, with a larger wavelength coverage and a higher signal-to-noise ratio, might be sensitive enough to infer the hot-spot geometry directly. Section 4 provides an example where the hot-spot shift is directly captured in the retrieval of Ariel data.

In our baseline model, we represented the planet without including clouds on the hot spot (the cloud-top pressure was fixed to 10⁶ Pa). This choice is justified by theoretical predictions from Lee et al. (2016), Parmentier et al. (2016), Lines et al. (2018), and Helling (2019), which suggests that the hot day side of irradiated exoplanets might be cleared up due to strong stellar irradiation from the host stars.

When clouds are included on the hot spot, we recovered a clear hot-spot solution, but the posterior solution also includes a cloudy hot-spot scenario with clouds up to 10³ Pa (see Appendix L). As a result, the corresponding temperature profile in the cloudy hot-spot solution presents larger uncertainties for higher pressures, which are in this case enabling a wider range of solutions without necessarily stratospheric thermal inversions. While opaque clouds are less likely present on the day side of WASP-43 b, recent theoretical models of this planet from Helling et al. (2020) suggest that this region could host some mineral clouds with a large particle size as well as photochemically driven hydrocarbon hazes. With a much stronger day-side emission, the information contained in the phase curve data is greater for these regions, thus potentially allowing us to detect the first pieces of evidence confirming these predictions.

It is known that combining data from multiple instruments can lead to strong degeneracies (Yip et al. 2020, 2021; Irwin et al. 2020; Changeat et al. 2020; Feng et al. 2020; Pluriel et al. 2020a) in retrieval results. While we believe the use of HST along with Spitzer provides a large wavelength coverage, necessary to simultaneously infer atmospheric properties in emission spectroscopy (temperature/clouds/chemistry), the stability of the solution is investigated by reproducing our retrieval on the HST data only. Appendix I provides the spectra centered around the HST wavelength range for our HST+Spitzer retrieval (colored spectra) and our HST-only case (gray). The runs are also extended to the Spitzer region in Figure 29 of Appendix I. Neither fit (with and without the Spitzer photometric points) shows major differences in the HST wavelength range and is well fitted to the observed data. Outside the HST wavelength coverage, however, we see that the inclusion of the Spitzer points is leading to large differences in the retrieval predictions.

For the thermal structure, we found similar results (see Appendix I), with larger errors on the retrieved profiles for the HST-only run. This behavior is expected from the decrease in information content because HST only covers a short wavelength range. The addition of the Spitzer data leads to observational constraints on the carbon-bearing species: CH₄ at 3.6 μm; CO and CO₂ at 4.5 μm. The chemistry, in the HST-only case, pushed toward very high metallicities (log(m) > 2.0 for the posterior distribution in Appendix I). We believe the addition of the Spitzer points is required to simultaneously retrieve the temperature and chemistry profiles in emission. The addition of the 3.6 and 4.5 μm Spitzer points provide valuable information to constrain the carbon-bearing species, thus impacting the retrieved chemistry greatly. We highlight, however, that the Spitzer data points are sensitive to the reduction employed. For the particular case of WASP-43 b, there are now five independent studies (Stevenson et al. 2017; Mendonça et al. 2018a; Morello et al. 2019; Bell et al. 2021; May & Stevenson 2020) that obtained very different reductions, potentially introducing biases to our results. In addition, most models use simple sine functions to reduce phase-curve data, and alternative physically motivated approaches might provide different results (Louden & Kreidberg2018).

In this work, we assumed a particular geometry when computing the phase curve emission. The solutions found with this new model contrast with previous analyses (Stevenson et al. 2017; Irwin et al. 2020; Feng et al. 2020) and highlight new possibilities for the atmosphere of WASP-43 b. The choices made here aimed to achieve a balance between our current understanding of the physical properties of irradiated planets and the complexity of the analyzed data. However, as of today, there is no obvious way to assess the amount of information that can be extracted from an exoplanet spectrum or a set of spectra. For example, we demonstrated that the analysis of phase curve data with unified phase curve models allows for the extraction of a more complete and complex picture than standard techniques. However, there is no guarantee that our model is complex enough or best represents the planet's geometry. In other terms, when fitting exoplanet data, one always recovers a model-dependent solution and in our case, other three-dimensional geometries might be more accurate or relevant than what we presented. In the case of phase curves, we believe complementary work (testing more complex/simple geometries, comparing with simulated data) must be performed to completely understand the model-dependent errors introduced by our choices. We highlight, however, that the technique presented in this work is very general and can be adapted to other three-dimensional geometries. We intend to evolve this model to reflect the advancements in our understanding of three-dimensional effects.

Similarly, this comment also applies to physical and chemical assumptions, which might lead to biases, too. In our most complex run, we assumed the planet chemistry is in thermochemical equilibrium and that the chemistry is coupled between the different regions. While suggested as a good approximation for hot-Jupiter planets in the temperature ranges of WASP-43 b, this may not be the case and previous studies already demonstrated the dangers of assuming equilibrium chemistry for exoplanets undergoing disequilibrium processes (Line & Yung 2013; Rocchetto et al. 2016; Mendonça et al. 2018b; Blumenthal et al. 2018; Changeat et al. 2019; A. F. Al-Refaie et al. 2020, in preparation; Changeat et al. 2020; Steinrueck et al. 2019; Drummond et al. 2020; Venot et al. 2020). In particular, the thermal inversion observed in the full model for pressures lower than 1 Pa is driven by the need to increase the temperature and thermally dissociate all the molecules, which makes the atmosphere transparent there. In fact, if other mechanisms are able to remove the molecules at those altitudes, the retrieval would most likely not require such strong thermal inversions, which we believe is linked to the equilibrium chemistry assumption. Such processes could, for example, come from photodissociation. In addition to this, the current equilibrium scheme we use does not include exotic molecules such as TiO, VO, AlO, or H^–. This might also lead to biases and have an important impact on the predicted molecular abundances, especially as the results from this analysis suggests the presence of thermal inversions on the day side of this planet. The definitive detection of such optical absorbers would provide further evidence to verify the thermal inversion hypothesis. Similarly, clouds were assumed as fully opaque gray opacities, which is a simplification that will have to be challenged when considering higher-quality spectra from next-generation telescopes (Bohren & Huffman 2008; Heng & Demory 2013; Lee et al. 2013; Powell et al. 2018, 2019; Barstow 2020; Helling et al. 2020). For WASP-43\,b, reflected light, which is not included here, is shown to have an important impact on the energy budget of this planet (Keating & Cowan 2017). Finally, our exploration of the parameter space was performed using the MultiNest routine (Feroz et al. 2009), which is well established in the exoplanet community. However, the recovered solutions might also depend on the sampling algorithms used, for example, alternative nested sampling algorithms are described in Handley et al. (2015), Speagle (2020), and Buchner (2021); the chosen settings; and the convergence criteria. We tested our full scenario with an increased number of live points (2500) to verify the stability of our solution to this parameter. The posterior distribution, shown in Appendix M, does not indicate different results compared to the 500 live points run for this example, but we note that retrieval abilities to separate multimodal solutions might greatly depend on this parameter.