A dormant overmassive black hole in the early Universe

Data description and reduction

All data used herein have been obtained from the JADES survey, the full description of which is available in ref. ³⁸. The spectroscopic survey consists of several tiers, characterized by their depth, ‘Medium’ or ‘Deep’, photometry from which targets were selected, HST or JWST, and the observed field, GOODS-S or GOODS-N. Here we make use of data from the Medium/HST in GOODS-N tier, which consists of a single NIRSpec spectrum in prism and R1000 gratings. We also make use of the accompanying NIRCam wide-band imaging data.

NIRSpec

The full description of the NIRSpec data used is available in ref. ². A summary is provided here for completeness. The observations in the Medium/HST tier in GOODS-N consisted of three medium-resolution gratings (G140M/F070LP, G235M/F170LP and G395M/F290LP) and low-resolution prism. The exposure time was 1.7 h per source in the prism and 0.8 h per source in the medium gratings. The data were processed according to the procedures laid out in ref. ³⁹ and other similar JADES papers, such as ref. ⁴⁰. A full description of the data reduction procedure will be presented in Carniani, S. et al. (in preparation). Here we note that the spectral data were reduced using the pipeline developed by the NIRSpec GTO team and the ESA NIRSpec Science Operations Team. As the primary interest of this study was the properties of the central, unresolved region containing the AGN, we use the one-dimensional spectra extracted from the central 3 pixels, corresponding to 0.3′, of each two-dimensional spectra. Path-loss corrections were calculated for each observation, taking into account the intra-shutter position, assuming a point-source geometry and a 5-pixel extraction box. Owing to the compact nature of the object, to maximize the signal-to-noise ratio (S/N), we used 3-pixel extractions. Although a 3-pixel box is not the extraction box we optimized the path-loss corrections for, we compared directly the two spectra and found no systematic difference within the uncertainties.

The full prism spectrum obtained is shown in Extended Data Fig. 1.

NIRCam

The imaging data consisted of seven wide (F090W, F115W, F150W, F200W, F277W, F356W and F444W) and one medium (F410M) filter bands of the NIRCam instrument in the GOODS-N field. We also used imaging in the F182M and F210M medium bands. The photometric data reduction procedure is presented in refs. ^41,42,43 with a full description to be made available in Tacchella et al. (manuscript in preparation). In summary, we used v.1.9.2 of the JWST calibration pipeline⁴⁴ together with the CRDS pipeline mapping context 1039. Stages 1 and 2 of the pipeline were run with our own sky-flat provided for the flat-fielding, otherwise keeping to the default parameters. After stage 2, custom procedures were performed to account for 1/f noise and subtract scattered light artefacts, ‘wisps’, along with the large-scale background. Astrometric alignment was performed using a custom version of JWST TweakReg, with corrections derived from HST F814W and F160W mosaics along with GAIA Early Data Release 3 astrometry. The images of individual exposures were then stacked in stage 3 of the pipeline with the final pixel scale being 0.03′ per pixel.

Spectral fitting and further spectral analysis

To identify the broad component in the Hα line, we used a Bayesian method to model it with two components models—one containing only narrow emission in the Hα line, [Nii]λ λ6548,6583 and [Sii]λ λ6716,6731 doublets, the broad-line model included a broad component in the Hα line. Narrow-line widths were constrained to be the same for every line, the ratio of the [NII] doublet fluxes was fixed to 3 (as from their Einstein coefficients ratios), the [Sii] doublet fluxes remained independent but constrained to be within the flux ratio 6,716/6,731 range expected in the low- and high-density regimes (0.45–1.45; ref. ⁴⁵). The priors on the peak widths were uniform with the fitted full width at half-maximum (FWHM) ranging between 700 and 1,500 km s⁻¹ for the narrow, and between 1,500 and 11,500 km s⁻¹ for the broad component, with the lower bound set by instrument resolution. The posterior is estimated with a Markov Chain Monte Carlo integrator⁴⁶. Redshifts for the narrow peaks and the BLR were fit independently with priors being set to narrow Gaussians centred on the overall redshift obtained through visual inspection and widths inferred from the pixel scale in redshift space. Line peak heights used log-uniform priors.

The performance of the two models was quantified using the Bayesian information criterion (BIC), defined as

