High-precision Atmospheric Characterization of a Y Dwarf with JWST NIRSpec G395H Spectroscopy: Isotopologue, C/O Ratio, Metallicity, and the Abundances of Six Molecular Species

The launch of the James Webb Space Telescope (JWST) marks a pivotal moment for precise atmospheric characterization of Y dwarfs, the coldest brown dwarf spectral type. In this study, we leverage moderate spectral resolution observations (R ∼ 2700) with the G395H grating of the Near-Infrared Spectrograph (NIRSpec) on board JWST to characterize the nearby (9.9 pc) Y dwarf WISEPA J182831.08+265037.8. With the NIRSpec G395H 2.88–5.12 μm spectrum, we measure the abundances of CO, CO2, CH4, H2S, NH3, and H2O, which are the major carbon-, nitrogen-, oxygen-, and sulfur-bearing species in the atmosphere. Based on the retrieved volume mixing ratios with the atmospheric retrieval framework CHIMERA, we report that the C/O ratio is 0.45 ± 0.01, close to the solar C/O value of 0.458, and the metallicity is +0.30 ± 0.02 dex. Comparison between the retrieval results and the forward modeling results suggests that the model bias for C/O and metallicity could be as high as 0.03 and 0.97 dex, respectively. We also report a lower limit of the 12CO/13CO ratio of >40, being consistent with the nominal solar value of 90. Our results highlight the potential for JWST to measure the C/O ratios down to percent-level precision and characterize isotopologues of cold planetary atmospheres similar to WISE 1828.


INTRODUCTION
Y dwarfs, the coldest brown dwarf spectral type with temperatures below 500K (Cushing et al. 2011;Kirk-patrick et al. 2011), offer a unique opportunity for studying the rich atmospheric chemistry and physics in cold giant planets.The combination of low temperature and molecular opacity of CH 4 , H 2 O, and NH 3 causes the Ydwarf spectrum to be dim in the near-infrared (1-2.5 µm) while bright in the mid-infrared (3-5 µm), coinciding with the opacity window of water vapor and methane (Marley & Leggett 2009).Mid-infrared observations are thus critical for mapping the spectral energy distribution and characterizing Y-dwarf atmospheric properties.The ease of observations compared to cold giant planets orbiting bright stars makes Y dwarfs ideal exoplanet analogs.By studying Y dwarf atmospheres, we can bridge the knowledge gap regarding how giant planets' atmospheric properties evolve across different temperatures, from the Jupiter-like temperature (∼150K) to the warm (> 500 K) temperature that is found among Tdwarfs and exoplanets like 51 Eri-b.Furthermore, comparative studies of Y-dwarf and other giant planet atmospheres will shed light on potential connections between atmospheric properties and the formation and evolution histories of cold giant planets.
To date, approximately 50 Y-dwarfs have been discovered, primarily through the Wide-field Infrared Survey Explorer (WISE) and Spitzer photometric observations (e.g., Kirkpatrick et al. 2019).Studies by Morley et al. (2012); Leggett et al. (2013Leggett et al. ( , 2017) ) reveal that the colormagnitude relation (e.g., H-W2 vs. M W 2 ) of Y-dwarfs deviate from that of the chemical equilibrium cloudless models.To better explain the Y-dwarf observed colors and absolute magnitudes, numerous atmospheric models have been developed to explore the possible chemical and physical properties, including water and sulfide clouds (Morley et al. 2012;Lacy & Burrows 2023), nonequilibrium chemistry (Phillips et al. 2020;Mukherjee et al. 2022;Lacy & Burrows 2023), non-adiabatic thermal structure (Leggett et al. 2021), non-solar metallicity (Marley et al. 2021), and a non-solar carbon-to-oxygen (C/O) ratio (Cushing et al. 2021).Despite the significant advance in atmospheric modeling of Y dwarfs, the existing discrepancy between observation and models underscores that our existing understanding of Y-dwarf atmospheres remains, as of yet, incomplete.
Spectroscopic data of Y dwarfs (e.g., Cushing et al. 2011;Tinney et al. 2012;Kirkpatrick et al. 2013;Leggett et al. 2014;Leggett & Tremblin 2023) are essential for characterizing their atmospheres and testing atmospheric models of cold giant planets.However, spectroscopic observations of Y dwarfs are challenging due to their faintness in the near-infrared, the high thermal background, and the effects of telluric absorption in ground-based mid-infrared observations (e.g., Miles et al. 2020).Space-based infrared spectroscopy is therefore vital for acquiring high-precision infrared spectra of Y dwarfs.One of the largest collections of Y-dwarf nearinfrared spectra was obtained using the HST/WFC3 instrument by Schneider et al. (2015).Based on the dataset, Zalesky et al. (2019) retrieved the temperaturepressure profiles, atmospheric composition, metallicity, and C/O ratios of 22 late-T and Y dwarfs.A recent JWST spectroscopic study of a Y dwarf by Beiler et al. (2023) identifies the ν 3 ammonia absorption features at 3 µm.These studies demonstrate the invaluable insights into Y-dwarf atmospheres gained from the spectroscopic data.
In this study, we utilize JWST mid-infrared 3-5 µm moderate-resolution (R ∼ 2700) spectroscopy to characterize the atmospheric properties of the Y dwarf WISEPA J182831.08+265037.8.With the unprecedented sensitivity and mid-infrared spectral resolution of JWST, we aim to answer the following questions: • What are the abundances of the detected gas species in the Y-dwarf atmosphere?
• How does the retrieved Y-dwarf atmospheric thermal structure compare to the radiative-convective equilibrium atmospheric models?
• Are there any spectral features of trace molecules and isotopologues present in the Y-dwarf spectrum?
Atmospheric modeling studies (e.g., Beichman et al. 2013;Leggett et al. 2013;Cushing et al. 2021) find that it is challenging to simultaneously fit models to the near-infrared YJHK and the mid-infrared (e.g., WISE W1 and W2) broadband photometry.Given the unusual colors and overluminosity, Beichman et al. (2013); Leggett et al. (2017) suggested WISE 1828 could be a tight binary system.Cushing et al. (2021) present the Hubble Space Telescope/Wide Field Camera 3 (WFC3) 1.1-1.7 µm spectrum of WISE 1828.They found that an atmospheric model with non-solar metallicity and C/O ratio provides a better fit to the near-infrared spectrum and mid-infrared photometry.Leggett et al. (2021) show that a binary system with a modified temperaturepressure atmospheric profile provides a decent fit to the observed photometry and spectrum.Recently, De Furio et al. (2023) reported no evidence of binarity with a separation beyond 0.5 au based on JWST NIRCam photometric observations.

