On the Compatibility of Ground-based and Space-based Data: WASP-96 b, an Example^*

The HST data of WASP-96 b were acquired by proposal 15469 led by Nikolay Nikolov and were taken in December 2018. We obtained the raw spatially scanned spectroscopic images from the Mikulski Archive for Space Telescopes ³ and used Iraclis, ⁴ a specialized, open-source software for the analysis of WFC3 scanning observations (Tsiaras et al. 2016b). The reduction process included the following steps: zero-read subtraction, reference pixels correction, nonlinearity correction, dark current subtraction, gain conversion, sky background subtraction, calibration, flat-field correction, and corrections for bad pixels and cosmic rays. For a detailed description of these steps, we refer the reader to Tsiaras et al. (2016b).

The reduced spatially scanned spectroscopic images were then used to extract the white and spectral light curves. We then discarded the first orbit of the visit as it presents stronger wavelength-dependent ramps. For the fitting of the white light curves, the only free parameters were the midtransit time and planet-to-star ratio. We did not fit for the inclination or reduced semimajor axis, as ingress and egress were not observed in each data set. The limb-darkening coefficients were selected from using the models of Claret et al. (2012, 2013) and using the stellar parameters from Hellier et al. (2014). The fitted white and spectral light curves for the G102 and G141 transmission observations are shown in Figures 1 and 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Additionally, two transits of WASP-96 b had been observed with Spitzer IRAC (program ID 14255). We used the Transit Light Curve Detrending Long Short-Term Memory (TLCD-LSTM) pipeline from Morvan et al. (2020) to detrend and fit the Spitzer data from the outer transit and centroids movements. The architecture is the same as in the aforementioned study except for the number of hidden units reduced to 64 and dropout rate set to 0.1 in order to prevent over fitting on the out-of-transit. Figure 3 shows the detrended light curves and the best-fit model to the data. For both data sets, the only free transit parameters were the planet-to-star radius ratio and the transit midtime, with the other model parameters fixed to those in Table 1.

Zoom In Zoom Out Reset image size

We also fitted the detrended light curve while allowing the planet semimajor axis to star radius ratio (a/R_s) and inclination (i) to vary. For both channels, the retrieved values for the transit depth and epoch remain very close to the ones found with the two orbital parameters fixed to the values from Nikolov et al. (2018). Furthermore, values retrieved for the semimajor axis and inclination shown in Table 2 are compatible with those from Nikolov et al. (2018) to 1σ.

Table 2. Comparison of Orbital Parameters between Different Data Sets

Parameter	Nikolov et al. (2018)	TESS	IRAC C1	IRAC C2
a/Rs	8.84 ± 0.1	8.85 ${}_{-0.10}^{+0.62}$	8.66 ${}_{-0.12}^{+0.14}$	8.78 ${}_{-0.06}^{+0.06}$
i [deg]	85.14 ± 0.2	85.55 ${}_{-0.34}^{+0.39}$	85.36 ${}_{-0.13}^{+0.14}$	85.21 ${}_{-0.06}^{+0.07}$

Download table as:  ASCIITypeset image

Keeping exoplanet transit ephemeris fresh is crucial for allowing further atmospheric characterization. Here, we use our Hubble observations along with data from Spitzer and the Transiting Exoplanet Survey Satellite (TESS, Ricker et al. (2014)) to update the orbital period and transit epoch of WASP-96 b. TESS data are publicly available through the MAST archive and we follow the procedure from Edwards et al. (2020b) to download, clean, and fit the 2 minute cadence Pre-search Data Conditioning (PDC) light curves (Smith et al. 2012; Stumpe et al. 2012, 2014). WASP-96 b was observed during Sector 2, and after excluding bad data, seven transits were recovered. These were individually fit and are shown in Figure 4. In our main analysis, the only free parameters were the transit midtime and the planet-to-star radius ratio. However, we also performed a fitting of all the TESS data where further orbital parameters, namely the inclination (i) and planet semimajor axis to star radius ratio (a/R_s) were allowed to vary. The results are summarized in Table 2 along with similar fits for the Spitzer data.

Zoom In Zoom Out Reset image size