where k is the number of free parameters and n is the number of data points fitted. Following the criteria in ref. ², we require ΔBIC = BIC_Narrow − BIC_Narrow+Broad to be above 5 for robust detection. Moreover, we require the fitted broad component to have a significance of at least 5σ. The significance of the broad component of GN-1001830 was found to be 9σ and ΔBIC was 31. The full summary of the fit is shown in a corner plot in Supplementary Fig. 2. As shown in this figure, the data are quite constraining on the different components with no significant degeneracies between them.

We note that the spectrum shown in Fig. 1 seems to have a feature at λ ≈ 4.94 μm, the origin of which is currently unclear. In the two-dimensional spectrum, this feature is offset by 1–2 pixels, suggesting that it might be Hα line of a foreground galaxy at z ≈ 6.53. However, one of the four exposures has an outlier at this location (Extended Data Fig. 2); although this has been masked, the feature might be a residual artefact. For these reasons, we choose to mask this feature in the final fit, although its inclusion would not significantly affect the results. We also note that the Hβ line in Extended Data Fig. 1 seems to contain a similar artefact, which manifests as an apparent broad wing. However, the inferred 2σ significance of the feature leads us to conclude that its origin is noise.

To check for evidence of outflows, we fit the Hβ line together with the [Oiii] doublet in the medium-resolution data. For this purpose, we fit these lines first with single components constrained to have the same width, with the ratio of [Oiii] doublet peaks fixed at 3, then introduce a broader outflow component into each line, and finally, we fit a broad component to Hβ. We find that single narrow component fits are preferred for each line, showing no evidence for outflows or a broad component in Hβ, as can be seen in Extended Data Fig. 3. The measured FWHM of the narrow lines in R1000 was \(22{5}_{-11}^{+11}\,{\rm{km}}\,{{\rm{s}}}^{-1}\), when corrected for instrumental broadening, using the point-source line spread function (LSF) models in ref. ⁴⁷. This broadening is 180 km s⁻¹ in the wavelength range considered.

We also use HST imaging of the source carried out in 2018 to check its variability and thus the presence of supernovae that could potentially produce broad Hα emission while in their nebular phase. Subtracting HST and JWST images taken in equivalent filters shows that the flux of the source varied by no more than 5% over in the 4-year period in the observed frame, corresponding to <5% variability over a 6-month period in the frame of the object. This conclusively shows that our observed broad component cannot be the product of a recent supernova explosion in the galaxy. We also rule out the possibility that the broad-line component may be attributed to the effect of multiple supernovae given the low SFR of the galaxy.

Furthermore, we leverage the higher depth of the prism to fit weaker emission lines, in particular, Hγ, [Oiii]λ4363, [Ne iii]λ3869 and [Oii]λ3727, as the [Oiii] auroral line may be used to constrain metallicity, whereas the remaining lines are candidate diagnostics for type 2 AGN. Each line has been fit with a single Gaussian profile together with Hβ and the [Oiii] doublet with redshift and FWHM fixed by the latter. We note that, although the [Oii]λ3727 line is part of a doublet, this doublet is completely blended in the prism spectrum and the overall detection significance ends up being marginal (Extended Data Table 1). All narrow emission lines fitted are summarized in Extended Data Table 1. It should be noted that we carry out our fitting using spectra extracted from the central 3 pixels of the source to enhance the S/N of lines as the region of emission is compact because of the AGN nature of the source. However, the flux corrections applied by the reduction pipeline are geared towards the 5-pixel extracted spectra; we thus redo our fits using the 5-pixel spectra and find that the derived line fluxes differ by less than 1σ.

