Stellar Properties of Observed Stars Stripped in Binaries in the Magellanic Clouds

Massive stars ( ∼ 8 – 25 M e ) stripped of their hydrogen-rich envelopes via binary interaction are thought to be the main progenitors for merging neutron stars and stripped-envelope supernovae. We recently presented the discovery of the ﬁ rst set of such stripped stars in a companion paper. Here, we ﬁ t the spectra of 10 stars with new atmosphere models in order to constrain their stellar properties precisely. We ﬁ nd that the stellar properties align well with the theoretical expectations from binary evolution models for helium-core burning envelope-stripped stars. The ﬁ ts con ﬁ rm that the stars have high effective temperatures ( T eff ∼ 50 – 100 kK ) , high surface gravities ( g log ~ 5 ) , and hydrogen-poor / helium-rich surfaces ( X H,surf ∼ 0 – 0.4 ) while showing for the ﬁ rst time a range of bolometric luminosities ( 10 3 – 10 5 L e ) , small radii ( ∼ 0.5 – 1 R e ) , and low Eddington factors ( Γ e ∼ 0.006 – 0.4 ) . Using these properties, we derive intermediate current masses ( ∼ 1 – 8 M e ) , which suggest that their progenitors were massive stars ( ∼ 5 – 25 M e ) and that a subset will reach core-collapse, leaving behind neutron stars or black holes. Using the model ﬁ ts, we also estimate the emission rates of ionizing photons for these stars, which agree well with previous model expectations. Further, by computing models for a range of mass-loss rates, we ﬁ nd that the stellar winds are weaker than predicted by any existing scheme ( M 10 wind 9  -  M e yr − 1 ) . The properties of this ﬁ rst sample of intermediate-mass helium stars suggest they both contain progenitors of type Ib and IIb supernovae, and provide important benchmarks for binary evolution and population synthesis models.


