Awesome SOSS: Atmospheric Characterisation of WASP-96 b using the JWST Early Release Observations

The newly operational JWST offers the potential to study the atmospheres of distant worlds with precision that has not been achieved before. One of the first exoplanets observed by JWST in the summer of 2022 was WASP-96 b, a hot-Saturn orbiting a G8 star. As part of the Early Release Observations program, one transit of WASP-96 b was observed with NIRISS/SOSS to capture its transmission spectrum from 0.6-2.85 microns. In this work, we utilise four retrieval frameworks to report precise and robust measurements of WASP-96 b's atmospheric composition. We constrain the logarithmic volume mixing ratios of multiple chemical species in its atmosphere, including: H$_2$O = $-3.59 ^{+ 0.35 }_{- 0.35 }$, CO$_2$ = $-4.38 ^{+ 0.47 }_{- 0.57 }$ and K = $-8.04 ^{+ 1.22 }_{- 1.71 }$. Notably, our results offer a first abundance constraint on potassium in WASP-96 b's atmosphere, and important inferences on carbon-bearing species such as CO$_2$ and CO. Our short wavelength NIRISS/SOSS data are best explained by the presence of an enhanced Rayleigh scattering slope, despite previous inferences of a clear atmosphere - although we find no evidence for a grey cloud deck. Finally, we explore the data resolution required to appropriately interpret observations using NIRISS/SOSS. We find that our inferences are robust against different binning schemes. That is, from low $R = 125$ to the native resolution of the instrument, the bulk atmospheric properties of the planet are consistent. Our systematic analysis of these exquisite observations demonstrates the power of NIRISS/SOSS to detect and constrain multiple molecular and atomic species in the atmospheres of hot giant planets.


INTRODUCTION
After launch in December 2021, a careful journey to L2, and a successful commissioning period, JWST finally began its long-awaited science operations on July 12, 2022. It is a credit to how far the field of exoplanet astronomy has progressed in the past couple of ★ E-mail: jake.taylor@physics.ox.ac.uk † NASA Sagan Fellow ‡ NASA Sagan Fellow decades that some of the very first observations with this revolutionary new observatory were of transiting exoplanets. JWST vastly extends the wavelength range with which exoplanet atmospheres can be probed from space. Previous state-of-the-art observations with the Hubble Space Telescope (HST) probed UV, optical, and near infrared wavelengths out to 1.7 µm. The Spitzer Space Telescope enabled predominantly photometric measurements further into the infrared, although only the bluest bandpasses at 3.6 and 4.5 µm remained in operation after the coolant ran out in 2009. In combination, the four instruments on board the JWST allows for the study of exoplanet atmospheres from 0.6-28 µm; enabling the characterization of their atmospheres at wavelengths never seen before, broadening our discovery space into uncharted territories. This is evident from the first results from the Early Release Science (ERS) observations of the hot-Jupiter WASP-39 b which yielded the first ever detections of CO 2 and SO 2 (The JWST Transiting Exoplanet Community Early Release Science Team et al. 2022;Alderson et al. 2022;Rustamkulov et al. 2022;Tsai et al. 2022). The Early Release Observations (ERO) program was designed to provide the astronomical community with publicly available data, touching on many of the key science objectives of JWST, immediately after the end of the commissioning period (Pontoppidan et al. 2022). For the exoplanet portion of the ERO program, transits of two hot-Saturns, WASP-96 b (Hellier et al. 2014) and HAT-P-18 b (Hartman et al. 2011;Fu et al. 2022) were observed with the Single Object Slitless Spectroscopy (SOSS) mode (Albert et al. submitted) of the Near Infrared Imager and Slitless Spectrograph (NIRISS) instrument (Doyon et al. submitted).
This work focuses on WASP-96 b, an inflated hot-Saturn exoplanet with a mass of 0.498± 0.03 M J , a radius of 1.2±0.06 R J , and an equilibrium temperature of ∼1300 K. It orbits a G8 star in the constellation of Phoenix with an orbital period of 3.4 d. The planet's short orbital period, combined with its low density makes it an ideal candidate for atmospheric spectroscopy. Indeed, there have been multiple previous atmospheric studies of this planet, with the first being a groundbased spectrum using VLT/FORS2 by Nikolov et al. (2018). This observation covered a spectral range of 0.36-0.82 µm using spectroscopic bins with widths of 0.016 µm. This spectroscopic precision allowed for the measurement of the pressure-broadened sodium D line with wings reported to cover six atmospheric pressure scale heights (Nikolov et al. 2018). The visibility of the sodium wings suggested that there are no clouds or hazes obscuring them in the atmosphere at the pressure ranges probed in transmission (Fortney 2005). This was supported by their atmospheric modeling which finds no evidence for additional opacity due to clouds. Nikolov et al. (2018) further conclude that the abundance of Na is consistent with the measured stellar value.
Prior to the commissioning period of JWST, Nikolov et al. (2022) published the transmission spectrum of WASP-96 b using HST and Spitzer, providing the first look at the infrared spectrum of the planet using space-based instrumentation. To explore the atmospheric constituents in detail they couple the HST and Spitzer observations with the previous VLT observations. They find an offset between the space and ground-based data, consistent with what was found by Yip et al. (2021) who explored the impact of combining space-based and ground-based observations. Together, the HST and Spitzer observations confirm previous findings that the transmission spectrum is consistent with a cloud-free atmosphere. They are able to put a constraint on the absolute sodium and oxygen abundance and find them to be 21 +27 −14 × and 7 +11 −4 × solar values, respectively. Later, McGruder et al. (2022) published a study of WASP-96 b adding to the existing ensemble of transit measurements using the ground-based telescope IMACS/Magellan as part of the ACCESS project (Jordán et al. 2013;Rackham et al. 2017). Their transmission spectrum covers the spectral range of 0.44-0.9 µm, which overlaps with the VLT/FORS2 observations of Nikolov et al. (2018), enabling an independent confirmation of the sodium feature with its pressure-broadened wings. They combine their two transits with the published VLT/FORS2 and HST data and perform spectral retrievals on the combined spectrum (0.4-1.644 µm). Their results indicate solar-to-super solar abundances of Na and H 2 O, with log-mixing ra-tios of −5.4 +2.0 −1.9 and −4.5 +2.0 −2.0 respectively, in rough agreement with Nikolov et al. (2022) and with previous suggestions of super-solar alkali abundances .
This study is the second paper in a two-part series providing an in-depth treatment of the ERO observations of WASP-96 b. The companion paper, Radica et al. (submitted) focuses on the reduction and extraction of the planet's transmission spectrum from the SOSS time series observations (TSO), as well as provides some initial insights into the composition of WASP-96 b's atmosphere through comparisons of the planet's transmission spectrum with grids of self-consistent atmospheric models. They conclude that the NIRISS/SOSS observations of WASP-96 b are best explained by a cloud-free atmosphere, with a solar-to-super-solar metallicity atmosphere and solar carbon-to-oxygen ratio (C/O). In this study, we perform a detailed atmospheric characterization of WASP-96 b using the spectrum presented in Radica et al. (submitted). NIRISS/SOSS provides spectral coverage from 0.6-2.8 µm, covering the red wing of the 0.59 µm sodium doublet as well as multiple water bands. Hence we can assess the robustness of previous observations (Nikolov et al. 2022;McGruder et al. 2022) and independently confirm them with an instrument specifically built to study exoplanet atmospheres.
In the era of HST, exoplanet atmosphere observations were generally binned to a resolution which achieved a sufficient signal to noise ratio to obtain spectral information. The choice of binning, though, has remained somewhat arbitrary; for example, with HST/WFC3 G141 (1.1-1.7 microns), Line et al. (2016) use 10 spectral bins for the secondary eclipse of HD 209458b and Kreidberg et al. (2014) use 22 spectral bins for the transmission spectrum and 15 spectral bins for the secondary eclipse of WASP-43b. Hence, for the same instrument, different spectral binning was used throughout the literature, and as of yet, there has been no exploration of whether different spectral bin choices impact the inferred atmospheric properties.
We therefore aim to answer the following questions: (i) What are the chemical species present in WASP-96 b's atmosphere are what are their abundances? Is the atmosphere indeed cloud free?
(ii) How robust against the choice of framework and model assumptions are the retrieved chemical abundances?
(iii) Does binning the data from native to lower resolution produce different inferred abundances?
This work is organized as follows: we provide a brief overview of the observations and data reduction in Section 2. We outline the different modelling strategies, both through inference (e.g., retrievals) and forward models, in Section 3, and present the modelling results in Section 4. Section 5 contains a brief discussion of these results, and we summarize our work in Section 6.

