# A mature quasar at cosmic dawn revealed by JWST rest-frame infrared spectroscopy

Sarah E. I. Bosman, Javier Álvarez-Márquez, Luis Colina, Fabian Walter, Almudena Alonso-Herrero, Martin J. Ward, Göran Å-Stlin, Thomas R. Greve, Gillian Wright, Arjan Bik, et al.

► To cite this version:

Sarah E. I. Bosman, Javier Álvarez-Márquez, Luis Colina, Fabian Walter, Almudena Alonso-Herrero, et al.. A mature quasar at cosmic dawn revealed by JWST rest-frame infrared spectroscopy. Nature Astronomy, 2024, 8, pp.1054-1065. 10.1038/s41550-024-02273-0. insu-04726421

HAL Id: insu-04726421

https://insu.hal.science/insu-04726421v1

Submitted on 8 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A mature quasar at cosmic dawn revealed by JWST rest-frame infrared spectroscopy

Received: 27 July 2023

Accepted: 15 April 2024

Published online: 17 June 2024

<figure><img src="image_2.png" /></figure>

Check for updates

Sarah E. I. Bosman*<sup>1,2</sup>✉, Javier Álvarez-Márquez<sup>3</sup>, Luis Colina<sup>3</sup>, Fabian Walter<sup>2</sup>, Almudena Alonso-Herrero<sup>4</sup>, Martin J. Ward<sup>5</sup>, Göran Östlin<sup>6</sup>, Thomas R. Greve<sup>7,8,9</sup>, Gillian Wright<sup>10</sup>, Arjan Bik<sup>6</sup>, Leindert Boogaard<sup>2</sup>, Karina Caputi<sup>11</sup>, Luca Costantin<sup>3</sup>, Andreas Eckart<sup>12</sup>, Macarena García-Marín<sup>13</sup>, Steven Gillman<sup>7,8</sup>, Jens Hjorth<sup>14</sup>, Edoardo Iani<sup>11</sup>, Olivier Ilbert<sup>15</sup>, Iris Jermann<sup>7,8</sup>, Alvaro Labiano<sup>16</sup>, Danial Langeroodi<sup>14</sup>, Florian Peißker<sup>12</sup>, Pierluigi Rinaldi<sup>11</sup>, Martin Topinka<sup>17</sup>, Paul van der Werf<sup>18</sup>, Manuel Güdel<sup>2,19,20</sup>, Thomas Henning<sup>2</sup>, Pierre-Olivier Lagage<sup>21</sup>, Tom P. Ray<sup>17</sup>, Ewine F. van Dishoeck<sup>18</sup> & Bart Vandenbussche<sup>22</sup>

The rapid assembly of the first supermassive black holes is an enduring mystery. Until now, it was not known whether quasar ‘feeding’ structures (the ‘hot torus’) could assemble as fast as the smaller-scale quasar structures. We present JWST/MRS (rest-frame infrared) spectroscopic observations of the quasar J1120+0641 at z = 7.0848 (well within the epoch of reionization). The hot torus dust was clearly detected at $\lambda_{rest} \simeq 1.3$ µm, with a black-body temperature of $T_{dust} = 1,413\,5_{-7.4}^{+5.7}$ K, slightly elevated compared to similarly luminous quasars at lower redshifts. Importantly, the supermassive black hole mass of J1120+0641 based on the Hα line (accessible only with JWST), $M_{BH} = 1.52 \pm 0.17 \times 10^9 M_\odot$, is in good agreement with previous ground-based rest-frame ultraviolet Mg II measurements. Comparing the ratios of the Hα, Paα and Paβ emission lines to predictions from a simple one-phase Cloudy model, we find that they are consistent with originating from a common broad-line region with physical parameters that are consistent with lower-redshift quasars. Together, this implies that J1120+0641’s accretion structures must have assembled very quickly, as they appear fully ‘mature’ less than 760 Myr after the Big Bang.

The origin of the first supermassive black holes (SMBHs) is an enduring puzzle. Measurements of the masses of SMBHs based on their rest-frame ultraviolet spectra have revealed the presence of a population of black holes (BHs) with masses up to $M = 1.5 \times 10^9 M_\odot$ at cosmological redshift $z > 7$, less than 800 Myr after the Big Bang<sup>1-6</sup>. The existence of such massive SMBHs requires either very massive BH seeds of $\ge 10^4 M_\odot$ at $z = 30$ or continuous growth rates in excess of the Eddington limit by a factor of a few<sup>7,8</sup>; (see ref. 9 for a review). Both scenarios pose substantial problems for current models of SMBH origins. The first accreting SMBHs, that is quasars, have, therefore, been

the target of extensive observational campaigns looking for clues to their growth mechanisms.

Observations of $z > 7$ quasars with optical and near-infrared instruments (corresponding to the quasars’ rest-frame ultraviolet) have revealed a complex picture. On one hand, early quasars appear to be accreting matter at below or only slightly above the Eddington rate<sup>5,10,11</sup>. The rest-frame ultraviolet broad emission lines of quasars at $5.7 < z < 7.5$ are identical in their luminosities and flux ratios to $2 < z < 5$ quasars, indicating unchanging levels of metal enrichment of quasar broad-line regions (BLRs) across 2.6 billion years of cosmic time<sup>12</sup>.
<figure><img src="image_3.png" /><figcaption>**Fig. 1 | JWST spectrum of the quasar.** Top, MRS spectrum of J1120+0641 at z = 7.0848. The detected emission lines are labelled. The last subarray of channel 4 is omitted, and the rest of the spectrum has been binned across four bins for display. The grey line indicates the 1σ uncertainty ΔF, calculated from 12</figcaption></figure>

randomly placed background apertures at each pixel. Coloured lines show the power-law (green) and black-body (blue) continuum emission, emission lines (red) and total model (orange). Bottom, residuals of the spectral fit in standard deviations, $F(\lambda) - \text{Fit}(\lambda)/\Delta F(\lambda)$.

On the other hand, quasar outflow properties do show signs of evolution at fixed quasar luminosity and mass. The highly ionized C IV emission line displays outflows of ≤1,000 km s⁻¹ on average in 1 < z < 5.5 quasars, which more than double to ~2,000–2,500 km s⁻¹ at z > 6.5 (refs. 10,13–16). Concurrently, the occurrence rate of quasars with broad absorption line outflow features (broad absorption line quasars) among the bright quasar population more than doubles at z > 5.8 compared to 2 < z < 4 (refs. 5,17,18).

These indications of changing quasar outflows have not yet led to much progress in understanding early SMBH growth mechanisms. If SMBHs grow from massive initial seeds by accretion processes identical to late-time quasars, then presumably, no changes at all should be seen in the earliest quasars. The main alternative to massive seeds is super-Eddington accretion in the earliest SMBHs, but changes in outflows are not a prediction of such models. In addition, no z > 6 quasars were observed to be sufficiently super-Eddington. To circumvent this issue, models have often invoked highly obscured (episodic) super-Eddington growth (for example, refs. 19–21). Episodic obscuration models reconcile the short apparent lifetimes of z > 6 quasars, some of which appear to be observed less than 0.1 Myr after the onset of the continuing luminous phase²²⁻²⁴.

Such short obscuration-shifting timescales might (hopefully) leave imprints in the rest-frame infrared properties of the first quasars. In particular, remnants of obscuration could be visible (1) in the hot dusty torus, which is most tightly linked to accretion processes, (2) in the obscuration of quasar BLR lines, which, if present, would affect the rest-frame ultraviolet measurements of the black hole mass ($M_{BH}$) employed to assess Eddington rates, or (3) in the density and ionization state of the BLR, as the mode of accretion shifts quickly from super-Eddington to sub-Eddington. The rest-frame infrared is, therefore, an ideal ground to search for signs of unusual quasar growth.

In this context, we used the recently commissioned mid-infrared instrument (MIRI)²⁵⁻²⁷ onboard the James Webb Space Telescope (JWST) to observe the rest-frame infrared spectrum of an early quasar.

The bright quasar J1120+0641 was the first to be discovered at z > 7 (ref. 1). It has, therefore, been the target of extensive multi-wavelength observations, which justify its inclusion as a priority target for guaranteed time observations in cycle 1. At the time of writing, nine quasars are known at z > 7 (for example, refs. 6,28,29; see ref. 30 for an up-to-date list).

Observations of J1120+0641 with NOEMA (ref. 31) and later ALMA (ref. 32) detected the submillimetric continuum and [C II] fine-structure line emission from its host galaxy at a systemic redshift of z = 7.0848 ± 0.0004 (refs. 33,34). A deep spectrum of J1120+0641 was obtained with the X-Shooter instrument³⁵ in the optical and at near-infrared wavelengths³⁶. This spectrum revealed a fast outflowing C IV emission line with a velocity $v = 2,136_{-27}^{+32}$ km s⁻¹ with respect to systemic¹⁶, a value that is significantly elevated compared to z < 5 quasars of comparable brightness³⁷ but average among z > 6 quasars.

The BH mass of J1120+0641 is 1.35 ± 0.04 × 10⁹ $M_{\odot}$ based on the Mg II line⁵, and its Eddington ratio is ~0.40 (ref. 11). J1120+0641 is radio-quiet³⁸,³⁹, and its X-ray properties are indistinguishable from other bright quasars at 1 < z < 5 with similar Eddington rates⁴⁰⁻⁴².