INTRODUCTION
Helium stars with masses intermediate between subdwarfs and Wolf-Rayet (WR) stars (∼ 2-8 M ) have been predicted to be created through mass transfer or common envelope Corresponding author: Y. Götberg, M. R. Drout ygoetberg@carnegiescience.edu maria.drout@utoronto.ca* Hubble Fellow ejection in binary stars with initial primary star masses of ∼ 8-25 M (e.g., Kippenhahn & Weigert 1967;Paczyński 1967;Ivanova 2011).These envelope-stripped stars should be common (Götberg et al. 2019;Shao & Li 2021), because a large fraction of massive binaries go through envelopestripping (∼30%, Sana et al. 2012), and the long-lasting helium-core burning phase usually remains after envelopestripping (e.g., Pfahl et al. 2002;de Mink et al. 2008, see however also Klencki et al. 2022).Because of their ubiquity, stripped stars have been proposed as the main progenitors of stripped-envelope supernovae (Smith et al. 2011b;Yoon et al. 2017;Sravan et al. 2019), which also matches with their low ejecta masses (Drout et al. 2011;Lyman et al. 2016).Envelope-stripping is also considered necessary for the creation of merging compact objects (Kalogera et al. 2007).For example, the evolutionary channel to merging binary neutron stars includes two stripped stars (Tauris et al. 2017;Vigna-Gómez et al. 2020;Ye et al. 2020).In addition, stripped stars are also so small that they can emit low-frequency gravitational waves detectable with the Laser Interferometer Space Antenna (LISA), when stripped by a compact object (Nelemans et al. 2004;Wu et al. 2018Wu et al. , 2020;;Götberg et al. 2020b;Kupfer et al. 2020;Liu et al. 2022) Furthermore, with their high effective temperatures (T eff ∼50-100kK), stripped stars should emit most of their radiation in the ionizing regime, thus providing a boost of ionizing emission several tens of millions of years after a starburst (Stanway et al. 2016;Götberg et al. 2019Götberg et al. , 2020a)).However, although "intermediate mass" stripped stars have many interesting implications, an observed sample of them was missing until recently.
Previous efforts have been made in the search for stripped helium stars, resulting in discoveries on the low-and highmass ends.In an impressive search for hot companions orbiting Galactic Be stars using ultraviolet (UV) spectroscopy, a set of hot subdwarf companions have been revealed (Wang et al. 2017(Wang et al. , 2018(Wang et al. , 2021)).With flux contributions of only up to ∼10% in the UV, the subdwarfs likely have low masses of ∼ 0.5-1.5 M (Klement et al. 2022a,b), which suggests that the bright and early type Be-star companions became more massive and more luminous after they gained significant mass from the donor star during conservative mass transfer.Subdwarfs that instead orbit faint companions have been studied for example by Schaffenroth et al. (2022).Also, during the recent searches for black holes, a number of inflated, low-mass (∼ 0.5 M ) stripped stars were unveiled instead (e.g., Irrgang et al. 2020;Bodensteiner et al. 2020;El-Badry et al. 2022).In addition, the star υ Sag, which was thought to be a ∼3 M intermediate mass helium giant (Dudley & Jeffery 1990), has recently been determined to have <1 M (Gilkis & Shenar 2022).In the higher mass range, searches for companions to WR stars that may have been responsible for the envelope-stripping (Vanbeveren et al. 1998) has been done (Shara et al. 2017(Shara et al. , 2020;;Shenar et al. 2019).In particular, the WR X-ray binary Cyg X-3 likely evolved via binary interaction, indicated from its short orbital period (van den Heuvel & De Loore 1973;van Kerkwijk et al. 1992).
While the above described studies are important for our understanding of interacting binaries, none of them included helium stars of intermediate mass.In fact, the only previously known intermediate mass stripped star is the ∼4 M quasi Wolf-Rayet (WR) star in the binary system HD 45166, however, even this star has recently been observed to have lower mass than previously thought (∼2 M , T. Shenar, pri-vate communication).However, in Drout & Götberg et al., under review, we presented a new sample of 25 stars in the Magellanic Clouds.Originally identified has having excess UV radiation in comparison to the main-sequence (Götberg et al. 2018), we demonstrate that they have colors, brightnesses, and optical spectra consistent with expectations for binary systems containing intermediate mass helium stars.In particular, their spectral morphologies fall into three broad categories, as expected for systems with a range of mass ratios: (1) those consistent with a stripped helium star dominating the optical flux of the system, (ii) those consistent with both a stripped star and a main-sequence companion contributing to the optical flux, and (iii) those consistent with a main sequence companion dominating the optical flux of the system.By comparing the measured equivalent widths of several diagnostic lines for the stars in Class 1, we were able to obtain rough estimates for their physical properties, demonstrating that they have hot temperatures (T eff 70kK), high surface gravities (log(g) ∼ 5), and depleted surface compositions (X H,surf 0.3), further solidifying their nature as intermediate mass helium stars.
Full characterization of the stripped star binary sample of Drout & Götberg et al., under review will deepen our understanding of binary interaction significantly, as it would produce direct constraints for binary evolution and population models.While the approximate effective temperatures, surface gravities and surface compositions presented in Drout & Götberg et al., under review were sufficient to establish their nature as intermediate mass stripped helium stars, more precise measurements and additional properties are needed to serve as benchmarks for detailed evolutionary models.In particular, obtaining bolometric luminosities would allow placement on the Hertzsprung-Russell diagram, stellar radii can inform their current evolutionary stage, and constraints on the stellar winds of stripped stars are important for understanding both the evolutionary past and future.Historically, envelope-stripping of massive stars were predominantly considered via strong stellar winds, but recent measurements of the mass-loss rates the suggested previous evolutionary stage, the red supergiants, are surprisingly low (Beasor et al. 2020).Low mass-loss rates of helium stars would further strengthen the binary-stripping scenario (Beasor & Smith 2022).For the future evolution, the stripped star winds directly affect the amount of hydrogen leftover from interaction and thus the supernova type (Gilkis et al. 2019).They also determine the orbital widening of short-period stripped star + compact object binaries and therefore also their ability to merge in gravitational wave events (Broekgaarden et al. 2022;Stevenson & Clarke 2022).
While full characterization of these stripped star binaries will ultimately require orbital solutions and ultraviolet spectroscopy, here we initiate the effort.We present a de-tailed analysis of the stellar properties of ten stripped stars that dominate over their companion stars even in their optical spectra using atmosphere modelling and spectral fitting.We provide precise measurements of their surface hydrogen and helium content, effective temperatures, surface gravities, stellar radii bolometric luminosities.We further estimate their stellar masses, emission rates of hydrogen-and heliumionizing photons, calculate their Eddington parameters, and estimate rough mass-loss rates via stellar winds.The paper is structured as follows.In Sect.2, we describe the specific sample of stars that we perform spectral fitting of in close detail, while in Sect.3, we describe how the spectra and photometry for this sample were obtained.Section 4 is dedicated to describing a newly computed spectral model grid and the methodology we use to fit the spectra and obtain stellar parameters for the observed stars.We summarize the best-fit properties with associated for the stellar parameters of the stars in Sect.5, while the full spectral fits for the individual stars are presented in Appendix A. In Sect.6, we motivate what evolutionary stage we believe the stars to be in.In Sect.7, we present a rough analysis for obtaining stellar wind mass-loss rate estimates, and in Sect.8 we present estimates for the emission rates of ionizing photons.In Sect.9, we discuss implications of the derived stellar parameters for massive binary evolution, and in Sect. 10 we summarize and conclude our findings.

STELLAR SAMPLE
The full sample of 25 stars presented in Drout & Götberg et al., under review was divided into three spectral groups.Specifically, they were divided based on a comparison of the equivalent widths of He ii λ5411 and Hη/He ii λ3835 lines (chosen to probe the presence of a hot helium star and a Btype MS star, respectively) for the observed stars to a model grid of helium star plus MS star binaries.We found that (i) 8 stars have significant He ii absorption and minimal shortwavelength Balmer lines, consistent with models where the stripped star contributes 80-100% of the optical flux (ii) 8 stars exhibit both He ii absorption and non-negligible shortwavelength Balmer lines, consistent with models where the stripped star contributes 20-80% of the optical flux, and (iii) 9 stars have strong Balmer lines an lack an detectable He ii absorption, only possible in the model grid if any stripped star component contributes < 20% of the optical flux.
In Drout & Götberg et al., under review these were designated Class 1: "Helium-star type", Class 2: "Composite type", and Class 3: "B-type", respectively.Members of all three classes with multiple epochs of spectroscopy showed evidence of radial velocity shifts, indicative of binary motion.While orbital solutions/spectral disentangling will ultimately allow for characterization of the spectral properties of both binary components in the full sample, here we de-scribe the motivation for the subset of 10 objects that we present detailed spectral fits for in this manuscript (Sect.2.1) and review the basic spectral features present in these stars (Sect.2.2).

Sample Selection
Our goal in this first follow-up manuscript is to provide detailed stellar properties for a set of intermediate mass helium stars.We therefore begin by selecting a set of 10 stars where we believe that the stripped star dominates the optical flux and the companion contributes minimally.For this sample, we can therefore adopt a simplified analysis and model the optical spectrum as a single star.Specifically, in this manuscript we will analyze: • The 8 stars of Class 1 from Drout & Götberg et al.,.Of these, stars 1-4 are located in the SMC and 5-8 in the LMC.
• A single object from Class 2 (star 16; located in the LMC).
• An additional star that was originally identified in search for stripped helium stars described by Drout & Götberg et al., under review, but rejected from their final sample based on its kinematics (star 26; likely a foreground halo object).
The optical spectra of these ten stars are displayed in Fig. 1.
We have used the full set of information available to us in assessing that the optical spectrum of a specific star is likely dominated by the flux of a single object.Here we elaborate on each item above.
The Class 1 stars from Drout & Götberg et al., under review all had spectral morphologies consistent with models for "isolated" helium stars, we are able to achieve a good spectral fit assuming contributions from a single star (see Sect. 5).In addition, while they all show radial velocity shifts, they appear as single-lined spectroscopic binaries.This requires that the companion stars are optically faint: either compact objects or low mass main-sequence stars (M 3 M ).However, in Drout & Götberg et al., under review we found that a MS companion star could potentially contribute up to 20% of the optical (V-band) flux and still be classified as a "helium-star type" spectrum.Therefore, in Appendix B we present a set of tests on how the presence of a MS companion may impact the results of our spectral fitting, concluding only minor effects could arise.
While star 16 was placed in the "Composite-type" class by Drout & Götberg et al., under review due to a combination of short-wavelength Balmer lines and He ii absorption, it is most likely an inflated stripped star.When inflated, the surface temperature and surface gravity of a stripped star will decrease, leading to stronger Balmer absorption if any hydrogen remains on the surface.This interpretation is strengthened by the good spectral fit (see Sect. 5 and Fig. 21) and the analysis of its evolutionary stage in Sect.6.It also exhibits radial velocity shifts indicative of a single-lined spectroscopic binary.This is in stark contrast to the other Class 2 objects from Drout & Götberg et al., under review, which (i) we were unable to achieve a reasonable spectral fit for assuming contributions from a single star and (ii) show indications of anti-correlated motion in their He ii and Balmer absorption lines, suggestive of double-lined spectroscopic binaries.
Finally, we address star 26.This object shows a significant UV excess in it spectral energy distribution and has an optical spectrum that would be grouped with the "Helium-star type" class from Drout & Götberg et al., under review due to strong He ii absorption and weak short-wavelength Balmer lines.However, it has a mean radial velocity and proper motions from Gaia DR3 that are sufficiently in-consistent with the bulk of stars in the LMC that we consider it a likely foreground, halo star (see Appendix C for a detailed kinematic assessment).Gaia does not detect a parallax at the 3σ level and we place a lower limit on its distance of ∼3.5 kpc (approximated by taking three times the parallax error provided by Gaia).Analysis presented in Sects.5 and 6 suggest that at a distance of 10 kpc, the properties of star 26 would be consistent with a subdwarf nature.In the rest of the paper, we therefore predominantly adopt the 10 kpc distance for star 26, but also present the stellar properties for the star assuming it is located in the LMC (completeness and for comparison).

Spectral Morphology
The optical spectra for the 10 stars are shown in Fig. 1.All objects show strong He ii absorption, indicative of high temperatures.Stars 1-8 and 26 all show weak short-wavelength Balmer/He ii -blends while star 16 shows stronger features in this regime, consistent with their classifications in Class 1 and 2, respectively, in Drout & Götberg et al., under review.He i lines are present in the spectra of stars 5, 7, 8, 16, and 26, while they are not present in the spectra of stars 1-4 and 6.Stars 1-4 and 6 all show N v lines in emission and/or absorption.Stars 5 and 6 display N iv λ4057 in emission, while in star 26 it appears in absorption.In the case of star 16, N iii lines are visible.Stars 7 and 8 have too poor signal-to-noise ratio spectra for these weak N-features to be detectable.In the case of stars 2, 5, 7, 8, and 26, carbon lines are visible.We will not discuss these further here, but address it in a future study on the CNO abundances of stripped stars.Finally, the Ca ii H & K doublet visible in several of the spectra at 3935 and 3970 Å is interstellar.

OBSERVATIONS
In order to derive detailed stellar properties for the stars described above, we utilize both the moderate-resolution optical spectra and UV-optical photometry.Data acquisition and reduction are described in detail in Drout & Götberg et al., under review.Here we briefly review key details of our methods.

Spectroscopy
We obtained multiple epochs of medium-resolution (R ∼ 4100) optical spectra (λ ∼ 3700 − 7000Å) for the stars detailed in Table 1 using the The Magellan Echellette (MagE) spectrograph on the Magellan/Baade 6.5m telescope at Las Campanas Observatory (Marshall et al. 2008).Spectra were taken during 22 dark/grey nights between December 2018 and February 2022 (PI: Götberg & Drout).Observations were typically taken at the parallactic angle, but on some occasions a rotation was applied to exclude other nearby stars from the slit.This can result in slightly lower signal-to-noise in the blue portion of the spectra (e.g., Star 7; Figure 1).
Initial data reduction was performed using the CarPy python-based pipeline 1 (Kelson et al. 2000;Kelson 2003).The pipeline performs bias/flatfield correction, sky subtraction, 1D spectral extraction, and wavelength calibration.Individual echelle orders were normalized by fitting low order polynomials to the continuum after performing 2.5σ clipping to reject contributions from absorption lines.Orders were then stitched together after normalization.We manually clip artifacts caused by both cosmic rays and by imperfect sky subtraction in cases where stars are located in bright/clumpy H ii regions (e.g., Star 6).Finding the true continuum is challenging, especially for the upper Balmer series (λ 3900Å), and we therefore carefully flatten each spectrum manually and exclude members of the Balmer series above Hδ in our analysis.We note that artifacts could be present in our final spectra that relate to slight variation of the continuum in the wings of broad lines or averaged spectra where the orders overlap.However, we do not consider that these artifacts are sufficiently large to significantly impact our results.
Finally, to produce the highest signal-to-noise ratio (SNR) spectra of each stars, we stack together observations taken on different occasions.However, all the stars considered here display radial velocity shifts and appear as single-lined spectroscopic binaries.We therefore must correct for binary motion when stacking spectra obtained days to months apart.This process is discussed in detail in Drout & Götberg et al., under review.The SNR is then calculated per pixel within the wavelength ranges, 4230-4300, 4400-4430, 4730-4830, and 5030-5250Å, and then averaged, resulting in final SNRs of our combined spectra ranging from ∼30-120 (see Table 1).These combined spectra are show in in Fig. 1 and will be made publicly available upon publication of this manuscript.Stars 1-8 and 16 were originally published in Drout & Götberg et al., under review, and we have now made Star 26 available as well.

Photometry
In this manuscript, we utilize photometry of the stars in our sample in 3 UV and 4 optical photometric bands:UVW2, UV M2, UVW1, U, B, V, and I. Specifically, this data is used to estimate the bolometric luminosity and extinction of each star by fitting magnitudes computed for the best-fit spectral models to the observed photometry.The optical photometry for all sources comes from the Magellanic Cloud Photometric Survey (Zaritsky et al. 2002(Zaritsky et al. , 2004)).Originally, these data are presented in the Vega magnitude system.We calculate zeropoint offsets to convert these to AB magnitudes by performing synthetic Vega and AB photometry on a subset of the stripped star models in our synthetic grid (described below) in order to minimize systematics due to the underlying spectral shape of the star.For a range of stripped star models, the resulting zeropoints vary by significantly less than the catalog magnitude uncertainties (<0.001 mag).The UV photometry was performed on images from the Swift-UVOT Magellanic Cloud Survey (Siegel et al. 2015b;Hagen et al. 2017) as described in Drout & Götberg et al., under review and Ludwig et al. in prep.In particular, to mitigate the effects of crowding in the Swift images, we performed forced point-spread-function photometry at the positions of the optical sources using the forward-modelling code The Tractor (Lang et al. 2016).Final magnitude calibration was then performed using standard HEASARC routines, and multiple observations the same source were averaged.
All photometric data that are used in this study are presented in Table 2. UV photometry for stars 1-8 and 16 were originally published in Drout & Götberg et al., under review, and we have now added magnitudes computed via the same method for star 26.

SPECTRAL FITTING
To obtain stellar properties for the stars in the spectroscopic sample, we compute a grid of spectral models and adopt a χ 2 minimization technique to identify the best-fit model and associated errors.Below, we describe these steps in detail.

Spectral model grid
We used the publicly available 1D non-LTE radiative transfer code CMFGEN (Hillier 1990;Hillier & Miller 1998) to compute a grid of stellar atmosphere models that we can use for spectral fitting and obtain properties of the stars in our spectroscopic sample.A subset of the models described here were used in Drout & Götberg et al., under review to estimate the effective temperature, surface gravity, and surface hydrogen mass fraction of stars 1-8 via a set of equivalent width diagnostics.We have now expanded this grid to cover a larger parameter space to aid in our spectral fitting.Below we describe the grid and computation method in detail.The signal-to-noise ratios are calculated per pixel for the stacked spectra, rounded to the nearest multiple of ten and then averaged over the wavelength ranges 4230-4300, 4400-4430, 4730-4820, 5030-5250Å.
Table 2. Photometric data with 1σ errors obtained from the UBVI survey at the Swope telescope (Zaritsky et al. 2002(Zaritsky et al. , 2004) ) and photometry perfomed on the Swift/UVOT images of the Magellanic Clouds (Siegel et al. 2014(Siegel et al. , 2015a(Siegel et al. , 2019) )  These spectral models are based on those presented in Götberg et al. (2018), which in turn stem from Groh et al. (2008) and the openly available O-star grid on the CMFGEN website 2 .For these models, we include the elements H, He, C, N, O, Si, and Fe.We compute the model spectra between 50 and 50,000 Å. Depending on the density of the wind, we adopt a suitable extent of the atmosphere, which is between 6 and 1000 times the surface radius.We use a minimum number of mesh points of 40, but up to more than 100, together with 15 core rays.
The CNO abundances originate from layers that once were part of the convective main-sequence core, and thus have experienced complete CNO processing.In the structure models of Götberg et al. (2018), the nitrogen and oxygen abun-dances have a rough constant level from the surface to the convective helium-burning core, while the carbon abundance increases by roughly a factor of three from the surface in to the hydrogen-free layer.This larger change in carbon is balanced by oxygen.However, because oxygen is more abundant, the fractional abundance change of oxygen is not prominent.Here, we refrain from a detailed analysis of possible variations of CNO abundances, which will be the topic of a future study.We note that none of the metal lines are used in our spectral fitting process (see Sect. 4.2).
To create a spectral model grid that easily can be scaled to the desired radius or luminosity, we fix the stellar radius in the models to 0.5 R .While this radius is typical for the expectations of envelope-stripped stars (see Tables 1 and B.1-B.3 of Götberg et al. 2018), we note that we scale the spectral models during the fitting procedure so that the radius is a free parameter.We let the luminosity adapt to the assumed radius and temperature, resulting in L bol ∼ 200 − 22, 000 L .We set the code to match the input temperature, radius and surface gravity at an optical depth of τ = 20 (quantities denoted by ), following Groh et al. (2008).We note that these properties are very similar at the photospheric optical depth of τ = 2/3 (quantities denoted by eff), but not exactly the same.Differences between the quantities at τ = 20 and τ = 2/3 are somewhat larger for models closer to the Eddington limit (see below).
Because the stars in the spectroscopic sample lack the typical emission lines originating from stellar winds, we adopt weak, fast, and relatively smooth stellar winds for the models in our primary grid.To do this, we assume mass-loss rates of 10 −9 M yr −1 , terminal wind speeds of 2500 km s −1 , which corresponds to one to several times the surface escape speed as has been measured for massive stars (Lamers et al. 1995), and modest clumping by assuming a volume filling factor, f vol , of 0.5.For the wind velocity profile, we assume a β-law (v(r) = v ∞ (1 − R /r) β ), setting β = 1.In section 7 we will vary these parameters to obtain rough estimates for the mass-loss rates for the stars in our sample.We adopt a turbulent velocity of 20 km s −1 , in common with Magellanic Cloud O-type stars (Ramírez-Agudelo et al. 2017).The impact of turbulence and thermal broadening is negligible for the diagnostic He ii and Hi lines, which are dominated by (Stark) pressure-broadening.There is no evidence for rotational broadening contributing significantly to the Pickering-Balmer lines, although we defer an investigation of rotation rates using metal lines to a future investigation.Before using the models for spectral fitting, we also degrade them to the spectral resolution of MagE using a Gaussian kernel.
The resulting spectral model grid covers most of the intended parameter space, as shown in Fig. 2. The figure shows that the difference between the temperature and surface gravity evaluated at τ = 20 and τ = 2/3 is negligible for most of the models, and at maximum the temperature (T and T eff ) and surface gravity (log 10 g and log 10 g eff ) differ by 10% and 5%, respectively.We encountered numerical convergence issues when high temperatures and low surface gravity are combined, because these combinations approach the Eddington limit (Γ e = 1, see Sect.4.2.4 and the dotted lines in Fig. 2).We note that the Eddington factor is independent of the assumed radius, mass, and bolometric luminosity.Because the spectral morphology changes significantly between 30 and 40 kK, we introduce the 35 kK models for all surface hydrogen mass fractions, and also the 33 and 37 kK models for the surface hydrogen mass fraction X H,surf = 0.01.In total, the grid contains 441 models and we make the full grid publicly available on Zenodo under a Creative Commons Public Domain license: doi:10.5281/zenodo.7976200.Please cite both the present article and the Zenodo dataset when reusing these model grids (Götberg et al. 2023a).

Fitting routine
We employ the χ 2 minimization technique to obtain the best-fit spectral model and models allowed within 1σ deviation for each star.This gives rise to measurements for their effective temperatures, effective surface gravity, hydrogen and helium surface mass fractions, and flux-weighted gravity.We then match the spectral models to the observed photometry to obtain extinction and luminosity, which in turn can be used to calculate the effective radius, spectroscopic mass, and Eddington factor.Finally, we use a set of evolutionary models and our derived bolometric luminosities to estimate evolutionary masses under the assumption the stars are central helium-burning (this assumption is investigated in Sect.6).
Using χ 2 minimization in our rather finely spaced and interpolated grid ensures that the model with the truly smallest χ 2 is found.Because all models within the chosen parameter space are included, the best-fit model will represent the true minimum and not a local minimum.Concerning the errors, artefacts related to the data reduction (Sect.3.1) and implementation of physical processes in CMFGEN (Massey et al. 2013) could mean that the formal 1σ errors we obtain in the χ 2 analysis are slightly underestimated.Below, we describe the details of the adopted fitting procedure3

Treatment of spectral lines
When fitting spectral models to the data, we choose to fit only to certain spectral lines (this is, for example, also done in the fitting procedure in the IACOB survey, see Simón-Díaz et al. 2011).The choice of what lines to fit to is important, because they are affected differently by parameter variations.This is demonstrated in Fig. 3, where we show the effect of varying the surface hydrogen to helium content, the temperature, and the surface gravity, on the four spectral lines He ii λ4100/Hδ, He ii λ4542, He i λ5876, and N v λ4604.For this figure, we start from the parameters T = 70 kK, log 10 g = 5.0 and X H,surf = 0.3 and vary each parameter.
The left panels of Fig. 3 show that the surface mass fraction of helium and hydrogen affect the central wavelength of the He ii λ4100/Hδ line blend along with the strength of He ii λ4100/Hδ, He ii λ4542 and He i λ5876.The effect on the nitrogen line is negligible.The central panels show that effective temperature significantly affects the strength of He i λ5876 and He ii λ4542 for T 70 kK, but these lines are minimally affected for variations at higher temperature.In fact, He i λ5876, the most temperature sensitive He i line in the spectral range, disappears for T >70-80 kK (see also Fig. 4).The nitrogen ionization balance is also sensitive to temperature variations for T > 70 kK.In Fig. 3 N v λ4604 is not present for T ≤ 50 kK, appears in absorption for T = 70 kK, and in emission for T = 100 kK.However, to fully trace these variations of the nitrogen features we would require both higher signal-to-noise spectra and to expand the model grid to vary the surface nitrogen mass fraction.Finally, the right panels of Fig. 3 show that variations in surface gravity affect both the strength and shape of the hydrogenic line transitions of He ii λ4100/Hδ and He ii λ4542.The effect of surface gravity on He i λ5876 and N v λ4604 is moderate.
Summarizing, to probe the parameters of the model grid when fitting the observed spectra, it is important to include (1) both pure He ii and He i lines when possible, since it gives the most accurate temperature determination and (2) a combination of pure He ii and H/He ii blended lines to trace surface hydrogen to helium content.This set will thus also include lines that are affected by Stark broadening and trace surface gravity.In choosing the final set of lines to fit, we avoid fitting to the α lines Hα/He ii λ6560 and He ii λ4686 because of their sensitivity to stellar wind and nebular contamination.This choice differs from analysis of the more luminous Wolf-Rayet and WN3/O3 stars where the α-lines often are used as primary diagnostic lines (Crowther et al. 1995;Neugent et al. 2017).The final set of lines used to fit the spectrum for each star are listed in Table 4 in Appendix A.
We renormalize the continuum for each spectral line individually before fitting with models.This is done by fitting a horizontal line to the ∼10 Å regions on both sides of each line.We then hold the continuum fixed in our χ 2 minimization.We select wavelength range that will be fit for each line by finding where the wings of the observed line first increase above the continuum level of 1 (due to noise fluctuations) on both sides of the central wavelength.
When computing χ 2 for one model, we compute the χ 2 for each line individually and then sum these together, meaning that all lines are weighted equally.Because some lines are narrower than others, this means that these will carry somewhat less importance to the fit compared to broader lines, which are composed of more data points.However, in tests with higher weighted narrow lines, we did not find significant improvements of the fits and therefore choose to not include different line weights.

Interpolating and constraining the spectral model grid
To obtain better fits and finer resolution in the measured parameters, we interpolate the spectral model grid.The interpolation is linear in T , log 10 g and X H,surf .We choose to sample T every 2 kK between 30 and 150 kK, log 10 g every 0.1 steps between 4.0 and 6.0, and X H,surf in steps of 0.05 between 0.05 and 0.7 (in addition to the computed models at 0.01).We do not extrapolate the grid, meaning that the high temperature and low surface gravity corner still is not populated with models (cf. Figure 2).
In addition, we use the presence of various nitrogen and He i lines to help constrain the the temperature range to from the full model grid to consider when fitting each individual star.While He ii and H lines are present throughout the entire model grid, the same is not the case for nitrogen and He i .Specifically, although the strength and detailed line profile of the nitrogen features are dependent on the abundance of nitrogen (which we do not vary in our grid) their presence can provide a sensitive temperature diagnostic at T > 60 kK.As demonstrated in the top panel of  As an example, we show here the models with surface hydrogen mass fraction X H,surf = 0.3, spread out in the temperature -surface gravity plane.Triangular markers show the presence of nitrogen or helium lines that are used to constrain the grid (see also Table 4).Gray circles indicate models in which none of the lines specified by the legend are present.
60 kK (pink triangles).We note that Fig. 4 only shows the part of the grid with surface hydrogen mass fraction X H,surf = 0.3 for illustration purposes, but the line presence only varies slightly with different surface hydrogen mass fraction.
When one or more of the described lines are present in an observed spectrum, we use it to constrain the model grid used in our fitting procedure.The constraints we use for each star are given in Table 4 in Appendix A. We do not use the absence of lines to constrain the model grid since poor signalto-noise ratio or exact nitrogen abundance can affect whether the line is visible.4.2.3.Spectral fitting to obtain T , T eff , log 10 g , log 10 g eff , X H,surf , X He,surf , and L For each star, we calculate the χ 2 of all models in the interpolated and constrained grid and determine which one is the best-fit model by finding the one with smallest χ 2 (desig-nated χ 2 min ).The models with χ 2 < χ 2 min +∆χ 2 are regarded as acceptable models and their properties are used to determine the errors on the fitted parameters.We determine ∆χ 2 by calculating the 68.27% confidence interval based on the number of degrees of freedom.The calculation of ∆χ 2 is done using the python function scipy.stats.chi2.ppf(see however Press et al. 1992).
We use the temperature and surface gravity at τ = 20 and τ = 2/3, along with the surface hydrogen mass fraction of the best-fit model as the best-fit values for these parameters (T , T eff , log 10 g , log 10 g eff , and X H,surf ).For the 1σ errors on these parameters, we use the maximum and minimum values among the models that fulfil χ 2 < χ 2 min + ∆χ 2 .Two more stellar parameters can be derived directly from these model fits.First, the surface helium mass fraction, which is simply X He,surf = 1 − X H,surf − Z.This corresponds to X He,surf = 0.98547, 0.89547, 0.69547, 0.49547, and 0.29547 for X H,surf = 0.01, 0.1, 0.3, 0.5, and 0.7 (see Sect. 4.1 for information on Z).Second, the inverse of the flux-weighted gravity, L ≡ T 4 eff /g (Kudritzki et al. 2003;Langer & Kudritzki 2014), can be calculated for each model and thus also determined using the χ 2 method outlined above.We present L in solar units, L , calculated assuming T eff, = 5, 777 K and g = 27, 400 cm s −2 .Note that the inverse of the fluxweighted gravity is very sensitive to uncertainties in the effective temperature, due to the fourth power in its definition.4.2.4.Obtaining L bol , A V , R eff , M spec , and Γ e In order to determine bolometric luminosities, we fit the spectral energy distributions (SEDs) of the acceptable models to the observed photometry of each star, including extinction as a free parameter.For each spectral model, we scale the spectrum to produce a range of bolometric luminosities between roughly 1 and 106 L .We then apply a range of extinction values between A V = 0 − 1.5 mag separated in steps of 0.01 mag, adopting the extinction curves from Gordon et al. (2003)  4 .For simplicity, we only adopt the average extinction curve for each of the Magellanic clouds, and do not explicitly include a separate Milky Way foreground component in the fitting.While the LMC and Milky Way extinction curves are comparable in the wavelength regions of interest, we discuss any impact of differences in the shape of the SMC and Milky Way curves in the ultraviolet in Sect. 5.The exception for this approach is star 26 evaluated at 10 kpc distance, where we only adopt the Milky Way extinction curve.We calculate the AB magnitudes of each resulting model in the Swift UVW2, UVM2, UVW1, and optical UBVI bands using the filter functions from the SVO filter 4 We employ the functions averages.G03 LMCAvg and averages.G03 SMCBar of the python package dust extinction for this calculation (https://dust-extinction.readthedocs.io/en/stable/).service5 (Rodrigo et al. 2012;Rodrigo & Solano 2020).We then calculate the chi-square statistic for the resulting modeled magnitudes compared to the observed photometric data, adopting distances of 50 kpc to the LMC (Pietrzyński et al. 2013) and 62 kpc to the SMC (Graczyk et al. 2020)  6 .Because extinction has larger influence in the UV compared to the optical, we prefer to use the described method fitting to photometry, rather than for example assessing flux calibrated optical spectra, which furthermore often have larger systematic uncertainties in absolute calibration.
We apply the above procedure to all models that fall within the χ 2 < χ 2 min + ∆χ 2 threshold from the spectral fitting (Section 4.2.3),resulting in a range of L bol and A V values for each star.(Because the photometric errors are small, we simply find a single best-fit value of these parameters for each spectral model.)For each star, we adopt the L bol and A V found for the best-fit spectral model from Section 4.2.3 as our baseline values.Errors are determined based on the minimum and maximum values found from fitting the larger sample of models accepted within 1σ from the spectral fitting.
For each model, we compute the effective radius using the bolometric luminosity and effective temperature following the Stefan-Boltzmann's law (L bol = 4πR 2 eff σT 4 eff ) and the spectroscopic mass by combining the surface gravity and effective radius (g eff = GM spec /R 2 eff ).As with extinction and bolometric luminosity, for each star we adopt the effective radius and spectroscopic mass found from the best-fit spectral model as our baseline values.Quoted errors similarly correspond to minimum and maximum values found from all models within 1σ based on the spectroscopic fit.
With the bolometric luminosity and spectroscopic mass, we can also estimate the Eddington factor for Thomson scattering, Γ e , which describes how close the star is to the Eddington limit (Gräfener et al. 2011).The Eddington factor is defined as follows where c is the speed of light, G is the gravitational constant, and κ e is the electron scattering opacity, defined as 4.2.5.Estimating the evolutionary mass, M evol Finally, we estimate the evolutionary masses for the stars in our sample using the relation between mass and luminosity for stripped stars that have reached half-way through central helium burning, defined as when X He,center = 0.5.To find this relation, we use the evolutionary models of Götberg et al. (2) This mass-luminosity relation should be a decent approximation for the mass-luminosity relation throughout heliumcore burning, since it does not significantly change during this phase.This is demonstrated in Fig. 5 where we use shaded background to show the variation in these parameters for central helium mass fractions between 0.8 and 0.1.However, we note that this definition of the evolutionary mass assumes that the stripped stars are in the phase of helium-core burning and are not currently contracting or expanding (cf.Laplace et al. 2020).We emphasize that using this relation to estimate the mass for stripped stars that are inflated may lead to an overestimated evolutionary mass.We will directly assess this for stars in our sample (e.g., for star 16) in Sect.6.
The models of Götberg et al. (2018) reach stripped star masses of ∼ 7.2 M and bolometric luminosities up to ∼ 10 5 L .As, in particular, star 1 could reach higher values, we allow for extrapolation of the mass-luminosity relation.

STELLAR PROPERTIES
In Table 3, we present the stellar properties that we obtain following the method described in Sect. 4. We show the fit for star 1 as an example in Fig. 6, while the fits for the other stars are presented in Appendix A. The top left panels of the figure shows the spectral lines used for the spectral fit.The observed spectrum is shown in black, while the thick colored lines indicate the best-fit spectral model.Other models acceptable within 1σ are shown as thin colored lines.1.17The top right panels show χ 2 as function of the effective temperature, surface gravity, and surface hydrogen mass fraction.The best-fit model (with the minimum χ 2 ) is shown as a big colored circle, while the models acceptable within 1σ are marked with smaller colored circles below the black line labeled 1σ.The models marked with gray dots are not acceptable within 1σ.As seen in these panels, none of the stars exhibit any ambiguity regarding where the true minimum and thus best-fit model lies.
The two middle panels show the normalized observed spectrum in black and the best-fit model overplotted in a thick colored line.The spectral lines used for the spectral fit are marked by shaded background.The bottom left panel shows, in black, the observed photometric data in AB magnitudes and centered on the central wavelengths of each filter.The best-fit model is shown in a thick colored line and large colored circles, while the models allowed within 1σ are plotted with thin lines.
Finally, the derived best-fit effective temperature and bolometric luminosity with associated errors are plotted using color in a Hertzsprung-Russell diagram at the bottom right.The models allowed within 1σ are shown using black dots.For reference, we also plot evolutionary tracks for a sequence of stripped star models from Götberg et al. (2018) using gray lines.These evolutionary models are for stripped stars with masses 1.5, 1.9, 2.5, 3.4, 4.5, 5.9, and 7.3 M , corresponding to initial masses of 5.5, 6.7, 8.2, 10, 12.2, 14.9, and 18.2 M .
In the remainder of this section, we summarize and discuss the stellar parameters found for the 10 stars in our spectroscopic sample.In several instances, we compare with the evolutionary models from Götberg et al. (2018).Work presented in this manuscript suggests that the observed wind mass loss rate (see Sect. 7) is lower compared to what we assumed for the evolutionary models.However, although winds are important for the spectral morphology and future evolution of stripped stars, winds only mildly affect their broad surface properties (Gilkis et al. 2019).
Effective temperature -We measure effective temperatures above 50 kK for all but one star.The best-fit effective temperatures are in the range 50 − 95 kK for stars 1, 2, 3, 4, 5, 6, 7, 8, and 26.Star 16 is somewhat cooler, with about 35 kK.The tightest constraints on the effective temperature can be made when both He i and He ii lines can be included in the spectral fit (see Sect. 4.2).However, for the hottest star (star 1) that does not display He i lines, the effective temperature can be well-constrained using the H and He ii lines alone, because of the high signal-to-noise ratio.In other cases where He i lines are not present (stars 2, 3, and 6) and/or when the signal-tonoise ratio is lower (stars 3, 4, 7, and 8), we obtain large, sometimes asymmetrical errors for the effective temperature.This occurs because the He ii lines have poor constraining power at high temperatures.
Surface gravity -We find typical surface gravities of log 10 g eff ∼ 57 -well above those of regular main-sequence stars, which are log 10 g eff ∼ 3.5 − 4.5, but below values for white dwarfs (log 10 g eff ∼ 6 − 9).Stars 5 and 16 have somewhat lower surface gravities, with log 10 g eff of about 4.5 and 4.2 respectively.The derived surface gravities for stars 3 and 26 are somewhat higher, with log 10 g eff of 5.4 and 5.7 respectively.We note that our obtained errors for surface gravity may be somewhat underestimated since it is challenging to identify the precise continuum adjacent to the broad Balmer and Pickering lines With constraints on effective temperature and surface gravity, the stars can be placed in Kiel diagrams, as shown in panels a) and b) of Fig. 7. Comparing to the Kiel diagram presented in Drout & Götberg et al., under review based on estimates of effective temperature and surface gravity using equivalent width diagnostics, this updated version is similar, illustrating the power of equivalent width analysis.In all panels of Fig. 7, we show the evolutionary tracks of donor stars in binary systems presented by Götberg et al. (2018).These models have initial masses of 4.5, 7.4, 9.0, 12.2, and 18.2 M , which results in masses of the stripped stars of 1.1(1.2),2.0(2.2),2.7(2.9),4.1(4.5),and 7.2(7.3)M for the LMC(SMC).We use the models with Z = 0.006 and Z = 0.002 to represent the LMC and SMC, respectively.We display the stars in the LMC using circles and the stars in the SMC with squares.Star 26 is displayed using a diamond.The figures show that stars 1-8 and 26 agree well with being helium-core burning stars stripped of their hydrogen-rich envelopes through mass transfer in binary systems.This can be seen by comparing their locations in the Kiel diagram to the binary evolution tracks that we have displayed for reference.Star 16 appears to be more inflated than typical helium-core burning stripped stars.
Inverse of flux-weighted gravity -For the inverse of the fluxweighted gravity, we obtain values of log 10 (L/L ) ∼ 2.5 − 4.5.Since the inverse of the flux-weighted gravity behaves as a luminosity, we create spectroscopic Hertzsprung-Russell diagrams in panels c) and d) of Fig. 7 using this quantity and the effective temperature.In this diagram, we see that all stars agree well with being donor stars stripped of their hydrogen-rich envelopes since they overlap with the expected location for stripped stars from the evolutionary models.Also in the spectroscopic Hertzsprung-Russell diagrams, the stars agree well with being central-helium burning stars, apart from star 16, which appears to be somewhat cooler than typical helium-core burning stripped stars.
Figure 7.The derived properties with associated errors for the spectroscopic sample shown with numbered markers plotted together with binary evolutionary models for donor stars in binary systems (Götberg et al. 2018).Stars in the LMC are marked using circles, stars in the SMC with squares, and the foreground object with a diamond.From top to bottom we show the Kiel diagram, the spectroscopic Hertzsprung-Russell diagram, and effective radius as function of effective temperature.The left panels are for the Large Magellanic Cloud and the right panels for the Small Magellanic Cloud.The evolutionary models are for stars with initial masses of 4.5, 7.4, 9.0, 12.2, and 18.2 M , with corresponding stripped star masses of 1.1(1.2),2.0(2.2),2.7(2.9),4.1(4.5)and 7.2(7.3)M for the Large(Small) Magellanic Cloud.The central helium burning is marked with a thicker and darker line and the evolutionary tracks are cut at central helium depletion.
Surface hydrogen and helium mass fraction -The best-fit surface mass fraction of hydrogen is well below what is expected for stars with hydrogen-rich envelopes, such as main-sequence stars.Five stars (star 1, 2, 4, 6, and 16) have surface hydrogen mass fractions between 0.3 and 0.4, while the remaining five stars (star 3, 5, 7, 8, and 26) have surface hydrogen mass fractions between 0 and 0.1.Conversely, the surface helium mass fraction for these two groups correspond roughly to between 0.6 and 0.7 and between 0.9 and 1.It is likely that three Extinction -We find small values for the extinction, between A V = 0.1 and 0.7 mag.Generally, we find lower extinction values for the stars located in the SMC (A V ∼ 0.1 − 0.4 mag) compared to those located in the LMC (A V ∼ 0.2 − 0.7 mag).These values agree with the low end of the distributions found for stars in the Magellanic Clouds by Zaritsky et al. (2002Zaritsky et al. ( , 2004)).This is expected since the stars were identified through their UV excess, meaning that our spectroscopic sample would be biased against stars whose sightlines are strongly affected by dust extinction.
Indeed, for a few stars (e.g.star 4 and star 8) the extinction values values are consistent with the expectation for foreground Milky Way extinction (Schlafly & Finkbeiner 2011), implying negligible internal extinction in the SMC/LMC, respectively.On this point, we note that the extinction curves we employ (Gordon et al. 2003) are averages over the Magellanic Clouds.They do well in representing the extinction curves for our observed sample as seen from the photometric fits, although the foreground should be better represented by a Milky Way average extinction curve.While the LMC and Milky Way extinction curves are similar over the wavelength regions we consider (Gordon et al. 2003), differences exist in the UV for the SMC.To ensure that the stellar parameters that depend on the extinction estimate are robustly estimated, we run the spectral fitting routine on the SMC star 4 using an average extinction curve for the Milky Way (Gordon et al. 2009), which, in contrary to the SMC curve, contains the bump around 2175Å.Despite this significant difference, we obtain estimates for the stellar parameters that are negligibly different from those obtained when using the SMC extinction curve.
Bolometric luminosity -The bolometric luminosities that we infer from the model fits are between 10 3 and 10 5 L .This range is typical, for example, for main-sequence stars with masses between ∼5 and ∼30 M (Georgy et al. 2013).The bolometric luminosity determination is sensitive to how well the effective temperature is determined since the peak of the spectral energy distribution is located in the un-observable ionizing regime and needs to be inferred from the shape of the modeled spectral energy distribution.This dependency is reflected in the larger errors on bolometric luminosity when the effective temperature also has larger errors (for example, see star 4, Figure 16).The bolometric luminosity is also dependent on the distance.This is not an issue for stars 1-8 and 16, which are members of the Magellanic Clouds, but affects star 26, which has a more uncertain distance.When placed in the Hertzsprung-Russell diagram in Fig. 8, it is again clear that the stars in our spectroscopic sample are poorly matched with main-sequence stars.Instead, they overlap with the helium main-sequence.The exception is again star 16, which instead appears to overlap with an inflated phase.The assumed 10 kpc distance of star 26 as displayed in Fig. 8  Effective radius -The effective radii we derive are well constrained and all close to 1R , spanning a range from 0.3 R to 1.4 R .Within the uncertainties, none of the stars exceed 1.6 R , suggesting that they are indeed much smaller than typical main-sequence stars with the same temperatures -the massive O-stars having radii 10 R .The measured radii agree well with predictions from binary stellar evolution models (0.6−1.4 R for stripped stars with masses between 2 and 7.2 M , Götberg et al. 2018).This can also be seen from panels e) and f) of Fig. 7.As shown in Table 3, star 26 has an estimated radius of 1.4 R when assumed to reside in the LMC, compared to 0.3 R when assumed at a distance of 10 kpc.Given its high surface gravity, the smaller size is more compelling, and in agreement with the star being located in the foreground.
Spectroscopic mass -We find spectroscopic mass estimates between 0.8 and 6.9 M for stars 1-8 and 16.For stars where we have very good model fits, such as for star 1, the errors in the spectroscopic mass are only ∼ 20%.For fits with larger uncertainties, such as for star 8, the errors are very large, reaching a factor of 10.Star 26 has an estimated spectroscopic mass of 38 M when assumed to reside in the LMC, but instead the more realistic 1.5 M when placed at 10 kpc distance.
Evolutionary mass -The evolutionary mass provides an additional handle on the stellar mass.On average, we find somewhat higher evolutionary masses than spectroscopic masses, stretching from 1.2 to 8.4 M .Among the sample, all but Figure 9.Comparison of the spectroscopic and evolutionary masses for the stars in the spectroscopic sample.The lines at 2.5 M are meant as approximations for the limit for stripped stars that reach core collapse vs evolve to white dwarfs.stars 8, 16 and 26 have evolutionary masses above 2.5 M , which can be used as an approximation for the boundary for what stars will undergo core collapse (Tauris et al. 2015).
We plot the evolutionary mass versus the spectroscopic mass found from our analysis in Fig. 9.The figure shows that the best constrained spectroscopic masses belong to stars with either high SNR (star 1) or spectra with both He i and He ii lines present (stars 5, 16, and 26, however not stars 7, or 8, likely because of their low SNR).We note that star 16 appears inflated (see above) and its mass may be poorly represented by the mass-luminosity relation we adopt when calcuating evolutionary mass (see Sect. 4.2.5).Dynamically inferred masses would be ideal to use for resolving what the true stellar masses are.
Eddington factor -We estimate that the stars in the spectroscopic sample have bolometric luminosities that mostly are far from their Eddington limits.Star 1 and star 5 are the closest to their Eddington limits, with Eddington factors of ∼ 0.4 and ∼ 0.25, respectively.The other stars all have Eddington factors of Γ e ∼ 0.006-0.15.The Eddington factors we find are quite similar to those of O-type stars (Lamers & Leitherer 1993).

