© The Authors 2026

1 Introduction

The atmospheric physics and the oxygen, nitrogen, and carbon chemistry of cold worlds beyond our Solar System contextualize the corresponding properties of Jupiter. For a long time, we have not been able to characterize such cold worlds in detail due to their faint fluxes in the accessible wavelengths. With the advent of the James Webb Space Telescope (JWST; Gardner et al. 2023; Rigby et al. 2023), we are now able to obtain their spectra in the mid-infrared (MIR) in a unique spectral resolving power and signal-to-noise (S/N), providing valuable information about the composition. Moreover, this ability allows us to unravel not just the features of the major absorbing species but also their faint isotopologs.

Understanding the formation of gas giants similar to Jupiter requires profound knowledge of the properties of such objects (e.g., Burrows et al. 1997; Knierim & Helled 2025). Recent discoveries of cold gas giant companions with direct imaging by the Mid-InfraRed Instrument on board JWST (JWST/MIRI, Wright et al. 2015; Wright et al. 2023) have populated a previously inaccessible region of the parameter space, including Eps Ind Ab, with an estimated effective temperature of ~275 K (Matthews et al. 2024); the similarly cold companion 14 Her c (Bardalez Gagliuffi et al. 2025); and a sub-Jovian planet in the young TWA 7 debris disk (Lagrange et al. 2025). More detailed characterizations of these cold worlds are also possible, as shown by the clear detection of ammonia in the atmosphere of GJ504b (Mâlin et al. 2025).

All of these discoveries provide promising insights into the formation, evolution, and dynamics of Jupiter-like planets, but in practice, the process of extracting high S/N spectra from direct observations suffers from the contrast to their host star. Thus, objects orbiting far away from their host stars are promising and easier targets for characterization of their atmospheric composition and structure. Further, free-floating brown dwarfs can also serve as laboratories to study the atmospheric properties of cold objects, such as their compositions, temperatures, and the governing physical processes that are thought to be comparable with cold gas giants (e.g., Burrows et al. 1997; Coulter et al. 2022).

While the spectra of brown dwarfs at the L/T transition have been studied intensively (e.g., Miles et al. 2020; Vos et al. 2023; Biller et al. 2024; Mollière et al. 2025; Nasedkin et al. 2025), the T/Y transition is still rather unexplored. With JWST we are now able to probe this cold end of the brown dwarf population (Barrado et al. 2023; Beiler et al. 2024a,b; Matthews et al. 2025; Kühnle et al. 2025; Vasist et al. 2025; Lueber et al. 2026). The chemistry of brown dwarfs at the T/Y transition with temperatures of ~500 K (Cushing et al. 2011) is dominated by water, ammonia, and methane, and it shows strong absorption lines in the MIR. To date, only one T-dwarf spectrum measured with the Medium Resolution Spectrometer (MRS) has been published. The first T dwarf observed with MIRI/MRS showed the unexpected presence of HCN and C₂H₂ (Matthews et al. 2025) and has a similar temperature compared to COCONUTS-2b.

A newly emerging avenue to probe the formation paths of cold brown dwarf and gas giant atmospheres has been enabled by studying isotopologs (Barrado et al. 2023; Gandhi et al. 2023; Kühnle et al. 2025; Ruffio et al. 2026). The unprecedented resolving power of MIRI/MRS (Wells et al. 2015; Wright et al. 2015; Argyriou et al. 2023) allows one to do so in the MIR. Isotopic ratios have been proposed as a formation tracer for cold gas giants and brown dwarfs, as it may allow us to link the atmospheric composition to the disk composition during formation (Mollière & Snellen 2019; Zhang et al. 2021a,b; Barrado et al. 2023; Nomura et al. 2023). Oxygen isotopes are abundant in the interstellar medium (ISM) and the Solar System, and they have been measured extensively (e.g., Wilson 1999; McKeegan et al. 2011). The local region in which our Solar System is located has been found to be enriched in ¹⁸O in comparison to the average ISM from ejecta from type II supernovae (Young et al. 2011). Gandhi et al. (2023) detected the oxygen isotopologs of CO in the planetary mass object VHS1256b, and Ruffio et al. (2026) also detected the same in the HR8799 planets c, d, and e.

Typically, these compositional results are derived using atmospheric retrievals, which have been widely used as a tool for the characterization of exoplanets and brown dwarfs (e.g., Madhusudhan & Seager 2009; Mollière et al. 2020; Barrado et al. 2023; Vos et al. 2023; Faherty et al. 2024; Kothari et al. 2024; Nasedkin et al. 2024a; Rowland et al. 2024; Matthews et al. 2025; Vasist et al. 2025; Nasedkin et al. 2025). Here, we present a retrieval analysis of the MIRI/MRS spectrum of the widely separated companion COCONUTS-2b using the publicly available retrieval code petitRADTRANS (Mollière et al. 2019; Blain et al. 2024; Nasedkin et al. 2024b).

COCONUTS-2b is one of the planetary-mass objects located furthest from its host star, and thus it is a prime target to study the atmosphere of a wide-separation planet (Zhang et al. 2025b). COCONUTS-2b was discovered as one of the first WISE targets by Kirkpatrick et al. (2011) as a field brown dwarf. Zhang et al. (2021c) showed using the GAIA early data release 3 and as part of the COol Companions ON Ultrawide OrbiTS survey that it is bound to the M3V dwarf L34-26 (also known as COCONUTS-2A), with a projected separation of ~6471 au or 594", and it is at a distance of 10.888 ± 0.002 pc (Bailer-Jones et al. 2021). The study of Zhang et al. (2025b) presented a Gemini/FLAMINGOS-2 on the Gemini-South telescope spectrum ranging from 1-2.5 μm, with an average resolving power of R~900. They used thermal evolution models, COCONUTS-2b’s estimated age of 150-800 Myr, and bolometric luminosity (log(L_bol/L_⊙) = −6.18 ± 0.25 dex) to constrain its bulk parameters, including an effective temperature of T_eff = ${483}_{-53}^{+44}$K, a surface gravity of log(g) = ${4.19}_{-0.13}^{+0.18}$ dex, a radius of $\mathrm{R}={1.11}_{-0.04}^{+0.03}{\mathrm{R}}_{\mathrm{J}}$, and a mass of M = 8 ± 2M_J (Zhang et al. 2025b). They also compared the Gemini/FLAMINGOS-2 spectrum of COCONUTS-2b to 16 grids of atmospheric models to estimate its subsolar to near-solar metallicity and C/O ratio. The latter classifies COCONUTS-2b as a planetary-mass object. Based on its J-, H-, and K-band flux, it was classified as a T9.5 dwarf and is thus spectroscopically right at the transition between T and Y dwarfs. In addition, CθCONUTS-2b shows signs of disequilibrium chemistry and clouds (Zhang et al. 2025b). Recently, Kiman et al. (2026) found that the COCONUTS-2 system is likely a member of the Corona of Ursa Major moving group and thus has an age of about 414±23Myr (Kiman et al. 2026). Using this age and an updated bolometric luminosity, they found a more precise mass for the companion of 7.5 ±0.4 M_J. The composition and bulk properties of the host star have not been widely studied. Kiman et al. (2026) constrained the host stars mass to be 0.4${}_{-0.02}^{+0.01}$ M_⊙, and Hojjatpanah et al. (2019) find the metallicity to be consistent with the solar value (with Z = 0.00±0.08). The lack of information on the host stars composition makes comparisons challenging, especially with respect to isotopic ratios. However, precise measurements of the ISM and the Sun as well as recent discoveries in other brown dwarfs and gas giants provide crucial references.

Both the atmospheric composition and the large separation of COCONUTS-2b from its host star raise the question of the formation of the companion. Possible scenarios include the following: (1) The companion formed in a manner similar to a star through gravitational collapse (Boss 1997). (2) It formed inside of a disk of the host star through gas accretion (Pollack et al. 1996) and migrated outward to its position observed today. (3) It evolved completely independent of the host star and was captured during a flyby (Zhang et al. 2021c; Marocco et al. 2024). Free-floating planets might be much more frequent than expected in a star-forming region (Miret-Roig et al. 2022), lending credence to the possibility of such a capture event. However, the COCONUTS system belongs to a looser moving group, making a flyby capture rarer. Further, the bulk parameter estimates are consistent with the age of the star and compatible with a gravitational collapse or a disk formation scenario. They are in line with both hot- and cold-start models, and thus we are not able to rule out either scenario (Zhang et al. 2025b).

In this paper, we present the MIRI/MRS spectrum of COCONUTS-2b reaching from 4.9 to 18 μm with a resolving power of up to R~3750. This dataset provides the highest resolving power spectrum in this wavelength range observed of COCONUTS-2b to date, and it allowed us to characterize the chemical composition and atmospheric structure of the companion using full resolution atmospheric retrievals. We introduce the observation and data reduction in Section 2. The retrieval setups are presented in Section 3. We present our results of the analysis in Section 4 and discuss them in 5. In Section 6 we conclude this work with a summary and provide an outlook for future work.

Fig. 1 Observation and cube detector images of channel 1A to 4A. Channel 4A was neglected from further analysis (see text). The red cross corresponds to the position of the target and the blue circle to the aperture from which the flux has been extracted from. The two negative spots on the detector correspond to two nods originating from the background subtraction. The yellow circle shows the position of a background object.

2 Observation and data reduction

COCONUTS-2b was observed on 05 September 2024 as part of the GO proposal number 6463 (PI: P. Patapis) using the MIRI/MRS instrument onboard of JWST. The observations consisted of 200 groups per integration, one integration per exposure, and one exposure. In total, 3.7 hours of telescope time were allocated for the observations and all three channels (SHORT, MEDIUM, LONG) were used to obtain a complete spectrum. To improve pointing accuracy, target acquisition was performed on the object. The science observation was performed in a negative 4-pt dither pattern in the FASTR1 readout. The observed spectrum with the corresponding detector images at 5, 10, 15 μm, and the summed channel 3B from 15.5 to 17.9 μm are presented in Fig. 1.

The data reduction has been carried out using the JWST pipeline (Bushouse et al. 2025) version 1.18.0, CRDS context files jwst_1364.pmap (version 12.1.6). Stage 1 was reduced from the _uncal files using the standard steps (including the ramp fitting, dark current subtraction, saturation, linearity, and jump corrections). Here, data gets processed from raw digital numbers (DN) to _rate files in DN/s. Stage 2 includes the flat field, stray light, fringe, and photometric correction. We added the “clean showers” option to the stray light correction, which improved significantly remaining striping artifacts due to cosmic rays that were not mitigated by the jump correction in stage 1. We use the nod subtraction method on the detector images to remove residual background. Here we subtract the dither pairs with the largest intersection to minimize PSF overlap from each other and replace the original image with the subtracted one. Together with the clean showers step mentioned above, this results in an exceptionally well subtracted background for this point source as seen in the cube images in Fig. 1. We included the standard fringe correction, but not the residual fringe correction. We tested for the latter, and the corrections were mainly within the pipeline error bars. Due to the small effect on the spectrum and as it may as well remove actual spectral features (Gasman et al. 2023), we neglect it in this reduction. Finally, stage 3 creates a cube from the two-dimensional detector image using the “drizzle” algorithm. To extract the spectrum we use the position in the cube determined by the extract_1d.ifu_autocen() function for all channels. As the S/N dropped significantly in channel 4 we decided to use only up until band 3C for this analysis. As COCONUTS-2b is a quasi-free-floating object and the host star is well outside the field of view, there is no residual light that needs to be subtracted. Thus, we set a circular aperture around the source with an radius of 1.5 × FWHM of the Point Spread Function (PSF). This number was chosen small enough to not include too much background systematics and large enough to limit sampling artifacts (Law et al. 2023). The cubes show an additional source on the top right compared to the point source marked in a yellow circle. This likely is a background galaxy that is bright in the MIR wavelengths, but not in the near infrared.

