Probing the dawn of galaxies: star formation and feedback in the JWST era through the GAEA model

TL;DR

This study tests the ability of the GAEA semi-analytic model, run on Planck-Millennium merger trees, to reproduce the high-redshift galaxy population in the JWST era. It evaluates the GSMF, UVLF (including AGN contributions and dust effects), and the MZR up to $z\sim13$, identifying tensions at $z\gtrsim10$ and exploring physical variants that could reconcile theory with observations. The authors show that the fiducial model matches many observables up to $z\lesssim10$ but underpredicts the number of the brightest galaxies at $z>10$; two variants—feedback-free starbursts at $z>10$ and suppressed high-$z$ stellar feedback—can raise counts to align with data, though each leaves distinct imprints on metallicity and star formation histories. The work highlights the potential for discriminating between high-$z$ scenarios using SFR–$M_*$ relations and the MZR, and underscores the need for improved modeling of dust, Pop III contributions, IMF evolution, and mini-halos to fully characterize early galaxy formation.

Abstract

The James Webb Space Telecope (JWST) opened a new window for the study of the highest redshift ($z>7$) Universe. This work presents a theoretical investigation of the very-high redshift Universe using the state-of-the-art GALaxy Evolution and Assembly (GAEA) model, run on merger trees from the Planck-Millennium $N$-body simulation. We show that GAEA successfully reproduces a wide range of high-$z$ observational estimates including: the galaxy stellar mass function up to $z\sim13$ and the total (galaxies and AGN) UV luminosity function (LF) up to $z\sim10$. We find that the AGN UV emission represents an important contribution at the bright end of the UVLF up to $z\sim8$, but it is negligible at higher redshift. Our model reproduces well the observed mass-metallicity relation at $z\leq4$, while it slightly overestimates the normalization of the relation at earlier cosmic epochs. At $z\geq11$, current UVLF estimates are at least one order of magnitude larger than model predictions. We investigate the impact of different physical mechanisms, such as an enhanced star formation efficiency coupled with a reduced stellar feedback or a negligible stellar feedback at $z>10$. In the framework of our model, both the galaxy stellar mass and UV luminosity functions at $z\geq10$ can be explained by assuming feedback-free starbursts in high-density molecular clouds. However, we show that this model variant leads to a slight increase of the normalization of the $z\geq10$ mass-metallicity relation, strengthening the tension with available data. A model with negligible stellar feedback at $z>10$ also predicts larger numbers of massive and bright galaxies aligning well with observations, but it also overestimates the metallicity of the interstellar medium. We show that these model variants can in principle be discriminated using the relation between the star formation rate and galaxy stellar mass.

Figures (15)

• Figure 1: Rest-frame galaxy UV Luminosity Function at $0\leq z \leq 14$ from GAEA, compared with observational estimates finkelstein2015hagen2015parsa2016bouwens2021apage2021castellano2023perez-gonzalez2023harikane2022aharikane2023aharikane2024adonnan2024mcleod2024franco2025whitler2025. The solid and dotted lines represent model predictions with and without dust attenuation, respectively. The grey horizontal lines show the space density corresponding to 10 sources in the PMS volume. The grey and orange shaded areas show the effect of cosmic variance on the dust-attenuated and the unattenuated UVLF, respectively (see text for details). Different symbols refer to different techniques for computing the LF: circle for the $1 / V_\mathrm{max}$ method, plus sign for the effective volume approach, cross sign for the STY method, diamond for the Bayesian approach, pentagon for the SWML method, square for the binned estimation of the LF. Arrows in downward direction are referring to upper limits, while arrows in upward direction to lower limits.
• Figure 2: Total (galaxies + AGN) rest-frame UV luminosity function at $3 \leq z \leq 9$ as predicted by GAEA, compared to observational data adams2023bcastellano2023perez-gonzalez2023donnan2024harikane2022aharikane2023aharikane2024a. Symbols and lines are as in \ref{['fig:GAEA_galUVLF_fiducial']}. The sub-panels show the ratio between the AGN UVLF and the total UVLF at each redshift, with (solid) and without (dashed) dust attenuation. The orange shaded areas show the effect of cosmic variance on the dust-unattenuated UVLF (see text for details).
• Figure 3: The galaxy stellar mass function at $4 \leq z \leq 14$ from the GAEA semi-analytical model used in this study. Symbols with associated error bars correspond to the observational measurements from harvey2025 for $7 < z < 13$, stefanon2015 for $4 < z < 7$, stefanon2021 for $6 < z < 10$, weaver2023 for $4 < z < 7$, muzzin2013 for $z = 4$, and shuntov2025 for $4 \leq z \leq 12$. Solid lines represent the intrinsic model predictions, while dashed and dotted lines are obtained assuming an uncertainty on stellar mass of 0.25 and 0.5, respectively. The green and grey shaded areas show the impact of cosmic variance on the intrinsic model predictions when survey sky areas from harvey2025 and the COSMOS2025 field are considered, respectively. The dotted horizontal grey lines mark the number densities corresponding to $10$ galaxies within the simulated volume.
• Figure 4: Stellar-to-Halo Mass Relation at $0 \leq z \leq 10$ from GAEA, compared with observational estimates stefanon2021shuntov2025paquereau2025 and theoretical work (lines with shaded areas). We include both empirical abundance matching models (MEAM from moster2013; UniverseMachine from behroozi2019) and the Santa-Cruz SAM yung2024a. The solid lines represent the median relations (the black ones refer to GAEA), while the shaded regions show the $1 \sigma$ confidence regions.
• Figure 5: Mass-metallicity Relation (MRZ) up to $z = 10$. GAEA model predictions (lines) are compared with observational data (as indicated in the legends within each sub-panel). The black solid line is the median value of the relation, while the shaded areas show the boundaries between the $16$-th and the $84$-th percentiles (dark green) and the $1$-st and $99$-th percentiles (light green). The predicted MZR has been downshifted by 0.2 (see text). The dot-dashed lines are the predicted MZR at fixed SFR. The black dashed lines show predictions from the model published in delucia2024, with the only difference being the chemical enrichment in low-mass haloes (see text for details).