EVOLUTIONARY STAGE: CONTRACTING, HELIUM-CORE BURNING, OR EXPANDING?
Stripped stars burn helium in their centers during the large majority of the remaining stellar lifetimes after envelopestripping is complete.Unlike the central hydrogen burning during the main-sequence, the radii of stripped stars only moderately change during the central helium burning phase (e.g., Götberg et al. 2019).There are, however, two shorterlasting inflated stages predicted for stripped stars.First, the contraction phase after envelope-stripping is complete, and, second, the expansion phase initiated after helium-core depletion (Laplace et al. 2020).
We show these evolutionary phases in Fig. 10, using the binary evolution models of Götberg et al. (2018).In the figure, we plot the radii of models of stripped stars with masses ∼ 1 − 7 M (corresponding to initial masses ∼ 4.5 − 18.2 M ) as function of their bolometric luminosity.The models are represented by solid black lines and arrows that demonstrate the evolutionary direction.In the top panel, we plot the contraction phase followed by the helium-core burning phase until the star reaches its minimum radius, while in the bottom panel we show the expansion phase during helium-shell burning, from the point where the star has reached its minimum radius, until death or the model evolves off the plot.We use dark gray background for the tracks to mark the central helium burning, which here is defined as when the central mass fraction of helium is between 0.9 and 0.01.The blue and red shading is used to show what fraction of the temporal duration of the stripped star phase has passed.Comparing the color shading with the dark gray background of the tracks, it is clear that central helium burning indeed coincides with the majority of the stripped star duration, while contraction and expansion correspond to about 10% and 1-5% of the stripped star phase, respectively.Thus, we expect that most stripped stars should be helium-core burning.
Figure 10 also shows that the radius change during central helium burning is somewhat mass dependent, with a larger change for the more luminous, higher-mass stripped stars.For example, we expect that a 7 M stripped star with L bol ∼ 10 5 L can have radii between ∼0.7 and 5 R during central helium burning, while a 3 M stripped star, with L bol ∼ 10 4 L , should be limited to radii between ∼0.6 and 1.5 R in the same evolutionary phase.The reason is twofold: first because more massive stars ignite helium in their cores earlier during the evolution, and second because of wind mass-loss, which allows deeper, more compact layers of the stellar models to be revealed (cf.Gilkis et al. 2019).We note that the binary evolution models we use were created for stars stripped via stable mass transfer, which leaves a layer containing hydrogen on the stellar surface (Götberg et al. 2017;Laplace et al. 2020).Stripped stars with no hydrogen layer are expected to be more compact and smaller than stripped stars that retain hydrogen (Yoon et al. 2017).
We overplot the stars in our spectroscopic sample in both panels of Fig. 10.All stars overlap with expectations for the central helium burning stage, apart from star 16.While it is possible that the stars are during the early stages of expansion, the different timescales make the helium-core burning stage more likely.More precise measurements for the stellar  2018), labeled by stripped star mass.We show the fraction of the stripped star duration using blue and pink shades and the central helium burning phase when 0.9 > X He,c > 0.01 using dark gray background for the evolutionary tracks.The stars in the spectroscopic sample are plotted using their effective radii and bolometric luminosities with numbered markers (see Table 3).The top panel shows that contraction lasts ∼ 10% of the stripped star duration, while the bottom panel shows the expansion phase lasts ∼ 1 − 5%.All stars but star 16 agree with the helium-core burning phase and the expansion phase, while star 16 could either be contracting or expanding.
masses than what we currently have could be used to determine the evolutionary stage more accurately.As an example, according to the models displayed in Fig. 10, star 1 could either match a helium-core burning star with mass ∼ 8 M or a ∼ 5 M expanding stripped star.Similarly, star 5, for example, matches either a ∼ 4 M helium-core burning stripped star or a ∼ 3 M expanding stripped star.
Star 16 is about twice as large compared to what is expected for helium-core burning stripped stars with its determined bolometric luminosity.We, therefore, consider that star 16 likely is experiencing an inflated stage (cf.Schootemeijer et al. 2018), which agrees with its lower surface gravity and lower effective temperature compared to the other stars in the sample (see Fig. 7 and Sect.5).Whether the star is in the contraction or expansion phase is not evident from current data: contraction stages should be slower and thus more common, but expansion phases should be brighter, favoring their detection (see Schootemeijer et al. 2018).Again, more precise mass measurements will provide insight in what evolutionary stage star 16 is in.
Even though we do not know the distance to star 26 very accurately, Figure 10 suggests that the star is likely a heliumcore burning subdwarf with mass of ∼ 1 M , demonstrated by the closeness to that evolutionary track.Especially its effective temperature also matches such a massive subdwarf scenario better than either that of a typical subdwarf Bstar or a helium-core burning stripped star in the LMC (cf.Götberg et al. 2018).If star 26 would have been located in the LMC (which would also require that it was a runaway star; Appendix C), it would overlap with an inflated stage (see Table 3), which does not match well with its high surface gravity.The 10 kpc distance we adopt here gives rise to a bolometric luminosity, stellar radius and spectroscopic mass that roughly match the expectations for a helium-core burning stripped star with the effective temperature of star 26 (Götberg et al. 2018), also accounting for the complete loss of hydrogen, which likely results in the slightly higher surface gravity and effective temperature.It is worth to note that star 26 has a significantly higher temperature (T eff > 50kK) than typical subdwarf B type stars (T eff ∼ 25kK), and is in fact much more similar to the ∼ 1.5 M subdwarf in the Galactic binary HD 49798 (Mereghetti et al. 2009;Brooks et al. 2017).

