An early dark matter-dominated phase in the assembly history of Milky Way-mass galaxies suggested by the TNG50 simulation and JWST observations

,


INTRODUCTION
In the ΛCDM framework, galaxy formation proceeds through the collapse of dark matter haloes, after which baryons can cool and form stars (e.g.White & Rees 1978).Galaxies grow through the accretion of (cool) gas and subsequent star formation, as well as by merging with nearby systems (e.g., Dekel et al. 2009;Oser et al. 2010).However, the details of these processes, their relative importance as a function of cosmic time, and the interplay between the dark and baryonic matter are still unclear.
Constraints on the mass budget, i.e. the relative contribution of dark matter, gas and stellar mass, can offer insight into the mass assembly histories of galaxies.In the local Universe, moderately massive galaxies degraaff@mpia.de(M * ∼ 10 10−11 M ⊙ ) are baryon-dominated near their centers (∼ 1 kpc), but contain substantial dark matter mass within the effective radius (f DM (< r e ) ∼ 40 − 70% for r e ∼ 5 − 10 kpc; Martinsson et al. 2013).Studies of ionized gas kinematics at z ∼ 1 − 3 paint a very different picture, as these find an increase in the baryonic mass fraction within r e toward higher redshift for galaxies of M * ≳ 10 9.5 M ⊙ , with galaxies appearing almost entirely baryon-dominated at z ∼ 2 (Wuyts et al. 2016;Price et al. 2016Price et al. , 2020;;Genzel et al. 2017Genzel et al. , 2020;;Nestor Shachar et al. 2023), although recent work suggests this may not be the case for all galaxies at this epoch (Sharma et al. 2023).For very massive galaxies at z ∼ 2 Genzel et al. (2017Genzel et al. ( , 2020) ) propose that such high baryon fractions may imply that the dark matter distributions are cored instead of cuspy, which may be caused by a very rapid formation process coupled with stellar and black hole feedback.
The James Webb Space Telescope (JWST) enables the measurement of galaxy kinematics at even higher redshifts and also lower masses, thereby probing the possible progenitors of more massive galaxies at cosmic noon.In de Graaff et al. (2024) we used JWST/NIRSpec observations to measure the ionized gas kinematics for six small (r e ∼ 1 kpc) galaxies at z ∼ 6 − 7 of low stellar mass, M * ∼ 10 7−9 M ⊙ .Interestingly, the dynamical masses inferred from the ionized gas kinematics exceed the stellar masses by a factor ∼ 30, and the baryonic masses by a factor ∼ 3, which may imply surprisingly high dark matter fractions of f DM (< r e ) ≈ 0.7.However, these JWST results are still subject to large observational and systematic uncertainties, raising the question whether these implied dark matter fractions are realistic, and whether these are compatible with observations at intermediate redshifts.
In this Letter, we use the cosmological magnetohydrodynamical galaxy simulation known as TNG50 (Nelson et al. 2019a;Pillepich et al. 2019), and part of the IllustrisTNG project 1 , to extract theoretical expectations for the baryonic mass fraction of galaxies on small spatial scales as a function of redshift and stellar mass.For relatively massive galaxies in the lowerredshift Universe (z ≲ 2 − 3), previous work within Il-lustrisTNG has provided an overall quantification of the dark matter fraction (Lovell et al. 2018) and a direct comparison to observations ( Übler et al. 2021).Here, we use TNG50 to focus on low stellar mass galaxies at z = 6, in order to compare the simulation expectations to the observational results from JWST at z ∼ 6 − 7. By tracing the descendant population in the simulation, we check whether the JWST results can be linked to, and are consistent with, ground-based measurements at cosmic noon.We use physical lengths when specifying distances and sizes throughout, unless otherwise stated.
1 www.tng-project.orgWith a baryonic mass resolution of m baryon = 8.5 × 10 4 M ⊙ , a dark matter particle resolution of m DM = 4.5 × 10 5 M ⊙ and a Plummer-equivalent gravitational softening length of ϵ z=0 DM, * = 0.29 kpc at z = 0, TNG50 has both the large volume and the high resolution needed to robustly evaluate the average mass profiles of galaxies on scales of ∼ 1 kpc across a wide range in redshift and stellar mass.To test the convergence of our results, we also use the TNG100 simulation (Nelson et al. 2019b), which employs the same galaxy formation model, but for a larger volume of ∼ 100 3 comoving Mpc 3 and lower resolution (m baryon = 1.4 × 10 6 M ⊙ ; m DM = 7.5 × 10 6 M ⊙ ; ϵ z=0 DM, * = 0.74 kpc).Haloes and subhaloes are identified using the Friendsof-Friends (FoF) (Davis et al. 1985) and subfind (Springel et al. 2001) algorithms, respectively.We define the effective radius as the 3D stellar half-mass radius of the subhalo: r e ≡ r 1/2, * .Because we trace a wide range in mass and redshift (and hence galaxy size), we define the total stellar mass as the stellar mass enclosed within a spherical aperture of twice the half-mass radius, M * ≡ M * (< 2r e ), instead of using a fixed aperture to compute stellar masses.
We define our primary sample by initially selecting central galaxies identified by subfind at z = 6 (snapshot 13) with a stellar mass M * > 10 7 M ⊙ , which are resolved by > 100 stellar particles.From this sample of 3926 systems we further select the subset of 381 galaxies with 10 8 M ⊙ < M * < 10 9 M ⊙ (i.e.resolved by > 1000 stellar particles), which we use to trace the descendant population and to compare with observations across different redshifts.