In the analysis presented, we include the archival Gemini/FLAMINGOS-2 observation presented in Zhang et al. (2025b) to complete the spectrum in the near-infrared (NIR). We exclude regions in the spectrum with low S/N, namely areas between 1.12-1.18 μm, 1.35-1.45 μm, 1.66-1.95 μm, and 2.3-4.95 μm in the retrieval analysis in line with Zhang et al. (2025b). More information on the data reduction can be found in Zhang et al. (2025b). Thus, we obtain a medium resolution spectrum from 1 to 18 μm with gaps between 2.3 and 4.9 μm and in the NIR.

3 Methodology

In this study, we used the atmospheric retrieval framework petitRADTRANS (Mollière et al. 2019; Blain et al. 2024; Nasedkin et al. 2024b, version 3.1.2) to characterize the atmosphere of COCONUTS-2b in terms of chemistry and physical properties based on the MIRI and Gemini data. In our initial step, we present two retrievals on a binned spectrum of spectral resolution up to $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$ using c-k opacity. The two retrievals vary by their pressure-temperature (PT) parameterization: first, we retrieved the temperature using ten pressure nodes, equally spaced in logarithmic pressures. The temperature at each node is parameterized as a free parameter, a factor a of the temperature at the node below (see Barrado et al. 2023; Kühnle et al. 2025). Note, we allowed inversions by choosing the prior range for each temperature node to be between [0.2, 1.2]. We interpolated between the nodes using a spline function. We did not regularize the profile, so no penalization on the change in slope is implemented. This is denoted as ‘free’ in the following. Secondly, we used another PT parameterization incorporating knowledge from self-consistent models to the retrievals as shown in Zhang et al. (2023). The retrieved values correspond to a deep temperature at the highest pressure node and ten nodes with values for dlogT/dlogP. We derived the nodes from the Sonora Elf Owl grid models (Mukherjee et al. 2024) with T_eff in [275, 800] K, logg in [3.25, 4.5] dex, [M/H] in [−1,+ 1] dex, CO in [0.5, 1, 2.5] times solar, and logarithmic K_zz in [2,4,7,8,9] according to the procedure presented in Zhang et al. (2023) for pressure levels of 10³−10⁻³ and has been extended for lower pressures of 10⁻⁶ bar in Zhang et al. (2025a). We present the priors in Table 1 and the retrieval using this parameterization is denoted as ‘constrained’ in the following. We vary the PT parameterization here to further analyze the effect of the change in PT structure on retrieval runs at full resolution.

In both retrievals, we included the treatment of error inflation to capture uncertainties in the model and a potential underestimation of the errors on the data. We retrieved for a value b, which is added to the measured error: ${\sigma }^2={\sigma }_{\mathrm{d}\mathrm{a}\mathrm{t}\mathrm{a}}^2+{10}^b$ as presented in Line et al. (2015) and has been used in various previous setups (e.g., Barrado et al. 2023; Matthews et al. 2025). We individually retrieved a value of b for each of the three channels in the MIRI dataset (1, 2, and 3) and the Gemini data values separately, adding four more parameters to the retrieval. To test whether this treatment is needed, we performed retrievals without this inflation and could see a significant decrease in the posterior width. When neglecting the error inflation, the retrievals become overconfident due to the small error bars. Therefore, by incorporating an error inflation we are able to mitigate this effect to some extent.

The chosen priors for the presented retrievals are shown in Table 1. For all retrievals, we chose priors derived from evolutionary models on the surface gravity and the radius based on Table 4 in (Zhang et al. 2025b), corresponding to a 3σ range of retrieved values from an extensive grid-model comparison.

In a second step, we used the retrieved outputs of the free and constrained retrieval to enable running the retrievals on the full resolution of the datasets applying line-by-line (lbl) opacities.

This allowed us to test for trace molecules and their isotopologs using the full resolving power of MIRI/MRS while still including the Gemini dataset. As running retrievals on full resolution is much more computationally expensive compared to running them on a binned spectrum, we needed to make some simplifications. As we want to probe for trace molecules, we are interested in the chemical composition. Thus, we only retrieved for molecular abundances and fix the PT structure, the radius, and surface gravity values, according to the outputs of the lower resolution retrievals. In a statistical sense, we are using the data twice, and thus we expect to reach more confident results in the second compared to the first set of retrievals. As the MIRI MRS instrument resolution varies, we calculated the forward model at a resolution of λ/Δλ = 15 000 (about 5 times the instrument resolution), which then got convolved to the mean instrument resolution per channel (Argyriou et al. 2023) and binned to the same spacing as the data.

We performed a leave-one-out (LOO) analysis on the latter, running one retrieval including all molecules we tested for, the so-called “base” retrieval, and eight retrievals where we excluded each one molecule. This enabled us to do a Bayes factor comparison and thus determine whether the corresponding molecule is significantly improving the fit or not. We ran two base retrievals for both the constrained and the free PT parameterization and we again retrieved for the error inflation parameter to fine-tune their values using the full resolution. For all the retrievals presented in the LOO analysis we fixed the error inflation values to the ones retrieved in the base retrieval. Tests showed that the values remained within one sigma for the b values in the retrievals with one molecule removed compared to the base, and thus fixing the error inflation parameter saved valuable computation time.

For the retrievals we included the following opacities in c-k resolution (R~1000) on the binned spectrum: H₂O (Rothman et al. 2010), CH₄ (Hargreaves et al. 2020), CO (Rothman et al. 2010), CO₂ (Yurchenko et al. 2020), NH₃ (Coles et al. 2019), H₂S (Azzam et al. 2016), PH₃ (Sousa-Silva et al. 2015), K (Mollière et al. 2019; Allard et al. 2016), and Na (Allard et al. 2016). We included as well collisional induced absorption line lists accounting for H₂-H₂ (Borysow et al. 2001; Borysow 2002), and H₂-He (Borysow et al. 1988, 1989; Borysow & Frommhold 1989) collisions. We neglected scattering in the presented retrievals as we deal with a clear atmosphere as we discuss in Section 4.2 and 5.8. We compare the retrieved abundances to a precomputed grid of chemical equilibrium calculations implemented in petitRADTRANS, easyChem (Mollière et al. 2017; Lei & Mollière 2025) in Section 4.3.

For the full resolution retrievals, we used the lbl line lists of R~10⁶, convolved and binned down to the corresponding resolving power of the instrument. We used the following opacities: H₂O (Polyansky et al. 2018), CH₄ (Hargreaves et al. 2020), NH₃ (Rothman et al. 2013), H₂S (Azzam et al. 2016), PH₃ (Sousa-Silva et al. 2015), K (Mollière et al. 2019; Allard et al. 2016), ¹⁵NH₃ (Rothman et al. 2013), H₂¹⁷O (Rothman et al. 2013), H₂¹⁸O (Rothman et al. 2013), ¹³CH₄ (Gordon et al. 2022), CH₃D (Rothman et al. 2013), HDO (Rothman et al. 2013), C₂H₂ (Rothman et al. 2013), and HCN (Harris et al. 2006). We included the same collisional-induced absorption lines and Rayleigh scattering as for the lower resolution retrievals. We retrieved for the abundances of H₂O, H₂¹⁸O, H₂¹⁷O, HDO, CH₄, ¹³CH₄, CH₃D, CO₂, CO, NH₃, ¹⁵NH₃, H₂S, PH₃, C₂H₂, HCN, and K. We chose the list of isotopologs based on the three most abundant trace gases: H₂O, CH₄ and NH₃. In addition, we tested for C₂H₂ and HCN as they have been detected unexpectedly in a T dwarf of similar temperature (Matthews et al. 2025). We neglected Na in the full resolution retrievals as it was unconstrained for the lower resolution retrieval and only has very low opacity in the wavelength range we probe with MIRI/MRS as shown in Fig. B.1 in the appendix.

To sample the prior space we used the nested sampling algorithm MultiNest (Skilling 2004; Feroz & Hobson 2008) with the corresponding python implementation given by pyMultinest (Buchner et al. 2014). We ran all retrievals using 1000 live points with a sampling efficiency of 0.05. All full resolution retrievals were run in nonconstant efficiency mode to probe the evidence precisely enabling a valid Bayes factor comparison, whereas the lower resolution retrievals were run in constant efficiency mode, as we were here not interested in the exact values of the evidence.

To decide whether a molecule is detected or not we use the Bayes factor (e.g., Benneke & Seager 2013; Thorngren et al. 2026). We applied the method as presented in Konrad et al. (2024) computed depending on the evidences of the two best-fits of the base Z_b and the one excluding the molecule in question ${\mathcal{Z}}_m:\mathrm{l}\mathrm{o}\mathrm{g}(\mathcal{B})=\frac{\mathrm{l}\mathrm{n}({\mathcal{Z}}_{\mathrm{b}})-\mathrm{l}\mathrm{n}({\mathcal{Z}}_{\mathrm{m}})}{\mathrm{l}\mathrm{n}(10)}$. Large values in log(B) correspond to the molecule significantly improving the fit and thus being detected. Values for log(B) below 0.5 show no significant improvement and negative values even favor the model that excludes the molecule (Thorngren et al. 2026). In addition, we calculated the BPICS value (Ando 2011; Thorngren et al. 2026) as a different metric for the LOO analysis. Here, smaller values indicate a preference of the model. The BPICS is a simplified version of the Bayesian Predictive Information Criterion (BPIC, Ando 2007) and based on Thorngren et al. (2026) computed like: BPICS = −2E_p[ln(L)] + 2N_p, where N_p is the number of parameter used for the model (ranging between 16 and 15 for the base and the LOO retrieval respectively) and E_p[ln(L)] the median value of the natural logarithm of the likelihood L (assuming a Gaussian distribution of the latter).

4 Results

In this section, we present the results of our retrieval analysis on the combined MIRI/MRS and Gemini/FLAMINGOS-2 data. We start by discussing the bulk parameters and PT structure, then the chemical composition and how it compares to the chemical equilibrium. Further, we discuss the full resolution retrieval, with emphasis on isotopologs.

4.1 Bulk parameters