CONSTRAINTS ON STELLAR WIND MASS-LOSS
In contrast to the original spectral models created for stripped stars by Götberg et al. (2018), the stars in our spec-troscopic sample do not show any strong/broad emission lines indicative of mass loss through stellar winds.However, it is possible that some wind is driven off the surfaces, for example through metal line driving and radiation pressure.The somewhat higher Eddington factors for stars 1 and 5 (see Table 3), for example, suggest some contribution from radiation pressure to the wind driving, and these stars could therefore perhaps have somewhat higher wind mass loss rates than the other stars.While ultraviolet spectroscopic will ultimately provide the most precise measurements of the wind properties from these stars, here we investigate what rough constraints can be placed from the optical spectra alone.
As seen in Fig. 1, the optical spectra contain only absorption features with the exception of weak N iv and N v emission lines.While these nitrogen lines may occur in emission, they are, in these cases, not signs of a stellar wind, instead the result of photospheric level inversion (cf. Rivero González et al. 2011, 2012).This is also clear from their narrow widths, which are not expected for the fast speed that is necessary for stellar winds to escape the surface of the compact stripped stars ( 1, 000 km s −1 ).In fact, for example, when the N v λλ 4604/20 doublet appears in emission, it is most likely because of high surface temperature causing the upper level to be pumped ( 90kK, see Figs. 3 and 4).
The lines that are most sensitive to wind mass-loss in the optical spectrum are Hα and He ii λ4686, since they are both α-lines (cf.e.g., the WN3/O3 stars discovered by Massey et al. 2014;Neugent et al. 2017, which show moderate wind mass-loss).Because Hα is very sensitive to contributions from surrounding H ii regions, we choose to focus on the effect of winds on He ii λ4686 to very roughly estimate the wind mass-loss rate of the observed sample of stars.
To estimate wind mass-loss rates, we take the best-fit spectral models for each star following the parameters presented in Table 3, and then compute new versions of these models assuming a range of wind mass-loss rates ( Ṁwind = 10 −10 , 10 −9 , 10 −8 , 10 −7 , and 10 −6 M yr −1 ), while fixing the terminal wind speed (v ∞ = 2500 km s −1 ), the amount of wind clumping ( f vol = 0.5), and the wind velocity profile (β = 1).While the wind speed is uncertain, we adopt 2500 km s −1 because it matches reasonably well with the ratio between terminal wind speed and surface escape speed, v esc , for massive O-stars, which is v ∞ /v esc ∼ 2.5 (Lamers et al. 1995).This ratio also matches reasonably well with the expectations for subdwarfs that was computed by Krtička et al. (2016) and the computed values for a range of helium star masses of Vink (2017).We estimate the surface escape speeds for the stars using the derived parameters (v esc = 2GM spec /R eff ) and present the values in Table 3.
After computing the spectral models with varying wind mass-loss rate, we find the upper limit for wind mass-loss rate acceptable for each star by identifying, by eye, the model with the highest wind mass-loss rate that still matches the line shape of He ii λ4686.This comparison is plotted in Fig. 11, where we show the observed spectra in black and the models with mass-loss rates 10 −10 , 10 −9 , 10 −8 , 10 −7 , and 10 −6 M yr −1 in yellow, green, blue, purple, and red, respectively.The left panels show a zoomed-out version displaying the development of wind emission, while the right panels show the detailed comparison between the models and the data.All wind mass-loss rates were not computed for all models.The 10 −10 M yr −1 models exist for stars 7 and 16, and the 10 −6 M yr −1 model exists for star 1.The reason is that the lowest wind mass-loss rate models are cumbersome to converge numerically and the highest wind mass-loss rate model was not necessary for other stars than star 1.
We find that stars 1 and 5 have some in-filling in He ii λ4686, suggesting there could be a stellar wind affecting the optical spectra.This aligns well with their somewhat higher Eddington factors of Γ e ∼ 0.38 and ∼ 0.26, respectively (see Table 3).The model with mass-loss rate 10 −7 M yr −1 and 10 −8 M yr −1 match best the He ii λ4686 line for star 1 and star 5, respectively.We, therefore, adopt these values as a rough mass-loss rate estimate for stars 1 and 5.For the remaining stars, no line-infilling is evident and all spectral line shapes are well-matched by the wind mass-loss rate models with Ṁwind = 10 −9 M yr −1 .We therefore adopt 10 −9 M yr −1 as the upper limit for the wind mass-loss rate for the remaining stars.In the case of star 7, it appears that the 10 −10 M yr −1 model produces a too deep spectral feature, therefore we do not consider the 10 −9 M yr −1 an upper limit for star 7, but a rough estimate.These low massloss rates match well given the lower Eddington factors of Γ e ∼ 0.04 − 0.15 for stars 2, 3, 4, 6, 7, 8, and 16, suggesting that wind driving from radiation pressure is small.Star 26 may be an exception, because we cannot distinguish between the 10 −9 and 10 −8 M yr −1 models and therefore adopt 10 −8 M yr −1 as an upper limit.However, we note that for this analysis, we adopted the stellar properties that correspond to membership of the LMC for star 26.We provide these rough estimates for the wind mass-loss rates in Table 3.We emphasize that the method we employ is approximate since the fixed wind parameters also influence the line shapes, although perhaps less than the wind mass-loss rates, within reasonable ranges.
The wind mass-loss rate of stripped stars is thought not only to change the spectral morphology, but primarily to affect the properties and future evolution of the stripped star (Yoon et al. 2017;Götberg et al. 2017;Gilkis et al. 2019;Laplace et al. 2020).Because of the lack of observed stripped stars, it has been difficult to construct a suitable wind massloss prescription.From the analysis of the Galactic quasi Wolf-Rayet star in HD 45166 (Groh et al. 2008), it previously appeared as if an extension of the empirical Wolf-Rayet Figure 12.Rough estimates for the mass-loss rate upper limits (and tentative number in the case of stars 1, 5 and 7) plotted as function of bolometric luminosity for the stars in the sample using colored and numbered symbols (because the symbols for stars 2 and 4 are behind other markers, we label them above).We also plot the massloss rate prescriptions from Nugis & Lamers (2000), Krtička et al. (2016), Vink (2017), andSander &Vink (2020) in beige, brown, light gray, and dark gray.We do not extrapolate the Krtička et al. ( 2016) scheme above 10 4 L since these models were created for subdwarfs.For the Sander & Vink (2020) scheme, we only show Z = 0.006 since the lower metallicity predictions are beyond the parameter space of the plot.wind mass-loss scheme of Nugis & Lamers (2000) was appropriate.However, a weaker wind prescription, for example, the one made for subdwarfs by Krtička et al. (2016) could also be accurate.Recently, efforts have been made to improve our understanding of wind mass loss from helium stars, in particular with the single-temperature models from Vink (2017) and the high-mass helium star models from Sander & Vink (2020).Interestingly, these studies predict lower wind mass-loss rates than what is expected from extrapolated Wolf-Rayet wind mass-loss schemes.Anticipating the results from these teams' ongoing theoretical efforts, we hope to provide a tentative, yet useful, comparison.
For radiation driven winds, mass-loss rate prescriptions are often described as luminosity dependent (see for example the review by Smith 2014).We, therefore, plot the estimates for wind mass-loss rates as function of the bolometric luminosity for the observed sample in Fig. 12.To compare, we also display the predictions from Nugis & Lamers (2000), Krtička et al. (2016), Vink (2017), and Sander & Vink (2020).For these, we adopt, when possible, surface helium mass fractions between 0.4 and 1, metallicity between 0.002 and 0.006, and effective temperature between 50 and 100 kK.These ranges result in the broad, colored bands that we display in Fig. 12.
Figure 12 shows that the mass-loss rate estimates from our observations are low compared to most schemes.None of the stars match the extrapolation of the Wolf-Rayet scheme from Nugis & Lamers (2000), and the massive helium star scheme from Sander & Vink (2020) does, understandably, not extend to sufficiently low luminosities.Stars 1, 5, 8, and 16 appear to agree with the predictions from the Vink (2017) scheme, but stars 2, 3, 4, 6, and 7 appear to have significantly lower mass-loss rates, resulting in a poor match.The flattening of the subdwarf prescription from Krtička et al. (2016) appears to better represent the low mass-loss rates of stars 2, 3, 4, 6, 7, 8, and 16, but it could be that the actual wind massloss rates are even lower than the expectations from this prescription.We also note that the prescription of Krtička et al. (2016) was fitted to data with L bol < 10 4 L and their models were tailored for cooler stars (T eff ∼ 15 − 55kK).We emphasize that, to obtain an accurate comparison, it is necessary to also allow other wind parameters than mass-loss rate to vary.If, for example, the winds were faster than the fixed v ∞ = 2500 km s −1 , higher mass-loss rates compared to our estimates would be allowed.
We note that the optical spectral lines that are sensitive to circumstellar gas cannot be used to determine the exact origin of this moving material.While stellar winds are expected for hot and helium-rich stars, these stars are binaries and gas could originate from disks, outflows, or ejecta (e.g., Gies et al. 1998;Smith et al. 2011a;Mauerhan et al. 2015).Such gas could, potentially, have an impact on these optical spectral lines that could be confused with stellar winds.To measure direction, speed, and better constrain the amount of circumstellar material -thus also its origin -UV spectroscopy is needed.This is the focus of an upcoming study in our series (HST/COS cycle 29 PI: Drout, HST/COS cycle 30 PI: Götberg).