MASS BUDGET AT Z ∼ 6
We start assessing the central mass budget of highredshift galaxies by looking at our sample of TNG50 central galaxies at z = 6, which is the median redshift of the JWST measurements by de Graaff et al. 2024, and measuring the baryonic mass fraction in a (3D) spherical aperture of the stellar half-mass radius: , (1) where M gas (M DM ) is the mass of all gas cells (dark matter particles) within the aperture.
Figure 1 shows f baryon (< r e ) versus the total stellar mass in TNG50, color-coded by r e : there is a strong dependence of f baryon on stellar mass.Galaxies of low stellar mass (M * ≲ 10 8.5 M ⊙ ) have baryon fractions in the range f baryon (< r e ) ∼ 0.1−0.3 and thus are dark matterdominated.Above M * ∼ 10 8.5 M ⊙ , f baryon (< r e ) increases rapidly in TNG50, and the most massive sys-tems (M * ∼ 10 10 M ⊙ ) are baryon-dominated.At fixed stellar mass, the scatter in f baryon (< r e ) among individual galaxies is modest (σ(f baryon |M * = 10 8 M ⊙ ) ≈ 0.05), and covariant with r e , so that larger effective radii imply lower values of f baryon .
These simulation results are in remarkably good agreement with the measurements from de Graaff et al. (2024), inferred from JWST observations of low-mass galaxies (M * ∼ 10 7−9 M ⊙ ) at z ∼ 6 − 7. The observational estimates of the baryonic mass were obtained from spectral energy distribution (SED) modeling to low-resolution spectroscopy, while the total mass estimates were inferred dynamically from measurements of ionized gas kinematics.We emphasize that these observed estimates differ conceptually in important ways from those for the simulated baryon fractions.For example, the observed effective radius represents the projected half-light radius of the emission line instead of a halfmass radius, although we find that these are ∼ 1 kpc in both simulated and observed galaxies.The stellar masses were inferred from SED modeling and the coldgas mass is empirically inferred from the star formation rate (SFR); both observational estimates carry significant systematic uncertainties that are not reflected in the error bars shown.
For a complete discussion of the systematic uncertainties in the different mass estimates we refer the reader to de Graaff et al. 2024, but we provide a brief summary here.The stellar mass in particular can be underestimated by as much as a factor 10 in extreme scenarios where only the integrated light of the galaxy is considered in the SED fitting and a luminous young stellar population outshines the underlying population of older stars in the galaxy (e.g.Giménez-Arteaga et al. 2023).However, the spatially-resolved analysis of one of the galaxies shown in Figure 1, with an SED that is representative of the full z ∼ 6 sample, yields an integrated stellar mass that is consistent with the mass inferred from SED fitting to the integrated light (Baker et al. 2023).The inferred gas masses are also uncertain, as these were obtained from the conversion of a star formation surface density under the assumption of a constant star formation efficiency (Kennicutt 1998; for further discussion see also Section 5), whereas we may expect this efficiency to increase toward higher redshift (Tacconi et al. 2020).The star formation efficiency in the low-mass, high-redshift regime is still poorly constrained observationally, but we note that Price et al. (2020) found that for galaxies at cosmic noon the gas masses estimated using the empirical relations by Kennicutt (1998)  ) of galaxies vs. total stellar mass.Throughout, by baryonic fraction we mean the ratio between the mass in gas and stars and the total mass of the galaxy (see Eq. 1).The TNG50 simulation predicts low baryon fractions for M * < 10 9 M⊙ and a sharp increase in f baryon toward higher masses.The half-mass radii are small (∼ 1 kpc in physical units) and, at fixed stellar mass, are anti-correlated with the baryon fraction.With the exception of one source, which may have high f baryon due to an underestimated dynamical mass, the baryon fractions measured at z ∼ 6 with JWST are in good agreement with the simulated galaxies.
in the inferred f baryon ) is still consistent within the scatter seen among the individual galaxies in TNG50.Finally, we note that observations are subject to projection effects that can lead to an underestimated circular velocity.This may for example explain the data point for which the baryonic mass exceeds the total dynamical mass (M baryon > M dyn ).In light of these different issues, the broad agreement found between the observations and the TNG50 data shown in Figure 1 is striking.

