On the difficulties of obtaining absolute transit depths with HST WFC3: KELT-11 b, an example

ABSTRACT

The study of exoplanetary atmospheres with low-resolution transmission spectroscopy relies on measuring minute changes in the transit depth with wavelength and a number of ground- and space-based instruments have been used to characterize exoplanets in different spectral bands. For the last decade, these instruments have each only probed a narrow spectral region, which has motivated the community to combine observations from different instruments in order to achieve a broader wavelength coverage. By analysing Hubble Space Telescope (HST) Wide Field Camera 3 (WFC3) data of KELT-11 b, we once again show the risks of following this now conventional approach. We demonstrate that changes in the reduction or analysis method can lead to drastic differences in the mean transit depth and that combining this with additional data can lead to discrepant interpretations of the atmospheric composition. With the launch of JWST, and its many available instruments and modes, observers may be tempted to combine data sets at longer wavelengths (e.g. NIRSpec – Near Infrared Spectrometer) with those from HST STIS (Space Telescope Imaging Spectrograph) or WFC3 without the consideration of offsets or other incompatibilities. Given the obvious potential issues, we caution against such an approach and encourage the community to thoroughly address the issue of data incompatibility instead of adhering to a de facto assumption of compatibility.

1. INTRODUCTION

Transmission and emission spectroscopy with the Wide Field Camera 3 (WFC3) of the Hubble Space Telescope (HST) have formed the cornerstone of exoplanet atmospheric characterization, enabling the discovery of water and clouds in many exo-atmospheres (e.g. Stevenson et al. 2014; Tsiaras et al. 2018; Skaf et al. 2020; Cubillos & Blecic 2021; Mansfield et al. 2021; Roudier et al. 2021; Changeat et al. 2022; Edwards et al. 2023b). While the detection of water is now routine for the hot-Jupiter class of planets, other molecules such as carbon-bearing or refractory species have been more challenging to detect because their spectral features occur at the edge of, or outside of, the spectral range of WFC3 (e.g. Changeat et al. 2020). To combat this issue, studies have often combined data from different instruments to facilitate the detection of additional species (e.g. Sing et al. 2016; Wakeford et al. 2017; Benneke et al. 2019; Pinhas et al. 2019; Welbanks et al. 2019; Alam et al. 2020; Colón et al. 2020; Pluriel et al. 2020; Changeat & Edwards 2021; Kawashima & Min 2021; Brande et al. 2022; Kreidberg et al. 2022; Mikal-Evans et al. 2023).

HST was never designed to undertake transit observations and, due to this and some basic characteristics of the observatory (e.g. its low-Earth orbit), spectroscopic time-series observations with WFC3 experience significant systematics (e.g. Kreidberg et al. 2014; Stevenson et al. 2014; Tsiaras et al. 2016a). While the correction of these with different pipelines usually leads to a robust recovery of the spectral features, offsets in the transit or eclipse depth are common (e.g. Anisman et al. 2020; Changeat et al. 2020; Colón et al. 2020; Mansfield et al. 2021; Mugnai et al. 2021; Changeat et al. 2022; Libby-Roberts et al. 2022; Edwards et al. 2023b). The choice of detrending method, limb-darkening coefficients, and orbital parameters can also cause differences in the mean transit depth (e.g. Stevenson et al. 2014; Tsiaras et al. 2018; Guo et al. 2020; Yip et al. 2020) and different observations of the same planet can also lead to differing white light-curve depths (e.g. Edwards et al. 2023b). Finally, offsets between HST WFC3 data and other instruments are also common (e.g. Luque et al. 2020; Murgas et al. 2020; Yip et al. 2021; McGruder et al. 2022; Constantinou, Madhusudhan & Gandhi 2023; Jiang et al. 2023).

When combining spectra from different instruments, offsets between these could cause differences in the retrieved atmospheric parameters. In some cases, the offsets make it evident that the data sets are completely incompatible in their native form (e.g. Yip et al. 2021), but this may not always be the case. While an offset can be fitted for in the retrieval (e.g. Luque et al. 2020; Murgas et al. 2020; Yip et al. 2021; McGruder et al. 2022; Edwards et al. 2023a), one cannot be sure of the compatibility of the spectra, particularly when there is no spectral overlap and/or a photometric point is used [e.g. from Kepler, Spitzer Infrared Camera Array (IRAC), or Transiting Exoplanet Survey Satellite (TESS)].

Until the launch of JWST (Gardner et al. 2023), only HST WFC3 G102 offered space-based spectral overlap with HST WFC3 G141, but this grism has rarely been used. Combining both the WFC3 IR grisms with data from HST STIS can be done with some confidence due to this spectral overlap (e.g. Wakeford et al. 2017), although the data are not taken in the same epoch so temporal changes may be an issue (e.g. Cho et al. 2003; Skinner & Cho 2021; Saba et al. 2022; Changeat et al. 2024).

The planet KELT-11 b was discovered in 2016 orbiting a bright G star (Kmag = 6.122), with an orbital period of 4.736 d (Pepper et al. 2017). Due to its very low density (0.093 g cm|$^{-3}$|⁠), it was immediately associated with a very large-scale height and was predicted to become one of the ‘benchmark’ planets for atmospheric characterization. Further observations from the ground and the Spitzer Space Telescope refined the orbital and stellar parameters (Beatty et al. 2017). Žák et al. (2019) analysed high-resolution data from HARPS, searching for sodium, hydrogen, and lithium in the atmosphere of KELT-11 b. They saw no evidence for these species, and they attributed the non-detection to the possible presence of high-altitude clouds or hazes. They also reported a low stellar activity of the host star, a result obtained by monitoring the Ca i and Mg i lines.

In 2018, HST observed a transit of KELT-11 b using the G141 grism of the WFC3 instrument and the data were published by two independent teams, Colón et al. (2020) and Changeat et al. (2020). Both studies noted the strange nature of the spectrum, which had a strong absorption feature at wavelengths longer than the usual 1.4 |$\mu$|m water band, a feature not seen in other WFC3 spectra. The spectra obtained by these teams were highly similar in shape, but the mean transit depths of the spectra disagreed by |$\sim$|260 ppm. In Colón et al. (2020), several pipelines were utilized to reduce and analyse the data, each finding similar spectral shapes but with differing mean depths of up to 360 ppm. Nevertheless, they choose to analyse their prime HST WFC3 spectrum jointly with TESS (0.6–1 |$\mu$|m) and Spitzer (3.2–3.9 |$\mu$|m) photometric points, conducting retrievals on these which indicated the presence of AlO as the TESS transit depth was much deeper than the HST WFC3 G141 data. They also tried including vertical offsets in the retrievals, but only considered offsets of up to 80 ppm, seemingly ignoring the magnitude of the variability highlighted by their different reduction and analysis methods. They stated that these offset retrievals still preferred the presence of AlO. However, conducting joint retrievals on the spectrum from Changeat et al. (2020), as well as other spectra derived by Colón et al. (2020), would likely have led to different results because of the relative change in depth between the HST WFC3 and TESS/Spitzer data. Due to the lack of wavelength overlap, it was not possible to tell which, if any, spectrum had the correct HST depth and thus any analysis using the TESS, or Spitzer, data are unreliable.