EMISSION RATES OF IONIZING PHOTONS
The emission rates of ionizing photons cannot be directly measured.But, they can be inferred from the shapes of the modeled spectral energy distributions.We estimate the emission rates of H, He, and He + ionizing photons, referred to as Q 0 , Q 1 , and Q 2 , by integrating the spectral energy distributions of the best-fit model and the models within 1σ error, following: where we integrate from 50Å, which is the shortest wavelength included in the spectral models, until λ lim , which is the ionization edge for the given atom or ion (912Å, 504Å, and 228Å for H, He, and He + , respectively) and thus sets whether c is the speed of light, λ is the wavelength, and L λ is the wavelength dependent luminosity.We also do not account for the effect of wind mass loss when estimating the ionizing emission rates.However, within the expected regime of weak winds (see Sect. 7), we do not expect large variations in either of the ionizing emission rates (cf.Schmutz et al. 1992).
We present the emission rates of ionizing photons in Table 3 and plot them in Fig. 13.The figure shows hardness diagrams, where we plot Q 1 as function of Q 0 in the left panel, and Q 2 as function of Q 0 in the right panel.The dotted lines show the ratio between the helium to hydrogen ionizing emission rates as labeled.The figures show that, while roughly half of the hydrogen-ionizing photons are also helium-ionizing photons (for all stars but star 16), only a small fraction of them are also He + -ionizing (typically ∼ 0.001 − 0.1%).
We expect that stars 1-8 have Q 0 ∼ 10 47.5 − 10 49 s −1 , Q 1 ∼ 10 47 − 10 49 s −1 , and Q 2 ∼ 10 43 − 10 47 s −1 .We compare these to the expected emission rates of ionizing photons from models of stripped stars with Z = 0.006 (Götberg et al. 2018) and models of OB main-sequence stars and WN-type WR stars from the 0.4 Z models from Smith et al. (2002) in Fig. 13.As the figure shows, the H-ionizing emission rates of stars 1-8 are similar to mid-late O-type main sequence stars, but lower by a factor of a few compared to WN-stars.Compared to OB-stars, stars 1-8 and 26 have harder ionizing emission, with typically more than an order of magnitude higher He 0 -ionizing emission rates compared to OB stars of the same Q 0 .Main-sequence stars with similar Q 0 as stars 2-8 are expected to emit many orders of magnitude lower rates of Q 2 .In fact, WN stars with similar temperatures as stars 2-8 also are expected to emit He + -ionizing photons at substantially lower rates, because of their opaque stellar winds.
Figure 13 demonstrates the important role the effective temperature plays for the emission rate of ionizing photons.Star 1 is the hottest star in the sample, and also the star with the hardest ionizing spectrum, where more than 1% of the hydrogen-ionizing photons also are He + -ionizing.In fact, star 1 is expected to have a similar emission rate of hydrogenionizing photons as an O7V-type star, but a three orders of magnitude higher emission rate of He + -ionizing photons (Smith et al. 2002).Götberg et al. (2018) predicted that stripped stars with masses ∼ 3−4 M should have Q 0 ∼ 10 48 s −1 , Q 1 ∼ 10 47.5 s −1 and Q 2 ∼ 10 44 − 10 45 s −1 .As seen from Table 3 and Fig. 13, stars 2-7 agree well with these predictions.We note that large variations in Q 2 were already predicted by Götberg et al. (2018) (see also Götberg et al. 2017) as a result of both metallicity variations and wind mass-loss rates.While the right panel of Fig. 13 exhibits an apparently smooth trend for Q 2 with Q 0 , we note that further observational explorations are needed to accurately determine the emission rates of ionizing Inferred emission rates of H-, He-, and He + -ionizing photons (Q 0 , Q 1 , and Q 2 , respectively), plotted against each other to explore ionizing hardness for the stars in the spectroscopic sample and using numbered colored symbols.A large fraction (∼ 50%) of the H-ionizing photons are He-ionizing, but only a small fraction (∼ 0.001 − 1%) are He + -ionizing.This shape of the spectral energy distribution is expected for stars with temperatures ∼ 50 − 100kK, but remains to be observationally confirmed.For comparison, we also display models with Z = 0.006 for stripped stars by Götberg et al. (2018) using pale blue and labeled with the stripped star mass, along with models with Z = 0.4Z from Smith et al. (2002) for OB-type main-sequence stars in dark gray, labeled by spectral types, and for WN-type WR stars in light gray, labeled by temperature in kK.
photons from stripped stars.Such observational explorations could include for example nebular ionization studies.