In Figure 2, we show the best-fit spectra of the retrievals using the free and the constrained PT profiles compared to the MIRI/MRS data in black and the Gemini/FLAMINGOS-2 data in gray. Panel a) shows the spectrum from 1 to 18 μm (with gaps between 1.12-1.18 μm, 1.35-1.45 μm, 1.66-1.95 μm, and 2.34.95 μm to exclude low S/N areas in the spectrum), panel b) a subset on the Gemini data, and c) a subset on the ammonia feature in the MIRI/MRS dataset. We see both fits explain the general spectrum well. The free retrieval finds a better agreement of the 2.2 μm peak, which however is a part of the spectrum with larger uncertainty on the data. The constrained retrieval fits the depth of the ammonia feature in the region between 8 and 12 μm slightly better. The largest differences between the two best-fit models are between 2.5 and 5 μm due to the lack of data. Especially, the retrieved abundance of CO is slightly different as the CO feature is deeper for the constrained case compared to the free case. Larger differences can be seen at 56 μm, probably originating from the CO abundance, not being matched precisely. The rest of the spectrum is explained well by both models and the residuals are on the order of 1 σ from the data, when including the error inflation.

The obtained bulk parameters for the free and constrained retrievals are presented in Fig. 3 and in Table 2. We obtain an effective temperature estimate using the built-in function in petitRADTRANS (compute_effective_temperature) integrating model spectra from 1 to 18 μm, which were sampled from the posterior. The effective temperature estimates fall well inside the one sigma range presented from the grid models (Zhang et al. 2025b) shown in black dashed and dotted lines (mean and error, respectively). The free retrieval results in an effective temperature of 471.3${}_{-1.4}^{+1.6}$ K compared to the value of 483${}_{-53}^{+44}$ K presented in Zhang et al. (2025b). This is compatible within three sigma with the value returned from the constrained retrieval of 468.1${}_{-1.1}^{+1.1}$ K. Kiman et al. (2026) present a slightly warmer estimate with a value of 493${}_{-9}^{+9}$ K. In general, the values found in this work and from Kiman et al. (2026) fall in the estimated range by Zhang et al. (2025b).

The metallicity is calculated based on the retrieved abundances using the built-in function volume_mixing_ratios2metallicity in petitRADTRANS. We find a solar to subsolar metallicity of 0.01${}_{-0.02}^{+0.02}$ dex and −0.12${}_{-0.02}^{+0.01}$ dex in the free and constrained case, respectively, which is in line with the findings of Zhang et al. (2025b) finding a solar to subsolar metallicity. Their best-fitting model gave a metallicity of −0.152${}_{-0.010}^{+0.012}$ dex. The difference in the metallicity is likely connected to changes in the PT structure, which we discuss in Section 4.2. The free retrieval is consistent with the host stars solar metallicity of 0.00±0.08 Hojjatpanah et al. (2019).

The C/O ratio is calculated from all carbon bearing species (namely CH₄, CO₂, and CO, based on the estimates of the low resolution retrievals) divided by the oxygen bearing molecules (H₂O, CO, and CO₂) accounting for the multiplicity of the respective atoms. For colder brown dwarfs, Calamari et al. (2024) determined an amount of 17.8${}_{-2.3}^{+1.7}$% of the oxygen budget being locked up in low lying condensing clouds using stoichiometric and mass balance calculations. Thus, we multiply our obtained C/O ratio by a factor of 1.22 to account for this oxygen sink according to the procedure presented in Rowland et al. (2024). The free and constrained retrievals constrain a slightly subsolar C/O ratio with a value of 0.43${}_{-0.01}^{+0.02}$ and 0.45${}_{-0.02}^{+0.02}$, respectively. Zhang et al. (2025b) finds a tendency to favor a subsolar C/O ratio, which is in line with our findings.