OBSERVATIONS
One transit of the hot-Saturn WASP-96 b was observed with NIRISS/SOSS on June 21, 2022 as part of the JWST Early Release Observations (ERO) program (Pontoppidan et al. 2022). The total duration of the time series observation (TSO) was 6.4 hr. The SUBSTRIP256 subarray configuration was used to capture the first three diffraction orders of the target star on the detector (Albert et al. submitted), which provides access to the full 0.6-2.8 µm wavelength range of the SOSS mode. Order 1 covers the 0.85-2.8 µm wavelength range, with the 0.6-0.85 µm information provided by order 2. The third order is generally too faint to be extracted, and does not provide any unique wavelength coverage (Albert et al. submitted).
The reduction of these data, using the supreme-SPOON pipeline 1 (Feinstein et al. 2023;Coulombe et al. 2023, Radica et al. submitted) are treated in depth in Radica et al. (submitted). In that work, they present a walkthrough of the critical reduction steps, including correction of the zodiacal background light and 1/ noise. Their final transmission spectrum, which we make use of in this work, was extracted with the ATOCA algorithm (Darveau-Bernier et al. 2022;Radica et al. 2022) to explicitly model the self-contamination of the first and second diffraction orders on the detector. The transit depths were additionally fit at the pixel level (that is, one transit depth per pixel column on the detector), and were post-processed to correct for contamination from background field stars, which can occur due to the slitless nature of the SOSS mode. This spectrum is shown at the pixel-level, as well as binned several lower resolutions in Figure 1. Their spectrophotometric light curves reach an average precision of 1.2× and 1.4× the photon noise for order 1 and 2 respectively, resulting in an average pixel-level transit depth precision of 522 ppm and 534 ppm respectively.

SPECTRAL ANALYSIS
We use two different modelling approaches to thoroughly explore WASP-96 b's atmosphere. The first is comparisons with forward models computed using 3D General Circulation Models (GCMs). These allow us to explore the potential formation of different species of cloud condensates in the atmosphere of WASP-96 b (e.g., Samra et al. 2023), and also consider chemical kinetics (Zamyatina et al. 2023). The second is a spectral retrieval analysis that allows us to infer the atmospheric properties of WASP-96 b, such as its chemical composition and temperature structure, directly from the Radica et al. (submitted) (Cubillos & Blecic 2021) to ensure that the atmospheric composition we retrieve is robust against the choice of retrieval code. The model setups for each of the two approaches are detailed below.