REDSHIFT EVOLUTION OF THE BARYON AND DARK MATTER FRACTION IN TNG50
Next, we use these same simulations to trace the mass assembly and mass budgets of z = 6 galaxies to lower redshifts.To do so, we select 381 galaxies within the stellar mass range 10 8 M ⊙ < M * < 10 9 M ⊙ (at z = 6) Figure 2. Mass assembly histories of low-mass TNG50 galaxies, selected to have 10 8 M⊙ < M * < 10 9 M⊙ at z = 6.Although all z = 6 galaxies are selected to be centrals in their dark matter haloes, many systems merge to form massive centrals or massive satellites in groups and clusters at lower redshifts.To construct the median mass history (solid black lines), we consider only unique subhaloes in each snapshot.The median stellar mass of the z = 6 descendant population increases 100-fold by z = 2 ; at z = 0 the median stellar mass of the population equals that of the Milky Way, but with large scatter (0.6 dex).In orange we show the mass assembly histories (median and 5-95 percentiles) for the population of TNG50 MW and M31 analogs from Pillepich et al. ( 2023), which were selected as analogs not only by stellar mass, but also by environment and stellar morphology.We find that 50 out of 213 galaxies of the z = 6 descendant population are indeed MW or M31 analogs at z = 0.
and use the merger trees identified with sublink (i.e., the DM-only merger trees from Rodriguez-Gomez et al. 2015) to link progenitors and descendants.