JWST observations⁴³ of J1120+0641 were obtained on 12 January 2023 (proposal ID 1263, P.I. Luis Colina) using the mid-infrared resolution spectrograph (MRS; refs. 44,45) of MIRI. The spectrograph provides integral field spectroscopy covering the 4.9 to 27.9 µm range with intermediate spectral resolution⁴⁴,⁴⁶,⁴⁷. The observations using the three MRS configurations (short, medium and long) covered the entire MRS spectral range with an integration time of 3,153 s per configuration.

We employed a Planck cosmology throughout the analysis⁴⁸. All quoted rest-frame wavelengths correspond to their values in a vacuum.

## Results

The spectrum of J1120+0641 is shown in Fig. 1. The reduction process and spectral modelling are described in Methods. We found that the residuals of our total spectral fit to the continuum emission and emission lines are adequately described by a Gaussian distribution over all channels, with no signs of additional broad emission lines or slowly varying deviations.
<figure><img src="image_4.png" /><figcaption>**Fig. 2 | Hot torus dust temperature across redshift.** Measurements of hot torus dust temperature across redshift. Black triangles indicate upper limits from ref. 63. Our spectroscopic dust temperature measurement of J1120+0641 has substantially smaller uncertainties. The occurrence rate of quasars with $T_{\text{dust}} > 1,400$ K might be increasing towards higher redshifts, as 4/16 quasars at</figcaption></figure>

$z > 6$ have dust hotter than 1,400 K compared to 1/137 at $z < 6$. Bounds containing 68% and 95% of temperatures are derived from the $N = 84$ quasars at $z < 5$ from ref. 64 which have spectroscopic measurements. All measurements are shown with $1\sigma$ uncertainties.

## Continuum emission and torus dust temperature

Quasars have long been known to display a prominent bump in emission around rest-frame wavelengths $1 \le \lambda_{\text{rest}} \le 20$ µm, which has been associated with emission from dust grains$^{49–51}$. Based on reverberation-mapping studies, the size of the dusty structure is known to be subparsec and to scale with luminosity$^{52,53}$, leading to the picture of a ‘hot torus’ located near the dust sublimation radius, with temperatures $T > 1,000$ K. The dusty torus has been modelled using both smooth dust distributions (for example, refs. 54–56) as well as more complex clumpy dust models (for example, refs. 57,58). Despite the potential complexity of the hot torus structure, observations have consistently found that the hot dust bump can be adequately fitted with a single-temperature black-body form across all quasar types and luminosities (for example, refs. 59–61). Some suggestions have been made that quasars at $z \ge 6$ are more often hot-dust-free than their later-time counterparts$^{62}$, but these tentative results have not been confirmed$^{63}$. In contrast, the hot dust upturn is clearly visible in our observations of J1120+0641.

The slope of the accretion disk continuum of J1120+0641 is best fitted as $\alpha = -0.372_{-0.035}^{+0.031}$, which indicates that there is a typical ‘blue’ continuum with no sign of extinction of the accretion disk. The emission from the accretion disk becomes equal to the brightness of the hot torus dust around observed-frame wavelengths $\lambda_{\text{obs}} = 10.5$ µm (Fig. 1). A hot dust temperature of $T = 1,417.5_{-7.2}^{+5.6}$ K provides the best fit to the continuum under the assumption of black-body emission.

The black-body temperature is in relatively good agreement with torus dust observed for quasars at $5 < z < 6.5$ (Fig. 2). Reference$^{63}$ estimated the uncertainties in photometry-based measurements to be around ±50 K. In contrast, the spectroscopic measurement enabled by the MIRI-MRS is accurate to better than 10 K (0.5%). The temperature we measured is slightly elevated compared to past literature, with only three quasars out of the sample of 68 in ref. 63 having higher temperatures (4% of the sample). J1120+0641 and the three other quasars with $T > 1,400$ K are all found at $z > 6$, where they constitute 25% of observed quasars (4 of 16). In contrast, no quasars with $T > 1,400$ K are found at $5 < z < 6$ (0 of 53). This suggests a possible evolution in dust temperature. Measurements of the hot dust temperature of 84 quasars at $z < 5$ are provided by ref. 64 based on low-resolution spectroscopy ($R \approx 20$) from Spitzer and AKARI. Combining the samples, the temperature trend becomes even stronger: 4 out of 16 quasars studied at $z > 6$ have dust hotter than 1,400 K, and only 1/137 for $0 < z < 6$ are over that temperature, for similar quasar bolometric luminosities. Using Poisson occurrence statistics$^{65}$, we determined that the evolution of

the quasar dust temperature was statistically significant at the $2\sigma$ level but did not reach $3\sigma$.

A potential caveat is the lack of correction for contamination by broad emission lines in past studies. The Paschen $\alpha$ and Paschen $\beta$ lines, in particular, contribute flux to the same observed wavelengths as the hot dust ‘upturn’, which is used for spectral energy distribution fitting of the temperature in some past photometric studies (for example, ref. 63). The effect of contamination by emission lines would be a systematic underestimation of the dust temperature, which may explain the shift between J1120+0641 and some of the past literature. We plan to quantify this bias and explore further trends of the dust temperature with, for example BH mass, in future work. More spectroscopic observations of hot dust in $z > 6$ quasars, but also in $z < 6$ quasars, are also clearly needed to confirm this potential temperature trend.

The potential cause(s) of an evolution towards hotter torii at earlier times is unclear. J1120+0641’s Eddington rate, luminosity, BH mass and X-ray properties are normal compared to the lower-redshift comparison samples. We speculate that a change in dust composition could create such an effect in principle. Both graphite and silicate grains have been suggested to make up quasar tori, since they can withstand the high measured temperatures without sublimating. Small dust grains could in principle lead to a higher torus density and, therefore, a higher temperature at the innermost radius, but small grains also sublimate more readily. Further modelling is required to quantify the potential impact of dust make-up on measured torus temperature.

## Emission-line BH mass estimates

The hydrogen Paschen-series emission lines are another crucial feature in the rest-frame infrared spectra of quasars$^{66}$. The Pa$\alpha$ and Pa$\beta$ lines are far less affected by extinction than hydrogen lines at shorter wavelengths, which makes them more reliable tracers of BH mass in highly extincted quasars$^{67}$. The hydrogen Balmer and Paschen lines, in contrast with other broad lines in the rest-frame ultraviolet, optical and infrared, are non-resonant and give access to a common BLR through several transitions. As a consequence, the flux ratios of the lines are good tracers of potential dust extinction in the BLR, as well the strength of ionization and density of the hydrogen BLR (refs. 67,68). The relevant BH mass scaling relations are presented and discussed in Methods.

The optimal fit to the H$\alpha$ line is shown in Fig. 3. The broadest component, with a full-width at half-maximum (FWHM) of $11,720_{-710}^{+560}$ km s$^{-1}$, shows signs of a mild outflow with speed $360 \pm 70$ km s$^{-1}$. In contrast, the narrower broad component (BC1) shows an inflow with respect to the quasar’s systemic redshift of $-315 \pm 37$ km s$^{-1}$.
<figure><img src="image_5.png" /><figcaption>Fig. 3 | Broad emission lines of the quasar. Spectra of the broad lines of J1120+0641 (flux in black and uncertainty in grey) with the optimal line model fit (orange), split into the continuum constituting of a power law and a black</figcaption></figure>

body (green) and the broad emission lines (red). The top right panel shows the emission-line complex of Paβ, He I and Pay in order of increasing wavelength. The latter two lines are blended.

We obtain a BH mass based on the FWHM and luminosity of the Hα line of:

$$M_{BH,H\alpha} = 1.55 \pm 0.22 \times 10^9 M_{\odot}. \qquad (1)$$

For comparison, the BH mass of J1120+0641 derived based on the rest-frame ultraviolet Mg II emission line is⁵

$$M_{BH,Mg II} = 1.35 \pm 0.04 \times 10^9 M_{\odot}. \qquad (2)$$

The BH mass based on the rest-frame ultraviolet C IV emission line was calculated in ref. 11 using the scaling relations of ref. 69 with the corrections of ref. 70:

$$M_{BH,C IV} = 2.40^{+0.06}_{-0.05} \times 10^9 M_{\odot}. \qquad (3)$$

These measurements are remarkably compatible, especially considering the scaling relations on which they rest contain intrinsic scatter at the level of 0.3–0.5 dex (that is a factor of 2–3). A similar level agreement is seen for two other quasars at z > 6 for which the Hα line was observed with NIRSpec (refs. 71,72).

A different BH mass estimate can be made based on the C IV emission line using the scaling relations of ref. 69 but without the corrections of ref. 70:

$$M_{BH,C IV, no corr} = 4.79^{+0.12}_{-0.11} \times 10^9 M_{\odot}. \qquad (4)$$

However, the corrections of ref. 70 are most relevant for quasars presenting an asymmetric C IV emission line, of which J1120+0641 is a good example.

The Eddington rate from the Hα BH mass is $f = 0.9 \pm 0.1$ with systematic uncertainties from the BH mass scaling relations of a factor of 2. Our measurement of J1120+0641’s Eddington rate is higher but within scatter compared to previous estimates¹¹: measure $f = 0.57 \pm 0.01$ based on the C IV emission line, with systematic uncertainties on the C IV BH mass of ~0.3–0.5 dex.