We retrieved the surface gravity and the radius and calculated a mass from these parameters. The radius in the free retrieval is constrained as 1.133${}_{-0.006}^{+0.005}$ R and in the constrained case as 1.160${}_{-0.003}^{+0.002}$ R_J. Both radii are slightly enlarged compared to the value of 1.07${}_{-0.04}^{+0.03}$ R_J found by Zhang et al. (2025b). In the retrieval, we added an evolutionary prior to the radius from 1.04 to 1.19 R_J. The logarithmic surface gravity of the free retrieval was estimated as 3.90${}_{-0.04}^{+0.04}$ dex and in the constrained case as 3.88${}_{-0.05}^{+0.04}$ dex. Again, we applied prior boundaries for the surface gravity between 3.7 and 4.4 dex based on the three sigma adopted parameter in Zhang et al. (2025b). Both of the gravity values are slightly smaller compared to previous values (4.19${}_{-0.13}^{+0.18}$ in Zhang et al. (2025b) and 4.17±0.02 in Kiman et al. (2026). The mass estimate was calculated from the radius and gravity values, and as our gravity estimates are smaller, we receive a smaller estimate of 4.1${}_{-0.3}^{+0.4}$ and 4.1${}_{-0.5}^{+0.4}$ M for the free and the constrained retrieval, respectively, compared to the 8 ±2 M_J (Zhang et al. 2025b) and 7.5±0.4 M_J Kiman et al. (2026). However, this mass estimate was obtained purely through the radiative transfer model and thus may be influenced due to degeneracies with other parameters. Specifically, it has been shown that surface gravity is difficult to estimate from atmospheric models alone (Mader et al. 2026). Thus, measuring the mass through evolutionary models that include a better age and luminosity estimate, similar to what was done in Kiman et al. (2026), would result in a more realistic mass measurement. To assess the degeneracies in the retrieved bulk parameters, we compare them in the corner plot in Fig. B.9 in the appendix. In general, we find a solar to subsolar metallic-ity, a subsolar C/O ratio compatible with Zhang et al. (2025b) and an inflated radius, a smaller gravity and a smaller mass in comparison to Zhang et al. (2025b) and Kiman et al. (2026).

Fig. 2 Best-fit spectra of the free and constrained atmospheric retrievals, respectively in green and violet, compared to the data in black. Panel a shows the full observation. In panel b, we show a subset at the Gemini wavelengths from 1 to 2.5 μm as indicated by the black box and in panel c we show another subset for the strong NH3 feature between 8 and 12 μm. The spectra are shown in lower resolution $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$.

Fig. 3 Bulk parameter for the free and the constrained retrievals. We show in panel a, the effective temperature estimate, in panel b the retrieved metallicity, in panel c the C/O ratio accounting for oxygen sequestration, in panel d the radius, in panel e the surface gravity, and in panel f the resulting mass based on the radius and gravity for the free and constrained retrieval (in green and violet, respectively). For panels a, b, and d-f, we show the mean values with dashed lines and the one sigma estimate with dotted lines, for the values obtained by Zhang et al. (2025b, Z25) in black and Kiman et al. (2026, K25) in blue.

Fig. 4 PT profiles of the free and the constrained retrieval (in green and violet, respectively). In dark blue we show the PT structure of Jupiter. The H2O, NH4SH, and NH3 condensation lines were taken from (Lodders & Fegley 2002), and the rest is from petitRADTRANS. Thicker lines in deeper atmospheric layers indicate the convective part of the atmosphere calculated based on Eq. (1).

4.2 Atmospheric structure

In Fig. 4, we show the retrieved PT structure for the free and constrained retrievals in green and violet respectively. The free retrieval provides a not strictly monotonic PT structure, with areas where the temperature gradient is not constant at around 50 bar, 0.5 bar, and 0.01 bar. As we have not used any form of regularization here, this might be due to the spline interpolation between the nodes resulting in such non-monotonies. Further, retrievals tend to overfit the data and such small changes in the PT structure might locally explain the data better. By comparing the free to the constrained PT profile, we see the strong effect of the constraint on the profile. The constrained retrieval results in a monotonically increasing temperature with pressure, as we would expect it for gas giants without any external heat sources. The overall shape is compatible with the free retrieval. Both retrievals become nearly isothermal above about 0.01 bar.

In this area, however, the contribution goes close to zero, and thus the changes seem not to affect the spectrum. In dark blue we add the measured PT profile of Jupiter (Seiff et al. 1998) as a reference for an irradiated gas giant. As expected, the profiles are shifted towards higher temperatures by about 420K. Also, the slope of the retrieved temperatures increase faster with pressure in comparison with the latter. We also plot the condensation lines of the most common species, namely NH₃, NH₄SH, H₂O, KCl, Na₂S, Fe, and silicates. The mean contribution function across all wavelengths showing where in the atmosphere the flux is emitted from is presented in dashed lines for both cases. Most of the flux is emitted from pressures between 10 and 0.1 bar, where no condensation lines are crossed by the PT profiles. Thus, we do not expect to see strong features from condensing clouds in the atmosphere, and indeed we can explain the spectrum well with a clear atmosphere. We discuss more about clouds in Section 5.8. The emission contribution of the free retrieval is shifted slightly to lower pressures compared to the constrained case. At lower pressures spectral features are less pressure broadened. Thus, to obtain the same spectrum in both cases, the abundances of the trace gases need to be increased (Mollière et al. 2015). Therefore, we find a higher overall metallicity for the free compared to the constrained case as presented in Fig. 3.

We approximated the boundary between the convective and the radiative part using the Schwartzschild criterion for convective stability, assuming the ideal gas and hydrostatic equilibrium, and we found the following expression for the convective part: $$\frac{\mathrm{\partial }\mathrm{T}}{\mathrm{\partial }\mathrm{P}}<\frac{\mathrm{R}\mathrm{T}}{{\mathrm{c}}_{\mathrm{p}}\mu \mathrm{P}},$$(1)

with T as the temperature, P as the pressure, R as the molecular gas constant (R=8.314 J/mol K), μ as the mean molecular weight (here μ ~2.33), and c_p as the specific heat capacity of a H and He dominated gas (here c_p ~27.5 J/mol K). Thus, we estimated the transition between the convective and radiative part to be at around 2-3 bars. We followed the approach presented in Heng (2017).

4.3 Chemical composition

The strongest absorption features of cold atmospheres below 500 K are due to H₂O, CH₄, and NH₃. In our low resolution retrieval analysis we find the following composition presented in Fig. 5. In green we show the values for the free retrieval and in violet the constrained values. We report all abundances in logarithmic volume mixing ratios. The exact numbers are given in Table B.1 in the appendix and we compare the values to the ones from the full resolution retrieval in Fig. B.8. We compare the retrieved abundances to the chemical equilibrium predictions at 0.28 bar, which corresponds to the height where the maximum emission emerges from in the case of the free PT parameterization. The chemical equilibrium is evaluated based on the retrieved PT structure and the corresponding retrieved metallicity and C/O ratio. In the panel b of Fig. 5, we show how the chemical equilibrium changes with pressure for CO, CO₂, CH₄, and H₂O for the free and the constrained case in black and gray respectively. When comparing the retrieved abundances to the equilibrium prediction at the height of maximum contribution (in the free retrieval), we see that they agree well for H₂O and CH₄. For CO they agree when taking into account quenching at the highlighted quench pressure in magenta. Further, this implies that choosing vertically constant abundances is a valid approximation for the pressure range covered by the emission contribution.

CH₄: for CH₄ we find a logarithmic abundance of −3.28${}_{-0.03}^{+0.03}$ and −3.42${}_{-0.02}^{+0.02}$ compared to a slightly lower equilibrium chemistry abundance of −3.33 and −3.46, for the free and constrained case, respectively. The full resolution retrievals obtain a value of −3.247${}_{-0.009}^{+0.01}$ and −3.353${}_{-0.009}^{+0.009}$ for the free and constrained retrieval respectively, in the same order of magnitude as the low resolution retrieval fits.

H₂O: we find slightly larger values compared to the equilibrium for H₂O. Here, we constrain a logarithmic volume mixing ratio of −2.96${}_{-0.03}^{+0.03}$ and −3.09${}_{-0.02}^{+0.02}$ compared to the chemical equilibrium of −3.04 and −3.18 for the free and constrained case respectively. For the full resolution free and constrained retrievals we find −2.713${}_{-0.007}^{+0.006}$ and −2.884${}_{-0.006}^{+0.006}$ respectively, compatible with the lower resolution results.

NH₃: larger discrepancies are found for NH₃, where we find values of −4.91${}_{-0.03}^{+0.03}$ and −5.09${}_{-0.02}^{+0.02}$ in the free and constrained case, which are significantly lower compared to the chemical equilibrium −3.94 and −4.06. The full resolution retrievals constrain values for the free and constrained retrievals of −4.925${}_{-0.005}^{+0.005}$ and −5.064${}_{-0.005}^{+0.005}$ respectively, so similar to the lower resolution fits.

CO and CO₂ : further, we find constraints on CO and CO₂, which are −3.95${}_{-0.05}^{+0.05}$ and −3.94${}_{-0.05}^{+0.04}$ for CO and −9.5${}_{-1.0}^{+0.9}$ and 9.6${}_{-0.9}^{+0.9}$ for CO₂ in the free and constrained case respectively. For CO the full resolution retrievals obtain values for the free and constrained cases −4.09${}_{-0.03}^{+0.03}$ and −4.06${}_{-0.03}^{+0.03}$ and for CO₂, −8.4${}_{-0.6}^{+0.3}$ and −8.2${}_{-0.3}^{+0.2}$ respectively. In chemical equilibrium, abundances of CO and CO₂ are negligible at the observed pressure, as temperatures are too low to chemically produce CO and CO₂ and rather CH₄ is favored. The fact that we constrain these two molecules with values much larger compared to the chemical equilibrium hints that other processes, such as vertical mixing must be present to create the observed disequilibrium chemistry. The implications of this is further discussed in Section 5.3. From the leave-one out comparison we find that CO₂ is preferred based on a Bayes factor comparison by 1.07. However, due to the only small effect of the molecule on the spectrum as shown in Fig. B.3 in the appendix, we conclude that this is a tentative detection. However, another larger band is visible in the 4.2 μm region visible in the NIRSpec data (Kiman et al. 2026), where it could indeed be detected more easily.

From the retrieved values we can calculate the amount of vertical mixing in the atmosphere parametrized by the K_zz value. We use the formulas for the chemical timescale where CO and CH₄ are in chemical equilibrium and mixing timescales t_mix = H²/K_zz with H the scale height based on Zahnle & Marley (2014) and rearranged as presented in Nasedkin et al. (2024a). We obtain a logarithmic K_zz of 3.3 cm²/s for both cases at a quench pressure of 43 and 53 bar for the free and the constrained case, respectively.

K and Na: we constrain the alkali species K, but not Na. This might have to do with the very small wavelength range where Na would be visible, which is close to 1 μm and thus at the edge of the Gemini dataset and in a noisy area of the spectrum (see Fig. B.1). Both species would not be expected in chemical equilibrium with logarithmic volume mixing ratios much smaller than −10 as they would only be stable in hotter atmospheres. We do not include any rainout of the alkali species and assume constant profiles throughout the atmosphere. As the absorption features of K are only visible between 1-2 μm (see Fig. B.1 in the appendix), this only influences the fit of the Gemini dataset. For the free and the constrained case we find an abundance of K of −10.0${}_{-0.7}^{+0.7}$ and −8.0${}_{-0.1}^{+0.1}$. In the full resolution case we constrain values for K of −10.2${}_{-0.7}^{+0.8}$ and −8.1${}_{-0.1}^{+0.1}$ for the free and constrained case respectively.

H₂S: for H₂S we find a higher logarithmic volume mixing ratio of −3.84${}_{-0.05}^{+0.05}$ in the free case and a smaller value of −4.25${}_{-0.08}^{+0.08}$ in the constrained case. In the free case the chemical equilibrium prediction is −4.66, and for the constrained case −4.79, which is smaller compared to our retrieved values. This is in line with what we find from the full resolution retrievals, with values of −3.83${}_{-0.03}^{+0.03}$ and −4.19${}_{-0.04}^{+0.04}$ for the free and constrained case, respectively. Including the molecule H₂S results in a clear preference according to the logarithmic (base 10) Bayes factor of 61.46 and a clear influence on the flux at around 1.6 μm is visible as shown in Fig. B.3 in the appendix. The abundance of H₂S is slightly larger compared to the chemical equilibrium prediction and varies dependent on the PT choice. As some of the lines that inform the fit are located in the FLAMINGOS-2 dataset, which is noisier, this could explain the differences. Also in the NIR, potentially unaccounted for clouds could vary the retrieved abundance when included in the fit. We discuss more on clouds in Sect. 5.8.

PH₃: Further, we constrain PH₃ in the free case with a value of −7.3${}_{-0.3}^{+0.2}$ and in the constrained case a value of and −7.3${}_{-0.2}^{+0.1}$. The full resolution retrievals find −7.10${}_{-0.07}^{+0.06}$ and −7.47${}_{-0.10}^{+0.09}$ for the free and constrained case, being compatible with the lower resolution results. In the chemical equilibrium we would expect a value of −6.46 and −7.36 for the free and constrained case. In equilibrium it drops significantly below pressures of about ~0.3 bar. We find that the full resolution retrieval including PH₃ is favored with a logarithmic Bayes factor of 7.03. The value we retrieve is of the order of magnitude we would expect from chemical equilibrium. However, we do not see a clear absorption feature at 10.1 μm as shown in Fig. B.3 in the appendix. At this wavelength we would expect to see the largest difference between the model including and the one excluding PH₃. Thus, we conclude we only tentatively detect this molecule. The opacity of PH₃ might contribute to the continuum and could potentially be influenced by other unaccounted effects, such as clouds. As the presence of this molecule in brown dwarfs is debated, we further discuss this in Section 5.3.

C₂H₂ and HCN: As presented in Table 3, we obtain a Bayes factor of 0.19 for C₂H₂ improving the fit very slightly compared to the base. However, as presented in Fig. B.2, we do not detect clear absorption features. Therefore, we do not detect this molecule. The logarithmic volume mixing ratios we find show a larger tail towards lower values with logarithmic abundances of −9.1${}_{-0.6}^{+0.3}$ and −8.9${}_{-0.5}^{+0.2}$ for the free and constrained case. HCN was not constrained in the full resolution retrievals with resulting values of −10.1${}_{-0.7}^{+0.8}$ and −10.1${}_{-0.7}^{+0.8}$ for the free and constrained case, respectively. Using the Bayes factor comparison in the LOO analysis, we find a value of −0.72 and thus conclude that HCN does not significantly improve the fit when included, and we can state a non-detection of HCN. This is confirmed by visual inspection of the fits presented in Fig. B.4 in the appendix. We discuss the absence of the hydrocarbons C₂H₂ and HCN in Section 5.7.

Fig. 5 Retrieved chemical abundances with the free and constrained retrievals in the spectral resolution $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$ (in green and violet, respectively) compared to the chemical equilibrium predictions for the free and the constrained retrievals (in dark and light gray, respectively). In panel a the abundances were taken at 0.28 bar, corresponding to the maximum contribution in the free case. In panel b the variations with pressure of the abundances involved in the chemical equilibrium, CO; CH4 , H2O, and CO2 are presented. We additionally indicate the pressure of the maximum contribution (0.28 bar) in blue and the quench pressure (43-53 bar) in magenta.

4.4 Isotopologs

To search for isotopologs in the MIRI/MRS spectrum, we perform retrievals on the full resolution of up to R~3750. We test for isotopologs of the most absorbing species in the infrared (ammonia, methane, and water). We included H₂¹⁸O, H₂¹⁷O, ¹⁵NH₃, ¹³CH₄, HDO, and CH₃D and compare the retrievals to the base fit in a LOO analysis.

In Figs. 6, 7, and 8, as well as the corresponding ones in the appendix we show the best-fit spectra for the base retrieval with the constrained PT parameterization without the molecule in question in red and including it in blue compared to the data in black. We also plot the difference between the data and the respective best-fits in the middle plots as well as the difference between the two best-fits in the lowest plots. For the first time in a cold brown dwarf, using MIRI/MRS we observe absorption features of the water isotopolog H₂¹⁸O and H₂¹⁷O in the 5-7 μm area, as can be seen in Figs. 6 and 7. In addition, including ¹⁵NH₃ helps explain several features in the 9-10 μm wavelength range, as presented in Fig. 8. More wavelength areas are covered in the appendix.

We detect three isotopes confidently with log(B) > 20 and a fourth isotope tentatively. In Table 4 we show the logarithmic evidence of the full resolution retrievals, comparing the base retrieval with the constrained PT profile to the LOO retrievals. H₂¹⁸O is favored with a difference of 72.486 compared to the base retrieval, which is the highest value in this comparison. Smaller differences are shown for H₂¹⁷O with a difference of 25.81 and for ¹⁵NH₃ with a difference of 24.68. Thus, these three isotopologs seem to be significantly improving the fit. In comparison, with a Bayes factor of 4.3 the presence of ¹³CH₄ is only slightly favored in the Bayes analysis. However, as shown in Fig. B.2, we do not detect clear features of ¹³CH₄ and therefore declare this detection as being tentative. HDO and CH₃D are not preferred in the Bayes comparison with negative logarithmic Bayes factors. Also, as shown in Figs. B.4 and B.8 in the appendix, we do not constrain HDO or CH₃D, but receive an upper limit, which is further discussed in Section 5.6. As mentioned previously, since we are statistically using the data twice by fixing the PT profile in a second retrieval, the detection significances may be overconfident. Thus, we require to distinguish the absorption features in the spectrum for the detected molecules, which we do for the oxygen and nitrogen isotopologs. As highlighted by Kipping & Benneke (2025, community commentary) and Thorngren et al. (2026), we refrain from discussing the σ values in this case.

In Fig. 9, we show the isotope ratios ¹⁴N/¹⁵N, ¹⁶O/¹⁷O, ¹⁶O/¹⁸O, and ¹⁸O/¹⁷O derived from the full resolution base retrievals including all molecules using the constrained PT profile. We compare the obtained values to literature values for the Sun and the ISM. For nitrogen we additionally compare the retrieved value to the values found from WISE 1828 (Barrado et al. 2023) and WISE 0855 (Kühnle et al. 2025) and for the oxygen molecules to the values found in VHS1256b (Gandhi et al. 2023) (from C¹⁸O). In addition, Ruffio et al. (2026) measure C¹⁸O in three of the four HR8799 gas giants. For ¹⁶O/¹⁸O they find larger ratios compared to our measurements with significantly larger error bars. The retrieved values are summarized in Table 5.

We obtain for ¹⁴N/¹⁵N values of 337${}_{-38}^{+47}$ and 324${}_{-40}^{+46}$ in the free and the constrained case, respectively. This is compatible within two sigma with the value found for the ISM of 274${}_{-18}^{+18}$ (Ritchey et al. 2015) and lower compared to the Solar value of 458.7${}_{-4.2}^{+4.2}$ (Marty et al. 2011). We can compare this value also to two Y-dwarfs where these molecules have been found: In WISE1828 a value of 670${}_{-211}^{+390}$ (Barrado et al. 2023) and 349${}_{-41}^{+53}$ (Kühnle et al. 2025) have been constrained (for a clear retrieval). Thus, the values found for COCONUTS-2b are consistent within one and two sigma with the values found in WISE0855 and WISE1828, respectively.

For the ratio ¹⁶O/¹⁸O we find a value of 247${}_{-24}^{+27}$ and 256${}_{-25}^{+29}$ for the free and constrained retrieval. These values seem significantly enriched in the heavier isotoplog when comparing to the ISM value of 557${}_{-30}^{+30}$ (Wilson 1999) and the solar value of 529.7${}_{-1.7}^{+1.7}$ (McKeegan et al. 2011). This ratio was also constrained for the planetary mass companion VHS1256b, with a value of 425${}_{-28}^{+33}$ which is significantly larger compared to our findings for COCONUTS-2b. We constrain the ratio of ¹⁶O/¹⁷O as 662${}_{-83}^{+98}$ and 934${}_{-139}^{+174}$ for the free and constrained retrieval compared to the ISM value of 2005.2${}_{-30}^{+30}$ (Wilson 1999) and the Solar value of 2776${}_{-14}^{+14}$ (McKeegan et al. 2011). The value presented in Gandhi et al. (2023) for VHS1256b is compatible with our free case with a value of 1010${}_{-100}^{+120}$. We find the ratio between the two heavier isotopes ¹⁸O/¹⁷O to be 2.7${}_{-0.4}^{+0.6}$ and 3.7${}_{-0.7}^{+0.8}$ They are within one to two sigma compatible with the values found for the ISM 3.6${}_{-0.2}^{+0.2}$ (Wilson 1999). The values are smaller compared to a solar value of 5.5 (Wilson 1999). Compared to VHS1256b's value of 2.4${}_{-0.3}^{+0.3}$ (Gandhi et al. 2023), we find compatible values within one to two sigma for the constrained and free retrieval. In summary, we find ¹⁴N/¹⁵N to be compatible with previous observations and the ISM value, ¹⁶O/¹⁷O and ¹⁶O/¹⁸O enriched compared to the solar and ISM value and ¹⁸O/¹⁷O in line with previous observations and the ISM value.

Fig. 6 Best-fit retrievals on full MIRI/MRS resolution data with and without including H218O (in blue and red, respectively) compared to the data in black. In panels a-c, one can see the three different parts in the spectrum where the difference in the best-fits is visible. The residuals between the data and the corresponding best-fits are shown in panels d-f and the differences in the fits in plots g-i. Here we present the wavelength ranges with the fourth, fifth, and sixth largest differences in the models compared to showing the ranges with the three largest difference, as there the models do not explain the data well overall. The latter and more wavelength areas are presented in the appendix in Fig. B.6.

Fig. 7 Best-fit retrievals on full MIRI/MRS resolution data with and without including H217O (in blue and red, respectively) compared to the data in black with the same panel structure as in Fig. 6. We present the wavelength ranges with the largest, second largest, and sixth largest differences in the models. We present more wavelength areas in the appendix in Fig. B.7.

Fig. 8 Best-fit retrievals on full MIRI/MRS resolution data with and without including 15NH3 (in blue and red, respectively) compared to the data in black with the same panel structure as in Fig. 6. Here, we present the wavelength ranges with the largest, second largest, and third largest differences in the models. More wavelength areas are provided in the appendix in Fig. B.5.

Fig. 9 Ratios of a) 14N/15N, b) 16O/18O, c) 16O/17O, and d) 18O/17O based on the full resolution retrievals. We compare the ratios to the values for the ISM (oxygen: Wilson (1999), 14N/15N: Ritchey et al. (2015) and solar values (16O/17O and 16O/18O: McKeegan et al. (2011), 18O/17O: Wilson (1999) and nitrogen: Marty et al. (2011) as well to the values of WISE 1828 (Barrado et al. 2023), WISE 0855 (Kühnle et al. 2025) for 14N/15N, and VHS1256b for the oxygen ratios (Gandhi et al. 2023).

5 Discussion

5.1 Full resolution retrievals: Advantages and shortcomings

We present a full resolution retrieval analysis and derive several elemental abundance ratios from it. The reported ratios, however, might be influenced by our assumption to fix the PT profile in the setup. Substantially reducing the number of parameters allows the retrievals to converge in a reasonable amount of time (on the order of days to weeks) and thus is a crucial step to enable this analysis. Generally, the PT structure is determined by the overall flux distribution which is well captured by the lower resolution retrieval. Still, the abundances of the isotopologs might be influenced by small unaccounted changes in the PT structure. To assess the variation in the model we run the same setup with the fixed value for the PT profile from the free and constrained retrievals and compare the derived values. As presented in Fig. 9 it has an effect on the ratio we retrieve for ¹⁶O/¹⁷O. However, the ratios of ¹⁴N/¹⁵N and ¹⁶O/¹⁸O are compatible with each other within one sigma. To investigate how the observed absorption features could be influenced due to changes in the PT structure, we calculated the models with and without the isotopologs for a hotter (positive shift of the interior temperature by 500 K), a colder (negative shift by 500 K), a steeper (increase in 10% per node), and a shallower (decrease in 10% per node) PT profile based on the one retrieved in the constrained retrieval. In Fig. B.10 in the appendix we show how the flux varies with respect to the changes. We find that indeed there are differences in the models based on the PT structure, however affecting only the overall shape of the spectrum and not the absorption features of the isotopologs ¹⁵NH₃, H₂¹⁸O, and H₂¹⁷O, as they stay visible for each test case. Thus, we conclude that the PT structure alone could not result in the detected spectral features and thus the assumption of fixing the PT structure is a valid approach. In future studies, more computational power is needed to enable retrieving the PT structure simultaneously with the abundances, which is currently beyond the scope of this setup.

Further, we fixed the radius and the surface gravity to the values of the lower resolution retrieval. We did this to only retrieve for the abundances simultaneously, as the goal of the LOO analysis was to test for the importance of each molecule. A change in radius only leads to a difference in the absolute flux and thus is not sensitive to the higher spectral resolution.

The amount of change in the abundances between the low and the full spectral resolution retrievals is assessed in Fig. B.8 in the appendix. Larger differences are seen for H₂O, CO, and CO₂. We retrieve more H₂O, less CO and a more contained upper estimate for CO₂ when going from low to the full resolution. It is visible that the full resolution retrievals result in narrower posteriors compared to the low resolution results. On the one hand, by fixing several parameter in the full resolution retrievals to values based on the low resolution retrievals, we statistically use the data twice in the setup. Thus, we expect narrower posteriors in the full resolution retrievals compared to the low resolution ones. On the other hand, providing more information from the data in the full resolution case, we receive a more precise estimate of the abundance.

In this analysis, we demonstrate that full resolution retrievals are possible using the medium resolution of the MIRI/MRS. Similarly, Hood et al. (2024, subm. to Nature Astronomy) and Ruffio et al. (2026) performed medium resolution retrievals on JWST/NIRSpec data. The advantage of running retrievals in full resolution is clearly to detect and measure the abundances of molecules and isotopologs with faint signatures, such as the ones here presented. In Fig. 10, we compare the two best-fits using the constrained PT profile in the area of 11.0-11.1 μm demonstrating the amount of detail gained in the full resolution fits. For this dataset, the small spectral signature of ¹⁵NH₃ located at 11.075 μm would have not been detectable using the lower resolution retrievals. Analyzing the full resolution allows us to study new molecular abundances in current and future mediumresolution observations with JWST. Especially, targets as cold as COCONUTS-2b with high S/N spectra in the MIR will profit from the higher resolution retrievals on MIRI/MRS data. This will help us study compositional differences between host stars and wide separation companions and potentially lead to observational constraints on the formation scenarios of the latter by measuring various elemental abundance ratios in addition to C/O (Öberg et al. 2011).

Fig. 10 Comparison of low and full resolution retrievals compared to the data. Panel a shows the best-fit spectra for the low resolution in green and the full resolution retrievals with and without 15NH3 (in blue and red, respectively). In panel b we show the corresponding residuals. Only with the improvement of the resolution were we able to detect the isotopolog.

5.2 Small error estimates

Previous versions of the JWST pipeline underestimated the uncertainties for MIRI/MRS data (e.g., Barrado et al. 2023; Kühnle et al. 2025) and the pipeline provided a too small estimate of the error bars. This has now been improved, and the pipeline version used here (v1.18.0) now propagates the errors properly. Small error bars on the data result in difficulties for the retrieval to fit the data as it tends to overfit minor differences. Therefore, the posterior distributions become very confident in the retrieved values. Even though the error bars from the pipeline have increased, we retrieve for an error inflation factor 10^b (Line et al. 2015) to account for model uncertainties in the 1D radiative transfer model and a potentially insufficient increase in the errors by the pipeline.

We retrieve a smaller surface gravity and a slightly increased radius compared to the estimate of Zhang et al. (2025b) and Kiman et al. (2026) based on self-consistent grid-models, while applying evolutionary model bounds to the priors. This then results in a significantly smaller mass than the previous estimates. This should highlight that obtained mass uncertainties in retrievals are not trustworthy. Again, the small posterior widths for the radius and the gravity lead to this effect.

The obtained abundances might well be overconfident, and thus the isotope ratios we retrieve are likely overly constrained. For the oxygen ratios, we would expect slightly larger error estimates due to the lower S/N in the part of the spectrum, where the molecules are detected. However, their error bars are consistent with previous measurements (e.g., Gandhi et al. 2023).

5.3 Disequilibrium chemistry

Based on Fig. 5, we observe CO clearly out of chemical equilibrium in COCONUTS-2b with abundances orders of magnitude larger than predicted by equilibrium chemistry. Thus, additional processes must be acting to modulate the chemical regime such that CO is present in the relatively cold atmospheric layers, for example, vertical mixing (Zahnle & Marley 2014). The resulting vertical mixing strength we calculate (K_zz of ~10^3.3 cm²/s) is compatible with quenching the equilibrium chemistry abundance of CO at a pressure of ~40-50 bar. Using the Sonora Elf Owl models and the Gemini/FLAMINGOS-2 data Zhang et al. (2025b) receive a stronger vertical mixing of log(K_zz) = 8.98${}_{-0.03}^{+0.01}$. Mukherjee et al. (2022) discuss the quench K_zz dependent on the effective temperature of various targets and find for T_eff ≈ 500K a typical K_zz of ~10⁶ cm²/s. This is larger compared to what we find; however, their sample has larger surface gravities than COCONUTS-2b. Smaller gravities lead to a decrease in Kzz as presented in Miles et al. (2020). Mukherjee et al. (2022) show that reaching the here observed values for CO requires small surface gravities or large vertical mixing. Our retrieved low surface gravity of about 3.9 dex and a rather small vertical mixing of K_zz =3.3cm²/s for both retrievals is compatible with this low-gravity scenario, even though the probed parameter space by Mukherjee et al. (2022) does not extend to such low gravities.

In contrast to CO, we find a smaller NH₃ abundance compared to the chemical equilibrium in both the free and constrained retrievals. Zahnle & Marley (2014) propose NH₃ as a potential tracer for surface gravity, as its abundances is largely insensitive to vertical mixing, in contrast to species such as CO and CO₂. The smaller NH₃ abundance is in line again with the relatively low surface gravity that we have retrieved.

Some studies have shown degeneracies between metallic-ity and surface gravity as the latter impacts the PT structure (Mollière et al. 2015; Zhang et al. 2021d; de Regt et al. 2025). Abundance ratios, such as the nitrogen isotope ratios presented here, seem to be more robust against this degeneracy. Even though the overall metallicities for the free and the constrained retrieval differ, the ¹⁴N/¹⁵N and C/O ratios for example stay compatible as shown in Figs. 3 and 9. Future analysis on the NIRSpec/G395H dataset (Copeland et al. in prep.) will provide additional insights on the abundances of CO and CO₂, since their strongest absorption lines are visible at around 4.2 and 4.4 μm and will thus enable tighter constraints on the metallicity and C/O ratio.

Further, our retrievals on the binned spectrum place constraints on the abundances of PH₃, despite the absence of clearly identifiable absorption features. Using the full resolution retrievals in Fig. B.3 in the appendix, we compare the best-fit spectra when including or removing one molecule from the line list. For PH₃, a weaker effect on the continuum improves the fit, while the strong feature expected at 10.1 μm is not visible. Although the logarithmic Bayes factor for including PH₃ (7.03) suggests statistical support and the retrieved abundance is broadly consistent with equilibrium chemistry, we consider the detection of of PH₃ to be tentative due to the lack of a distinct absorption feature. PH₃ is a molecule that has received significant attention with JWST. It has a significant abundance in Jupiter (Larson et al. 1977) (due to its about 3 times higher metallicity, Mahaffy et al. (2000), and was expected to be easily visible in brown dwarf spectra (Miles et al. 2020). However, early JWST observations did not detect this feature (e.g., Beiler et al. 2023; Matthews et al. 2025; Vasist et al. 2025), and in particular (Beiler et al. 2024b) highlighted that chemical models tend to over-predict PH₃ abundance relative to observations of T- and Y-dwarfs with temperatures below 500K. To date, PH₃ has only been observed in WISE0855, in which it was constrained to have an abundance of one part per billion (Rowland et al. 2024) and recently at a larger abundance in Wolf1130C (Burgasser et al. 2025). The retrieved statistics favoring PH₃ in COCONUTS-2b, alongside the absence of a clear feature, motivate further work on this mysterious molecule.

5.4 Ammonia isotope ratio consistent with ISM value

¹⁵NH₃ is only observable in the atmospheres of cold objects, and we are able to measure it using medium to high resolution data. As shown in Fig. 10, the small features of the isotopolog are only discernible when using the full spectral resolution. This is in line with what Matthews et al. (2025) found by comparing the features of NH₃ and ¹⁵NH₃ for different temperature regimes: the warmer the targets the less evident the features, with a transition occurring at around 500K. Thus, COCONUTS-2b is just cold enough for a successful ¹⁵NH₃ detection. So far, COCONUTS-2b is the warmest object where we can detect ¹⁵NH₃, and thus we demonstrate that we are able to uncover even such faint features using the full resolution retrievals. For the colder objects, WISE0855 and WISE1828, ¹⁵NH₃ was detected already in binned spectra of spectral resolution $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$. We find the isotope ratio of ¹⁴N/¹⁵N for COCONUTS-2b to be compatible with the ISM value within one to two sigma for the constrained and free retrievals. Given that COCONUTS-2b is a relatively young system, an isotope ratio close to the ISM value is consistent with formation from a molecular cloud with an average composition of the ISM. However, the retrieved ratio is also compatible within one to two sigma with the only two existing measurements of Y dwarfs (Barrado et al. 2023; Kühnle et al. 2025), which are assumed to be either older due to their low effective temperatures or have lower masses as they would cool faster. This shows the need to measure more isotope ratios in a larger sample of brown dwarfs in the future to identify possible trends with age, temperature or other bulk parameters.

5.5 Heavy oxygen enrichment?

We detect the oxygen isotopes ¹⁸O and ¹⁷O in water for the first time in a substellar companion’s atmosphere. The retrieved ratios between the most common isotope ¹⁶O and the rarer (and heavier) isotopes indicate a significant enrichment in the rarer isotopes. Various processes have been discussed in the literature potentially leading to such an enrichment.

In case that COCONUTS-2b formed in a disk, a possible process would be photochemical dissociation in the disk. The more abundant isotopolog gets dissociated less likely due to selfshielding in the presence of incoming radiation of the host star. This leads to an enrichment of the less abundant, heavier isotopologs in the gas phase, which could subsequently be inherited by a gas giant forming in the disk. This process has been proposed for the planetary-mass object VHS-1256b (Gandhi et al. 2023). Such a mechanism requires a disk-formation pathway followed by outward migration. However, Ciesla et al. (2026) suggests that this process (isotope-selective photodissociation) in a disk is unlikely to produce oxygen isotope fractionation at levels that would be detected in the atmospheric composition of a companion forming from the disk.

Processes inside the atmosphere are unlikely explanations for the enrichment. On Earth, fractionation of heavier water isotopologs occurs in the hydrological cycle due to differentiation through condensation and evaporation (Gat 1996). While condensation of water may occur in the cool atmosphere of COCONUTS-2b, the values of enrichment would be several orders of magnitude smaller than the values we received.

Alternatively, COCONUTS-2b might have formed from a molecular cloud that was itself enriched in heavy oxygen isotopes. In this scenario, both the companion and its host star would be expected to share similar isotope ratios. Future high-precision observations of the host star’s oxygen isotope composition will therefore be crucial for testing this hypothesis and constraining the formation pathway of the companion.

Finally, it is possible that our retrieved values underestimate the oxygen isotope ratios in the atmosphere. Unaccounted atmospheric properties, such as clouds or the presence of additional molecules, could introduce degeneracies in the retrieved ratios. Atmospheric inhomogeneity, such as patchy clouds or hot spots, could further update the retrieved parameters. In fact, models that account for inhomogeneities predict different parameters compared to homogeneous models in a study on another planetary-mass brown dwarf (Zhang et al. 2025a). Further investigations on the model improvements are beyond the scope of this paper but will be important to address in future studies. Nevertheless, the enrichment appears robust across different PT parameterizations and the values for ¹⁸O/¹⁷O are consistent with the ISM and previous observations. We are only now able to measure oxygen ratios in cold gas giant and brown dwarf atmospheres offering a new avenue for models to interpret enrichment. It is crucial to provide more observational measurements in the future to benchmark them.

5.6 Sensitivity analysis on D/H

Deuterium is especially interesting when discussing the mass estimate of a target. Brown dwarfs and gas giants are often distinguished by the deuterium burning limit at around 13 M_J (e.g., Burrows et al. 1997; Spiegel et al. 2011). Below this mass gas giants are not able to burn deuterium to He and thus deuterium-bearing isotoplogs would be present. Above this limit, objects are massive enough to burn deuterium and, depending on their mass, burn a large percentage of deuterium in the first hundreds of Myr (Spiegel et al. 2011). Thus, we have tested for the presence of deuterium-bearing species HDO and CH₃D but were unable to constrain the abundances of either of them. We nonetheless find for both upper limits beyond which the retrieval rules out the presence. For HDO the posterior distribution peaks before it falls off at the upper limit, hinting tentatively towards the presence of this molecule in the COCONUTS-2b atmosphere. CH₃D in contrast, falls off directly as seen in Fig. B.11 in the appendix.

CH₃D has been proposed to be present and detectable in brown dwarf atmospheres (Morley et al. 2019). This isotopolog was found with a proto-solar D/H ratio in the atmosphere of the coldest brown dwarf WISE0855 (Rowland et al. 2024), which was interpreted as evidence of a mass below the deuterium burning limit. The reported value of Spiegel et al. (2011) for the ratio between deuterium and hydrogen D/H for the ISM is (2.3±0.2)×10⁻⁵. Although our data do not allow us to directly constrain D/H in COCONUTS-2b, we calculated the upper bound of the D/H ratio with the derived upper limits on HDO and CH₃D according to $$\begin{gathered} (\mathrm{D}/\mathrm{H}{)}_{\mathrm{C}{\mathrm{H}}_4}=\frac{\mathrm{V}\mathrm{M}\mathrm{R}(\mathrm{C}{\mathrm{H}}_3\mathrm{D})}{4\times \mathrm{V}\mathrm{M}\mathrm{R}(\mathrm{C}{\mathrm{H}}_4)}, \\ (\mathrm{D}/\mathrm{H}{)}_{{\mathrm{H}}_2\mathrm{O}}=\frac{\mathrm{V}\mathrm{M}\mathrm{R}(\mathrm{H}\mathrm{D}\mathrm{O})}{2\times \mathrm{V}\mathrm{M}\mathrm{R}({\mathrm{H}}_2\mathrm{O})}. \end{gathered}$$

The calculations enabled us to test whether this upper limit could be consistent with the ISM value. With the equations above, our calculated 3σ (99.7% percentile) upper limits for the (D/H)_CH4 ratio are 1.42×10⁻⁵ for the constrained retrieval and 6.29×10⁻⁶ for the free case. For (D/H)_H2O we obtain an upper limit of 6.76×10⁻⁵ for the constrained retrieval and 1.34×10⁻⁴ for the free case. In Fig. B.11 in the appendix, we show the ratio in comparison with the ISM value. Interestingly, the free retrieval constrains the (D/H)_H2O value around the ISM value, and thus the upper limit we calculated would be consistent with our expectation of being consistent with the ISM value. However, the constrained case for (D/H)_H2O gives values smaller than the ISM value and falls off just right before the ISM value. As the latter PT profile is more constrained, it may not fit the small features that the free profile can, and thus it may not reach an upper upper limit for HDO. CH₃D seems to be even more difficult to constrain as both PT parameterizations result in unconstrained values for (D/H)_CH4 clearly excluding values smaller than ~10⁻⁵. The lower two panels in Fig. B.11 in the appendix show the spectral features from the retrieved case compared to a model excluding HDO or CH₃D. The spectral differences are subtle, suggesting that the expected features of those molecules lie close to the detection threshold with the MIRI/MRS data. While the current dataset does not have the S/N or spectral resolution to conclusively detect deuterium-bearing isotopes, these results provide a tantalizing hint that deuterium may be present in the atmosphere. If present, this would confirm previous studies that have determined the mass to be in the planetary regime, below the deuterium-burning limit. Therefore, future observations of COCONUTS-2b reaching a higher S/N or spectral resolution are needed to conclude on the presence of the deuterium-bearing molecules.

Fig. 11 Comparison of the spectrum of COCONUTS-2b to WISE 045 8 (Matthews et al. 2025). In panel a both spectra are plotted on top of each other, with the spectrum of WISE 045 8 scaled by a factor of 0.5. In panel b the difference between COCONUTS-2b and WISE 045 8 is shown. The purple and orange boxes indicate where in the spectrum the differences in absorption of C2H2 and HCN are located. We do not see the absorption of these molecules in COCONUTS-2b in spite of the presence in WISE 045 8 (Matthews et al. 2025) and their similar spectral type.

5.7 Non-detection of hydrocarbons

WISE J045853.90+643451.9 (hereafter WISE0458) is a T8.5-T9.0 brown dwarf binary and was observed with MIRI/MRS with an unexpected detection of HCN and C₂H₂ (Matthews et al. 2025). Due to the similar spectral type and effective temperature of ~500K, this brown dwarf binary is comparable to COCONUTS-2b; in fact up to now the only other T-dwarf observed with MIRI/MRS making COCONUTS-2b the first opportunity to study the production of these molecules in more detail. In Fig. 11, we overlay both acquired spectra, where we scale the flux of WISE0458 with a factor of 0.5. The lower panel shows the difference between the two spectra. As WISE0458 is slightly warmer in effective temperature compared to COCONUTS-2b, we see a slope in the residuals with negative residuals for shorter and positive for larger wavelengths, demonstrating that the compositions are generally similar. In WISE0458, clear absorption features of C₂H₂ and HCN are visible in the purple and orange boxes, respectively. We do not detect these molecules in COCONUTS-2b, which is indicated by the offset in the difference plot and confirmed by a corresponding retrieval analysis, where we do not detect the two molecules with a logarithmic Bayes factor of 0.19 and −0.72 for C₂H₂ and HCN, respectively.

The presence of HCN and C₂H₂ has been discussed as a general property of T-dwarf atmospheres (Matthews et al. 2025). Models predict HCN in cold atmospheres to be connected to strong vertical mixing and high surface gravities (Zahnle & Marley 2014). As we find a rather small value for the vertical mixing parameter (K_zz = 3.3) and gravity (logg≈3.9 dex), a non-detection is plausible. Explaining the presence of C₂H₂ in WISE0458 indicates an incomplete chemical network and/or more complex and unaccounted atmospheric processes (Matthews et al. 2025). The non-detection in COCONUTS-2b being in line with current models adds up to the complexity and shows the diversity in compositions for T-type objects. Thus, we can rule out their presence being strictly spectral type dependent. Instead other characteristics, such as the surface gravity, need to be taken into account. This is another indication that our current understanding of the processes leading to the disequilibrium chemistry in such cold atmospheres is incomplete. More data in the near and MIR for late T- and Y-dwarfs will help improve our understanding of carbon and nitrogen chemistry in these cold atmospheres.

5.8 A note on clouds

Salt clouds such as Na₂S and KCl have been proposed in objects at the T-/Y- transition (Morley et al. 2012; Manjavacas et al. 2022). However, in contrast to the strong silicate and iron cloud features in hotter L/T- transition objects (e.g., Suárez & Metchev 2022; Miles et al. 2023; Mollière et al. 2025), there are no distinct spectral features of salt clouds in the observed wavelength range with MIRI/MRS. Low lying silicate clouds have been proposed to change the available chemical budget through rain-out processes (Calamari et al. 2024). This has been shown to change the retrieved elemental abundances, such as C/O. To account for this, we apply the proposed correction factor to our stated C/O values. In contrast to pure absorption features, scattering from cloud particles may still affect the retrieved spectrum. However, in our analysis, we can already explain the overall flux distribution well without including any scattering effects from clouds. Especially, we are able to fit well both the Gemini/FLAMINGOS-2 and the MIRI/MRS parts together with a clear atmosphere model. We have explored various setups using cloudy retrievals including Na₂S, KCl and water clouds; however, they were not conclusive. The strongest variations in the resulting fits were located between 3 and 5 μm, where we do not include any data in this analysis. Future studies, ideally ranging over the full spectral energy distribution, may focus on whether condensing, potentially patchy, clouds still might be present.

5.9 Compositional comparison to the host star

Few observations have been taken so far of the host star L34-26. Its metallicity is estimated to be solar, with a value of [M/H] = 0.00 ± 0.08 (Hojjatpanah et al. 2019), while no measurements of its C/O ratio have been reported. Future observations would be needed to constrain the elemental and isotope abundances of the M dwarf, enabling comparative studies between the host star and the companion. In principle, similar elemental and isotopic ratios would indicate a potentially binary-system-like formation scenario, in which both the host and companion formed from the same molecular cloud. Conversely, significant differences in elemental or isotopic compositions could point to alternative formation pathway scenarios, such as disk-related processes or gravitational capture (Marocco et al. 2024). Future observations of cold T- and Y- brown dwarfs using MIRI/MRS and NIRspec, as planned in JWST programs GO 3647 (PI: P. Patapis), GO 5765 (PI: E. Matthews) and GO 8441 (PI: J. Faherty) will provide crucial comparative reference measurements for elemental and isotopic ratios in substellar atmospheres. Extending similar analyses to additional cold systems, such as Eps Ind, TWA 7, and Her 14c, will further improve our understanding of the diversity of formation pathways of cold companions. Finally, ground-based telescopes, such as the future ELT, will provide even higher resolution measurements on these faint companions along with their host stars on smaller wavelength coverages of the spectrum.

6 Summary and outlook

We have presented the MIRI/MRS observations of the cold (T_eff ≈ 480K) far-out planetary-mass companion COCONUTS-2b (Program ID: 6463, PI: P. Patapis), with an exceptional S/N of up to 40 at 11.8 μm. To study its chemical composition, we performed a comprehensive atmospheric retrieval analysis of the full resolution of the spectrum. So far, retrievals of MIRI/MRS data have been performed on binned spectra using a spectral resolution of R~1000. The higher resolution allows nitrogen and oxygen isotopologs to be detected and allows for a more detailed analysis of the composition of this atmosphere than possible in previous studies. In our setup, we chose two different PT param-eterizations, either a free or a constrained one. We summarize the results of both in the following:

• We detected the molecules H₂¹⁸O and H₂¹⁷O and constrained the ratios of ¹⁶O/¹⁸O to be 247${}_{-24}^{+27}$ and 256${}_{-25}^{+29}$, ¹⁶O/¹⁷O to be 662${}_{-83}^{+98}$ and 934${}_{-139}^{+174}$, and ¹⁸O/¹⁷O 2.7${}_{-0.4}^{+0.6}$ and 3.7${}_{-0.7}^{+0.8}$ for the free and constrained retrieval, respectively, indicating an enrichment in the heavier oxygen isotopes;

• We detected the nitrogen isotopolog ¹⁵NH₃ and constrained the ratio to be 337${}_{-38}^{+47}$ and 324${}_{-40}^{+46}$ for the free and constrained profile, respectively. The values are within one to two sigma of the value we found for the ISM and previous Y dwarf observations and are therefore compatible;

• We did not detect HCN and C₂H₂, and we only detected PH₃, CO₂, and ¹³CH₄ tentatively using the Bayes factor comparison;

• As we did not constrain HDO and CH₃D but obtained only upper limits for both species, we likely are at the sensitivity limit to detect deuterium-bearing molecules;

• We found a subsolar to solar metallicity and a subsolar C/O ratio compatible with previous measurements.

The measured isotopolog abundances will be important when making comparisons to the composition of host stars in order to study the formation of COCONUTS-2b. In addition, future data with an even higher S/N or higher spectral resolution will potentially result in more accurate ratios and might allow us to detect deuterium-bearing species - if present. This may provide an additional constraint on the mass estimate for this wide-separation companion, whose dynamical mass will be difficult to obtain. An atmospheric retrieval analysis using the NIRSpec dataset will potentially lead to a better understanding of the vertical mixing due to better constraints on CO and CO₂. Further comparisons using self-consistent grid models, such as the study done on MIRI/LRS data, will be crucial for a broader characterization of COCONUTS-2b (Ravet et al., in review). Combining the NIRSpec, MIRI/LRS, and MIRI/MRS datasets will ultimately allow more complex features to be included in our models, such as clouds or vertically changing abundance profiles. Thus, the atmosphere of cold COCONUTS-2b still holds many mysteries waiting to be solved.

Data availability

The raw data of the MIRI/MRS dataset used in this analysis is available on MAST (https://mast.stsci.edu/) under the PID: 6463 or directly through the DOI: https://doi.org/10.17909/12mr-hg42. The reduced data will be made available under the following link: https://doi.org/10.5281/zenodo.19887581.

Acknowledgements

We thank the anonymous referee for the comments that improved the quality of this paper. We want to thank Céline Nussbaumer and Jérémie Pierre for their insights on the retrieval setups through their semester projects and Jean Hayoz and Henrik Knierim for fruitful discussions. Further, we would like to thank Fabian Grübel for valuable additions. PP thanks the Swiss National Science Foundation (SNSF) for financial support under grant number 200020_200399. DB is supported by Spanish MCIN/AEI/10.13039/50H000H033 grant PID2023-150468NB-I00 and No. MDM-2017-0737. NW acknowledges funding from NSF award #223 8468, #1909776, and NASA Award #80NSSC22K0142. NW acknowledges support for program #06463 was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-03127. MB acknowledges support in France from the French National Research Agency (ANR) through project grants ANR-20-CE31-0012 and ANR-23-CE31-0006 and the Action Thématique Exosystèmes as part of the CNRS/INSU.

References

Appendix A Observation

We present the S/N for the resulted observation form 5 to 21 μm in Fig. A.1. The maximum median S/N across the channels was reached in channel 3A with 22.1 and reaching a maximum S/N of up to ~40. The lowest was reached in channel 4A with a value of 2.7. Thus, we do not include the latter in the analysis. In Fig. A.2 we present the errors on the dataset in mJy. The largest errors are found in channel 4A and 3C after that. The lowest can be found in channel 3A.

Appendix B Retrievals

In this appendix we provide several supporting plots for the retrieval analysis. We present in Fig. B.1 the opacities used in ck resolution. We also highlight the Gemini/FLAMINGOS-2 and MIRI/MRS covered wavelength ranges. We present the median, 16th, and 84th percentiles of the posteriors of the low resolution retrievals in Tab. B.1. The full resolution retrieval comparisons for species not clearly detected are shown for ¹³CH₄ and C₂H₂ in Fig. B.2, for H₂S, PH₃, and CO₂ in Fig. B.3, and HCN, HDO, and CH₃D in Fig. B.4. Here, we present the areas in the spectrum where the difference in the model is the strongest. In Fig. B.5, B.6, and B.7 we present more wavelength areas in the spectrum compared to the ones presented in Fig. 6,7, and 8, where the absorption of ¹⁵NH₃, H₂¹⁸O and H₂¹⁷O are visible, respectively. The comparison of all retrieved molecular abundances including the ones from the low and full resolution retrievals are summarized in Fig. B.8. In Fig B.9 we present the correlations between some of the bulk parameter retrieved and calculated from the low resolution retrieval results. Fig. B.10 shows how the variations in the PT profile vary the spectrum either including the isotopologs or not. We present four cases, where we either increase or decrease the bottom temperature or vary the slopes to either a stepper or shallower profile. In Fig. B.11, we present the (D/H)_H2O and (D/H)_CH4 ratio for the free and constrained case and compare them to the ISM value.

Fig. A.1 Signal-to-noise ratio per channel for the MIRI/MRS observation with median values highlighted in dashed lines and presented in the legend.

Fig. A.2 Error from the pipeline across the wavelengths.

Fig. B.1 Opacities used for the retrievals at $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$.

Fig. B.2 Best-fit retrievals on full MIRI/MRS resolution with and without including 13CH4 and C2H2 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.

Fig. B.3 Best-fit retrievals on full MIRI/MRS resolution with and without including H2S, PH3 and CO2 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.

Fig. B.4 Best-fit retrievals on full MIRI/MRS resolution with and without including HCN, HDO and CH3D in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.

Fig. B.5 Best-fit retrievals on full MIRI/MRS resolution with and without including 15NH3 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the fourth to eight largest differences between the best-fits.

Fig. B.6 Best-fit retrievals on full MIRI/MRS resolution with and without including H218O in blue and red respectively compared to the data in black with the same panel structure as in Fig. 6. Here we show the part of the spectrum with the largest, second, third, seventh and eight largest differences between the best-fits.

Fig. B.7 Best-fit retrievals on full MIRI/MRS resolution with and without including H217O in blue and red respectively compared to the data in black with the same panel structure as in Fig. 7. Here we show the part of the spectrum with the fourth to eight largest differences between the best-fits.

Fig. B.8 Abundances in comparison between the full and the low resolution as well as the base and the LOO retrievals.

Fig. B.9 Corner plot of the retrieved and calculated bulk parameters for the free and the constrained case in green and violet, respectively: the radius, logg, mass, and metallicity. These results are based on the low resolution retrievals

Fig. B.10 Effect on the model when varying the PT profile. Panel a) shows the variation of the PT profile for either a shift or a change in the slope of the PT profile. The effect of on the spectrum either including (’iso’) or neglecting the isotopolog is shown in panels b)-e). In b) we shift the PT profile by 500K towards lower temperatures, c) a shift by 500K towards higher temperatures, d) a steeper profile by multiplying each node by 10%, and e) a shallower profile by 10%.