Here, we present a thorough analysis of KELT-11 b’s atmosphere via optical and near-infrared transit data. We assiduously explore the effects of different extraction methods, limb-darkening coefficients, and long-term trends on the WFC3 G141 data. We show how combining these data with photometric data from TESS without considering an offset can lead to extreme differences in the interpretation of the atmosphere and how, because of a lack of wavelength overlap, one cannot be sure which solution is correct. Furthermore, we present the analysis of another transit of KELT-11 b from HST WFC3, this time using the G102 grism. We use the wavelength overlap of the two WFC3 IR grisms to conduct retrievals on the combined spectrum, comparing these results with those from the WFC3 G141 and TESS + WFC3 G141 retrieval results. We also perform retrievals on the HST WFC3 spectrum, and photometric data, from Colón et al. (2020) to highlight the impact of restricting the offsets between transit depth measurements from different instruments. Finally, we again strongly caution against combining transit depth measurements which do not have wavelength overlap without at least including an offset parameter, hoping that future analyses with HST, and JWST, will consider that the data sets may not be automatically compatible and act accordingly.

2. HST WFC3 TRANSIT OBSERVATIONS

Here, we analyse two data sets from HST WFC3. The first was taken with G141 (1.088–1.68 |$\mu$|m) in April 2018 (PN: 15225, PI: Knicole Colón, Colón 2017) and has previously been analysed by two independent works (Changeat et al. 2020; Colón et al. 2020). The second data set also consists of a single transit of KELT-11 b, although this time with the G102 grism (0.80–1.15 |$\mu$|m), and it was taken in December 2020 (PN: 15926, PI: Knicole Colón, Colón et al. 2019). Both HST visits consisted of nine orbits with the data being taken using the spatial scanning mode. During an exposure using the spatial scanning mode, the instrument slews along the cross-dispersion direction, allowing for longer exposure times and increased signal-to-noise ratio, without the risk of saturation. Both forward (increasing row number) and reverse (decreasing row number) scanning were used for these observations to increase the duty cycle. In both cases, the detector settings were SUBTYPE = SQ512SUB, SAMP|$_{\rm SEQ}$| = SPARS25, NSAMP = 4, and APERTURE = GRISM512, and the scanning speed was 0.96 arcsec s|$^{-1}$|⁠. The final images had a total exposure time of 46.695 518 s and a total scanning length of 51.312 arcsec. These led to a maximum signal level of |$\sim$|36 000 electrons per pixel for the G102 observations and |$\sim$|38 000 electrons per pixel for G141 observations. For calibration purposes, a 2.55 908 s non-dispersed (direct) image of the target was taken at the beginning of each orbit, using the F130N filter and the following settings: SUBTYPE = SQ512SUB, SAMP|$_{\rm SEQ}$| = RAPID, NSAMP = 4, and APERTURE = GRISM512.

We obtain the publicly available data from the HST MAST archive.¹ We carry out the analysis of both HST transits using Iraclis, our highly specialized, open-source² software for processing WFC3 spatially scanned spectroscopic images (Tsiaras et al. 2016a, c, 2018, 2019). The reduction process includes the following steps: zeroth read subtraction, reference pixel correction, non-linearity correction, dark current subtraction, gain conversion, sky background subtraction, calibration, flat-field correction, and bad pixel/cosmic-ray correction. Then, we extracted the white (0.800–1.150 |$\mu$|m and 1.088–1.680 |$\mu$|m for G102 and G141, respectively) and the spectral light curves from the reduced images, taking into account the geometric distortions caused by the tilted detector of the WFC3 infrared channel. For a more detailed account of the procedures followed in Iraclis, we refer the reader to Tsiaras et al. (2016c).

In the extraction of the flux from the calibrated images, Iraclis offers two techniques. The first extracts the flux from the image as it is at the end of the scan. The second splits the data into the up-the-ramp reads (i.e. the non-destructive reads), subtracting each one sequentially from the next read and extracting the flux from each resulting image. The flux is then extracted from each of the differential images, as shown in the Appendix, before being summed. We refer to these as ‘full-scan’ and ‘splitting’ extractions, respectively. The splitting extraction is useful for ensuring that no background stars have overlapping contributions to the target star’s spectrum. Additionally, differences between the flux extracted with each method can also occur due to persistence. While the level of persistence is dependent upon the accumulated charge, it is also correlated with the time the charge has spent on the detector as well as the time since the last read (Anderson et al. 2014). In the context of HST WFC3, this means the persistence is dependent upon the brightness of the host star, the scanning rate, and the readout scheme employed.

To observe KELT-11 b (Kmag = 6.122), a rapid scan rate was required to avoid saturating the detector. The set-up is such that if one takes a raw image for either visit, and subtracts the last non-destructive read from the destructive read, the persistence effect can clearly be seen. In Fig. 1, we show this effect for the second forward scan of the second orbit (the first image used in the analysis) for each visit. We also looked at the last forward scan of the last orbit (the penultimate image used in the analysis), but did not find any noticeable difference with time. However, differences are clear when comparing the G102 data sets to those taken with the G141 grism: the persistence effect is greater in the latter case, likely due to the higher flux levels as the scan rate and readout schemes were identical. Such an effect can also be seen in other HST WFC3 spatial scanning data sets (e.g. Edwards et al. 2023a).

Therefore, due to the persistence effect, the extracted flux using the full-scan and splitting methods are different. The light curves for these are given in the Appendix and one also notices that the flux extracted from the forward and reverse scans is different. The upstream/downstream effect is caused by the forward and reverse scans having slightly different exposure times (McCullough & MacKenty 2012). For the fits undertaken in this work, the splitting extraction has been utilized unless explicitly stated otherwise.

2.1 Baseline white light-curve fitting process

We fit the light curves using our transit model package PyLightcurve, with the transit parameters from Beatty et al. (2017) and limb-darkening coefficients calculated based on the PHOENIX (Allard, Homeier & Freytag 2012) model (see tables in Appendix  A). We use the claret non-linear formula (Claret 2000) and take the stellar parameters from Beatty et al. (2017). These were computed using ExoTETHyS (Morello et al. 2020). During our fitting of the white light curve, the planet-to-star radius ratio and the mid-transit time were the only free parameters, along with a model for the systematics (Stevenson et al. 2014; Tsiaras et al. 2016c).

It is common for WFC3 exoplanet observations to be affected by two kinds of time-dependent systematics: the long- and short-term ‘ramps’. The first affects each HST observation sequence and is usually represented by a linear behaviour, while the second affects each HST orbit and is modelled as having an exponential behaviour. The formula we used for the white light-curve systematics (R_w) was the following:

where t is time, |$n^{\rm scan}_{\rm w}$| is a normalization factor, |$T_{\rm 0}$| is the mid-transit time, |$t_{\rm o}$| is the time when each HST orbit starts, |$r_{\rm a}$| is the slope of a linear systematic trend along each HST visit, and (⁠|$r_{b_1},r_{b_2}$|⁠) are the coefficients of an exponential systematic trend along each HST orbit. The normalization factor we used (⁠|$n^{\rm scan}_{\rm w}$|⁠) was changed to |$n^{\rm for}_{\rm w}$| for upward scanning directions (forward scanning) and to |$n^{\rm rev}_{\rm w}$| for downward scanning directions (reverse scanning). The reason for using different normalization factors is the slightly different effective exposure time due to the known upstream/downstream effect (McCullough & MacKenty 2012).

We did not use the first orbit of each visit in the analysis due to the stronger ramps seen in this data set. Furthermore, we removed the first forward and reverse scans from each visit. We then use emcee (Foreman-Mackey et al. 2013) to explore the parameter space and we utilize 200 walkers with 1000 000 iterations including a burn-in of 700 000. We fit the white light curves using the formulae above and the uncertainties per pixel, as propagated through the data reduction process. However, it is common in HST WFC3 data to have additional scatter that cannot be explained by the ramp model. For this reason, we scale up the uncertainties in the individual data points, for their median to match the standard deviation of the residuals, and repeat the fitting (Tsiaras et al. 2018). We show the white light-curve fits resulting from the splitting extraction in Fig. 2. We note that, in both cases, the residuals of the fit are clearly not Gaussian. However, due to the divide-by-white method, the residuals on the spectral light curves are far better, as discussed in Section 2.3.

2.2 Alternative white light-curve fits

Assumptions made in the analysis process described in Section 2.1 can have an impact on the final spectrum derived. Hence we performed a number of additional fits to the data where we adjusted the assumptions to explore the effect on the transmission spectrum. Table 1 provides a guide to all the different fits conducted, with additional detail given in Sections 2.2.1–2.2.5. For brevity, the plots for these additional white light-curve fits are given in Appendix  A.

2.2.1 Long-term trend

One of the main systematic effects seen in HST observation sequences is the visit-long trend. As discussed in Section 2.1, this is usually assumed to be a linear effect with time. However, other models for this systematic have been suggested or attempted. From Fig. 2, it is clear that the systematics present in the data set have not been fully accounted for in the linear model. While the divide-by-white method means that the spectral light curves are well-fitted, this still provides the motivation to explore different long-term trends.

Guo et al. (2020) noted that the linear detrending model did not appear to fully account for the systematics seen in the HST WFC3 G141 transit light curves of HD 97658 b. They showed that the choice of detrending model (linear, quadratic, logarithmic, or exponential) affected the recovered white light-curve depth but did not overly affect the spectral features seen. Another literature example of using a quadratic detrending model is Line et al. (2016) who fitted the HST WFC3 eclipse observations of HD 209458 b under this assumption. Changeat et al. (2022) then showed that a much deeper eclipse depth was recovered if a linear trend was used, further suggesting that the choice of long-term trend can have a significant impact on the light-curve fitting. Similar results were found for the HST WFC3 transit observation of WASP-12 b (Stevenson et al. 2014).

The two previous analyses of the HST WFC3 G141 data of KELT-11 b used different long-term trend models: Changeat et al. (2020) assumed a linear long-term trend while Colón et al. (2020) utilized a quadratic model. The mean transit depth recovered by each study was very different (⁠|$\sim$|260 ppm), but both had similar spectral features.

Therefore, here we explore the effect of using a quadratic long-term trend, conducting fits using it for all the variations mentioned ahead. In this case, the model for the systematics (equation 1) becomes

where r|$_{\rm a2}$| is the additional fitted parameter.

2.2.2 Extraction method

As discussed in Section 2, and shown in the Appendix, the choice of extraction method has an effect on the extracted flux. Therefore, we conduct a fitting on the data that did not split the data by the up-the-ramp reads (i.e. the full-scan extraction) to explore whether the lower flux levels extracted would affect the subsequent transit spectrum. In this case, the Claret limb-darkening law (Claret 2000) was always used and we fitted the light curve using both the linear and quadratic long-term trends.

2.2.3 Limb-darkening coefficients

The emission from the visible disc of the exoplanet host star is not homogeneous, with the limbs contributing less flux. Therefore, over the course of a transit, the planet does not block out an equal amount of light at each point in time. The effect is particularly noticeable during ingress and egress, giving visible and near-infrared transits their curved shape rather than a more box-like profile. Many limb-darkening laws have been utilized to model this effect and the choice of coefficients can affect the recovered transit depth (e.g. Tsiaras et al. 2018). To explore this effect, we conduct an extra fitting of the HST light curves using linear limb-darkening coefficients. These have been frequently used in the literature for fitting HST data (e.g. Kreidberg et al. 2014) and these were again calculated using ExoTETHyS.

2.2.4 Orbital parameters

In our standard fit of the HST WFC3 light curves, we fix all transit light-curve parameters except for the transit mid-time and the planet-to-star radius ratio. The approach is common in the analysis of HST data as there significant gaps in the transit coverage due to HST’s low Earth orbit which results in the Earth obscuring the target for portions of each HST orbit. These gaps often mean that ingress and egress are poorly sampled, limiting the information content of the light curve. Nevertheless, some studies have fitted for other parameters, such as the inclination of the planet’s orbit with respect to the observer, and Alexoudi et al. (2018) found that using the different inclinations could lead to varying slopes in optical data. In this case, the HST WFC3 G141 data sampled both ingress and egress, so we also attempt a fit with the inclination as a free parameter. We note that the previous studies of HST WFC3 G141 data took different approaches: Colón et al. (2020) fitted for the inclination while Changeat et al. (2020) did not. Although the G102 data samples neither ingress nor egress, we also performed a similar fit to these data for completeness.

2.2.5 Changing the number of baseline orbits used

The out-of-transit orbits of HST data are used to provide a baseline flux to determine the transit depth. The first orbit of each observation is often removed due to the stronger ramp seen, although we note that Zhou et al. (2017) have proposed a method of accounting for this charge-trapping effect. For the majority of transit observations with HST, removing this first orbit leaves one with two out-of-transit orbits (one before transit and one after) as well as one or two in-transit orbits. However, in the case of KELT-11 b, the transit duration is far longer than that of most planets studied and, therefore, more baseline orbits have also been taken. In the case of the G102 data there are four baseline orbits while for the G141 visit there are three.

As previously noted, in Fig. 2, is it clear that the long-term systematics have not been well fitted. In Changeat et al. (2020), additional baseline orbits were removed in an attempt to achieve a better white light-curve fit while still using only a linear long-term trend. Whether such an approach is advisable is, of course, debatable: more baseline should give us a more accurate and precise constraint on the transit depth, but only if our model is capable of reproducing the systematics in the data. Therefore, in this study we also attempted fits in which we removed additional out-of-transit orbits to see the impact on the final spectrum, both in terms of the depths recovered but also the precision. In all cases, we made sure there was at least one orbit both before and after the transit.

2.3 Spectral light-curve fitting