Our tentative (3.6σ) detection of [Oiii]λ4363, together with the [Oiii]λ λ5007,4959 doublet, can be used to estimate electron temperature based on their ratio. The electron density cannot be reliably estimated from the spectrum because of the lack of an [Sii]λλ6716,6731 detection. However, we note that the inferred electron temperature is relatively insensitive to a density between 100 cm⁻³ and 10,000 cm⁻³, typical of narrow-line regions of AGN⁴⁸. We thus assume electron density of the order of 1,000 cm⁻³ and calculate the electron temperature to be \({T}_{{\rm{e}}}\approx {25,\,000}_{-4,800}^{+3,200}\,{\rm{K}}\). Assuming the main ionization mechanism to be AGN activity, we follow the methods in ref. ⁴⁹ to derive the contributions to metallicity from different ionic species of oxygen—\(12+\log ({{\rm{O}}}^{++}/{\rm{H}})=7.2{8}_{-0.1}^{+0.16}\) and \(12+\log ({{\rm{O}}}^{+}/{\rm{H}}) < 6.13\), the latter value presenting an upper limit due to the low detection significance of [Oii]λ3727. The contribution of higher ionization oxygen species is probably negligible because of the lack of [Oiv] and Heii line detections. Heii lines have similar ionization potential to [Oiv], but helium is much more abundant; thus, the lack of detection of Heii lines implies that the radiation is not sufficiently hard to produce significant amounts of highly ionized oxygen. The final oxygen abundance ratio estimate is thus \(12+\log ({\rm{O}}/{\rm{H}})=7.3{2}_{-0.10}^{+0.16}\), corresponding to Z ≈ 0.04Z_⊙. However, the [Oiii]λ4363 line is heavily blended with Hγ in prism. Therefore, the simple two-Gaussian fit may underestimate the relevant uncertainties. As a conservative estimate, we obtain a lower limit on the [Oiii]λ4363 flux by fitting the blended lines with a single Gaussian profile and subtracting from its flux the Hγ flux obtained in a fit without [Oiii]λ4363 included. This method resulted in \({F}_{[{\rm{OIII}}]\lambda 4363}\ge 9.1\times 1{0}^{-20}\,{\rm{erg}}\,{{\rm{s}}}^{-1}\,{{\rm{cm}}}^{-2}\). Repeating the former analysis gives T_e ≥ 15,000 K and 12 + log(O/H) ≤ 7.72. A lower limit on metallicity was derived by fitting the blended feature with [Oiii]λ4363 only, which results in T_e ≤ 34,000 K and \(12+\log ({\rm{O}}/{\rm{H}})\ge 7.08\). These limits are consistent with the best-fit value and place our source below the mass–metallicity relation^50,51 at similar redshifts (Supplementary Fig. 3).

As a final check on these results, we fit the blended [Oiii]λ4363 and Hγ feature by fixing the ratio of Hβ and Hγ to the appropriate Balmer decrement and find \(12+\log ({\rm{O}}/{\rm{H}})=7.2{5}_{-0.17}^{+0.20}\), which is completely consistent with the above estimates.

Morphological and photometric fitting and stellar population properties

To constrain the stellar mass and SFR of the host galaxy, we use fractional spectral energy distribution (SED) fitting. To do this, we use the tool ForcePho (Johnson, B., manuscript in preparation), which enables us to forward model the light distribution using a combination of Sérsic profiles. We perform spatially resolved photometry with ForcePho following the methodology detailed in refs. ^42,43. In short, we model the AGN and host galaxy as a central point-source component and underlying host galaxy component, respectively, and fit the light distribution in the individual exposures of all 10 NIRCam bands simultaneously. This enables us to obtain accurate spatially resolved fluxes and morphological parameters for the galaxy. This approach has been used previously in refs. ^2,42,43.

Extended Data Fig. 5 shows the data, residual, model and point-source-subtracted host galaxy for the ForcePho fit. We can see that the galaxy + point-source model has fit the data well without leaving significant residuals and the galaxy component appears bright enough for reliable photometry, with S/N ranging from 6 to 40 across our filters. The resulting best-fit morphological parameters are reported in Extended Data Table 2, which shows that the host galaxy is compact (R_e = 137 ± 8 pc) with a disk-like profile (Sérsic index n ≈ 1). The quoted statistical-only error on R_e is rather small, considering the marginally resolved nature of the source.

The PSF model that is approximated by ForcePho is based on WebbPSF—incorporating forward modelling of the optics of the telescope, with additional calibration provided by field stars. To investigate the uncertainty coming from the PSF approximation used by ForcePho, we also re-fit the data with a different PSF approximation, which includes charge transfer effects, and obtain a 16% smaller radius. We thus adopt a 16% systematic error floor, which results in a final estimate of R_e = 137 ± 23 pc. The fluxes for the point-source component and host galaxy can be seen in Extended Data Fig. 6. We find that the galaxy component dominated (at the 90% level) in all filters except for the F444W, for which the contribution from the AGN broad Hα dominates.