Fig. B.11 The calculated (D/H)H2O ratio from the upper limits for the free and constrained retrieval in green and violet and the ratio (D/H)CH4 for the free and constrained case in blue and pink compared to the ISM value (Spiegel et al. 2011) in black in panel a). Panel b) presents the modeled spectrum without HDO and CH3D in yellow, the retrieved best-fit in violet including the two molecules compared to the data. The residuals to plot b) are shown in panel c).

All Tables

All Figures

Fig. 1 Observation and cube detector images of channel 1A to 4A. Channel 4A was neglected from further analysis (see text). The red cross corresponds to the position of the target and the blue circle to the aperture from which the flux has been extracted from. The two negative spots on the detector correspond to two nods originating from the background subtraction. The yellow circle shows the position of a background object.

In the text

Fig. 2 Best-fit spectra of the free and constrained atmospheric retrievals, respectively in green and violet, compared to the data in black. Panel a shows the full observation. In panel b, we show a subset at the Gemini wavelengths from 1 to 2.5 μm as indicated by the black box and in panel c we show another subset for the strong NH3 feature between 8 and 12 μm. The spectra are shown in lower resolution $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$.

In the text

Fig. 3 Bulk parameter for the free and the constrained retrievals. We show in panel a, the effective temperature estimate, in panel b the retrieved metallicity, in panel c the C/O ratio accounting for oxygen sequestration, in panel d the radius, in panel e the surface gravity, and in panel f the resulting mass based on the radius and gravity for the free and constrained retrieval (in green and violet, respectively). For panels a, b, and d-f, we show the mean values with dashed lines and the one sigma estimate with dotted lines, for the values obtained by Zhang et al. (2025b, Z25) in black and Kiman et al. (2026, K25) in blue.

