Gravitational Wave Measurement of the Mbh-Mbulge Intrinsic Scatter at High Redshift

Abstract

The observed GWB spectrum is higher in amplitude than model predictions by a factor of 2-3. Using a semi-analytic model, we evaluate the effect of a high-scatter supermassive black hole (SMBH) scaling relation (Mbh-Mbulge) on models of the nanoHertz gravitational wave background (GWB). By implementing an intrinsic scatter of the Mbh-Mbulge relation, which is larger at higher redshift, but matches local observations, we find that the amplitude of GWB models increases to be consistent with the low-frequency end of the GWB spectrum. This amplitude increase is not uniform across frequencies, a strongly evolving scatter preferentially increases the number density of the most massive SMBHs which, in the GWB spectrum, minimizes the strength of the low-frequency turnover. Our models with positively evolving intrinsic scatter can reproduce the electromagnetically observed overmassive SMBHs at 4 < z < 6 without changing the Mbh-Mbulge normalization though we find that including moderate normalization evolution improves fits to the GWB data. We conclude that the Mbh-Mbulge relation which best describes the available GWB and electromagnetic data sets has intrinsic scatter which evolves as epsilon(z) = epsilon_0 + (0.56 +/- 0.4) log10(1 + z) and normalization which evolves as alpha(z) = alpha_0 (1 + z)^(0.84 +/- 0.35). The results of this work imply that the Mbh-Mbulge relation we see today is not universal throughout cosmic time and that a diversity of seeding models and growth mechanisms may be at play in the early stages of SMBH-galaxy evolution.

Figures (9)

• Figure 1: Here we show the separate effects of an evolving $M_\mathrm{BH}$--$M_\mathrm{bulge}$ scatter (top) and normalization (bottom) on the $z = 2$ BHMF (left) and the GWB spectrum (right). We keep all parameters in the models constant and change only the power law for evolving scatter, $\varepsilon_z$ (or that for evolving normalization, $\alpha_z$) where $\varepsilon_z, \alpha_z = 0$ indicates no redshift evolution and $\varepsilon_z, \alpha_z > 0$ indicates positive redshift evolution. Both types of evolution increase the number density of the most massive SMBHs where the amount of the increase is correlated with the strength of evolution, but only an evolving normalization significantly increases the number density of SMBH masses $M_\mathrm{BH} < 10^{9} M_\odot$. Gravitational wave frequency is inversely proportional to SMBH chirp mass. Therefore the relative changes to number densities across the mass range have corresponding impacts on the GWB spectrum. Both types of evolution increase the amplitude at the lowest frequencies (high masses), but only $\alpha_z$ affects the high frequency (low mass) end of the spectrum.
• Figure 2: Left: The posterior distributions for $\varepsilon_z$ from model M$_{\varepsilon}$ (light blue, solid) and model M$_{\varepsilon +}$ (dark blue, dashed). The vertical lines represent the median of the posteriors, for model M$_{\alpha}$ (red, dotted) $\varepsilon_z$ was fixed to a value of 0.5 so there is no distribution associated with the median line. The median value of model M$_{\varepsilon}$ is higher than that of model M$_{\varepsilon +}$ because it had no additional free parameters to increase the amplitude and so $\varepsilon_z$ had to compensate more. Right: The GWB spectra associated with our three models. We see that all models have similar amplitudes and shapes which are all within the violins of the data. The low-frequency turnover is present in all three models though it is not strong. Notably, model M$_{\varepsilon}$ lies below the other two models and the bulk of the data only truly aligning with the lowest frequency bin. Interestingly, the spectra of model M$_{\varepsilon +}$ and model M$_{\alpha}$ are quite similar despite model M$_{\varepsilon +}$ having no evolution in the $M_\mathrm{BH}$--$M_\mathrm{bulge}$ normalization.
• Figure 3: Here we show how the GWB amplitude of a model correlates with the sampled value of $\varepsilon_z$ for model M$_{\varepsilon}$ and model M$_{\varepsilon +}$ measured at the lowest frequency bin $f_\mathrm{gw} \sim 2$ nHz. For this plot we measure the amplitude in the lowest frequency bin since that is where the impact of $\varepsilon_z$ is greatest. The dark blue points with large spread are from model M$_{\varepsilon +}$ where the blue squares represent their median binned along the x-axis. The light blue line is from model M$_{\varepsilon}$, which only sampled $\varepsilon_z$. On average, there is a clear positive correlation between the GWB amplitude and the value of $\varepsilon_z$. We allowed many parameters to vary in model M$_{\varepsilon +}$ and so the correlation between GWB amplitude and $\varepsilon_z$ is less direct. It is still present (blue squares), but there is a larger spread. Other parameters in the model have a strong correlation with the GWB amplitude (e.g. local $M_\mathrm{BH}$--$M_\mathrm{bulge}$ normalization, $\alpha_0$) so the model does not need to rely solely on $\varepsilon_z$ to reproduce the GWB amplitude. This effect explains why the posterior distribution for model M$_{\varepsilon +}$ favors lower values with a larger spread than that of model M$_{\varepsilon}$.
• Figure 4: A comparison between our result predictions and observed SMBH and galaxy masses. In gray, we show the 1-, 2-, and 3-$\sigma$ regions of scatter around the $M_\mathrm{BH}$--$M_\mathrm{bulge}$ relation from Kormendy_Ho_2013 with three different options for evolution. We add 0.4 dex in quadrature to the intrinsic scatter to approximate measurement error. In all panels the solid black line represents the local $M_\mathrm{BH}$--$M_\mathrm{bulge}$ relation and the dashed line is the location of the evolving relation. The pink circles in each panel represent EM-observed SMBH--galaxy measurements from the literature. From left to right we include data from Kormendy_Ho_2013, Ding_2020, Zhang_2023, and the final panel includes data from Sun_2025_corrected, Li_2025_normal, Jones_2025, and Brooks_2025. Top row: The intrinsic scatter of the relation evolves according to the median posterior value from model M$_{\varepsilon +}$, $\varepsilon_z = 0.56$. The normalization is unchanging. Middle row: Both the intrinsic scatter and normalization are evolving according to the results from model M$_{\alpha}$, $\varepsilon_z = 0.5$, $\alpha_z = 0.84$. Bottom row: Only the normalization is evolving here. This represents the best-fit model from Matt_2026. All three cases predict a population of overmassive SMBHs consistent with high-redshift JWST observations without violating local constraints. The two options with evolving scatter also predict an undermassive population, which can reproduce the "normal mass" SMBHs found by Li_2025_normal. This lower-mass population is not predicted by the high normalization, low scatter model.
• Figure 5: The value of scatter in the $M_\mathrm{BH}$--$M_\mathrm{bulge}$ relation versus redshift for the data in \ref{['fig:mmb_obs_comp']}. The lines represent the predicted evolution of $\varepsilon(z)$ using the median of the posterior for $\varepsilon_z$ in model M$_{\varepsilon}$ (upper solid line) and the fixed value in model M$_{\alpha}$ (lower dotted line) with their 68% confidence intervals indicated by the shaded regions. We estimate the intrinsic scatter of the EM data in our four redshift bins by calculating the standard deviation of the distance from the local Kormendy_Ho_2013 relation and subtracting the measurement error in quadrature (pink points). These two models bracket the estimated intrinsic scatter evolution in the EM data which lie closest to the line from model M$_{\varepsilon}$.