The next stage is to fit this SED to obtain the stellar population properties of the host galaxy. To do this, we use the Bayesian SED fitting code Prospector²⁰, which uses Flexible Stellar Population Synthesis⁵² (FSPS) with MIST isochrones⁵³, nebular line and continuum emission⁵⁴ and a Chabrier⁵⁵ initial mass function (IMF). We also use the Bayesian SED fitting code Bagpipes¹⁹, which uses stellar population models in ref. ⁵⁶, alongside the nebular line and continuum emission⁵⁴ with a Kroupa⁵⁷ IMF. For both codes, we include a flexible two-component dust model following ref. ⁵⁸ consisting of a birth cloud component (affecting only light from the birth clouds themselves, for example, stars younger than 10 Myr) and a separate diffuse component (affecting light from all sources). We assume a flexible star formation history (SFH) with a continuity prior⁵⁹, in which we fit for the ratio between the six SFH bins. Finally, we exclude the F410M and F356W filters from the fit as those contain the [Oiii]λλ5007,4959 doublet and thus may have been contaminated by AGN ionization. This can be seen in the SED in Extended Data Fig. 6, in which we see the host galaxy has increased flux in the F410M and F356W bands compared with the point source although these bands contain [Oiii]λλ5007,4959. This suggests that emission line flux from the AGN may probably still be contributing to the recovered host galaxy SED in these bands justifying our exclusion of these bands from the SED modelling.

Extended Data Fig. 7 shows the resulting Prospector fit (black) to the observed photometry (yellow) with the χ values below. These fits yield a stellar mass of \(\log ({M}_{* }/{M}_{\odot })=9.0{0}_{-0.25}^{+0.27}\) and \(\log ({M}_{* }/{M}_{\odot })=8.7{1}_{-0.52}^{+0.24}\) and instantaneous (within the past 10  Myr) \({\rm{SFR}}=1.4{8}_{-0.42}^{+0.95}\,{{\rm{M}}}_{\odot }/{\rm{yr}}\) and \({\rm{SFR}}=1.2{5}_{-0.89}^{+1.13}\,{M}_{\odot }\,{{\rm{yr}}}^{-1}\) for BAGPIPES and Prospector, respectively. These results are consistent within 1σ; thus, we combine the chains given by each code for the final estimate—\(\log ({M}_{* }/{M}_{\odot })=8.9{2}_{-0.31}^{+0.30}\) and \({\rm{SFR}}=1.3{8}_{-0.45}^{+0.92}\,{M}_{\odot }\,{{\rm{yr}}}^{-1}\).

We also test the effect of fitting the SED of the combined photometry, that is host galaxy + point source, to estimate how much we would overestimate the stellar mass of the galaxy by not decomposing the AGN and host galaxy components. We find that the results are altered by less than 1σ.

To further test the validity of our decompositions, we perform morphological analysis on the stacked NIRCam images of the source, comparing them with a model PSF. This comparison is carried out in two NIRCam bands—the F115W and F277W. The former was chosen because its PSF is the smallest among the bands not contaminated by the Lyman break, which prevents the use of F090W. The F277W band was chosen as it smears the source over more pixels, better sampling the PSF, whereas our ForcePho fits indicate that the PSF component is still sub-dominant even in this filter.

When analysing the F277W band, we perform an isophotal fit for both the source and the model PSF using the Photutils package⁶⁰. The results of this fit are shown in the top row of Supplementary Fig. 4. As shown in this figure, the best-fit ellipsoidal isophotes for GN-1001830 are slightly elongated, suggesting some extended morphology, whereas the isophotes of the PSF are circular. The radial profiles (Supplementary Fig. 4, rightmost column) show that our object is more extended than the PSF.

Radial profiles in the F115W band were estimated by placing concentric circular apertures on both our source and the model PSF as the flux was dispersed over too few pixels, undersampling the PSF and making isophotal fits non-viable. Nevertheless, as shown in Supplementary Fig. 4 (bottom right), our source, although compact, is still significantly more extended than the PSF.

All numeric properties of the host galaxy and its black hole obtained from the spectral and SED fitting are summarized in Extended Data Table 2.