The Paα broad line is described well by a two-Gaussian model (Fig. 3). However, some low-frequency spectroscopic fringing is present around the location of the line, which necessitated fixing the centre of the broadest component (BC2) to be at the same redshift as the broadest component of the Hα line to achieve convergence.

With this constraint in place, the best-fitting redshift of the narrower broad component of the Paα line (BC1) implies an inflow of

~415 ± 150 km s⁻¹, a value consistent with the core component of the Hα line within uncertainties. It is, therefore, highly probable that (the core part of) the two emission lines originate from the same part of the BLR, even though the FWHM of the core components is not fully consistent.

The FWHM of the broadest component of Paα is 11,900±2100 km s⁻¹, fully compatible with the FWHM of the broadest component of Hα. The flux ratio of the core components (BC1) to the broadest components (BC2) is ~0.8 and ~0.6 for Hα and Paα, respectively. However, the fluxes in the BC1 and BC2 components are degenerate in the Paα line, such that the small difference is unlikely to be statistically significant. The measurement of the total flux of the Paα emission line is not affected by the degeneracy. We obtained a BH mass using the relation from ref. 67:

$$M_{BH,Paα} = 9.8^{+2.8}_{-2.3} \times 10^8 M_{\odot}. \qquad (5)$$

The resulting mass is consistent with the estimates based on the Hα and Mg II emission lines within the intrinsic scatter of the respective scaling relations. The mass is, however, inconsistent with the estimate based on C IV.

The Paβ emission line was more challenging to fit with two Gaussian components. The core component (BC1) has a relatively high signal-to-noise ratio (SNR), and its FWHM and redshift are consistent with the properties of the core components of Hα and Paα. This again suggests an origin from a common BLR. However, we encountered strong difficulties in constraining Paβ’s broader component, if it is detected at all. As can be seen in Fig. 3, the spectrum surrounding the Paβ emission line contains many fringing residuals giving rise to spurious and fairly narrow features. Despite masking a feature that otherwise causes the fit to fail to converge, the constraints we obtain are still rather poor. The optimal FWHM of BC2, 11,040±2,200 km s⁻¹, reflects this uncertainty while still being fully consistent with the widths of Hα’s and Paα’s broadest components. Moreover, the shift of the broadest component is almost certainly impacted by the presence of spurious spectral features and cannot be used to draw conclusions on the BLR kinematics. In conclusion, the lower SNR of the line or the presence of spectroscopic fringing do not enable us to confidently characterize the line’s broadest component. The line’s total flux is relatively little impacted by these considerations as the underlying continuum is very well constrained. We obtained the BH mass estimate based on Paβ using the scaling relation of ref. 73:

$$M_{BH,Paβ} = 8.7^{+2.1}_{-1.7} \times 10^8 M_{\odot}. \qquad (6)$$
Again, the mass estimate based on Paβ is consistent with all other measurements obtained so far within the relation's intrinsic scatter, with the exception of the estimate based on the C IV line. Our observations, therefore, confirm suspicions that C IV is one of the least reliable tracers of BH mass due to its highly ionized and scattered nature (for example, ref. 70).

## Other emission lines

The He I and Pay lines are blended, which strongly hinders the determination of their fluxes and kinematics. We, therefore, modelled them with one Gaussian component each, nevertheless finding many of their fitting parameters to be degenerate. Based on our best-fitting decomposition, the He I lines appear to contain about twice the flux of the Pay line. Both lines are quite broad, with FWHMs of 4,050±310 km s⁻¹ and 3,520±390 km s⁻¹, respectively. That Pay's FWHM exceeds the width of the core components of the previous three hydrogen transitions could be due to the presence of an additional broader component, which we do not have sufficient SNR to separate. Alternatively, the blended He I line could also contain broad wings that contaminate the fit over Pay.

Paδ, in contrast, was not contaminated by blending, but its lower SNR of ≤2 per pixel warrants only a single Gaussian component (leftmost line in Fig. 3, top right panel). The width of the line is captured well and fully consistent with the widths of the core components of Hα, Paα and Paβ. This further suggests that the narrower components of all detected hydrogen transitions originate from the same BLR and that the anomalously large measurement in Pay is due to blending contamination.

As the Paδ line provides a well-constrained FWHM without contamination by fringing, we opted to use it in addition to Paα and Paβ in our nominal Paschen-series-based BH mass estimate with the infrared-only estimate method of ref. 61 (Methods). We first calculated the continuum luminosity at 1 µm. Using only the continuum model (accretion disk power law and hot dust black body), we obtained $V L_{1µm} = (1.39 \pm 0.25) \times 10^{46}$ erg s⁻¹. To obtain the optimal combined FWHM of the hydrogen emission lines, FWHM<sub>H</sub>, we combined the total line widths of Paα, Paβ and Paδ with inverse-variance weighting and obtained FWHM<sub>H</sub> = 2,590 ± 100 km s⁻¹. Finally, we obtained:

$$M_{BH,IR} = 8.9 \pm 1.4 \times 10^8 M_{\odot}. \qquad (7)$$

This BH mass is again fully consistent with the previous ones with the exception of C IV. In summary, all other mass tracers seem to be functioning no differently than at later times. An overview of the BH mass measurements is shown in Fig. 4.

The O I emission line is the only non-hydrogen emission line detected that is not blended. We found the FWHM of O I, 2,900±410 km s⁻¹, to be consistent with that of the core components of the hydrogen lines, which is sensible given the ionization charge equilibrium between the O I ion and neutral hydrogen.

Combining several transitions of O I is a test of the ionization strength in the joint O I + Fe II BLR, using an approach like that we employed for the broad hydrogen lines⁷⁴⁻⁷⁶. However, we have only one other line available to us, the rest-frame ultraviolet O I 1,304 Å line³⁶. We calculated the line fluxes of both transitions taking into account the slit losses in the X-Shooter spectrum and obtained a ratio for these lines $f_{O_{1,0.8446µm}/f_{O_{1,1304Å}}} = 0.36 \pm 0.08$. Reference ⁷⁶ reports higher values for six quasars at $z < 1.5$, ranging from 1.5 to 4.7. This deviation would suggest a slightly lower ionization parameter U or a lower O I BLR density in J1120+0641 compared to the six quasars at $z < 1.5$. However, the O I 1,304 Å line is blended with a Si II line from which it cannot be separated in the X-Shooter spectrum. The comparison sample is too small to draw conclusions on the possible evolution of O I.

<figure><img src="image_6.png" /><figcaption>Fig. 4 | Measurements of the BH mass of J1120+0641 using different tracers.</figcaption></figure>

Measurements of the BH mass of J1120+0641 using different tracers (N = 1). The blue boxes represent 1σ measurement uncertainties. The longer error bars (whiskers) represent the 1σ intrinsic scatter in the BH mass scaling relations. Mg II and C IV-based masses are from ref. 5 and ref. 11, respectively, whereas the other measurements originate in this paper. IR, infrared.

## BLR density and ionization parameter

We modelled the hydrogen BLR using a single-slab approximation with Cloudy (Methods). The results are shown in Fig. 5, together with measurements from J1120+0641 and comparison samples at low redshift. Based on interpolating between outputs, we find that the ratios of Hα, Paα and Paβ are like $0 < z < 4$ luminous quasars, and the lines are modelled well with a BLR density of $log N_H/cm^{-3} \ge 12$ and an ionization parameter $-7 < log U < -4$. Very little extinction of the BLR, of the order $E(B-V) \le 0.1$, is required to explain the observed emission-line ratios. Compared to the sample of 12 quasars of ref. 77 at $0 < z < 0.2$ for which all three emission lines are detected, J1120+0641 seems to have a somewhat low ionization parameter U (although not the lowest) as indicated by a slightly high value of Paα/Paβ = 1.54 ± 0.34 but a normal value of Paα/Hα = 0.106 ± 0.010. Unlike many $z < 4$ quasars⁶⁶, J1120+0641's hydrogen broad emission lines are marginally compatible with pure case B recombination (at $2\sigma$), which predicts Paα/Paβ = 2 and Paα/Hα = 0.1 (ref. 78, cyan point in Fig. 5). Overall, the ratios of the Paschen lines in J1120+0641 are very much like those of quasars at $0 < z < 5$.

## Discussion

### Structure of the BLR

Our observations reveal the kinematics of the Hα and first four Paschen-series emission lines of hydrogen in unprecedented detail in $a = 7$ quasar. Perhaps surprisingly, we found that the hydrogen broad lines are fully consistent with a simple model of the hydrogen BLR. The 'core' components of the Hα and Paβ lines, as well as the single visible component of Paδ, are consistent in their widths of ~2,880 km s⁻¹, well within their ~5% uncertainties. The broadest component, when it is detected in Hα and Paα and marginally detected in Paβ, is similarly extremely consistent across lines in its FWHM of ~11,000 km s⁻¹. The velocity shift of the core components is similarly highly consistent with its value of ~320 km s⁻¹ (inflowing) across Hα, Paα and Paβ. The broadest component's shift is less well attested and can be measured with confidence only for Hα, being ~360 ± 70 km s⁻¹ (outflowing). Deviations can be accounted for by the presence of residual fringing (broadest components of Paα and Paβ), line blending (Pay) and low SNR (Paδ).