In the text

	Fig. 4 PT profiles of the free and the constrained retrieval (in green and violet, respectively). In dark blue we show the PT structure of Jupiter. The H2O, NH4SH, and NH3 condensation lines were taken from (Lodders & Fegley 2002), and the rest is from petitRADTRANS. Thicker lines in deeper atmospheric layers indicate the convective part of the atmosphere calculated based on Eq. (1).
In the text

Fig. 5 Retrieved chemical abundances with the free and constrained retrievals in the spectral resolution $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$ (in green and violet, respectively) compared to the chemical equilibrium predictions for the free and the constrained retrievals (in dark and light gray, respectively). In panel a the abundances were taken at 0.28 bar, corresponding to the maximum contribution in the free case. In panel b the variations with pressure of the abundances involved in the chemical equilibrium, CO; CH4 , H2O, and CO2 are presented. We additionally indicate the pressure of the maximum contribution (0.28 bar) in blue and the quench pressure (43-53 bar) in magenta.

In the text

Fig. 6 Best-fit retrievals on full MIRI/MRS resolution data with and without including H218O (in blue and red, respectively) compared to the data in black. In panels a-c, one can see the three different parts in the spectrum where the difference in the best-fits is visible. The residuals between the data and the corresponding best-fits are shown in panels d-f and the differences in the fits in plots g-i. Here we present the wavelength ranges with the fourth, fifth, and sixth largest differences in the models compared to showing the ranges with the three largest difference, as there the models do not explain the data well overall. The latter and more wavelength areas are presented in the appendix in Fig. B.6.