OBSERVATION AND DATA REDUCTION
The JWST/NIRSpec observations of J1828 were executed on 2022 July 28, as part of the GTO Program PID 1189.At the beginning of the observations, a NIR-Spec Wide Aperture Target Acquisition (WATA) image with a 3.6s exposure time and a clear filter was obtained to place the target at the slit center for fixedslit spectroscopy.The NIRSpec F290LP/G395H grating was used to obtain a J1828 spectrum ranging from 2.8 to 5.2 µm under the NRSIRS2RAPID readout mode.Three-point dithering along the slit S200A1 was performed during the observation.At each dither pointing, there were two integrations which were averaged over 26 groups each.The total exposure time for the G395H observation amounted to 2363.4 seconds.For G395H spectra obtained with the S200A1 slit, there is a gap between detectors NRS1 and NRS2.The detector gap introduces a wavelength gap between 3.685 and 3.789 µm in the 2.880-5.142µm spectra.
The data reduction was performed by mostly following the standard STScI pipeline with modified background subtraction, cosmic ray removal, and spectral extraction steps, described below.The data were processed through the JWST pipeline version with a CRDS version of 11.16.21under the context of jwst 1089.pmap.The data reduction steps are mainly done in three stages.In Stage 1, the uncal.fitsfiles were processed through the superbias subtraction, reference pixel correction, non-linearity correction, dark current correction, cosmic rays detection, ramp fitting step, and gain scale correction that eventually returned count-rate images with pixel values in the units of electrons per second.In Stage 2, the count rate images were assigned to the world coordinate system before proceeding to the customized background subtraction step.For each exposure, a columnaveraged background was estimated by taking the median of the pixel values at each wavelength in the crossdispersion direction among the source-free regions.The variances of the median values over each column were treated as the uncertainty of the one-dimensional sky background and added to the Poisson noise per pixel.After the background subtraction, the data was processed through wavelength correction, flat-fielding, and pathloss corrections, and was converted from an electron count rate to Janskys per steradian.The dispersed spec-tra were then resampled and rectified so that the x-axis of the spectral images were equivalent to the spectral dispersion axis.To identify any cosmic rays or bad pixels missed in the pipeline, we measured the centroid of each image column over a cross-dispersion aperture size of six pixels.A three-sigma clipping of the measured centroid values was used to flag the column numbers of spectral images that were possibly affected by cosmic rays or bad pixels.Finally, an aperture size of six pixels was used in the spectral extraction step of the Stage 2 pipeline.The centers of the spectral traces for the three dithering positions were estimated as y= 9, 20, and 28 in the configuration file of spectral extraction step extract 1d.

Estimation of flat-field uncertainty
The wavelength-dependent flat-field uncertainty of NIRSpec was not yet available during our data reduction stage and is regularly being updated with more in-flight calibration data.Based on the in-flight observation of a white dwarf, the NIRSpec G395H spectrum flux has an 0.91 ± 1.97% RMS residual compared to the model flux template (Böker et al. 2023).The estimate is likely a lower limit because the reduction uses the flat-field reference file from the same dataset and does not include slitloss and systematic errors.By comparing the spectra at different dithering positions, we found that the portion of spectra with estimated signal-to-noise ratios above 30 differ from each other systematically by around the 6% level.The flat-field uncertainty at around 6% is higher than the data uncertainty that is as low as 2% based on the read noise and Poisson noise.We therefore include an additional data uncertainty parameter as a nuisance parameter in our spectral fitting in Section 3.1.

RESULTS
The median spectra averaged over three dithering positions are shown in Figure 1.Based on the estimated noise that includes read noise, photon noise, and the background subtraction uncertainty, the peak of the spectra in the 4-4.5 µm exceeds a signal-to-noise ratio (SNR) of 60 while the 3.00-3.10µm spectrum has SNRs of around 10.The spectrum shares overlapping wavelength coverage with the Spitzer [4.5] band.We integrate the G395H spectrum and derive a Spitzer [4.5]band magnitude of 14.272 ± 0.017, which is 0.048 mag, or 1.8σ higher than the measured Spitzer [4.5] mag of 14.32 ± 0.02 (Kirkpatrick et al. 2019).
We then utilize two complementary atmospheric modeling tools to fit the G395H spectrum and to characterize the atmospheric thermal structure, compositions, C/O ratio, and metallicity in the following subsections.

Description of atmospheric modeling tools
Self-consistent radiative-convective equilibrium (RCE) models and atmospheric retrieval methods provide complementary insights into the atmospheric physics and chemistry of Y-dwarf atmospheres.In RCE models, the thermal and chemical structures are constrained by known physical processes.On the other hand, atmospheric retrievals are free to explore the parameter space which provides the best fit to the data.Therefore, cross-checking the retrieval results against the best-fit RCE models is useful for examining the possible physical and chemical processes and insuring that the retrieval results are not non-physical.We leverage both the RCE models (Sonora Elf-Owl Grid, Murkherjee et al. in press) and an atmospheric retrieval framework (CHIMERA, Line et al. 2014a,b) to fit the spectrum and characterize the temperature, gravity, C/O ratio, metallicity, and atmospheric composition.