A slightly redshifted 'core' component of the broad hydrogen emission line is also seen in Hβ in a few quasars at $z > 6.5$ (ref. 79), matching our observations of an ~300 km s⁻¹ 'inflow'. The physical origin of this
spectral cube⁸⁷, was run using the background-subtracted and fully calibrated detector images and skipping the master background subtraction and the outlier detection steps. This process produced 12 three-dimensional spectral cubes, one for each band (short, medium and long) of the MRS channels, with spatial and spectral sampling of 0.13' × 0.13' × 0.8 nm for channel one, 0.17' × 0.17' × 1.3 nm for channel two, 0.20' × 0.20' × 2.5 nm for channel three and 0.35' × 0.35' × 6 nm for channel four. The resolving power of the MRS observations ranged from 4,000 in channel one to 1,500 in channel four, corresponding to a FWHM of 75 to 200 km s⁻¹, respectively⁴⁶,⁴⁷. Third-level (pixel-level) defringing was applied, but its impact was minimal.

We performed one-dimensional spectral extractions individually for each of the MRS cubes using a circular aperture of radius equal to 1 × FWHM(λ), where FWHM(λ) = 0.3 arcsec for λ < 8 µm and FWHM (λ) = 0.31 × λ[µm]/8 arcsec for λ > 8 µm. The selected FWHM(λ) values follow the MRS point spread function FWHM. We used the MRS point spread function models⁸⁸ to correct for the aperture losses in the one-dimensional extracted spectra. The percentage of the flux that was lost from the selected aperture was ~30% for all MRS channels. Finally, the channels were stitched together by interpolating the flux over the overlapping wavelength regions and averaging without weighting.

## Uncertainties

Due to the complexity of the reduction process, it was not possible to obtain the spectroscopic uncertainties from first principles concurrently with the spectral extraction. Instead, we empirically determined the spectroscopic uncertainties by distributing 12 apertures randomly in the field surrounding the quasar at each slice (spectral pixel) of the datacubes. The 1σ uncertainties were then given by the standard deviation of the fluxes in the apertures obtained using the same procedure as for the quasar.

We checked that the consistency of the uncertainties obtained through this process agree with the observed residuals of the fitting procedure in each subchannel. After a first pass of the fitting procedures (see below), we assumed that the fitting residuals were described well by a random Gaussian distribution. For each MRS subchannel, we fitted a Gaussian distribution to the residuals and used its FWHM as an estimate of the uncertainty over the whole subchannel. We found that the uncertainties estimated from random apertures agree well with the empirical uncertainties from spectral fitting. This has the advantage of capturing the slight increase of variance with wavelength within each subarray. We, therefore, adopted the uncertainties from random apertures throughout the analysis. Finally, we also tested the effect of varying the empirically estimated uncertainties on our fitted parameters by arbitrarily rescaling the uncertainties by factors up to 3, finding no statistically significant effect on the recovered optimal parameter values. Only the parameter uncertainties were affected by this rescaling. The spectral uncertainties were averaged and divided by √2 over the overlapping edges of the subchannels to reflect the improved SNR in those regions.

Although we found no evidence of interchannel flux offsets, the overall MRS flux calibration has an accuracy of only 10%. This systematic effect was taken into account in all measurements of fluxes and luminosities.

## Spectral analysis

We aimed to derive the BH mass, the dust temperature in the torus and the physical properties of the BLR of J1120+0641. To achieve this, we needed a model-fitting method flexible enough to fit the continuum and the complex emission-line profiles simultaneously, and we needed to derive covariant uncertainties between all parameters. We, therefore, employed the Sculptor software, which provides an open-source framework for viewing and analysing quasar spectroscopy⁸⁹.

We added custom functions to Sculptor to model the rest-frame infrared continuum emission and used built-in functions to fit the

emission lines. Sculptor provides a convenient way of specifying bounds for all model parameters and for tying parameters together, as described in the following sections. We used the built-in Monte Carlo Markov chain sampler with 200 walkers and 50 rounds of burn-in to generate parameter constraints and covariance matrices. All model parameters were fitted concurrently.

## Continuum emission

Quasar continuum emission at 0.09 ≲ λ<sub>rest</sub> ≲ 1 µm is dominated by radiation from the accretion disk, which is empirically described well by a power law:

$$F_{\lambda} \propto \lambda^{\alpha}. \qquad (8)$$

Over this wavelength range, several spectral breaks are known to be present, most prominently at λ<sub>rest</sub> ≈ 0.12, 0.3 or 0.6 µm (refs. 64,90). As our MRS spectrum covers the range λ<sub>rest</sub> > 0.61 µm, a single unbroken power law is a sufficient description of the accretion disk emission. We set a broad flat prior on the power-law slope of α ∈ [-2.5, 2.5] and no prior on the power-law amplitude. The power-law parameters were fitted together with all other model parameters over the entire wavelength range except for the minor masking of fringing features (see later).

Based on observations of a strong emission bump in the rest-frame infrared, hot dust is known to be present in the accretion structure surrounding the accretion disk, that is, the torus. The torus dust is traditionally assumed to be highly optically thick, and therefore, it is usually modelled as a pure black body:

$$F_{\nu} \propto \frac{M_{\text{dust}}}{M_{\odot}} \lambda^{-3} \left( \exp\left(\frac{hc}{k\lambda T}\right) - 1 \right)^{-1}. \qquad (9)$$

Hot torus dust is generally thought to consist of graphitic dust due to the lack of strong dust absorption features seen in obscured quasars and because silicate dust sublimes at T ≈ 1,300–1,500 K (refs. 91,92), although hot dust temperatures of T > 1,700 K are seen in Seyfert galaxies (for example, ref. 93). Graphitic hot torus dust instead sublimes at temperatures T ≈ 2,000 K, meaning that it can start to dominate the quasar spectrum at λ ≥ 0.9 µm. Moreover, black-body-emitting dust colder than T ≈ 300 K will not noticeably emit light in the range covered by MIRI-MRS channels 1, 2 or 3 unless the (surface) mass of dust is unphysically large (M<sub>dust</sub> > 10<sup>5</sup> M<sub>☉</sub>). Consequently, we set a prior on T ∈ [300, 3,000] K.

Under the assumption of pure black-body radiation, the mass of dust cannot be constrained without a model for its physical distribution, as only the radiating surface of the dust torus is assumed to be visible. The factor M<sub>dust</sub> in equation (9) is, therefore, better understood as a dust surface area, but we retain the notation for consistency with older works and with the optically thin regime (see below). Because the dependence of equation (9) on T is exponentially stronger than on M<sub>dust</sub>, we set a logarithmically flat prior on log (M<sub>dust</sub>/M<sub>☉</sub>) ∈ [-5, 5].

**More complicated dust models.** It is not known a priori whether hot torus dust emission is indeed optically thick. In principle, another modified black-body term (1 − exp(−τ(v))) can be added as a prefactor to equation (9), where τ(v) is a frequency-dependent optical depth proportional to v. Past studies were generally limited in fitting spectral energy distributions to a small number of photometric points, and both the optically thin and optically thick spectral energy density forms have been employed in the literature to constrain the temperature. The optically thick formula shifts the location of the emission’s peak to longer wavelengths, resulting in larger temperature estimates than the optically thin assumption. The much higher resolution and sensitivity of MIRI-MRS compared to past far-infrared instruments could in principle distinguish between the optically thick and thin
solutions using the shape of the continuum upturn, motivating us to conduct a separate fitting attempt using the optically thin formula.

The proportionality factor between $F_v$ and $M_{dust}$ then depends on the physical cross section of the emitting object, which is unknown for a quasar's hot dust torus. The dust mass, therefore, still cannot be constrained unless we make the assumption of fully optically thin dust ($\tau \to 0$), in which case the spectrum becomes

$$F_v \propto \frac{M_{\text{dust}}}{M_{\odot}} \lambda^{-(3+\beta)} \left( \exp\left(\frac{hc}{\lambda\bar{T}}\right) - 1 \right)^{-1}, \quad (10)$$

where $\beta$ is the modified black-body slope, usually taken to be $1.5 \le \beta \le 2.0$. We fitted the data using both a fixed $\beta = 2$ and a free $\beta \in [1.5, 2]$.

The exponential dependence of the emission on temperature causes the rest-frame infrared spectrum to be dominated by a relatively much smaller amount of the hottest dust that is present. We, therefore, began by considering dust at a single temperature, corresponding to the hottest dust present at the inner edge of the accretion torus. This single-temperature assumption is natural in the black-body radiation case, as only radiation originating from the surface of a very dense torus is assumed to be visible. However, interferometric observations of low-z active galactic nuclei, for example with GRAVITY (ref. 94), have revealed that (somewhat colder) torus dust can be entrained along the polar axis$^{95,96}$. The presence of a small amount of optically thin hot dust is, therefore, plausible in principle for luminous quasars (but see refs. 97,98 for the lack of such an effect in fainter active galactic nuclei, and see also results from the Matisse instrument revealing the more complicated dust distribution in ref. 99).

We present the results of this analysis in Supplementary Fig. 1. As we found no statistically significant deviation of the spectrum from the one-temperature continuum, we did not consider more complicated dust models in which the emission arises from a distribution of dust masses and temperatures. Such models would entail a large increase in the number of model parameters, which we would not be able to meaningfully constrain, and their interpretation would be, anyway, hindered by our lack of knowledge of $\tau(v)$.