Supplementary Fig. 5 shows the SFR of the host galaxy against lookback time (and redshift) presenting the SFH as derived by both Prospector and Bagpipes. Both SFHs are consistent with a flat SFH although the two codes show systematic offsets, probably resulting from different assumptions in the codes, in particular, the different stellar populations used. This suggests that the galaxy experienced almost constant star formation over the past few 100 Myr. A comparison of GN-1001830 with the star-forming main sequence at similar redshifts is shown in Extended Data Fig. 4. Our source lies below the star-forming locus of similarly massive galaxies by a factor of about3 and would take about 1 Gyr (that is, about the age of the Universe at z ~ 6) to double its mass with the current SFR. This indicates that the galaxy is currently fairly quiescent and may have been so for quite some time, although uncertainties on the SFH are large. The presence of AGN activity might offer a possible explanation for this state of the host, suggesting that AGN negative feedback might be responsible for suppressing star formation.

We also fit the extracted point source (AGN) component with a reddened power law with a fixed slope of β = −1.55 (F_λ ∝ λ^β), corresponding to an average slope of type 1, unobscured high-z quasars⁶¹. Assuming the SMC extinction curve, the resulting A_V is 2.68 ± 1.00, which, although poorly constrained, is fully consistent with the value derived from spectroscopy. The bolometric luminosity inferred from the fitted λL_λ at 5,100 Å, with the bolometric correction from ref. ⁶², and corrected for absorption, is \({4.6}_{-3.2}^{+7.5}\times 1{0}^{44}\,{\rm{erg}}\,{{\rm{s}}}^{-1}\), which is highly uncertain but consistent with the value obtained from the broad component of Hα.

Estimation of the black hole mass and accretion rate

As mentioned in the main text, the broad component of Hα can be used to infer the black hole mass by assuming the local scaling relations are valid at high redshift and, specifically, from the equation^17,18:

The best fitting values for the broad component of Hα (FWHM of \({5,700}_{-1,100}^{+1,700}\,{\rm{km}}\,{{\rm{s}}}^{-1}\) and flux of \({27.3}_{-4.0}^{+4.1}\times 1{0}^{-19}\,{\rm{erg}}\,{{\rm{s}}}^{-1}\,{{\rm{cm}}}^{-2}\)) give \(\log {M}_{{\rm{BH}}}/{M}_{\odot }\,=\) \(8.2{3}_{-0.36}^{+0.38}\), with the scatter on equation (2) included in the uncertainties. Coupled with the bolometric luminosity of 2 × 10⁴⁴ erg s⁻¹, estimated from the broad component of Hα (following the scaling relation given by ref. ⁶³), gives an Eddington ratio \({\lambda }_{{\rm{Edd}}}=0.00{9}_{-0.003}^{+0.005}\), with a systematic scatter of 0.5 dex.

We note that a recent measurement of the black hole mass in a super-Eddington accreting quasar at z ~ 2 has cast doubts on the validity of the ultraviolet (UV) virial relations for single-epoch black hole measurements³⁰. The same work, however, points out that when using Hα to measure the black hole mass, the discrepancy is only a factor of 2.5. Moreover, the discrepancy has been ascribed to deviations in the BLR size in the super-Eddington regime, which is certainly not the case for GN-1001830. In summary, the black hole mass measurement inferred from the broad Hα in GN-1001830 is reasonably solid.