Mass assembly histories
Although all galaxies were selected to be the central galaxy in their FoF halo at z = 6 , at lower redshift these systems may fall into larger haloes and hence become satellite galaxies.Moreover, the z = 6 galaxies are not necessarily the main progenitor of their z = 0 descendants.Therefore, we forward trace the merger history for a subhalo selected at z = 6 by iteratively searching for its next descendant; this is straightforward if the subhalo at a given snapshot is the main progenitor of the descendant in the following snapshot.If this is not the case (usually due to a merger with a more massive subhalo, either from the initial sample or a massive system not in the initial selection), we trace the subsequent evolution of the product of the merger.This implies that from the 381 galaxies at z = 6, only 213 descendants remain by z = 0.When calculating population statistics, we only consider unique descendant galaxies, i.e. we do not double count objects that have merged to form a single system at z < 6.
In Figure 2 we illustrate the mass assembly histories of the z = 6 descendant population through the redshift evolution of the stellar and the host halo mass (M 200 ).Thin gray lines show individual mass histo-ries, and the black line shows the evolution of the ensemble median.In TNG50, the low-mass galaxies at z = 6 typically grow to 10 10 M ⊙ by z = 2, and are therefore comparable in stellar mass to the objects that are typically studied kinematically at cosmic noon (e.g.Wuyts et al. 2016;Price et al. 2016Price et al. , 2020;;Wisnioski et al. 2015Wisnioski et al. , 2019)).By z = 0 the median stellar mass is 10 10.5 M ⊙ -comparable to the stellar mass of the Milky Way (e.g.see Licquia & Newman 2015; Bland-Hawthorn & Gerhard 2016) -although with an ensemble scatter of 0.6 dex.In contrast, the typical descendant halo mass (∼ 10 13 M ⊙ ) far exceeds the virial halo mass of the Milky Way (∼ 10 12 M ⊙ ).
For comparison, we show in Figure 2 the sample of Milky Way and Andromeda (hereafter MW and M31) analogs from Pillepich et al. (2023), who selected 198 TNG50 galaxies at z = 0 with stellar masses, morphologies and 1 Mpc scale environment similar to the MW and M31.At fixed redshift, these analogs typically have slightly lower stellar masses and significantly lower host halo masses than the sample of z = 6 descendants selected in this work.However, we note that there is overlap between the two populations, as 50 galaxies of the 213 descendants in our selection are also in the MW/M31 analog sample of Pillepich et al. (2023).As mentioned in Section 3, the observations and simulations differ significantly in their measurement methodology for all components of the overall mass budget: stars, gas and dark matter.Notably, the effective radii represent different physical quantities, one being based on stellar mass and the other on emission line extent.Moreover, although the mass-size relation in TNG is broadly in agreement with observations at z ≤ 4 and M * > 10 9 M ⊙ (Pillepich et al. 2019), it is still unclear whether this holds at higher redshifts and lower stellar mass.A detailed forward modeling of the simulations to create and analyze realistic mock images and mock spectra would be needed to perform an apples-to-apples comparison with the observations at z ∼ 6.However, this is beyond the scope of the current work.
Instead, to mitigate uncertainties in the sizes of simulated galaxies and their comparison to observed galaxies, in what follows we measure the mass budget within a fixed physical aperture.Note that for simulations this is still a spherical radius, and for the observations a projected one.We set the aperture radius equal to the average half-light radius found in the z ∼ 6 JWST observations: r aper = 1 kpc.We hence measure the stellar, gas and dark matter aperture masses for the z = 6 descendant sample across 0 < z < 6.Separately, for each snapshot we also measure the aperture masses for all central galaxies in the snapshot, in order to map the relative mass components as a function of redshift in narrow bins of stellar mass (analogous to the work by Lovell et al. 2018).
We show the redshift evolution in f baryon (< 1 kpc) in the left panel of Figure 3 for the z = 6 descendant population (black line, shaded regions).Dashed lines show the redshift evolution at fixed stellar mass, for a (colorcoded) range of stellar masses.Unsurprisingly, TNG50 predicts that the baryon fraction depends strongly on stellar mass.But remarkably, Figure 3 shows that f baryon (< 1 kpc | M * ) is nearly constant from z ∼ 6 to the present across a wide range of (fixed) stellar masses.Dash-dotted lines show the analogous estimates for the lower-resolution TNG100 simulation, implying that resolution has only a weak effect on this result.In fact, we have checked with the three different resolution levels of the same TNG50 volume that the simulation outcomes are converging in most regimes.
For comparison, we also extract the mass components of central galaxies in the EAGLE cosmological hydrodynamical simulation (Schaye et al. 2015;Crain et al. 2015), specifically the reference model run for a volume of 100 3 comoving Mpc 3 (Ref-L100N1504), using halo and subhalo catalogs that were constructed with the exact same methods as for the TNG simulations2 .
The EAGLE simulation has similar baryonic and dark matter particle mass resolution to the TNG100 simulation and assumes a similar cosmology.In the context of this paper, the most important difference between the two simulations is the subgrid galaxy formation model, as the two use different prescriptions for, e.g., star formation, stellar feedback and black hole feedback.It is therefore interesting that -just like the TNG simulations -EAGLE predicts f baryon (< 1 kpc) depends strongly on stellar mass at fixed redshift, while f baryon (< 1 kpc | M * ) is only weakly dependent on redshift, noting that the apparent decrease from z = 1 to z = 0 likely reflects the fact that the mass profiles of z ∼ 0 EAGLE galaxies are affected by the resolution of the simulation for r ≲ 2 kpc (see Schaller et al. 2015;de Graaff et al. 2022).
Turning to the descendants of the selected z = 6 galaxies in TNG50, we find that for any individual descendant f baryon (< 1 kpc) increases rapidly toward lower redshift: while the central region is still dark matter-dominated at z = 6 with f baryon (< 1 kpc) ∼ 0.25, it is highly baryon-dominated at z ∼ 2 (f baryon (< 1 kpc) ∼ 0.7).Considering the strong dependence of f baryon on M * , this evolution can be understood by the rapid growth of total stellar mass over time found in Figure 2.This M *dependence also explains the redshift-dependent ensemble scatter (i.e. the width of the contours): the sample at z = 6 was selected in a narrow range in stellar mass, where the scatter in f baryon is small (see also Figure 1).As the stellar mass distribution broadens toward lower redshift (Figure 2), the scatter in f baryon (< 1 kpc) at fixed redshift also increases.
We also want to put the stellar-to-gas mass ratio of de Graaff et al. (2024) in the context of the TNG50 simulation.We therefore examine the stellar-to-gas mass ratio within the same 1 kpc aperture.To plot the mass ratio M * /M gas (< 1 kpc) in the case where M gas (< 1 kpc) = 0, we set a floor of M gas = 0.01 × M * .The right panel of Figure 3 shows that also this mass ratio depends strongly on stellar mass at fixed redshift.In contrast with f baryon (< 1 kpc), the TNG50 simulation shows a strong redshift dependence of M * /M gas at fixed stellar mass: the ratio is an order of magnitude higher at z = 0 than it was at z = 6.This shows that also M * /M gas (< 1 kpc) ∼ 0.25 for low-stellar-mass galaxies in TNG50 at z ∼ 6 is broadly consistent with the JWST results of de Graaff et al. (2024), who found M * /M gas (≲ 1 kpc) ∼ 0.1.We note that, although the TNG100 simulation shows qualitatively similar trends to the high-resolution TNG50 simulation, at lower stellar masses and higher redshifts the values of M * /M gas (< 1 kpc) can differ by up to a fac-tor 3 between the two TNG simulations.This indicates that the precise value of M * /M gas (< 1 kpc) depends on resolution.Therefore, it seems difficult to make a robust quantitative comparison between the observed and simulated values.
Toward lower redshift the stellar-to-gas mass ratio increases sharply for the z = 6 descendants, and this central build-up of stellar mass coincides with the rapid growth of the total stellar mass (Figure 2), reminiscent of the compaction phase found in zoom-in simulations of massive galaxies at z ∼ 2 − 4 (e.g.Zolotov et al. 2015;Lapiner et al. 2023).At z ≲ 2 we find that M gas (< 1 kpc) ∼ 0, and the mass ratio therefore plateaus at M * /M gas (≲ 1 kpc) ∼ 100 due to our imposed gas mass floor.This likely reflects the kinetic feedback mode from active galactic nuclei (AGN) in the Illus-trisTNG model (Pillepich et al. 2018b), in which black hole-driven winds remove gas from their surroundings.
Comparing with the EAGLE simulation, we again find qualitative agreement between the EAGLE and TNG simulations.At z = 6 and low stellar mass the value of M * /M gas (< 1 kpc) agrees well with that of TNG50 and therefore also agrees with the JWST observations within a factor ∼ 2. At higher masses and lower redshift the stellar-to-gas mass ratio differs increasingly strongly from the TNG simulations, reflecting the differences in the EAGLE and IllustrisTNG galaxy formation models.
Putting the results of Figure 3 together, we can gain insight into the physical processes behind the redshiftindependence of f baryon at fixed stellar mass.The strong dependence of f baryon (< 1 kpc) on stellar mass implies there is a M * -dependent balance between the cooling of baryons on the one hand, and feedback from stellar evolution and/or AGN on the other hand.Indeed, Nelson et al. (2019a) showed that this feedback in TNG50 depends strongly on stellar mass: at z > 4, where feedback is driven primarily by stellar winds and supernovae, the mass loading factor is a factor ≈ 30 higher for lowmass (∼ 10 8 M ⊙ ) galaxies than for the highest masses (∼ 10 11 M ⊙ ), and explains why at these high redshifts the most massive galaxies are able to rapidly form a stellar mass-dominated center.
Toward lower redshifts the mass loading factor increases, especially for more massive galaxies, and the M * -dependence for the feedback is thus reduced.Yet, despite this redshift evolution in the feedback efficiency, which implies that gas can be expelled more easily from the centers of galaxies at all masses at lower redshift, f baryon (< 1 kpc | M * ) remains constant with redshift.This can be partially explained by the fact that the cooling of baryons is also redshift-dependent.However, perhaps more importantly, the ratio of M gas /M * changes with redshift as well: at z < 3 a substantial fraction of the baryonic mass is in the form of stars (even at low stellar mass), which are not affected explicitly by feedback.We note that AGN feedback, which becomes prominent at z ≲ 2 and high stellar masses in TNG50, therefore also plays little role in the evolution of f baryon , as by this time the stellar mass-dominated center has already been formed.
It is thus likely the interplay between cooling, feedback, and central stellar mass growth that conspires to the redshift-independence of f baryon (< 1 kpc | M * ) in the TNG simulations.The fact that the EAGLE simulation with a very different subgrid model shows a qualitatively similar trend is intriguing.At z > 5 there appears to be a slight upturn in f baryon (< 1 kpc | M * ) in the EAGLE simulation, which may reflect physical differences in the stellar feedback or cooling prescriptions with respect to the TNG model at high redshifts.