3D Forward Modelling using GCMs
Although all previous observational studies of this planet have concluded a cloud-free upper atmosphere for WASP-96 b (Nikolov et al. 2018(Nikolov et al. , 2022McGruder et al. 2022), the idea that clouds could be present in the atmosphere of this planet was recently explored by Samra et al. (2023). Their 3D GCM model considers a kinetic, non-equilibrium formation model for mixed material cloud particles. Their GCM models show that clouds could indeed be ubiquitous in the low-pressure, terminator regions of WASP-96 b's atmosphere, with silicate and metal oxide clouds being the most prominent condensate species. They conclude that the Nikolov et al. (2022) transmission spectrum can also be fit with cloudy models. However, whether the clouds predicted by the kinetic model actually form in WASP-96 b depends on whether they are cold-trapped below the photosphere Powell et al. 2018), a mechanism that cannot currently be resolved with the kinetics models.
We perform our own GCM modelling to investigate the plausibility that WASP-96 b could host clouds. For our first GCM analysis, we use the non-grey SPARC/MITgcm (Showman et al. 2009). Specifically, we make use of the large grid of models generated by Roth et al. (in prep). The GCM setup is very similar to that described in Parmentier et al. (2018Parmentier et al. ( , 2021, but does not consider cloud condensation. The model grid spans a wide range of equilibrium temperatures, atmospheric metallicity, orbital period, and surface gravities which is then interpolated into the specific WASP-96 b parameters. However, other parameters such as the planet radii are fixed and the models have an infinite drag timescale. The resulting thermal profiles are then interpolated to the system parameters of WASP-96 b. The thermal profiles are read into CHIMERA (see Section 3.2.1 for more details) to produce a transmission spectrum of an atmosphere that has a solar C/O and metallicity of 1× and 10× solar. The thermal profiles and transmission spectra are shown in Figure 2. To capture the PT profile parameter space spanned by our range of considered metallicities, we denote the 1× solar profile in a solid line, and the edge of the shaded area denotes the 10× solar profile. It can be seen that the PT structures cross the condensation curves for various cloud species. Specifically, the SPARC/MITgcm predicts that the morning limb favours high-altitude Na 2 S clouds with deeper MnS clouds, whereas in the evening limb only MnS clouds could condense in the pressure regions probed by transmission. Silicate clouds should form only in the 1× solar metallicity case, and only in the deep layers of the atmosphere (∼1 bar). Depending on whether vertical mixing is large enough, they could be efficiently mixed up to the pressure levels probed by the observations or remain trapped in the deep layers of the atmosphere (Powell et al. 2018).
Our second GCM analysis utilizes the Met Office Unified Model (UM) run specifically for WASP-96 b. We used the same basic model setup as in Drummond et al. (2020) and Zamyatina et al. (2023) with the following changes: (a) the PHOENIX BT-Settl stellar spectrum (Rajpurohit et al. 2013) closely matching WASP-96, (b) WASP-96 b parameters from Hellier et al. (2014), and (c) the UM version 11.6, initialised with (d) WASP-96 b dayside-average pressure-temperature profile obtained with the 1D radiative-convective-chemistry model ATMO (Drummond et al. 2016) assuming chemical equilibrium for the chemical species present in the Venot et al. (2012) chemical network. We further assume the atmosphere to be cloud/haze free and have a solar metallicity and C/O ratio based on the initial modelling of Radica et al. (submitted).
Within the UM framework, we ran two simulations, each with a different chemical scheme: one assuming chemical equilibrium, and the other a chemical kinetics scheme, which computes the production and loss of the chemical species present in the Venot et al. (2019) reduced chemical network. We will refer to these simulations as "UM 1× solar equilibrium" and "UM 1× solar kinetics", respectively, with the latter simulation accounting for the opacity changes not only due to changes in pressure and temperature, but also due to the transport of chemical species in the atmosphere. The left panel of Figure 2 shows that both UM simulations predict similar limb-average PT profiles (weighted over all latitudes and ±20 • longitude), with the morning limb being colder than the evening limb at pressures <1 bar. Contrary to the SPARC/MITgcm, pressure-temperature profiles from the UM suggest that only MnS clouds could form on WASP-96 b's limbs. This is because the UM predicts a shallower temperature gradient at pressures <10 −2 bar causing the UM to have temperatures 100-200 K higher than those predicted at comparable pressures by the SPARC/MITgcm. Both GCMs predict similar positions for the MnS cloud decks on both limbs, when the assumed metallicity is 1× solar. However, given that the MnS nucleation rate is relatively low (Gao et al. 2020), these clouds might not form quickly enough for their opacity to be relevant for WASP-96 b.   The right panel of Figure 2 shows that both the UM and the SPARC/MITgcm simulations produce transmission spectra that agree well with WASP-96 b's JWST NIRISS/SOSS transmission spectrum in the range 1.3-2.15 µm. Blueward of 1.3 µm, however, K and H 2 O features are muted relative to those predicted by the haze and cloud free GCM models, suggesting the presence of a scattering opacity source. Redward of 2.15 µm, the observations and the models broadly agree, but the observed transit depths vary highly with wavelength. Of particular note is the region between 2.15-2.5 µm, where the "UM 1× solar kinetics" simulation predicts a higher transit depth than the "UM 1× solar equilibrium" simulation. This difference is caused by an enhancement of the abundance of CH 4 due to transportinduced quenching, which is captured only in the UM kinetics simulation. However, we are not able to robustly distinguish between these two cases with the current data. Another difference is that, because the solar-composition SPARC/MITgcm predicts cooler limb temperatures than the UM simulation, the spectral features of the SPARC/MITgcm are shallower than the observations. However, the SPARC/MITgcm 10× solar metallicity leads to a hotter thermal profile and thus to a better match to the data. This difference highlights the intrinsic dependence of the observables to the modelling framework when using complex, 3D, GCMs (Showman et al. 2020).
Overall, both clear-sky GCMs used in this study provide good agreement with our JWST NIRISS/SOSS transmission spectrum of WASP-96 b. However, we are not able to robustly distinguish between 1× and 10× solar metallicity models with the current data, and both models struggle to reproduce the observations blueward of 1.3 µm. This further motivates an in-depth investigation using atmospheric retrievals.

Atmospheric Retrieval
Atmospheric retrievals are a powerful tool to extract information about an exoplanet atmosphere directly from the data (Madhusudhan & Seager 2009). We explore the data in a hierarchical way, from simple (e.g., cloud free, free abundances, isothermal) to complex models (e.g., inclusion of hazes and clouds, chemical equilibrium, non-isothermal), with multiple retrieval codes. The first set of retrievals we perform are 'free chemistry' retrievals, which directly infer the volume mixing ratios (VMR) for a set of chemical species assumed to be present in the atmosphere (the VMRs are assumed constant with altitude). Each retrieval framework assumed that the atmosphere is dominated by H 2 -expected for objects that have physical properties similar to Saturn, and included the same molecules as opacity sources. All frameworks use the WASP-96 system parameters reported in Hellier et al. (2014). The second set of retrievals are performed assuming that the vertical abundances of the chemical species are in thermochemical equilibrium.
To robustly interpret our WASP-96 b observations, we employ four different retrieval frameworks: CHIMERA (Line et al. 2013), Aurora (Welbanks & Madhusudhan 2021), POSEIDON (MacDonald & Madhusudhan 2017;MacDonald 2023) and PyratBay (Cubillos & Blecic 2021). A multiple-retrieval approach allows us to compare our results in the regime of high-precision data (Barstow et al. 2020(Barstow et al. , 2022, thereby quantifying the stability of our atmospheric inferences to model implementations. The common molecules to each code are: H 2 O, CO, CO 2 , CH 4 , NH 3 , HCN, Na and K, these have a prior U(−12,−1) 2 for all VMRs. The set-up of each code is explained in the following subsections. Furthermore, CHIMERA is also used to run a chemical equilibrium retrieval, as an additional test.

CHIMERA
We use CHIMERA 3 to perform both free and chemically consistent spectral retrievals. CHIMERA is the only framework in this study that uses the correlated-approach (Lacis & Oinas 1991) when computing transmission through the atmosphere. The -tables are computed at a resolution of R=3000; the line-by-line data used to calculate the -tables are from the following sources: H 2 O (Polyansky et al. 2018;Freedman et al. 2014), CO 2 (Freedman et al. 2014), CO (Rothman et al. 2010a), CH 4 (Rothman et al. 2010a), HCN (Barber et al. 2014), Na (Kramida et al. 2018;Allard et al. 2019), and K (Kramida et al. 2018;Allard et al. 2016a), and were computed following the methods described in Gharib-Nezhad et al. (2021); Grimm et al. (2021). We assume the atmosphere is dominated by H 2 , with a He/H 2 ratio of 0.1764; therefore, we also model the H 2 -H 2 and H 2 -He collisioninduced absorption (CIA) (Richard et al. 2012).
To compute the thermal structure, we use the parameterisation described in Madhusudhan & Seager (2009). This approach splits the atmosphere into three layers: the upper atmosphere, where no inversion can take place, a middle region, where an inversion is possible, and a deep layer, where the thermal structure is isothermal. We also consider a scenario in which the temperature structure is isothermal, and find that the abundances do not dependent on the thermal structure parameterisations. We also consider an atmosphere that is just parameterised by an isothermal model, it can be seen that all retrievals tend towards an isothermal temperature structure 4.
Our chemically consistent retrievals aim to explore the impact of physical coupling between the atmospheric composition and temperature structure. Specifically, the molecular and atomic vertical abundances are assumed to be in thermochemical equilibrium. The equilibrium abundances are computed using the NASA CEA (Chemical Equilibrium with Applications) model (Gordon & McBride 1994) for a given C/O, metallicity, and temperature structure. Thus, the C/O ratio and metallicity are free parameters for these retrievals instead of the chemical abundances themselves.
We model hazes following the prescription of Lecavelier Des Etangs et al. (2008), which treats hazes as enhanced H 2 Rayleigh scattering with a free power law slope. This parameterisation expresses the opacity as Hazes = 0 ( / 0 ) , where is the Rayleigh enhancement factor and is the scattering slope (equal to −4 for H 2 Rayleigh scattering). 0 is the H 2 Rayleigh cross section at 0 , given by 2.3×10 −27 cm 2 and 430 nm respectively. Alongside the haze calculation, we fit for a constant-in-wavelength grey cloud with opacity cloud . Hence, we term this model "Simple Haze + Cloud model". To explore the parameter space, we coupled our parametric forward model with the Bayesian Nested Sampling algorithm PyMultiNest (Feroz et al. 2009;Buchner et al. 2014).

Aurora
We complement our atmospheric analysis by inferring the atmospheric properties of WASP-96 b using Aurora Welbanks & Madhusudhan (2021), a Bayesian atmospheric retrieval framework for the interpretation of ground-and space-based observations of transiting exoplanets. Our atmospheric model setup generally follows a similar approach to previous atmospheric studies (e.g., Welbanks & Madhusudhan 2019) with the same priors for WASP-96 b as in the analysis of the existing VLT observations (Nikolov et al. 2018) presented in . Our atmospheric model computes line-by-line radiative transfer in transmission geometry in a plane-parallel atmosphere. The pressure structure of the atmosphere assumes hydrostatic equilibrium for a varying-with-height gravity, in a grid of 100 layers uniformly distributed in log-pressure from 10 −7 to 100 bar. The Bayesian inference is performed using the framework MultiNest (Feroz et al. 2009) through its Python implementation Py-MultiNest (Buchner et al. 2014) using 2000 live points.
We explore a series of atmospheric model scenarios with Aurora including the possibility of multidimensional clouds and hazes (e.g., , terminator inhomogeneities (e.g., Welbanks & Madhusudhan 2022), and other modelling assumptions regarding the number of free parameters in our retrievals (e.g., . Through this exploration of models we determined a fiducial 19 parameter model for our 'free retrieval' analysis using Aurora and other similar frameworks. This model setup considers a non-isothermal pressure-temperature structure parameterized using the six parameter prescription of Madhusudhan & Seager (2009). Eight sources of opacity are considered in our models. These species, expected to be the main absorbers for hot gas giants (e.g., Madhusudhan 2019), are parameterized by their logarithmic volume mixing ratios assumed to be constant with height. We consider the presence of clouds and hazes in our atmospheric models using the modeling strategy for inhomogeneous terminator cover presented in Line & Parmentier (2016). We consider the presence of scattering hazes as deviations from the Rayleigh scattering in the models by following the parameterization of Lecavelier Des Etangs et al. (2008) as described above. The spectroscopic effect of clouds is included by considering the presence optically thick cloud decks at a specific pressure level. The combination of inhomogeneous clouds and hazes is implemented following the single-sector prescription as explained in  using four additional free parameters. Finally, we use one free parameter to infer the reference pressure corresponding to the assumed planetary radius. To compare our high-resolution (R∼30,000) spectra to the NIRISS/SOSS observations we follow the model binning strategy presented in Pinhas et al. (2018).

POSEIDON
The third atmospheric retrieval code we employ is POSEIDON (Mac-Donald & Madhusudhan 2017; MacDonald 2023). POSEIDON is a well-established atmospheric modelling and spectral retrieval code that was recently released as an open-source 4 Python package (Mac-Donald 2023). The radiative transfer technique underlying POSEI-DON's transmission spectrum forward model is described in Mac-Donald & Lewis (2022). Our POSEIDON retrieval samples the parameter space using the Bayesian nested sampling algorithm Multi-Nest, deployed via its Python wrapper PyMultiNest (Feroz et al. 2009;Buchner et al. 2014).
Our WASP-96 b POSEIDON retrieval analysis employs a 19-parameter model accounting for non-isothermal pressuretemperature profiles, inhomogeneous clouds and hazes, and the eight common chemical species described above. One parameter encodes the planetary radius at a 10 mbar reference radius. The five-parameter PT profile follows the prescription in Madhusudhan & Seager (2009), 4 POSEIDON is available here: https://github.com/ MartianColonist/POSEIDON modified to place the reference temperature parameter at 10 mbar. The four-parameter inhomogeneous aerosol model follows MacDonald & Madhusudhan (2017). Finally, eight parameters specify the constant-in-altitude free abundances of H 2 O, CO, CO 2 , CH 4 , HCN, NH 3 , Na, and K. The model constructs an atmosphere ranging from 10 −8 -100 bar, with 100 layers uniformly distributed in log-pressure, and assumes a H 2 + He-dominated background atmosphere with He/H 2 = 0.17. The Bayesian retrieval of this 19-parameter space used 1,000 PyMultiNest live points.
At each location in the parameter space, POSEIDON computed WASP-96 b transmission spectra at a resolution of = 20,000 from 0.55-2.9 µm. The radiative transfer uses opacity sampling of highresolution pre-computed cross sections ( ∼ 10 6 ) from the following  (Ryabchikova et al. 2015), and K (Ryabchikova et al. 2015). We additionally include continuum opacity from H 2 and He CIA (Karman et al. 2019) and H 2 Rayleigh scattering (Hohm 1994). We convolve each = 20,000 model spectrum with the instrument Point Spread Function (PSF), before binning down to the resolution of the observations (here, = 125) to compute the likelihood of each parameter combination. We treat NIRISS/SOSS orders 1 and 2 separately during the convolution and binning procedure, accounting for their different intrinsic PSFs and instrument transmission functions.

PyratBay
Lastly, we also employed PyratBay, the Python Radiative-transfer in a Bayesian framework. PyratBay 5 is an open-source framework for exoplanet atmospheric modelling, spectral synthesis, and Bayesian retrieval. It utilizes the most up to date line-by-line opacity sources from ExoMol (Tennyson et al. 2016), HITEMP (Rothman et al. 2010a), atomic species Na, K (Burrows et al. 2000) and collision-induced opacities of H 2 -H 2 (Borysow et al. 2001;Borysow 2002) and H 2 -He pairs (Borysow et al. 1988. For effective use in retrieval, we compress these large databases (while retaining information from the dominating line transitions), using the available package (Cubillos 2017). To model the vertical temperature structure, we implement three parameterization schemes: isothermal, and Line et al. (2013) and Madhusudhan & Seager (2009) prescriptions. This retrieval framework also implements a self-consistent 1D radiative-convective equilibrium scheme (Malik et al. 2017), the classic "power law+gray" prescription, a "single-particle-size" haze profile, a "patchy cloud" prescription for transmission geometry (Line & Parmentier 2016), and two complex Mie-scattering cloud models. The first is a fully self-consistent microphysical kinetic cloud model of Helling & Woitke (2006) which follows the formation of seed particles, growth of various solid materials, evaporation, gravitational settling, elemental depletion and replenishment (Blecic et al. 2023). The other is a parametrized Miescattering thermal stability cloud model (Kilpatrick et al. 2018;Venot et al. 2020).
In this work, we assumed the atmosphere of WASP-96 b to be hydrogen dominated (He/H 2 = 0.17) and include collision-induced absorption of H 2 -H 2 and H 2 -He. We included molecular opacity sources of H 2 O (Polyansky et al. 2018), CH 4 (Hargreaves et al. 2020), NH 3 (Yurchenko et al. 2011;Yurchenko 2015), HCN (Harris et al. 2006(Harris et al. , 2008, CO (Li et al. 2015), and CO 2 (Rothman et al. 2010b) and resonant-line cross-sections of Na and K. In addition, we account for the Rayleigh-scattering cross-section of H 2 (Dalgarno & Williams 1962) and an unknown haze particulate, by applying a power law prescription of Lecavelier Des Etangs et al. (2008). Our radiative transfer routine uses opacity sampling of high-resolution precomputed cross sections tables generated at a resolution of R∼4×10 7 , calculates the transmission spectra at = 20,000, and computes the likelihood of each model by binning it down to resolution of = 125. We generated the atmosphere between 10 −9 -100 bar, with 81 layers uniformly distributed in log-pressure, retrieving in addition to the constant-with-altitude molecular and alkali volume mixing ratios listed above, also the Lecavelier Des Etangs et al. (2008) haze parameters and the planetary radius at the reference pressure of 0.1 bar. To find the best modelling setup we tested our available temperature parametrizations and the full range of the cloud models from simple to complex Mie-scattering clouds, assuming species expected to be seen on this temperature regimes. We compared these models using the Bayesian Information criteria (BIC Liddle 2007). We found the lowest BIC for the model assuming Madhusudhan & Seager (2009) temperature prescription with patchy opaque cloud deck and hazes, accounting for both Lecavelier Des Etangs et al. (2008) and Dalgarno & Williams (1962) haze particles and opacities from H 2 O, CO 2 , CO, Na and K. To explore the phase space of these parameters, we have coupled our atmospheric model with the Bayesian Nested Sampling algorithm PyMultiNest (Feroz et al. 2009;Buchner et al. 2014) and the Multi-core Markov-chain Monte Carlo code MC3 (Cubillos et al. 2016). Both algorithms returned the same constraints.

RESULTS
In this section we present the results from our retrieval analysis. We also discuss the impact the resolution of the data has on our inferred abundances.

Retrievals
Using the frameworks described above, we infer the atmospheric properties of WASP-96 b using the NIRISS/SOSS observations binned to four different constant resolutions (R = 125, 250, 500, and pixel level). As discussed below (see Section 4.2), we find our inferences robust regardless of the resolution of the binned observations. Therefore, we present our results using the R = 125 binned observations for clarity. Our first consideration is the possible presence of clouds and hazes. As described above, our atmospheric frameworks compute scenarios representative of cloud-free atmospheres, hazy atmospheres, cloudy atmospheres, and atmospheres with inhomogeneous cloud and haze cover. Comparing these atmospheric scenarios using their Bayesian evidence and comparing them to a 'sigma' scale (e.g., Benneke & Seager 2013; Welbanks & Madhusudhan 2021), we find a 6 model preference for inhomogeneous clouds and hazes over simple cloud-free atmospheres. However, we note that it is primarily a Rayleigh scattering slope which we detect as opposed to any opacity from a grey cloud deck (see Section 4.1.2). We thus limit our discussion to the "inhomogeneous clouds and hazes" model runs moving forwards.
The use of more complex prescriptions separating the spectroscopic effects of clouds from those of hazes across inhomogeneous terminators (e.g., Welbanks & Madhusudhan 2021), may result in lower model preferences but consistent inferred atmospheric properties. The full retrieved posterior distributions for Aurora, POSEI-DON, PyratBay and CHIMERA can be found in Figures A2, A3, A4 and A5 respectively.

Retrieved Abundances
The results for the inhomogeneous haze and cloud model runs for all four frameworks are presented in Figure 4, where we present the best-fit transmission spectra, thermal structure and posteriors for H 2 O, CO, CO 2 , Na and K. We do not present the posteriors for NH 3 , HCN or CH 4 as they remain mostly unconstrained given existing observations. The retrieved abundances from all codes are summarized in Table 1, and generally remain consistent within 1-, demonstrating that the retrieved atmospheric properties are robust against different model implementations. They are also largely consistent with a solar metallicity atmosphere, in agreement with the interpretation of Radica et al. (submitted) using self-consistent radiative thermochemical equilibrium models.
Using the Aurora framework, we then assess the detection significance of each molecule. This is done by computing the Bayesian evidence for a model without each molecule and comparing to the original model with all species included. We present the breakdown in Table 2.
As a final test, we perform a chemically consistent retrieval on the same data using CHIMERA in order to directly retrieve the atmosphere log(C/O) and log(Met). Like the free retrieval, we fit for our Simple Haze + Cloud model. We  Figure A6. Therefore we find that the data is consistent with a model that has a solar C/O ratio within 1 and a solar metallicity within 1 . These results are consistent with the modelling work presented in Radica et al. (submitted). We further demonstrate the consistency with Radica et al., (submitted) in Figure 3, the left panel compares the free retrieved results compared to the volume mixing ratios obtained from the best fitting model in Radica et al., (submitted), it can be seen that the abundances obtained in our free retrieval are consistent with these profiles. The outlier is the abundance of CO 2 , which we find to be consistent with a volume mixing ratio of 10x solar. The right panel shows the retrieved volume mixing ratios for the chemical equilibrium framework, these are again consistent with the free retrieval and the models of Radica et al., (submitted).

Retrieved Cloud Parameters
We describe in more detail the model preference for inhomogeneous clouds and hazes over the cloud free model described above. The models considering the presence of inhomogeneous clouds and hazes suggest a large fraction (≳70% i.e., = 0.88 +0.  PyratBay) suggesting that the spectroscopic impact of these gray clouds is minimal. Similarly, the low cloud opacity (e.g., log( cloud ) = -32.66 +1.62 −1.48 ) retrieved by our CHIMERA analysis suggests low impact due to clouds.
On the other hand, our inferred haze scattering properties suggest they make a significant contribution in our WASP-96 b observations. While the scattering slope is retrieved to be largely Rayleigh-like The best-fit retrieved vertical volume mixing ratios obtained from the chemical equilibrium retrieval, the shading representing the 1 uncertainty. These are compared to retrieved results from the free retrieval from Aurora. We note that the vertical location of the retrieved free abundances are arbitrary and do not represent the region probed.  Table 1. Retrieved abundances and their accompanying error for our baseline model. We present the results from all of the retrieval frameworks. In the final row we present the abundance of each species calculated at solar metalicity and C/O. The elemental abundances were obtained from Lodders & Fegley (2002). We present the vertical volume mixing ratio for these species and compare them to the retrieved values in Figure 3.

2.02
Clouds and Hazes -6.69 Table 2. Retrieved abundances from Aurora and their accompanying detection significance. We also present the detection significance of using our cloud+haze model compared to a cloud free model. one order of magnitude (log 10 ( ) = 1.85 +0.73 −0.47 Aurora; 1.70 +0.60 −0.41 POSEIDON; 2.49 +0.95 −0.77 PyratBay). The inferences from the chemical equilibrium retrievals with CHIMERA remain largely in agreement and suggestive of spectroscopic signatures of Rayleigh scattering rather than clouds (e.g., log 10 ( cloud ) = −33.21 +1.20 −1.11 , = −3.31 +0.43 −0.50 , log 10 ( ) = −1.39 +0.26 −0.27 and f = 0.93 +0.05 −0.10 ). All models thus tell the story of an atmosphere with small aerosol particles that produce a Rayleigh scattering slope at short wavelengths, but no evidence for a grey cloud deck, which, as is also the case with our chemical inferences above, is consistent with the interpretation of Radica et al. (submitted).

Resolution Testing
Atmospheric retrievals can be computationally demanding, and the spectral resolution of the forward model is a large factor in determining the speed of the calculation. To thoroughly study the spectrum of an exoplanet atmosphere, one needs to perform multiple retrieval studies, each study requiring on the order of 10 4 to 10 5 model calculations; which can become unfeasible at the native R∼700 resolution of NIRISS/SOSS. In this section, we seek to answer the question: Do we infer the same abundances if we bin the native resolution data to lower resolutions?
To answer this we perform a retrieval analysis on three different transmission spectrum resolutions: R = 125, R = 250, and R = 500, shown in Figure 1. We use the same parameterised model presented in Figure 4 and correlated-tables calculated at R=3000; hence the model has a resolution six times greater than the maximum data resolution. We find that the retrieved abundances for data with a resolution of R=125 are the same as with a resolution of R=500. Hence, no information is lost when binning the data. We present the posteriors of H 2 O, CO 2 , and K in Figure 5. The colours correspond to those in Figure 1. We binned our observations to R=500, R=250, and R=125 to explore if this binning down of the data causes a loss of information. Posteriors for R=125 are shown in red, blue is R=250 and purple is R=500. The points and error bars show the median retrieved value, and 1 credible interval for each test. A retrieval on each resolution yields consistent abundances to well within 1 , allowing us to conclude that no information is lost when binning our WASP-96 b NIRISS/SOSS transmission observations. held the unique privilege of being one of the few "cloud-free" exoplanets known. Subsequent studies (Nikolov et al. 2022;McGruder et al. 2022) added HST/WFC3 transit depths, as well as additional ground-based transmission observations from Magellan/IMACS, however the conclusion of the cloud-free nature of WASP-96 b's upper atmosphere remained unchanged. The GCM models of Samra et al. (2023), though, found that the terminator region of WASP-96 b should be entirely covered in clouds given the temperature structure of the planet. Moreover, they show that cloudy transmission spectra can provide an equally good fit to the ensemble of transmission data analyzed in Nikolov et al. (2022). Our two independent GCM models also predict that clouds should be able to form at the terminator of WASP-96 b in the pressure regions probed by transmission spectroscopy (see Figure 2). These models predict that the atmosphere is likely dominated by MnS and Na 2 S clouds. MgSiO 3 clouds should form in the deep layers of the atmosphere and would be observable only if the vertical mixing was extremely large to easily replenish the upper atmosphere in cloudforming material, an assumption that is inherent to the Samra et al. (2023) calculation.

Since
One solution to this discrepancy could be that smaller particles than predicted by Samra et al. (2023) form in larger quantities at low pressures in WASP-96b's atmosphere. These could be composed of 2 or KCl which would naturally form at much lower pressures than the silicate clouds that dominate the cloud composition in the 100 to 10 mbar range. However, the detection of sodium and potassium in WASP-96b's atmosphere seems to rule out this possibility. MnS is another candidate for forming clouds at low pressures Morley et al. 2012), however, Gao et al. (2020) predicts that the nucleation rates for MnS are so low that they should hardly form. Another option would be the formation of a high altitude haze layer formed of photochemically produced particles. Photochemistry is known to naturally form small particles at low pressures that can produce strong scattering slopes (Lavvas & Koskinen 2017;Kawashima & Ikoma 2019;Helling et al. 2020;Steinrueck et al. 2021). Additional information about the cloud composition could be gathered by targeting the resonant features of the cloud-forming material in the JWST/MIRI LRS bandpass.
We further note that our detection of a strong scattering slope in the optical is partially degenerate with the abundance of gaseous sodium in the atmosphere. Indeed when a scattering slope is not included in the retrievals, we obtain an unphysical alkali abundances (e.g., log 10 (Na) = −2.54 0.28 −0.34 with CHIMERA). However, including enhanced Rayleigh scattering, the Na abundances drops to slightly super-solar to solar values, in agreement with Nikolov et al. (2022). Our inferred abundances of Na and of the presence of a scattering slope therefore needs to be carefully interpreted because of this degeneracy, driven by the fact that the NIRISS/SOSS bandpass cuts off at 0.6 µm, and is therefore only able to probe the red wing of the Na feature. Without fully resolving the Na feature peak, it is difficult to differentiate between a slope caused by a Rayleigh scattering haze, or the red wing of a broadened Na feature. More work needs to be conducted to further understand this degeneracy in the context of observations with NIRISS/SOSS.

Comparison to Radica et al., (submitted)
A suite of forward models were compared to the data in our companion paper (Radica et al., submitted). Three different grids of models were used: PICASO, ATMO and ScCHIMERA, producing a picture of an atmosphere that has a metallicity of 1 -5× solar and a solar C/O. Our free retrieval results demonstrate that we are obtaining an abundance of H 2 O that is consistent with solar values and a CO 2 abundance that is super solar, this demonstrates that our results are consistent with Radica et al., submitted. Similarly to Radica et al., submitted, we need to invoke enhanced Rayleigh scattering slope to match the observations at the shortest wavelengths, but find no spectroscopic impact from a grey cloud deck. We compare the vertical volume mixing ratios obtained from the best fit ScCHIMERA model with our retrieved results in Figure 3 which shows we are obtaining a consistent picture of the atmosphere.

CONCLUSIONS
In this paper we have performed a detailed atmospheric characterisation of WASP-96 b using the transmission spectrum obtained with NIRISS/SOSS as part of the Early Release Observations, and first presented in Radica et al. (submitted).
We ran GCM simulations in order to model the planet's atmosphere using the SPARC MIT/gcm and the UM. These clear-sky models are able to well fit the spectrum redward of 1.3 µm, and favour an atmosphere with solar metallicity. However, blueward of 1.3 µm, the GCMs under predict the observed transit depths, likely indicating missing opacities such as a scattering haze.
We then performed a suite of retrievals using four different modelling frameworks: CHIMERA, Aurora, PyratBay, and POSEI-DON. We find that a model with patchy clouds and hazes best describe the data and that each framework produces results which are consistent within 1 . We report the retrieved abundances from Aurora as log 10 (H 2 O) = −3.59 +0.35 −0.35 , log 10 (K) = −8.04 +1.22 −1.71 , log 10 (CO) = −3.25, and log 10 (CO 2 ) = −4.38 +0.47 −0.57 . We find a large tail in the posterior CO, so we describe this abundance as an upper limit. Further transmission observations with JWST, particularly with NIRSpec G395H are necessary to more accurately constrain the abundance of CO.
The retrieved abundance of H 2 O is consistent with Yip et al. (2021) andMcGruder et al. (2022). Our precision is ∼10× better than McGruder et al. (2022) and ∼4× better Yip et al. (2021). Our range of retrieved abundances of Na is consistent with Nikolov et al. (2022), McGruder et al. (2022, Yip et al. (2021) and , however given that NIRISS' wavelength coverage does not capture the complete Na feature, this results in a degeneracy between the abundance of Na and a Rayleigh scattering slope. This is also reflected in the extremely low detection significance for Na (1.24 ). We therefore caution against any strong interpretations of this Na abundance. We also report a constrained abundance of potassium, although with only a marginal detection significance (∼2 ), in the atmosphere of WASP-96 b which was not found in previous studies due to the lower resolution of the optical data. The strong potassium constraint in the atmosphere of WASP-39 b from NIRISS/SOSS (Feinstein et al. 2023), and the tentative detection here, demonstrates how powerful this instrument is to study alkali metals, and opens the door for a new tracer of formation history, the K/O ratio (Feinstein et al. 2023).
Our chemically consistent retrievals favour an atmosphere that has a solar C/O ratio within 1 and solar metallicity within 1 . We find the log(C/O) = −0.30 +0.17 −0.37 and log(Met) = −0.63 +0.64 −0.44 , where solar values are log(C/O) = −0.26 and log(Met) = 0. This is consistent with the GCM models and the grid models in Radica et al. (submitted) which favour an atmosphere that is 1× solar C/O and 1 -5× solar M/H. We explore the appropriate resolution to study observations obtained with NIRISS/SOSS. We find that binning the data from native to R=125 does not impact the inferred abundances. This is useful, given that retrievals at native resolution are computational demanding. In the era of JWST, we need to explore more complex models, which are computational demanding in themselves, therefore we should trade data resolution for model complexity.
Finally, it is critical to note that the previous studies retrieved on a transmission spectrum created through the combination of multiple instruments, with six transits required to construct the Nikolov et al. (2022)

DATA AVAILABILITY
All data used in this study is publicly available from the Barbara A. Mikulski Archive for Space Telescopes 6 . The models generated in this paper can be made available on request.