## Emission lines

Broad quasar emission lines from the BLR are known to contain several Gaussian (or Lorentzian) components of different widths, with the broadest components typically showing the strongest outflows with respect to the systemic redshift. Accordingly, we tried to model the broad lines with two Gaussian components when the data quality allowed it. More details on individual broad lines are given below. The broad-line parameters were fitted along with the continuum emission parameters to take potential covariances into account in the posteriors. All best-fitting line parameters and their uncertainties are given in Supplementary Table 1.

We calculated emission-line fluxes $f_X$ by integrating the (double) Gaussian best-fitting models and the line luminosities with $L_X = 4\pi D_X^2 f_X$. The bolometric luminosity of J1120+0641 is $L_{bol} = 1.705 \pm 0.005 \times 10^{47}$ erg s$^{-1}$, as calculated by ref. 11 based on the quasar's luminosity at 3,000 Å using the empirical relation of ref. 100 as updated by ref. 101:

$$L_{\text{bol}} = 5.15 \times (3,000 \AA L_{\lambda,3,000 \AA}). \quad (11)$$

Finally, the Eddington ratio for accretion for J1120+0641 was obtained with

$$\lambda_{\text{Edd}} = \frac{L_{\text{bol}}}{L_{\text{Edd}}} = \frac{L_{\text{bol}}}{1.257 \times 10^{38}} \left( \frac{M_{\text{BH}}}{M_{\odot}} \right)^{-1}, \quad (12)$$

where $L_{Edd}$ is the Eddington accretion rate for a BH with mass $M_{BH}$. The intrinsic scatter in the $L_{bol}$ scaling relation has a factor of ~2, so that

it dominates the error budget compared to measurement uncertainties. However, by using the same value of $L_{bol}$ as past literature, we ensured that comparisons of $\lambda_{Edd}$ across several emission lines are meaningful.

In our fitting analysis, we neglected the instrument broadening of the MRS, which introduces a wavelength-dependent FWHM of 75–200 km s$^{-1}$. All the quasar emission features that we analysed are substantially broader, such that the effect of resolution is at most 10 times weaker than our fitting uncertainties for the FWHM of emission lines.

**Hα.** The broad Hα line was modelled with two Gaussian components with free centroids, amplitudes and widths, amounting to six free parameters. We initiated the fit with a reasonable first guess obtained by eye in the Sculptor interface. Both broad Gaussian components were given broad flat priors for the centroid in [5.2, 5.4] µm and the width in [300, 20,000] km s$^{-1}$ and an unconstrained amplitude. We found that all parameters were constrained well with no need for an additional narrow-line component of Hα and no signs of [N II] emission. The lack of stellar emission lines is consistent with the relative faintness of the host galaxy in submillimetric [C II] (ref. 34). The host galaxy has an estimated mass in the range $6 \times 10^9$ to $4 \times 10^{10} M_{\odot}$, only five to 20 times more massive than the SMBH (refs. 102,103). We, therefore, did not include an additional narrow component for the Hα line.

BH masses based on the broad Hβ line have been calibrated well from reverberation-mapping studies. Since the widths of the Hα and Hβ lines follow a tight relationship, Hα can also be used to calculate the BH virial mass on its own$^{101,104}$:

$$\log\left(\frac{M_{\text{BH,Hα}}}{M_{\odot}}\right) = 0.379 + 0.43 \log\left(\frac{L_{\text{Hα}}}{10^{42} \text{ erg s}^{-1}}\right) + 2.1 \log\left(\frac{\text{FWHM}_{\text{Hα}}}{\text{km s}^{-1}}\right). \quad (13)$$

This scaling relation has an intrinsic systematic scatter of a factor of ~2, which is, therefore, the dominating uncertainty in our measurement, far above the uncertainties of the parameters $L_{Hα}$ and FWHM$_{Hα}$, which were calculated from fitting emission lines.

FWHM$_{Hα}$ is the width of the total line model where its value is equal to half of its maximum. For all emission lines, we preferred the BH mass scaling relations based on FWHM instead of $\sigma$ since the high SNRs of our observations mean that a single Gaussian fit is a poor model of the lines, leading to high values of $\xi^2$, which ultimately broaden the uncertainties on emission-line parameters.

**Paα.** Like Hα, we attempted to model the Paα (1.8756 µm) line with two broad Gaussian components. Although the central component was constrained well, the broad component's centroid did not converge. We suspect that this was caused by low-frequency fringing in the spectrum, which is visible at $15.5 < \lambda < 16.5$ µm and may extend over the neighbouring Paα emission line (Fig. 3). To surmount this difficulty, we fixed the redshift of the broadest line component to be the same as for the Hα emission line. The centroid of the other component was given a flat prior in [15.0, 15.5] µm, with the same width and amplitude priors as the Hα emission line.

A BH mass estimate can be obtained from the Paα directly$^{67}$:

$$\frac{M_{\text{BH,Paα}}}{M_{\odot}} = 10^{7.31} \left( \frac{L_{\text{Paα}}}{10^{42} \text{ erg s}^{-1}} \right)^{0.48} \left( \frac{\text{FWHM}_{\text{Paα}}}{10^3 \text{ km s}^{-1}} \right)^{1.68}, \quad (14)$$

where $L_{Paα}$ and FWHM$_{Paα}$ are the luminosity and the FWHM of the line, respectively. The intrinsic scatter of the Paα-based BH mass estimate is potentially very low (~0.2 dex; ref. 67), but it has been tested on fewer BHs with known mass from reverberation mapping compared to scalings based on Hα.
Paβ. The Paβ (1.2822 µm) line was similarly modelled with two Gaussian components with the same amplitude and width constraints as the previously discussed lines and a prior on the Gaussian centroids in [10.3, 10.5] µm. The presence of a spurious fringing feature (or unidentified narrow emission line) at λ = 10.223 µm biased the fitting, and therefore, we masked this feature. The fit is shown in Fig. 3.

A recent BH mass estimator based on the Paβ line is given by$^{73}$:

$$ \frac{M_{\text{BH,Paβ}}}{M_{\odot}} = 10^{7.04} \left( \frac{L_{\text{Paβ}}}{10^{42} \text{ erg s}^{-1}} \right)^{0.48} \left( \frac{\text{FWHM}_{\text{Paβ}}}{10^{3} \text{ km s}^{-1}} \right)^{2} . \quad (15) $$

Like the estimate based on the Paα line, the Paβ estimator has a low intrinsic scatter of ~0.2 dex, but has been tested on smaller samples of quasars with reverberation-mapping measurements of the BH mass$^{67}$.

**Other emission lines.** We modelled the Pay (1.0941 µm) + He I (1.0830 µm) complex at λ ≈ 8.82 µm with two Gaussians (Fig. 3). To distinguish the two components, we constrained the centroid of the He I component to [8.6, 8.77] µm and the centroid of Pay to [8.77, 9.0] µm. Due to the relatively low SNR of ~4.5 per pixel, the fit posterior did not completely drop at the boundary between the centroids, that is, the components were not fully separated at the 3σ level. However, we were still able to obtain posteriors on the line widths at 2σ.

The remaining detected emission lines, O I (0.8446 µm) and Paδ (1.0052 µm), were modelled with single Gaussians and were constrained well (Fig. 3).

It is possible to derive a BH mass directly from the rest-frame infrared Paschen-series lines and the rest-frame infrared continuum. Following the method of ref. 61, the Paα and Paβ lines can be calibrated to reverberation-mapping measurements, which yields the relation:

$$ \log M_{\text{BH,IR}} = (0.88 \pm 0.04) \times \left( 2 \log \frac{\text{FWHM}_{\text{H}\alpha}}{\text{km s}^{-1}} + 0.5 \log \frac{\text{vL}_{\text{µm}}}{\text{erg s}^{-1}} \right) - (17.39 \pm 1.02), \quad (16) $$

where FWHM$_{H\alpha}$ is the FWHM of one or several Paschen-series hydrogen lines (under the assumption that they originate in a common BLR). Although this relation has 1 dex of intrinsic scatter, we will be able to check its consistency with other BH mass tracers.

## Cloudy modelling

To extract information about the hydrogen BLR from the ratio of the Hα and Paschen series lines, we implemented a simple model of the BLR with the software Cloudy (ref. 105). We used a simplified version of the model proposed by ref. 106, in which the frequency-dependent emission F(v) is incident on a single BLR cloud with density N$_{H}$I. The form of F(v) consists of an ultraviolet bump, power-law emission over the X-rays and a high-energy cutoff:

$$ F(v) \propto v^{-0.2} e^{-hv/kT_{\text{cut}}} e^{-0.136 \text{eV}/hv} + \alpha v^{-1.8}, \quad (17) $$

where T$_{cut}$ = 150,000 K and α is normalized such that F(2 keV)/F(2,500 Å) = 2.25 × 10⁻⁴. The values of the ultraviolet slope and X-ray cutoff are consistent with observations of J1120+0641 within uncertainties. We parametrized the strength of the incident emission with the parameter U (the 'ionization parameter') such that U = Φ/cN$_{H}$I where Φ(s⁻¹cm⁻²) is the photon flux. Emission-line ratios remained nearly constant with U. We refer interested readers to ref. 106 for more details of the model's physical derivation and to ref. 68 for an application to the hydrogen BLR. Unlike those authors, we did not include metals in the BLR cloud, and we assumed no turbulence of the gas. We recorded the output of the model for three values of log N$_{H}$I(cm⁻³) = 10, 12, 14 and values of log U = -7, -6, ..., 0.