IMPLICATIONS FOR BINARY EVOLUTION
With the parameter determinations described in this paper, there are several topics interesting to discuss in the context of interacting massive binary stars.We choose a subset here.

Resulting surface composition from envelope-stripping
The stripped stars in our sample have a range of surface hydrogen mass fractions, from about 0.4 down to negligible amounts (see Sect. 5 and Table 3;and also Appendix B).This suggests that envelope-stripping results in both hydrogenpoor and hydrogen-free stars.Because leftover hydrogen can affect both the effective temperature, ionizing emission rates, future expansion and thus binary interaction, and supernova type, this result suggests that approximating stripped stars with pure helium stars may lead to a poor representation.
A range of surface hydrogen mass fractions has been predicted from models (e.g., Yoon et al. 2017) and is thought to arise from how deeply the stars are stripped into the chemical gradient that results from the receding main-sequence core.The depth of stripping could depend on how large the Roche lobe was at detachment (for the case of stable mass transfer), the metallicity and thus opacity of the stellar envelope (e.g., Sravan et al. 2019), and perhaps also whether the envelope was stripped via common envelope ejection or stable mass transfer (e.g., Ivanova 2011).Given the weak stellar winds, we consider it unlikely that wind mass loss after envelope-stripping significantly affects the surface hydrogen content of these stars.Because, with a typical wind mass-loss rate of 10 −9 M yr −1 and typical stripped star durations of 1 Myr, only about 0.001 M of material can be removed during the stripped star phase.The total mass of hydrogen expected for stripped stars with surface hydrogen mass fraction of 0.3 and stellar masses 2-7 M is 0.03-0.06M (Götberg et al. 2018).
To establish the relation between the amount of leftover hydrogen and the envelope-stripping mechanism, orbital monitoring is needed.If stripped stars with hydrogendepleted surfaces predominantly have short ( 1 day) orbital periods, this would suggest that common envelope ejection removes more hydrogen.The surface hydrogen content could thus provide an easy way to determine the envelopestripping mechanism and identify different types of binary systems.