For each of our white light-curve fits, we also perform a fit to the spectral light curves. In this case, the only free parameter within the transit model is the planet-to-star radius ratio. We also fit a model for the systematics (⁠|$R_\lambda$|⁠) that includes the white light curve (divide-by-white method, Stevenson et al. 2014) and a wavelength-dependent, visit-long slope (Tsiaras et al. 2016c). The model is given by

where |$\chi _\lambda$| is the slope of a wavelength-dependent long-term systematic trend along each HST visit, |$LC_{\rm w}$| is the white light curve, and |$M_{\rm w}$| is the best-fit model for the white light curve. Again, the normalization factor we use (⁠|$n^{\rm scan}_\lambda$|⁠) was changed to (⁠|$n^{\rm for}_\lambda$|⁠) for upward scanning directions (forward scanning) and to (⁠|$n^{\rm rev}_\lambda$|⁠) for downward scanning directions (reverse scanning). Also, in the same way as for the white light curves, we perform an initial fit using the pipeline uncertainties and then refit while scaling these uncertainties up, for their median to match the standard deviation of the residuals. For each spectral light-curve fit, the Markov chain Monte Carlo (MCMC) has 100 walkers and 300 000 iterations with a burn-in of 200 000.

The best-fit models for the spectral light curves for the standard fits are shown in Fig. 3. From Fig. 3, it is immediately clear that the first two spectral bins of the G102 data are not well fitted and the same is found for all fits of this data set.

The inability of Iraclis to provide an adequate modelling of the systematics for these spectral channels stems from the methodology used and the change in the flux when the stellar light is reduced by the transmission of each element of the HST WFC3 G102 optical chain. As noted earlier, the model for the white light-curve systematics is applied to the spectral light curves with an additional, wavelength-dependent normalization factor (⁠|$n^{{\rm scan}}_\lambda$|⁠). The ramps seen in each orbit of HST data are due to charge-trapping and are heavily linked to the flux levels. In most observations taken by HST thus far, the change in flux levels is relatively small with wavelength and thus the ramps have a similar behaviour in each spectral channel: the normalization factor is enough to deal with any differences and thus lead to a good light-curve fit. However, in this case, the flux levels are very different between the bluest and reddest spectral bin of the G102 data, leading to different ramp behaviours.

Correcting for this effect would require changing the fundamental approach that Iraclis takes to analysing HST data. Therefore, we did not attempt to correct for it here and instead choose simply to not include the two bluest spectral bins in any of our retrievals. However, efforts should clearly be conducted in the future to overcome this issue.

3. ADDITIONAL TRANSIT OBSERVATIONS

Due to the narrow wavelength coverage of HST WFC3 G141, it has become common practice combining these data with observations from other instruments and facilities. Generally, these observations do not have any spectral overlap and it is assumed that the data sets are compatible (i.e. that both facilities have correctly recovered the absolute transit depth). We add-in additional data to explore the impact of this assumption.

3.1 TESS

We use the 2 min cadence Pre-search Data Conditioning light curves (Smith et al. 2012; Stumpe et al. 2012, 2014). KELT-11 b was studied in Sectors 9 and 62. The data are clipped, with only data within |$\pm$|3 transit durations (21.3 h) of the predicted mid-times of the nine full transits being utilized. These are fitted simultaneously with the planet-to-star radius ratio |$R_{\rm p}/R_{\rm s}$| and transit mid time (⁠|$T_0$|⁠) as free parameters. As with the HST fits, the other parameters are fixed to those from Beatty et al. (2017). For detrending systematics, we fit a linear visit long trend for each transit. We perform the fits using Claret limb-darkening coefficients and, as with the HST data, these are calculated using ExoTETHyS (Morello et al. 2020). We use the emcee (Foreman-Mackey et al. 2013) package to sample the parameter space, utilizing 100 walkers and 150 000 iterations with a burn-in of 50 000.

3.2 Spitzer

A transit observation was taken with Spitzer’s IRAC at 3.6 |$\mu{\rm m}$| (GO-12096, PI: Thomas Beatty). The observation was previously analysed in Beatty et al. (2017) and Colón et al. (2020). However, we decided not to fit these data in this analysis because the transit occurred earlier than expected, meaning the pre-ingress part of the light curve is missing (Beatty et al. 2017). As noted in Changeat et al. (2020), the ‘ramp’ effect is especially pronounced during the settling of IRAC observations (see e.g. Agol et al. 2010) so the processing of this light curve by any of the standard detrending techniques (e.g. May & Stevenson 2020; Morvan et al. 2020) could lead to larger uncertainties in the retrieved transit depth and potentially lead to a poor accuracy as well. We note that in Section 5.3 we include the Spitzer observation as fitted by Colón et al. (2020).

4. RETRIEVAL SET-UP

To explore how the choice of detrending technique affected our understanding of the atmospheric composition of KELT-11 b, we conduct a number of retrievals. Each spectrum is analysed using our Bayesian retrieval framework TauREx 3 (Waldmann et al. 2015a, b; Al-Refaie et al. 2021, 2022). We utilize the absorption cross-sections from the ExoMol data base (Tennyson & Yurchenko 2012; Tennyson et al. 2016; Chubb et al. 2021) and explore the parameter space with the algorithm MultiNest (Feroz, Hobson & Bridges 2009; Buchner et al. 2014) with 750 live points and an evidence tolerance of 0.5. We adopt uniform priors for all the free parameters.

The atmosphere is modelled between |$p \in [10^{-4}, 10^6]$| Pa using 100 plane-parallel layers uniformly partitioned in log-space. For the trace gases we consider the molecules: H|$_2$|O (Barton et al. 2017; Polyansky et al. 2018), CH|$_4$| (Hill, Yurchenko & Tennyson 2013; Yurchenko & Tennyson 2014), CO (Li et al. 2015), CO|$_2$| (Rothman et al. 2010), NH|$_3$| (Coles, Yurchenko & Tennyson 2019), HCN (Harris et al. 2006; Barber et al. 2013), TiO (McKemmish et al. 2019), VO (McKemmish, Yurchenko & Tennyson 2016), FeH (Bernath 2020), and AlO (Patrascu, Yurchenko & Tennyson 2015).

We consider only free chemistry for all retrievals and explore the effect of different data sets and detrending techniques. In these free chemistry retrievals, we fit each of the molecule abundances in volume mixing ratios with priors of log|$_{\rm 10}$|(VMR) |$\in [-15, -1]$|⁠.

Based on the expected chemical equilibrium profiles (Changeat et al. 2020), chemical variations with altitude can occur for H|$_2$|O and CO|$_2$|⁠, so we employed a two-layer model as described in Changeat et al. (2019). The rest of the atmosphere is composed of H|$_2$| and He for which the ratio is fixed to the solar value (He/H|$_2$| = 0.1765). On top of the molecular absorption, we include absorption from the H- ion (John 1988) using the methodology described in Edwards et al. (2020), Rayleigh scattering (Cox 2015), and Collision Induced Absorption processes from H|$_2$|–H|$_2$| (Abel et al. 2011) and H|$_2$|–He (Abel et al. 2012) pairs.