Comparison with observations at cosmic noon
Despite the caveat that there are substantial differences in the observed and simulated baryonic mass fractions, it is instructive to draw a qualitative picture of the redshift evolution of the mass budget within the effective radius across cosmic time.
In Figure 4 we show the dark matter fraction, f DM (< r e ) ≡ 1−f baryon (< r e ), as a function of redshift.For the TNG50 z = 6 descendant sample the median evolution in f DM (solid line) is color-coded by the median stellar mass at the given redshift, and contours indicate the variation within the population.Similar to Figure 3, the ensemble median and scatter in f DM (r e ) at fixed redshift depend strongly on the stellar mass distribution, which broadens in width and shifts to higher stellar masses toward lower redshift.The evolution in r e also plays an important role, leading to a slight upturn in f DM (< r e ) at z < 1, and further increasing the ensemble scatter (see also Figure 1).We also show results at z ∼ 1 − 3 from various ground-based studies that obtained dynamical masses from ionized gas kinematics (Wuyts et al. 2016;Genzel et al. 2020;Price et al. 2020;Puglisi et al. 2023).All observed data points are color-coded by the total stellar mass, typically inferred from SED modeling to photometric data.
The effective radius for the simulated galaxies is the 3D stellar half-mass radius, which differs from the observed half-light radii, and therefore affects the comparison in f DM (< r e ).Specifically, based on the analysis of mock JWST images of TNG50 galaxies by Costantin et al. (2023), we estimate that the half-light radii of low-mass TNG50 galaxies at z ∼ 6 are ≈ 0.3 dex smaller than the half-mass radii.On the other hand, at z = 3 and M * ∼ 10 10 M ⊙ the half-light radii are ≈ 0.2 dex larger than the half mass radii.Pillepich et al. (2019) and Rodriguez-Gomez et al. (2019) presented half-light radii at lower redshifts, from which we also estimate a difference of ≈ 0.2 dex between the halflight and half-mass radii of galaxies with stellar masses ∼ 10 10−11 M ⊙ .Therefore, we also show the median evolution for f DM (< 0.5r e ) and f DM (< 1.5r e ) in Figure 4: the effect of the different apertures is to shift the median evolution in f DM down (up) by a factor ∼ 1.4 (∼ 1.2).
Taking into account the large scatter in both the observations and simulations, the two appear to agree well both at intermediate and high redshifts.In detail, however, there may be systematic discrepancies, due to both the precise galaxy sample selection and the method used to infer the dark matter fraction observationally.For example, Übler et al. (2021) constructed and analyzed mock integral field spectroscopic observations for a sample of 7 massive (M * > 10 10.6 ) and highly star-forming galaxies at z = 2 from the TNG50 simulations.By exploring mock observations at multiple viewing angles per object, they showed that there is a slight systematic bias in the inferred dark matter fraction (f DM (< r e ) ≈ 0.34 in TNG50 vs. f DM (< r e ) ≈ 0.20 for real galaxies) and large scatter per object depending on the viewing angle (f DM,max (< r e ) − f DM,min (< r e ) ∼ 0.05 − 0.5).
Furthermore, the dependence of f DM on stellar mass suggests that the observations at z ∼ 2 and z ∼ 6 may be entirely consistent with each other, if, as suggested by the TNG50 merger histories, the z ∼ 6 galaxies are progenitors of the more massive baryon-dominated systems observed at z ∼ 2. If the observed systems at z ∼ 6 can be confirmed to be dark matter-dominated (discussed in Section 5), these objects would form the first direct detection of the progenitors of z ∼ 2 massive galaxies.