## Data availability

All raw data used in this work can be retrieved from the public Mikulski Archive for Space Telescopes, filed under Program ID 1263 (P.I. Colina).

## Code availability

The tools used in the data analysis are publicly available: the JWST MRS data-reduction pipeline⁸¹ (with improvements as described and referenced in Methods) and the fitting and visualization tool Sculptor⁸⁹.

## References

1. Mortlock, D. J. et al. A luminous quasar at a redshift of z=7.085. Nature 474, 616–619 (2011).

2. Wang, F. et al. The discovery of a luminous broad absorption line quasar at a redshift of 7.02. Astrophys. J. Lett. 869, L9 (2018).

3. Yang, J. et al. Exploring reionization-era quasars. IV. Discovery of six new z≥6.5 quasars with DES, VHS, and unWISE photometry. Astron. J. 157, 236 (2019).

4. Yang, J. et al. Pōniuā'ena: a luminous z=7.5 quasar hosting a 1.5 billion solar mass black hole. Astrophys. J. Lett. 897, L14 (2020).

5. Yang, J. et al. Probing early supermassive black hole growth and quasar evolution with near-infrared spectroscopy of 37 reionization-era quasars at 6.3<z≤7.64. Astrophys. J. 923, 262 (2021).

6. Wang, F. et al. A luminous quasar at redshift 7.642. Astrophys. J. Lett. 907, L1 (2021).

7. Volonteri, M. The formation and evolution of massive black holes. Science 337, 544 (2012).

8. Pacucci, F., Natarajan, P., Volonteri, M., Cappelluti, N. & Urry, C. M. Conditions for optimal growth of black hole seeds. Astrophys. J. Lett. 850, L42 (2017).

9. Inayoshi, K., Visbal, E. & Haiman, Z. The assembly of the first massive black holes. Annu. Rev. Astron. Astrophys. 58, 27–97 (2020).

10. Shen, Y. et al. Gemini GNIRS near-infrared spectroscopy of 50 quasars at z≥5.7. Astrophys. J. 873, 35 (2019).

11. Farina, E. P. et al. The X-shooter/ALMA sample of quasars in the epoch of reionization. II. Black hole masses, Eddington ratios, and the formation of the first quasars. Astrophys. J. 941, 106 (2022).

12. Lai, S. et al. Chemical abundance of z-6 quasar broad-line regions in the XQR-30 sample. Mon. Not. R. Astron. Soc. 513, 1801–1819 (2022).

13. De Rosa, G. et al. Evidence for non-evolving Fe II/Mg II ratios in rapidly accreting z-6 QSOs. Astrophys. J. 739, 56 (2011).

14. Mazzucchelli, C. et al. Physical properties of 15 quasars at z>6.5. Astrophys. J. 849, 91 (2017).

15. Meyer, R. A., Bosman, S. E. I. & Ellis, R. S. New constraints on quasar evolution: broad-line velocity shifts over 1.5≤z≤7.5. Mon. Not. R. Astron. Soc. 487, 3305–3323 (2019).

16. Schindler, J.-T. et al. The X-SHOOTER/ALMA sample of quasars in the epoch of reionization. I. NIR spectral modeling, iron enrichment, and broad emission line properties. Astrophys. J. 905, 51 (2020).

17. Bischetti, M. et al. Suppression of black-hole growth by strong outflows at redshifts 5.8–6.6. Nature 605, 244–247 (2022).

18. Bischetti, M. et al. The fraction and kinematics of broad absorption line quasars across cosmic time. Astrophys. J. 952, 44 (2023).

19. Vito, F. et al. High-redshift AGN in the Chandra deep fields: the obscured fraction and space density of the sub-L population. Mon. Not. R. Astron. Soc. 473, 2378–2406 (2018).

20. Trebitsch, M., Volonteri, M. & Dubois, Y. Black hole obscuration and duty-cycles mediated by AGN feedback in high-redshift galaxies. Mon. Not. R. Astron. Soc. 487, 819–831 (2019).

21. Davies, F. B., Hennawi, J. F. & Eilers, A.-C. Evidence for low radiative efficiency or highly obscured growth of z>7 quasars. Astrophys. J. Lett. 884, L19 (2019).
22. Eilers, A.-C., Davies, F. B. & Hennawi, J. F. The opacity of the intergalactic medium measured along quasar sightlines at z~6. *Astrophys. J.* **864**, 53 (2018).

23. Eilers, A.-C. et al. Detecting and characterizing young quasars. I. Systemic redshifts and proximity zone measurements. *Astrophys. J.* **900**, 37 (2020).

24. Eilers, A.-C., Hennawi, J. F., Davies, F. B. & Simcoe, R. A. Detecting and characterizing young quasars. II. Four quasars at z~6 with lifetimes <10⁴ yr. *Astrophys. J.* **917**, 38 (2021).

25. Rieke, G. H. et al. The mid-infrared instrument for the James Webb Space Telescope. I. Introduction. *Publ. Astron. Soc. Pac.* **127**, 584 (2015).

26. Wright, G. S. et al. The mid-infrared instrument for the James Webb Space Telescope. II. Design and build. *Publ. Astron. Soc. Pac.* **127**, 595 (2015).

27. Wright, G. S. et al. The mid-infrared instrument for JWST and its in-flight performance. *Publ. Astron. Soc. Pac.* **135**, 048003 (2023).

28. Bañados, E. et al. An 800-million-solar-mass black hole in a significantly neutral Universe at a redshift of 7.5. *Nature* **553**, 473–476 (2018).

29. Matsuoka, Y. et al. Discovery of the first low-luminosity quasar at z>7. *Astrophys. J. Lett.* **872**, L2 (2019).

30. Bosman, S. All z>5.7 quasars currently known. Zenodo https://doi.org/10.5281/zenodo.3634964 (2020).

31. Guilloteau, S. et al. The IRAM interferometer on Plateau de Bure. *Astron. Astrophys.* **262**, 624 (1992).

32. Wootten, A. & Thompson, A. R. The Atacama Large Millimeter/ Submillimeter Array. Proc. IEEE **97**, 1463–1471 (2009).

33. Venemans, B. P. et al. Detection of atomic carbon [C II] 158 µm and dust emission from a z=7.1 quasar host galaxy. *Astrophys. J. Lett.* **751**, L25 (2012).

34. Venemans, B. P. et al. Kiloparsec-scale ALMA imaging of [C II] and dust continuum emission of 27 quasar host galaxies at z~6. *Astrophys. J.* **904**, 130 (2020).

35. Vernet, J. et al. X-shooter, the new wide band intermediate resolution spectrograph at the ESO Very Large Telescope. *Astron. Astrophys.* **536**, A105 (2011).

36. Bosman, S. E. I. et al. A deep search for metals near redshift 7: the line of sight towards ULAS J1120+0641. *Mon. Not. R. Astron. Soc.* **470**, 1919 (2017).

37. Bosman, S. E. I. & Becker, G. D. Re-examining the case for neutral gas near the redshift 7 quasar ULAS J1120+0641. *Mon. Not. R. Astron. Soc.* **452**, 1105–1111 (2015).

38. Momjian, E., Carilli, C. L., Walter, F. & Venemans, B. The highest redshift quasar at z=7.085: a radio-quiet source. *Astron. J.* **147**, 6 (2014).

39. Bañados, E. et al. Constraining the radio-loud fraction of quasars at z>5.5. *Astrophys. J.* **804**, 118 (2015).

40. Page, M. J. et al. X-rays from the redshift 7.1 quasar ULAS J1120+0641. *Mon. Not. R. Astron. Soc.* **440**, L91–L95 (2014).

41. Moretti, A. et al. X-ray observation of ULAS J1120+0641, the most distant quasar at z=7.08. *Astron. Astrophys.* **563**, A46 (2014).

42. Wang, F. et al. A significantly neutral intergalactic medium around the luminous z=7 quasar J0252-0503. *Astrophys. J.* **896**, 23 (2020).

43. Gardner, J. P. et al. The James Webb Space Telescope mission. *Publ. Astron. Soc. Pac.* **135**, 068001 (2023).

44. Wells, M. et al. The mid-infrared instrument for the James Webb Space Telescope. VI. The medium resolution spectrometer. *Publ. Astron. Soc. Pac.* **127**, 646 (2015).

45. Argyriou, I. et al. JWST MIRI flight performance: the medium-resolution spectrometer. *Astron. Astrophys.* **675**, A111 (2023).

46. Labiano, A. et al. Wavelength calibration and resolving power of the JWST MIRI medium resolution spectrometer. *Astron. Astrophys.* **656**, A57 (2021).

47. Jones, O. C. et al. Observations of the planetary nebula SMP LMC 058 with the JWST MIRI medium resolution spectrometer. *Mon. Not. R. Astron. Soc.* **523**, 2519–2529 (2023).

48. Planck Collaboration, et al. Planck 2018 results. VI. Cosmological parameters. *Astron. Astrophys.* **641**, A6 (2020).

49. Barvainis, R. Hot dust and the near-infrared bump in the continuum spectra of quasars and active galactic nuclei. *Astrophys. J.* **320**, 537 (1987).