Changeat et al. (2020) noted that, when analysing the HST WFC3 G141 data of KELT-11 b, the recovered chemistry is highly dependent upon the inputs to the retrieval. The purpose of this investigation is not to discover the definitive atmospheric model for KELT-11 b per se, but to understand the impact on the retrieved atmospheric chemistry when applying different detrending techniques and combining different instruments. Hence, we do not further explore this dependency, instead keeping a constant set-up which, we believe, included all relevant absorbers. Therefore, any changes found by the retrievals are due to changes in the information content of the input data.

As mentioned, different assumptions during the fitting of data can lead to offsets in the mean transit depth. While this generally has little impact if analysing a single data set, when combining data from different instruments this offset could lead to discrepant results. Therefore, to attempt to overcome potential offsets between data sets, several studies have introduced an offset parameter within the atmospheric retrievals (e.g. Luque et al. 2020; Murgas et al. 2020; Yip et al. 2021; McGruder et al. 2022; Edwards et al. 2023a). Hence, for several retrievals in this work, we add an offset parameter (⁠|$\Delta$|⁠) to explore its ability to correct for inconsistencies between data sets. In the majority of cases, we set the bounds of this to be |$\Delta \in [-1000, +1000]$| ppm to ensure it was not overly constricted. The exception is for one retrieval on the data from Colón et al. (2020), where we replicated their methodology of using bounds of |$\Delta \in [-80, +80]$| ppm.

5. RESULTS

We fit the two HST WFC3 data sets using a variety of common assumptions in the fitting and detrending process. We find that these can cause large differences in the white light-curve depth and thus in mean transit depth across the spectral data. In Fig. 4, we show the R|$_{\rm p}$|/R|$_{\rm s}$| from each of the fits for both the G102 and G141 data sets. We find a large variation in the white light-curve depth, particularly between the fits with different long-term trends. The spectrum from each fitting are shown in Figs 5 and 6, and we find generally consistent spectral features. We find an offset between the standard fits to the HST WFC3 G102 and G141 data sets, which is obvious due to the spectral overlap.

In Section 5.1, we quantitatively discuss the differences that are seen. Then, in Section 5.2, we explore the impacts on the atmospheric chemistry caused by analysing these different spectra, firstly alone and then when combining these different fits to the HST WFC3 G141 data with photometric observations from TESS. We then discuss the results of fitting for an offset between the TESS and WFC3 data before showing the impact of the G102 spectra. Finally, we discuss the retrievals on the data from Colón et al. (2020) and show the impact of overly constraining the offset parameters in the retrieval.

5.1 Variations in the transit depth due to fitting assumptions

5.1.1 Effect of different long-term trends

We find that the choice of long-term trend has a serious impact on the recovered transit depth. Throughout all the fits, the quadratic long-term trend leads to shallower transit depths. In the case of our standard fitting, the quadratic trend results in a spectrum which had a mean transit depth which was 280 ppm shallower than with the linear trend for the G141 data and 250 ppm shallower in the G102 case. For comparison, the uncertainties on these white light-curve depths is |$\sim$|16 ppm. Fig. 7 clearly shows that the additional parameter added for the quadratic fit, |$r_{\rm a2}$|⁠, is correlated with the transit depth (R|$_{\rm p}$|/R|$_{\rm s}$|⁠). In the quadratic fit, the normalization factors (⁠|$n^{\rm for}_{\rm w}$|⁠, |$n^{\rm rev}_{\rm w}$|⁠) are correlated with the transit depth but also with |$r_{\rm a2}$|⁠. However, for the linear fit, there are no strong correlations between these normalization factors and the transit depth, suggesting the additional term added for the quadratic fit is driving this correlation. Nevertheless, as noted in the Appendix, the quadratic trend often leads to white light-curve fits where the residuals were more Gaussian than in the cases where a linear long-term trend was used. Therefore, it seems that the quadratic trend may be better at reproducing the systematics present in the data but with a loss of confidence in the absolute white light-curve depth. We note that our fits with a quadratic trend give white light-curve-similar depths to the value report by Colón et al. (2020).

5.1.2 Effect of different extraction methods

We find that utilizing different extraction methods leads to small offsets in the HST WFC3 data. As discussed in Section 2, the splitting extraction leads to a higher flux level compared with an extraction using only the destructive reads due to flux lost to charge traps during the scan. Therefore, the G141 depth recovered from the splitting extraction is shallower by 70 ppm in the case where a linear long-term trend is used for both fits. As the persistence effect was less pronounced, the G102 case was less extreme, differing only by 8 ppm which is less than the 1|$\sigma$| uncertainties on the white light-curve depth.

5.1.3 Effect of different limb-darkening coefficients

We find that utilizing different limb-darkening coefficients leads to small offsets in both the HST WFC3 data. Intriguingly, the offset induced by changing the limb-darkening coefficients was not equal across the two grisms: the linear coefficients led to a G102 depth which was 15 ppm deeper while the shift in the G141 data was over 70 ppm. These values are quoted for the linear long-term trend while the quadratic trend leads to differences of 12 and 23 ppm for G102 and G141, respectively.

5.1.4 Effect of removing baseline orbits

We find that the number of out-of-transit orbits utilized also changes the white light-curve depth. The impact that removing a baseline orbit had was generally dependent upon the long-term trend that was being used to fit the data. For the linear trend, removing baseline orbits decreased the white light-curve depth while in the quadratic case the depth generally increased. The absolute difference in the G141 white light-curve depth from the standard fit varied from 24 to 175 ppm. For G102, the variation was between 45 and 145 ppm.

5.2 Effect on the recovered atmospheric composition

As discussed, the mean depth of the HST WFC3 G102 and G141 spectra change depending on the detrending assumptions. However, these spectra have highly similar spectral features. We now explore the changes in the composition seen when performing atmospheric retrievals on each of these. Firstly, we investigate the impact when only data from HST WFC3 G141 is used before also adding additional data. While there are now data with the G102 grism, we initially consider combining the TESS and HST WFC3 G141 data to highlight potential issues that arise from combining data sets without spectral overlap, a common practice within the field. Finally, we consider the G102 data, fitting for the offset seen between it and the G141 data. With the exception of adding in offsets between data sets as a free parameter, the retrieval set-ups are always identical.

For brevity, we do not preform retrievals on all of our reductions. Instead, we discuss the differences caused by the linear and quadratic long-term trends as these have among the highest impact and the choice of this long-term trend has been noted to affect the mean transit depth of other HST data sets (e.g. Stevenson et al. 2014; Guo et al. 2020; Changeat et al. 2022). Additionally, we explore the impact of using the splitting extraction or full-scan extraction as this also had a significant impact and both are common extraction choices.

5.2.1 HST WFC3 G141 only