The retrieval of the transmission spectrum was performed using the publicly available retrieval suite TauREx 3 (Al-Refaie et al. 2019). ⁵ For the star parameters and the planet mass, we used the values from Hellier et al. (2014). In our runs we assumed that WASP-96 b possesses a primary atmosphere with a ratio He/H₂ = 0.17 (i.e., solar abundance). To this, we added trace gases and included the molecular opacities from the ExoMol (Tennyson et al. 2016), HITRAN (Gordon et al. 2016), and HITEMP (Rothman & Gordon 2014) databases for the following molecules: H₂O (Polyansky et al. 2018), CH₄ (Yurchenko et al. 2017), CO (Li et al. 2015), CO₂ (Rothman et al. 2010) NH₃ (Coles et al. 2019), K, and Na (Kramida et al. 2013). The line-broadened profiles for the resonance doublets of Na and K are computed using Allard et al. (2016) and Allard et al. (2019). For the clouds, we use gray opaque clouds and the Mie cloud model from Lee et al. (2013). On top of this, we also included collision-induced absorption (CIA) from H₂ to H₂ (Abel et al. 2011; Fletcher et al. 2018) and H₂-He (Abel et al. 2012) as well as Rayleigh scattering for all molecules. We assumed isothermal and isochemical profiles throughout all our retrievals.

In our retrieval analysis, we used log uniform priors for all parameters as described in Table 3. Finally, we explored the parameter space using the nested sampling algorithm Multinest (Feroz et al. 2009) with 750 live points and an evidence tolerance of 0.5.

The top panel of Figure 5 shows the best-fit spectrum from our retrieval based on G102 and G141 data only. The best-fit model contains modulations that correspond to multiple absorption features of water, clearly indicating the presence of the molecule in the atmosphere. Our free Bayesian retrieval analysis recovered a water abundance of log(H₂O) = −3.08${}_{-1.81}^{+1.08}$, consistent with predictions from equilibrium chemistry models (Agúndez et al. 2012; Stock et al. 2018; Woitke et al. 2018). Assuming an isothermal temperature profile, our retrieval determines a temperature T = ${609}_{-120}^{+173}$ K, which is lower than the equilibrium temperature of the planet and much lower than the temperature obtained by Nikolov et al. (2018). Simulations have demonstrated that the complexity of the terminator limb, which includes three-dimensional asymmetries in the chemical and thermal structure of the terminator region, often affects the absolute retrieved temperature and can have biases (Caldas et al. 2019; MacDonald et al. 2020; Pluriel et al. 2020a). A lower than expected temperature is often retrieved from HST WFC3 data (e.g., Skaf et al. 2020). Similarly, narrow wavelength coverage can result in wrong estimations of the atmospheric temperatures, which can be unstable if a single molecular band is probed (e.g., Rocchetto et al. 2016; Tsiaras et al. 2018; Pinhas et al. 2019). On the other hand, our analysis also shows a gray cloud top pressure at log(P_clouds) = ${4.39}_{-1.21}^{+0.95}$, which seems to suggest a relatively clear atmosphere, confirming the findings of Nikolov et al. (2018). We did not fit for more complicated scattering models, as the HST wavelength range does not cover a sufficiently large wavelength range (in particular shorter wavelength) to constrain atmospheric scattering model (see Appendix A for results from fitting a more complicated scattering model). Due to the weak absorptions, especially compared to water, of other molecules in the wavelengths considered, we were not able to determine the molecular abundance of NH₃, CO, and CH₄. For NH₃, we were able to extract a 1σ upper bound of log(NH₃)_upper = −6.51. The posterior distributions for this retrieval are shown in Figure 6.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The top panel of Figure 5 shows the HST spectra (G102 and G141) analyzed here. We also display the raw spectrum obtained with the VLT in Nikolov et al. (2018). With no correction, we can immediately observe that the two sets of spectra are not compatible. In particular, a significant offset between the the ground-based (orange) and space-based data (blue and dark blue) at around 0.8 μm. There are a number of potential sources for the differences seen here: variations in the stellar properties, instrument systematics, differences in the reduction pipelines, telluric corrections and the use of different orbital parameters or limb-darkening coefficients.

An imperfect correction of instrument systematics certainly has the potential to significantly alter the recovered transit depth, the best-fit models of the systematics, namely the orbital and long period ramps, are shown in Figure 7. The corrections applied to some exposures are greater than the offset seen between the HST and VLT data, which may explain the offset observed. However, the ramps seen in the G102 and G141 data are very different yet they have both been fitted such that the final data products are seemingly in good agreement. The observations had different exposure times (179.05 s; 156.70 s), scan lengths (241; 360) and scan rates (0013 s⁻¹; 0022 s⁻¹), so different systematics are to be expected.