Sonora Elf Owl models
The Sonora Elf Owl model grid (Murkherjee et al. in press) is a RCE cloudless atmospheric model grid that includes a self-consistent treatment of non-equilibrium chemistry of NH 3 , CO, H 2 O, CH 4 , CO 2 , N 2 , PH 3 though 1D vertical mixing as described in Mukherjee et al. (2022).The model grid was computed with the PICASO atmospheric model (Batalha et al. 2019;Mukherjee et al. 2023) and spans five parameters including effective temperature (T ef f = [275,2400] K), gravity (log g cms −2 = [3.25,5.5]),atmospheric metallicity ([M/H] =[-1,+1]), carbon-to-oxygen ratio (C/O = [0.22,1.14]), and vertical eddy diffusion coefficient (log(K zz (cgs))= [2,9]).The Elf Owl model grid assumes a constant K zz throughout the atmosphere.The high resolution spectra from the atmospheric model grid was computed using the resampled opacity grid at a wavelength resolution of 60,000 from Batalha et al. (2020).For this work, the relevant line lists that were used to compute the opacities are CO (Li et al. 2015) There is one subtlety to the Elf Owl grid in how the PH 3 abundance is computed.The Elf Owl grid computed the PH 3 abundance assuming quenching of the species for all but PH 3 .For PH 3 using the quenching methodology increases the abundance beyond chemical equilibrium.Because there is no observational evidence of strong PH 3 absorption, thePH 3 abundance is kept at chemical equilibrium.We linearly interpolated the model spectra for parameter points in between the Elf Owl model grid points.We then performed instrumental broadening, rotational broadening v sin i, and wavelength shifting, or equivalently radial velocity shifting, on the model spectrum.The instrumental broadening kernel is a wavelength-dependent function that is the same as the NIRSpec G395H spectral dispersion profile (1 ) and was scaled to match the unresolved emission line width of the calibration dataset (PID 1125).The rotational broadening was done using PyAstronomy (Czesla et al. 2019).Finally, we scaled the model spectrum with a scaling factor of r 2 /d 2 , where r is the radius of the brown dwarf and d is the distance of 9.92 pc (Kirkpatrick et al. 2019).We use the "MLFriends" Nested Sampling Algorithm (Buchner 2016(Buchner , 2019) ) implemented in the open source code UltraNest (Buchner 2021) and find best-fit solution by maximizing the log-likelihood function below: where the ϵ is the flux uncertainty and log(f ) is the logarithmic scaling factor for flux uncertainty.In total the PICASO/Elf Owl model fit has nine free parameters (five from the grid, and then the scaling factors for the flux uncertainty, radius, wavelength shift, and v sin i).