It is difficult to estimate dust attenuation in this object. The lack of broad Hβ does not provide strong constraints: even assuming a standard case B recombination ratio of 2.8, the broad Hα line flux implies that the broad Hβ is not detectable. We, therefore, assume, as found for other AGN at high z, that the bulk of the obscuration towards the BLR also affects the narrow components⁶⁴. From the observed ratio of the narrow components of Hα and Hβ (\(5.5{1}_{-0.69}^{+0.86}\)), and assuming an SMC extinction law (appropriate for high-z AGN^65,66), we infer \({A}_{V}=2.0{0}_{-0.41}^{+0.44}\,{\rm{mag}}\). We also repeat the estimate using ratios of Hα and Hβ lines to Hγ. This yields \({A}_{V}=2.3{1}_{-0.61}^{+0.81}\,{\rm{mag}}\) and \({2.6}_{-2.1}^{+2.6}\,{\rm{mag}}\), respectively. These values are quite uncertain because of the lower brightness of Hγ and its proximity to the Hβ line but remain consistent with the previous estimate. The extinction-corrected black hole mass, bolometric luminosity and Eddington ratio are \(\log \,{M}_{{\rm{BH}}}/{M}_{\odot }=8.6{1}_{-0.37}^{+0.38}\) (with the uncertainty including intrinsic scatter on the virial relation), L_bol = 10⁴⁵ erg s⁻¹ and \({\lambda }_{{\rm{Edd}}}=0.02{4}_{-0.008}^{+0.011}\), respectively, with the same intrinsic scatter of 0.5 dex. The extinction correction is uncertain because of the use of the narrow lines. However, the extinction-corrected values of black hole mass and Eddington ratio are still consistent with the uncorrected ones within 2σ. Hence, our conclusions are not significantly altered by the presence of dust. Furthermore, in the next section, we show that fitting the nuclear source detected by NIRCam with a dust-reddened AGN slope results in an extinction consistent with that inferred from the narrow lines. In the previous section, we also infer the bolometric luminosity from the fitting of the nuclear SED and, although with large uncertainties, independently obtain a value consistent with that obtained from the broad component of Hα. As an additional check, we infer the bolometric AGN luminosity from the luminosity of the narrow Hβ and [Oiii]λ5007 lines using the scaling relations in ref. ⁶⁷. We obtain \({L}_{{\rm{bol}},{\rm{H}}\beta ({\rm{N}})}={4.7}_{-0.7}^{+0.7}\times 1{0}^{44}\,{\rm{erg}}\,{{\rm{s}}}^{-1}\) and \({L}_{{\rm{bol}},[{\rm{OIII}}]}={4.9}_{-0.1}^{+0.1}\times 1{0}^{44}\,{\rm{erg}}\,{{\rm{s}}}^{-1}\), which, although lower than the broad Hα estimate, are still consistent with it once the 0.3–0.4 dex scatter on the calibrations is taken into account. Moreover, we use the 5,100 Å luminosity from the previous section to independently infer the black hole mass, which results in \(\log \,{M}_{{\rm{BH}}}/{M}_{\odot }=8.1\pm 0.8\), with the large error coming from the uncertainties in source decomposition (being the AGN light sub-dominant), reddening and intrinsic scatter on the virial relations. This value is consistent with estimates using broad Hα, but the uncertainty makes it rather unconstraining.

We note that even without correcting for extinction, the resulting lower limit on the black hole mass would still imply a black hole to stellar mass ratio several 100 times above the local relation.

Black hole scaling relations with σ and dynamical mass

As discussed in ref. ², although high-z AGNs are offset on the black hole–stellar mass plane relative to the local relation, they are much closer to the local relation between black hole mass and stellar velocity dispersion σ_* relation and to the local relation between black hole mass and host galaxy dynamical mass. Here we explore the location of GN-1001830 on the latter two scaling relations.

We use the width of the [Oiii] line as a proxy of the velocity dispersion, FWHM = 255 ± 38 km s⁻¹, as measured from the medium-resolution grating spectrum, deconvolved for the LSF for compact sources. The uncertainty in the FWHM value for the narrow lines includes the 20% systematic uncertainty in the LSF broadening, which is not significant for the broad component. We then derive the stellar velocity dispersion by correcting gaseous velocity dispersion by 0.1 dex to obtain the stellar velocity dispersion, following ref. ⁶⁸, giving \({\sigma }_{* }=12{1}_{-16}^{+16}\,{\rm{km}}\,{{\rm{s}}}^{-1}\), with 0.3 dex intrinsic scatter for \(\log {\sigma }_{* }=2.08\pm 0.32\). The resulting location of GN-1001830 in Extended Data Fig. 8 (magenta circle) shows that our object is close to the local relation (solid line, with dispersion shown with a shaded region), as are other high-z AGNs previously measured by JWST.

We also use the host galaxy parameters to estimate its dynamical mass using the same approach as in refs. ^2,4, which makes use of the equation