DISCUSSION & CONCLUSIONS
We have used the cosmological simulation TNG50 to extract theoretical expectations for the stellar, gas and dark matter mass components in the central parts of high-redshift galaxies.We have explored these components' fractional contributions as a function of stellar mass and redshift.In particular, we have shown that TNG50 predicts that the central baryon fraction of galaxies is independent of redshift at fixed stellar mass: the centers of low stellar mass galaxies have always been dark matter-dominated.The central baryon fraction, is however, a strong function of stellar mass at all redshifts, and so greatly increases along the evolutionary paths of individual galaxies.The central stellar-to-gas mass ratio, on the other hand, grows strongly towards lower redshifts, even at a fixed stellar mass.(5-95) percentiles of the evolution of TNG50 galaxies selected at z = 6 with 10 8 M⊙ < M * < 10 9 M⊙.The thick solid line shows the median evolution of fDM(< re) of the TNG50 galaxies; to bracket possible systematic differences in the comparison to observations, we also show the median evolution of fDM(< 0.5re) (dotted) and fDM(< 1.5re) (dashed).The median evolutionary tracks and data points from high-redshift observations are colorcoded by the total stellar mass.Qualitatively, this suggests a link between observations at different mass and redshift epochs, and agreement between the TNG50 simulation and observations.
Turning to high redshifts, we compare these results to recent measurements of the mass budget at z ∼ 6 from JWST observations of galaxies of low stellar mass (M * ∼ 10 7−9 M ⊙ ; de Graaff et al. 2024).Overall, these observational results appear in good agreement with TNG50 galaxies at z = 6: the baryonic mass fractions are approximately equally low at low stellar masses f baryon (< 1 kpc; M * = 10 8−9 M ⊙ ) ∼ 0.2 − 0.4.The relative mass fraction among the baryonic components (stars and gas) also agrees at least qualitatively.We find M * (< 1 kpc) ≪ M gas (< 1kpc) for both the simulated and observed galaxies at z ∼ 6 and M * < 10 9 M ⊙ .Quantitatively, on the other hand, we find that -at face value -the stellar-to-gas mass ratio is a factor two higher in the TNG50 simulated systems compared to the observations.However, we deem this comparison inconclusive as this mass ratio in TNG depends on the resolution of the simulation (at the 20 − 40% level across a factor of 16 in mass resolution at z ∼ 5 − 6 and M * ∼ 10 8−9 M ⊙ ).
Understanding the cold-gas content in these galaxies is of broader interest, as it speaks to the efficiency of star formation (M gas /SFR) at high redshift.Indeed, de Graaff et al. (2024) had estimated the gas content by assuming present-epoch scaling relations for the star formation efficiency, and could therefore be underestimated by a factor 3. However, this would counter the general observed trend that the depletion time decreases toward higher redshifts (Tacconi et al. 2020).It would also imply a discrepancy in the central gas masses at the level of 1 dex with respect to TNG50.
This leaves only two plausible scenarios for the observations of the mass budget at z ∼ 6. (i) The inferred low baryon fractions imply that these high-redshift galaxies -likely progenitors of massive galaxies at z < 2 -are dark matter-dominated even within 1 kpc.(ii) Alternatively, the ionized gas kinematics, on which the observational mass budget estimates rely, could be dominated by non-gravitational motions due to outflows from feedback, most probably of stellar origin.If the resulting velocity gradients and dispersions are biased high, the inferred dynamical mass will be substantially overestimated, and hence the baryon fraction will be underestimated.
The fact that the TNG50 model is in agreement with the first scenario of a dark matter-dominated phase early in the galaxy mass assembly history shows that this first interpretation is at least plausible.Furthermore, the simulation offers insight into the evolution of the descen-dants of z = 6 galaxies with M * = 10 8−9 M ⊙ , and shows that these dark matter-dominated systems are indeed progenitors of more massive (M * ∼ 10 10 M ⊙ ) galaxies at z ∼ 2. In addition to the rapid stellar mass growth, and consistent with observations at cosmic noon, these systems are baryon-dominated in the center and have higher central stellar masses than gas masses, indicating a rapid period of stellar mass assembly in the center through star formation or accretion.Moreover, we find that these galaxies typically grow to MW-mass systems at z ∼ 0 (M * ∼ 10 10.5 M ⊙ ).In qualitative agreement with observations (e.g.Martinsson et al. 2013), the central kpc is dominated by the stellar component with minimal dark matter (f DM (< 1 kpc) ∼ 0.2), whereas the dark matter fraction within the half-mass radius is substantially higher (f DM (< r e ) ∼ 0.4 − 0.5).
However, we cannot rule out at present the second scenario in which the dynamics are biased by nongravitation motions.Further observations using integral field spectroscopy to robustly resolve the velocity fields of these high-redshift systems will be critical to distinguish between the two scenarios.
The comparison in the mass budget between the highredshift observations and simulations demonstrates the great potential constraining power of such observations.The TNG50 simulation makes clear predictions for the central baryon fraction and its dependence on stellar mass, and (lack thereof) on redshift.This is similarly true for the stellar-to-gas mass ratio, which also appears more sensitive to the details of the simulation (resolution, subgrid model).Therefore, systematically mapping these mass ratios as a function of stellar mass and redshift through observations would provide direct insight into galaxy mass assembly.Current surveys with JWST that simultaneously probe the stellar population, ISM and kinematic properties of galaxies at z > 1 across a wide range in stellar mass provide an ideal dataset for such a study (Maseda et al. 2024;Eisenstein et al. 2023).Combined with a careful modeling analysis of observational biases by the use of mock observations from simulations, these observations will be able to place strong constraints on theoretical models and constrain the fundamental question: how efficient was star formation in the early universe?