In the text

	Fig. 7 Best-fit retrievals on full MIRI/MRS resolution data with and without including H217O (in blue and red, respectively) compared to the data in black with the same panel structure as in Fig. 6. We present the wavelength ranges with the largest, second largest, and sixth largest differences in the models. We present more wavelength areas in the appendix in Fig. B.7.
In the text

	Fig. 8 Best-fit retrievals on full MIRI/MRS resolution data with and without including 15NH3 (in blue and red, respectively) compared to the data in black with the same panel structure as in Fig. 6. Here, we present the wavelength ranges with the largest, second largest, and third largest differences in the models. More wavelength areas are provided in the appendix in Fig. B.5.
In the text

Fig. 9 Ratios of a) 14N/15N, b) 16O/18O, c) 16O/17O, and d) 18O/17O based on the full resolution retrievals. We compare the ratios to the values for the ISM (oxygen: Wilson (1999), 14N/15N: Ritchey et al. (2015) and solar values (16O/17O and 16O/18O: McKeegan et al. (2011), 18O/17O: Wilson (1999) and nitrogen: Marty et al. (2011) as well to the values of WISE 1828 (Barrado et al. 2023), WISE 0855 (Kühnle et al. 2025) for 14N/15N, and VHS1256b for the oxygen ratios (Gandhi et al. 2023).

In the text

	Fig. 10 Comparison of low and full resolution retrievals compared to the data. Panel a shows the best-fit spectra for the low resolution in green and the full resolution retrievals with and without 15NH3 (in blue and red, respectively). In panel b we show the corresponding residuals. Only with the improvement of the resolution were we able to detect the isotopolog.
In the text

Fig. 11 Comparison of the spectrum of COCONUTS-2b to WISE 045 8 (Matthews et al. 2025). In panel a both spectra are plotted on top of each other, with the spectrum of WISE 045 8 scaled by a factor of 0.5. In panel b the difference between COCONUTS-2b and WISE 045 8 is shown. The purple and orange boxes indicate where in the spectrum the differences in absorption of C2H2 and HCN are located. We do not see the absorption of these molecules in COCONUTS-2b in spite of the presence in WISE 045 8 (Matthews et al. 2025) and their similar spectral type.

In the text

	Fig. A.1 Signal-to-noise ratio per channel for the MIRI/MRS observation with median values highlighted in dashed lines and presented in the legend.
In the text

	Fig. A.2 Error from the pipeline across the wavelengths.
In the text

	Fig. B.1 Opacities used for the retrievals at $\frac{\lambda }{\mathrm{\Delta }\lambda }=1000$.
In the text

	Fig. B.2 Best-fit retrievals on full MIRI/MRS resolution with and without including 13CH4 and C2H2 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.
In the text

	Fig. B.3 Best-fit retrievals on full MIRI/MRS resolution with and without including H2S, PH3 and CO2 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.
In the text

	Fig. B.4 Best-fit retrievals on full MIRI/MRS resolution with and without including HCN, HDO and CH3D in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the normalized largest difference between the best-fits.
In the text

	Fig. B.5 Best-fit retrievals on full MIRI/MRS resolution with and without including 15NH3 in blue and red respectively compared to the data in black with the same panel structure as in Fig. 8. Here we show the part of the spectrum with the fourth to eight largest differences between the best-fits.
In the text

	Fig. B.6 Best-fit retrievals on full MIRI/MRS resolution with and without including H218O in blue and red respectively compared to the data in black with the same panel structure as in Fig. 6. Here we show the part of the spectrum with the largest, second, third, seventh and eight largest differences between the best-fits.
In the text

	Fig. B.7 Best-fit retrievals on full MIRI/MRS resolution with and without including H217O in blue and red respectively compared to the data in black with the same panel structure as in Fig. 7. Here we show the part of the spectrum with the fourth to eight largest differences between the best-fits.
In the text

	Fig. B.8 Abundances in comparison between the full and the low resolution as well as the base and the LOO retrievals.
In the text

	Fig. B.9 Corner plot of the retrieved and calculated bulk parameters for the free and the constrained case in green and violet, respectively: the radius, logg, mass, and metallicity. These results are based on the low resolution retrievals
In the text

Fig. B.10 Effect on the model when varying the PT profile. Panel a) shows the variation of the PT profile for either a shift or a change in the slope of the PT profile. The effect of on the spectrum either including (’iso’) or neglecting the isotopolog is shown in panels b)-e). In b) we shift the PT profile by 500K towards lower temperatures, c) a shift by 500K towards higher temperatures, d) a steeper profile by multiplying each node by 10%, and e) a shallower profile by 10%.

In the text

Fig. B.11 The calculated (D/H)H2O ratio from the upper limits for the free and constrained retrieval in green and violet and the ratio (D/H)CH4 for the free and constrained case in blue and pink compared to the ISM value (Spiegel et al. 2011) in black in panel a). Panel b) presents the modeled spectrum without HDO and CH3D in yellow, the retrieved best-fit in violet including the two molecules compared to the data. The residuals to plot b) are shown in panel c).

In the text