CHIMERA
CHIMERA (Line et al. 2013(Line et al. , 2014a,b,b) is an atmospheric retrieval framework that has been tested on brown dwarf and exoplanet spectra (e.g., Line et al. 2015;Zalesky et al. 2019;Hood et al. 2023).CHIMERA consists of 33 free parameters, including 15 pressure-temperature knots, gravity, radius, rotational broadening v sin i, radial velocity, cloud vertical sedimentation efficiency (f sed ), cloud volume mixing ratio, cloud base pressure, and uncertainty inflation parameters, two hyperparameters for thermal structures (γ and log(β)), and constant volume mixing ratios of H 2 O, 12 CO, 13 CO, CO 2 , CH 4 , H 2 S, NH 3 , and Na.CHIMERA adopts the Markov chain Monte Carlo (MCMC) method to explore the parameter space and estimate the posterior probability distribution of the parameters.We refer to the details of the model setup to Line et al. (2015).

Model Fitting Results
The best-fit model spectra are plotted in Figure 2. The best-fit CHIMERA and Elf-Owl model spectra give chi-squared values of 17180 and 42564 over 3233 data points respectively, corresponding to reduced chisquares of 5.40 and 13.2.We note that the apparent high reduced chi-squared is partially caused by the underestimated data uncertainty, which includes only photon, readout, and background subtraction noise, but not flatfield uncertainty.Including an additional 6% flat-field uncertainty (see Section 2.1) will lower the reduced chisquared for the best-fit CHIMERA and Elf-Owl models to 1.28 and 1.65, respectively.
In Figure 2, the NRS2 part of the spectra manifests larger residuals than in NRS1.We observe that the residuals of the CHIMERA retrieval not only gives lower chi-squared but also appear to be less wavelength correlated than that of Elf-Owl model's residual.This is not surprising because the CHIMERA model has over three times more parameters (33 vs. 10) than that of the Elf-Owl model.The fitted cloud volume mixing ratio (VMR) by CHIMERA is low, with a log(V M R cloud ) of −12.23 ± 5. Therefore, we interpret that there is no significant cloud opacity contribution to the 3-5 µm spectrum.The fitted v sin i values of CHIMERA and Elf-Owl are lower than the resolution limit of NIR-Spec/G395H (∼110 km/s for R∼2700) and should be interpreted as a free parameter that accounts for both instrumental and rotational broadening.
In the following subsections, we focus on the fitted thermal structure, atmospheric composition, the C/O ratio, metallicity, and isotopologue constraints of the best-fit CHIMERA and Elf-Owl models.We then inspect the residuals and search for potential minor species in Section 3.5.

Thermal structure and Contribution Function
We plot the retrieved temperature-pressure (TP) profile in Figure 3 with 68% and 95% confidence intervals.As indicated in Figure 3, the TP profile is tightly constrained in the 1-20 bar region.By integrating the spectra from 1 to 200 µm based on the CHIMERA retrieval results, we estimate that the effective temperature is 534 +8 −25 K.We note that our best-fit model indicates that the G395H spectrum covers about 47% of the bolometric luminosity.
Comparison of the TP profiles fitted with different modeling approaches provides a qualitative assessment on the potential systematic uncertainty due to the model assumptions.In the left panel of Figure 3, we also plot the TP profile of the best-fit Elf-Owl model.The best-fit Elf-Owl model has an effective temperature of 425 +0.46 −0.33 K and a log(g) of 4.3 ± 0.01.The uncertainty of the Elf-Owl models are much smaller than the grid spacing and likely underestimated.Given a fixed pressure, the temperature of the CHIMERA retrieval is cooler than that of the Elf-Owl model.Qualitatively, the two TP profiles are relatively similar in the 1-5 bar region and increasingly differ at higher pressures.Therefore, caution should be exercised when interpreting the fitted TP profiles at a pressure of five bars or higher.
To illustrate the pressures from which the flux at different wavelengths is emitted, we use PICASO to calculate the contribution function based on the fitted atmospheric composition and thermal structure by CHIMERA and Elf-Owl models.In Figure 3, we show the relative contribution function of the CHIMERA retrieval, which is normalized by the flux density.Based on the calculated contribution function, the G395H spectrum probes 0.6-19 bar region based on the bestfit Elf-Owl model (not shown in Figure 3) and probes around the 0.8-27 bar region based on the CHIMERA retrieval results.Therefore, the fitted Elf-Owl and CHIMERA model spectra probe similar pressure ranges.The contribution function results are also consistent with the tight constraints of the CHIMERA-retrieved TP profile in Figure 3.

Atmospheric Composition
In Figure 4, we plot the best-fit molecular abundances of the CHIMERA and Elf-Owl models.We also plot the pressure at which the optical depth of individual molecule reaches 0.01 in Figure 5 to visualize the contribution from each molecule on the observed absorption feature.The corner plot and marginalized distribution of the key parameters are shown in Appendix 5.In Table 1, we list the median values and 1-σ confidence intervals of the fitted molecular abundances.Guided by the marginalized distribution of molecules that show bounded lower and upper limits, we report the detection of H 2 O, CH 4 , 12 CO, NH 3 , H 2 S, and CO 2 in the G395H spectra.
To our knowledge, the only previous detections of H 2 S in ultracool brown dwarfs are Tannock et al. (2021); Hood et al. (2023) with ground-based medium-to-high resolution near-infrared spectroscopy.In Figure 6, we plot the increase of the residuals with the removal the H 2 S from the CHIMERA retrieval.Our study presents one of the first bounded constraints on the H 2 S abundance in a brown dwarf atmosphere.On the other hand,  the posterior distribution of 13 CO/ 12 CO ratio of the CHIMERA retrieval does not show a clear lower bound, so we report no detection for 13 CO in the data.For the best-fit Elf-Owl model that includes PH 3 , we notice that the absorption by PH 3 that is in chemical equilibrium causes the model to underestimate the flux density in the 4.1-4.3µm region.Meanwhile, the residuals of the CHIMERA retrieval, which does not include PH 3 , does not exhibit spectral features similar to that expected from PH 3 absorption at 4.3 µm.Therefore, we report no evidence of PH 3 absorption in the WISE 1828 spectra.
Overall, the CHIMERA retrieval gives higher ( > 1σ) molecular abundances of H 2 O, CH 4 , 12 CO, NH 3 , H 2 S, and CO 2 when compared to the best-fit Elf-Owl model grid.In particular, the CO 2 abundance, which is sensitive to the metallicity (e.g., Zahnle & Marley 2014), of the best-fit Elf-Owl model is almost five orders of magnitude lower than that of the CHIMERA retrieval (VMR(CO 2 of 10 −8.8 vs. 10 −12 ).Upon inspection of the residuals in the 4.2-4.3µm at which CO 2 absorption occurs, we find that CHIMERA retrieval residuals are lower than the Elf-Owl model.Furthermore, the residuals of the Elf-Owl model in the 4.1-4.3µm is dominated by the excess PH 3 absorption that overlaps with the CO 2 absorption in the 4.2-4.3µm region.Therefore, we argue that the fitted CO 2 abundance by the CHIMERA retrieval is more reliable than that by the Elf-Owl model.
In Table 2, we derive the relative molecular abundances, which is likely to be less sensitive to the potential degeneracy between the gravity and molecular abundance.Because the molecular abundance changes with altitude in the Elf-Owl models, we listed the molecular abundances at five bars, which is the averaged photospheric pressures probed by the 3-5 µm spectrum.The differences in the relative abundances are smaller than that in the absolute abundances between the two mod- els.Therefore, the relative molecular abundances are more robust against the model assumptions to reproduce the observed spectral features.

Elemental Abundance
We can derive atmospheric elemental abundances based on the fitted molecular abundances from the CHIMERA retrieval, while making certain assumptions about the background undetected gas species.We do not do this for the fitted Elf-Owl models because the elemental abundances of forward models are set by the assumed metallicity and C/O ratio; For the CHIMERA retrieval framework, we can derive the elemental abundance [C/H],[N/H], [O/H], and [S/H] based on the detected molecular abundances CH 4 , CO 2 , CO, NH 3 , H 2 S, and H 2 O.We calculate the hydrogen,carbon,nitrogen,oxygen abundance by the following equations: We then normalize the abundances by the corresponding solar abundance (e.In addition to the detected molecules, there are other non-detected molecules that are also important as the potential reservoir of C,N,S,and O elements.For example, the derived N/H is probably underestimated because of the missing N 2 abundance.Based on the chemical abundance of the best-fit Elf-Owl models, the N 2 abundance is about 60% of NH 3 in the 0.1-20 bar region.Therefore, the contribution from N 2 could increase the estimated [N/H] to +0.04, which is near solar [N/H] value.We refer to Section 4.6 for further discus-

C/O ratio and metallicity
The best-fit Elf-Owl model has a metallicity of [M/H]= −0.57± 0.01 and a C/O ratio of 1.01×Solar C/O value (0.458) (Lodders et al. 2009).The corner plot of the grid-model fitting shows that the metallicity is negatively correlated to the C/O ratio, as shown in Figure 7. Based on the elemental abundances of the CHIMERA retrieval, we calculate the C/O ratio to be 0.43 ± 0.01, similar to the solar C/O value.We also calculate the metallicity by summing the elemental abun-dance heavier than hydrogen.The derived metallicity [M/H] is equal to 0.30 ± 0.01 dex, or about two times solar metallicity.Therefore, the two models share similar C/O ratios, but the Elf-Owl model has a metallicity 0.87 dex lower than the CHIMERA retrieval.
Similar to the estimation of oxygen abundance, our estimation of C/O and metallicity did not account for the potential silicate condensate cloud formation.Therefore, the actual oxygen bulk abundance and the corresponding inferred metallicity should be higher than +0.30dex while the C/O ratio should be lower than 0.43.We note that the retrieved molecular abundances and metallicity by CHIMERA are correlated with gravity.We further discuss the interpretation and implications of the atmospheric composition, derived metallicity, and C/O ratio in Section 4.3 and 4.4.
3.4.3. 13CO/ 12 CO ratio The CHIMERA retrieval gives a log[ 13 CO/ 12 CO ( 13 CO/ 12 CO) solar ] of −1.61 +0.29 −0.39 , which corresponds to log( 13 CO) of -8.58.However, the marginalized posterior distribution of log( 13 CO/ 12 CO) does not manifest a clear lower boundary.To further study the impact of 13 CO opacity on the modeling results, we use PICASO to model the spectra with the best-fit 13 CO abundances (log( 13 CO)=-8.58)and with 5-σ 13 CO abundances (log( 13 CO)=-7.61) in Figure 8.In Figure 8, we plot the modeled spectra with two different 13 CO abundances.The residuals of the CHIMERA retrieval are within 7 σ in the 4.6-4.96µm region.Upon the manual inspection of the elevated residuals, we conclude that we can rule out a 13 CO enriched atmosphere with a log( 13 CO) of -7.61 because the additional 13 CO opacity increases the maximum residuals to 10 σ and increases the chi-squares by ∼1200.In conclusion, we do not detect significant 13 CO opacity in our data and place an upper limit of log( 13 CO) of less than -7.96, or 12 CO/ 13 CO higher than 40.We observe that there are numerous residuals in the model fit that exceed 5 σ levels for both CHIMERA and Elf-Owl/PICASO models (e.g., Figure 2).In several instances, the residuals exhibit wavelength-correlated behavior pointing to potential trace species that were not included in the model fit.We, therefore, examined the detected molecules in Jupiter and Saturn's atmosphere, which are about 200-300 K colder than WISE 1828, for potential molecules present in the atmosphere of WISE 1828.The 3-5 µm infrared spectra of Jupiter and Saturn show absorption features of a variety of trace molecular species, including PH 3 (4.6-4.8 µm) (Bezard et al. 1989), GeH 4 (4.737µm) (Bezard et al. 1989;Giles et al. 2017), AsH 4 (4.704µm,4.728µm) (Bezard et al. 1989;Noll et al. 1989;Giles et al. 2017), , and C 2 H 2 (∼ 3.7 µm, ∼ 4.6 µm) (Ridgway 1974;Rothman et al. 2005).However, we do not detect strong deviations at the line positions of those possible molecular species based on the residuals of the best-fit models.We find that the Elf-Owl model fitting residuals (data−model) have a ∼ 20σ positive deviation, which is higher than the averaged residuals (see Figure 2) at 4.55 µm.We attribute the positive residual to the inclusion of CH 3 D in the opacities used in PICASO.Specifically, the HAPI code (Roman et al. 2015) was used to compute the opacity of CH 4 , which automatically weights the isotopologues by their Earth abundances.For CH 3 D, this is 1% the abundance of CH 4 2 .On the other hand, although CHIMERA uses identical CH 4 line lists as PICASO (Hargreaves et al. 2020), it only includes the main isotopologue.CH 3 D is known to have a distinct absorption feature at 4.55 µm (e.g., Morley et al. 2018).Therefore, the strong (∼20 σ) positive residual of the Elf-Owl best-fit model at 4.55 µm suggests that WISE 1828's atmosphere is likely to have less than 1% CH 3 D/CH 4 .This result is consistent with CHIMERA retrieval's residuals that do not exhibit high residuals at 4.55 µm because CHIMERA does not include CH 3 D in the CH 4 opacity calculation.The molecular abundances of the Elf-Owl grid models are not free parameters and are altitude dependent.We listed the Elf-Owl model abundances at five bars, which is the averaged peak locations of the contribution function in Figure 3.The molecular abundances of CHIMERA are constant with altitude.See the text in Section 3.1 for the description of the estimated 68% confidence intervals (C.I.) for the Elf-Owl model parameters with a bootstrapping method.The presented 3-5 µm spectrum of WISE 1828 with JWST has demonstrated the value of R∼3000 midinfrared spectroscopy for atmospheric characterization of brown dwarfs.In the cool Y-dwarf atmosphere, the major carbon-bearing species is CH 4 .The 3-5 µm spectrum simultaneously constrains CH 4 , CO 2 , and CO which is useful for understanding the non-equilibrium chemistry in the Y-dwarf atmosphere.The best-fit Elf-Owl models suggest a moderate vertical mixing strength with a log(K zz (cm 2 s −1 )) of 4. Our fitted vertical mixing is comparable with the fitted vertical mixing strengths in the legacy Spitzer study of L and T dwarfs (Stephens et al. 2009).More atmospheric studies of Y dwarfs will be necessary to probe the variations of vertical mixing strengths and the corresponding non-equilibrium chemical composition in these cool, methane-dominated atmospheres.

Search for possible trace species
The acceptable fit of the non-equilibrium model fitting results by the Elf-Owl models also validates the Mukherjee et al. ( 2022)'s non-equilibrium chemical model in the cool Y-dwarf atmospheric regime.The lack of PH 3 absorption features in the 4.2-4.3µm region, where PH 3 is predicted by non-equilibrium chemistry (Visscher et al. 2006), suggests that our understanding of phosphorus chemistry is incomplete.

Comparison to evolutionary models and potential binarity
Brown dwarf evolutionary models (e.g., Phillips et al. 2020;Marley et al. 2021) simulate the evolution of brown dwarf atmospheres with different metallicities and masses, and predict their radius and temperature changes with age.The Sonora Bobcat evolutionary models (Marley et al. 2021) suggest that WISE 1828 has a mass of 9.982 M Jup and an age of 1.4Gyr based on the fitted Elf-Owl temperatures and gravity.The CHIMERA results suggest that the best-fit radius and gravity are 1.23 ±0.01 R Jup and log(g) of 5.21±0.01respectively.The fitted gravity by CHIMERA is higher than the expected gravity at 10 Gyr by the evolutionary models (see Figure 9).Therefore, the fitted Elf-Owl grid-model solution, which has a lower metallicity than the CHIMERA retrieval, is more consistent with the evolutionary model prediction.We defer to Section 4.3 for a more detailed discussion on the impact of metallicity on the results.The inconsistent gravity between atmospheric retrieval frameworks and evolutionary models is not uncommon, and has been reported in the atmospheric study of other brown dwarfs (e.g., Zalesky et al. 2019;Kitzmann et al. 2020).
Based on the Bobcat brown dwarf evolutionary models, a brown dwarf with a T ef f of 533 K, solar or 2x solar metallicity, and a log(g) of around 5.08 has a radius of around 0.87 R Jup and 33M Jup .The fitted radius is therefore about 40% higher than the evolutionary model prediction.One possible explanation for the discrepancy between our results and the evolutionary model radius is that WISE 1828 is an unresolved equal-mass binary.In the binary scenario, the radius of WISE 1828's single component will be equal to the square root of the bestfit radius, or 0.87 R Jup .A radius of 0.87 R Jup is then roughly consistent with the evolutionary model radius at 10Gyr.We also inspect if the fitted radial velocity is consistent with a potential binary scenario.We estimate the semi-major axis based on the following assumptions as to the inclination, mass, and semi-amplitude of the radial velocity of the potential binary.We first assume that the possible inclinations of the binary system follows an isotropic distribution.We assume the mass of each binary component is 33M Jup as indicated by the CHIMERA best-fit T ef f and gravity.Finally, we assume that the binary system has a semi-amplitude of 32.0 ± 0.2 km/s, which is the same as the CHIMERA best-fit radial velocity.Based on these assumptions, our Monte Carlo estimation results suggest that such a potential binary system would have a semi-major axis of 0.0098 +0.002 −0.006 au, or 20 +4 −12 R Jup .This derived semi-major axis is smaller than the 0.5 au upper limit estimated by De Furio et al. ( 2023).However, this order-of-magnitude estimation does not account for the radial velocity of the potential binary system and the uncertainty of the absolute wavelength calibration, which is below 14km/s3 .Future multi-epoch observations will be essential to detect potential radial velocity variations of WISE 1828.

Caveats of atmospheric retrieval results
Our CHIMERA retrieval results provide a decent fit to the 3-5 µm spectrum with a reduced chi-squared of 5.4.The retrieved gravity and radius are higher than the evolutionary model prediction.Here, we discuss the impact of potential bias in the retrieved gravity and radius on the results.
Our retrieved molecular abundances are positively correlated with gravity (see the figure in Appendix 5).The correlation between gravity and abundance is because the photospheric pressure (P(τ λ )=1) is proportional to the gravity and opacity ratio under the assumption of hydrostatic equilibrium: Therefore, the photospheric pressure remains the same when both the volume mixing ratio and gravity increase by the same fractional amount.Indeed, the molecular abundances except CO 2 (see Section 3.4 for the discussion of CO 2 abundance) of the Elf-Owl model are about 0.5-1.0dex lower than that of CHIMERA, while the gravity of Elf-Owl models is also about 0.8 dex lower than the CHIMERA retrieval (see Table 1), similar to the expected correlation based on the analytical equation above.The interpretation of the effective temperature, metallicity, and gravity is less certain because of the degeneracy between gravity and molecular abundance.The best-fit Elf-Owl results present a low gravity, low temperature, and sub-solar metallicity scenario while the CHIMERA retrieval points to a high gravity, high temperature, and above-solar metallicity scenario.The H 2 -H 2 and H 2 -He collision-induced absorption (CIA) (e.g., Saumon et al. 2012), whose absorption features fall outside of the data wavelength coverage in this study, are sensitive to the photospheric pressure and gravity.Therefore, further spectral analysis that extends to a wider wavelength coverage may help better constrain the gravity and distinguish the two scenarios.
The C/O ratio is less correlated to the fitted gravity because the abundance of carbon-and oxygen-bearing species are both positively correlated with gravity, therefore the C/O ratios should remain roughly the same even if the gravity changes.Indeed, the derived metallicity and C/O ratios from the MCMC samples do not strongly correlate with each other, as shown in Figure 10.For the same reason, the relative abundances between two molecules in our retrieval results are also less correlated to the retrieved gravity.

C/O ratio and metallicity in the solar neighborhood
The C/O ratio and metallicity measured from planetary atmospheres are proposed as potential tracers of planet formation scenarios ( Öberg et al. 2011;Madhusudhan et al. 2011).However, connecting atmospheric abundances to the complex and under-constrained planet formation scenarios is challenging (e.g., Mollière et al. 2022).As a brown dwarf's mass range overlaps with the low-mass end of star formation and the highmass end of planet formation, measurements of C/O and metallicity of a large sample of brown dwarfs will test if planet formation and star formation pathways could lead to different atmospheric properties.
We compare the derived C/O ratio and metallicity with the retrieval results in Zalesky et al. (2019).In Zalesky et al. (2019), they use CHIMERA to retrieve the C/O ratios and metallicity of 14 late-T and Y dwarfs based on the HST spectra (Schneider et al. 2015).The retrieved C/O ratios range from 0.25 to 0.67 with a metallicity range from around -0.4 to +0.6 dex.Therefore, our retrieved metallicity and C/O ratios of WISE 1828 are consistent with that from Zalesky et al. (2019), as well as that of solar neighborhood stars of Hinkel et al. (2014) (see Figure 5 in Zalesky et al. 2019).The retrieval results in Zalesky et al. (2019) have C/O uncertainties of around 0.1 or around 20% based on the HST 1.1-1.7 µm spectra, and our results with the same retrieval framework has a C/O ratio uncertainty of 0.01 or 2% level, which is an order-of-magnitude lower than the previous results.Our results illustrate that we can constrain the C/O ratios and metallicity to percent-level precision with the JWST spectroscopy.4.5.Implication of 12 CO/ 13 CO ratio Recent 13 CO detections in young giant planets and brown dwarfs (Zhang et al. 2021b,a;Gandhi et al. 2023) present an exciting opportunity to study the potential value of isotopologues in connecting the atmospheric properties to planet and star formation pathways.A brown dwarf formed through gravitational collapse formation pathways should share similar isotopologues with the ISM values, which are set through processes including ion-exchange reactions (Langer et al. 1984), isotope-selective photodissociation (Bally & Langer 1982;van Dishoeck & Black 1988;Visser et al. 2009), and gas-to-ice isotopologue partitioning (Acharyya et al. 2007;Smith et al. 2015).If a giant planet is formed in a disk, then additional processes, such as the dynamical CO depletion that varies across disk radius (Zhang et al. 2019), could also be important and drive the isotopologue ratios away from the ISM values (Yoshida et al. 2022).Furthermore, the isotopologue ratios of a planet may also change with atmospheric evolution (for reviews, see Nomura et al. 2022).
Our results place a lower limit of the 12 CO/ 13 CO of 40.The isotopologue ratio lower limit of 40 is about 1σ higher than that of the young giant planet TYC 8998-760-1 b (31 +17 −10 Zhang et al. 2021b).The lower limit is consistent with the isotopologue ratio of the Sun (86.8 ± 3.9 Scott et al. 2006;Asplund et al. 2009), the brown dwarf 2MASS J03552337+1133437 (97 +25 −18 Zhang et al. 2021a), and with the local ISM value within 1kpc of around 30-180 (see Figure 5 in Smith et al. 2015).More modeling and observational studies of nearby stars and planets (e.g., Crossfield et al. 2019;Zhang et al. 2021b,a) are required to identify the dominating processes that shape the isotopologue ratio.Our JWST results demonstrate the potential of measuring the isotopologue ratio of cold giant planets for identifying potential trends and building up the sample size of planetary-mass objects with detected isotopologue.4.6.Comparison to Barrado et al. 2023 results Barrado et al. (2023) (hereafter B23) presented the JWST Mid-Infrared Instrument (MIRI) Medium-Resolution Spectrometer (MRS) 5-20 µm spectra of WISE 1828 and reported the detection of 15 NH 3 .Therefore, both 15 NH 3 and 13 CO isotopologue constraints are available for the WISE 1828 atmosphere.B23 also reported the atmospheric properties derived from five different atmospheric retrieval frameworks.. Four of the five retrievals include both the MIRI spectra and the Hubble Space Telescope 1.1-1.7 µm spectra.We summarize the key results in B23 and compare them with our results in Table 3.The derived gravity, metallicity, molecular abundances in B23 are consistent with our results within two-sigma levels.B23's retrieval results prefer a lower effective temperature and a larger radius compared to our results.The implied bolometric luminosity of B23, which is proportional to temperature T 4 and radius r 2 , is similar to that of the Elf-Owl model within around 10% level.The temperature-pressure profile of B23 is similar to our results in the 1-10 bar region in Figure 3. Therefore, the B23 results show similar chemical composition to both CHIMERA and Elf-Owl model results.The best-fit effective temperature of B23 is lower than that of both CHIMERA and Elf-Owl models, but the bolometric luminosity is similar with that of the Elf-Owl models.In this section, we will discuss our findings within the context of prior modeling studies concerning WISE 1828 as reported by Cushing et al. (2021) and Leggett et al. (2021).
In their study, Leggett et al. (2021) assumed that WISE 1828 is a binary system.They employed modified T-P profile models developed by (Tremblin et al. 2015(Tremblin et al. , 2019) ) to fit the Spectral Energy Distribution (SED) of .5.When considering the possibility of an unresolved binary system, their best-fit Sonora Bobcat binary model indicated a T ef f of 300 K and 350 K, log(g) of 4.0, metallicity of -0.5 and -0.5 dex, and C/O ratios of 0.5 and 0.25 times the solar value, respectively.One significant factor contributing to the disparities between our results and those of previous studies likely arises from differences in the wavelength coverage of the dataset.Previous modeling efforts primarily centered on fitting models to the near-infrared spectrum and midinfrared photometry.In contrast, our analysis focused on fitting the models to the 3-5 µm spectrum.Given that the 3-5 µm spectrum is particularly sensitive to the abundance of CO 2 , CH 4 , CO, and H 2 O, our findings are expected to provide a more precise estimate of the C/O ratio compared to prior studies.
It's important to note that the near-infrared spectrum probes pressure levels higher than those probed by the mid-infrared spectrum, as illustrated in Figure 3. Consequently, the T-P profile derived from the previous modeling studies should offer a more accurate estimate of the temperature at pressures greater than around twenty bars.A combined analysis incorporating archival near-infrared and mid-infrared spectra would provide a more comprehensive and less biased view of the temperature-pressure profile of WISE 1828.

CONCLUSIONS
Our study marks one of the first in-depth atmospheric investigations of a Y-dwarf using JWST NIRSpec spectroscopy.Our rich results provide a glimpse of the powerful capabilities of JWST in revolutionizing our understanding of cold giant planet atmospheres.We summarize our key findings as follows: 1. We present the JWST NIRSpec/G395H 2.880-5.142µm spectrum of WISE 1828.Our moderate spectral resolution (R∼ 2700) spectra of WISE 1828 cover ∼47 % of the bolometric luminosity of Y-dwarf WISE 1828 with a peak signalto-noise ratio exceeding 90.
2. We utilize two complementary atmospheric modeling tools, the CHIMERA atmospheric retrieval framework and the Elf-Owl radiative-convective equilibrium grid models, to characterize the atmospheric structure and composition of WISE 1828.
3. We find that the best-fit Elf-Owl grid model has an effective temperature of 425 percent-level precision of C/O ratio and metallicity.We caution that the metallicity is highly model dependent while the fitted C/O derived from the CHIMERA retrieval and Elf-Owl models agree well to each other within 7%.
6.We attribute the difference in the forward modeling and the atmospheric retrieval results to the degeneracy between gravity and molecular opacity.We derived relative abundances of WISE 1828 that are less dependent on the fitted gravity parameter in the atmospheric models.Future observations with a wider wavelength coverage will be useful to disentangle the degeneracy and provide a more accurate estimate of the effective temperature, elemental abundances, metallicity, and gravity of WISE 1828.
7. The retrieved radius and gravity are both higher than that of the brown dwarf evolutionary models.One possible explanation is that WISE 1828 is a tightly bound binary pair that share similar temperatures.We discuss the implication of the potential binarity based on the fitted results and evolutionary models.Our crude estimation of the semi-major axis is consistent with previous studies.

Figure 1 .
Figure 1.The reduced spectrum of WISE 1828.The parts of the spectra covered by the NRS1 and NRS2 detectors are shown in the top and bottom panels, respectively.The one-sigma uncertainties are plotted as grey shaded regions.

Figure 2 .
Figure 2. Top panel: the best-fit Elf-Owl model (semi-transparent green lines) and the CHIMERA atmospheric retrieval (semi-transparent orange lines) spectra.For clarity, the Elf-Owl model spectra and a copy of data spectra are shifted by −8 × 10 −11 Wm −3 in the y-axis.The original and shifted data spectra are plotted as semi-transparent blue lines for comparison with the CHIMERA and Elf-Owl models respectively.Middle and bottom panels: The residual for the CHIMERA and Elf-Owl models plotted in orange and green lines respectively.The sharp positive residual at 4.55 µm of the Elf-Owl model is because of the lower-than-solar CH3D abundance (See Section 3.5).

Figure 3 .
Figure 3. Left panel: the temperature-pressure profile of CHIMERA (T ef f = 534 +8 −25 K, log(g) = 5.02±0.01), the best-fit Elf-Owl model (T ef f = 425 +0.5 −0.3 K, log(g) = 4.38 ± 0.01), and of the petitRADTRANS retrieval in Barrado et al. (2023).The colored dots indicate the pressures at which the local temperatures equal to the effective temperatures of the models.See Section 4.6 for the discussion of the thermal structure.Right panel: the contribution function of CHIMERA.The 3-5µm spectrum probes around 0.8-27 bar.The probed pressure range corresponds to the region in which the temperature-pressure profiles have small uncertainties.

Figure 4 .
Figure 4.The fitted chemical abundances of the forward model Elf-Owl and the CHIMERA retrieval framework.The dashed and solid lines show the abundance of the best-fit model spectrum from the Elf-Owl models and CHIMERA, respectively.

Figure 5 .
Figure 5.The spectrum is plotted in dark-grey lines along with the P(τ =0.01) of the detected molecules.The left y-axis shows the flux density of the spectrum while the right y-axis shows the photospheric pressure at which the molecular optical depth reach 0.01.The three most abundant molecules CH4, H2O, and NH3 are plotted in thin light blue, green, and orange lines respectively.The three less abundant molecules H2S, CO2, and CO are plotted in thick red, purple, and blue lines respectively.The P(τ =0.01) lines of molecules are binned down into resolution of R=1350 for clarity.

Figure 6 .
Figure 6.Top panel: The best-fit CHIMERA model spectrum with and without H2S opacity are plotted in blue and orange lines, respectively.The data is plotted in black line.Middle panel: The fractional flux density change of the best-fit spectrum after the removal of H2S opacity.Bottom panel: the residuals relative to the data uncertainties for the CHIMERA retrieval results with and without H2S opacity are plotted in blue and orange lines.Exclusion of H2S opacity leads to an increase of χ 2 by over 900.sion about ammonia abundances and nitrogen isotopologue.The [O/H] is also likely underestimated because of silicate cloud (e.g., Mg 2 SiO 4 ) formation, which could lower the detected [O/H] by ∼30% and increase the C/O ratio by ∼ 40% (e.g., Zalesky et al. 2019; Calamari et al. 2022).After considering the effect of other nondetected molecules, our results suggest that WISE 1828 has above-solar [O/H], [S/H], and [C/H] abundances and likely near-solar [N/H] abundance.

Figure 7 .
Figure 7.The corner plot of the Elf-Owl model grid fitting results.The nine model parameters are effective temperatures, gravity (log(g)), vertical mixing log(Kzz), metallicity log(M/H), C/O ratio, radius (RJup), wavelength shift in unit of Å, and log(v sini km/s).

Figure 8 .
Figure 8. Top panel: The retrieved spectrum with the best-fit and the 5σ upper limit 13 CO abundance; Middle panel: the change of the spectra after removing 13 CO and 5% 12 CO opacity; Bottom panel: the spectral difference between the two model spectra relative to the observational data uncertainty.

Table 2 .
The relative chemical abundances of CHIMERA and the Elf-Owl models CHIMERA Elf-Owl (5 bar) log(CH4/H2O) −0.36 ± 0Importance of detected molecular abundances for atmospheric physics and chemistry of cool giant planets

Figure 9 .
Figure9.Based on the Sonora Bobcat evolutionary models, the fitted CHIMERA temperature and gravity implies an age of older than 10Gyr and mass of around 33 MJ , or the fitted gravity is likely unphysically high (seeSaumon & Marley 2008).The fitted Elf-Owl temperature and gravity indicates an age of 1.4 Gyr and 10 MJ .The colored contour lines show the age contour of a brown dwarf while the grey contour lines show the mass at an effective temperature and gravity.The solid and dashed lines represent the evolutionary model predictions for half-solar and solar metallicity objects respectively.

Figure 10 .
Figure 10.The metallicity and C/O ratios posterior distributions do not strongly correlate to each other.

Figure 11 .
Figure11.The corner plot of CHIMERA retrieval results.The parameters for TP profiles are excluded for clarity.See text for the explanation of the listed parameters.

Table 1 .
The fitted key model parameters of Elf-Owl grid models and the CHIMERA retrieval framework

Table 3 .
Barrado et al. (2023)eling results toBarrado et al. (2023) times the solar C/O ratio, a T ef f of 350 K, and a log(g) of 4.0.The best-fit cloudy model, on the other hand, exhibited a T ef f of 275 K and log(g) of 4 +0.46 −0.33 K, log(g) of 4.38 ± 0.01, log(K zz ) of 4.65 ± 0.04, C/O ratio of 0.46 ±0.01, and metallicity [M/H] of −0.57± 0.01.The CHIMERA retrieval has an effective temperature of 534 +8 −25 K, metallicity of +0.30 +0.01 −0.02 dex, and a C/O ratio of 0.43 ± 0.01.The reduced chisquared values of the best-fit CHIMERA and Elf-Owl models are 5.4 and 13.2 respectively.While the CHIMERA retrieval has a lower reduced chisquared value, the fitted gravity of the Elf-Owl model is more consistent with the Bobcat evolutionary model prediction.4. Based on the CHIMERA retrieval results, we reported the detection and bounded abundance constraints of CH 4 , H 2 O, CO, NH 3 , H 2 S, and CO 2 molecules.We report the first measurement of the H 2 S abundance in a Y-dwarf atmosphere.Our residual analysis finds no evidence of trace molecular species PH 3 , GeH 4 , C 2 H 2 , AsH 4 .We do not detect 13 CO in the G395H data with CHIMERA.We place a lower limit in the 12 CO/ 13 CO isotopologue ratios of 40.We derive the elemental abundances of [C/H], [N/H], [O/H], and [S/H].We find qualitative evidence of a lower-than-nominal CH 3 D abundance (i.e., CH 3 D/CH 4 < 1%) in WISE 1828 based on the Elf-Owl best-fit model's residuals at 4.55 µm. 5.The CHIMERA retrieval and Elf-Owl model fitting results demonstrate that the capability of JWST in characterizing Y-dwarf atmosphere to