We firstly explore the impact of different detrending techniques and assumptions on the atmospheric composition when using only data for HST WFC3 G141. We find that the retrievals generally prefer similar molecular abundances but, as seen in Fig. 8, there are some exceptions. The retrieval on the data set with a linear long-term trend and the full-scan extraction leads to a preference for the presence of HCN, though the abundance is suspiciously high and WFC3 data alone are known to not be reliable when inferring the presence of carbon-bearing species (Changeat et al. 2020). There is also some variation in the preferred abundance of CO|$_2$| (see Fig. 8), but the same point applies. However, by and large the detrending method has seemingly not led to any strong differences in our interpretation of the atmosphere, which is reassuring. All data sets suggest the presence of H|$_2$|O, CO, and TiO, but we note that the retrieved 10 Bar radius, cloud pressure, and temperature differ.

5.2.2 TESS and HST WFC3 G141

While there are now data with the G102 grism, we firstly consider only the TESS and HST WFC3 G141 spectra to highlight potential issues that arise from combining data sets without spectral overlap. From our retrievals, we see that the preferred atmospheric models when combining the photometric point from TESS with the HST WFC3 G141 spectrum recovered from the linear and quadratic trends are very different. One could have guessed this a priori from the vast difference in the spectral depths: in the linear case, a swift decrease in transit depth is needed at shorter wavelengths to match the TESS data point, while in the quadratic case a strong absorption is needed instead. Hence, the differences in the mean G141 depths cause discrepancies in the inferred atmospheric composition. As shown in Fig. 9, some retrievals prefer the presence of TiO while others do not. One retrieval finds evidence for CO while others prefer HCN.³ Even the H|$_2$|O abundance varies between models, which is the primary species that G141 data are sensitive to.

Given the enormous variation seen in the mean G141 transit depth in Section 5.1, and the variation in the inferred composition in these two cases, it is obvious that there is no way to tell which solution is closest to the true composition of the planet. While one might tempted to try to statistically determine which light-curve fit to the WFC3 data are preferred, and thus which atmospheric model is preferred, the results can hardly be trusted given that it is clear that the absolute transit depth of the WFC3 G141 cannot be known and so all spectra are likely to be offset from the true value.

Hence, despite being the standard methodology for characterizing exoplanet atmospheres at low-resolution, adding together data from instruments without spectral overlap and without accounting for potential offsets can clearly led to incorrect interpretations of the atmosphere.

5.2.3 Accounting for offsets without spectral overlap

To attempt to overcome potential offsets between data sets, several studies have introduced an offset parameter within the atmospheric retrievals (e.g. Luque et al. 2020; Murgas et al. 2020; Yip et al. 2021; McGruder et al. 2022; Edwards et al. 2023a; Jiang et al. 2023). However, if only one of the data sets is spectroscopic in nature, and if the data sets do not have any wavelength overlap, it is likely to be difficult for the retrieval to converge to a reliable value for this parameter. Nevertheless, we attempted to do this for the TESS and HST WFC3 G141 data for KELT-11 b. For this exercise, we again use the data sets fitted with Claret limb-darkening coefficients and the linear or quadratic long-term trend. As it is used in all four retrievals, we fix the TESS depth and instead allow the mean depth of the HST WFC3 G141 data to be shifted. To avoid over-constraining the model, the priors on the offset are set to |$\pm$|1000 ppm.

We find that the retrievals generally converge to a single solution, one containing HCN and TiO as well as H|$_2$|O and CO|$_2$|⁠. As displayed in Fig. 10, the HST WFC3 G141 are offset by each retrieval until they reach roughly the same depth, with the offsets reducing the mean depth by 300–700 ppm. Given that, as shown in Section 5.2.1, the retrievals on each different detrending of the WFC3 G141 data lead to similar solutions, and given that a broad uniform prior is used for the offset parameter, this result is expected. Hence, these retrievals appear to confidently suggest that H|$_2$|O, CO|$_2$|⁠, TiO, and HCN are present in the atmosphere of KELT-11 b. However, the need to force such a large offset suggests that the conclusions of these retrievals may not be robust and that the solutions are essentially driven by the WFC3 spectrum alone.

5.2.4 Impact of the HST WFC3 G102 data

Different methodologies for fitting the HST WFC3 G102 data lead to a spectra with similar features but very different mean depths. G102 data has spectral overlap with the data from G141 and an offset is always clearly apparent between the data sets, even when the same methodology is applied to both grisms. Hence, we now explore the impact on the recovered atmospheric chemistry. Here we again use only those fitted with linear and quadratic long-term trends, with the orbital parameters fixed, and both the splitting and full-scan extractions. In light of the clear discrepancies, we allow the mean depth of the G141 spectrum to be offset with respect to the G102 spectrum.

Our best-fit models are shown in Fig. 11. In these, we see that the chemistry is very similar across each of the retrievals. The exceptions are again some of the carbon-bearing species. For instance, the CO|$_2$| abundance profile retrieved in the linear splitting extraction case differs from the rest and lacks the preference for solutions which contain a significant (log|$_{10}$|(VMR)|$\gt $|−5) abundance of HCN.

We also conduct retrievals which also include the TESS photometric point. In these cases, the TESS photometric point is fixed and the mean depth of both the G102 and G141 spectra are allowed to vary. As shown in Fig. 12, we find that these four retrievals also converge to approximately the same solution but note that it is different to the solution without the TESS data in several ways. Namely, any evidence for HCN or AlO is removed.

Changeat et al. (2020) showed that, for the WFC3 G141 data alone, the retrieval results were very dependent upon the set-up. One might find similar results on this combined spectrum.

5.3 Retrievals on the data from Colon et al. 2021

Having noted that different pipelines led to offsets in the mean transit depth of the HST WFC3 G141 data of around 0.036 per cent (360 ppm), Colón et al. (2020) also conducted a retrieval with an offset parameter included. In this retrieval, they allowed the TESS and WFC3 data points to shift relative to the Spitzer IRAC photometry. However, they restricted this offset to be a maximum of 0.008 per cent (80 ppm) from the baseline, arguing that, based on the precision of the TESS and HST data, the offset was unlikely to be greater than this value. Therefore, they seemingly neglected the fact that their own reductions of the WFC3 G141 data had led to far larger offsets (up to 360 ppm). While no posterior distributions were shown, they stated that this retrieval led to a result that was compatible with their standard retrieval where no offset was considered.

Here also we conduct retrievals on their data sets. We perform fits to their HST WFC3 G141 data alone, as well as conducting several fits with the addition of the photometric data points from TESS and Spitzer. In the case where multiple instruments are used, we consider retrievals both with, and without, offsets between the data. We perform two retrievals with offsets: in the first we also limit the offset to be 80 ppm, which is far less the offsets found, while in the second we allow it to vary by 0.1 per cent (1000 ppm), far more than the offsets found.

We find that, in the restricted offset case, the retrieval is clearly overly constrained by the prior on the offset, as the peak in the probability distribution is close to the upper bound of the prior, as shown in Fig. 13. Furthermore, when the offset is allowed to be free, the abundances of the chemical species retrieved changes. Intriguingly, the retrievals with the data sets fixed and with a restricted offset both found evidence for AlO. However, evidence for this species disappeared once the offset between data sets was not restricted, with a higher abundance of TiO and HCN instead being used to fit the data, though we note that these abundances do not seem feasible when compared with the expected abundances in chemical equilibrium. Our retrieval on just the WFC3 data from Colón et al. (2020) also prefers using TiO to explain the data instead of AlO.