50. Kishimoto, M., Hönig, S. F., Beckert, T. & Weigelt, G. The innermost region of AGN tori: implications from the HST/NICMOS type 1 point sources and near-IR reverberation. *Astron. Astrophys.* **476**, 713–721 (2007).

51. Mor, R. & Netzer, H. Hot graphite dust and the infrared spectral energy distribution of active galactic nuclei. *Mon. Not. R. Astron. Soc.* **420**, 526–541 (2012).

52. Minezaki, T. et al. Inner size of a dust torus in the Seyfert 1 galaxy NGC 4151. *Astrophys. J. Lett.* **600**, L35–L38 (2004).

53. Suganuma, M. et al. Reverberation measurements of the inner radius of the dust torus in nearby Seyfert 1 galaxies. *Astrophys. J.* **639**, 46–63 (2006).

54. Efstathiou, A. & Rowan-Robinson, M. Dusty discs in active galactic nuclei. *Mon. Not. R. Astron. Soc.* **273**, 649–661 (1995).

55. Manske, V. & Henning, T. Two-dimensional radiative transfer with transiently heated particles: methods and applications. *Astron. Astrophys.* **337**, 85–95 (1998).

56. Schartmann, M., Meisenheimer, K., Camenzind, M., Wolf, S. & Henning, T. Towards a physical model of dust tori in active galactic nuclei. Radiative transfer calculations for a hydrostatic torus model. *Astron. Astrophys.* **437**, 861–881 (2005).

57. Schartmann, M. et al. Three-dimensional radiative transfer models of clumpy tori in Seyfert galaxies. *Astron. Astrophys.* **482**, 67–80 (2008).

58. Hönig, S. F. & Kishimoto, M. The dusty heart of nearby active galaxies. II. From clumpy torus models to physical properties of dust around AGN. *Astron. Astrophys.* **523**, A27 (2010).

59. Riffel, R. A., Storchi-Bergmann, T., Dors, O. L. & Winge, C. AGN-starburst connection in NGC7582: Gemini near-infrared spectrograph integral field unit observations. *Mon. Not. R. Astron. Soc.* **393**, 783–797 (2009).

60. Mor, R. & Trakhtenbrot, B. Hot-dust clouds with pure-graphite composition around type-I active galactic nuclei. *Astrophys. J. Lett.* **737**, L36 (2011).

61. Landt, H. et al. A near-infrared relationship for estimating black hole masses in active galactic nuclei. *Mon. Not. R. Astron. Soc.* **432**, 113–126 (2013).

62. Jiang, L. et al. Dust-free quasars in the early Universe. *Nature* **464**, 380–383 (2010).

63. Leipski, C. et al. Spectral energy distributions of QSOs at z>5: common active galactic nucleus-heated dust and occasionally strong star-formation. *Astrophys. J.* **785**, 154 (2014).

64. Hernán-Caballero, A., Hatziminaoglou, E., Alonso-Herrero, A. & Mateos, S. The near-to-mid infrared spectrum of quasars. *Mon. Not. R. Astron. Soc.* **463**, 2064–2078 (2016).

65. Gehrels, N. Confidence limits for small numbers of events in astrophysical data. *Astrophys. J.* **303**, 336–346 (1986).

66. Soifer, B. T. et al. Spitzer infrared spectrograph (IRS) observations of the redshift 3.91 quasar APM 08279+5255. *Astrophys. J. Suppl. Ser.* **154**, 151–154 (2004).

67. Kim, D., Im, M. & Kim, M. New estimators of black hole mass in active galactic nuclei with hydrogen Paschen lines. *Astrophys. J.* **724**, 386–399 (2010).

68. Oyabu, S. et al. AKARI near- and mid-infrared spectroscopy of APM 08279+5255 at z=3.91. *Astrophys. J.* **697**, 452–457 (2009).

69. Vestergaard, M. & Peterson, B. M. Determining central black hole masses in distant active galaxies and quasars. II. Improved optical and UV scaling relationships. *Astrophys. J.* **641**, 689–709 (2006).
70. Coatman, L. et al. Correcting C iv-based virial black hole masses. *Mon. Not. R. Astron. Soc.* **465**, 2120–2142 (2017).

71. Marshall, M. A. et al. GA-NIFS: black hole and host galaxy properties of two $z \approx 6.8$ quasars from the NIRSpec IFU. *Astrophys. J.* **678**, A191 (2023).

72. Loiacono, F. et al. A quasar-galaxy merger at $z \approx 6.2$: black hole mass and quasar properties from the NIRSpec spectrum. Preprint at https://doi.org/10.48550/arXiv.2402.13319 (2024).

73. Kim, Y. et al. discovery of a faint quasar at $z=6$ and implications for cosmic reionization. *Astrophys. J. Lett.* **813**, L35 (2015).

74. Kwan, J. & Krolik, J. H. The formation of emission lines in quasars and Seyfert nuclei. *Astrophys. J.* **250**, 478–507 (1981).

75. Grandi, S. A. Reddening indicators for quasars and Seyfert 1 galaxies. *Astrophys. J.* **268**, 591–601 (1983).

76. Matsuoka, Y., Oyabu, S., Tsuzuki, Y. & Kawara, K. Observations of O I and Ca II emission lines in quasars: implications for the site of Fe II line emission. *Astrophys. J.* **663**, 781–798 (2007).

77. Landt, H. et al. The near-infrared broad emission line region of active galactic nuclei. I. The observations. *Astrophys. J. Suppl. Ser.* **174**, 282–312 (2008).

78. Osterbrock, D. E. Astrophysics of Gaseous Nebulae and Active Galactic Nuclei (Univ. Science Books, 1989).

79. Yang, J. et al. A spectroscopic survey of biased halos in the reionization era (ASPIRE): a first look at the rest-frame optical spectra of $z>6.5$ quasars using JWST. *Astrophys. J. Lett.* **951**, L5 (2023).

80. Campanelli, M., Lousto, C., Zlochower, Y. & Merritt, D. Large merger recoils and spin flips from generic black hole binaries. *Astrophys. J. Lett.* **659**, L5–L8 (2007).

81. Labiano, A. et al. The MIRI medium resolution spectrometer calibration pipeline. In *Proc. SPIE Conference Series* Vol. 9910 (eds Peck, A. B. et al.) 99102W (SPIE, 2016).

82. Álvarez-Márquez, J. et al. MIRI/JWST observations reveal an extremely obscured starburst in the $z=6.9$ system SPT0311-58. *Astron. Astrophys.* **671**, A105 (2023).

83. Álvarez-Márquez, J. et al. Nuclear high-ionisation outflow in the Compton-thick AGN NGC 6552 as seen by the JWST mid-infrared instrument. *Astron. Astrophys.* **672**, A108 (2023).

84. Morrison, J. E. et al. JWST MIRI flight performance: detector effects and data reduction algorithms. *Publ. Astron. Soc. Pac.* https://doi.org/10.1088/1538-3873/acdea6 (2024).

85. Patapis, P. et al. Geometric distortion and astrometric calibration of the JWST MIRI medium resolution spectrometer. *Astron. Astrophys.* **682**, A53 (2024).

86. Gasman, D. et al. JWST MIRI/MRS in-flight absolute flux calibration and tailored fringe correction for unresolved sources. *Astron. Astrophys.* **673**, A102 (2023).

87. Law, D. R. et al. A 3D drizzle algorithm for JWST and practical application to the MIRI medium resolution spectrometer. *Astron. J.* **166**, 45 (2023).

88. Patapis, P. et al. Geometric distortion and astrometric calibration of the JWST MIRI Medium Resolution Spectrometer. *Astron. Astrophys.* https://doi.org/10.1051/0004-6361/202347339 (2024).

89. Schindler, J.-T. Sculptor: interactive modeling of astronomical spectra. Astrophysics Source Code Library, record ascl:2202.018 (ASCL, 2022).

90. Glikman, E., Helfand, D. J. & White, R. L. A near-infrared spectral template for quasars. *Astrophys. J.* **640**, 579–591 (2006).

91. Salpeter, E. E. Formation and destruction of dust grains. *Annu. Rev. Astron. Astrophys.* **15**, 267–293 (1977).

92. Lodders, K. Solar System abundances and condensation temperatures of the elements. *Astrophys. J.* **591**, 1220–1247 (2003).

93. Landt, H., Ward, M. J., Steenbrugge, K. C. & Ferland, G. J. Strong variability of the coronal line region in NGC 5548. *Mon. Not. R. Astron. Soc.* **454**, 3688–3696 (2015).

94. GRAVITY Collaboration First light for GRAVITY: phase referencing optical interferometry for the Very Large Telescope Interferometer. *Astron. Astrophys.* **602**, A94 (2017).

95. López-Gonzaga, N., Burtscher, L., Tristram, K. R. W., Meisenheimer, K. & Schartmann, M. Mid-infrared interferometry of 23 AGN tori: on the significance of polar-elongated emission. *Astron. Astrophys.* **591**, A47 (2016).

96. Hönig, S. F. & Kishimoto, M. Dusty winds in active galactic nuclei: reconciling observations with models. *Astrophys. J. Lett.* **838**, L20 (2017).

97. GRAVITY Collaboration The resolved size and structure of hot dust in the immediate vicinity of AGN. *Astron. Astrophys.* **635**, A92 (2020).

