Using Keck Planet Imager and Characterizer highresolution R  35000 spectroscopy from 2.29 to 2.49 m we present uniform atmospheric retrievals for eight young substellar companions with masses of 1030 M SUBJupSUB orbital separations spanning 50360 au and T SUBeffSUB between 1500 and 2600 K. We find that all companions have solar CO ratios and metallicities to within the 12 level with the measurements clustered around solar composition. Stars in the same stellar associations as our systems have nearsolar abundances so these results indicate that this population of companions is consistent with formation via direct gravitational collapse. Alternatively core accretion outside the CO snowline would be compatible with our measurements though the high mass ratios of most systems would require rapid core assembly and gas accretion in massive disks. On a population level our findings can be contrasted with abundance measurements for directly imaged planets with m lt 10 M SUBJupSUB which show tentative atmospheric metal enrichment compared to their host stars. In addition the atmospheric compositions of our sample of companions are distinct from those of hot Jupiters which most likely form via core accretion. For two companions with T SUBeffSUB  17002000 K  And b and GSC 6214210 b our bestfit models prefer a nongray cloud model with gt3 significance. The cloudy models yield 23 lower T SUBeffSUB for these companions though the CO and CH still agree between cloudy and clear models at the 1 level. Finally we constrain SUP12SUPCOSUP13SUPCO for three companions with the highest signaltonoise ratio data GQ Lup b HIP 79098b and DH Tau b and report inlineformula mmlmath overflowscrollmmlmivmmlmimmlmisinmmlmimmlmiimmlmimmlmath inlineformula and radial velocities for all companions.
['Exoplanet atmospheres', 'Brown dwarfs', 'High resolution spectroscopy', 'Exoplanet atmospheric composition', 'Atmospheric clouds', 'Exoplanet formation', 'Substellar companion stars', '487', '185', '2096', '2021', '2180', '492', '1648', 'Astrophysics - Earth and Planetary Astrophysics', 'Astrophysics - Solar and Stellar Astrophysics']
Are These Planets or Brown Dwarfs Broadly Solar Compositions from Highresolution Atmospheric Retrievals of 1030 M SUBJupSUB Companions
{'Are these planets or brown dwarfs? Broadly solar compositions from high-resolution atmospheric retrievals of ∼ 10-30 M Jup companions': "Jerry W. Xuan , 1 Chih-Chun Hsu , 2 Luke Finnerty , 3 Jason J. Wang , 4 Jean-Baptiste Ruffio , 5 Yapeng Zhang , 1 Heather A. Knutson , 6 Dimitri Mawet, 1, 7 Eric E. Mamajek , 7 Julie Inglis , 6 Nicole L. Wallack , 8 Marta L. Bryan , 9 Geoffrey A. Blake, 6 Paul Molli'ere, 10 Neda Hejazi , 11 Ashley Baker, 1 Randall Bartos, 7 Benjamin Calvin , 3 Sylvain Cetre, 12 Jacques-Robert Delorme , 12 Greg Doppmann, 12 Daniel Echeverri , 1 Michael P. Fitzgerald , 3 Nemanja Jovanovic , 1 Joshua Liberman , 1, 13 Ronald A. L'opez , 3 Evan Morris, 14 Jacklyn Pezzato, 1 Ben Sappey , 5 Tobias Schofield, 1 Andrew Skemer, 14 J. Kent Wallace , 7 Ji Wang , 15 Shubh Agrawal , 16 and  \nKatelyn Horstman 1  \n1 Department of Astronomy, California Institute of Technology, Pasadena, CA 91125, USA  \n2 Center for Interdisciplinary Exploration and Research in Astrophysics (CIERA), Northwestern University, 1800 Sherman, Evanston, IL, 60201, USA  \n3 Department of Physics & Astronomy, 430 Portola Plaza, University of California, Los Angeles, CA 90095, USA  \n4 Center for Interdisciplinary Exploration and Research in Astrophysics (CIERA) and Department of Physics and Astronomy, Northwestern University, Evanston, IL 60208, USA  \n5 Department of Astronomy & Astrophysics, University of California, San Diego, La Jolla, CA 92093, USA  \n6 Division of Geological & Planetary Sciences, California Institute of Technology, Pasadena, CA 91125, USA  \n- 7 Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Dr.,Pasadena, CA 91109, USA  \n8 Earth and Planets Laboratory, Carnegie Institution for Science, Washington, DC, 20015  \n9 David A. Dunlap Institute Department of Astronomy & Astrophysics, University of Toronto, 50 St. George Street, Toronto, ON M5S 3H4, Canada  \n10 Max-Planck-Institut fur Astronomie, Konigstuhl 17, 69117 Heidelberg, Germany  \n11 Department of Physics and Astronomy, University of Kansas, Lawrence, KS 66045, USA  \n12 W. M. Keck Observatory, 65-1120 Mamalahoa Hwy, Kamuela, HI, USA  \n13 James C. Wyant College of Optical Sciences, University of Arizona, Meinel Building 1630 E. University Blvd., Tucson, AZ. 85721  \n14 Department of Astronomy & Astrophysics, University of California, Santa Cruz, CA95064, USA  \n15 Department of Astronomy, The Ohio State University, 100 W 18th Ave, Columbus, OH 43210 USA  \n16 Department of Physics and Astronomy, University of Pennsylvania, Philadelphia PA 19104", 'ABSTRACT': 'Using Keck Planet Imager and Characterizer (KPIC) high-resolution ( R ∼ 35000) spectroscopy from 2 . 29 -2 . 49 µ m, we present uniform atmospheric retrievals for eight young substellar companions with masses of ∼ 10-30 M Jup , orbital separations spanning ∼ 50-360 au, and T eff between ∼ 1500-2600 K. We find that all companions have solar C/O ratios, and metallicities, to within the 1-2 σ level, with the measurements clustered around solar composition. Stars in the same stellar associations as our systems have near-solar abundances, so these results indicate that this population of companions is consistent with formation via direct gravitational collapse. Alternatively, core accretion outside the CO snowline would be compatible with our measurements, though the high mass ratios of most systems would require rapid core assembly and gas accretion in massive disks. On a population level, our findings can be contrasted with abundance measurements for directly imaged planets with m< 10 M Jup , which show tentative atmospheric metal enrichment compared to their host stars. In addition, the atmospheric compositions of our sample of companions are distinct from those of hot Jupiters, which most likely form via core accretion. For two companions with T eff ∼ 1700-2000 K ( κ Andb and GSC 6214-210 b), our best-fit models prefer a non-gray cloud model with > 3 σ significance. The cloudy models yield 2 -3 σ lower T eff for these companions, though the C/O and [C/H] still agree between cloudy and  \nclear models at the 1 σ level. Finally, we constrain 12 CO / 13 CO for three companions with the highest S/N data (GQ Lup b, HIP 79098 b, and DH Tau b), and report v sin i and radial velocities for all companions.', '1. INTRODUCTION': "High-contrast imaging surveys have revealed a population of substellar companions, generally classified as giant planets ( ∼ 2 -13 M Jup ) or brown dwarfs ( ∼ 13 -75 M Jup ), orbiting at large separations ( ∼ 3 -1000 au) from their host stars (see reviews by Bowler 2016; Currie et al. 2023). Between giant planets and brown dwarfs, there are also dozens of low-mass substellar companions ( m ∼ 10-30 M Jup ) at wide orbital separations (dozens to hundreds of au). These objects have often been termed 'planetary-mass companions' (e.g. Ireland et al. 2010; Kraus et al. 2013), though there is no conclusive evidence as to whether they form like planets. Insights into their formation processes would help provide more physically-based definitions for giant planets and brown dwarfs (Schlaufman 2018), with giant planets being the product of bottom-up core accretion (Pollack et al. 1996), and brown dwarfs the product of top-down gravitational collapse either in a disk or molecular cloud (e.g. Offner et al. 2010; Bate 2012; Kratter & Lodato 2016).  \nThese widely-separated 10-30 M Jup companions have occurrence rates of only a few percent (e.g. Nielsen et al. 2019). The rarity of these companions aligns with the difficulties that they pose to both planet-like and starlike formation processes. The current orbital locations of many of these companions are too far for either core accretion or disk instability to operate efficiently given low surface densities at large distances ( > 100 au) in the disk (e.g. Dodson-Robinson et al. 2009). On the other hand, cloud fragmentation has issues explaining the extreme mass ratios (a few percent) of these systems (e.g. Bate 2012). If these companions form via core accretion at closer distances followed by outward scattering, there should be close-in companions that served as the scatters, which have not yet been detected (Bryan et al. 2016; Pearce et al. 2019).  \nTo understand the nature of directly imaged companions, the field has focused on two complementary approaches. The first examines their orbital architectures as a function of companion mass. Such studies have found evidence for distinct distributions of semi-major axis, orbital eccentricity, and stellar obliquity around a dividing mass of ∼ 10 -20 M Jup (Nielsen et al. 2019; Bowler et al. 2020; Bowler et al. 2023; Nagpal et al. 2023), though the exact results can be sensitive to the specific dividing mass (Do ' O et al. 2023). This suggests a fuzzy boundary between giant planets and brown dwarfs and highlights the importance of further under-  \nnding the intermediate-mass companions with ∼ 1030 M Jup . The second approach relies on the analysis of spectro-photometry, which contains information about the physical processes and chemical inventory of their atmospheres. Indeed, the atmospheric abundances of substellar companions encode fossil information about their accretion histories, and could potentially inform different formation scenarios (e.g. Gravity Collaboration et al. 2020; Molli'ere et al. 2022).  \nEarly studies highlighted the carbon and oxygen abundances of the atmosphere as informative observables (e.g. Oberg et al. 2011; Madhusudhan 2012). To first order, the C/O ratios of solids in the disk are predicted to vary as a function of disk radius. Planets that form via core accretion in a protoplanetary disk, which is a relatively slow process occurring on Myr timescales, can incorporate varying quantities of gas and solids into their atmospheres, potentially resulting in a wide range of atmospheric metallicities and C/O ratios. On the other hand, companions that form rapidly on dynamical timescales via direct gravitational collapse are expected to inherit C/O and metallicities similar to those of their host stars, analogous to the case of binary star systems (Hawkins et al. 2020).  \nHowever, these predictions can be complicated by a range of effects. In particular, if widely-separated ( ≳ 100 au), ∼ 10-30 M Jup companions can form via core accretion outside the CO snowline, they are expected to have stellar C/O and metallicities as the solids at these locations are of stellar composition and the bulk of the metals is in the solid phase (Chachan et al. 2023). We may also see systematically lower atmospheric metallicities for objects that form via core accretion in the outer disk, as small grains could rapidly drain inward to the star in the absence of pressure gaps. On the other hand, if these companions form via disk instability, pressure bumps and spiral structures can lead to local enhancements or reductions in the disk metallicity that can be inherited by the companions (e.g. Boley et al. 2011). Therefore, mapping the measured composition of a single planet/companion to a specific formation pathway is far from a simple one-to-one process and requires sophisticated disk models to fully disentangle the intricacies (Molli'ere et al. 2022).  \nWhile there could be significant uncertainty in interpreting the composition for a single object, a dominant process would be more apparent as a trend in the population. In this regard, Wang (2023a) noted that several imaged planets with m ≈ 3 -13 M Jup have metallicities  \nhigher than their star's by ≈ 0.1-0.7 dex, with typical errors of 0.2 dex (Gravity Collaboration et al. 2020; Petrus et al. 2021; Brown-Sevilla et al. 2023; Wang et al. 2023; Molli'ere et al. 2020; Wang 2023a), suggesting that they may have formed via core accretion. Recently, Zhang et al. (2023) also reported a potential ∼ 30-170 × metal enrichment for AF Lep b ( ≈ 3 M Jup ) relative to its star. These measurements are not without caveats. For example, Zhang et al. (2023) could not reliably constrain the C/O of AF Lep b from their low-resolution data and Landman et al. (2023) showed that high-resolution retrievals of β Pic b lead to sub-stellar metallicity that disagree with the super-stellar metallicity found by Gravity Collaboration et al. (2020) using low-resolution data. Despite these caveats, there is a possible trend of superstellar metallicities for at least some directly imaged planets. In contrast, higher mass ( m ∼ 50 -70 M Jup ) brown dwarf companions generally exhibit both C/O and bulk metallicities consistent with their stellar values (e.g. Line et al. 2015; Xuan et al. 2022; Wang et al. 2022; Phillips et al. 2024), which is expected given their presumed binary-star like formation pathways. A few exceptions to this trend of chemical homogeneity between high mass BDs and their stars have been attributed to missing physics in the modeling, rather than real differences (e.g. Calamari et al. 2022; Balmer et al. 2023).  \nSystematically measuring the atmospheric compositions of ∼ 10-30 M Jup companions could help determine their nature. To date, only a few of these companions have reported abundance measurements (Hoch et al. 2020; Zhang et al. 2021; Hoch et al. 2022; Palma-Bifani et al. 2023; Demars et al. 2023; Inglis et al. 2023), with a trend of approximately solar C/O values (Hoch et al. 2023). With the exception of Inglis et al. (2023) however, all these studies employed medium-resolution spectroscopy ( R ∼ 4000) and used self-consistent grid models to estimate the companion's abundances. Grid models with low dimensionality can provide poor fits to data but yield unrealistically tight constraints on the model parameters (e.g. Ruffio et al. 2021), though this could be accounted for by inflating the uncertainties (Hoch et al. 2020). In addition, older grid models may contain outdated line lists, while the retrieval approach enables the incorporation of new line lists more easily. Retrievals also allow more flexibility in defining the cloud models and fitting for isotopic abundances. However, retrievals are not without caveats either. For example, retrieval studies often produce unphysically small radii (e.g. Gonzales et al. 2020; Burningham et al. 2021; Lueber et al. 2022; Hood et al. 2023), and overly isothermal profiles which might suggest inadequacies in the cloud models (e.g. Burningham et al. 2017; Brown-Sevilla et al.  \n2023). Ultimately, it is important to compare both approaches, for example, by using the information from self-consistent thermal profiles as informed priors in retrievals (Zhang et al. 2023; Xuan et al. 2024).  \nIn this paper, we present systematic atmospheric retrievals for a sample of eight young ( ∼ 1 -100 Myr) companions with m ∼ 10-30 M Jup using Keck Planet Imager and Characterizer (KPIC) high-resolution spectroscopy ( R ∼ 35,000, K band). To ensure physical solutions, our retrievals are informed by mass and radius priors from evolutionary models, and self-consistent thermal profiles following Xuan et al. (2024). We measure the C and O abundances of all companions and constrain the 12 C / 13 C for three companions with the highest signalto-noise (S/N) data. With the statistical leverage of our sample and uniform analysis framework, we aim to understand whether this population of objects is more akin to high-mass giant planets or low-mass brown dwarfs.  \nThis paper is organized as follows: § 2 overviews the known properties of our eight systems and uniformly estimates the stellar ages and companion bulk properties from evolutionary models. In § 3, we describe the KPIC observations and data reduction. § 4 lays out our spectral analysis framework, including the atmospheric retrieval setup. The results of our retrievals are summarized in § 5, and we discuss the implications of our measurements and analysis in § 6. Finally, we present our conclusions in § 7.", '2. SYSTEM PROPERTIES': "In Fig. 1, we place our sample in the context of directly imaged companions with both C/O and metallicity measurements. 1 Our sample consists of six late K to early M type stars (GQ Lup, DH Tau, ROXs 12, ROXs 42B, 2M0122, GSC 6214-610) and two B stars (HIP 79098, κ And), and are either confirmed or likely  \n1 Before this work, directly imaged companions with measured abundances from both retrievals and grid model fits are: AF Lep b (Zhang et al. 2023; Palma-Bifani et al. 2024), HR 8799 b, c, d, e (Lavie et al. 2017; Wang et al. 2020; Molli'ere et al. 2020; Wang et al. 2023; Konopacky et al. 2013; Barman et al. 2015; Ruffio et al. 2021), β Pic b (Gravity Collaboration et al. 2020; Landman et al. 2023), HIP 62426 b (Wang 2023b; Petrus et al. 2021), HD 206893 B (Kammerer et al. 2021), VHS J1256-1257 b (Gandhi et al. 2023; Petrus et al. 2023a,b). The objects with abundances from retrievals are 51 Eri b (Whiteford et al. 2023; Brown-Sevilla et al. 2023), YSES 1 b (Zhang et al. 2021), HD 4747 B (Xuan et al. 2022), HR 7276 B (Wang et al. 2022), Gl 229 B (Calamari et al. 2022; Howe et al. 2022), HD 72426 B (Balmer et al. 2023), HD 33632 B (Hsu et al. 2024a), HD 984 B (Costes et al. 2024). Finally, objects with abundances from grid models are PDS 70 b and c (Wang et al. 2021b), κ Andb (Hoch et al. 2020), HD 284149 b (Hoch et al. 2022), AB Pic b (PalmaBifani et al. 2023), GQ Lup b and GSC 6214-210 b (Demars et al. 2023).  \nFigure 1. Confirmed directly imaged substellar companions that have published C/O and metallicity values from retrievals (red) and/or grid model fits (blue). Objects that have abundance measurements from both grid models and retrievals are shown as red points with a blue outline. The eight companions studied in this paper are denoted with a black outline. As shown, some of our objects have previous abundance measurements, which we summarize in § 2.3.  \n<!-- image -->  \nmembers of various nearby star-forming regions and/or young moving groups. Indeed, five of the systems are located in the Scorpius-Centaurus association. Below, we summarize the properties of the systems, with a focus on parameters relevant to our retrieval study.", '2.1. Stellar Ages': "The stellar, and by extension, system ages can inform the evolutionary states of substellar companions, including their radius and mass. To estimate ages of our stars, we either perform isochrone fitting or adopt a literature age when isochrone fitting is complicated by factors such as unresolved binarity. For isochrone fitting, we use the Baraffe et al. (2015) models (BHAC15) and fit using literature bolometric luminosity and T eff measurements for the stars. We inflate the T eff error bars to 150 K to be conservative, which is larger or equal to the reported T eff uncertainties for our stars. For L bol , we apply a correction based on the stars' Gaia DR3 parallaxes, as several measurements were reported using a pre-Gaia distance. We implement a rejection-sampling method to interpolate the models in mass and age space following Dupuy & Liu (2017). In short, for each mass-age pair, we compare the interpolated T eff and L bol to the measured values and accept values based on their probability. The median and 68% credible interval of the resulting age posteriors are listed in Table 1, and as we  \ndiscuss in § 2.3, are all consistent with previous age estimates in the literature.", '2.2. Companion Mass and Radius': "Given the stellar ages and literature measurements of the companions' L bol , we can derive the expected companion mass and radius from substellar evolutionary models that have been shown to reasonably reproduce the bulk properties of benchmark substellar companions with dynamical masses (e.g. Dupuy & Liu 2017) - although the models can differ from each other. To derive these priors while accounting for model uncertainty, we consider four different evolutionary models: ATMO 2020 (Phillips et al. 2020; Chabrier et al. 2023), SM08 (Saumon & Marley 2008), AMES-Dusty (Allard et al. 2001), and BHAC15 (Baraffe et al. 2015). 2  \nWith the same rejection-sampling technique as described above, we use the stellar ages and companion L bol to derive posteriors for mass, radius, T eff , and log g . As described in § 4.3, we will use the evolutionary-model derived mass and radius as priors in our atmospheric retrievals. Since different models predict slightly different mass and radius, we visually determine mass and radius priors that encompass most of the posterior range for different models. These priors are listed in Table 1, and we show an example in Fig. 2. Plots for other companions are included in Appendix A. We note that our derived masses and radii are consistent with previous estimates in the literature for each companion, and we rederive them for the sake of uniformity.", '2.3. Overview of the Systems': 'Below, we summarize the properties of each system, including our derived stellar ages and companion properties. Note that most of our companions show little orbital motion since their discovery, so we quote their projected orbital separations at the epoch of the KPIC observations. Statistically, the most likely orbital semimajor axis is similar to the observed projected separation (Yelverton et al. 2019).', '2.3.1. ROXs 42B': "ROXs 42B is a resolved binary with K -band flux ratio ∼ 3 (Ratzka et al. 2005) and a member of the ρ Ophiuchus cloud ( ρ Oph). Our isochrone fits for the primary component yield an age of 2 . 2 +1 . 8 -1 . 0 Myr, consistent with the ρ Oph age from Miret-Roig et al. (2022),  \nFigure 2. The interpolated mass, radius, log g , and T eff from three different evolutionary models (ATMO 2020 in blue, Phillips et al. 2020; AMES-Dusty in red, Allard et al. 2001; BHAC15 in purple, Baraffe et al. 2015) for GQ Lup b. The dashed black lines indicate the adopted priors for mass and radius in the retrievals. Plots for the other companions are shown in Appendix A.  \n<!-- image -->  \nwho confirmed the star's ρ Oph membership from Gaia DR3 kinematics. A candidate companion was identified around the binary by Ratzka et al. (2005) and later confirmed by Kraus et al. (2013). The companion is located at a projected separation ρ ≈ 1 . 2 '' from the central binary, or about 170 au. 3 The companion has a spectral type of L 1 ± 1 (Bowler et al. 2014a). Using the companion's L bol from Currie et al. (2014), we estimate m = 13 ± 5 M Jup and r = 2 . 10 ± 0 . 35 R Jup for ROXs 42B b, consistent with previous studies (e.g. Kraus et al. 2013; Currie et al. 2014). An early retrieval study was performed by Daemgen et al. (2017), who used only 15 µ m photometry and did not provide constraints on the companion's chemical abundances. Recently, Inglis et al. (2023) performed retrievals on pre-upgrade Keck/NIRSPEC spectra ( R ∼ 25,000) of the companion, finding C / O = 0 . 50 ± 0 . 05, [Fe / H] = -0 . 67 ± 0 . 35, and v sin i = 10 . 5 ± 0 . 9 km s -1 . ROXs 42 b has a mid-IR excess from Spitzer indicative of a circumsubstellar disk (Martinez & Kraus 2021), although the companion does not show any accretion features.", '2.3.2. ROXs 12': "ROXs 12 is likely a member of ρ Oph (Miret-Roig et al. 2022). However, Luhman (2022) assign the star memberships of either ρ Oph or Upper Scorpius (Upper Sco). For consistency, we perform isochrone fitting to obtain 6 . 5 +3 . 8 -2 . 6 Myr, consistent with the age estimate from Kraus et al. (2013). A candidate companion to ROXs 12 was first noted by Ratzka et al. (2005) and later confirmed by Kraus et al. (2013). The com-  \nnion (ROXs 12 b) is located at a projected separation of ≈ 1.8 '' or about 250 au. Based on the companion's L bol estimated by Bowler et al. (2017), we find m =19 ± 5 M Jup and r =2.2 ± 0 . 35 R Jup . Bowler et al. (2017) perform a detailed characterization of the system, and we summarize the results. They determine a spectral type of L 0 ± 2 for the companion and find that the companion is likely on a misaligned orbit relative to the host star's spin axis. These authors also find evidence of an outer tertiary component in the system at 5000 au, which shares common proper motion and radial velocity as ROXs 12. A lack of Pa β emission indicates there is no evidence of a disk around ROXs 12 b. Bryan et al. (2020a) use pre-upgrade Keck/NIRSPEC to measure v sin i = 8 . 4 +2 . 1 -1 . 4 km s -1 for ROXs 12 b.", '2.3.3. DH Tau': "DH Tau is a member of Taurus. Given the large age scatter in Taurus (e.g. Luhman 2023), we perform isochrone fits to derive 0 . 7 +0 . 3 -0 . 1 Myr. The companion to DH Tau was discovered by Itoh et al. (2005) and Luhman et al. (2006), and orbits at a projected separation of ≈ 2.3 '' or ≈ 310 au from DH Tau, which is part of an ultra-wide binary system (2210 au) with DI Tau (Kraus & Hillenbrand 2009). Bonnefoy et al. (2014) determine a spectral type of M9 -9 . 5 for DH Tau b. Based on the companion's L bol from Luhman et al. (2006), we estimate m =12 ± 4 M Jup and r =2.6 ± 0 . 6 R Jup for the companion. Xuan et al. (2020) measured v sin i = 9 . 6 ± 0 . 7 km s -1 for DH Tau b using pre-upgrade Keck/NIRSPEC data and detected CO and H 2 O in its spectrum. DH Tau b is likely accreting via a circumsubstellar disk, as evidenced by the presence of the Hα line, excess optical continuum emission (Zhou et al. 2014), the Pa β line (Bonnefoy et al. 2014), detection of linear  \npolarization (Holstein et al. 2021), and mid-IR excess emission seen in Spitzer (Martinez & Kraus 2021).", '2.3.4. GQ Lup': "GQ Lup is an on-cloud member of Lupus 1. Galli et al. (2020) estimate an age of 1 . 2 -1 . 8 Myr for Lupus 1. For consistency, we carry out isochrone fits and obtain 2 . 5 +1 . 5 -0 . 9 Myr. The companion GQ Lup b was discovered by Neuhauser et al. (2005) and has a projected separation of ≈ 0 . 7 '' or ≈ 110 au. Alcal'a et al. (2020) found a wide ∼ 0 . 15 M ⊙ component at 2400 au, which they conclude to be most likely bound to GQ Lup A, making this a likely triple system. Like DH Tau b, GQ Lup b likely hosts a circumsubstellar disk, as indicated by Hα and Pa β emission lines and an elevated optical continuum (Seifahrt et al. 2007; Zhou et al. 2014; Demars et al. 2023). Stolker et al. (2021) fitted 0.65 µ m spectro-photometry of the companion and determined a spectral type of M9. They also found excess emission at 4-5 µ m that can be explained from a blackbody with T ≈ 460 K , which they attribute to a disk around GQ Lup b. Demars et al. (2023) used VLT/SINFONI medium-resolution data and grid models based on ATMO to estimate C/O and metallicity for GQ Lup b. Their values are broadly consistent with a solar composition, but discrepant between different observing epochs at the ∼ 60% level in C/O and > 0 . 4 dex in metallicity. From VLT/CRIRES spectroscopy of GQ Lup b, Schwarz et al. (2016) measured v sin i = 5 . 3 +0 . 9 -1 . 0 km s -1 and made detections of CO and H 2 O in the companion's atmosphere. Based on the companion's L bol from Stolker et al. (2021), we estimate m = 33 ± 10 M Jup and r = 3 . 7 ± 0 . 7 R Jup for GQ Lup b.", '2.3.5. GSC 06214-00210': "GSC 06214-00210 (hereafter GSC 6214-210) is an Upper Sco member according to Miret-Roig et al. (2022), who determine the star to be in the slightly older ' π Sco' sub-group. Our isochrone fits yield 22 . 2 +10 . 7 -8 . 0 Myr, consistent with the results from Pearce et al. (2019) who found 16 . 9 +1 . 9 -2 . 9 Myr. The companion was discovered by Ireland et al. (2010) and is separated by 2 . 2 '' on the sky or ≈ 240 au. Bowler et al. (2014b) determine a spectral type of M9 . 5 ± 1, and Bowler et al. (2011) detect Pa β line emission, indicating GSC 6214-210 b possesses a circumsubstellar disk. Demars et al. (2023) used VLT/SINFONI data and ATMO grid models to estimate C/O and [M/H] for GSC 6214-210 b. Their values are broadly consistent with solar, but discrepant between different epochs at the ∼ 70% level in C/O and ∼ 0 . 3 dex in metallicity. Bryan et al. (2018) use pre-upgrade Keck/NIRSPEC to measure v sin i = 6 . 1 +4 . 9 -3 . 8 km s -1 for this companion. Using the companion's L bol from  \nPearce et al. (2019), we estimate m = 21 ± 6 M Jup and r = 1 . 55 ± 0 . 25 R Jup for GSC 6214-210 b.", '2.3.6. 2MASS J01225093-2439505': "2MASS J01225093-2439505 (hereafter 2M0122) is a member of AB Dor (Malo et al. 2013). Our isochrone fitting yields 144 +105 -82 Myr, where the large error bars are due to the grid spacing. This age is consistent with the AB Dor age of 149 +51 -19 Myr from Bell et al. (2015). The companion 2M0122 b was detected by Bowler et al. (2013) at a projected separation of ≈ 1 . 4 '' or ≈ 50 au, who determine a spectral type of L 5 ± 1. Using the companion's L bol from Hinkley et al. (2015), we estimate m = 25 ± 12 M Jup and r = 1 . 2 ± 0 . 2 R Jup for 2M0122 b. 4 Bryan et al. (2020b) use pre-upgrade Keck/NIRSPEC to measure v sin i = 13 . 4 +1 . 4 -1 . 2 km s -1 for 2M0122 b, which enabled a measurement of the companion's obliquity when combined with its photometric rotation period of 6 . 0 +2 . 6 -1 . 0 hr from Hubble Space Telescope (Zhou et al. 2019).", '2.3.7. κ And': "κ And is a probable member of Columba (Zuckerman et al. 2011), and has a range of previous age measurements as summarized by Hoch et al. (2020). Most recently, isochrone fitting from Jones et al. (2016) aided by an interferometric radius measurement of the star yields an age of 47 +27 -40 Myr, broadly consistent with a Columba age. The BHAC15 grid does not go to high enough T eff for κ And ( T eff ≈ 11000 K), so we adopt a uniform age prior between 5-100 Myr for this star based on the Jones et al. (2016) result to estimate the companion mass and radius. The companion κ Andb was discovered by Carson et al. (2013), has spectral type of L 0 -1 Currie et al. (2018), and was at a projected separation of ≈ 0 . 8 '' at time of KPIC observations. Hoch et al. (2020) present a spectral analysis with R ∼ 4000 K band Keck/OSIRIS spectra and report C/O= 0 . 70 +0 . 09 -0 . 24 and a metallicity of -0 . 2 -0 . 0. They also carry out orbit fits to estimate a = 80 +50 -20 au. With the companion's L bol from Currie et al. (2018), we estimate m = 22 ± 9 M Jup and r = 1 . 35 ± 0 . 25 R Jup .", '2.3.8. HIP 79098': "HIP 79098 is likely an unresolved binary in Upper Sco (Janson et al. 2019; Luhman & Esplin 2020). MiretRoig et al. (2022) find HIP 79098 to be in the ' σ Sco'  \nsub-group of Upper Sco. Because the binary properties are unknown and the star exceeds the BHAC15 grid limits, we adopt the Upper Sco age of 10 ± 3 Myr for this system from Pecaut & Mamajek (2016). The companion HIP 79098 b was discovered by Janson et al. (2019) at a projected separation of ≈ 2 . 4 '' or ≈ 360 au. This is the least studied companion in our sample, with only J , H , and K band photometry. Unlike for previous companions, we use HIP 79098 b's absolute K magnitude to estimate its L bol with the empirical K -L bol relation for young brown dwarfs from Sanghi et al. (2023). From this, we estimate m = 28 ± 13 M Jup and r = 2 . 6 ± 0 . 6 R Jup for the companion.", '3.1. KPIC Observations': "We observed the companions in this study using the upgraded Keck/NIRSPEC (Martin et al. 2018; L'opez et al. 2020). The data were collected using both the first version of the KPIC fiber injection unit (FIU) (20192021; Delorme et al. 2021), and the upgraded phase 2 system (2022-2023; Echeverri et al. 2022). The FIU is located downstream of the Keck II adaptive optics system and is used to inject light from a selected target into one of the single-mode fibers connected to NIRSPEC. For all targets, we obtained R ∼ 35 , 000 spectra in K band, which is broken up into nine echelle orders from 1.942.49 µ m. The observing strategy is similar to that of Wang et al. (2021c), although in some datasets we 'nod' between two fibers to enable background subtraction between adjacent frames. The relative astrometry of each companion was computed using whereistheplanet.com (Wang et al. 2021a), based on literature orbital solutions and unpublished data for κ Andb (J. Wang, private communication). For calibration purposes, we also acquire spectra of the host stars before observing the companions, and spectra of a nearby telluric standard star (A0 or B9 spectral type) at similar airmass as the science target. The standard star is observed right before or after the associated science observations. Table 2 summarizes the observations reported in this paper.", '3.2. Data Reduction': 'We only briefly summarize the data reduction procedure in this paper and refer to Xuan et al. (2024) for additional details. For datasets using a single science fiber, we remove the thermal background from the raw images using combined instrument background frames taken before or after the night of observation. For datasets where we perform fiber-nodding, we apply nod-subtraction between adjacent frames, as the spectral traces of each fiber land on a different location in the detector. We  \nalso remove persistent bad pixels identified from the background frames. Then, we use data from a telluric standard star to fit the trace of each column in the four fibers and nine spectral orders, which gives us the position and standard deviation of the point spread function (PSF) in the spatial direction at each column. The trace positions and widths are smoothed using a cubic spline to mitigate random noise. We adopt the trace locations and widths as the line spread function (LSF) positions and widths in the spectral dispersion dimension.  \nFor every frame, we then extracted the 1D spectra in each column of each order. To remove residual background light, we subtracted the median of pixels that are at least 5 pixels away from every pixel in each column. Finally, we used optimal extraction to sum the fluxes using weights defined by the 1D Gaussian LSF profiles calculated from spectra of the telluric star.  \nFor our analysis, we use three spectral orders from 2.29-2.49 µ m, which contain strong CO and H 2 O absorption lines from the companions. The three spectral orders have gaps in between, and cover wavelengths of 2 . 29 -2 . 34 µ m (order 33), 2 . 36 -2 . 41 µ m (order 32), and 2 . 44 -2 . 49 µ m (order 31), respectively.', '4.1. Forward Model of KPIC High-resolution Spectra': "Our forward model for KPIC data follows the framework of previous KPIC papers (e.g. Wang et al. 2021c; Xuan et al. 2024). Here, we present a brief summary. First, we generate atmospheric templates with petitRADTRANS (Molli'ere et al. 2019, 2020), which are shifted in RV and rotationally broadened using the function from Carvalho & Johns-Krull (2023). Then, we convolve the RV-shifted and rotationally-broadened templates with the instrumental LSF determined from spectral trace widths in the spatial direction. 5  \nNext, the atmospheric template is multiplied by the telluric and instrumental response ( T in Eq. 1), which is determined by dividing the standard star spectra by a PHOENIX-ACES model (Husser et al. 2013) matching the standard star's T eff and log g . Since our standard stars have A0 or B9 spectral types, they have nearly no spectral lines in the wavelength region we use in our analysis (2.29-2.49 µ m), mitigating errors due to an imperfect stellar spectrum. For six of our companions, which have projected separations ∼ 1.2-2.5 '' and generally low companion-star contrasts, we find that the  \nTable 1. System Properties for Young Companion Sample  \nReferences -(1) Itoh et al. (2005), (2) Bowler et al. (2013), (3) Kraus et al. (2013), (4) Ireland et al. (2010), (5) Janson et al. (2019), (6) Neuhauser et al. (2005), (7) Carson et al. (2013), (8) Currie et al. (2018), (9) Stolker et al. (2021), (10) Pearce et al. (2019), (11) Bowler et al. (2014a), (12) Hinkley et al. (2015), (13) Bowler et al. (2017), (14) Donati et al. (2012), (15) Luhman et al. (2006), (16) Jones et al. (2016), (17) Yu et al. (2023), (18) Ratzka et al. (2005), (19) Sebastian et al. (2021), (20) Pecaut & Mamajek (2016), (21) Bowler et al. (2014b), (22) Bowler et al. (2011), (23) Gaia Collaboration (2022)  \nspeckle intensity is negligible from preliminary analysis; when allowing for a speckle contribution the fit quality does not improve. The two exceptions are GQ Lup b and κ Andb, which are separated by < 0.8 '' from their host stars (see Table 2). For these two datasets, we account for the significant speckle flux in the companion spectra using the on-axis observations of their host stars, taken immediately before the companion exposures.  \nThe last step is to flux-normalize the companion and/or stellar models and multiply them by flux scale factors, which are in units of NIRSPEC counts. After scaling, the companion and speckle models are added in the case of GQ Lup b and κ Andb, while for the other companions, we only consider the companion flux. To summarize, the forward model is:  \nFM b = α b TM b + α s D s (1)  \nwhere FM b denotes the forward model fitted to the data, α b and α s are the flux scaling factors of the companion and speckle, T is the combined telluric and instrumental response, M b is the companion template from petitRADTRANS , and D s is the observed KPIC spectra of the host star, which already has T factored in. Note that for the six other companions besides GQ Lup b and κ Andb, α s is taken to be zero.  \nLastly, to remove the smoothly varying continuum in the KPIC spectra, we apply high-pass filtering with a median filter of 100 pixels ( ∼ 0 . 002 µ m) on the data and forward model (FM b ) before computing the residuals. The choice of 100 pixels was determined from a series of injection-recovery tests by Xuan et al. (2022) as the optimal size for accurately retrieving molecular abundances in KPIC data. In Appendix B, we also show results from an alternative continuum treatment with a spline model (Ruffio et al. 2023a; Agrawal et al. 2023).", '4.2. Preliminary Analysis for Molecular Detection': "To confirm detection of the companion signal in our data, we fit the KPIC spectra of each companion using atmospheric models from the cloudless Sonora grid (Marley et al. 2021). We select Sonora models with T eff and log g that best match each companion's bulk properties, as estimated from the evolutionary models in § 2.2. Using the forward model framework described above, we estimate the maximum likelihood value for both the companion flux and speckle flux in the data as a function of RV shift, following Ruffio (2019) and Wang et al. (2021c). This forward model framework allows us to estimate the companion flux (in data counts) as a function of RV, which can be interpreted as a crosscorrelation function (CCF). In this paper, we refer to  \nTable 2. KPIC observations presented in this work. The throughput is end-to-end throughput measured from top of the atmosphere, and varies with wavelength due to differential atmospheric refraction and the instrumental blaze function. We report the 95% percentile throughput over the K band, averaged over all frames. We also report the median spectral SNR per pixel from 2.29-2.49 µ m.  \nthis as the CCF for simplicity, but note that it is not the same as the traditional CCF (e.g. as in eq.1 of Xuan et al. 2020). To estimate the CCF S/N, we perform the same fitting procedure using a spectral trace that contains only background flux, and take the standard deviation of this background CCF as noise. We calculate CCFs for templates with CO, H 2 O, and CO + H 2 O, as shown in Fig. 3. The companions are detected with CCF S/N between 7 to 80 when using the combined CO+H 2 O template, and both CO and H 2 O are individually detected with S / N > 3 for all companions. By detecting the major carbon and oxygen-bearing species in their atmospheres, we can constrain their atmospheric compositions with retrievals.", '4.3. Atmospheric Retrieval Setup and Inputs': 'We use the radiative-transfer code petitRADTRANS to generate synthetic companion templates for use in atmospheric retrievals. These synthetic templates represent M b in Eq. 1 of the forward model. We use the lineby-line opacity sampling mode, and down-sample the native R = 10 6 opacities by a factor of 3 to speed up the retrievals. In our retrievals, we fit for the chemical abundances ( § 4.3.2), cloud structure ( § 4.3.3), and temperature profile ( § 4.3.4). We impose mass and radius priors motivated by evolutionary models in the retrievals ( § 4.3.5). Other parameters such as RV, v sin i , and flux scales ( α b , α s ) are also fitted for in our forward model. As an example, we summarize the fitted parameters in Table 3 for κ Andb.', '4.3.1. Opacities': 'For the hottest companions in our sample, we find that there is contribution to the emission spectrum from  \nregimes with T > 3000 K (see Fig. 4), which exceeds T max = 3000 K of default petitRADTRANS opacity tables. Therefore, we adopt the opacities generated by Xuan et al. (2024), which go up to 4500 K . We include the line opacities of H 16 2 O (Polyansky et al. 2018), C 16 O, 13 CO (Rothman et al. 2010), OH (Brooke et al. 2016), FeH (Dulick et al. 2003; Bernath 2020), TiO (McKemmish et al. 2019), AlH (Yurchenko et al. 2018), VO (McKemmish et al. 2016), H 2 S (Azzam et al. 2016), and CH 4 (Hargreaves et al. 2020). In addition, we include atomic opacities from Na, K, Mg, Ca, Ti, Fe, Al, and Si (Kurucz 2011). For continuum opacities, we include the collision induced absorption (CIA) from H 2 -H 2 and H 2 -He, as well as the H- bound-free and free-free opacity.', '4.3.2. Chemistry': "We parameterize the chemical abundances with C/O and [C/H], where [C/H] is equivalent to the bulk metallicity. In other words, we assume [C/H]=[Fe/H]=[N/H] and so on, 6 while the oxygen abundance is determined by [C/H] and C/O. We denote the bulk metallicity as [C/H] since we are only sensitive to C-, O-bearing species (i.e. CO, H 2 O) in the companion atmospheres.  \nAs in Xuan et al. (2022), we use an equilibrium chemistry grid in our retrievals computed with easyCHEM , a  \nFigure 3. Cross-correlation functions (CCFs) with CO (red, dashdot lines), H 2 O (blue, dotted lines), and CO+H 2 O (purple, solid lines) molecular templates from Marley et al. (2021). Each panel is for a different companion, but shares the same legend. These CCFs are computed using 3 spectral orders from 2 . 29 -2 . 49 µ m. The gray lines are CCFs of the CO+H 2 O templates with background flux in the slit, and we use the standard deviations of the gray lines to estimate the CCF noise. The effect of rotational broadening is clearly visible in the κ Andb CCF.  \n<!-- image -->  \nGibbs free energy minimizer described in Molli'ere et al. (2017). In this paper, we update the chemical grid to use the updated solar elemental abundances from Asplund et al. (2021). Our chemical grid stores the massmixing ratios of both gas-phase species and condensates (which we use for the cloud model in the next section). We tested the option of including a quench pressure ( P quench ) to allow for carbon disequilibrium chemistry, which fixes the abundances of H 2 O, CO, and CH 4 where P < P quench using the equilibrium values found at P quench (Zahnle & Marley 2014). However, most of our companions are too hot ( T eff ≳ 1800 K ) for CH 4 to be detectable, and simultaneous detection of CO and CH 4 is necessary to constrain carbon quenching (Xuan et al. 2022). Since we find that P quench is unconstrained for all companions from preliminary tests, we do not include it in the reported retrievals.", '4.3.3. Clouds': "Condensate cloud opacity is expected to gradually decrease with increasing temperature from L to M spectral types, as important cloud particles such as MgSiO 3 , Fe, and Al 2 O 3 start to evaporate between ∼ 1600-1900 K at a pressure of 0.1 bars (Molli'ere et al. 2022). Since our companions have spectral types ranging from L to early M spectral types, we consider both clear and cloudy models in order to explore the sensitivity of our retrieved abundances to the assumed cloud models.  \nFor the cloudy models, we use 1) a simple gray opacity model which adds a constant opacity to the atmosphere, and the 2) EddySed model (Ackerman & Marley 2001) as implemented in petitRADTRANS , which includes the effect of scattering from clouds (Molli'ere et al. 2020). We now describe fitted parameters of the Ed-  \nTable 3. Fitted Parameters and Priors for κ Andb Retrievals  \nNote -U stands for a uniform distribution, with two numbers representing the lower and upper boundaries. N stands for a Gaussian distribution, with numbers representing the mean and standard deviation. The P -T parameters are described in § 4.3.4 and the cloud parameters are described in § 4.3.3.  \ndySed model. First, log( ˜ X Mg S iO 3 ) is the scaling factor for the cloud mass fraction, so that log( ˜ X Mg S iO 3 ) = 0 is equal to the equilibrium mass fraction. The equilibrium mass fraction is determined by the chemical grid (see § 4.3.2). For each cloud condensate, this scaling factor along with f sed , K zz , and σ g set the cloud mass fraction as a function of pressure and the cloud particle sizes (for details, see Ackerman & Marley 2001; Molli'ere et al. 2020). Here, f sed is the sedimentation efficiency, K zz is the eddy diffusion coefficient, and σ g is the width of the lognormal cloud particle size distribution. Following Zhang et al. (2023), we fit a different f sed for each cloud species but a global K zz and σ g . Therefore, when including two different cloud species, there are a total of six cloud parameters.  \nGiven the range in T eff of our objects, we consider models with MgSiO 3 + Fe for the colder objects ( T eff ≲ 2000 K) and Fe + Al 2 O 3 for the hotter objects. Since the cross sections of these cloud species have similar slopes over the small wavelength range (2 . 29 -2 . 49 µ m) we are modeling (Wakeford & Sing 2015), our choice of the cloud species primarily serves to set the cloud base locations in the EddySed model. For each companion, we choose two cloud species which intersect their P -T profile at the deepest pressures (i.e.  \ncloser to the photosphere), as these clouds would more meaningfully impact the emission spectra.", '4.3.4. Temperature Structure': "We adopt the pressure-temperature ( P -T ) parameterization from Xuan et al. (2024), which is motivated by Piette & Madhusudhan (2020). Our profile is parameterized by seven ∆ T/ ∆ P values between eight pressure points and the temperature at one of these pressures, T anchor . Because the photosphere for each companion is located at slightly different locations, we manually customize the pressure extent for each companion's retrieval to optimally encompass the companion's emission contribution. Specifically, we set at least four pressure points in the region where 90% of the flux originates, which we determine by computing the wavelength-weighted emission contribution function. We choose the other points to be approximately equally spaced in log pressure. The selected pressure points are labeled in Fig. 4 for a few companions, and listed in Table 3 and Appendix C for all companions. For the radiative transfer, the eight P -T points from our profile are interpolated onto a finer grid of 100 P -T points using a monotonic cubic interpolation as recommended by Piette & Madhusudhan (2020). We do not  \napply smoothing to our profiles as Rowland et al. (2023) showed that smoothing can bias retrieval results.  \nUsing this P -T prescription, we are able to set meaningful priors on the fitted ∆ T/ ∆ P values by considering atmospheric profiles from self-consistent models. After defining the pressure extent and pressure points used for a given companion, we fit SPHINX (Iyer et al. 2023) and Sonora (Marley et al. 2021) profiles with T eff and log g similar to the companion's T eff and log g using our P -T function. The expected T eff and log g for each companion are determined from evolutionary models ( § 2.2). We save the best-fit P -T values of each self-consistent profile, noting that the ∆ T/ ∆ P values are very similar between profiles with different T eff and log g as the slope is set by radiative-convective equilibrium. In the retrieval, we then set uniform P -T priors that bracket the mean of the best-fit values with a wide prior range of ± 50 -100% of the mean value. While the self-consistent models we considered are cloudless, the wide uniform priors we impose allow sufficient flexibility in the profile shapes, as we discuss in § 5.2. We note that our method of setting physically-motivated priors on the P -T profile is similar in spirit to the approach from Zhang et al. (2023), who imposed Gaussian priors on the temperature gradient d ln T/d ln P informed by self-consistent models in their retrievals. Both approaches leverage information from self-consistent P -T profiles to prevent the retrievals from returning overly isothermal or unphysical profiles.", '4.3.5. Mass and Radius Priors from Evolutionary Models': "The KPIC spectra are not flux-calibrated and cover a very small wavelength range. Hence, these data provide little information on a companion's radius and mass. We find that preliminary retrievals occasionally yield unphysical radii and log g for our companions (e.g. log g > 5 for a young, low-gravity companion). Therefore, we adopt mass and radius priors from evolutionary models in our retrievals. These priors are described in § 2.2 and listed in Table 1.", '4.3.6. Model Fitting with Nested Sampling': 'We use nested sampling as implemented by dynesty (Speagle 2020) to find the posterior distributions for the model parameters. Specifically, we use between 500-800 live points and adopt the stopping criterion that the estimated contribution of the remaining prior volume to the total evidence is less than 1%.  \nOne advantage of adopting nested sampling is that we can use the Bayesian evidence from each fit to calculate the Bayes factor B , which quantifies the relative probability of model M 2 compared to M 1 . We use the Bayes  \nfactor to compare different models throughout this paper to determine whether a given M 2 is justified over M 1 .', '5. ATMOSPHERIC RETRIEVAL RESULTS': 'We summarize the retrieval results in this section. Given the large number of objects, we focus on representative examples here. In Fig. 5, we plot the data and best-fit model for GQ Lup b, our highest S/N dataset. Plots for other companions are shown in Appendix D. Key parameters from the retrievals for all companions are listed in Table 4. In § 5.1, we show the insensitivity of our results to the mass and radius priors used in the retrievals. We discuss the retrieved P -T profiles and compare the retrieved T eff with predicted T eff from evolutionary models in § 5.2. In § 5.3, we discuss whether cloudy models are preferred over clear models for our companions, and whether clouds impact the retrieved abundances. Finally, we discuss our C/O, [C/H] measurements in § 5.4 and 12 C / 13 C constraints for three companions in § 5.5.', '5.1. Effect of Mass and Radius Priors': 'We adopt mass and radius priors in our retrievals to prevent the fits from yielding unphysical radii and log g . In principle, jointly fitting photometry or lowresolution data with the KPIC HRS could provide the flux information needed to anchor the L bol and radius. For example, Stolker et al. (2021) fit the 0.65 µ m spectro-photometry for GQ Lup b accounting for dust extinction and emission from a circumplanetary disk and find 3 . 77 ± 0 . 10 R Jup , which is consistent with our priors derived from evolutionary models. However, as Stolker et al. (2021) illustrates, for broadband spectro-photometry it is necessary to model both extinction from dust and potential circumplanetary disk emission for the youngest companions. Several companions in our study likely possess circumplanetary disks as well (e.g. Bowler et al. 2011; Holstein et al. 2021; Stolker et al. 2021). Doing such modeling for all the companions in our sample is beyond the scope of this paper, which is focused on measuring atmospheric abundances from high-resolution spectra. Hence, we choose to adopt mass and radius priors to achieve the same purpose. We verified that these mass and radius priors do not affect the retrieved abundances by comparing against retrievals without the mass and radius priors, and found that the measured abundances vary by much less than 1 σ when the priors are imposed.', '5.2. P -T Profile and Effective Temperature': "Our retrieved P -T profiles mostly follow selfconsistent (cloudless) profiles in the deeper atmosphere  \nTable 4. Results of Spectral Retrievals for Eight Substellar Companions  \nNote -Selected atmospheric parameters and their central 68% credible interval with equal probability above and below the median are listed. These values only account for statistical error. For clarity, we list the metallicity as C/H ( × solar) instead of the fitted [C/H]. The retr. T eff is the retrieved T eff , while the evol. T eff refers to the expected T eff derived from a set of evolutionary models ( § 2.2). The rightmost column lists the log Bayes factor ln( B ) for each retrieval. We compute ln( B ) with respect to the clear model for each companion, i.e. the model with ln( B )=0. For each companion, the adopted model is in bold, which is by default the clear model unless one of the cloudy models is preferred over the clear model by > 3 σ , which corresponds to > 3 in ln( B ). Negative values of ln( B ) for an alternative model means the clear model is preferred over this model. For GQ Lup b, HIP 79098 b, and DH Tau b, we ran retrievals without 13 CO opacities to quantify the detection significance of 13 CO. The companion RVs from our retrievals are provided in Appendix F.  \na Based on Trotta (2008), ln( B ) of 3 and 11 correspond to 3 σ and 5 σ significance, respectively.  \n(see Fig. 4 and Appendix E), but they can be hotter and more isothermal than the corresponding self-consistent models above the cloud bases. This could be due to several reasons. First, the cloudless model P -T profiles do not provide the best point of comparison for cloudy atmospheres, and cloud formation is expected to heat up the atmosphere above the cloud base resulting in a detached radiative zone (Tsuji 2002; Burrows et al. 2006). This effect can be seen as a kink the P -  \nT profile, which is most obvious for κ Andb at pressures lower than ∼ 1 bar (see Fig. 4). Upcoming selfconsistent models that include cloud radiative feedback (Morley et al. 2024) would provide a better comparison. Second, our retrievals could be showing evidence of the well-known cloud and isothermal P -T correlation (e.g. Burningham et al. 2017; Molli'ere et al. 2020), whereby a cloudier atmosphere can be re-produced by one with a more isothermal temperature gradient. While we at-  \nFigure 4. Retrieved P -T profiles and emission contribution functions for three example companions (one per row). Left panels: We plot random draws of the retrieved P -T profiles in blue. The gray lines show Sonora cloudless models (Marley et al. 2021) with similar bulk properties as the companions. Different cloud condensation curves are plotted as colored dashed lines. The horizontal blue lines mark pressure points between which we fit ∆ T values in our P -T parameterization. Right panels: The emission contribution functions of the best-fit models. For a given wavelength, darker colors mean that a larger fraction of emission originates from that pressure level. The darker, bar-like structures seen for κ Andb and GSC 6214-210 b coincide with the cloud base locations. In our best-fit EddySed models, we find that the total cloud optical depth ( τ ) is about 1.4/0.6 for κ Andb/GSC 6214-210 b at their retrieved cloud bases, so a non-negligible part of the flux beneath the clouds can propagate through the cloud base. These two companions prefer the EddySed cloud model over the clear model to > 3 σ significance (see § 5.3). Plots for the other companions are included in Appendix E.  \n<!-- image -->  \npt to address this behavior with our informed P -T  \nprofile priors, our priors are uniform and wide enough  \nFigure 5. KPIC data of for GQ Lup b from one science fiber are shown in black, with error bars in gray. This represents half the data, since we nodded between two fibers. Each spectral order used in the retrieval is broken into two panels. The full model is shown in red dashed lines and consists of the RV-shifted and broadened companion model in blue ( M b in Eq. 1), the stellar model in purple ( M s ) to model the speckle contribution, and the telluric and instrumental response ( T ). The residuals are shown as gray points. For clarity, an offset of +50 counts was added to the companion model. Plots for the other companions can be found in Appendix D.  \n<!-- image -->  \nthat isothermal behavior is still allowed. A stricter, Gaussian prior such as implemented by Zhang et al. (2023) could prevent this behavior at the expense of limiting the parameter space for the retrieval to explore.  \nTo assess whether our priors on the P -T profile affect the retrieved parameters, we repeat the baseline HIP 79098 b retrieval with wide priors from 0 -2000 K for each ∆ T value. This second retrieval yielded the  \nsame median value and error bars for the abundances, RV, and v sin i as our baseline retrieval, indicating that our default P -T priors are conservative enough and not biasing the results.  \nFrom our retrievals, we compute T eff by sampling the posteriors to generate low-resolution models from 0.1 to 30 µ m, calculating the integrated flux, and applying the Stefan-Boltzmann law with the retrieved radii.  \nTo see how our results compare with evolutionary models, we compare our retrieved T eff posteriors with the evolutionary-model predicted T eff that were estimated in § 2.2. As with the mass and radius priors, we estimate a range of evolutionary T eff that encompass all the models, which are listed in Table 4.  \nFor the two companions that showed > 3 σ preference for the EddySed cloud model, GSC 6214-210 b and κ Andb, our retrieved T eff from the EddySed model is lower compared to those from the clear and gray opacity models (see Fig. 6). In addition, the AMES-Dusty and SM08 models which include clouds generally predict lower T eff for the same object compared to the cloudless ATMO 2020 and AMES-COND models (Appendix A). The EddySed-retrieved T eff for these two companions are closer to the predictions of the cloudy evolutionary models. For example, the EddySed-retrieved T eff is 1680 +60 -50 K for κ Andb. For this companion, the cloudy SM08 evolutionary model predicts T eff ≈ 1760 K while the cloudless ATMO 2020 evolutionary model predicts T eff ≈ 1860 K (see Fig. A1). Thus, the overlap between the EddySed-retrieved T eff and predicted T eff from cloudy evolutionary models for GSC 6214-210 b and κ Andbis consistent with the fact that we find evidence of clouds in these objects (discussed further in § 6.3).  \nNext, we discuss the six companions that did not show strong preference for cloudy models; these companions have similar retrieved T eff posteriors between different cloud models. First, for ROXs 12 b, ROXs 42 Bb, and 2M0122 b, we find good agreement ( < 1 σ ) between their retrieved and evolutionary T eff . For the remaining three, GQ Lup b, DH Tau b, and HIP 79098 b, which are among our hottest companions with late M spectral types, the retrieved T eff from each cloud model is lower by ∼ 300 -400 K compared to the evolutionary model predictions (a ∼ 2 σ discrepancy). We find that the retrieved log( L bol /L ⊙ ) of these companions are also slightly lower than those predicted from the models by ∼ 0 . 1 -0 . 2 dex (or by 1 -2 σ ). On the other hand, the median-retrieved radius is slightly higher by ∼ 5 -10% compared to the median of the radius prior used in the retrieval.  \nXuan et al. (2024) also noted T eff and radius discrepancies from retrievals of a M7.5 stellar companion with similar T eff , where the T eff /radius was slightly lower/higher than evolutionary models. One potential explanation is that high-resolution spectra after highpass filtering are more sensitive to the slope of the P -T profile than the absolute temperature value (Landman et al. 2023), and therefore may not provide accurate T eff . In our case, the retrieved T eff is also influenced by the  \nmass and radius priors we placed. However, this alone does not explain why only these companions show disagreement with the evolutionary models. Another reason for the T eff discrepancy could be shortcomings in the treatment of clouds in evolutionary models, especially at the M to L transition. As shown by a systematic study of brown dwarfs in Sanghi et al. (2023), clouds may be highly inaccurate at the M to L transition. This means that the evolutionary T eff are not necessarily correct for late M objects.  \nFor purposes of this paper, we simply check that our retrieved abundances are not affected by the retrieved T eff by repeating retrievals with the P -T profile fixed to a self-consistent profile that has higher T eff . Specifically, for HIP 79098 b, we tried fixing the P -T profile to a SPHINX model profile (Iyer et al. 2023) with T eff = 2600 K and log g = 4 . 0 (overplotted in Fig. 4) and repeated the retrieval. We find that all retrieved parameters shift by < 1 σ compared to the baseline retrieval. Therefore, we conclude that discrepancies in T eff and radius have negligible impact on the retrieved abundances from our high-resolution data.", '5.3. Impact of Clouds': "For most of our companions, the posteriors from our clear and cloudy retrievals overlap significantly; the median or best-fit retrieved abundance parameters are identical within ≪ 1 σ . We compare the cloudy models and clear models statistically with the log Bayes factor, or ln( B ), which are listed in Table 4. The cloudy models are compared against the clear model, so a positive ln( B ) for a cloudy model means it is preferred over the clear model. While all companions prefer cloudy models, the preference is only significant for κ Andb and GSC 6214-210 b, which show 3 . 7 σ and 3 . 4 σ preferences for the EddySed cloud model. For these two companions, our median retrieved abundances from the EddySed models can differ by ≈ 1 σ compared to the clear model, as we show in Fig. 6. For all other companions, the cloudy models are preferred by 0 -2 . 6 σ , and the addition of clouds negligibly impact their retrieved abundances and other parameters.  \nWe illustrate the impact of clouds on the retrieved spectra of κ Andb and GSC 6214-210 b in Fig 6, which shows that both the continuum and line depths are varying in a distinctly non-gray behavior. This explains why the spectra could not be fit as well with the gray opacity model. As the data have relatively low S/N ( ∼ 5 and ∼ 2 per pixel for κ Andb and GSC 6214-210 b, respectively), the residual noise makes it difficult to clearly see the companion lines. Thus, we try to visualize the effect of clouds in the data after binning down the companion  \n+  \n0  \n.  \n93  \nFigure 6. Top left: The best-fit Eddysed cloud model (blue) and clear model (red dashdot) for κ Andb. These models are normalized companion templates from petitRADTRANS , and show non-gray variations in line depths and continuum location. We also overplot the normalized companion spectra (data - best-fit speckle model), which has been binned down to R = 5000 to better visualize the companion lines. Top right: the joint posterior distributions of v sin i , C/O, [C/H], and T eff from the different cloud models. There are ≈ 1 σ variations in the abundances and v sin i , while the T eff is ≈ 350 K lower for the EddySed model compared to the clear model. Bottom panel: same but for GSC 6214-210 b. Note the different x-axis scales for the spectral plots, which are chosen for visualization purposes. The retrieved T eff for GSC 6214-210 b is ≈ 500 K lower in the EddySed model compared to the clear model. Since the continuum-to-line contrast increases for lower T eff atmospheres, this could cause absorption lines in the cloudy model to be deeper than those in the clear model despite the added cloud opacity, which is evident in the case of GSC 6214-210 b.  \n<!-- image -->  \nflux contribution of in data. Since the κ Andb spectra have a non-zero speckle component (see Eq. 1), we subtract the best-fit speckle model from its KPIC spectra and then bin them down to R = 5000. As shown in Fig 6, we see that the companion line depths do match better with the cloudy, EddySed model compared to the clear model for κ Andb. However, the same effect is not as clear when we plot the GSC 6214-210 b companion flux, which is more noisy due to its lower S/N. In addi-  \non, the relative difference between the clear and cloudy model is smaller for this object. Thus, we conclude that a clear confirmation of the effect of clouds in highresolution spectra would require higher S/N data. This could be achieved with longer integration times on our targets, or existing high S/N spectra of isolated brown dwarfs.  \nWhy do only κ Andb and GSC 6214-210 b prefer EddySed cloud models in our retrievals? Based on Table 4,  \nFigure 7. Left panel: The CCF between the 13 CO-only template and (data - best-fit reduced model) in blue for GQ Lup b. The reduced model retrieval does not include 13 CO opacities. The CCF between the 13 CO-only template and (data - best-fit full model) is shown in red. The full model retrieval includes 13 CO opacities. The fact that the blue CCF shows a peak at the companion's rest frame (gray solid line) indicates a real 13 CO detection. For comparison, the gray dotted line is the telluric rest frame. In the red CCF, we do not expect a peak since 13 CO is fitted for in this model, so the residuals should be free of 13 CO. Middle and right panels: same but for HIP 79098 b and DH Tau b. The 13 CO detection for DH Tau b is tentative.  \n<!-- image -->  \nwe see that they are among the colder companions. Except for 2M0122 b, all other companions are predicted to be hotter. For κ Andb, we retrieve T eff = 1680 +60 -50 K from the EddySed model and T eff = 2050 +130 -140 K for a clear atmosphere. For GSC 6214-210 b, we retrieve T eff = 1860 +170 -110 K and T eff = 2420 ± 150 K from the EddySed and clear retrievals, respectively. As noted in § 5.2, the EddySed-retrieved T eff for these two companions matches predictions from cloudy substellar evolutionary models (SM08, AMES-dusty) better than the T eff retrieved from the clear model, which are overly high. The fact that we see evidence of clouds for only κ Andb and GSC 6214-210 b suggests that the impact of clouds on the spectra is smaller for hotter companions and therefore harder to detect despite their higher S/N data (e.g. GQ Lup b has S/N per pixel 12 v.s. 5 for κ Andb). When examining the P -T profiles for κ Andb and GSC 6214-210 b, we see that their cloud bases are predicted to be closer to the K band continuum due to their lower T eff (see Fig. 4). In contrast, for hotter companions, the cloud bases are at much lower pressures than the continuum. Physically, this means that there is more condensable cloud material for κ Andb and GSC 6214-210 b at pressures closer to their photospheres, thereby causing them to appear cloudier.  \nTherefore, our results match the trend of decreasing cloud opacity with rising temperature that is observed across the L to M transition, for instance from the decreasing silicate cloud absorption strength at ≈ 8-11 µ m seen in mid-IR Spitzer/IRS spectra of brown dwarfs (e.g. Cushing et al. 2006; Su'arez & Metchev 2022). We note that the coldest of our companions, 2M0122 b, shows a weak preference (2 σ ) for the EddySed model, but its S/N is the lowest among the sample ( ∼ 1 per pixel), which may prevent a stronger constraint of its clouds  \nproperties. Alternatively, viewing geometry also affects the strength of the silicate cloud feature of brown dwarfs, whose equators tend to be cloudier than their poles (Vos et al. 2017; Su'arez et al. 2023). Given the measured companion spin axis inclination of 33 +17 -9 deg (Bryan et al. 2020b), 2M0122 b is close to a pole-on geometry, which may reduce the effect of clouds on its measured emission spectrum relative to companions with more equatorial viewing geometries. Finally, because spectral features become deeper compared to the continuum at lower T eff for brown dwarfs (e.g. Ruffio et al. 2023b), this causes CO lines to be deeper in the EddySed model for GSC 6214-210 b, despite the impact of clouds (Fig. 6). This suggests a potential trade-off between T eff (set by the PT profile) and clouds. Future work using broader wavelength coverage or flux-calibrated data should be able to constrain T eff better in order to confirm our findings. We further discuss the effect of clouds on high-resolution spectroscopy in § 6.3.", '5.4. Carbon and Oxygen Abundances': 'Our retrieved C/O and [C/H] values for the companions are listed in Table 4. As with previous highresolution studies (e.g. Finnerty et al. 2023; Xuan et al. 2024), our constraints on the relative abundances of different molecular species as indicated by the ratio of their relative line depths, which are used to calculate C/O, are tighter than our constraints on the absolute abundances of individual molecular species, which are used to determine the overall atmospheric [C/H]. Specifically, our uncertainties range between 0 . 02 -0 . 08 for C/O and 0 . 2 -0 . 5 dex for [C/H] depending on the object and S/N. Both the C/O and [C/H] of our companions are consistent at the 1-2 σ level with the solar composition from Asplund et al. (2021). The implications of our mea-  \nured abundances for formation pathways are discussed in § 6.4.', '5.5. 12 CO / 13 CO Measurements': "We obtain bounded constraints on 12 CO / 13 CO for GQ Lup b, HIP 79098 b, and DH Tau b, which are the three companions with the highest S/N detections as shown by the CCFs (Fig. 3). To quantify the detection significance of 13 CO, we run additional retrievals where we leave out 13 CO opacities, and compute the Bayes factor between these 'reduced models' without 13 CO and the 'full model' that includes 13 CO and fits for 12 CO / 13 CO. The resulting log Bayes factors are listed in Table 4, and correspond to 5 . 7 σ , 4 . 3 σ , and 3 . 4 σ detections of the 13 CO isotopologue for GQ Lup b, HIP 79098 b, and DH Tau b, respectively. We note that ROXs 12 b also has a 12 CO / 13 CO posterior which peaks at ∼ 100, but the corresponding detection significance is < 3 σ so we do not consider it in the following discussion. The other four companions show unbounded 12 CO / 13 CO posteriors from our tests, so their baseline models do not include 12 CO / 13 CO.  \nFollowing Xuan et al. (2022), we perform a crosscorrelation analysis to obtain a complementary perspective on the robustness of these 13 CO detections. The goal of this analysis is to assess whether the full models prefer 13 CO independent of the Bayes factor calculation. First, we compute the CCF between a 13 CO-only model and the (data - model without 13 CO). The latter is equivalent to the residuals of the reduced model, and will contain residual 13 CO lines if the data contain 13 CO. Then, we compute the CCF between the 13 COonly model and the (data - model with 13 CO). This second CCF should not show a detection, as 13 CO is already fitted for in the full model. We generate 13 COonly models by manually zeroing the opacities of all other line species except 13 CO when computing the full model. The CCF calculations follow the framework described in § 4.2.  \nFrom the CCF analysis, we find that the CCFs for GQ Lup b and HIP 79098 b show convincing peaks at the companion's rest frame that indicate a strong 13 CO detection (see Fig. 7). However, the 13 CO CCF for DH Tau b is much noisier, showing residual structure comparable to the peak. Therefore, we consider the 13 CO detection in DH Tau b to be tentative, since it is possible that 13 CO is being used to fit for systematics that are unaccounted for in the model. This is consistent with the larger uncertainties in the 12 C / 13 C of DH Tau b, and is expected given the lower S/N of the DH Tau b data compared to the data for the other two  \ncompanions. We further discuss our 12 C / 13 C measurements in the context of previous results in § 6.5.", '6.1. Comparison with Previous Work': "As noted in § 2, many our companions have previous medium- or high-resolution studies which reported on their v sin i and/or atmospheric abundances. Table 5 compares our baseline KPIC measurements with previous results.  \nHoch et al. (2020) used Keck/OSIRIS mediumresolution spectra ( R ∼ 4000, 2 . 22 -2 . 4 µ m) to measure the abundances of κ Andb using a custom grid of PHOENIX models (Barman et al. 2011; Barman et al. 2015) that vary in metallicity and C/O. They report C / O = 0 . 70 +0 . 09 -0 . 24 and [M / H] = 0 . 0 -0 . 2, which are consistent with our measured C / O = 0 . 58 +0 . 05 -0 . 04 and [C / H] = -0 . 12 +0 . 28 -0 . 19 . Our measured v sin i = 39 . 4 ± 0 . 9 km s -1 for κ Andb also agrees with the recent study by Morris et al. 2024, accepted, who find 38 . 42 ± 0 . 05 km s -1 using KPIC data from a different observing night.  \nSchwarz et al. (2016) used VLT/CRIRES spectra ( R ∼ 100 , 000, 2 . 30 -2 . 33 µ m) to measure v sin i = 5 . 3 +0 . 9 -1 . 0 km s -1 for GQ Lup b, which is consistent with our v sin i = 6 . 4 +0 . 3 -0 . 4 km s -1 to within ≈ 1 . 1 σ . This agreement is reassuring given the ∼ 3 × lower spectral resolution of Keck/KPIC compared to VLT/CRIRES, and confirms the ability of KPIC to measure v sin i values below the resolution limit of ∼ 7 . 5 km s -1 , as demonstrated by Xuan et al. (2024) as well.  \nDemars et al. (2023) used VLT/SINFONI mediumresolution spectra and the self-consistent model ATMO (Tremblin et al. 2015, 2016) to fit the abundances of GQ Lup b and GSC 6214-210 b. They adopted C / O = 0 . 44 +0 . 13 -0 . 11 and [M / H] = 0 . 23 ± 0 . 06 for GQ Lup b and C / O = 0 . 48 +0 . 16 -0 . 12 and [M / H] = -0 . 16 ± 0 . 17 for GSC 6214-210 b. Overall, our measured abundances agree with their results at the 1 -2 σ level. However, we note that the Demars et al. (2023) results are discrepant between different observing epochs at the ∼ 60% level in C/O and > 0 . 4 dex in [M/H], and sometimes hit the limits of their grid. In contrast, as we show in Appendix B, our KPIC results between different observing epochs of GQ Lup b are consistent at the < 3% level in C/O and ∼ 0 . 1 dex in [C/H]. This demonstrates the reliability of high-resolution spectra and our forward model + retrieval approach in measuring atmospheric abundances.  \nFive of our targets (ROXs 42 Bb, ROXs 12 b, DH Tau b, 2M0122 b, GSC 6214-210 b) have reported v sin i measurements with pre-upgrade Keck/NIRSPEC spectra ( R ∼ 25 , 000) as part of a spin survey for lowmass substellar companions (Bryan et al. 2018; Xuan  \nTable 5. Comparison of Our Measurements with Previous Work  \nNote -TP refers to this paper.  \nReferences -(1) Schwarz et al. (2016), (2) Bryan et al. (2018), (3) Bryan et al. (2020b), (4) Bryan et al. (2020a), (5) Xuan et al. (2020) (6) Hoch et al. (2020), (7) Demars et al. (2023), (8) Inglis et al. (2023), (9) Morris et al. 2024, accepted  \net al. 2020; Bryan et al. 2020a). We list the literature values along with our new values in Table 5, and discuss them below.  \nFor the faster rotators GSC-6214-610 b and 2M0122 b, our measured v sin i are consistent with the pre-upgrade NIRSPEC values from Bryan et al. (2018) and Bryan et al. (2020b) at the 1 σ and 2 σ levels, respectively. For the slower rotators (ROXs 12 b, ROXs 42 Bb, and DH Tau b), we find ≳ 3 σ discrepancies in v sin i , with a trend that KPIC measured v sin i values are lower than those from pre-upgrade NIRSPEC studies. However, we note that for ROXs 12 b and ROXs 42 Bb our KPIC data only prefer a non-zero spin with ≈ 1 σ significance, so future data is required to definitively measure their v sin i . For ROXs 42 Bb, Inglis et al. (2023) perform atmospheric retrievals with pre-upgrade NIRSPEC data and report v sin i = 10 . 5 ± 0 . 9 km s -1 , which agrees with the v sin i = 9 . 5 +2 . 1 -2 . 3 km s -1 found by Bryan et al. (2018) for the same data. However, we find v sin i = 4 . 2 +1 . 8 -2 . 2 km s -1 for ROXs 42 Bb. Despite these differences in v sin i , our retrieved abundances for ROXs 42 Bb agree with those from Inglis et al. (2023). Finally, in Xuan et al. (2020) we found v sin i = 9 . 6 ± 0 . 7 km s -1 for DH Tau b using pre-upgrade NIRSPEC data, while we measure v sin i = 5 . 7 +0 . 8 -1 . 0 km s -1 with KPIC.  \nLow S/N spectra can preclude confident spin measurements, especially for very slow rotators. However, in all three cases above, the KPIC data have higher CCF S/N compared to the pre-upgrade NIRSPEC data. For ROXs 12 b, ROXs 42 Bb, and DH Tau b, we obtain ∼ 25 σ , ∼ 15 σ , and ∼ 20 σ detections (Fig. 3), respectively, while the previous NIRSPEC detections have CCF S/N ∼ 10 σ (Xuan et al. 2020; Inglis et al. 2023). If the companions are in fact rotating faster, our KPIC data should have been able to reveal this.  \nDespite the higher CCF S/N of the KPIC data, the uncertainties of our new v sin i measurements are comparable with the older measurements. This suggests that there may be additional sources of uncertainty that the older NIRSPEC studies did not consider. For example, the retrieval approach allows us to vary the P-T profile and therefore adjust the pressure-broadened line shapes in each iteration. In contrast, Xuan et al. (2020) and Bryan et al. (2020a) used fixed atmospheric templates to compute v sin i , and may under-estimate the uncertainty introduced by the uncertain atmospheric properties of the companion.  \nWe note that for a high S/N dataset of the isolated brown dwarf 2MASS J03552337+1133437 (hereafter 2M0355), Zhang et al. (2021) found similar discrep-  \nancies in v sin i compared to Bryan et al. (2018). Using the same NIRSPEC data for the brown dwarf, Zhang et al. (2021) found an upper limit of 4 km s -1 for v sin i , consistent with a non-detection of spin, but Bryan et al. (2018) reported v sin i = 14 . 7 ± 1 km s -1 . Using earlier R ∼ 25 , 000 NIRSPEC data, Blake et al. (2010) also found a similarly large spin for this brown dwarf which agrees with the Bryan et al. (2018) value. Recently, however, Zhang et al. (2022) used much higher resolution VLT/CRIRES+ spectra ( R ∼ 80,000) to confirm the slow spin of 2M0355, finding v sin i = 2 . 5 -3 . 0 km s -1 . The fact that Zhang et al. (2021) and Bryan et al. (2018) used the same NIRSPEC data and found different v sin i for 2M0355 point to details in the data reduction and spectral extraction process as possible sources of discrepancy.  \nFor example, the higher spin values from older NIRSPEC data may be the result of systematic underestimation of the instrumental line spread function, which is degenerate with v sin i . For KPIC, the LSF is estimated using the spectral trace widths in the spatial direction, and conservatively allowed to vary between 1.0 to 1.2 times this width following Wang et al. (2021c), who found that the LSF is ≈ 1 . 1 times broader in the spectral direction than the spatial direction (see § 4.1). A major advantage of the single-mode fiber injection used by KPIC is the stability of the LSF, whose shape does not vary with the adaptive optics (AO) correction quality (only the intensity varies). Furthermore, we confirm that the KPIC LSF is constant over the course of an observing night by re-fitting the LSF every 20 minutes over a 3.5-hour period. The LSF weights in each refitting vary by < 5%, with most of the variation in the line wings, suggesting we are limited by the signal-to-noise of each substack and that the flux-weighted LSF is even more stable. In addition, the LSF in KPIC is measured to vary as a function of wavelength by ∼ 10 -20% within and across different spectral orders. In regular NIRSPEC, the resolution should also vary with wavelength at a similar level, but this wavelength-dependence is neglected in pre-upgrade NIRSPEC spin studies since it is hard to measure (Xuan et al. 2020).  \nIn contrast, for the pre-upgrade NIRSPEC analysis, the LSF was determined by fitting telluric lines in the host star spectra. For example, Xuan et al. (2020) found R = 24 , 800 ± 1000 from the DH Tau dataset, which was adopted as a Gaussian instrumental broadening kernel and assumed to be constant with wavelength. We note that NIRSPAO (NIRSPEC in AO mode) was used for these observations, and the size of the PSF was comparable or smaller than the slit width of 0.041 '' (Xuan et al. 2020). In this case, variations in the AO correc-  \nquality cause varying PSF sizes, and when the PSF becomes smaller than the slit this would cause a narrower LSF. However, the effect of a varying LSF should be captured by measuring the instrumental resolution from telluric lines in the data, which are varying in the same way. However, as shown in Xuan et al. (2020), the approach of determining the LSF width from telluric lines can sometimes produce poor fits and discrepant results between different spectral orders, which may bias the adopted LSF estimate.  \nWe conclude that future work is required to explain the v sin i discrepancies in more detail. Given the higher spectral resolution, higher S/N, and more stable LSF of KPIC, and the fact that our KPIC-measured spin of GQ Lup b agrees with R ∼ 100,000 CRIRES data (Schwarz et al. 2016), we conclude that our measurements are likely more reliable than pre-upgrade NIRSPEC observations of the same objects. Despite some differences, our new v sin i results are in good agreement with the findings from Bryan et al. (2020a) that low-mass substellar companions rotate much slower than their breakup velocities, pointing to mechanisms such as magnetic breaking that can efficiently reduce the companion's angular momentum on short timescales (Bryan et al. 2018; Batygin 2018). Finally, our measured abundances for κ Andb and ROXs 42B b both agree with previous work using medium- and high-resolution spectra (Hoch et al. 2020; Inglis et al. 2023), which is a validation of the power of these data to constrain chemical abundances in substellar companions and exoplanets. Except for ROXs 42B b, our KPIC measurements of the other seven companions are the first high-resolution retrievals of these objects.", '6.2. Rotation Rates and RVs': "Here, we place our measured v sin i into context with literature spin measurements for low-mass companions and exoplanets. Given the discussion above, we update the v sin i 's for DH Tau b, ROXs 42B b, ROXs 12 b, 2M0122 b, and GSC 6214-210 b. For bound companions, we compile literature measurements from Bryan et al. (2020c), and recent KPIC studies from Wang et al. (2021c); Xuan et al. (2022); Wang et al. (2022); Xuan et al. (2024); Hsu et al. (2024a). For field brown dwarfs, we compile results from Crossfield et al. (2014); Tannock et al. (2021); Hsu et al. (2021, 2024b); Vos et al. (2022). The spin measurements come from both v sin i and rotational period measurements. For the v sin i values, we assume isotropically distributed inclinations to remove the unknown i .  \nWe compute how close the objects' final rotation speeds compare to their final break-up velocities by  \nevolving their ages, radii, and spin velocities to 5 Gyr assuming constant angular momentum evolution, following the methodology detailed in Wang et al. (2022) and Hsu et al. (2024a). We adopt the evolutionary model from Baraffe et al. (2003) for this. Evolving the rotation to the same age allows us to remove age-dependent effects from literature measurements of v sin i or photometric rotational periods. Fig. 8 shows that the spin measurements in our sample fall within the overall trend compared to the literature measurements, and the bestfit rotational trend from Wang et al. (2021c), who identified a tentative anti-correlation between fractional rotation speed and object mass. Five companions in our sample have ages ≲ 10 Myr (ROXs 12 b, ROXs 42 Bb, DH Tau b, GQ Lup b, HIP 79098 b). These companions also have slower v sin i values below 10 km s -1 and are likely still undergoing gravitational contraction and gradually spinning-up. This fits the findings of Bryan et al. (2020a), who showed that younger substellar objects generally have lower rotation speeds than older objects, and that their rotation speeds increase as their radii contract with age following constant angular momentum evolution. On the other hand, the fractional rotation velocities of isolated brown dwarfs (squares) are much more scattered, and these field objects do not exhibit a clear trend for rotation. However, we note that field brown dwarfs with masses ≲ 30 M Jup fall into the Y dwarf regime, for which rotation rates are extremely challenging to measure due to their faintness. To confirm the tentative trend between rotation rate and mass, it would be useful to extend v sin i or rotational period measurements to young, directly imaged planets with masses ∼ 1 -10 M Jup . A detailed analysis of all KPIC rotation measurements will be presented in a future study.  \nIn addition to measuring v sin i , we also measure the RV of the companions at the observed epochs. We provide the companion RVs in Appendix F, which could be used to refine their orbits (Ruffio et al. 2019; Do ' O et al. 2023; Xuan et al. 2024).", '6.3. The Effect of Clouds on Narrowband High-Resolution Spectroscopy': "Arather unexpected finding of this paper that narrowband high-resolution spectra can be sensitive to clouds for T eff ∼ 1700 -2000 K objects, even after high-pass filtering the data and models. We illustrate the impact of clouds in Fig. 6 for κ Andb and GSC 6214-210 b, the two companions which show > 3 σ preferences for the EddySed cloud model (Ackerman & Marley 2001). In our retrievals with petitRADTRANS , the EddySed model accounts for both absorption and scattering from clouds,  \n0  \nFigure 8. Top: The fractional final rotational velocities versus companion masses for sources in our sample (labeled) compared to companions with spin measurements from the literature, as well as Jupiter and Saturn. The literature values from v sin i and rotational periods are shown in blue and gray circles, respectively. The best-fit rotational trends from Wang et al. (2021c) are plotted as light blue lines. The blue and gray dashed lines indicate typical measurement limits of v sin i (5.0 km s -1 ) and period (10 hr) from the instruments used (Hsu et al. 2024a). Bottom: Same as the top panel but with the addition of field brown dwarfs as squares.  \n<!-- image -->  \nwhose opacities are computed from optical constants of real condensates (Molli'ere et al. 2019). We find that a simple, gray opacity cloud model is not sufficient to explain their spectra. The sensitivity to clouds appears to arise from line depth variations with wavelength, which we can visualize in the κ Andb spectra after binning it down to average out the residual noise. While our  \n<!-- image -->  \nFigure 9. Left: C/O as a function of projected separation or semi-major axis for the eight companions in this study. The points are sized by companion mass, using the mass legend in the right panel. The blue line and region show the solar C/O value from Asplund et al. (2021). Right: In color: 1 σ contours of C/O and [C/H] for the eight companions, which are all consistent with solar composition (intersection of dashed gray lines) to the 2 σ level. In § 6.4.1, we justify that their host stars likely have solar composition as well. We overplot hot Jupiters which have abundance measurements for at least one oxygen-bearing and one carbon-bearing species from JWST or high-resolution spectroscopy as black points (Alderson et al. 2023; August et al. 2023; Bean et al. 2023; Bell et al. 2023; Brogi et al. 2023; Finnerty et al. 2024; Xue et al. 2024), which are also sized by mass. The hot Jupiters show considerable scatter and super-solar (and super-stellar) metallicities except for one planet, unlike our companions which cluster around solar (and stellar) composition.  \n<!-- image -->  \ndata does not contain absolute flux information, the line depth variations are preserved after the HPF continuum removal procedure.  \nDespite their preference for clouds, we note that the retrieved C/O and metallicity for κ Andb and GSC 6214-210 b are only weakly affected by clouds ( ≲ 1 σ level shifts). However, compared to the clear models, the EddySed cloudy models retrieve a much lower T eff . Specifically, the T eff from the EddySed models are more consistent with cloudy evolutionary models (see § 5.2). The fact that κ Andb and GSC 6214-210 b prefer clouds is consistent with their T eff of 1700 -2000 K and early L spectral types. Specifically, silicate clouds cause a broad absorption feature around ≈ 8 -11 µ m which has been directly observed in L dwarfs using mid-IR spectra from Spitzer and JWST (e.g. Cushing et al. 2006; Su'arez & Metchev 2022; Miles et al. 2023). For brown dwarfs with T eff ≳ 2000 K , the silicate absorption band starts to disappear (Su'arez & Metchev 2022), which is consistent with the fact that the hotter T eff > 2200 K companions in our sample show zero or weak preference clouds (0 -2 . 6 σ ).  \n2M0122 b, with T eff ∼ 1500 -1700 K , presents an exception to this pattern. This companion should also be cloudy, but its low S/N data or nearly pole-on viewing angle (see § 5.3) may prevent a stronger constraint on clouds. In addition, due to our limited S/N, it is difficult to directly visualize the effect of clouds in the GSC 6214-  \n210 b spectra. We emphasize that future work with flux-calibrated data and wider wavelength coverage is required to confirm our findings of clouds in κ Andb and GSC 6214-210 b. Retrieval studies of high S/N, highresolution spectra for cloudy L dwarfs would also provide a good test of our findings.  \nFinally, we note that our results represent an update to Xuan et al. (2022), who found that the K band spectrum of brown dwarf companion HD 4747 B ( T eff ≈ 1400 K and log g ≈ 5.3) is insensitive to clouds because the MgSiO 3 and Fe cloud bases are expected to lie below its K band photosphere. Our new findings indicate that when the cloud bases do intersect with the P -T profile near the photosphere, as is the case for T eff ∼ 1700-2000 K companions, K band high-resolution spectroscopy could be sensitive to cloud opacities despite continuum removal. There remain significant challenges in retrieving accurate abundances in the presence of clouds from low-resolution data (e.g. Burningham et al. 2017; Xuan et al. 2022; Lueber et al. 2022; Inglis et al. 2023), so high-resolution retrievals of L dwarfs may provide a complementary way forward.", '6.4. Towards Atmospheric Abundance Trends for Directly Imaged Companions': "6.4.1. Host Star Abundances  \nKnowledge of the host star abundances is important for formation inferences, since they represent proxies for  \nthe natal elemental abundances in the protoplanetary disk or molecular cloud. For the young stars ( ∼ 1 -100 Myr) in our sample, however, it is difficult to measure their C and O abundances due to complicating factors such as rotation, magnetic fields, veiling, and stellar activity. Even for field stars, C and O abundance calibrations are only recently being worked out for M and K spectral types (e.g. Souto et al. 2022; Hejazi et al. 2023); six of the eight stars in our sample fall into this category. The remaining two B9 stars, HIP 79098 and κ And, have extremely rapid rotation ( v sin i > 100 km s -1 ) that leads to significant spectral line broadening. This causes neighboring individual lines to blend with each other, most of which cannot then be distinguished even with high-resolution spectroscopy.  \nTo circumvent these challenges, studies have focused on more favorable targets: early K to late F young stars with lower rotation rates and no apparent veiling. The C and O abundances of these more solar-like stars can be estimated from atomic C and O lines in the visible using established abundance calibrations (Reggiani et al. 2024). Alternatively, we can utilize abundance measurements for other species, including Ca, Mg, Si, and Fe, which also have strong, isolated lines in the visible (e.g. Santos et al. 2008; Biazzo et al. 2012). Importantly, Ca, Mg, and Si are alpha elements along with C and O. Using the elemental abundances of ≈ 6000 stars from the Hypatia catalog (Hinkel et al. 2014), we confirm that the abundances of Ca, Mg, and Si scale together with those of C and O to the ≲ 0 . 2 dex level. Below, we summarize previous abundance measurements for stars in the same star-forming associations as our host stars. Since stars in open clusters are found to be chemically homogeneous at the < 0 . 03 dex level (e.g. De Silva et al. 2006; Bovy 2016; Ting et al. 2012; Poovelil et al. 2020), stars in the same star-forming associations should also be chemically homogeneous at a similar level (Reggiani et al. 2024).  \nIn terms of C and O, Reggiani et al. (2024) showed that the F7 star HD 181327 in the β Pic moving group has [C/H]= -0 . 08 ± 0 . 06, [O/H]= -0 . 10 ± 0 . 06, and C/O= 0 . 62 ± 0 . 08, which is very close to solar. Kinematic studies have shown that the β Pic moving group likely originated from the Scorpius-Centaurus (Sco-Cen) association (Mamajek & Feigelson 2001; Ortega et al. 2002). This means that to first order, stars in the ScoCen association should also have roughly solar C and O abundances. Five of our stars belong to sub-regions within Sco-Cen (HIP 79098, GSC 6214-610, GQ Lup, ROXS 12, and ROXs 42B). Wang et al. (2020) measured [C/H]= 0 . 04 ± 0 . 12, [O/H]= 0 . 08 ± 0 . 14, and C/O= 0 . 54 +0 . 12 -0 . 09 for HR 8799, which is an early F λ Boo  \nstar. The kinematic association of HR 8799 is unclear, but a recent study posits that HR 8799 may have formed either alone or in a since-dispersed small stellar group within a larger star-forming complex that gave birth to the Columba and Carina groups (Faramaz et al. 2021). If κ And is a member of Columba as suggested by the age measurement from Jones et al. (2016), then it likely has a solar composition as well.  \nIn terms of other elements, Santos et al. (2008) measured Fe and Si abundances for six star-forming regions including Lupus, ρ Oph, and Taurus and found solar [Fe/H] and [Si/H] values across their sample. Four of our stars belong to these regions (GQ Lup, ROXS 12, ROXs 42B, DH Tau). More recently, Biazzo et al. (2017) analyzed the spectra of six other stars in Lupus and also found [Fe / H] ≈ 0 . 03 on average. We note that the Lupus and ρ Oph regions are both embedded in the larger ScoCen association. In Taurus, D'Orazi et al. (2011) found an average [Fe / H] = -0 . 01 ± 0 . 05 and [Si / H] consistent with zero from seven stars. In the older AB Dor moving group, where 2M0122 is, Biazzo et al. (2012) find average values of [Fe / H] = 0 . 10 ± 0 . 03, [Mg / Fe] = -0 . 03 ± 0 . 03, and [Ca / Fe] = -0 . 01 ± 0 . 05, again consistent with solar abundances to within 0.1 dex. Together, these studies paint a broad picture of solar abundances in young starforming regions (e.g. see also Biazzo et al. 2011; Spina et al. 2014).  \nTherefore, to first order, all our host stars are expected to have nearly solar C and O abundances to the ∼ 0 . 1 -0 . 2 dex level. We therefore proceed under the assumption that our host stars have solar compositions, but emphasize that future work should attempt to perform precise C and O abundance measurements for young stars in various star-forming regions in order to provide context for abundance measurements of directly imaged planets and brown dwarfs, a majority of which are young. Specifically, high-resolution, nearinfrared spectroscopy has proven vitally important in characterizing the fundamental parameters and measuring the abundances of lower-mass host stars with spectral types later than mid-K (e.g. L'opez-Valdivia et al. 2021; Souto et al. 2022; Hejazi et al. 2023; Cristofari et al. 2023). In these latter type stars, the atomic lines of C and O in the optical are too weak for abundance measurements, but instead, C- and O-bearing molecules such as OH and CO can be utilized to infer the C and O abundances. As mentioned, performing abundance measurements for young, late-type stars will involve a detailed accounting of the effects of magnetic fields, line veiling, and rotation, which we are working on including in a follow-up study in order to directly determine the stellar abundances.", '6.4.2. Formation Pathways of Widely-separated, 10-30 M Jup Companions': "The C/O and [C/H] measurements for our eight companions are summarized in Fig. 9. The companions all have C/O consistent to solar within 1 σ level with the exception of our lowest S/N target, 2M0122 b, which has a lower C/O of 0 . 37 ± 0 . 08 but is still consistent with solar at the 2 σ level. Their [C/H] values are likewise all consistent with solar at the ≈ 1 -2 σ level. Since we established that their host stars most likely possess solar C and O abundances, our measurements suggest chemical homogeneity between these low-mass companions and their stars. The trend of solar C/O ratios for widely-separated companions with m ≈ 10 -30 M Jup was first noted by Hoch et al. (2023), who compiled literature measurements and their own Keck/OSIRIS results. Our only overlap with the Hoch et al. (2023) sample is κ Andb, as discussed in § 6.1. While Hoch et al. (2023) only studied C/O, our measurements of solar [C/H] provides a complementary piece of information. High-resolution studies of high-mass ( m ≈ 60 -70 M Jup ) brown dwarf companions indicate that their atmospheric abundances are consistent with those of their host stars at the ≲ 2 σ level (Wang et al. 2022; Xuan et al. 2022). Therefore, our results, along with Hoch et al. (2023), indicate that widely-separated ( ≳ 50 au) 10-30 M Jup companions, which have masses in between those of directly imaged planets and high-mass brown dwarfs, have abundance pattern that more closely resemble those of brown dwarfs.  \nOn a population level, we can also compare our abundance measurements with those made for close-in hot Jupiters (HJs; m ∼ 0 . 4 -10 M Jup , a ∼ 0 . 01 -0 . 1 au), which are found to have a large scatter in C/O and generally super-stellar atmospheric metallicities (e.g. Brogi et al. 2023; Alderson et al. 2023; August et al. 2023; Bell et al. 2023; Bean et al. 2023; Finnerty et al. 2024; Xue et al. 2024). These hot Jupiters almost certainly form via core accretion (within the CO snowline), though a variety of post-formation evolution processes are possible (see review by Fortney et al. 2021). The strong correlation between total heavy element mass and the total planet mass suggests these planets have solid cores consistent with core accretion (Thorngren et al. 2016), and the correlation between host star metallicity and their occurrence (Petigura et al. 2018; Osborn & Bayliss 2020) also indicates that metals are needed for their formation. As shown in Fig. 9, our widely separated, 10-30 M Jup companions occupy a region of the C/O and metallicity space that is clustered around solar (and likely stellar) composition, which is distinct than the the HJs. We carry out a weighted two-sample t-test to quantify the  \nstatistical likelihood that the compositions of HJs and those of our sample are drawn from the same population. We find p-values of 4 . 2 × 10 -4 for C/O (3.5 σ ) and 7 . 5 × 10 -6 for [C/H] (4 . 5 σ ), indicating that our measurements are inconsistent with the null hypothesis that the atmospheric compositions of HJs and our sample belong to the same population. Therefore, this serves as empirical evidence that our companions likely did not form via core accretion inside the CO snowline.  \nOur ∼ 10-30 M Jup companions orbit a diverse group of host stars, with stellar masses spanning an order of magnitude between ≈ 0 . 35 -0 . 45 M ⊙ for DH Tau and 2M0122 (Sebastian et al. 2021; Yu et al. 2023) to ≈ 3 -4 M ⊙ for kap And b and HIP 79098 (Janson et al. 2019; Gaia Collaboration 2022). This translates to a range of mass ratios ( q ) from ≈ 0 . 007 to 0 . 07 (with ∼ 30 -60% uncertainties on q mostly due to the uncertain companion masses). Previous studies find that multiple star formation operates in a largely scale-invariant manner, with mass ratio being a more informative quantity than masses of either component (e.g. Goodwin 2013; Duchˆene et al. 2023). For example, Duchˆene et al. (2023) found a dearth of lowmass stellar companions around intermediate-mass stars ( M = 1 . 75 -4 . 5 M ⊙ ) with q ≈ 0 . 05 -0 . 10, which matches the 'brown dwarf desert' around solar-mass stars (e.g. Grether & Lineweaver 2006; Metchev & Hillenbrand 2009; Sahlmann et al. 2011). They postulate that this mass ratio desert delineates two distinct formation regimes with systems with q ≲ 0 . 02 forming in disks (either disk instability or core accretion), and those with q ≳ 0 . 07 forming through cloud fragmentation. Seven of our companions have q ≲ 0 . 03, while 2M0122 b has q ≈ 0 . 07, so a mix of formation mechanisms are possible.  \nHow do our results fare in terms of the different formation mechanisms? Perhaps most directly, our findings are consistent with birth via direct gravitational collapse in a massive protostellar disk or cloud fragmentation, which should produce broadly stellar compositions as is seen for stellar binaries (e.g. Hawkins et al. 2020). While some studies predict that low-mass objects which form via disk instability can be subsequently enriched by the accretion of solids post-formation (Boley et al. 2011), it would be difficult to significantly alter the metallicity of the 10-30 M Jup objects studied here (see also Inglis et al. 2023). For example, a 10 × increase in metallicity for a 20 M Jup object would require the addition of ∼ 2 M Jup of solids. Assuming a solid-to-gas mass ratio of 0.01 for the disk, accretion of ∼ 2 M Jup of solids would require a disk-to-star ratio exceeding 0.2 for a solar-mass star and that all the solids in the disk accrete onto the companion. Therefore, significant metal enrichment for  \nour companions is challenging, especially as most of our stars are less massive the Sun. By a similar logic, significantly sub-solar metallicities are also unexpected for 10-30 M Jup companions that form via disk instability.  \nAlternatively, our composition measurements could also be compatible with core accretion outside the CO snowline (i.e. in situ). Outside the CO snowline, nearly all the metals are condensed into solids as there is very little gas. As a result, the solids inherit a stellar composition. This means that unless the planet accreted more gas than solids and became extremely metal poor, accretion of solids at these locations would typically produce stellar C/O and metallicities (e.g. Chachan et al. 2023). Since our companions have nearly solar and stellar metallicities, they are chemically consistent with forming via core accretion outside the CO snowline. We note that planetesimal accretion is less efficient at wider orbital distances, but pebble accretion may be able to produce a sufficiently massive core on a reasonable timescale (Lambrechts & Johansen 2012). For example, Gurrutxaga et al. (2023) demonstrated pebble accretion pathways to form the giant planets PDS 70 b ( m ≈ 2 -6 M Jup , a ≈ 21 au) and c ( m ≈ 3 -10 M Jup , a ≈ 34 au; planet values from Wang et al. 2021b). However, it remains to be seen whether realistic pebble accretion scenarios can produce our higher mass and wider separation (up to 360 au) companions. For example, half of our companions have q ≳ 0 . 03, which means that companion formation must begin early in a disk with disk-to-star mass ratio > 0 . 03 and finish before the bulk of the gas reservoir accretes onto the star. By the protoplanetary disk phase, typical disk-to-star mass ratios of ∼ 0 . 01 have insufficient mass to produce > 10 M Jup companions around solar-mass stars even if all the disk mass goes into the companion (let alone around the M dwarfs in our sample). Finally, early, massive circumstellar disks can be prone to gravitational instability (Kratter & Lodato 2016), so an additional challenge of forming companions via core accretion at large distances is the need to prevent disk fragmentation from occurring.", '6.5. The 12 C / 13 C Isotopic Ratio': "Minor isotopologues of carbon, oxygen, and nitrogen have recently been detected in exoplanets and brown dwarfs (Zhang et al. 2021; Zhang et al. 2021; Line et al. 2021; Barrado et al. 2023; Gandhi et al. 2023), allowing constraints on isotopic ratios such as 12 C / 13 C, 12 O / 18 O and 14 N / 15 N. Xuan et al. (2024) also demonstrated carbon and oxygen isotopic homogeneity between a late-M dwarf companion and its late-K primary star with Keck/KPIC spectroscopy. In this paper, we  \nFigure 10. 12 C / 13 C measured from CO as a function of projected separation or semi-major axis. The points are color-coded by companion-star mass ratio, and the three open squares are from this paper. The dashed and dashdot lines denote the solar value (Lyons et al. 2018) and ISM value (Wilson 1999), respectively. The other points are from Line et al. (2021) for WASP-77 Ab, Zhang et al. (2021) for YSES 1 b, Xuan et al. (2024) for HIP 55507 B, Costes et al. (2024) for HD 984 B, and Gandhi et al. (2023) for VHS 1256 b.  \n<!-- image -->  \nfind bounded constraints on 12 CO / 13 CO for three companions, GQ Lup b, HIP 79098 b, and DH Tau b, though the 13 CO detection for DH Tau b is tentative.  \n12 C / 13 C is proposed to be a formation diagnostic as CO ice outside the CO snowline is expected to be enriched in 13 C due to isotopic selective fractionation processes such as CO self-shielding (Zhang et al. 2021). Planets forming via core/pebble accretion outside the CO snowline may therefore show a lower 12 C / 13 C due to this effect. Our three companions show a range of 12 C / 13 C values between ≈ 50-150. The values for DH Tau b (53 +50 -24 ) and HIP 79098 b (52 +22 -17 ) are broadly consistent with the local ISM value of 12 C / 13 C = 69 ± 6 (Wilson 1999), while GQ Lup b shows a higher 12 CO / 13 CO of 153 +43 -31 . In Appendix B, we compare retrieval results from several different nights of KPIC data for GQ Lup b, and find ≈ 0.1 dex differences in the retrieved log( 12 CO / 13 CO) between different nights, which is consistent with the findings of Xuan et al. (2024). When including the night-to-night scatter of 0 . 10 dex as a systematic error, GQ Lup b's measurement is consistent with the solar value of 12 C / 13 C = 93 . 5 ± 3 . 1 (Lyons et al. 2018) at the 1 σ level. Overall, these findings do not provide evidence that our companions accreted a significant amount of ice with low 12 C / 13 C, consistent with their likely star-like formation pathways.  \nIn Fig. 10, we plot our new measurements along with previous constraints on the 12 C / 13 C for exoplanets and low-mass companions. There is considerable scatter in the measurements, and no clear trend with orbital dis-  \ntance or mass ratio. Some of these results should be revisited with better quality or higher-resolution data to reduce the error bars. For directly imaged companions, the exquisite S/N from JWST is poised to provide more precise measurements of isotopologue ratios, as demonstrated by Gandhi et al. (2023). Another direction is to obtain higher S/N data to detect additional isotopologues such as C 18 O and H 18 2 O (Zhang et al. 2022; Xuan et al. 2024) Finally, to better interpret isotopic ratio data, more work is required on the modeling side to understand the details of isotopic variability and fractionation chemistry in disks ( Oberg et al. 2023).", '6.6. System and Orbital Architectures': "The formation of substellar companions should also be considered in the context of the host stellar system. As noted in § 2, five out of eight of our systems are either confirmed or likely multiple-star systems. ROXs 42B is a resolved binary with projected separation ∼ 12 au (Kraus et al. 2013) while HIP 79098 is a suspected binary (Janson et al. 2019), making the b components around them 'circumbinary.' One way to explain these configurations is via disk fragmentation, which generally produces multiple fragments that interact gravitationally (e.g. Kratter & Lodato 2016). For example, the protostellar disk around ROXs 42B could first fragment and form the secondary star, which continues to grow given high infall rates from the surrounding nebula. At a later stage when the disk mass is lower, but still gravitationally unstable, a second fragmentation could produce the ∼ 13 M Jup companion. Dynamical interactions may cause the substellar companion to move away from the star while the secondary star moves inward.  \nOn the other hand, DH Tau, ROXs 12, and GQ Lup have wide stellar companions at thousands of au (Kraus & Hillenbrand 2009; Bowler et al. 2017; Alcal'a et al. 2020). Alcal'a et al. (2020) hypothesized that GQ Lup A and its wide ( ≈ 2400 au) stellar companion may have formed via turbulent fragmentation of a molecular cloud core, while GQ Lup b fragmented out of the circumprimary disk around GQ Lup A. Similar formation pathways may be responsible for the DH Tau and ROXs 12 systems as well. While disk fragmentation simulations do produce low-mass substellar companions (e.g. Stamatellos & Whitworth 2009), it is challenging to prevent these companions from further accretion to become low-mass stars or be tidally destroyed during a rapid inward migration phase powered by infall (Zhu et al. 2012; Forgan & Rice 2013). Indeed, to produce the 1030 M Jup companions we observe today, fragmentation must occur in a narrow window after the period of high infall rates, but before the disk mass becomes too small  \nfor Toomre instability (Kratter et al. 2010). Therefore, if disk fragmentation is a common formation mechanism, we should observe more low-mass stars around the same stellar types compared to substellar or planetary companions. Published observations are in qualitative agreement with this prediction. For example, Duchˆene et al. (2023) find a significantly higher number of lowmass stellar companions around intermediate-mass stars compared to substellar companions.  \nOrbital architectures provide another piece of the puzzle. One informative probe is the relative alignment between different angular momentum vectors in the system. Recently, Bowler et al. (2023) presented a summary of stellar obliquity constraints for directly imaged substellar companions, finding that misalignment between the companion orbits and stellar spin axes are common. Specifically, GQ Lup b may be on a nearly perpendicular orbit with respect to the circumstellar disk around GQ Lup A, and is also misaligned with the star's spin axis (Stolker et al. 2021). Similarly, the orbit of ROXs 12 b is also likely misaligned compared to the stellar spin axis (Bowler et al. 2017; Bowler et al. 2023). Bryan et al. (2020b) provided the first measurement of companion obliquity for 2M0122 b, finding that the companion's spin axis is tentatively misaligned with respect to the stellar spin. 7 Since stellar binaries that form via turbulent fragmentation tend to have misaligned spin and orbital orientations (e.g. Jensen & Akeson 2014; Lee et al. 2016; Offner et al. 2016), random spin-orbit distributions should also be prevalent for systems with substellar companions if they form in this manner, which is consistent with current findings. Alternatively, wide, misaligned stellar companions in these systems could torque the disk around the primary, causing the companions that fragment in the disk to inherit such misalignments (e.g. Bowler et al. 2017). Finally, disk turbulence itself could result in random spin-orbit misalignments (Jennings & Chiang 2021). Regardless of the detailed mechanism, the system and orbital architectures of our systems are broadly consistent with star-like formation. This is consistent with the trend of chemical homogeneity we infer between these companions and their stars.", '7. CONCLUSION': "7 We note that our higher KPIC v sin i measurement of 19 . 6 +3 . 0 -2 . 5 km s -1 for 2M0122 b (see Table 5) would increase the line-of-sight inclination of the companion's spin axis, making it more aligned with the stellar spin axis. Due to the low S/N of our 2M0122 b spectra however, we emphasize that future data will be required to confirm the exact v sin i of this companion.  \nWe have carried out a uniform atmospheric retrieval study of eight widely-separated substellar companions ( ∼ 50-360 au, ∼ 10-30 M Jup ) with Keck/KPIC highresolution K band spectroscopy. From these retrievals, we measure the companion's C/O, metallicity (denoted [C/H]), isotopic abundances, in addition to their cloud properties, spins, radial velocities, and temperature profiles. To complement the continuum-removed highresolution data, we adopt mass and radius priors from evolutionary models in the retrievals.  \nFirst, we find that these companions have broadly solar composition (to within 2 σ level), and likely stellar composition given the trend of solar abundances seen for stars in the same star-forming associations. Their abundance pattern is similar to systems with high-mass brown dwarf companions and stellar binaries, which show chemically homogeneity (e.g. Hawkins et al. 2020; Wang et al. 2022; Xuan et al. 2022). On the other hand, their abundances are distinct from hot Jupiters, which show a range in C/O and generally super-stellar metallicities (see Fig. 9), and directly imaged giant planets with m ∼ 3 -10 M Jup , which show tentative metal enrichment (Wang 2023a). Thus, the population of lowmass substellar companions from direct imaging likely traces the tail-end of star formation processes such as gravitational disk instability and cloud fragmentation, making them low-mass brown dwarfs instead of 'superJupiter' planets or 'planetary-mass companions.' Alternatively, we note that our composition measurements are also consistent with core accretion outside the CO snowline where these companions are observed today, since accretion at these locations would also yield stellar C/O and metallicities. However, such a scenario requires core accretion to proceed early and rapidly in a massive, protostellar disk in order to explain the accretion of 10-30 M Jup of material, especially for the systems with lower mass host stars.  \nSecond, we find evidence of clouds in two of the colder companions (Fig. 6), kap And b ( T eff = 1640 +220 -170 K ) and GSC 6214-210 b ( T eff = 1860 +170 -110 K ), with the EddySed cloud model being preferred by ≳ 3 σ compared to both the clear and gray cloud models. This indicates that narrow-band, high-resolution spectra can be sensitive to non-gray, scattering clouds when the cloud opacity is high near the photosphere, as is the case for these companions. This result highlights the potential of high-resolution spectroscopy in constraining both abundances and clouds for brown dwarfs and exoplanets.  \nThird, we present three new measurements of 12 CO / 13 CO for GQ Lup b (153 +43 -31 ), HIP 79098 b (52 +22 -17 ), and DH Tau b (53 +50 -24 ). A cross-correlation analysis shows solid detections for the first two compan-  \nions, and a tentative detection of 13 CO for DH Tau b. From retrievals of six independent KPIC datasets of GQLup b data (Appendix B), we find ≈ 0 . 1 dex systematic errors in 12 CO / 13 CO, which is consistent with the findings of Xuan et al. (2024). After accounting for systematics, our measurements agree with either the ISM or solar 12 C / 13 C to within 1 σ . We place these measurements in context of previous work and do not identify any clear trends between 12 C / 13 C and mass, mass ratio, or orbital distance at this stage (Fig. 10). More precise and accurate measurements of 12 C / 13 C and complementary measurements of the isotopic ratio in the host stars are necessary to further interpret these results.  \nFinally, we present radial velocity and spin measurements for the companions. We find some discrepancies between our v sin i values and previous studies using pre-upgrade Keck/NIRSPEC data (Bryan et al. 2020a), while our v sin i for GQ Lup b agrees with previous VLT/CRIRES measurements at a higher resolution of R ∼ 100 , 000 (Schwarz et al. 2016). Due to the higher spectral resolution ( R ∼ 35 , 000), higher S/N, and more stable line spread function of KPIC (which uses postupgrade NIRSPEC; Martin et al. 2018) compared to pre-upgrade NIRSPEC ( R ∼ 25 , 000), we adopt our new spin measurements in this paper. Seven out of eight of our companions have relatively slow v sin i ≈ 4 -20 km s -1 , with kap And b being an outlier with v sin i = 39 . 4 ± 0 . 9 km s -1 (see also Morris et al. 2024, accepted). Our objects follow the overall trend of literature spin measurements, and display a large scatter in their fractional rotational velocities, v/v breakup , as shown in Fig. 8.  \nLooking forward, it would be useful to obtain higher S/N spectra for some of the companions in order to obtain higher precision measurements. The results for 2M0122 b in particular should be re-visited, as this is our lowest S/N dataset and we also retrieve a much lower C/O for this companion than the rest of the sample. Higher S/N data would also help increase the quantity and quality of 12 CO / 13 CO measurements. For example, Gandhi et al. (2023) showed that JWST/NIRSpec enables isotopologue measurements with unprecedented precision for widely-separated substellar companions. In addition, detailed forward modeling of the residual starlight could also open the door to studying closein, high-contrast companions with the JWST/NIRSpec integral field unit (Ruffio et al. 2023a). The wide wavelength coverage of JWST would be useful for constraining cloud properties and measuring more elemental abundances besides C and O in order to break degeneracies in formation inferences as well (e.g. Cridland  \net al. 2020; Turrini et al. 2021; Ohno & Fortney 2023; Chachan et al. 2023).  \nMany of the companions presented in this work have relatively large separations and low star-to-companion contrasts, making them accessible with traditional NIRSPEC in AO mode (Bryan et al. 2018). A natural next step would be to extend our measurements to companions with lower masses and smaller orbital distances. Keck/KPIC and similar instruments like VLT/HiRISE (Vigan et al. 2024) are ideal for these targets, as shown by Wang et al. (2021c). By targeting more directly imaged planets with high-resolution spectroscopy, we could further test the tentative trend of metal enrichment suggested by Wang (2023a) for these planets. Such measurements are especially important given literature discrepancies in the retrieved abundances for the same planet, either due to using data with different spectral resolution (e.g. β Pic b; Gravity Collaboration et al. 2020; Landman et al. 2023), or using different retrieval methods on similar low-resolution data (e.g. 51 Eri b; Brown-Sevilla et al. 2023; Whiteford et al. 2023). As demonstrated by this paper and previous work, atmospheric abundance measurements from medium-to-high resolution spectroscopy are providing results that are reliable to uncertain assumptions about clouds (Zhang et al. 2021; Xuan et al. 2022; Inglis et al. 2023), and consistent across different observing nights (Ruffio et al. 2021; Xuan et al. 2024; Landman et al. 2023). Therefore, medium-to-high resolution spectroscopy is poised to improve our understanding of the nature and formation of directly imaged companions in the near future, bringing more clarity in the delineation between giant planets and their brown dwarf counterparts.", 'ACKNOWLEDGMENTS': "J.X. thanks Melanie Rowland, Yayaati Chachan, Michael Liu, Kazumasa Ohno, Douglas Lin, Ricardo L'opez-Valdivia, Henrique Reggiani, and Jingwen Zhang for helpful discussions. J.X. is supported by the NASA Future Investigators in NASA Earth and Space Science and Technology (FINESST) award #80NSSC23K1434. J.X. also acknowledges support from the Keck Visiting Scholars Program (KVSP) to commission KPIC Phase II capabilities. D.E. is supported by the NASA FINESST award #80NSSC19K1423. D.E. also acknowledges support from the Keck Visiting Scholars Program (KVSP) to install the Phase II upgrades. Funding for KPIC has been provided by the California Institute of Technology, the Jet Propulsion Laboratory, the Heising-Simons Foundation (grants #2015-129, #2017318, #2019-1312, #2023-4598), the Simons Foundation, and the NSF under grant AST-1611623. The computations presented here were conducted in the Resnick High Performance Center, a facility supported by Resnick Sustainability Institute at the California Institute of Technology. W. M. Keck Observatory access was supported by Northwestern University and the Center for Interdisciplinary Exploration and Research in Astrophysics (CIERA). The data presented herein were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California and the National Aeronautics and Space Administration. The Observatory was made possible by the generous financial support of the W. M. Keck Foundation. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to conduct observations from this mountain. This work benefited from the 2023 Exoplanet Summer Program in the Other Worlds Laboratory (OWL) at the University of California, Santa Cruz, a program funded by the Heising-Simons Foundation and NASA. The research here acknowledges use of the Hypatia Catalog Database, an online compilation of stellar abundance data as described in Hinkel et al. (2014), which was supported by NASA's Nexus for Exoplanet System Science (NExSS) research coordination network and the Vanderbilt Initiative in Data-Intensive Astrophysics (VIDA).  \nFacilities: Keck (KPIC)  \nSoftware: petitRADTRANS (Molli`ere et al. 2019), (Speagle 2020)  \ndynesty", 'A. INTERPOLATED PROPERTIES FROM EVOLUTIONARY MODELS': 'Here, we show the evolutionary-model predicted mass, radius, log g , and T eff distributions for the seven other companions. Different colors indicate different evolutionary models, with ATMO 2020 in blue (Phillips et al. 2020; Chabrier et al. 2023), AMES-Dusty in red (Allard et al. 2001), BHAC15 in purple (Baraffe et al. 2015), and SM08 in gray (Saumon & Marley 2008). The distributions for GQ Lup b are shown in Fig. 2.  \n<!-- image -->  \nFigure A1. Same as Fig. 2, but for the other seven companions.  \n<!-- image -->', 'B. COMPARING CONTINUUM REMOVAL METHODS WITH MULTIPLE DATASETS OF GQ LUP B': 'We compare the retrieved abundances and v sin i between different observing epochs for GQ Lup b, our highest S/N target. In addition, we compare two methods of removing the continuum: the high-pass filter (HPF) method used in the paper and a spline model. As described in § 4.1, a median-filter with a size of 100 pixels ( ≈ 0 . 002 µ m) is adopted as the HPF. The choice of the filter size is motivated by injection-recovery tests in Xuan et al. (2022), while the median-filter is found to perform better compared to Fourier-based filters or Gaussian filters (Wang et al. 2021c).  \nFor the alternative, spline model, we follow the method in Ruffio et al. (2023a) using the open-source package breads (Agrawal et al. 2023). 8 Since the spline model allows us to modulate the companion continuum and the speckle continuum separately, we experimented with a few different choices. Our default spline model modulates the continuum of both speckle and companion with a third-order spline model. Specifically, we use three spline nodes per spectral order. We note that such a spline operates over a much larger scale ( ≈ 1000 pixels) compared to our 100-pixel median filter, so is not an exact comparison to the HPF method. When we instead used 10 nodes per spectral order for the speckle continuum, the results were consistent with the 3-node model but the continuum oscillations appear more stochastic and less smooth, leading us to disfavor these higher number of nodes. Finally, we also experimented with a simpler model of not modulating the companion continuum and only modulating the speckle continuum. We find that this did not change the results significantly either.  \nTo compare the two continuum removal methods, we perform retrievals on six independent datasets for GQ Lup b across three different nights (two fibers were used each night), UT 2023 June 23, June 24, and June 29. The same retrieval setup and priors are used. For the HPF, we jointly fit for flux scales for the companion and the speckle contribution, α b and α s (see Eq. 1). Similarly, the spline model performs a linear optimization to determine the spline parameters for each proposed model.  \nIn Fig. B1, we compare the retrieved C/O, [C/H], and log( 12 CO / 13 CO) between the different datasets for each continuum removal method. First, we note that the retrieved values between the two methods agree well; taking the median across the results from six datasets, we obtain C / O = 0 . 65, [C / H] = 0 . 17, log( 12 CO / 13 CO) = 2 . 07 from the HPF method, and C / O = 0 . 67, [C / H] = 0 . 21, log( 12 CO / 13 CO) = 2 . 15 from the spline method. Next, we quantify the relative agreement between different datasets with the standard deviation ( s ) of the median retrieved values. We find that the HPF yields a slightly lower spread in C/O by ≈ 35% compared to the spline method. The spread in [C/H] and log( 12 CO / 13 CO) is comparable between the two methods. Therefore, we conclude that the performances of the two continuum removal methods are comparable for our data. We adopt the HPF as the continuum removal method in this paper to be consistent with previous KPIC retrieval papers (e.g. Wang et al. 2021c; Xuan et al. 2022, 2024), but note that our results will not be affected if we chose the spline model instead.  \nBy comparing results from different datasets, we can also assess the level of systematic error in the retrieved parameters. Here, we focus on the HPF results. Approximating the standard deviation as the systematic error, we find night-to-night systematics of 0 . 022, 0 . 14 dex, and 0 . 10 dex for C/O, [C/H], and log( 12 CO / 13 CO). These systematic errors are comparable in size to those found in Xuan et al. (2024) by comparing two independent KPIC datasets for a M7.5 companion.  \n5  \n<!-- image -->  \nFigure B1. Left: posterior distributions of three key parameters from using the spline model. Six different datasets are compared (3 nights × 2 fibers per night). Right: the same but using the HPF (median-filter of 100 pixels). The titles on each histogram show the standard deviation ( s ) of the median between the six different posteriors. The legend is shared between both plots.  \n<!-- image -->', 'C. PRIORS FOR RETRIEVALS': 'Here, we list the fitted parameters and adopted priors for retrievals of the seven other companions. The corresponding information for κ Andb is shown in Table 3.  \nTable C1. Fitted Parameters and Priors for Retrievals  \nNote -For definitions and table notes see Table 3. In the table entries for companions after GQ Lup b, we omit a few common parameters that use the same priors. The different cloud species considered for each companion are indicated. Note that we only fit for speckle flux in the GQ Lup b and kap And b data, as explained in § 4.1.  \nFigure D1. Same as Fig. 5, but for HIP 79098 b.  \n<!-- image -->', 'D. BEST-FIT MODELS FROM RETRIEVALS': 'Here, we plot the KPIC spectra and best-fit models for the seven other companions. The corresponding plot for GQ Lup b is shown in Fig. 5.  \nFigure D2. Same as Fig. D1, but for DH Tau b.  \n<!-- image -->  \nFigure D3. Same as Fig. D1, but for kap And b.  \n<!-- image -->  \nFigure D4. Same as Fig. D1, but for GSC 6214-210 b.  \n<!-- image -->  \nFigure D5. Same as Fig. D1, but for ROXs 12 b.  \n<!-- image -->  \n(counts)  \nFlux  \n(counts)  \nFlux  \nFigure D6. Same as Fig. D1, but for ROXs 42B b.  \n<!-- image -->  \nFigure D7. Same as Fig. D1, but for 2M0122 b. Visually, and from the CCF analysis (Fig. 3), this is our lowest S/N dataset. Retrievals for 2M0122 b should be re-visited with improved S/N data.  \n<!-- image -->', 'E. PRESSURE-TEMPERATURE PROFILES AND EMISSION CONTRIBUTION FUNCTIONS': 'Here, we plot the retrieved temperature profiles and emission contribution functions for the five other companions. The corresponding plots for κ Andb, GSC 6214-210 b, and HIP 79098 b are shown in Fig. 4. We note that a superadiabatic region is visible in the P-T profile of GQ Lup b around 5 × 10 -1 bars. Deviations from self-consistent 1D models such as SPHINX and Sonora are not unexpected due to 3D effects such as rotation-induced horizontal transport (Tan & Showman 2021). However, our data alone is insufficient to show that this super-adiabatic region is real. For purposes of this paper, we check that the retrieved P-T profile of GQ Lup b does not bias its retrieved parameters by running a separate retrieval with a fixed P-T profile matching the T eff = 2600 K, log g = 4 . 0 SPHINX profile. The resulting posteriors are consistent to within 1 σ for all parameters, demonstrating that the results are not sensitive to the exact shape of the P-T profile.  \n<!-- image -->  \n(bar)  \nPressure  \n10  \n-  \n10  \n4  \n3  \n-  \n10  \n-  \n10  \n-  \n2  \n1  \n0  \n10  \n1  \n0  \n10  \n2M 0122 b  \n1000  \n2000  \nTemperature (  \nMgSiO3  \nFe  \nT  \neff  \nT  \neff  \nK  \n)  \nFigure E1. Same as Fig. 4, but for the other five companions.  \n<!-- image -->  \n2.30  \n2.30  \n2.35  \nWavelength (um)  \nWavelength (um)  \n=  \n1600  \nK  \n, log(  \ng  \n)=4.5  \n=  \n1500  \nK  \n, log(  \ng  \n)=4.5  \n3000  \n4000  \n10-3  \nL  \n10-1  \n100  \n10-5  \n10-4  \nL  \n10-2  \n10-1  \n100', 'F. RADIAL VELOCITIES FOR THE SUBSTELLAR COMPANIONS': 'Here, we list the measured RVs of the eight companions, which have been corrected for barycentric motion.  \nTable F1. KPIC radial velocity measurements for the eight substellar companions studied in this work. We have applied the barycentric correction to the RVs, so their reference is the solar system barycenter.', 'REFERENCES': "Ackerman, A. S., & Marley, M. S. 2001, The Astrophysical Journal, 556, 872, doi: 10.1086/321540  \nAgrawal, S., Ruffio, J.-B., Konopacky, Q. M., et al. 2023,  \nAJ, 166, 15, doi: 10.3847/1538-3881/acd6a3  \nAlcal'a, J. M., Majidi, F. Z., Desidera, S., et al. 2020,  \nAstronomy and Astrophysics, 635,  \ndoi: 10.1051/0004-6361/201937309  \nAlderson, L., Wakeford, H. R., Alam, M. K., et al. 2023, Nature, 614, 664, doi: 10.1038/s41586-022-05591-3  \nAllard, F., Hauschildt, P. H., Alexander, D. R., Tamanai, A., & Schweitzer, A. 2001, The Astrophysical Journal, 556, 357, doi: 10.1086/321547  \nAsplund, M., Amarsi, A. M., & Grevesse, N. 2021,  \nAstronomy and Astrophysics, 653,  \ndoi: 10.1051/0004-6361/202140445  \nAugust, P. C., Bean, J. L., Zhang, M., et al. 2023, The Astrophysical Journal Letters, 953, L24,  \ndoi: 10.3847/2041-8213/ace828  \nAzzam, A. A., Tennyson, J., Yurchenko, S. N., & Naumenko, O. V. 2016, Monthly Notices of the Royal Astronomical Society, 460, 4063,  \ndoi: 10.1093/mnras/stw1133  \nBalmer, W. O., Pueyo, L., Stolker, T., et al. 2023, ApJ, 956, 99, doi: 10.3847/1538-4357/acf761  \nBaraffe, I., Chabrier, G., Barman, T. S., Allard, F., & Hauschildt, P. H. 2003, A&A, 402, 701,  \ndoi: 10.1051/0004-6361:20030252  \nBaraffe, I., Homeier, D., Allard, F., & Chabrier, G. 2015, Astronomy and Astrophysics, 577,  \ndoi: 10.1051/0004-6361/201425481 Barman, T. S., Konopacky, Q. M., Macintosh, B., & Marois, C. 2015, ApJ, 804, 61, doi: 10.1088/0004-637X/804/1/61 Barman, T. S., Macintosh, B., Konopacky, Q. M., & Marois, C. 2011, ApJ, 733, 65, doi: 10.1088/0004-637X/733/1/65 Barrado, D., Molli`ere, P., Patapis, P., et al. 2023, arXiv e-prints, arXiv:2311.08054, doi: 10.48550/arXiv.2311.08054 Bate, M. R. 2012, MNRAS, 419, 3115, doi: 10.1111/j.1365-2966.2011.19955.x Batygin, K. 2018, AJ, 155, 178, doi: 10.3847/1538-3881/aab54e Bean, J. L., Xue, Q., August, P. C., et al. 2023, Nature, 618, 43, doi: 10.1038/s41586-023-05984-y Bell, C. P. M., Mamajek, E. E., & Naylor, T. 2015, MNRAS, 454, 593, doi: 10.1093/mnras/stv1981 Bell, T. J., Welbanks, L., Schlawin, E., et al. 2023, Nature, 623, 709, doi: 10.1038/s41586-023-06687-0 Bernath, P. F. 2020, Journal of Quantitative Spectroscopy and Radiative Transfer, 240, doi: 10.1016/j.jqsrt.2019.106687  \nBiazzo, K., D'Orazi, V., Desidera, S., et al. 2012, Monthly Notices of the Royal Astronomical Society, 427, 2905, doi: 10.1111/j.1365-2966.2012.22132.x  \nWhiteford, N., Glasse, A., Chubb, K. L., et al. 2023, Monthly Notices of the Royal Astronomical Society, 525, 1375, doi: 10.1093/mnras/stad670  \nWilson, T. L. 1999, Rep. Prog. Phys, 62, 143  \nXuan, J. W., Bryan, M. L., Knutson, H. A., et al. 2020,  \nThe Astronomical Journal, 159, 97,  \ndoi: 10.3847/1538-3881/ab67c4  \nXuan, J. W., Wang, J., Ruffio, J.-B., et al. 2022, The Astrophysical Journal, 937, 54, doi: 10.3847/1538-4357/ac8673  \nXuan, J. W., Wang, J., Finnerty, L., et al. 2024, The  \nAstrophysical Journal, 962, 10,  \ndoi: 10.3847/1538-4357/ad1243  \nXue, Q., Bean, J. L., Zhang, M., et al. 2024, ApJL, 963, L5, doi: 10.3847/2041-8213/ad2682  \nYelverton, B., Kennedy, G. M., Su, K. Y. L., & Wyatt,  \nM. C. 2019, MNRAS, 488, 3588,  \ndoi: 10.1093/mnras/stz1927  \nYu, J., Khanna, S., Themessl, N., et al. 2023, The  \nAstrophysical Journal Supplement Series, 264, 41,  \ndoi: 10.3847/1538-4365/acabc8"}