Zoom In Zoom Out Reset image size

For both HST observations we fitted a linear long-term trend, in line with many previous studies. Guo et al. (2020) suspected, from visual inspection, that the trends seen in the WFC3 data of HD 97658b deviated from a linear fit. Thus, they experimented with a number of different trends (quadratic, exponential, logarithmic) finding that, while the light-curve depths recovered from each observation were roughly consistent with one another, the chosen trend affected the transit depth recovered from the white light curves. They noted that the quadratic model produced the most uncertain depth as well as being the most discrepant between visits. Agol et al. (2010) also saw a bias when fitting Spitzer data with a quadratic trend. Here, visual examination of the raw light curves does not indicate that long-term trends are nonlinear, thus we do not explore different detrending modules.

In Nikolov et al. (2018), the authors noted that the derived depths from their two visits had marginal disagreement (1.4σ). Said difference was 720 ppm but the authors noted that this level of variation is consistent with the photometric variability of the star, which is associated with active regions on its surface, of 920 ppm. The variation in the stellar flux could therefore be the cause of much of the offset seen here. We note that this could also be causing differences between the G102 and G141 observations though the derived depths are within 1σ. These were taken 10 days apart while the VLT data sets were acquired with a gap of 24 days.

The effect of combining different instruments is studied in the literature. Alexoudi et al. (2018) highlighted the importance of the using the correct orbital parameters for data covering visible wavelengths: otherwise slopes can be induced, or removed. Yip et al. (2020) highlight the danger of combining data from HST and the Spitzer Space Telescope, which are not sharing a common baseline and where the information redundancy for carbon-based species is limited. Despite this, studies often combine these data (e.g., Sing et al. 2016; Pinhas et al. 2019). The combination of individual analyses from two different sets of observations may not necessarily agree with each other and should be approached with care.

The case we presented here is an obvious example of when the offset can be visually inspected but it is not the first study to introduce an offset to a data set. Bruno et al. (2020) showed that, for the active star WASP-52, stellar spots could create incoherent observations. The team combined transit data from HST/STIS, HST/WFC3, and Spitzer/IRAC of WASP-52 b and corrected for the offset in HST/WFC3 by accounting for the effect of stellar spots before retrieving the planet's atmospheric composition based on the corrected observations.

Kirk et al. (2019) explored the effect of stellar activity on the atmospheric retrieval of WASP-39 b when combining different data sets. In that study, the infrared data came from HST WFC3 G102 and G141, while the optical data sets were from the ground-based ACAM instrument on the 4.2 m William Herschel Telescope or from HST STIS, providing continuous coverage from 0.4 to 1.6 μm. They found no noticeable offset and little difference in the retrieved parameters when accounting for stellar activity, suggesting that, for WASP-39 b, the data sets could be compatible.

Other studies have applied, or fitted for, offsets without wavelength overlap. In a study of WASP-74 b with photometry from ground-based instruments along with HST/WFC3 and Spitzer/IRAC, Luque et al. (2020) fitted for an offset for the HST data, discovering a best-fit value of 434 or 615 ppm, depending upon the additional data used. Wilson et al. (2020) obtained FORS2 data of WASP-103 b and combined it with other ground-based data as well as data from HST and Spitzer. They applied a small offset to match the GMOS and FORS data before including an offset parameter in their retrievals to account for any further discrepancies. Meanwhile, in their study of HAT-P-12 b, Yan et al. (2020) attempted to fit for two offsets to improve their fitting of HST WFC3, HST STIS, and LBT data but found it did not lead to solutions that were statistically more valid. Finally, when fitting for an offset in the case of WASP-69 b between HST WFC3 and OSIRIS observations, a value of 479 or 618 ppm was recovered by Murgas et al. (2020) with/without also fitting for a spot correction. They noted the cause could be biases due to instrument systematics, but that largest semiamplitude of the WASP photometry, 13 mmag, would result in a flux variation of 2.4%, which could account for the 2.2% increase in flux between the observations.