If we had restricted the offset size in Section 5.2.3, it is likely we would have found that different reduction methods lead to different compositions. Therefore, accounting for an offset clearly affects the chemistry in this case as well. Additionally, the priors utilized for the offset also impact the final conclusions.

6. DISCUSSION

In this study, we clearly demonstrate that different spectroscopic data sets of transiting exoplanets can be incompatible due to offsets, a result which is neither new nor surprising, and so once again show that multiple data sets should not be carelessly combined. The literature is littered with examples of how offsets in the mean transit depth can be easily be incurred (e.g. Stevenson et al. 2014; Chachan et al. 2019; Anisman et al. 2020; Changeat et al. 2020, 2022; Colón et al. 2020; Guo et al. 2020; Mansfield et al. 2021; Mugnai et al. 2021; Yip et al. 2021; Edwards et al. 2023b) and yet the vast majority of studies (e.g. Sing et al. 2016; Wakeford et al. 2017; Benneke et al. 2019; Pinhas et al. 2019; Welbanks et al. 2019; Alam et al. 2020; Colón et al. 2020; Brande et al. 2022; Kreidberg et al. 2022; Mikal-Evans et al. 2023) ignore these potential discrepancies in pursuit of a wider wavelength coverage and thus a seemingly better understanding of the chemistry of the planet(s) studied. By highlighting the case of KELT-11 b, it is hoped that future works will, at the very least, consider the possibility of that data sets can be incompatible due to offsets.

Our analysis of the HST WFC3 data of KELT-11 b shows that the assumptions made during the reduction and fitting of data can change the mean transit depth. While the impact on the atmospheric chemistry is minimal if only a single instrument is used, adding in additional data sets without accounting for a potential offset leads to disparate solutions. Given that it is unclear which detrending and fitting techniques are ‘best’, it becomes difficult to robustly infer the atmospheric chemistry. Here, fitting for an offset between the TESS and HST WFC3 G141 data for various WFC3 reductions leads to roughly consistent results. However, once the G102 data are analysed it seems that this solution is incorrect as the preferred models predict a steep increase in transit depth which is not seen in the G102 data (see Figs 10 and 11). Therefore, fitting for an offset between data sets with no spectral overlap cannot be considered completely reliable. On the other hand, our retrievals on the combined spectra from various reductions of G102 and G141 data led to consistent solutions, though they too may be unreliable if the assumptions made in the retrievals are not suitable or if the planet’s transmission spectrum varies temporally.

The issue of data incompatibility when combining data sets is not going to disappear any time soon. Even with the wider wavelength coverage offered by the instruments onboard JWST, they still do not generally probe the entirety of the spectral range desired for constraining the majority of the key molecules. Therefore, data sets from different JWST instruments will undoubtedly be combined, with data from other instrument such as HST STIS or WFC3 also being thrown into the mix. Indeed, studies have already combined data from these facilities and offsets have been seen (e.g. Constantinou et al. 2023; Dyrek et al. 2024). Furthermore, Holmberg & Madhusudhan (2023) showed that different treatments of JWST NIRISS data (e.g. limb-darkening) also affected the transmission spectrum. Even with instruments which provide a wide wavelength coverage in a single observation [e.g. JWST NIRSpec PRISM (Birkmann et al. 2022), Twinkle (Edwards et al. 2018; Stotesbury et al. 2022), and Ariel (Tinetti et al. 2018, 2021)], there may still be the temptation to add further data sets to enhance them such as including data at ultraviolet wavelengths from HST WFC3 UVIS (Ultra-violet and Visible light channel) or at the longer wavelengths covered by JWST MIRI (Mid-Infrared Instrument) LRS (Low Resolution Spectrometer) (Kendrew et al. 2015).

For observations with a higher resolution or signal-to-noise ratio (SNR), offsets should be more readily apparent as retrievals may be unable to fit the distinct spectral features. However, for data with a lower resolution or SNR, offsets may yet wreak havoc. These lower SNR observations are likely to be of smaller planets, including rocky worlds within the habitable zone of their star, and so any conclusions drawn from these studies are likely to receive wide-spread attention. As such it is imperative that the possibility of incompatibilities are seriously considered in any such study, and that the measures listed below are taken in an attempt to nullify concerns that the results are biased by issues of data set compatibility.

If one is to combine instruments, we lay out several key steps and checks that can be undertaken in an attempt to maximize the likelihood that they will be compatible or to check the degree to which the results depend upon the assumption of compatibility. Some of these are HST WFC3 specific, while others are more general. These recommendations are:

• Extraction of WFC3 scanning data: Choosing whether or not to utilize the non-destructive reads from WFC3 scanning data may affect the extracted flux. Given that charge-trapping can clearly be seen here in the raw data when subtracting the first read from the last, utilizing the non-destructive reads is the prudent choice, especially as using the non-destructive reads also helps to avoid including background sources contributing to the extracted spectrum.

• Fitting parameters: The same orbital parameters and limb-darkening law should be used to analyse all data sets considered in the study as changing them can cause offsets in the mean transit depth. However, it should not be assumed that the literature values for the orbital parameters guarantee a correct spectrum.

• Detrending technique: Different assumptions for the removal of systematics, including using different pipelines, can lead to offsets in the mean transit depth recovered. The corner plots of these fits should be checked for correlations between the transit depth and the parameters used to model the systematics. For HST WFC3 data, it is likely that the linear trend will usually give the least correlation with the transit depth, but of course this does not mean it statistically provides the best-fit nor that it allows for the recovery of the absolute transit depth.

• The impact of new data: Different data sets are generally combined to provide sensitivity to additionally species. Given the results of this study and others, this approach should generally be discouraged unless extreme caution is used. If one is to combine data sets, exploring how the addition of extra data impacts the retrieved atmospheric abundances is a prudent strategy to undertake. Doing so allows one to specific which data sets are responsible for the detection of different species and to caveat detections based on compatibility. For instance, Pluriel et al. (2020) showed that carbon-bearing species are not detected when only analysing the WFC3 G141 transmission spectrum of KELT-7 b, but are seemingly constrained with the addition of a photometric point from Spitzer. In this case, the detection hinges on the ability to reliably co-add the data sets. Additionally, the effect of using different detrending or fitting assumptions, or reduction pipelines, on the retrieved atmospheric composition should also be considered.

• Fitting for offsets: Including a parameter for an offset between data sets can be a useful approach to understanding the effect of potential incompatibilities. However, the offset parameter should not be overly constrained, particularly not to a level which is smaller than the offsets seen between different reductions or analyses of the data. Furthermore, the reliability of fitting for an offset is reduced if there is no spectral overlap, particularly if photometric data are being used. Moreover, the preferred value for the offset parameter can be driven by the data set with the most data points, and thus may have a preference to simply converge to the solution of this driving data set. Therefore, it cannot be assumed to be a perfect solution to the problem of data set incompatibility and far from guarantees uncovering the absolute transit or eclipse depth.

