GoKawiilIn this paper, uncertainties are given as 1σ or 68% confidence intervals. Upper limits are indicated at the 2σ level unless otherwise stated. We adopt cosmological parameters measured in ref. 51, that is, a Λ cold dark matter (ΛCDM) model with total matter density in units of the critical density, Ω m = 0.310, and Hubble constant, H 0 = 67.7 km s−1 Mpc−1.
Spectroscopic sample
This study makes use of the public JWST data collected as part of several observational programmes with the NIRSpec spectrograph52 with PIDs: 1345 (CEERS)53, 1181 (JADES)54, 1210 (JADES)55,56, 2674 (PI A. Haro)57, 3215 (JADES Origins Field)58, 4106 (PI E. Nelson)59, 4233 (RUBIES)18, 2565 (PI K. Glazebrook)60, 2750 (PI A. Haro)61 and 6541 (PI E. Egami)62. These observations have been uniformly reduced and published as part of the DAWN JWST Archive (https://dawn-cph.github.io/dja) (DJA)18,19. Using the v.3 reductions in DJA, we selected all galaxies observed in the medium-resolution grating spectra with a broad Hα component and a spectroscopic redshift produced by msaexp (https://github.com/gBrammer/msaexp)63. To this, we added JWST broad-line objects reported in the literature with publicly available data and processed in DJA. We selected objects with a full width at half-maximum (FWHM) linewidth greater than about 1,000 km s−1 from the objects in the archive. We then selected spectra with high SNR (median SNR >5 per 10 Å for the continuum-subtracted region ±2,000 km s−1 around the Hα line) and also included objects for the stacked spectrum using broad-line objects with lower SNR (5 > SNR/10 Å > 1) to ensure that we are not biased by our SNR selection.
We note that our sample spans a range of colours when using the existing selection criteria (Extended Data Fig. 2). Although only three objects of 12 are classified as LRDs (AEH), more objects have a similar inflection point around the Balmer series limit (CDGH, as seen in the PRISM spectra in Extended Data Fig. 1) which is not picked up by the selection criteria due to redshift effects and the contribution of the strong optical emission lines making the colour gradient flatter. However, some objects are bluer than photometrically selected LRDs (CFI). Therefore, although most of our sample has classical LRD-like features, some objects probably have a wider range of properties, possibly produced by differing extinctions or contributions from the host galaxy. This difference can possibly be explained, in part, by the incompleteness of current photometric selection criteria in bluer F277W–F444W colours and fainter rest-UV magnitudes64. Despite this range of colours, the presence and magnitude of electron scattering by the ionized gas cocoons does not depend on the location in this colour space (as suggested by the optical depth in Extended Data Table 1 or exponential width in Extended Data Table 3).
Emission line models
All best-fit results in this paper were produced using the Monte Carlo Markov Chain (MCMC) NUTS sampler as part of the package PyMC v. 5.17.0 (ref. 65), except objects A and D, which were fitted using the Ensemble sampler emcee v. 3 (ref. 66) (due to an incompatibility of the P Cygni model used here with the tensor formalism in PyMC). We sampled the posterior distributions with 4k walkers (where k is the number of free parameters) and 105 samples per walker. We use the mode values of the posterior parameter distributions as the best values and the 68% highest-density interval as the range of uncertainty. Finally, we find that the resolution of NIRSpec gratings is higher than the nominal value. Using the resolved widths of narrow [O iii] lines in the high-resolution G235H grating of object A, we estimate that the medium-resolution grating G235M has about 1.7 times higher R than the nominal value. This scaling factor on R has been assumed for G395M and also for the G395H grating, which agrees with the results of forward modelling of the NIRSpec instrument response for point sources67.
In modelling the broad Hα profile, we assumed a broadening mechanism: either a Doppler velocity broadening or a Compton scattering broadening. The former is modelled using a Gaussian function \(f(\lambda ;A,\mu ,\sigma )=A\,\exp ({(\lambda -\mu )}^{2}/(2{\sigma }^{2}))\) with amplitude A, line centre μ and velocity dispersion σ. For the Compton-scattered profile we use a symmetric exponential \(g(\lambda ;B,{\lambda }_{0},W)=B{{\rm{e}}}^{-| \lambda -{\lambda }_{0}| /W}\) with amplitude B, line centre λ 0 , and e-folding scale W.
First, we tested both broadening mechanisms by fitting two sets of models. Extended Data Fig. 3 shows the comparison between the broadening models. To fit the data reasonably well, the models also included (narrow) Gaussians for the host galaxy Hα and [N ii]doublet with fixed centroids and velocity dispersions tied to the same value and limited to <1,000 km s−1. The ratio of the two [N ii] lines is set to 0.33 (ref. 68). In some cases, an additional Gaussian absorption component (object E) or a P Cygni profile (objects A and D) are required to accurately model the broad Hα component. We note here that the [N ii] lines are only required in the fits in objects B and G, which fitted better than redshifted absorption. Whether this is an artefact related to the complex spectral shape due to possible self-reversal in the line, or whether [N ii] really is observed in these cases, is unclear and would require higher resolution spectra and more sophisticated models to establish. Finally, wavelength regions around the emission lines [O i] λ6,302, He i λ6,678 and [S ii] λλ6,717, 6,731 are excluded from the fit.
Although most Hα lines in the sample are predominantly exponential with very high statistical significance (Extended Data Fig. 3), to reconstruct the intrinsic Doppler widths, we model the lines with a Gaussian convolved with an exponential, instead of a pure exponential. These models are convolved with the instrumental resolution of the relevant gratings (which were taken from the JWST JDox website for NIRSpec) at the Hα peak (we assume the actual resolution is about 1.7 times better than the nominal; see the description above). To alleviate the complexity of some of our models and more accurately constrain the narrow Hα components, we use the velocity widths of the optical [O iii] lines as a Gaussian prior on Hα and [N ii] widths, in which relevant spectral coverage is available. The best-fit profiles and intrinsic Doppler components with their widths are presented in Fig. 3 and Extended Data Table 3.
We also test whether, for example, gas turbulence or Raman scattering could be responsible for line broadening by comparing a basic Lorentzian27,69 and an exponential line shapes. The former is defined as a symmetric profile with FWHM 2γ centred at λ 0 : \(h(\lambda ;C,{\lambda }_{0},\gamma )=C\frac{\gamma }{{(\lambda -{\lambda }_{0})}^{2}-{\gamma }^{2}}\). The exponential is a significantly better fit in most objects or an equivalent fit in objects H and L (Extended Data Table 2 and Extended Data Fig. 5). This indicates that any potential contribution from turbulence broadening is not significant, we do not assume a more physically motivated profile of a Gaussian convolved with Lorentzian (that is, a Voigt profile). Another reason to exclude this model is that the ratios of the line widths between Hα, Hβ and Paβ, for example, are expected to differ by about a factor of 2–3 in velocity27, something which is not generally observed in these types of objects. Furthermore, turbulent broadening may sometimes result in enhanced red-wing profiles69, whereas all objects here have symmetric Hα wings (Extended Data Fig. 3). However, although it has been argued that the lines cannot be Lorentzian on this basis29, optical depth effects, which would be quite different for the different lines, could affect the relative line widths and more careful non-LTE radiative transfer analysis would help explain this issue.
... continue reading