where K(n) = 8.87 − 0.831n + 0.0241n², with Sérsic index n, \(K(q)\,=\) \({[0.87+0.38{{\rm{e}}}^{-3.71(1-q)}]}^{2}\), where q is the axis ratio⁶⁹, R_e is the estimated effective radius and σ_* is the stellar velocity dispersion. Using equation (3) along with values for R_e, n and q given by ForcePho (summarized in Extended Data Table 2) gives \(\log {M}_{{\rm{dyn}}}/{M}_{\odot }=9.5{0}_{-0.39}^{+0.39}\), with the uncertainty dominated by the intrinsic scatter on equation (3). We warn that the errors on this estimate are probably underestimated because of the absence of high-resolution spectral observations for GN-1001830. Moreover, the M_dyn value may be underestimated as our object is not fully centred in the slit (Supplementary Fig. 1), which may cut off part of the rotation curve in case of source rotation. The position of our object with respect to other JWST sources and the local scaling relation on the M_BH–M_dyn and M_BH–\(\log {\sigma }_{* }\) star plane is shown in Extended Data Fig. 8. The source remains above the local relation; however, the difference is not as severe as in the black hole–stellar mass relation.

Gas fraction and SFE

Following our dynamical and stellar mass estimates, we can obtain an estimate for the gas mass in the host galaxy, using M_gas = M_dyn − M_*, assuming little contribution from dark matter, especially at such early epochs, within the central few 100 pc. This gives \(\log {M}_{{\rm{gas}}}/{M}_{\odot }=9.3{7}_{-0.39}^{+0.39}\) and a gas fraction f_gas of \(0.7{4}_{-0.30}^{+0.18}\). The depletion time for our object can thus be estimated as \(\frac{{M}_{{\rm{gas}}}}{{\rm{SFR}}}=1.5{9}_{-0.85}^{+0.91}\,{\rm{Gyr}}\). The depletion time given by the scaling relation from ref. ⁷⁰, assuming a star-forming main sequence of the form given in ref. ⁷¹, evaluates to 0.66 Gyr. Although these values are consistent within the intrinsic scatter on the relation, as shown in Supplementary Fig. 6, this is still suggestive of star formation being inhibited relative to the population of normal star-forming galaxies at this epoch.

Completeness analysis

To quantify the selection bias affecting the identification of AGN with broad emission lines, we run the fitting procedure described in 1.2 on simulated broad Hα profiles. These profiles were simulated by using the error extension of the prism data for GN-1001830 to simulate Gaussian noise and adding Gaussian line profiles along with a power-law continuum on top of it. The FWHM and luminosity for the narrow Hα component were uniformly sampled from FWHM ∈ [200, 600] km s⁻¹ and \(\log [L({\rm{erg}}\,{{\rm{s}}}^{-1})]\in [42,43]\) ranges, respectively. Redshifts of the simulated sources were uniformly drawn from between 6 and 7. Continuum normalization was set to be 100 times smaller than the Hα narrow peak height and slopes were sampled from a uniform distribution (−1 < α < 1). The model grid for the broad component of Hα was computed by varying log λ_Edd between −3.0 and 0.5 in steps of 0.25 and the \(\log {M}_{{\rm{BH}}}/{M}_{\odot }\) between 5.5 and 9 in steps of 0.25. For bins with \(\log {M}_{{\rm{BH}}}/{M}_{\odot } < 7\), we performed the simulation assuming the grating spectrum R1000 error extension to assuage the effects of the lower resolution of the prism. Equation (2) and scaling relations from ref. ⁶³ were used to convert the Eddington ratio and black hole mass values to the luminosity and width of the simulated broad component. This yielded a 15 × 15 grid, each point of which contained 100 spectra simulated according to the above recipe.

The fitting of the simulated data was carried out with the same parameters as those described in 1.2, with a bounded least-squares procedure being used to make the fitting computationally tractable. The completeness of each grid point was calculated as the ratio of the sources recovered according to our criteria (in terms of both broad-line significance and ΔBIC) to the sources inserted. The final completeness function is presented in Fig. 3 and shows that we are inherently biased against highly sub-Eddington black holes predicted by simulations, as most of those lie in regions of low completeness. We caution that the procedure described here provides only completeness with respect to the final step of AGN selection—the fitting of the BLR—and covers only the medium-depth tier of the JADES survey. The full selection function for the JADES survey is more complex as it spans multiple survey tiers, telescope instruments and source selection methods. Its full treatment is thus beyond the scope of the paper, and our results here should be regarded as closer to an upper limit.