• Spectral overlap: If the data sets that are being combined have spectral overlap then one can be far more certain as to whether the data appears to be compatible or not. However, we note that temporal variations, either in the planet or star, cannot be ruled out in this case and so data sets may still not be compatible even if the recovered depths are similar.

• Recovering absolute transit depths: To perfectly recover the transit depth, one must entirely remove additional flux contributions (background, dark-current, etc.) as well as flawlessly modelling instrument effects which alter the flux levels (e.g. persistence). If one does a good job at reducing and fitting the data, the true transit depth should be within the 1|$\sigma$| uncertainties on the measured transit depth to 68 per cent confidence. However, in this work we have shown that different fitting methods lead to variations in the mean depth that are far greater than the 1|$\sigma$| uncertainties. Therefore, when analysing HST WFC3 data one cannot be sure that the measured depths are in fact close to the true value. Instead, one can only have confidence in the relative depth between different HST WFC3 spectral bins (i.e. the spectral features). As different JWST pipelines or analysis techniques are also finding a similar issue (e.g. Holmberg & Madhusudhan 2023), it is clear more work is needed before we can trust the absolute measurements of transit depth from either facility.

• Retrieval set-up: Finally, we reiterate the advice of Changeat et al. (2020) as it is clear with retrievals that what you get out often depends upon what you put in. Fully exploring the robustness of any ‘detection’ is clearly the most prudent approach. Furthermore, as mentioned earlier for the fitting of an offset, the choice of priors in these retrievals is also crucial.

7. CONCLUSIONS

As differences in the detrending or analysis technique change the mean transit depth, confidently recovering absolute transit depths with data from HST WFC3 is clearly not possible. Yet if one is to combine together data from different instruments, then it is imperative that these data sets recover the absolute transit depth. We have shown that there is a wide range of uncertainty on the transit depth of KELT-11 b from HST data that precludes a simple co-addition of spectra from different instruments. Therefore, combining spectra from different instruments requires us to consider offsets between them. However, our work shows that combining WFC3 with a photometric point from TESS or Spitzer can lead to incorrect interpretations of the atmosphere, even with an offset parameter in the retrieval. Hence, to robustly combine data sets, each must independently resolve at least one spectral feature, with spectral overlap between the data sets also helping to robustly recover the atmospheric chemistry.

With such a myriad of potential issues to overcome, the transiting exoplanet community faces a tough challenge in the coming years. Unless we are thorough in our analyses, as well as honest about the potential biases and flaws in our work, we will not be able to uncover the true nature of exoplanetary atmospheres.

DATA AVAILABILITY

Our work is based on observations made with the NASA/ESA Hubble Space Telescope. These publicly available observations were taken as part of proposal GO-15225 and GO-15926, both led by Knicole Colón (Colón 2017; Colón et al. 2019). These were obtained from the Hubble Archive which is part of the Barbara A. Mikulski Archive for Space Telescopes (MAST). Additionally, this paper includes data collected by the TESS mission which is funded by the NASA Explorer Program. TESS data are also publicly available via the Barbara A. MAST. All these data sets can be accessed via DOI: 10.17909/x8ap-3987. We are thankful to those who operate this archive, the public nature of which increases scientific productivity and accessibility (Peek et al. 2019). The project also received funding from the Science and Technology Funding Council (STFC) grants ST/S002634/1 and ST/T001836/1.

Computing: We acknowledge the availability and support from the High-Performance Computing platforms (HPC) from DIRAC, and OzSTAR, which provided the computing resources necessary to perform this work. This work utilized the Cambridge Service for Data-Driven Discovery (CSD3), part of which is operated by the University of Cambridge Research Computing on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). The DiRAC component of CSD3 was funded by BEIS capital funding via STFC capital grants ST/P002307/1 and ST/R002452/1 and STFC operations grant ST/R00689X/1. DiRAC is part of the National e-Infrastructure. This work utilized the OzSTAR national facility at Swinburne University of Technology. The OzSTAR program receives funding in part from the Astronomy National Collaborative Research Infrastructure Strategy (NCRIS) allocation provided by the Australian Government.

Software:Iraclis (Tsiaras et al. 2016c), TauREx3 (Al-Refaie et al. 2021), pylightcurve (Tsiaras et al. 2016b), ExoTETHyS (Morello et al. 2020), Astropy (Astropy Collaboration 2018), h5py (Collette 2013), emcee (Foreman-Mackey et al. 2013), Matplotlib (Hunter 2007), Multinest (Feroz et al. 2009), Pandas (McKinney 2011), Numpy (Oliphant 2006), SciPy (Virtanen et al. 2020), and corner (Foreman-Mackey 2016).

Facilities: Hubble Space Telescope (HST), Transiting Exoplanet Survey Satellite (TESS), and Spitzer Space Telescope.

ACKNOWLEDGEMENTS

QC is funded by the European Space Agency (ESA) under the 2022 ESA Research Fellowship Program. BE and AT received travel support for this work under the ESA Science Faculty Funds, funding reference ESA-SCI-SC-LE117.

Footnotes

References

APPENDIX A: ADDITIONAL FIGURES, TABLES OF PARAMETERS, AND WHITE LIGHT-CURVE FITS

In this Appendix, we present all the white light-curve fits from Section 2.2 as well as additional information about the observations.

A raw spatial scan from each observation sequence is shown in Fig. A1. The scans have the same spatial length but different spectral width. The colour scale is consistent for both images and the drop in flux at shorter wavelengths is notable for the G102 scan. It is this strong change in flux that seemingly means that the divide-by-white method cannot accurately reproduce the systematics at these wavelengths.

As noted in the main text, the white light-curve fits often had non-Gaussian residuals. During the G102 visit in particular, there were shifts in the position of the spatial scan on the detector, as shown in Fig. A2. While these were not so large as to restrict the possibility of analysing this data, they could be inducing additional systematics and/or reducing the precision achievable (Stevenson & Fowler 2019). Whatever the cause, the model in our standard fit is clearly not capable of fully removing systematics.

As noted in the main text, Iraclis offers two methods from extracting the flux. The process of each of these techniques is shown in Fig. A3. These lead to light curves with different flux levels, as shown in Fig. A4.

The parameters used for our light-curve fits, and in the retrievals, are given in Table A1. Figs A5 and A6 show all the white light-curve fits for the G102 and G141 observations, respectively. Despite the larger shifts in the position of the spectrum on the detector, the residuals of the G102 fits are generally more Gaussian. Additionally, using the quadratic fit usually lead to a fit where the AC function was lower (i.e. the residuals were more Gaussian), potentially indicating a better fit. However, it should be noted that, as discussed in Section 5.1.1, the quadratic fit leads to a far more variable white light-curve depth when other parameters were changed. Furthermore, the parameters for modelling the systematics are generally correlated with the transit depth (see Fig. 7).

We report the spectra obtain from our standard fits with the linear and quadratic long-term trends in Table A2, along with the limb-darkening coefficients used.