98. GRAVITY Collaboration An image of the dust sublimation region in the nucleus of NGC 1068. *Astron. Astrophys.* **634**, A1 (2020).

99. Gámez Rosas, V. et al. Thermal imaging of dust hiding the black hole in NGC 1068. *Nature* **602**, 403–407 (2022).

100. Richards, G. T. et al. Spectral energy distributions and multi-wavelength selection of type 1 quasars. *Astrophys. J. Suppl. Ser.* **166**, 470–497 (2006).

101. Shen, Y. et al. A catalog of quasar properties from Sloan Digital Sky Survey Data Release 7. *Astrophys. J. Suppl. Ser.* **194**, 45 (2011).

102. Venemans, B. P. et al. The compact, ~1kpc host galaxy of a quasar at a redshift of 7.1. *Astrophys. J.* **837**, 146 (2017).

103. Yue, M. et al. EIGER V. Characterizing the host galaxies of luminous quasars at $z \ge 6$. *Astrophys. J.* https://doi.org/10.3847/1538-4357/ad3914 (2024).

104. Greene, J. E. & Ho, L. C. Estimating black hole masses in active galaxies using the Hα emission line. *Astrophys. J.* **630**, 122–129 (2005).

105. Ferland, G. J. et al. The 2017 release Cloudy. *Rev. Mex. Astron. Astrofis.* **53**, 385–438 (2017).

106. Tsuzuki, Y. et al. Fe II emission in 14 low-redshift quasars. I. Observations. *Astrophys. J.* **650**, 57–79 (2006).

## Acknowledgements

S.E.I.B. is supported by the German Research Foundation (Emmy Noether Grant No. BO 5771/1-1). S.E.I.B., F.W. and L.B. acknowledge funding from the European Research Council under the European Union's Horizon 2020 research and innovation programme (Grant Agreement No. 740246 'Cosmic Gas'). L. Colina acknowledges financial support from the Community of Madrid (Talent Attraction Grant No. 2018-T2/TIC-11612). The Cosmic Dawn Center (DAWN) is funded by the Danish National Research Foundation (Grant No. 140). J.A.-M., L. Colina and A.L. acknowledge support from the Spanish Ministry of Science and Innovation (Project MCIN/AEI/10.13039/501100011033, Grant No. PIB2021-127718NB-100). A.A.-H. is supported by the Spanish Ministry of Science and Innovation/State Agency of Research (Project MCIN/AEI/10.13039/501100011033, Grant No. PID2021-124665NB-100) and by 'ERDF A way of making Europe'. A.B. and G.Ö. acknowledge support from the Swedish National Space Administration. O.I. acknowledges funding from the French National Agency for Research (Project iMAGE, Grant No. ANR-22-CE31-0007). J.H. and D.L. were supported by the Villum Foundation (Investigator Grant for Project No. 16599). K.I.C. acknowledges funding from the Netherlands Research School for Astronomy and the Dutch Research Council (Vici Grant No. VI.C.212.036). T.P.R. would like to acknowledge support from the European Research Council (Advanced Grant No. 743029 (EASY)). A.E. and F.P. acknowledge support through the German Space Agency (Grant Nos. DLR 500S1501 and DLR 500S2001 from 2015 to 2023). J.H. was supported by the Villum Foundation (Research Grant No. VIL54489). The work presented is the effort of the entire MIRI team and the enthusiasm within the MIRI partnership is a crucial factor
in its success. MIRI draws on the scientific and technical expertise of the following organizations: Ames Research Center, USA; Airbus Defence and Space, UK; CEA-Irfu, Saclay, France; Liège Space Center, Belgium; National Research Council, Spain; Carl Zeiss Optronics, Germany; Chalmers University of Technology, Sweden; Danish Space Research Institute, Denmark; Dublin Institute for Advanced Studies, Ireland; European Space Agency, Netherlands; ETCA, Belgium; ETH Zurich, Switzerland; the Goddard Space Flight Center, USA; Institute d'Astrophysique Spatiale, France; National Institute of Aerospace Technology, Spain; the Institute for Astronomy, Edinburgh, UK; Jet Propulsion Laboratory, USA; Laboratoire d'Astrophysique de Marseille, France; Leiden University, Netherlands; Lockheed Advanced Technology Center, USA; the Optical and Infrared Group of the Netherlands Research School for Astronomy at Dwingeloo, Netherlands; Northrop Grumman, USA; the Max Planck Institute for Astronomy, Heidelberg, Germany; Laboratoire d'Etudes Spatiales et d'Instrumentation en Astrophysique (LESIA), France; Paul Scherrer Institut, Switzerland; Raytheon Vision Systems, USA; RUAG Aerospace, Switzerland; Rutherford Appleton Laboratory (RAL Space), UK; the Space Telescope Science Institute, USA; Netherlands Organisation for Applied Scientific Research, Netherlands; UK Astronomy Technology Centre, UK; University College London, UK; University of Amsterdam, Netherlands; University of Arizona, USA; University of Cardiff, UK; University of Cologne, Germany; University of Ghent, Belgium; University of Groningen, Netherlands; University of Leicester, UK; University of Leuven, Belgium; University of Stockholm, Sweden; Utah State University, USA. A portion of this work was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with NASA. We would like to thank the following national and international funding agencies for their support of the MIRI development: NASA, European Space Agency, the Belgian Science Policy Office, the French National Centre for Space Studies, the Danish National Space Centre, the German Aerospace Centre, Enterprise Ireland, Ministerio De Economía y Competitividad, the Netherlands Research School for Astronomy, the Netherlands Organisation for Scientific Research, the Science and Technology Facilities Council, the Swiss Space Office, the Swedish National Space Board and the UK Space Agency. This work is based on observations made with the NASA/ESA/CSA JWST. Some data were obtained from the Mikulski Archive for Space Telescopes at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-03127 for JWST; and from the European JWST archive operated by the Science Data Centre of the European Space Astronomy Centre.

## Author contributions

S.E.I.B. performed the data analysis, model fitting, comparison with the literature and authored most of the manuscript. J.A.-M. performed the data reduction and authored the corresponding part of the manuscript. L. Colina designed the experiment, chose the target and led the observational campaign. F.W., A.A.-H., M.W., G.O., T.R.G. and G.W. provided detailed feedback, edited the manuscript, helped design the means of comparison of the data with past literature and contributed relevant references. A.B., L.B., K.C., L. Costantin, A.E., M.G.-M., S.G., J.H., E.I., O.I., I.J., A.L., D.L., F.P., P.R., M.T. and P.v.d.W. provided comments on the analysis, helped with the interpretation of the results and contributed their expertise with the MIRI instrument as members of the MIRI team. M.G., T.H., P.-O.L., T.R., E.F.v.D., B.V., L.C. and G.W. are the leaders of the MIRI European Consortium, on behalf of which this work was conducted.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary information** The online version contains supplementary material available at https://doi.org/10.1038/s41550-024-02273-0.

**Correspondence and requests for materials** should be addressed to Sarah E. I. Bosman.

**Peer review information** *Nature Astronomy* thanks Yuichi Harikane and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

© The Author(s), under exclusive licence to Springer Nature Limited 2024

<footnote n="1">¹Institute for Theoretical Physics, Heidelberg University, Heidelberg, Germany. ²Max-Planck-Institut für Astronomie, Heidelberg, Germany. ³Centro de Astrobiología (CAB), CSIC-INTA, Madrid, Spain. ⁴Centro de Astrobiología (CAB), CSIC-INTA, Villanueva de la Cañada, Spain. ⁵Centre for Extragalactic Astronomy, Department of Physics, Durham University, Durham, UK. ⁶Department of Astronomy, Oscar Klein Centre, Stockholm University, Stockholm, Sweden. ⁷Cosmic Dawn Center (DAWN), Copenhagen, Denmark. ⁸DTU Space, Technical University of Denmark, Kongens Lyngby, Denmark. ⁹Department of Physics and Astronomy, University College London, London, UK. ¹⁰UK Astronomy Technology Centre, Royal Observatory Edinburgh, Edinburgh, UK. ¹¹Kapteyn Astronomical Institute, University of Groningen, Groningen, The Netherlands. ¹²Physikalisches Institut der Universität zu Köln, Cologne, Germany. ¹³European Space Agency, Space Telescope Science Institute, Baltimore, MD, USA. ¹⁴DARK, Niels Bohr Institute, University of Copenhagen, Copenhagen, Denmark. ¹⁵LAM (Laboratoire d'Astrophysique de Marseille), Aix Marseille Université, CNRS, UMR 7326, Marseille, France. ¹⁶Telespazio UK for the European Space Agency, ESAC, Villanueva de la Cañada, Spain. ¹⁷Astronomy & Astrophysics Section, Dublin Institute for Advanced Studies, Dublin, Ireland. ¹⁸Leiden Observatory, Leiden University, Leiden, The Netherlands. ¹⁹Department of Astrophysics, University of Vienna, Vienna, Austria. ²⁰Institute of Particle Physics and Astrophysics, ETH Zurich, Zurich, Switzerland. ²¹AIM, CEA, CNRS, Université Paris-Saclay, Université Paris Diderot, Sorbonne Paris Cité, Gif-sur-Yvette, France. ²²Institute of Astronomy, KU Leuven, Leuven, Belgium. e-mail: bosman@thphys.uni-heidelberg.de</footnote>