The role of selection effects on the M
_BH–M
_* relation

As discussed in the main text, the finding that most of the black holes at high-z newly discovered by JWST are overmassive on the M_BH–M_* could be because of a selection effect, that is, the scatter of the relation is much larger at high redshift and more massive black holes tend to be preferentially selected at high redshift because they are, on average, more luminous. We have shown and confirmed that selection effects play an important part; however, our findings indicate that they cannot completely explain the offset relative to the local M_BH–M_* relation.

Specifically, the scenario proposed in ref. ²⁹ (which predicts a black hole–stellar relation similar to the local one but with an order of magnitude scatter) hardly reaches the black hole to stellar mass ratio observed in our object, and their observational bias scenario would require an observed luminosity of about 10⁴⁵⁻⁴⁶ erg s⁻¹, that is, 1–2 dex above the luminosity observed in GN-1001830 before dust obscuration correction (we note that they do not take into account dust extinction, whereas this and most JWST-discovered AGN are affected by extinction).

The Trinity simulation²⁸ can produce overmassive black holes as observed in GN-1001830, but they require even higher luminosities, in excess of 10⁴⁸ erg s⁻¹, hence totally incompatible with the luminosity of our object. Therefore, although selection effects are important (as we illustrate more thoroughly in section ‘Completeness analysis’), our finding suggests that the overmassive nature of high-z black holes is also associated with an intrinsic offset of the black hole–stellar mass relation, as also suggested by other studies¹³.

Finally, our finding that GN-1001830, as well as many other JWST-discovered AGN², is closer (or even consistent for many AGN) with the M_BH–σ and M_BH–M_* relations, indicates that the selection effects on black hole mass cannot play a major part, or else the same strong, orders of magnitude, offset should also be present on these relations.

These aspects are, however, outside the scope of this paper and will be discussed more extensively in a dedicated paper.

Comparison with the FABLE simulations

To further contextualize our findings, we also compare the properties of our source with predictions from the FABLE (Feedback Acting on Baryons in Large-scale Environments)⁷² simulations. These are carried out with the massively parallel AREPO code⁷³ with new comoving 100 h⁻¹ Mpc boxes. We consider both Eddington-limited and super-Eddington accretion, bounding our black hole Bondi–Hoyle–Lyttleton accretion rate in two different simulations by 1 and 10× Eddington. Details of other subgrid models, which are largely based on the Illustris galaxy formation models⁷⁴, can be found in ref. ⁷².

A comparison between the simulated sources from FABLE and our object is provided in Extended Data Fig. 9, with the same coding as in Fig. 2, that is, grey symbols show the Eddington-limited scenarios and the red contours the distribution of simulations in which super-Eddington accretion is allowed. As shown in this figure, the super-Eddington simulation can produce more massive black holes with respect to both luminosity and stellar mass in comparison with the Eddington-limited simulation, resulting in a better match with our observations. However, GN-1001830 still lies above our super-Eddington simulations in the right-hand M_BH–M_* plot, probably because of the FABLE simulations lacking in volume. To better sample a larger volume, we also include simulation data from magnified simulations of a massive protocluster (M_h > 10¹²M_⊙ at z = 6.66), previously used in ref. ¹⁴, which is taken from the larger Millennium simulation volume of 500 h⁻¹ Mpc (ref. ⁷⁵). These points occupy the high-mass region in both panels of Extended Data Fig. 9. The objects obtained in these simulations have stellar masses slightly larger than that inferred from GN-1001830 and, in runs including super-Eddington accretion up to 10× Eddington (and also earlier black hole seeding, unlike the other FABLE results shown here; ref. ¹⁴), have more massive black holes. Although a more quantitative match will be explored in future work, these results qualitatively show that super-Eddington bursts more readily explain the properties and presence of objects such as GN-1001830, particularly in increasing the black hole to stellar mass ratio, relative to Eddington-limited scenarios.

We also note that, as with CAT sources discussed in the main text, most of the simulated highly sub-Eddington sources reside in the low completeness region of the JADES survey (Extended Data Fig. 10) and are thus hard to detect at these redshifts even with current instruments.