Companion types
In this paper, we have chosen to analyze stripped stars whose flux dominates the optical spectrum and for which no evident sign of a bright companion is present (see also Appendix B).Despite this apparent lack of a companion star, the stripped stars exhibit radial velocity variations consistent with orbital motion.This suggests that optically faint companion stars are present.Such companions can only be lower-mass main sequence companions or compact objects.
In Drout & Götberg et al., under review, we found that stripped star + main-sequence star systems will appear as "Helium-star-type" if the main-sequence star is (1) 0.6 times as massive as the stripped star, and (2) early on its main-sequence evolution (which is expected from binary evolution if the companion is that much less massive).Assuming that stripped stars typically are about a third as massive as their progenitors, this critical mass ratio of q crit = 0.6 translates to a critical initial mass ratio of q crit,init = 0.6 × 1/3 = 0.2.If interaction is initiated in a system with q init < 0.4, it is thought that a common envelope should develop (Hurley et al. 2002).We have therefore reason to believe that the stripped stars of "Helium-star-type" are the result of common envelope ejection when orbiting MS stars or stable mass transfer/common envelope ejection when orbiting compact objects.
To better explore what kinds of objects have stripped these stars, orbital monitoring, lightcurve studies, and X-ray observations will be important.The "composite-type" and "Btype stars" with UV excess presented by Drout & Götberg et al., under review provide an opportunity to study companion stars and assess how they were affected by the previous envelope-stripping phase, which could have led to mass gain and spin-up for the accretor stars.To further explore the masses and types of accretor stars, methods such as those of Wang et al. (2018Wang et al. ( , 2021)), who used cross-correlation of spectra in the ultraviolet regime to search for subdwarf companions to rapidly rotating Be stars, could be of interest, since it successfully reaches the part of the population of stripped star systems that do not exhibit UV excess.9.3.Future evolution to supernovae and compact objects According to our evolutionary mass estimates, seven stars are more massive than 2.5 M , meaning that they most likely will reach core collapse (cf.Tauris et al. 2015), and thus explode as stripped-envelope supernovae (e.g., Drout et al. 2011;Lyman et al. 2016;Yoon et al. 2017).With some that have leftover hydrogen and others that are consistent with no leftover hydrogen (Sect.9.1), in conjunction with low wind mass-loss rates (Sect.7), these stars likely will result in both type Ib (hydrogen-free) and type IIb (hydrogen-poor) supernovae.
The structure models of stripped stars with mass > 2.5 M from Götberg et al. (2018) have surface hydrogen mass fractions of X H,surf ∼ 0.25 − 0.30 and corresponding total hydrogen masses of 0.04 − 0.06 M .According to computations from Hachinger et al. (2012), such hydrogen masses should result in type IIb supernovae.If the stellar structure of these models is representative of stripped stars, this should mean that stars 1, 2, 4, and 6 should result in IIb supernovae.Stars 3, 5, and 7 have substantially lower or negligible surface hydrogen mass fractions (see Table 3).The type of their resulting stripped-envelope supernovae is less evident, and they could result in either IIb (Dessart et al. 2011) or Ib (Hachinger et al. 2012).
It is possible (likely for short-period systems) that the stripped star will fill its Roche-lobe anew after central helium depletion, during helium-shell burning (Laplace et al. 2020).This interaction stage should remove some or all leftover hydrogen, depending on when the interaction is initiated and how much hydrogen is left.The helium can only be removed for extremely short period systems (P orb 0.5 days, cf.Tauris et al. 2013Tauris et al. , 2015)), thus limiting the evolutionary pathways leading to type Ic supernovae, unless any leftover helium remains hidden during the explosion (e.g.Piro & Morozova 2014).
Assuming core-collapse will lead to the creation of a 1.4 M neutron star, we expect that the stripped stars in our sample should produce ejecta masses of ∼ 1.5 − 2.7 M for all stars with masses > 2.5 M apart from star 1, which could have as much as ∼ 7 M ejecta.These numbers agree with the obsessionally constrained ejecta masses for most stripped-envelope supernovae (e.g.Drout et al. 2011;Lyman et al. 2016).
Because of its higher mass, it is possible that star 1 will create a black hole.While it is difficult to know what mass such a black hole would have, it could be similar to the mass of the carbon/oxygen core.Laplace et al. (2021) estimate the carbon/oxygen core mass to be 6.2 M for a 8.2 M heliumcore mass, which is similar to the evolutionary mass of star 1.In conjunction with its low metallicity, this could make star 1 a good calibrator for evolutionary pathways leading to merging black hole binaries.
Stars 8, 16, and 26 (assuming it is residing in the foreground) have lower predicted masses compared to the rest of the sample, sand should lead to white dwarf creation.Stars 8 and 16 likely have current masses above the Chandrasekhar limit and therefore should lose some material before white dwarf creation.Assuming the mass lost will be the outermost layers, they should lose all of the remaining hydrogen and could thus result in DB type white dwarfs.Given that stars 8, 16, and 26 most likely are, or will be, helium-burning objects, they should evolve into C/O white dwarfs.
Depending on the magnitudes of potential kicks present at compact object formation, the orbit of these binaries will be affected.Orbital solutions for the current systems will help constrain possible future evolutionary pathways, in some cases potentially leading to double compact object formation.

SUMMARY & CONCLUSIONS
We present a spectroscopic analysis to obtain the stellar properties for a set of 10 stars first presented in Drout & Götberg et al., under review that we argue are stripped of their hydrogen-rich envelopes via binary interaction.We measure directly from the spectral fitting, for all but one star, effective temperatures confidently above 50kK, surface gravities log g ∼ 5 and surface hydrogen (helium) mass fractions ∼0-0.4 (∼1-0.6).By fitting the spectral energy distribution of the models to UV and optical photometry, we obtain low extinction values (A V ∼ 0.1 − 0.65) and bolometric luminosities of ∼ 3 × 10 3 -10 5 L .Combined with effective temperature and surface gravity, we then estimate stellar radii ∼ 0.6 − 1.5R and spectroscopic masses ∼ 0.8 − 6.9M .Using a mass-luminosity relation from binary evolution models, we estimate the evolutionary masses to ∼ 1.2 − 8.4M .
These properties agree well with the expectations from detailed binary evolution models for helium-core burning stars that have been stripped of their hydrogen-rich envelopes in binaries.This confirms the prediction that the large majority of hydrogen-rich envelopes can be stripped off during binary interaction, leaving the helium core exposed with no or only a thin layer of hydrogen-polluted material left on the surface (Götberg et al. 2017).
Our analysis of the observed properties of stripped stars helps to strengthen several expectations about envelopestripping in binaries that have existed for several years, but which have remained untested: 1. Stars stripped in binaries can be sufficiently massive to reach core-collapse.Thus, they most likely can produce neutron stars and black holes.However, they can also be progenitors for white dwarfs.
2. Stars stripped in binaries can have some or no residual hydrogen left on their surfaces after envelopestripping.This suggests that binary-stripped stars are progenitors of both Ib and IIb supernovae.
3. Stars can be stripped by compact objects or low-mass stars.This must be true because the stripped stars we analyze here dominate the optical spectrum.
4. The stellar properties expected from binary evolution models where stars are stripped via stable mass transfer reflect the observed stellar properties reasonably well.
5. While detailed analysis of ultraviolet spectra is needed, the optical spectra indicate that the wind mass-loss rates from stripped stars are likely lower ( Ṁwind 10 −9 M yr −1 ) than expected from extrapolations of Wolf-Rayet wind mass-loss schemes, and possibly also single-temperature helium star schemes.These low mass-loss rates suggest that winds are unimportant in the removal of residual hydrogen or stripping of the helium layer, suggesting such removal only can happen through future binary interaction.
The derived stellar masses and general stellar properties of the stripped stars indicate that we have filled the gap in the helium-star mass range, creating a bridge between subdwarfs and Wolf-Rayet stars.This observed stellar sample offers opportunities to constrain uncertain physics, such as understanding wind mass loss from hot and helium-rich stars and the period evolution of interacting binaries.
To explore the full parameter space of stripped star binaries, studies reaching systems with massive and exotic companions, along with a Galactic sample, will be needed.A more complete coverage over the binary parameter space will provide better constraints for binary evolution and population synthesis models.Larger samples will also provide the opportunity to study the effect of metallicity on massive binary interaction, which could lead to a better understanding of the distant, young Universe when metallicity was low.The research field of massive stars, and especially stripped helium stars, is and will be even more dependent on incoming ultraviolet data from the Hubble Space Telescope.These data are crucial for studying stellar winds, but also likely the vast majority of stripped stars, which are thought to orbit brighter and more massive main-sequence stars (Wang et al. 2021).Conversely, identifying and studying the effects on the companion stars, affected by significant mass accretion and spinup due to binary interaction, will require UV spectroscopy.In this appendix, we show, for each star, the detailed fits that give rise to the properties that we present in this paper (see Table 3).A description for how these fits are performed, see Sect.4.2.
We use a set of the strongest and most robustly modeled spectral lines of hydrogen and helium for the spectral fitting.These usually include Hδ/He ii λ4100, He ii λ4200, Hγ/He ii λ4339, He ii λ4542, and Hβ/He ii λ4859, and when present we also include He i λ5876.We avoid to use Hα and He ii λ4686 for the fits since they are α lines and are therefore very sensitive to stellar wind and surrounding ionized gas, which can impact the determination of the stellar properties we focus on here (see Sect. 7).We also avoid using He ii λ5412 when possible because this spectral line sometimes has contributions from the outer parts of the stellar atmosphere, which is affected by the density and thus also the stellar wind.
Which exact spectral lines that we use for the different stars are presented in Table 4.In the case of stars 2 and 3, He ii λ4200 is affected by noise and we therefore chose to include also He ii λ5412 in the fits.When observing star 7, we needed to rotate the telescope out of the parallactic angle to avoid including nearby stars in the slit, which led to poor signal-to-noise ratio in the blue part of the spectrum and we therefore chose to exclude Hδ/He ii λ4100 and He ii λ4200.
In Figs. 6,14,15,16,17,18,19,20,21,22,and 23,we show the detailed fits to the spectroscopy and photometry of the stars.Each set of panels display the same things for each star and we describe them below.
The top left panels show zoom-in panels for the wavelength range of each spectral line that is used for the spectral fit.The black line with errorbars show the observed spectrum, the colored thick line shows the best-fit model, and the colored thin lines show the models allowed within 1σ errors.
The 1σ errors are determined using χ 2 (see Sect. 4.2), and we therefore display the χ 2 for each included model as function of the three parameters the model grid spans (effective temperature, surface gravity and surface hydrogen mass fraction) in the top right panels.The best-fit model, which has the minimum χ 2 , is marked with a large colored circle and the models allowed within 1σ are shown with colored circles located below the black line marked 1σ.Models that are not allowed within 1σ are shown as gray circles.The properties resulting directly from the spectral fit are written at the very top right.
To demonstrate that the best-fit model also matches other spectral features, we show a larger wavelength range together with the best-fit model in the two middle panels.For convenience, we mark the lines used for the spectral fit with colored background and we also give a rough estimate for the signal-to-noise ratio (SNR) of the observed spectrum.
We show the fit to the photometry in the bottom left panel.The panel shows the observations with associated errors from Swift (the three bluest datapoints) and Swope (the four reddest datapoints) in black and located at the mid-wavelength of the filter function (Rodrigo et al. 2012;Rodrigo & Solano 2020).All models allowed within 1σ from the spectroscopic fit are shifted to their respective best-fit magnitude and extinction and shown in color.The best-fit model from the spectroscopic fit is shown with large colored circles and a thick line.The resulting bolometric luminosity and extinction are written in the middle at the bottom together with the estimates for stellar radius and spectroscopic mass that follows (see Sect. 4.2).The evolutionary mass is estimated from the mass-luminosity relation described in Sect.4.2.5.
In addition, we also show the models allowed within 1σ in the Hertzsprung-Russell diagram and marked with black dots.The best-fit model is showed as a large colored circle and the errorbars indicate the extent of the models allowed within 1σ.For reference, we display detailed evolutionary models for donor stars in binary systems from Götberg et al. (2018) and for initial masses of 5.5, 6.7, 8.2, 10, 12.2, 14.9, and 18.2 M , which correspond to stripped star masses of 1.5, 1.9, 2.5, 3.4, 4.5, 5.9, and 7.3 M .The evolutionary models are monotonically brighter with mass.For stars in the LMC and SMC, we show models from the Z = 0.006 and Z = 0.002, respectively.