Figure 1 .
Figure1.Baryonic mass fraction within the stellar halfmass radius (TNG50) and the emission line half-light radius (JWST; de Graaff et al. 2024) of galaxies vs. total stellar mass.Throughout, by baryonic fraction we mean the ratio between the mass in gas and stars and the total mass of the galaxy (see Eq. 1).The TNG50 simulation predicts low baryon fractions for M * < 10 9 M⊙ and a sharp increase in f baryon toward higher masses.The half-mass radii are small (∼ 1 kpc in physical units) and, at fixed stellar mass, are anti-correlated with the baryon fraction.With the exception of one source, which may have high f baryon due to an underestimated dynamical mass, the baryon fractions measured at z ∼ 6 with JWST are in good agreement with the simulated galaxies.

4. 2 .Figure 3 .
Figure3.Evolution in the baryonic mass fraction (left) and stellar-to-gas mass ratio (right) within a fixed spherical aperture of radius 1 kpc according to the IllustrisTNG and EAGLE cosmological hydrodynamical galaxy simulations.Colored lines show the median f baryon and M * /Mgas ratio within narrow stellar mass ranges for central galaxies in the TNG50 (dashed) and TNG100 (dash-dotted) simulations and, for comparison, also the EAGLE Ref-L100N1504 simulation (dotted).Solid black lines show the median time evolution of the descendant population of low-mass (M * ∼ 10 8−9 M⊙) galaxies selected from TNG50 at z = 6, with shaded regions indicating the 16-84 and 5-95 percentiles.As the z = 6 galaxies grow in stellar mass, f baryon and M * /Mgas increase rapidly, resulting in stellar mass-dominated centers by z ∼ 2.

Figure 4 .
Figure4.Redshift evolution of the dark matter fraction within the stellar half-mass radius (TNG50) and the half-light radius (observations).The dark (light) shaded regions indicate the 16-84 (5-95) percentiles of the evolution of TNG50 galaxies selected at z = 6 with 10 8 M⊙ < M * < 10 9 M⊙.The thick solid line shows the median evolution of fDM(< re) of the TNG50 galaxies; to bracket possible systematic differences in the comparison to observations, we also show the median evolution of fDM(< 0.5re) (dotted) and fDM(< 1.5re) (dashed).The median evolutionary tracks and data points from high-redshift observations are colorcoded by the total stellar mass.Qualitatively, this suggests a link between observations at different mass and redshift epochs, and agreement between the TNG50 simulation and observations.