However, there have been other instances in the literature where ground-based and space-based observations are combined and analyzed without wavelength overlap nor offset correction. For example, Danielski et al. (2014) combined observations from NASA Infrared Telescope Facility (IRTF)/SpeX instrument and WFC3 observations on HD 189733b, many others followed similar trends (e.g., Bean et al. 2013; Mancini et al. 2013; Stevenson et al. 2016; Sotzen et al. 2020). While these data sets could be compatible, there is no guarantee and this should not be taken for granted.

In our case, without assigning a specific cause, we attempted to correct for the offset by fitting a single flat offset parameter in a combined retrieval. The parameter applies a shift to the entire VLT spectrum vertically to create a coherent observation with the HST data set. This choice of correction does not necessarily represent the complexity of the problem here and it is not guaranteed that the instrument systematics are wavelength-independent.

The best-fit spectrum is plotted in the bottom panel of Figure 5. Our combined retrieval unveiled an offset of 1198${}_{-102}^{+135}$ ppm between the two instruments, as shown in the posterior distributions in Figure 8. The increased wavelength coverage (continuous coverage from 0.4 to 1.6 μm) allows us to fit for more complicated cloud models and probe the presence of two species (H₂O and Na). The extension to the visible also helped to provide better constraint on the temperature on the terminator. The retrieved temperature (T = 954${}_{-195}^{+198}$ K) is close to the expected terminator temperature from a linear trend derived in Skaf et al. (2020), given the equilibrium temperature of the planet. An interesting potential explanation for the large temperatures difference between the 1700 K inferred from the VLT only data and the one retrieved from the HST only data could be that the signal from those two molecules comes from different regions of the terminator. Caldas et al. (2019), taking the example of H₂O and CO, predicted that some molecules could have a signal from the dayside part of the terminator region only or inversely the nightside part only.

Zoom In Zoom Out Reset image size

The retrieved abundances for H₂O and Na in our combined retrieval are within 2σ agreement with the individual analysis reported in our work and Nikolov et al. (2018; see Table 4 for a comparison between retrieval scenarios). The consistency in the water abundances indicates that the water feature may be stable enough to be retrieved accurately with HST, even in the case of the low terminator temperature recovered in Section 3.1. Such a result is expected given the strong features of H₂O in the WFC3 range and lack thereof in the visible.

Table 4. Comparison between Different Retrieved Quantities for Three Different Scenarios (HST Only, VLT Only, and HST+VLT)

Parameter	HST	(Nikolov et al.)	HST+VLT
log(H2O)	$-{3.08}_{-1.81}^{+1.08}$	N/A	$-{3.65}_{-0.94}^{+0.90}$
log(Na)	N/A	$-{5.1}_{-0.4}^{+0.6}$	$-{3.88}_{-0.82}^{+1.05}$
Rp [RJ ]	${1.22}_{-0.01}^{+0.01}$	N/A	${1.21}_{-0.01}^{+0.01}$
T [K]	${609}_{-120}^{+173}$	${1710}_{-200}^{+150}$	954 ${}_{-195}^{+198}$

Note. Results from VLT alone are reproduced from Nikolov et al. (2018). We have omitted quantities that could not be constrained by any of the scenarios.

Download table as:  ASCIITypeset image

Based on the retrieval result and the visually compatible observations in Figure 5, it may be tempting to conclude that the correction has been successful. However, we would like to emphasize here that this kind of correction is an ad hoc solution to the problem and does not contain any theoretical support. Any conclusion drawn from these kinds of combined observations should be treated lightly and comparisons to model fitting on single data sets should be made.

A retrieval study of 10 hot Jupiters by Pinhas et al. (2019) found that optical data played a significant role in ensuring that a reliable constraint on the abundances derived from infrared data could be placed. When combining optical and infrared data of HD 209458b from HST (STIS + G141), they found the constraints on water to be narrowed by a factor of 3. However, the abundances retrieved in each case were drastically different: log(H₂O)_WFC3 = −3.3${}_{-0.75}^{+0.80}$ and log(H₂O) _STIS+WFC3 = −4.66${}_{-0.30}^{+0.39}$.