B. IMPACT OF THE COMPANION STAR ON FIT
In this paper, we chose to fit the spectra of stars with "Helium-star-type" spectral morphology, approximating their spectra as single, although these stars exhibit binary motion.While these stars at maximum have a very minor contribution from a main-sequence companion, because their spectral morphologies do not show typical signs of main-sequence stars, it is valid to investigate whether a minor contribution can affect the derived stellar properties.
Here, we test the performance of the spectral fitting routine when (1) removing the contribution from a main-sequence companion from the spectrum of star 6, and (2) adding the contribution from a main-sequence companion to the spectrum of star 5.Because we expect that a main-sequence companion should contribute with hydrogen lines, we choose star 6 for the first experiment, since it has measured surface hydrogen content.This experiment is meant to explore whether we could have mistaken the contribution from a main-sequence companion for surface hydrogen content of the stripped star.If true, fitting the spectrum after subtracting a companion star should result in a good fit as well.Similarly, for the second experiment, we choose star 5, because it does not show any signs of surface hydrogen content.If a main-sequence companion could be mistaken by surface hydrogen content, fitting the composite spectrum should result in good fits, but higher derived surface hydrogen content for star 5.
For both tests, we use a spectral model of a late B-type star created using the modeled stellar properties from a 2.2 M evolutionary model, 20% through the main-sequence evolution (see supplementary material of Drout & Götberg et al., under review).We scale the contribution of the B-star such that it contributes both 10% and 20% of the total optical flux in the binary composite.The B-type model does not show any He i lines and its spectrum is dominated by Balmer lines in the optical.We do not simulate smearing of its spectral features that should occur by stacking after correcting for radial velocity shifts of the stripped star in stars 5 and 6.However, we expect that the effect from such smearing on the spectral features is small.We also do not adapt the B-type model for stellar rotation, since it is likely such systems are created through common envelope ejection.
We then fit the test spectra with the models as described in Sect.4.2.When removing the contribution from the B-type star from star 6's spectrum, we find poor spectral fits both when assuming 10% and 20% contribution, as visualized in Fig. 24.This illustrates that the Balmer lines from the B-type companions are so prominent that subtracting their contribution results in spectral features (in particular hydrogen lines) that are poorly fit by single stripped star models.
When instead adding the B-type contribution to the spectrum of star 5, we find a poor fit when assuming 20% contribution, but a realistic fit when assuming 10% contribution with only slightly deep Balmer lines, as evidenced in Fig. 25.This suggests that the presence of a B-type companion that contributes 20% of the flux should be detectable from the spectral morphology.It results in poor fits to the single stripped star models, requiring a fit to two components simultaneously.However, a 10% flux contribution could potentially be missed.The derived stellar properties for the fit with 10% contribution are very similar to those derived for star 5, but with a slightly higher hydrogen mass fraction (X H,surf = 0.05).
Deeper investigation of the binary companions is needed, but requires several additional analyses and will be addressed in a future study.However, from the analysis presented in this appendix, we conclude that the optical contribution from a companion star must be small for the spectral model fits to be good.Therefore, if any, we expect small influence from the companion star on the derived stellar properties.

C. KINEMATIC ASSESSMENT OF STAR 26
Here we carry out a detailed kinematic assessment of star 26 compared to the bulk of objects in the LMC, following the same methodology outlined in Drout & Götberg et al., under review. In  From this, we see that the mean radial velocity of 162 km s −1 is slightly low for the LMC.It overlaps with only the extreme tail of the full sample of OB stars listed on Simbad, and falls below the common threshold of 200 km s −1 often adopted for membership (see e.g.González-Fernández et al. 2015, Davies et al. 2018).In addition, the proper motion values of (µ α ,µ δ ) = (2.86,−4.71) mas yr −1 are significantly offset from the bulk of LMC stars, which have median values of (µ α ,µ δ ) = (1.83,0.30) mas yr −1 .Comparing these proper motion values with the distribution of likely LMC members, we find a χ 2 value of ∼165.This indicates that star 26 is located significantly outside the region that contains 99.7% of likely LMC members (designated by χ 2 < 11.6).In addition, Gaia DR3 lists zero excess noise and an astrometric goodness-of-fit close to zero (astrometric gof al = −0.28)for this object, indicating the the astrometric fit was high quality.
While it is possible for stripped helium stars can receive a kick upon the death of their companion stars, the proper motions observed for star 26, would imply a systematic velocity of ∼1200 km s −1 relative to the mean values for the LMC (assuming a distance of 50 kpc).These values are significantly larger the those predicted for runaway stripped stars of ∼100 km s −1 by Renzo et al. (2019).Thus, we consider it more likely that star 26 is a foreground halo object.This is supported by the fits presented above, which exhibit both a cooler temperature and higher surface gravity than other objects modeled here, consistent with a subdwarf interpretation.In Table 5 we provide the same kinematic information presented for all objects in the sample of Drout & Götberg et al., under

Figure 1 .
Figure 1.The normalized spectra of the stars in the observed spectroscopic sample, described in Sects. 2 and 3.These stars are selected from the sample of Drout & Götberg et al., under review and thought to be stars stripped of their hydrogen-rich envelopes via binary interaction.

Figure 2 .
Figure2.Coverage of the spectral model grid used as base for spectral fitting to obtain stellar properties of stripped stars.This visualization shows which models have reached convergence using a colored circle, where the small and large circles correspond to the temperature and surface gravity at optical depth τ = 20 and 2/3, respectively.Thin, dotted lines indicate where the Eddington factor is 0.1 (leftmost), 0.5 (middle), and 1 (rightmost).

Figure 3 .
Figure 3.Effect of varying surface hydrogen/helium mass fraction (left), temperature (center), and surface gravity (right) on the spectral lines Hδ/He ii λ4100 (first row), He ii λ4542 (second row), He i λ5876 (third row), and N v λ4604 (fourth row).We use the model with T eff = 70 kK, log 10 g = 5.0 and X H,surf = 0.3 as a base (black solid line) and vary each parameter according to the legends.
Fig. 4, N iv λ 4057 is present in emission roughly at T ∼ 60 − 80 kK (dark blue triangles), N v λ 4604 appears in absorption for T ∼ 60 − 90 kK (downward cyan triangles), while it flips into emission for T 100 kK (upward cyan triangles), and N v λ 4945 appears in emission for T 70 kK (teal triangles).On the low temperature end, He i can provide a similar discriminant.As the bottom panel of Fig. 4 shows, He i λ5876 is present for T 70 kK (purple triangles) and He i λ4471 for T

Figure 4 .
Figure4.We use the presence of a set of nitrogen lines (top) and helium lines (bottom) to constrain the model grid used when fitting the observed spectra.As an example, we show here the models with surface hydrogen mass fraction X H,surf = 0.3, spread out in the temperature -surface gravity plane.Triangular markers show the presence of nitrogen or helium lines that are used to constrain the grid (see also Table4).Gray circles indicate models in which none of the lines specified by the legend are present.

Figure 5 .
Figure5.The relation between mass and luminosity for stripped stars half-way through central helium burning (X He,center ≡ 0.5) is used to estimate the evolutionary mass.Here, we show the relation for stars stripped via stable mass transfer using evolutionary models fromGötberg et al. (2018) for Z = 0.006 (pink line) and Z = 0.002 (dark red line).The shaded background color demonstrates the variation throughout helium-core burning when the central helium mass fraction is between 0.8 and 0.1.
photosphere (τ = 2/3) apart from T and log 10 g , which we display for comparison and that correspond to the temperature and surface gravity at τ = 20.

Figure 6 .
Figure 6.Fit for star 1. See Appendix A for the fits of the remaining stars.See Table3for the derived temperature and surface gravity at optical depth τ = 20.
stars (stars 5, 7, and 26) are completely hydrogen free.These values are broadly consistent with the estimates presented in Drout & Götberg et al., under review based on equivalent width diagnostics.

Figure 8 .
Figure 8.The stars in our spectroscopic sample, shown with numbered markers, match well with models of stars stripped in binaries (gray lines) in the Hertzsprung-Russell diagram.The left panel shows the stars from the LMC plotted together with models of Z = 0.006 and the right panel shows stars in the SMC plotted together with models of Z = 0.002.Star 26, which likely is a foreground object, is plotted using an assumed distance of 10 kpc and diamond-shaped marker.We label the zero-age main-sequences and gray-shade the parts of the diagrams with cooler temperatures.Wolf-Rayet stars in each of the clouds are shown using purple circles and a shaded region, while the expected locations of bright subdwarfs are marked with a green-shaded ellipse.The weak-wind WN3/O3 stars in the LMC are indicated using lighter purple color.
matches well with the expected location for heliumcore burning, massive subdwarfs.Compared to the set of Wolf-Rayet stars (dark purple circles, Hainich et al. 2014, 2015; Shenar et al. 2016), WN3/O3 stars (lighter purple circles in the LMC plot, Neugent et al. 2017), and the expected location of subdwarfs in the two clouds (teal shaded regions, cf.Heber 2016), it is clear that the stars in our spectroscopic sample create a connecting bridge between faint subdwarfs and bright Wolf-Rayet stars.

Figure 10 .
Figure 10.Contraction (top) and expansion (bottom) phases for stripped stars demonstrated using the Z = 0.006 stripped star models of Götberg et al. (2018), labeled by stripped star mass.We show the fraction of the stripped star duration using blue and pink shades and the central helium burning phase when 0.9 > X He,c > 0.01 using dark gray background for the evolutionary tracks.The stars in the spectroscopic sample are plotted using their effective radii and bolometric luminosities with numbered markers (see Table3).The top panel shows that contraction lasts ∼ 10% of the stripped star duration, while the bottom panel shows the expansion phase lasts ∼ 1 − 5%.All stars but star 16 agree with the helium-core burning phase and the expansion phase, while star 16 could either be contracting or expanding.

Figure 11 .
Figure11.The shape of the He ii λ4686 spectral line is very sensitive to surrounding gas and we, therefore, use it to estimate wind mass-loss rates.The observed He ii λ4686 lines are shown from top to bottom for each star along with models for a range of wind mass-loss rates ( Ṁwind = 10 −6 -red, 10 −7 -purple, 10 −8 -blue, 10 −9 -green, and 10 −10 M yr −1 -yellow).Right panels show zoom-ins of the observed spectral lines, while left panels show zoom-outs that also include the expectations for wind emission.

Figure
Figure13.Inferred emission rates of H-, He-, and He + -ionizing photons (Q 0 , Q 1 , and Q 2 , respectively), plotted against each other to explore ionizing hardness for the stars in the spectroscopic sample and using numbered colored symbols.A large fraction (∼ 50%) of the H-ionizing photons are He-ionizing, but only a small fraction (∼ 0.001 − 1%) are He + -ionizing.This shape of the spectral energy distribution is expected for stars with temperatures ∼ 50 − 100kK, but remains to be observationally confirmed.For comparison, we also display models with Z = 0.006 for stripped stars byGötberg et al. (2018) using pale blue and labeled with the stripped star mass, along with models with Z = 0.4Z fromSmith et al. (2002) for OB-type main-sequence stars in dark gray, labeled by spectral types, and for WN-type WR stars in light gray, labeled by temperature in kK.

Figure 18 .
Figure 18.Fit for star 6.We have clipped out the line cores of He ii λ4339/Hγ, He ii λ4860/Hβ, He ii λ5412, and He ii λ6563/Hα.

Figure 22 .
Figure 22.Fit for star 26.Assuming a foreground distance of 10 kpc.

Figure 23 .
Figure 23.Fit for star 26.Assuming it is a member of the LMC.

Figure 24 .
Figure 24.Best-fit models for star 6 after having removed the contribution from a 2.2 M late B-type companion star, assuming it contributed 10% (top) and 20% (bottom) of the optical flux.In both examples, the fits are poor.

Figure 25 .
Figure 25.Best-fit models for star 5 after having added the contribution from a 2.2 M late B-type companion star, assuming it contributes 10% (top) and 20% (bottom) of the optical flux.While the 20% contribution results in a poor fit, the 10% contribution is acceptable and almost reproduces the effective temperature, surface gravity and surface hydrogen mass fraction derived for star 5.
Fig. 26 we show both the average radial velocity measured for star 26 (left panel; based on 10 epochs of observations between 2018 and 2022) and the proper motion in RA and DEC from Gaia EDR3 Gaia Collaboration et al. (2020).For comparison, we also show (i) the 16 LMC members presented in Drout & Götberg et al., under review (colored dots; both panels), (ii) a sample of OB stars pulled from Simbad that overlap with the LMC and have radial velocity measurements (grey dots; left panel), and (iii) a sample of bright likely LMC members pulled from Gaia EDR3 (grey dots; right panel; see Drout & Götberg et al., under review for details of sample selection). review.

Figure 26 .
Figure 26.Left: Comparison of the mean radial velocity observed for Star 26 to known OB stars in the LMC (grey histogram).Horizontal "errorbars" designate the range of radial velocities observed at different epochs.Right: Comparison of the Gaia proper motions measured for Star 26 to likely LMC members (grey points).Other LMC stars presented in Drout & Götberg et al., under review are shown as colored/numbered circles in both panels.Figures adapted from Drout & Götberg et al., under review.

Table 1 .
Observations to obtain optical spectra.
(seeDrout & Götberg et al., under review and  Ludwig et al., in preparation).These apparent magnitudes are presented in the AB system, where we have converted the optical data from Vega magnitudes following the description in Sect.3.2.

Table 3 .
Stellar properties of the stripped stars in our spectroscopic sample.

Table 5 .
Kinematic Information for Star 26.