In contrast, the water abundances recovered here, with and without optical data, are in good agreement with one another (well within 1σ). This may well be due to the addition of the G102 grism, which has not been used to observe HD 209458b, highlighting that it has the capability to provide excellent constraints on the water abundance when combined with G141. Additionally, while here we focus on a case where ground-based and space-based instruments demonstrate an offset, such a discrepancy could also occur in space-based data sets from different instruments. Combining HST STIS, HST WFC3 G141, and Spitzer IRAC has become commonplace in the field (e.g., Sing et al. 2016). Given that there is no wavelength overlap in these studies, there is a risk of offsets occurring that could bias the results of subsequent atmospheric retrievals. For instances, the Spitzer IRAC bands cover spectral regions where carbon-bearing molecules such as CH₄, CO, and CO₂ absorb, which could be biased when the instrument is not well-calibrated, leading to wrong estimates of the C/O ratio. The G102 grism remains an underutilized instrument for exoplanet spectroscopy but would offer extra confidence that the STIS and G141 data sets are compatible by providing wavelength overlap with both.

In addition to VLT observation, we have also explored the scenario when observations from Spitzer and TESS are added to the HST data. Figure 9 shows the TESS and Spitzer data plotted alongside the best-fit spectrum from the HST and corrected VLT retrieval from the previous section. The two Spitzer points (3.6 μm and 4.5 μm) are within 1σ of the best-fit solution. Hence, it is possible that Spitzer is consistent with the other instruments. However, we remain cautious, and as there is no wavelength overlap for Spitzer, we cannot use the methodology employed for the correction of the VLT data. Additionally, given the size of the error bars on the Spitzer data, little additional spectral information would be gained. The TESS point, on the other hand, is about 2σ away from the solution. Seven transits were observed with TESS and the depth of each of these is shown in Figure 10. While several of the individual observations are consistent with the HST + VLT best-fit (to 1σ), the weighted mean of these observations is larger than the model (orange bar). The average transit depth is larger than that obtained with HST, and the Spitzer points are seemingly consistent with the HST data. This provides further indications that the source of the offset seen between data sets may be caused by the reduction and analysis of the VLT observations. The transit depth data for all instruments is given in Table 6.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We found that the observed HST and TESS transits were just outside the 1σ literature ephemeris. Hence, we refined the period and reference midtransit time using the original ephemeris from Hellier et al. (2014), and the new data analyzed. We determined the ephemeris of WASP-96 b to be P = 3.42525650 ± 0.00000043 days and T₀ = 2457665.84332 ± 0.00014 BJD_TDB, where P is the planet's period, T₀ is the reference midtime of the transit and BJD_TDB is the barycentric Julian date in the barycentric dynamical frame. Our derived period is 0.32 s shorter than that from Hellier et al. (2014) and we improved the precision of the period by a factor of 6, thus reducing the current uncertainty on the transit time. The observed minus calculated plots are given in Figure 11 and all transit midtimes used for the fitting are listed in Table 5. We note that the Spitzer points are both seemingly poor fits to the trend, while the second HST observation (G141 grism) gives extremely tight bounds on the midtime despite the gaps within the light curve. We attempted an ephemeris fit without the HST observations and found little change in the period. TESS will soon reobserve WASP-96 b, which will allow for further refinement of its period. The midtimes have been uploaded to ExoClock ⁶ , an initiative to ensure transiting planets are regularly followed up, keeping their ephemeris up-to-date for the ESA Ariel mission (Tinetti et al. 2018; Edwards et al. 2019).

Zoom In Zoom Out Reset image size

Water appears to be ubiquitous in exoplanetary atmospheres. To understand the distribution in the abundance of this molecule, we compared the retrieved water (log) abundance (log(H₂O) of exoplanets with similar sizes (±0.5 R_J ) and masses (±0.2 M_J ) against their respective equilibrium temperature (see the top panel of Figure 12) using data from Tsiaras et al. (2018), Pinhas et al. (2019), and Skaf et al. (2020). While these trends are not statistically significant yet, at temperatures of 1200–1400 K, planets appear to generally have a high water abundance. These planets are also closest in terms of size and mass to WASP-96 b (lower panel of Figure 12). An HST WFC3 study of Kepler-51 b and d, which have radii of 0.61 R_J and 0.84 R_J , respectively, uncovered flat spectra for both, with no discernible atmospheric features (Libby-Roberts et al. 2020). However, these planets have noticeably cooler temperatures than WASP-96 b (400–550 K) and very low densities (<0.07 g cm⁻³). While the addition of our detection adds weight toward the tendency of having water rich atmosphere for this class of hot planets, more objects of a similar class are needed in order to statistically verify this claim.

Zoom In Zoom Out Reset image size