# Deciphering the Nature of Virgil: An Obscured Active Galactic Nucleus Lurking within an Apparently Normal Ly α Emitter during Cosmic Reionization

Pierluigi Rinaldi, Pablo G Pérez-González, George H Rieke, Jianwei Lyu,

Francesco D'eugenio, Zihao Wu, Stefano Carniani, Tobias J Looser, Irene

Shivaei, Leindert A Boogaard, et al.

► To cite this version:

Pierluigi Rinaldi, Pablo G Pérez-González, George H Rieke, Jianwei Lyu, Francesco D'eugenio, et al. Deciphering the Nature of Virgil: An Obscured Active Galactic Nucleus Lurking within an Apparently Normal Ly α Emitter during Cosmic Reionization. The Astrophysical Journal, 2025, 994 (1), pp.86. 10.3847/1538-4357/ae089c . hal-05470354

HAL Id: hal-05470354

https://hal.science/hal-05470354v1

Submitted on 21 Jan 2026

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 CC BY 4.0 - Attribution - International License
# Deciphering the Nature of Virgil: An Obscured Active Galactic Nucleus Lurking within an Apparently Normal Ly$\alpha$ Emitter during Cosmic Reionization

Pierluigi Rinaldi¹, Pablo G. Pérez-González², George H. Rieke¹, Jianwei Lyu¹, Francesco D'Eugenio³,⁴, Zihao Wu⁵, Stefano Carniani⁶, Tobias J. Looser³,⁴, Irene Shivaei², Leindert A. Boogaard⁷, Tanio Diaz-Santos⁸,⁹, Luis Colina², Göran Östlin¹⁰, Stacey Alberts¹, Javier Álvarez-Márquez², Marianna Annuziatella², Manuel Aravena¹¹, Rachana Bhatawdekar¹², Andrew J. Bunker¹³, Karina I. Caputi¹⁴,¹⁵, Stéphane Charlot¹⁶, Alejandro Crespo Gómez¹⁷, Mirko Curti¹⁸, Andreas Eckart¹⁹, Steven Gillman¹⁵,²⁰, Kevin Hainline¹, Nimisha Kumari²¹, Jens Hjorth²², Edoardo Iani²³, Hanae Inami²⁴, Zhiyuan Ji¹, Benjamin D. Johnson⁵, Gareth C. Jones³,⁴, Álvaro Labiano²⁵, Roberto Maiolino³,⁴,²⁶, Jens Melinder¹⁰, Thibaud Moutard¹²,²⁷, Florian Peissker¹⁹, Marcia Rieke¹, Brant Robertson²⁸, Jan Scholtz³,⁴, Sandro Tacchella³,⁴, Paul P. van der Werf⁷, Fabian Walter²⁹, Christina C. Williams³⁰, Chris Willott³¹, Joris Witstok¹⁵,³², Hannah Übler³³, and Yongda Zhu¹

¹Steward Observatory, University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA; prinaldi@arizona.edu

²Centro de Astrobiología (CAB), CSIC-INTA, Cra. de Ajalvir Km. 4, 28850- Torrejón de Ardoz, Madrid, Spain

³Kavli Institute for Cosmology, University of Cambridge, Madingley Road, Cambridge, CB3 0HA, UK

⁴Cavendish Laboratory, University of Cambridge, 19 JJ Thomson Avenue, Cambridge, CB3 0HE, UK

⁵Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge MA 02138, USA

⁶Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy

⁷Leiden Observatory, Leiden University, PO Box 9513, NL-2300 RA Leiden, The Netherlands

⁸Institute of Astrophysics, Foundation for Research and Technology-Hellas (FORTH), Heraklion, 70013, Greece

⁹School of Sciences, European University Cyprus, Diogenes Street, Engomi, 1516 Nicosia, Cyprus

¹⁰Department of Astronomy, Stockholm University, Oscar Klein Centre, AlbaNova University Centre, 106 91 Stockholm, Sweden

¹¹Instituto de Estudios Astrofísicos, Facultad de Ingeniería y Ciencias, Universidad Diego Portales, Av. Ejército 441, Santiago 8370191, Chile

¹²European Space Agency (ESA), European Space Astronomy Centre (ESAC), Camino Bajo del Castillo s/n, 28692 Villanueva de la Cañada, Madrid, Spain

¹³Department of Physics, University of Oxford, Denys Wilkinson Building, Keble Road, Oxford OX1 3RH, UK

¹⁴Kapteyn Astronomical Institute, University of Groningen, P.O. Box 800, 9700AV Groningen, The Netherlands

¹⁵Cosmic Dawn Center (DAWN), Copenhagen, Denmark

¹⁶Sorbonne Université, CNRS, UMR 7095, Institut d'Astrophysique de Paris, 98 bis bd Arago, 75014 Paris, France

¹⁷Space Telescope Science Institute (STScI), 3700 San Martin Drive, Baltimore, MD 21218, USA

¹⁸European Southern Observatory, Karl-Schwarzschild-Strasse 2, 85748 Garching, Germany

¹⁹Physikalisches Institut der Universität zu Köln, Zülpicher Str. 77, 50937 Köln, Germany

²⁰DTU Space, Technical University of Denmark, Elektrovej, Building 328, 2800, Kgs. Lyngby, Denmark

²¹AURA for European Space Agency, Space Telescope Science Institute, 3700 San Martin Drive. Baltimore, MD, 21210, USA

²²DARK, Niels Bohr Institute, University of Copenhagen, Jagtvej 155A, 2200 Copenhagen, Denmark

²³Institute of Science and Technology Austria (ISTA), Am Campus 1, 3400 Klosterneuburg, Austria

²⁴Hiroshima Astrophysical Science Center, Hiroshima University, 1-3-1 Kagamiyama, Higashi-Hiroshima, Hiroshima 739-8526, Japan

²⁵Telespazio UK for the European Space Agency, ESAC, Camino Bajo del Castillo s/n, 28692 Villanueva de la Cañada, Spain

²⁶Department of Physics and Astronomy, University College London, Gower Street, London WC1E 6BT, UK

²⁷Aix Marseille Univ, CNRS, CNES, LAM, Marseille, France

²⁸Department of Astronomy and Astrophysics, University of California, Santa Cruz, 1156 High Street, Santa Cruz, CA 95064, USA

²⁹Max Planck Institut für Astronomie, Königstuhl 17, D-69117, Heidelberg, Germany

³⁰NSF National Optical-Infrared Astronomy Research Laboratory, 950 North Cherry Avenue, Tucson, AZ 85719, USA

³¹NRC Herzberg, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada

³²Niels Bohr Institute, University of Copenhagen, Jagtvej 128, DK-2200, Copenhagen, Denmark

³³Max-Planck-Institut für extraterrestrische Physik (MPE), Gießenbachstraße 1, 85748 Garching, Germany

Received 2025 April 4; revised 2025 August 18; accepted 2025 September 3; published 2025 November 17

## Abstract

We present a comprehensive analysis of the MIRI Extremely Red Object Virgil, a Ly$\alpha$ emitter at $z_{spec} = 6.6379 \pm 0.0035$ with the photometric properties of a Little Red Dot. Leveraging new JWST/MIRI imaging from the MIDIS and PAHSPECS programs, we confirm Virgil's extraordinary nature among galaxies in JADES/GOODS-South, exhibiting a strikingly red NIRCam-to-MIRI color (F444W-F1500W = $2.84 \pm 0.04$ mag). Deep NIRSpec/PRISM spectroscopy from the OASIS program offers key insights into the host galaxy, revealing properties of an average star-forming galaxy during Cosmic Reionization, such as a subsolar metallicity, low-to-moderate dust content, and a relatively high ionization parameter and electron temperature. By estimating the star formation rate of Virgil from UV and H$\alpha$, we find evidence that the galaxy is either entering or fading out of a bursty episode. Although line-ratio diagnostics employed at high $z$ would classify Virgil as an active galactic nucleus (AGN), this classification becomes ambiguous once redshift evolution is considered. Nonetheless, Virgil occupies the same parameter space as recently confirmed AGNs at similar redshifts. The new deep MIRI data at $15 \ \mu m$ reinforce the AGN nature of Virgil, as inferred from multiple spectral energy distribution (SED) fitting

<figure><img src="image_2.png" /></figure>
codes. Virgil’s rising infrared SED and UV excess resemble those of Dust-Obscured Galaxies (DOGs) studied with Spitzer at Cosmic Noon, particularly blue-excess HotDOGs. Our results highlight the need for a multiwavelength approach incorporating MIRI to uncover such extreme sources at $z \gtrsim 6$ and to shed light on the interplay between galaxy evolution and early black hole growth during Cosmic Reionization.

Unified Astronomy Thesaurus concepts: Active galactic nuclei (16); High-redshift galaxies (734); Galaxy evolution (594); Near infrared astronomy (1093); AGN host galaxies (2017); Galaxy formation (595); Photoionization (2060); Spectral energy distribution (2129); Infrared astronomy (786); Galaxies (573); Infrared photometry (792)

# 1. Introduction

Infrared (IR) astronomy has advanced dramatically in recent decades. The Infrared Astronomical Satellite (IRAS; G. Neugebauer et al. 1984) and the Infrared Space Observatory (ISO; G. Helou & M. F. Kessler 1995) provided the first all-sky IR surveys, laying the foundation for extragalactic studies. These efforts paved the way for the Spitzer Space Telescope (R. D. Gehrz et al. 2007)—the last of the Great Observatories<fnref n="34" />—and later, the Herschel Space Observatory (G. L. Pilbratt et al. 2010). While Herschel extended coverage into the far-IR, Spitzer remained unmatched in mid-IR sensitivity and resolution throughout the pre-JWST era (J. P. Gardner et al. 2023).

With the Infrared Array Camera (IRAC; G. G. Fazio et al. 2004) and the Multiband Imaging Photometer (MIPS; G. H. Rieke et al. 2004), Spitzer allowed us to advance our understanding of galaxy evolution. It enabled robust estimates of stellar mass ($M_*$), age, and star formation rate (SFR) in high-redshift galaxies (e.g., C. Papovich et al. 2006; H. Yan et al. 2006; M. Bradač 2020), uncovered massive systems at $z > 6$, and observed some of the most distant sources known at the time—e.g., GN-z11 (P. A. Oesch et al. 2014, 2016) and MACS0647-JD (D. Lam et al. 2019; V. Strait et al. 2020). It also revealed strong nebular emission through photometric excesses (e.g., K.-H. Huang et al. 2016), Balmer breaks, and early dust content (C. Papovich et al. 2006). In parallel, it advanced studies of active galactic nuclei (AGNs; M. Lacy & A. Sajina 2020) and clarified the role of luminous and ultraluminous infrared galaxies (LIRGs/ULIRGs) in driving feedback between starbursts, black holes, and the interstellar medium (ISM; L. Armus et al. 2020). Many of these phenomena—such as early dust and strong emission lines—are now confirmed by JWST photometry and spectroscopy (P. Rinaldi et al. 2023; K. Boyett et al. 2024; A. Kuruvanthodi et al. 2024; D. Langeroodi et al. 2024).

In particular, Spitzer played a key role in advancing the study of Extremely Red Objects (EROs; R. Elston et al. 1988, 1989)—sources with very red mid-IR to optical colors—triggering extensive investigations into their nature (e.g., G. Wilson et al. 2004; H. Yan et al. 2004; D. Stern et al. 2006; G. Wilson et al. 2007). Initially believed to be $z \approx 10$ sources (e.g., B. Mobasher et al. 2005), EROs were later found to include Submillimeter Galaxies (SMGs; A. W. Blain et al. 2002) and Dust-Obscured Galaxies (DOGs; A. Dey et al. 2008). DOGs—categorized as power-law or “bump” types (with a spectral energy distribution (SED) flattening at long wavelengths)—likely trace a transition from dusty star formation to AGN activity during gas-rich mergers, with their red colors driven by heavy obscuration and mid-IR reemission. Additionally, EROs were found to overlap with high-redshift galaxies (R. J. Bouwens et al. 2011) selected via the Lyman break technique (C. C. Steidel et al. 1996), highlighting their heterogeneous nature.

Despite its breakthroughs, Spitzer was limited by its spatial resolution. This limitation has now been overcome by JWST, which combines unparalleled sensitivity, angular resolution, and wavelength coverage with both the Near-Infrared Camera (NIRCam, 0.7–4.8 µm; M. J. Rieke et al. 2023a) and the Mid-Infrared Instrument (MIRI, 5.6–25.5 µm<fnref n="35" />; G. S. Wright et al. 2023). The launch of JWST marks a significant leap forward compared to its predecessors, including the Hubble Space Telescope (HST), heralding a golden era of IR extragalactic astronomy.

Since its launch, JWST has pushed the boundaries of the redshift frontier, with groundbreaking discoveries at very high redshift, such as GS-z12-0, GS-z13-0, GS-z14-0, and GS-z14-1 (E. Curtis-Lake et al. 2023; A. J. Bunker et al. 2024; S. Carniani et al. 2024a) from the JWST Advanced Deep Extragalactic Survey (JADES; A. J. Bunker et al. 2020; D. J. Eisenstein et al. 2023a; M. J. Rieke et al. 2023b) team, as well as other very high-$z$ sources identified by independent teams (e.g., P. Arrabal Haro et al. 2023; B. Wang et al. 2023). These observations have offered unprecedented insights into the very early Universe, with the Near Infrared Spectrograph’s Micro-Shutter Assembly (NIRSpec/MSA; P. Ferruit et al. 2022; P. Jakobsen et al. 2022) playing a crucial role in spectroscopically confirming galaxies at $z > 10$.

With JWST, the ERO “industry” is back in the spotlight, unveiling even more compact and red sources (I. Labbé et al. 2023), now known as “Little Red Dots” (LRDs; J. Matthee et al. 2024). The study of these sources and their nature has triggered a huge amount of literature in a very short time (e.g., L. J. Furtak et al. 2023; M. Killi et al. 2024; V. Kokorev et al. 2023; H. Übler et al. 2023; H. B. Akins et al. 2025; G. Barro et al. 2024a; J. E. Greene et al. 2024; D. D. Kocevski et al. 2025; V. Kokorev et al. 2024; K. N. Hainline et al. 2025; P. G. Pérez-González et al. 2024a; P. Rinaldi et al. 2025a; C. C. Williams et al. 2024).

Recent studies highlighted MIRI’s critical role in uncovering EROs, revealing “exotic” objects undetectable even in the deepest NIRCam surveys, such as JADES (D. J. Eisenstein et al. 2023a; M. J. Rieke et al. 2023b), or showing unexpected SED upturns in objects that NIRCam would typically classify as star-forming galaxies (SFGs). In this context, the MIRI Deep Imaging Survey (MIDIS; G. Östlin et al. 2025) has been fundamental in identifying the first MIRI Extremely Red Objects (MEROs): Cerberus (P. G. Pérez-González et al. 2024b) and Virgil (E. Iani et al. 2025).

<footnote n="34">³⁴ As described in G. H. Rieke (2006).</footnote>

<footnote n="35">³⁵ Here we refer to the MIRI imager. However, the MIRI Medium Resolution Spectroscopy (MRS) extends up to ≈ 27 µm.</footnote>
While the nature of Cerberus—a source uniquely detected at 10 µm—remains entirely unknown because it is so faint, Virgil offers a great opportunity for detailed investigation thanks to its detection across multiple facilities, including HST, JWST (NIRCam and MIRI), and MUSE on the Very Large Telescope (VLT) (R. Bacon et al. 2023). As reported in E. Iani et al. (2025), Virgil is identified as a Lyα emitter (LAE) at $z_{spec} = 6.6312 \pm 0.0019$ (from VLT/MUSE) located in GOODS-South (GOODS-S; M. Dickinson et al. 2003; M. Giavalisco et al. 2004). It exhibits very red colors between NIRCam and MIRI bands (F444W–F1000W > 2 mag; E. Iani et al. 2025), thus revealing that the MIRI photometric information is crucial at $z \gtrsim 6$ to identify it and similar objects, as the strong upturn happens entirely outside the wavelength range of NIRCam.

E. Iani et al. (2025) concluded that Virgil could be either a dusty starburst—similar to what has been proposed for GN20 (e.g., L. Colina et al. 2023)—or an SFG hosting a dust-obscured AGN. Recent JWST/NIRSpec observations of GN20 reveal complex kinematics, extended Hα emission, and signatures of AGN-driven feedback (H. Übler et al. 2024a). However, unlike GN20, Virgil is both much less massive and more compact, making its extreme red colors and potential AGN activity even more remarkable. In this case, no AGN counterparts were identified in any of the existing multiwavelength catalogs of GOODS-S (P. Ranalli et al. 2013; B. Luo et al. 2017; J. Lyu et al. 2022; S. Gillman et al. 2025), including the most recent AGN catalog from J. Lyu et al. (2024) with MIRI data, which is based on the Systematic Mid-infrared Instrument Legacy Extragalactic Survey (SMILES; S. Alberts et al. 2024; G. H. Rieke et al. 2024). In the latter case, this could be due to the limited depth of the SMILES observations, which are not deep enough to clearly detect ($>3\sigma$) this object in the MIRI bands above 7.7 µm.

Moreover, E. Iani et al. (2025) demonstrated that Virgil's SED closely resembles that of a typical LRD (see their Figure 10). Although Virgil does not meet all photometric color criteria for LRDs (e.g., V. Kokorev et al. 2024), it satisfies the compactness criterion in the F444W band, possibly making it one of the first LRDs where the host galaxy is well detected. Interestingly, E. Iani et al. (2025) highlight the peculiar ultraviolet (UV) morphology of Virgil, consistent with the findings of P. Rinaldi et al. (2025a), where nearly 30% of the photometrically selected LRDs exhibit a complex UV morphology, suggesting that interactions may play a significant role in triggering AGNs in LRDs.

In this paper, we take a step forward by performing a detailed analysis of Virgil, leveraging the very recent data collected by the MIDIS (PID: 6511; PI: Östlin), PAHSPECS (PID: 5279; PIs: Shivaei, Díaz-Santos, and Boogaard), and OASIS (PID: 5997; PIs: Looser and D'Eugenio) teams. These include exceptionally deep MIRI imaging at 5.6, 7.7, 10, and 15 µm, and deep NIRSpec/PRISM data. The deep NIRSpec/PRISM data enable the most comprehensive analysis of an LRD to date, providing critical insights into Virgil's physical properties at UV and optical wavelengths, while the new MIRI data at 15 µm further confirm its extreme nature, suggesting a strong similarity to other EROs previously studied in the Spitzer era.

This paper is structured as follows. In Section 2, we describe the data set used in this study, including HST and JWST observations, with particular focus on the MIRI data reduction.

<table><caption>**Table 1**
Definitions of Line Ratios Adopted in this Work</caption><thead><tr><th>Index</th><th>Definition</th></tr></thead><tbody><tr><td>R3</td><td>log<sub>10</sub>([O III]λ5007/Hβ)</td></tr><tr><td>O32</td><td>log<sub>10</sub>([O III]λλ4959, 5007/[O II]λλ3727, 3729)</td></tr><tr><td>R23</td><td>log<sub>10</sub>([O III]λλ4959, 5007 + [O II]λλ3727, 3729/Hβ)</td></tr><tr><td>Ne3O2</td><td>log<sub>10</sub>([Ne III]λ3869/[O II]λλ3727, 3729)</td></tr><tr><td>Ne3O2Hd</td><td>log<sub>10</sub>([Ne III]λ3869 + [O II]λλ3727, 3729)/Hβ</td></tr><tr><td>O3Hg</td><td>log<sub>10</sub>([O III]λ4363/Hγ)</td></tr><tr><td>O33</td><td>log<sub>10</sub>([O III]λ5007/[O II]λλ3727, 3729)</td></tr></tbody></table>

**Note.** We use the O32 definition from A. Calabrò et al. (2024a), which has a fixed offset of $log_{10}(1.3)$ between our values and the more common definition of O32 (≡ $[O III]λ5007/[O II]λλ3727, 3729$).

Section 3 presents a detailed analysis of the spectral properties of Virgil, based on NIRSpec data obtained by the OASIS team. Sections 4 and 5 present updated photometry and SED fitting, incorporating the newly acquired MIRI data from the MIDIS and PAHSPECS teams. In Section 6, we discuss the nature of Virgil, and in Section 7, we summarize our main conclusions.

Throughout this paper, we consider a cosmology with $H_0 = 70 \text{ km s}^{-1} \text{Mpc}^{-1}$, $\Omega_M = 0.3$, and $\Omega_\Lambda = 0.7$. All magnitudes are total and refer to the AB system (J. B. Oke & J. E. Gunn 1983). A P. Kroupa (2001) initial mass function (IMF) is assumed (0.1–100 $M_\odot$). Moreover, in Table 1, we report the line ratios that will be adopted throughout this work.

## 2. Data Set

In this work, we used data from HST and JWST, which are available in GOODS-S. In particular, we will briefly describe the data reduction for the MIDIS data and the NIRSpec/PRISM data collected by OASIS.

### 2.1. HST

Regarding the HST data, we used the ACS/WFC and WFC3/IR data from the Hubble Legacy Field (HLF) observations that cover GOODS-S. The HLF provides deep imaging in nine HST bands covering a wide range of wavelengths (0.4–1.6 µm), from the optical (ACS/WFC F435W, F606W, F775W, F814W, and F850LP filters) to the near-infrared (WFC3/IR F105W, F125W, F140W, and F160W filters). We refer the reader to K. E. Whitaker et al. (2019) for a more detailed description of this data set.

### 2.2. NIRCam

We made use of NIRCam data from JADES Data Release 2 (JADES DR2—PIDs: 1180 and 1210; PIs: D. Eisenstein and N. Luetzgendorf; D. J. Eisenstein et al. 2023a, 2023b), which include observations from the JWST Extragalactic Medium-band Survey (JEMS—PID: 1963; PIs: C. C. Williams, S. Tacchella, and M. Maseda; C. C. Williams et al. 2023) and the First Reionization Epoch Spectroscopically Complete Observations (FRESCO—PID: 1895; PI: P. Oesch; P. A. Oesch et al. 2023).

The JADES/NIRCam data allowed us to cover a wide range in wavelengths ($\approx 1–5$ µm) with the following bands: F090W,
F115W, F150W, F182M, F200W, F210M, F277W, F335M, F356W, F430M, F444W, F460M, and F480M. The estimated $5\sigma$ depth ranges from 30.5 to 30.9 mag (depths are calculated using 0.2" diameter circular apertures assuming point-source morphologies; see K. N. Hainline et al. 2024 for details).

## 2.3. MIRI

We utilized MIRI data from the MIDIS (PID: 6511; PI: Östlin), PAHSPECS (PID: 5279; PIs: Shivaei, Díaz-Santos, and Boogaard), and SMILES (PID: 1207; PI: Rieke) programs, which allowed us to cover a wide range of wavelengths (5.6–25.5 µm). Specifically, MIDIS provides observations at 5.6, 7.7, and 10 µm, PAHSPECS enables photometry at 15 µm, and SMILES samples 12.8, 18, 21, and 25.5 µm. For SMILES, we relied on publicly available images released by the SMILES team. In contrast, we independently reduced the MIDIS data, while the PAHSPECS team provided a 5" × 5" cutout centered on Virgil. Details on the MIRI/F1500W data reduction will be presented in I. Shivaei et al. (2025 in preparation), based on the procedures described in P. G. Pérez-González et al. (2024a) and G. Östlin et al. (2025).

To process the MIDIS data, we employed our custom MIRI pipeline (already calibrated in P. Rinaldi et al. 2023 and E. Iani et al. 2024, and further improved), built on the latest *jwst* pipeline (1.17.1) with *pmap_1321*. We used the default *stage 1* and *stage 2* steps of the *jwst* pipeline to produce calibrated data (*_cal.fits*). Additionally, we implemented custom steps for creating a master background, performing a superbackground subtraction and homogenization, and alignment to a common world coordinate system (WCS) based on the JADES DR2 catalog (D. J. Eisenstein et al. 2023a).

Our custom pipeline relies on a multitiered approach. First, it combines the*_cal.fits* files into an initial mosaic using a custom drizzle algorithm based on the official package provided by STScI.<fnref n="36" /> This initial mosaic serves as the baseline for creating a source mask to exclude sources in the individual*_cal.fits* files. To construct this mask, our pipeline uses SOURCE EXTRACTOR (E. Bertin & S. Arnouts 1996). A preliminary *master_background* is then built by reprojecting this preliminary source mask onto each exposure, stacking all exposures into a 3D array, and computing the *biweight location* along the *z*-axis. This method minimizes the influence of outliers and provides a robust estimate of the 2D background (i.e., *master_background*).

For each exposure, our pipeline iteratively estimates a scale factor (A) to optimize the subtraction of the *master_background*. The scale factor is recalculated at each iteration using clipped data to ensure convergence and minimize residuals between the corrected and original images. Once the *master_background* is subtracted, the pipeline filters each exposure along rows and columns to remove residual patterns while masking out the sources. This procedure is similar to the standard approach used to mitigate wisp effects in NIRCam data, as described by M. B. Bagley et al. (2024).

Then, the pipeline drizzles together all the single exposures to produce a new mosaic, from which it creates a *master_source_mask*. At this level, the pipeline makes use of both SOURCE EXTRACTOR and NOISECHISEL (M. Akhlaghi & T. Ichikawa 2017) to create the *master_source_mask*: SOURCE EXTRACTOR effectively identifies compact sources.

In contrast NOISECHISEL detects faint sources and extended regions around bright objects. The *master_source_mask* is then reprojected onto each single exposure to improve the *master_background* and superbackground subtraction and homogenization, following a methodology similar to the one outlined in M. B. Bagley et al. (2024).

After this step, the exposures are aligned to the JADES DR2 reference catalog using a modified version of the *tweak-reg step* within the *jwst* pipeline.

The final steps include running the *outlier step* and *resampling step*, producing the final mosaic. This final mosaic is then used to refine the *master_source_mask* and iterate the entire procedure for further optimization. Similar methodologies have already been widely used when reducing MIRI data (e.g., S. Alberts et al. 2024; G. Östlin et al. 2025; P. G. Pérez-González et al. 2024a).

Finally, we estimated the $5\sigma$ depth in our reduced MIRI images, which is 28.6, 28.3, 27.1, and 24.9 mag at, respectively, 5.6, 7.7, 10, and 15 µm ($r = 0.3$). The $5\sigma$ depth for the SMILES images at 12.8, 18, 21, and 25.5 µm is 24.4, 23.3, 22.8, and 20.8 mag (S. Alberts et al. 2024; $r = 0.4–0.7$). Depths are estimated using circular apertures assuming point-source morphologies.

## 2.4. NIRSpec/PRISM

The OASIS program consists of two ultradeep NIRSpec/MSA pointings, each with a total exposure time of 28 hr, using the PRISM/CLEAR configuration to obtain low-resolution ($R \approx 30–300$) spectroscopy over the 0.6–5.3 µm range. The observations are executed in multiple visits, with each pointing incorporating three nodding positions within the MSA to improve background subtraction, and multiple dithered observations to mitigate detector artifacts. The MSA configuration is designed to target a sample of 223 galaxies at $5 < z < 8$. The NIRSpec observations are performed using the NRSRAID6 readout pattern with three groups per integration and a total of six integrations per exposure.

To reduce the NIRSpec/PRISM data, we made use of the pipeline developed by the ESA NIRSpec Science Operations Team (SOT) and the Guaranteed Time Observations (GTO) NIRSpec teams, with optimizations aimed at improving background subtraction, rectification, 1D extraction, and spectral combination. While broadly consistent with the STScI pipeline used for generating archive products, these enhancements ensure more accurate spectral calibration, particularly for compact sources affected by spatial offsets within the slit. A wavelength correction is applied to mitigate this bias, which arises when sources are not perfectly centered along the dispersion axis. For a more detailed description of the NIRSpec data reduction, we refer the reader to A. J. Bunker et al. (2024) and F. D'Eugenio et al. (2025). A more detailed description of this data set will be presented in T. Looser & F. D'Eugenio et al. (2025, in preparation).

# 3. On the Nature of Virgil: Insights from NIRSpec/PRISM

In this section, we delve into the nature of Virgil based on NIRSpec/PRISM data. Due to the extended nature of Virgil in the NIRSpec wavelength range, we corrected the flux following the method outlined in S. Carniani et al. (2024b). We applied a first-order polynomial ($\alpha_1 + \alpha_0\lambda$) to model the wavelength-dependent slit losses and calibrated the flux to

<footnote n="36">³⁶ https://github.com/spacetelescope/drizzle.</footnote>
<table><caption>**Table 2**
Emission Line Fluxes and Derived Properties for Virgil</caption><thead><tr><th>Line</th><th>Value</th><th>Error</th></tr></thead><tbody><tr><td>N IV] λ1487</td><td>47.3</td><td>19.0</td></tr><tr><td>C IV λ1549</td><td>28.3</td><td>18.7</td></tr><tr><td>C III] λ1906</td><td>48.4</td><td>16.4</td></tr><tr><td>[O II] λλ3727,3729</td><td>67.4</td><td>5.4</td></tr><tr><td>[Ne III] λ3867</td><td>34.6</td><td>5.5</td></tr><tr><td>He I λ3889</td><td>20.9</td><td>5.5</td></tr><tr><td>[Ne III] λ3968</td><td>22.9</td><td>4.5</td></tr><tr><td>Hδ λ4102</td><td>20.5</td><td>4.2</td></tr><tr><td>Hγ λ4340</td><td>32.4</td><td>4.9</td></tr><tr><td>[O III] λ4363</td><td>16.5</td><td>4.7</td></tr><tr><td>Hβ λ4861</td><td>92.5</td><td>4.1</td></tr><tr><td>[O III] λ4959</td><td>206.0</td><td>5.2</td></tr><tr><td>[O III] λ5007</td><td>587.9</td><td>5.2</td></tr><tr><td>He I λ5877</td><td>15.5</td><td>4.2</td></tr><tr><td>Hα λ6563<sup>a</sup></td><td>304.9</td><td>6.6</td></tr><tr><td>Hα (narrow) λ6563</td><td>219.9</td><td>12.5</td></tr><tr><td>Hα (broad) λ6563</td><td>168.8</td><td>17.9</td></tr><tr><td colspan="3"><strong>Galaxy Properties</strong></td></tr><tr><td>FWHM (broad Hα) [km s<sup>-1</sup>]</td><td>1750</td><td>575</td></tr><tr><td>E(B - V)<sub>gas</sub> [mag]</td><td>0.24</td><td>0.07</td></tr><tr><td>E(B - V)<sub>stars</sub> [mag]</td><td>0.10</td><td>0.03</td></tr><tr><td>A<sub>V</sub> (Hα/Hβ) [mag]</td><td>0.65</td><td>0.18</td></tr><tr><td>12+log(O/H)<sup>b</sup></td><td>7.78</td><td>0.18</td></tr><tr><td>T<sub>e</sub>(O<sup>+</sup>) [K]</td><td>14,800</td><td>600</td></tr><tr><td>T<sub>e</sub>(O<sup>2+</sup>) [K]</td><td>18,300</td><td>3200</td></tr><tr><td>log<sub>10</sub>(U)<sup>b</sup></td><td>-2.03</td><td>0.08</td></tr><tr><td>n<sub>e</sub> <sup>c</sup> [cm<sup>-3</sup>]</td><td>&lt;500</td><td>...</td></tr><tr><td>EW<sub>0</sub>(Hα) [Å]</td><td>422</td><td>20</td></tr><tr><td>EW<sub>0</sub>([O III]λλ4959.5007) [Å]</td><td>1514</td><td>18</td></tr><tr><td>SFR(Hα) [M<sub>☉</sub> yr<sup>-1</sup>]</td><td>6.10</td><td>0.93</td></tr><tr><td>SFR(UV) [M<sub>☉</sub> yr<sup>-1</sup>]</td><td>2.96</td><td>1.02</td></tr><tr><td>Σ<sub>SFR(Hα)</sub> [M<sub>☉</sub> yr<sup>-1</sup> kpc<sup>-2</sup>]</td><td>5.25</td><td>0.80</td></tr><tr><td>B</td><td>0.32</td><td>0.16</td></tr><tr><td>M<sub>UV</sub> [mag]</td><td>-18.70</td><td>0.17</td></tr><tr><td>β</td><td>-1.76</td><td>0.18</td></tr><tr><td>f<sub>esc,LyC</sub></td><td>0.02</td><td>0.01</td></tr><tr><td>ξ<sub>ion</sub> [log<sub>10</sub>(Hz erg<sup>-1</sup>)]</td><td>24.91</td><td>0.18</td></tr><tr><td>z<sub>spec</sub></td><td>6.6379</td><td>0.0035</td></tr></tbody></table>

**Notes.** Observed fluxes and errors (not corrected for dust) are given in units of $10^{-20}$ erg s$^{-1}$ cm$^{-2}$. We also show the Hα fluxes for the narrow and broad components obtained by MSAEXP. Derived galaxy properties with associated uncertainties are provided.

a Refers to Hα corrected for the contribution from [N II]λλ6548, 6583.

b Refers to the estimated median quantities.

c Refers to upper limits (2σ).

match the NIRCam photometry. That is, we focus on the rest-frame UV and optical spectral regions. All the line fluxes and properties estimated in this section are summarized in Table 2.

## 3.1. Redshift Assessment of Virgil

We used MSAEXP (G. Brammer 2023; K. E. Heintz et al. 2025) to estimate the spectroscopic redshift. MSAEXP fits the data using a combination of Gaussian continuum splines and emission-line templates. We set nsplines = 11 as an input parameter and searched for the minimum reduced $\chi^2$ value within the redshift range $z = 6-7$, based on the prior spectroscopic estimate from VLT/MUSE (R. Bacon et al. 2023). As shown

in the 2D spectrum (Figure 1), a faint Lyα feature—possibly blended with N V λ1240—is visible, consistent with the detection reported by R. Bacon et al. (2023). The line appears spatially offset, suggesting it may originate from a distinct kinematic component, potentially explaining the slight redshift discrepancy with the MUSE measurement. Several other emission lines are detected, yielding an unambiguous redshift of $z_{spec} = 6.6379 \pm 0.0035$, which is slightly higher than the VLT/MUSE estimate<fnref n="37" /> but still consistent to the second decimal place.

## 3.2. Virgil's Physical Properties from Its UV/Optical Spectrum

### 3.2.1. Balmer Decrement and Dust Content

We estimated the color excess, $E(B - V)_{gas}$, by analyzing the Balmer decrement between Hβ and Hα. We measured an observed Balmer decrement of $3.70 \pm 0.19$. However, since we cannot resolve [N II]λλ6548, 6583 individually, we applied the correction from P. Anders & U. Fritze-v. Alvensleben (2003) for low-metallicity objects (10%–20% Z$_{⊙}$; see next section), where the corrected Hα flux is ≈91% of the observed (Hα + [N II]λλ6548, 6583). After this correction, the Balmer decrement is $3.30 \pm 0.16$. Assuming Case B recombination (D. E. Osterbrock & G. J. Ferland 2006), and adopting an electron temperature of $1.8 \times 10^4$ K (as derived in the next section) and an electron density of $10^2$ cm$^{-3}$, we use PYNEB (V. Luridiana et al. 2015) to compute an intrinsic Hα/Hβ ratio of 2.76. The measured ratio corresponds to $E(B - V)_{gas} = 0.24 \pm 0.07$ mag, assuming the Small Magellanic Cloud (SMC) reddening law (K. D. Gordon et al. 2003), which has been shown to reproduce well the dust attenuation in high-redshift galaxies (e.g., N. A. Reddy et al. 2018) and reddened quasars (e.g., P. F. Hopkins et al. 2004). This value leads to an attenuation of $A_V = 0.65 \pm 0.18$ mag, adopting $R_V = 2.74 \pm 0.13$ (K. D. Gordon et al. 2003). It further implies $E(B - V)_{stars} = 0.10 \pm 0.03$ mag.<fnref n="38" />

To further assess the dust attenuation, we also compared the Balmer decrement for Hβ with respect to Hγ and Hδ (expected theoretical ratios of 0.474 and 0.263, respectively, under Case B recombination; D. E. Osterbrock & G. J. Ferland 2006). The resulting $E(B - V)_{gas}$ values were consistent with a low level or zero attenuation.

If Virgil instead hosted a Type 1 AGN (see Section 3.4), the broad component of Hα must be excluded when computing the Balmer decrement. In this case, considering only the narrow component of Hα and treating the detected broad Hβ as a narrow line due to its low significance, we would obtain a lower Balmer decrement of $2.38 \pm 0.18$, consistent with a dust-free scenario. Alternatively, in this scenario, the Balmer decrement may be intrinsically lower, as observed in some $z \approx 4-7$ galaxies with density-bounded regions or in low-mass extreme emitters (L. Sandles et al. 2024; W. McClymont et al. 2025; C. Scarlata et al. 2024).

### 3.2.2. Gas-phase Metallicity and Ionization State

Both MSAEXP and our internal routine successfully deblended Hγ and [O III]λ4363. Additionally, we detect the [O II]λλ3727, 3729 doublet. These lines, along with [O III]

<footnote n="37">³⁷ The redshift from VLT/MUSE is based on a single line identification.</footnote>

<footnote n="38">³⁸ Assuming the differential dust attenuation, we adopt $E(B - V)_{stars} = 0.44 \times E(B - V)_{gas}$ (Y. Koyama et al. 2019; I. Shivaei et al. 2020).</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Left panel: 2D MSA/PRISM spectrum produced by MSAEXP. We optimally scaled the trace to highlight all of the significant line detections. We show the data in azur, while the uncertainty is highlighted in gray. Assuming the best-fit MSAEXP $z_{\text{spec}} = 6.6379 \pm 0.0035$, we show the positions of the emission lines as dashed vertical lines. Right panel: JWST/NIRCam RGB image ($2.5' \times 2.5'$) along with the slit positions showed in white. We highlight Virgil and the foreground LAE at $z \approx 4.77$ already studied in J. Matthee et al. (2022) with VLT/MUSE data.</figcaption></figure>

$\lambda\lambda4959, 5007$, enabled us to estimate the electron temperatures ($T_e$) of $O^+$ and $O^{2+}$ and the gas-phase metallicity (e.g., Y. I. Izotov et al. 2006; M. Curti et al. 2024; R. L. Sanders et al. 2024). For this purpose, we employed a nonparametric approach based on D. Langeroodi & J. Hjorth (2024a), utilizing the GENESIS-METALLICITY tool,<fnref n="39" /> which is built on PYNEB (V. Luridiana et al. 2015). We derived $T_e(O^+)$ = (1.48 ± 0.06) × 10⁴ K and $T_e(O^{2+})$ = (1.83 ± 0.32) × 10⁴ K (adopting the direct method). We highlight that $T_e(O^{2+})$ is higher than typical values observed in galactic H II regions, but is consistent with the temperatures measured in the most metal-poor galaxies (e.g., E. D. Skillman et al. 2013).

Furthermore, the GENESIS-METALLICITY tool offers the option to estimate gas-phase metallicity using either the direct ($T_e$) method or strong-line calibrations, depending on the available line fluxes. Both approaches yielded consistently low metallicities for Virgil, with 12 + log(O/H) = 7.70 ± 0.14 for the direct method, and 12 + log(O/H) = 7.91 ± 0.04 for the strong-line method. Based on the commonly accepted solar oxygen abundance (12 + log(O/H) = 8.69; M. Asplund et al. 2009), this gives us $Z/Z_\odot \approx 0.11-0.18$.

In the recent literature, other emission lines have been proposed to trace the gas-phase metallicity (mostly from UV), such as C IV] $\lambda$1548 and C III] $\lambda$1907 (see M. Mingozzi et al. 2022). Based on their Equation (7), which exploits the equivalent width (EW) of C III] $\lambda$1907, 09, we estimated that 12 + log(O/H) = 7.78 ± 0.21 (in agreement with our previous measurements). As shown in the mass-metallicity relation (MZR; Figure 2), Virgil is in good agreement with the relation recently reported by M. Curti et al. (2024) for galaxies at $z \gtrsim 6$ by leveraging the JADES sample as well as the recent literature from JWST.

We further investigate the ionization state of Virgil's gas by examining the ionization parameter, $\log_{10}(\mathcal{U})$, a key ISM property that reflects the intensity and hardness of the ionizing radiation, gas density, and the geometry of gas relative to

ionizing sources. The O32 index, a well-established tracer of $\log_{10}(\mathcal{U})$, yields a value of -2.03 ± 0.08 based on the empirical calibration from C. Papovich et al. (2022; their Equation (4)), which was derived from galaxies at 1.1 < z < 2.3 in the CLEAR survey (R. C. Simons et al. 2023). Notably, this value is consistent with the estimate from UV-based tracers identified by M. Mingozzi et al. (2022), who calibrated log(C III]/O III)], log(C IV/C III)], log(EW(O III))), and log(EW(C IV)) against O32. Applying their Equations (8) and (10), we find an average $\log_{10}(\mathcal{U}) = -1.99 \pm 0.15$, confirming a high level of ionization in agreement with the O32-based result. Alternatively, using the relation between $\log_{10}(\mathcal{U})$ and the Ne3O2 index from J. Witstok et al. (2021; their Equation (3)), we derive $\log_{10}(\mathcal{U}) = -2.33 \pm 0.07$. While slightly lower than the previous estimates, this value remains indicative of a highly ionized ISM, which is in line with the average galaxy population during the Epoch of Reionization (EoR).

We were also able to fit the [S II]$\lambda\lambda6717, 6731$ doublet, although with a low signal-to-noise ratio (SNR, $\approx 2\sigma$). Following the theoretical framework of L. J. Kewley et al. (2019), we estimated an upper limit on the electron density, $n_e \approx 500 \text{ cm}^{-3}$, considering its weak dependence on $T_e$ (X. Zhang 2024). However, the low spectral resolution of NIRSpec/PRISM<fnref n="40" /> prevents deblending the [O II] $\lambda\lambda3727, 3729$ doublet, limiting further constraints on $n_e$. Recently, M. W. Topping et al. (2025) analyzed deep JWST/NIRSpec spectra of SFGs at $z \approx 1.4-10$ as part of the AURORA survey (A. E. Shapley et al. 2025), focusing on the [S II]$\lambda\lambda6717, 6731$ doublet (detected up to $z \approx 6.8$), its evolution over cosmic time, and its correlation with various parameters (e.g., line ratios and star formation rate). Notably, our upper limit aligns with the average $n_e$ expected at $z \approx 6.5$, and it is consistent with the estimate obtained using their relation with the star formation rate (SFR; their Equation (2)) and Ne3O2 (their Equation (3)).

<footnote n="39">³⁹ The tool is available at https://github.com/langeroodi/genesis_metallicity/tree/main.</footnote>

<footnote n="40">⁴⁰ $R \approx 30-420$, depending on the observed wavelength.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** The stellar mass–gas-phase metallicity relation. Virgil is shown as a red filled circle. The stellar mass is taken from E. Iani et al. (2025) (median value), with an error bar reflecting the uncertainty in the estimate, consistent with the results obtained in this paper (see Section 5). Filled squares represent recent JWST measurements at $z \approx 3-10$ (collectively labeled as JWST literature), including galaxies from JADES (M. Curti et al. 2024), CEERS (K. Nakajima et al. 2023), and the SMACS J0723.3-7327 field (M. Curti et al. 2023). For comparison, we show MZR determinations at lower redshifts from M. Curti et al. (2020; Sloan Digital Sky Survey; SDSS) and R. L. Sanders et al. (2021; MOSDEF at $z \approx 2-3$). Additionally, we include the best fit of the low-mass end of the MZR at $z \approx 3$ from M. Li et al. (2023), based on JWST/NIRISS slitless spectroscopy. We also show the MZR determinations from the comprehensive study of M. Curti et al. (2024) as dashed lines, spanning $z = 3-10$ and split into $z < 6$ and $z > 6$. Finally, we also show theoretical predictions at high redshift from ASTRAEUS (G. Ucci et al. 2023), FIRSTLIGHT (I. Langan et al. 2020), the First Light and Reionisation Epoch Simulations (FLARES; C. C. Lovell et al. 2021), and TNG50 simulations (D. Nelson et al. 2019).</figcaption></figure>

As a final note, we report the detection of the N IV]λ1487 emission line, while N III]λ1750 is detected at a significance well below $2\sigma$. Additionally, MSAEXP identifies a possible N Vλ1240 feature, although it may be blended with Ly$\alpha$ due to the low resolution of NIRSpec/PRISM. The presence of N IV] λ1487 has been proposed as a signature of low-luminosity quasars (E. Glikman et al. 2007; A. J. Bunker et al. 2023; A. J. Cameron et al. 2023). However, we note that NIRSpec spectra often exhibit spurious oscillations in this wavelength range, and higher-resolution, deeper observations will be necessary to confirm the reality of these lines. For this reason, we do not discuss this further.

### 3.2.3. UV β slope, Burstiness, and Ionizing Photon Production Efficiency

We modified the official MSAEXP code to enable estimation of the UV β slope (in the rest-frame wavelength range $\lambda \approx 1500-2500$ Å) by fitting a linear relation ($\log_{10}(f_{\lambda}) = \beta \log_{10}\lambda + q$) using nonlinear least squares optimization via `curve_fit` (SCIPY). We also employed a Monte Carlo

approach, performing 1000 linear fits to spectral realizations with randomly perturbed fluxes and adopting the median and standard deviation of the resulting β values. We retrieved $\beta = -1.76 \pm 0.18$.

To investigate the potential role of Virgil during the EoR, we estimated its ionizing photon production efficiency, $\xi_{\text{ion}}$. Following the method outlined in P. Rinaldi et al. (2024), we first derived the Lyman continuum escape fraction ($f_{\text{esc, LyC}}$) from the observed UV β slope using the prescription of J. Chisholm et al. (2022). We then computed $\xi_{\text{ion}}$ based on its dependence on $(1 - f_{\text{esc, LyC}})^{-1}$ (Equation (2)) in P. Rinaldi et al. 2024). We find $f_{\text{esc, LyC}} = 0.02 \pm 0.01$ and $\log_{10}(\xi_{\text{ion}}/\text{Hz erg}^{-1}) = 24.91 \pm 0.18$.41 The latter is below the canonical value ($25.2 \pm 0.1$; B. E. Robertson et al. 2013), but still in agreement with the average high-z population at $z \approx 6-7$ (J. Matthee et al. 2017; P. Rinaldi et al. 2024; C. Simmonds et al. 2024a). This result is not a surprise, as the rest-frame equivalent width (EW$_{0}$) for H$\alpha$ is about 422 ± 20 Å, and very high $\xi_{\text{ion}}$ is often associated with very strong line emitters (i.e., large EW$_{0}$(H$\alpha$); see Figure 5 from P. Rinaldi et al. 2024).

Following the work presented in H. Atek et al. (2022) and R. Navarro-Carrera et al. (2024), high-z galaxies are expected to grow through bursty episodes of star formation, driven by large temporal fluctuations in their SFR—a trend recently suggested also by D. Langeroodi & J. Hjorth (2024b).

The NIRSpec/PRISM spectrum enabled us to estimate the star formation rate (SFR)42 directly from emission lines, specifically H$\alpha$, using the calibration from R. C. Kennicutt & N. J. Evans (2012) based on a Kroupa IMF. In particular, we adopt:

$$ \log_{10}\left(\frac{\text{SFR}}{\text{M}_{\odot}\,\text{yr}^{-1}}\right) = \log_{10}\left(L_{\text{H}\alpha}/(\text{erg s}^{-1})\right) - 41.64, \quad (1) $$

where the coefficient 41.64 corresponds to the calibration at subsolar metallicity (10%–20% $Z_{\odot}$; R. L. Theios et al. 2019).

We find that $\log_{10}(\text{SFR(H}\alpha)/(\text{M}_{\odot}\,\text{yr}^{-1})) = 0.79 \pm 0.07$. Based on the estimate from E. Iani et al. (2025) for the effective radius ($r_{\text{eff}} \approx 0.43$ pkpc), we retrieve $\Sigma_{\text{SFR(H}\alpha)} = \text{SFR(H}\alpha)/(2\pi r_{\text{eff}}^2) = 5.25 \pm 0.80 \text{ M}_{\odot} \, \text{yr}^{-1} \, \text{kpc}^{-2}$, which is in line with recent results from A. Calabro et al. (2024b) at $z \approx 6.5$.

The detection of H$\alpha$ allowed us to estimate another key parameter: the burstiness, defined as $B = \log_{10}(\text{SFR(H}\alpha)/\text{SFR(UV)})$ (H. Atek et al. 2022). We derived SFR(UV) at 1500 Å using the calibration from R. C. Kennicutt & N. J. Evans (2012), adapted for subsolar metallicity following R. L. Theios et al. (2019), obtaining $\log_{10}(\text{SFR(UV)}/\text{M}_{\odot}\,\text{yr}^{-1}) = 0.47 \pm 0.15$. This yields a burstiness parameter of $B = 0.32 \pm 0.16$, indicating that the source is currently forming stars at a rate more than 2 times higher than its recent past average.

The timescales of H$\alpha$ and UV emission provide further insights. Following the discussion presented in A. L. Faisst et al. (2019), E. Iani et al. (2024), and R. Navarro-Carrera et al. (2024),43 after a star formation episode, H$\alpha$ is initially

<footnote n="41">41 Although the systematic uncertainty reported by J. Chisholm et al. (2022) is on the order of 0.05, the lower bound of our estimate makes $f_{\text{esc, LyC}}$ consistent with 0%. Nonetheless, the resulting $\xi_{\text{ion}}$ would not differ significantly from the reported value.</footnote>

<footnote n="42">42 Both SFR(H$\alpha$) and SFR(UV) are corrected for dust attenuation.</footnote>

<footnote n="43">43 See Figure 12 in A. L. Faisst et al. (2019) based on different SFHs, Figure 15 in E. Iani et al. (2024) based on STARBURST99, and Figure 7 in R. Navarro-Carrera et al. (2024) based on FIRSTLIGHT simulations.</footnote>
enhanced due to newly formed O and B stars, while the UV luminosity is still rising. This results in SFR(Hα)/SFR(UV) > 1, which has been suggested as a signature of a bursty star formation phase (H. Atek et al. 2022). As massive stars rapidly evolve and die, Hα emission declines, while longer-lived UV-emitting stars sustain the UV luminosity, leading to SFR(Hα)/SFR(UV) < 1. As these UV-bright stars also fade, the ratio eventually returns to unity.

Our findings suggest that Virgil likely experienced a previous burst of star formation well before the most recent one (i.e., the one currently producing Hα), which explains its SFR(UV) and likely placed its *B* below unity. We now measure *B* = 0.32 ± 0.16, indicating that the source is either entering a new bursty phase or fading out of one. Notably, this is consistent with its position in the SFR-*M* plane (P. Rinaldi et al. 2022, 2025b), where it overlaps with the separation line for the starburst cloud proposed by K. I. Caputi et al. (2017, 2021). This suggests we may be observing Virgil at a transitional stage in its star formation history.

### 3.3. The ISM Conditions of Virgil: Excitation versus Ionization

Since the NIRSpec/PRISM spectrum provides access to [O II]λλ3727,3729, [Ne III]λ3869, Hδ, Hβ, and [O III] λλ4959,5007 at $z \approx 6.64$, we investigated the interstellar medium (ISM) conditions of Virgil by leveraging the following line ratios: O32, R23, Ne3O2, and Ne3O2Hd. To fully explore the complexity of the ISM in this source, we closely followed the methodology described in A. Calabrò et al. (2023) and A. Calabrò et al. (2024a), running PYCLOUDY (v0.9.15) with CLOUDY version 17.01 (G. J. Ferland et al. 2017). We adopted their setup (see A. Calabrò et al. 2024a, Section 2.3) to generate line predictions for direct comparison with our results.

Briefly, we modeled SFGs with a spherically symmetric, radiation-bounded shell of gas surrounding a population of young (O- and B-type) stars, adopting the incident radiation field from BPASS stellar population models (with an IMF extending up to 100 $M_*$ and a continuous star formation in the past 30 Myr; J. J. Eldridge et al. 2017). Following A. Calabrò et al. (2024a), we considered the metallicity range from 0.05 to 1 times solar (i.e., 0.05, 0.1, 0.15, 0.2, 0.3, 0.4, 0.5, 0.7, and 1), with the solar reference consistent with M. Asplund et al. (2009). In all cases, we derived predictions for four different ionization parameters $log_{10}(U) = -3, -2.5, -2, \text{and } -1.5$ and for three gas density values ($10^2, 10^3$, and $10^4$ cm$^{-3}$).

Given the low metallicity we inferred for this object, we decided to follow the approach adopted in A. Calabrò et al. (2024a) also for the dust depletion, where the metals are depleted in the beginning of our CLOUDY calculations, and we considered that this depletion is metallicity dependent, as discussed in A. Calabrò et al. (2023; see their Table 2).

Finally, we explored the density-bounded scenario for the nebula by varying the stopping criterion in CLOUDY from a Lyman continuum (LyC) optical depth = 10 (fully ionization-bounded case) to 0.1 (fully density-bounded case). This would correspond to an escape fraction of ionizing photons going from 0% to 100%. A density-bounded scenario typically enhances the flux of high ionization species relative to low-ionization species. That is, this would increase the O32 and Ne3O2 indices, mimicking the effect of a high ionization parameter. However, as previously found, Virgil has a very low $f_{\text{esc, LyC}}$ (≈2%), so the increase in O32 and Ne3O2 would

be very minimal (<0.05 dex), thus not impacting our conclusions. We refer the reader to A. Calabrò et al. (2023) and A. Calabrò et al. (2024a) for a more detailed description of the models adopted in this work.

As reported by A. Calabrò et al. (2024a), the ionization of the ISM increases with redshift, as seen in the R23–O32 diagram (Figure 3, left panel) from the local Universe (SDSS) to high redshifts (see redshift bins from $z \approx 4$ to $z \approx 12$ in A. Calabrò et al. 2024a). For comparison, we include recent high-redshift ($z \approx 4–12$) NIRSpec results from A. J. Cameron et al. (2023), S. Mascia et al. (2023), K. Nakajima et al. (2023), K. Boyett et al. (2024), A. Calabrò et al. (2024a, 2024b), and G. Roberts-Borsani et al. (2024), along with the low-redshift sample from the SDSS (D. G. York et al. 2000; G. Kauffmann et al. 2003). Additionally, we include metal-poor and high-ionization systems selected in the local Universe to approximate the properties of high-redshift galaxies (F. Bian et al. 2018).

Virgil lies in the upper region of this diagram, but it is not as extreme as, for instance, GHZ2 at $z \approx 12$, analyzed by A. Calabrò et al. (2024a) and first reported by J. A. Zavala et al. (2025). Overall, its behavior aligns with other high-z galaxies observed during the EoR. Comparing its position with photoionization models for SFGs, Virgil overlaps with a wide range of models, from high-ionization systems ($log_{10}(U) = -1.5$) with nearly solar metallicity to lower-ionization ($log_{10}(U) = -2.5$) metal-poor environments, spanning a variety of electron densities ($10^2–10^4$ cm$^{-3}$). This degeneracy complicates its interpretation. However, independent measurements (see Section 3.2.2) suggest $log_{10}(U) \approx -2$ and a metal-poor composition (≈10%–20% $Z_\odot$), helping to constrain its nature. Taken together, these results suggest that Virgil shares properties with the average SFG population during the EoR.

We also investigated the Ne3O2Hd–Ne3O2 diagram (Figure 3, right panel). As with R23–O32, we include SDSS galaxies from the local Universe and high-z analogs from F. Bian et al. (2018). Virgil falls within the region occupied by other high-z galaxies (e.g., K. E. Heintz et al. 2025) and lies along the evolutionary sequence reported by A. Calabrò et al. (2024b) at $z \approx 4–12$. Compared to R23–O32, this diagram exhibits less degeneracy with photoionization models. However, model predictions in this case would suggest a lower $log_{10}(U)$ for metal-poor solutions, but still within the range inferred for this object when using different tracers.

Without the prior assumption that Virgil may host an AGN, as proposed in E. Iani et al. (2025), these two diagrams, together with the previously derived physical properties, indicate that Virgil has a hard ionizing radiation field and ionization properties consistent with the average SFG population observed during the EoR. However, its exceptionally high equivalent width for [O III]λλ4959, 5007 (1514 ± 18 Å), combined with its low metallicity (10%–20% $Z_\odot$) and relatively low $M_*$ ($log_{10}(M_*/M_\odot) \approx 8.5$ from E. Iani et al. 2025), suggests that Virgil is a relatively low-$M_*$ extreme line emitter, consistent with the sample studied by K. Boyett et al. (2024), which spans similar redshifts.

Interestingly, following the discussion in J. E. Rhoads et al. (2023) and the spectral properties analyzed so far, Virgil could resemble a typical Green Pea (GP; C. N. Cardamone et al. 2009; A. E. Jaskot & M. S. Oey 2013; A. Henry et al. 2015; Y. I. Izotov et al. 2016; S. W. Brunker et al. 2020, among
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Left panel: O32–R23. Right panel: Ne3O2Hd–Ne3O2. Virgil is shown as a red filled circle, with higher opacity representing the case without dust correction and lower opacity indicating the dust-corrected one. In both panels, we show the photoionization models from CLOUDY with the following setup: $n_e$ (10$^2$, 10$^3$, and 10$^4$ cm$^{-3}$; increasing thickness), $log_{10}(U)$ (–3, –2.5, –2, –1.5; different shades of blue), and different metallicities (0.05 to solar; color coded in viridis using “x” symbols). For comparison, we show low-z star-forming galaxies from SDSS (D. G. York et al. 2000; G. Kauffmann et al. 2003) and local high-z analogs from F. Bian et al. (2018). To put our results in context, we also show recent observations at high-z from A. J. Bunker et al. (2023), A. J. Cameron et al. (2023), S. Mascia et al. (2023), K. Nakajima et al. (2023), K. Boyett et al. (2024), A. Calabrò et al. (2024a, 2024b), G. Roberts-Borsani et al. (2024), and K. E. Heintz et al. (2025) at $z \approx 4–12$.</figcaption></figure>

others) or a Blueberry (BB; e.g., A. J. Cameron et al. 2023; D. Langeroodi et al. 2023), which are often referred to as the best nearby analogs to high-redshift galaxies.

### 3.4. Does the Spectrum Indicate that Virgil is an AGN?

E. Iani et al. (2025) modeled Virgil across all available photometric bands (from HST and JWST), performing SED fitting using multiple codes, and concluded that it contains an AGN. We will reexamine this approach with improved data in Section 4. To investigate this possibility further, we exploit NIRSpec/PRISM data to search for spectral signatures indicative of an AGN.

#### 3.4.1. Searching for Broadening in the Balmer Lines

In E. Iani et al. (2025), we found that Virgil’s SED resembles that of an LRD, which often exhibits broad Balmer lines (e.g., L. J. Furtak et al. 2023; V. Kokorev et al. 2023; J. Matthee et al. 2024; P. Rinaldi et al. 2025a). Using MSAEXP, we modeled the NIRSpec/PRISM data with and without a broad component for the Balmer lines ($H\alpha$ and $H\beta$), closely examining the residuals. Figure 4 presents the results in two panels. To take into account the effects of the line-spread function, we additionally convolve our model with a Gaussian of variable resolution (Y. Isobe et al. 2023; A. de Graaff et al. 2024). On the left, we show the fit to the $H\alpha + [NII]$ complex without a broad-$H\alpha$ component, while on the right, we show the fit with the broad-$H\alpha$ component. In the latter case, the broad component is detected at the level of approximately 20σ and exhibits a full width at half-maximum (FWHM) of 1750 ± 575 km s$^{-1}$.

To assess the quality of these fits, we computed the Bayesian Information Criterion (BIC; A. R. Liddle 2007) for

both scenarios. The model, including the broad component, achieved the lowest BIC, indicating the best balance between goodness-of-fit and model complexity. However, the $\Delta$BIC is not large enough to definitively conclude that Virgil is a Type 1 AGN ($\Delta$BIC < 2). We also observe a tentative detection of a broad component in $H\beta$, but its significance falls below 2σ, leaving its presence uncertain. We want to highlight that these results may be affected by the source’s position within the slit (Figure 1, right panel): the slit only partially covers the galaxy, with the red knot located near the edge—or possibly even in the gap between shutters—potentially limiting our sensitivity to broad-line emission.

Nonetheless, under the assumption that Virgil exhibits $H\alpha$ with a broad component, we estimated the central black hole mass ($M_*$). Assuming the gas in the broad-line region (BLR) is virialized, $M_*$ was derived from the spectral properties of the $H\alpha$ BLR region using the calibration proposed by A. E. Reines et al. (2013):

$$ \log_{10}\left(\frac{M_{\mathrm{BH}}}{M_{\odot}}\right) = \alpha + \log_{10}(\varepsilon) + \beta \log_{10}\left(\frac{L_{H\alpha, \mathrm{broad}}}{1 \times 10^{42} \mathrm{erg} \mathrm{~s}^{-1}}\right) + \gamma \log_{10}\left(\frac{\mathrm{FWHM}_{\mathrm{broad}}}{1 \times 10^3 \mathrm{km} \mathrm{~s}^{-1}}\right), \qquad (2) $$

where $\alpha = 6.57$, $\beta = 0.47$, and $\gamma = 2.06$, and $\epsilon$ is the scaling factor that depends on the structure, kinematics, and orientation of the BLR. Different studies report $\epsilon$ ranging from 0.75 to 1.4 (e.g., A. E. Reines et al. 2013). Following A. E. Reines & M. Volonteri (2015), for our $M_*$ we adopted $\epsilon = 1.075 \pm 0.325$. We retrieved $\log_{10}(M_*/M_\odot) = 6.98$. This $M_*$ estimate should be regarded as a lower limit, given the absence of a dust
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Left panel: fit of the $H\alpha$ + [N II] complex without including a broad component. Right panel: fit of the $H\alpha$ + [N II] complex with adding a broad component. In this case, the FWHM is 1750 ± 575 km s⁻¹. The residuals appear to behave better when a broad component is included, although the Bayesian Information Criterion shows only a marginal improvement with the added complexity.</figcaption></figure>

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

attenuation correction for the $H\alpha$ broad component. Nonetheless, it agrees with the lower-limit estimate of E. Iani et al. (2025), derived from the $L_{bol}$ obtained through SED fitting (see their Section 4.2.2). We estimated the bolometric luminosity ($L_{bol}$) by using Equation (6) in J. Stern & A. Laor (2012), finding that $log_{10}(L_{bol}/(erg\ s^{-1})) = 44.03 \pm 0.43$. We also adopted Equation (1) from H. Netzer (2009) from which we derive $log_{10}(L_{bol}/(erg\ s^{-1})) = 45.17 \pm 0.02$, which likely represents an upper limit. Also in this case, these estimates are in agreement with E. Iani et al. (2025). In Figure 5 we show our measurements of $M_*$ for Virgil in the context of the recent literature.

Based on the derived $M_*$ from the $H\alpha$ BLR, we also derived the Eddington luminosity. We found that $log_{10}(L_{Edd}/(erg\ s^{-1})) = 4\pi GM_*m_p c/\sigma_T = 45.07 \pm 0.35$, where $G$ is the gravitational constant, $m_p$ the proton mass, $c$ the speed of light, and $\sigma_T$ the Thomson scattering cross section. Finally, we estimated the Eddington ratio as $\lambda_{Edd} = L_{bol}/L_{Edd}$. Depending on the choice of $L_{bol}$, the Eddington ratio ($\lambda_{Edd}$) ranges from ≈0.1 to 1.23, consistent with the typical values found for LRDs (e.g., A. Nobori-guchi et al. 2023). We show our results in Figure 5 and summarize these measurements in Table 3.

### 3.4.2. Searching for AGN Signatures with Line Ratio Diagnostics

Disentangling SFGs from AGNs at high redshift is particularly challenging. Classical diagnostics, such as the Baldwin, Phillips & Terlevich (BPT; J. A. Baldwin et al. 1981) diagram, were calibrated for low-redshift galaxies and strongly depend on metallicity, making them less reliable in the early Universe. At high redshift, galaxies are typically more metal-poor (e.g., M. W. Topping et al. 2024), with younger stellar populations and higher ionization parameters. These conditions shift SFGs toward the AGN locus in BPT-like diagrams, while AGN Narrow Line Regions (NLRs) often overlap with the SFG sequence. As a result, traditional line ratio diagnostics

become increasingly ambiguous, as also predicted by theoretical models (K. Nakajima & R. Maiolino 2022).

Although we are not able to resolve [N II]$\lambda$6548 from $H\alpha$, and both [O I]$\lambda$6302 and the [S II]$\lambda$6717, 6731 doublet have SNR ≲2, it is important to stress that classical BPT-like diagrams may not reliably determine whether Virgil hosts an AGN. Even when these lines are available, recent studies have shown that AGNs can occupy the same region as SFGs in BPT diagrams. For example, the $z \approx 5.55$ galaxy GS_3073 shows clear AGN signatures through broad $H\alpha$ and $H\beta$ emission, yet it lies within the SFG region in the BPT diagram (H. Übler et al. 2023). Like Virgil, GS_3073 is extremely metal-poor (12 + log(O/H) = 8.02), underscoring how low metallicity can blur the distinction between star formation- and AGN-dominated galaxies at high redshift. Similar conclusions have been reported by J. Scholtz et al. (2025) and D. D. Kocevski et al. (2023), in agreement with photoionization models from K. Nakajima & R. Maiolino (2022), which predict significant overlap between low-metallicity SFGs and AGNs in classical diagnostics.

Since the launch of JWST, the limitations of classical diagnostics have led to renewed efforts to develop alternative line ratio diagrams optimized for high-redshift galaxies (e.g., M. Hirschmann et al. 2023; G. Mazzolari et al. 2024; A. E. Shapley et al. 2025; B. E. Backhaus et al. 2025). In the following, we employ a set of diagnostic diagrams that have been extensively used to distinguish AGNs from SFGs at high redshifts.

**The “OHNO” diagram.** One widely used high-redshift diagnostic is the “OHNO” diagram (L. Trouille et al. 2011; G. R. Zeimann et al. 2015; B. E. Backhaus et al. 2022, 2023; N. J. Cleri et al. 2023; J. R. Trump et al. 2023; L. M. Feuillet et al. 2024), which is based on Ne3O2 and R3 and has proven to be a robust ionization diagnostic for high-redshift galaxies. Ne3O2 is particularly useful as it involves emission lines with similar ionization energies and closely spaced wavelengths, minimizing dust attenuation effects. As with the R23–O32 and
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Left panel: the derived *M*, as a function of $L_{bol}$, assuming Virgil has a broad component in $H\alpha$. We also show other recent findings from the literature (with some of them selected as LRDs; P. Rinaldi et al. 2025a), Y. Matsuoka et al. (2019), L. J. Furtak et al. (2023), Y. Harikane et al. (2023), D. D. Kocevski et al. (2023), V. Kokorev et al. (2023), R. L. Larson et al. (2023), H. Übler et al. (2023), and R. Maiolino et al. (2024). The dashed lines represent bolometric luminosities corresponding to Eddington ratios of $L_{bol}/L_{Edd} = 0.01, 0.1$, and 1.0. Right panel: the derived *M*, as a function of $M_*$ (stellar mass for Virgil comes from E. Iani et al. 2025, and its error bar reflects the overall uncertainty in fitting this object). For comparison, we show the bold black dashed line indicating the best fit to $z \approx 0$ AGN samples (J. E. Greene et al. 2016). The trend at higher redshifts is based on the recent analysis presented in F. Pacucci et al. (2023).</figcaption></figure>

<table><caption>**Table 3**
Black Hole Properties</caption><thead><tr><th>Property</th><th>Value</th></tr></thead><tbody><tr><td>FWHM<sup>a</sup> [km s<sup>-1</sup>]</td><td>1750 ± 575</td></tr><tr><td>log<sub>10</sub>(M<sub>*</sub>/M<sub>☉</sub>)<sup>b</sup></td><td>6.98</td></tr><tr><td>log<sub>10</sub>(L<sub>bol</sub>/(erg s<sup>-1</sup>))<sup>c</sup></td><td>44.03 ± 0.43 – 45.17 ± 0.02</td></tr><tr><td>log<sub>10</sub>(L<sub>Edd</sub>/(erg s<sup>-1</sup>))</td><td>45.07 ± 0.35</td></tr><tr><td>λ<sub>Edd</sub></td><td>0.1–1.23</td></tr></tbody></table>

**Notes.** These measurements rely on the assumption that the broadening of $H\alpha$ is real. Nonetheless, these values are in agreement with previous estimates from photometry presented in E. Iani et al. (2025).

<sup>a</sup> Refers to the broad component of $H\alpha$.

<sup>b</sup> Indicates that this mass estimate should be considered a lower limit, as no reliable dust correction can be applied to the assumed broad component.

<sup>c</sup> Refers to the two values estimated by using both J. Stern & A. Laor (2012) and H. Netzer (2009) formulas.

Ne3O2Hd–Ne3O2 diagrams, we also included photoionization models, following the configuration adopted by A. Calabrò et al. (2024a). For AGN models, we adopted the default CLOUDY prescription, assuming a multicomponent power-law continuum with a “blue bump” temperature of $10^6$ K and spectral energy indices of $\alpha_{UV} = -0.5$, $\alpha_X = -1.35$, and $\alpha_{OX} = -1.4$ (B. A. Groves et al. 2004) for the UV, X-ray, and optical-to-X-ray ranges, respectively. We used the same grid of metallicities (0.05 to 1 solar) as for the models for SFGs. We refer the reader to A. Calabrò et al. (2023) and A. Calabrò et al. (2024a) for a more detailed discussion of these models.

It is worth noting that the OHNO diagram may be affected by its dependence on metallicity (R. Tripodi et al. 2024), where metal-poor objects, such as Virgil, can appear as outliers (J. Scholtz et al. 2025).

In Figure 6, we show the “OHNO” diagram with Virgil as well as the recent literature on other high-$z$ sources (M. Killi et al. 2024; D. D. Kocevski et al. 2023; V. Kokorev et al. 2023; R. L. Larson et al. 2023; J. R. Trump et al. 2023; J. Scholtz et al. 2025; H. Übler et al. 2024b; A. Calabrò et al. 2024a; P. Rinaldi et al. 2025a; with some of them identified as Type 1 and Type 2 AGNs) and, for comparison, the sample of SFGs and AGNs at $z \approx 0$ from D. G. York et al. (2000) and G. Kauffmann et al. (2003).

Virgil exhibits moderately high Ne3O2 and R3 ratios,<sup>44</sup> and it lies well above the separation line proposed by B. E. Backhaus et al. (2022) at $z \approx 1.5$. Thus, it is consistent with AGN model predictions at subsolar metallicities ($Z_{gas}/Z_{☉} \lesssim 0.2$), moderate ionization parameters ($\log_{10}(\mathcal{U}) \approx -3$), and $n_e \approx 10^3–10^4 \text{ cm}^{-3}$. However, Virgil is also consistent with star formation models, overlapping with the parameter space expected for H II regions with high $n_e (\approx 10^3–10^4 \text{ cm}^{-3})$, high ionization parameter ($\log_{10}(\mathcal{U}) \approx -2$), and a wide range of metallicities.

Interestingly, Virgil occupies a region of the parameter space consistently populated by LRDs (some of which exhibit clear broadening in their Balmer lines; P. Rinaldi et al. 2025a), as well as other recently discovered AGNs (including Type 2; e.g., D. D. Kocevski et al. 2023; V. Kokorev et al. 2023; R. L. Larson et al. 2023; J. Scholtz et al. 2025; H. Übler et al. 2024b). However, as shown in A. Calabrò et al. (2024a) and supported by CLOUDY model predictions, this region of parameter space is also consistently populated by SFGs with a different set of properties.

We then utilized the recent results of B. E. Backhaus et al. (2022, 2025) to extrapolate the separation line up to $z \approx 6.5$.

<footnote n="44"><sup>44</sup> We show our measurements before and after applying dust corrections.</footnote>
<figure><img src="image_9.png" /></figure>

**Figure 6.** The “OHNO” diagram, displaying the following ratios: [O III]λ5007/H$\dot{\beta}$ vs. [Ne III]λ3870/[O II]3727, 3728. Virgil is shown as a red filled circle, with higher opacity representing the case without dust correction and lower opacity indicating the dust-corrected one. We also show the recent literature at $z \gtrsim 4$ from M. Killi et al. (2024), D. D. Kocevski et al. (2023), V. Kokorev et al. (2023), R. L. Larson et al. (2023), J. R. Trump et al. (2023), H. Übler et al. (2024b), A. Calabrò et al. (2024a), and P. Rinaldi et al. (2025a). In particular, sources from P. Rinaldi et al. (2025a) have been photometrically selected as LRDs, with some of them clearly showing broad components in their Balmer lines from NIRSpec data. For comparison, we also show SDSS (SFGs and AGNs) at low $z$ ($z \approx 0$) from D. G. York et al. (2000) and G. Kauffmann et al. (2003). The separation line comes from B. E. Backhaus et al. (2022) at $z \approx 1.5$. We also show an extrapolation at $z \approx 6.5$ (increased thickness), based on B. E. Backhaus et al. (2024). For comparison, theoretical predictions (for both SFGs and AGNs; dashed-dotted and dashed, respectively) are shown with the following settings: $n_e$ (10$^2$, 10$^3$, and 10$^4$ cm$^{-3}$; increasing thickness), $\log_{10}(\mathcal{U})$ (–3, –2.5, –2, –1.5; different shades of blue), and different metallicities (0.05 to solar; color coded in viridis using “x” symbols).

(thick gray line), relying on the proposed redshift evolution of R3 and N3O2. In this scenario, Virgil, along with sources from M. Killi et al. (2024), J. Scholtz et al. (2025), and J. R. Trump et al. (2023), as well as LRDs from P. Rinaldi et al. (2025a), and the BL AGN from V. Kokorev et al. (2024), would fall below the extrapolated separation line, placing Virgil within the SFG region—or at least blurring the net distinction that would otherwise arise using the relation at $z \approx 1.5$.

We emphasize that this is only a test to explore the implications of shifting the separation line to higher redshifts, based on the evolutionary trends proposed by B. E. Backhaus et al. (2022, 2025), and should not be interpreted as a definitive redefinition. Nonetheless, applying diagnostic diagrams calibrated at lower redshifts complicates the classification of Virgil and could provide misleading results in general.

Considering the degeneracies revealed by CLOUDY models in the region of the parameter space occupied by Virgil and the possible redshift evolution of R3 and Ne3O2, both the SFG and AGN scenarios (or a mixed nature) remain plausible for this source.

*The [O III]λ4363 diagrams.* G. Mazzolari et al. (2024) proposed new separation lines for three diagnostic diagrams based on [O III]λ4363, which has been linked to AGN activity (J. Brinchmann 2023) with higher ISM temperatures driven by AGN ionizing radiation (H. Übler et al. 2024b):

1. O32 versus O3Hg

2. Ne3O2 versus O3Hg

3. O33 versus O3Hg

From a theoretical point of view, enhanced [O III]λ4363 emission could hint at the presence of an AGN, as the ionizing photons produced by an AGN are typically more energetic than those from star formation, thus leading to significantly more efficient gas heating. Notably, it has been proposed that strong [O III]λ4363 emitting regions may coincide with high-ionization nuclear emission-line regions (HINERs; L. Binette 1985), and that an enhanced [O III]λ4363/[O III]λ5007 ratio ($R_{[O III]}$) would imply a strong connection with AGN activity. For Virgil we find that $R_{[O III]} = 0.032 \pm 0.01$, which is similar to what has been observed in other Seyfert 2 objects:
ESO 138-G01, Mrk 1210, and NGC 4507, for which L. Binette et al. (2024) concluded that a high-density NLR component is likely present. In particular, the combination of such high $R_{[O III]}$ and relatively high R3 would imply that their plasma is much hotter (≈18,000–20,000 K), possibly as a result of fast shocks (e.g., L. Binette et al. 1985). Following this line, T. Nagao et al. (2001) found, for a sample of low-z objects, that a high $R_{[O III]}$ ratio correlates with objects exhibiting hotter mid-IR (MIR) colors.

Below, we briefly overview each diagnostic diagram (already presented in G. Mazzolari et al. 2024).

We want to highlight that dust reddening can significantly affect the line ratios, particularly for the O32–O3Hg diagram, due to the wavelength separation between the adopted lines. This effect can cause AGNs to shift toward the SFG+AGN locus when there is substantial reddening. Nonetheless, the inverse scenario does not hold, as SFGs would not populate the AGN locus under any circumstances. The Ne3O2–O3Hg diagram, instead, is very similar to the “OHNO” diagram, but it is more stable as the separation between SFG and AGN is mainly based on a different gas temperature.

For O32–O3Hg, G. Mazzolari et al. (2024) showed that high-z sources are generally distributed toward the upper region of the diagram compared to local samples (SDSS, contour area). Moreover, they point out that a nonnegligible fraction of high-z sources not classified as AGNs falls within the region occupied by local AGNs. In contrast, some high-z AGNs overlap with the area covered by their local analogs. Following the discussion in G. Mazzolari et al. (2024) and their photoionization models, SFGs should occupy a well-constrained region of the O32–O3Hg diagram, with their upper boundary closely matching the distribution of SDSS SFGs and local analogs. In contrast, AGNs are expected to populate both the SFG region and the area occupied by local AGN samples. They show that SFG models do not extend into the AGN-dominated region. However, AGN models can reach the SFG region for some local and high-z AGNs. They also conclude that these trends remain unchanged even when considering the entire grid of SFG models from A. Feltre et al. (2016), reinforcing the idea that AGN-driven photoionization can produce significantly higher $[O III]\lambda4363/H\gamma$ ratios than ionization from hot stars. This is likely due to the higher-energy ionizing photons in AGNs, which heat the gas more efficiently at a given ionization parameter.

As with the O32–O3Hg diagram, normal SFGs and local analogs primarily occupy the lower region of the Ne3O2–O3Hg parameter space, while AGNs are now distributed across a broader area, including the upper-left region, which is not reached by any of the SFG models that G. Mazzolari et al. (2024) tested or by observed SFG samples. Notably, they also find that sources classified as AGNs in the O32–O3Hg diagram remain above the SFG distribution in this diagnostic. Likewise, sources not explicitly identified as AGNs but found within the AGN region of the first diagnostic continue to occupy the same region in the Ne3O2–O3Hg diagram.

The third diagnostic compares O33 and O3Hg. The overlap between local SFGs and AGNs is significantly larger than in the previous diagnostics. However, a distinct region characterized by high $[O III]\lambda4363/H\gamma$ and high $[O III]\lambda5007/[O III]\lambda4363$ remains populated exclusively by local AGNs. As noted by G. Mazzolari et al. (2024), very few high-z galaxies fall within the AGN-only region, with most shifting toward lower

$[O III]\lambda5007/[O III]\lambda4363$ values. They also show that current photoionization models do not cover this region.

According to G. Mazzolari et al. (2024), the distribution of SFG models closely follows that of local galaxies, particularly at the boundary with the AGN-only region. In contrast, AGN models extend across both the SFG and AGN regions. They also noted that when considering the full parameter grid of A. Feltre et al. (2016), SFG and AGN models would overlap entirely.

In the case of Virgil,<fnref n="45" /> we find that it consistently lies above or on the separation lines proposed by G. Mazzolari et al. (2024) (Figure 7), overlapping with AGNs and LRDs identified in recent studies based on JWST (e.g., V. Kokorev et al. 2023; K. Nakajima et al. 2023; J. Scholtz et al. 2025; H. Übler et al. 2023; I. Juodžbalis et al. 2024; P. Rinaldi et al. 2025a).

Nonetheless, one should note that $[O III]\lambda4363$ could be contaminated by $[Fe II]\lambda4360$. However, M. Curti et al. (2017) showed that such contamination primarily occurs in galaxies with an enriched ISM, whereas Virgil is metal-poor. Additionally, G. Mazzolari et al. (2024), by using CLOUDY models over a wide range of metallicities ($-2 \le \log(Z/Z_\odot) \le 0.5$) and ionization parameters ($-4 \le \log_{10}(\mathcal{U}) \le -1$), found that the $[Fe II]\lambda4288/[Fe II]\lambda4360$ ratio is roughly constant (≈1.25). Therefore, they searched for $[Fe II]\lambda4288$ in their sample populating the “AGN-only” region and found no evidence of its presence, concluding that any possible contamination from $[Fe II]\lambda4360$ can be considered negligible.

As with the “OHNO” diagram, we also explored the potential redshift evolution of key line ratio diagnostics. Specifically, we examined the Ne3O2 versus O3Hg plane (Figure 7, middle panel), based on the empirical evolution recently proposed by B. E. Backhaus et al. (2025). We display their separation lines at $z \approx 0$ and $z \approx 6.5$ (gray lines, with increased thickness). Under this redshift evolution, Virgil would fall in the SFG/AGN region, suggesting a mixed nature.

We want to highlight that B. E. Backhaus et al. (2025), based on photoionization models, found that only AGN-dominated sources populate the region above the separation line proposed by G. Mazzolari et al. (2024) in this diagram. Indeed, taken at face value, these line ratio diagnostics would clearly place Virgil in the AGN regime along with other claimed AGNs at similar redshifts. A similar argument can be made for the O32 and O33 diagnostics shown in the other two panels.

Altogether, this again underscores the complexity of distinguishing AGNs from SFGs at high redshift—especially when accounting for redshift evolution—as illustrated by Virgil’s position in a region indicative of a mixed nature.

Being in the “AGN” zone in the adopted diagnostic diagrams is a conditio sine qua non for being an AGN, as shown in Figures 6 and 7, although some galaxies in these regions have substantial star formation. If interpreted with no redshift correction, the two sets of diagnostic diagrams analyzed in this section would indicate that Virgil may host an AGN: (1) its position in the “OHNO” diagram aligns with regions typically occupied by AGNs and LRDs; while (2) in the $[O III]\lambda4363$ based diagrams from G. Mazzolari et al. (2024), Virgil occupies a region predominantly populated by recently discovered AGNs at high $z$. However, this classification strongly depends on the assumption that the criteria to

<footnote n="45">⁴⁵ For Virgil, $[O III]\lambda4363$ and $H\gamma$ were successfully deblended, with both lines detected at SNR > 5.</footnote>
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Left panel: O32 vs. O3Hg. Middle panel: Ne3O2 vs. O3Hg. Right panel: O33 vs. O3Hg. Virgil is shown as a red filled circle, with higher opacity representing the case without dust correction and lower opacity indicating the dust-corrected one. The dashed black lines indicate the separation criteria proposed by G. Mazzolari et al. (2024) to distinguish AGNs from SFGs. The contour areas correspond to the SDSS sample in the local Universe (K. N. Abazajian et al. 2009) for comparison. Filled circles represent pre-JWST literature data for SFGs and AGNs (mostly at low redshift; C. K. Seyfert 1943; Y. I. Izotov et al. 2006; D. A. Berg et al. 2012; R. Amorín et al. 2015; M. Perna et al. 2017; H. Yang et al. 2017a, 2017b; Y. I. Izotov et al. 2018; O. L. Dors et al. 2020; M. Armah et al. 2021; S. A. Pustilnik et al. 2021). For comparison with recent JWST-based studies, we include results from Y. Harikane et al. (2023), K. Nakajima et al. (2023), J. Scholtz et al. (2025), H. Übler et al. (2024b), I. Juodžbalis et al. (2024), and P. Rinaldi et al. (2025a) (with some of them having both SFGs and AGNs). In particular, the sample from P. Rinaldi et al. (2025a) consists of photometrically selected LRDs in GOODS-S, some exhibiting clear broadening in the Balmer lines. The gray dashed line in the middle panel represents the separation criterion at $z \gtrsim 4$ recently proposed by B. E. Backhaus et al. (2025) to distinguish BLAGNs from non-BLAGNs.</figcaption></figure>

separate SFGs and AGNs remain unchanged with redshift. If the redshift evolution of the relevant line ratios is considered—particularly in light of recent findings (e.g., B. E. Backhaus et al. 2025)—the classification of Virgil becomes more ambiguous, potentially suggesting a mixed nature or even a galaxy dominated by star formation.

In summary, we have conducted a detailed analysis of Virgil based on its emission-line spectrum. At first glance, Virgil appears to be a typical galaxy for its redshift, with low-to-moderate attenuation in its emission lines and low but representative metallicity. While its emission-line ratios could suggest AGN activity, this classification becomes ambiguous when accounting for the redshift evolution of diagnostic diagrams. Had it not been identified through its pronounced spectral steepening in the MIRI bands (E. Iani et al. 2025), Virgil could easily be misclassified as an ordinary high-$z$ SFG, sharing many properties with the bulk of the galaxy population during the EoR.

# 4. Photometry and Spectral Energy Distribution of Virgil

The nature of Virgil remains ambiguous when relying solely on data below 5 µm, even when rest-frame optical emission lines are considered, but changes significantly with the addition of longer-wavelength observations. We leveraged these observations to perform updated photometry and SED fitting since deeper MIRI data are now available at 7.7, 10, and 15 µm from the MIDIS and PAHSPECS teams. We build on the SED fitting presented in E. Iani et al. (2025), now incorporating the latest MIRI data. As the derived physical

properties remain essentially unchanged, we refer the reader to E. Iani et al. (2025) for a detailed description of the modeling and results. Here, we focus on a more interpretative discussion in light of the spectral analysis presented in this work.

## 4.1. Forward Modeling and Photometry of Virgil

We performed a set of customized photometric measurements using the PHOTUTILS package (L. Bradley et al. 2022). Object centroids were computed using the “windowed positions” method implemented in Source Extractor (E. Bertin & S. Arnouts 1996). Depending on the filter, we employed a combination of Kron apertures (R. G. Kron 1980) with a Kron parameter of 2.5 and circular apertures. This approach was crucial as Virgil appears elongated, particularly at shorter wavelengths (NIRCam short channel), as noted in E. Iani et al. (2025). In Figure 8, we show Virgil’s color F444W–F1000W (2.50 ± 0.02 mag) with respect to the SMILES sample in GOODS-S as well as a sample of LRDs at $z \approx 4–8$ already reported in V. Kokorev et al. (2024), P. G. Pérez-González et al. (2024b), and P. Rinaldi et al. (2025a). Interestingly, at 15 µm, Virgil shows F444W–F1500W = 2.87 ± 0.04 mag, which highlights how extreme this source is at the MIRI wavelengths. An overview of Virgil’s appearance is presented in Figure 9.

For circular aperture photometry, we adopted filter-dependent radii, with apertures as large as 0.4–0.5 in the MIRI bands. In each case, we measured the background in a close region around Virgil and calculated the noise using non-adjacent pixels (5 pixels apart) to take into account noise
<figure><img src="image_11.png" /><figcaption>**Figure 8.** The color–magnitude diagram shows F444W along x-axis and F444W–F1000W along y-axis, including Virgil, data from SMILES and JADES (hexagons), and a sample of LRDs in GOODS-S at $z \approx 4–8$ previously reported by V. Kokorev et al. (2024), P. G. Pérez-González et al. (2024b), and P. Rinaldi et al. (2025a) with MIRI detection in the SMILES catalog in GOODS-S.</figcaption></figure>

correlation, as explained in P. G. Pérez-González et al. (2023). We corrected the photometry to total flux using aperture corrections from JADES (for HST and NIRCam) and SMILES (for MIRI). We validated our new measurements by comparing them with those reported in E. Iani et al. (2025), finding good agreement.

We used FORCEPHO to model and subtract the nearby LAE ($z_{\text{spec}} \approx 4.77$) to avoid light contamination during the photometric measurements. The LAE and Virgil were modeled using J. L. Sérsic (1968) profiles, with one Sérsic profile representing the extended flux of the LAE and four additional compact Sérsic profiles modeling its stellar clumps. We observed that Virgil consists of two knots, one of which likely drives the strong upturn in the SED at MIRI wavelengths, leading us to model these knots separately with two Sérsic profiles. We utilized FORCEPHO to sample all model parameters using NIRCam and HST images, employing the Hamiltonian Monte Carlo Markov Chain method. Each galaxy was then forward-modeled in the MIRI images by fitting the flux of each component using nonnegative least squares regression, incorporating effective point-spread-function (ePSF) models provided by M. Libralato et al. (2024). Figure 10 presents the LAE-subtracted image using the best-estimated model, with minimal residuals at the LAE’s location, demonstrating the quality of our modeling. Although we obtained a model flux of Virgil during the process, we do not adopt it as fiducial, as the complex morphology of Virgil is not well represented by two Sérsic profiles. Nevertheless, we find good agreement between the model flux and the photometry results from the LAE-subtracted images.

## 4.2. Virgil is a Little Red Dot

First, we examined the photometric criteria for LRDs. E. Iani et al. (2025) demonstrated that Virgil satisfies most of the commonly used color criteria for LRDs, except for a

relatively blue F277W–F444W (V. Kokorev et al. 2024; P. Rinaldi et al. 2025a), which explains why it was not selected in recent LRD studies in GOODS-S (V. Kokorev et al. 2024; P. G. Pérez-González et al. 2024a; P. Rinaldi et al. 2025a). To further characterize Virgil, we adopted the methodology of D. D. Kocevski et al. (2025) and estimated the continuum slope ($f_{\lambda} = \lambda^{-\beta}$) between two continuum bands to identify a V-shaped SED. Specifically, following R. Lin et al. (2025), we applied the criteria $\beta_{UV} < -0.37$ and $\beta_{\text{opt}} > 0$, where $\beta = -0.4(m_1 - m_2)/\log_{10}(\lambda_1/\lambda_2)$ $-2$. Using F430M and F444W, we obtained $\beta_{\text{opt}} \approx 0.87$, while for $\beta_{UV}$, we found $\beta_{UV} \approx -0.68$ by relying on F115W and F150W. This places Virgil within the locus of V-shaped SED galaxies, reinforcing its classification as an LRD. These findings underscore the complexity of LRD selection (K. N. Hainline et al. 2025) and suggest that refining color criteria, potentially incorporating MIRI photometry when available (e.g., G. Barro et al. 2024b), may be necessary to select Virgil-like sources during the EoR. Finally, as highlighted in E. Iani et al. (2025), we found that Virgil satisfies the compactness criterion, showing $c_{F444W}(\equiv F(0.5)/F(0.25)) \approx 1.5$.

# 5. Stellar Properties of Virgil

## 5.1. SED Fitting

We followed the approach adopted in E. Iani et al. (2025) and employed multiple SED-fitting codes to estimate the stellar properties of Virgil. While E. Iani et al. (2025) present a detailed SED-fitting analysis, we focused on deriving key stellar parameters to contextualize our spectral results and investigate the nature of this object. Specifically, we use BAGPIPES (A. C. Carnall et al. 2019), CIGALE (M. Boquien et al. 2019), PROSPECTOR (B. D. Johnson et al. 2021), and SYNTHESIZER-AGN (P. G. Pérez-González et al. 2003, 2008, 2024a). As demonstrated and largely discussed in E. Iani et al. (2025), Virgil is well described only when an AGN component is included, except for SYNTHESIZER-AGN, which tends to prefer a two-population solution. With respect to E. Iani et al. (2025), we expanded the two-population scenario also to BAGPIPES.

For BAGPIPES, we adopt the same setup as in P. Rinaldi et al. (2025a) for LRDs, utilizing synthetic templates from G. Bruzual & S. Charlot (2003) with a Kroupa IMF, a stellar mass cutoff of $100 M_{\odot}$, and nebular emission modeled with CLOUDY (G. J. Ferland et al. 2013). We employ a continuity nonparametric star formation history (SFH) model (J. Leja et al. 2019), defining age bin edges (in lookback time) based on Virgil’s spectroscopic redshift, following a logarithmic distribution from $z = 30$ to the cosmic time at that redshift. The same approach is applied to the age parameter.

We let stellar masses range from $10^5$ to $10^{13} M_{\odot}$ (uniform prior in log). We adopt a Calzetti reddening law (D. Calzetti et al. 2000), allowing $A_V$ to vary between 0 and 6, while metallicity ($Z/Z_{\odot}$) spans 0 to 2.5. The ionization parameter is constrained between -4 and 0.001. For the AGN contribution, we follow A. C. Carnall et al. (2023; see their Table 1), adopting broad, flat (uniform) priors for all parameters.

For the two-population scenario, we adopt a delayed SFH with $\tau$ ranging from 0.1 to 13 Gyr. Each population follows the same metallicity, age, and dust constraints as the AGN-included run.
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Cutouts (2.5 × 2.5) centered on Virgil, showing imaging from HST (from F775W only), NIRCam, and MIRI, each displayed at their native pixel scales (30 mas for HST and NIRCam, 60 mas for MIRI), and color coded in blue, orange, and red, respectively. Virgil is highlighted in F356W and F560W. Its elongated morphology, featuring two distinct knots, is clearly visible in the NIRCam bands thanks to their superior resolution. F560W, F770W, and F1000W are from MIDIS; F1500W is from PAHSPECS; while the remaining MIRI bands are from SMILES. Each postage stamp is displayed using its intensity scale centered around the median pixel value to enhance visual contrast across filters.</figcaption></figure>

For CIGALE, we assume a delayed exponentially declining SFH (delayed-$\tau$ model) with two stellar populations. We adopt G. Bruzual & S. Charlot (2003) stellar population models with both solar and subsolar metallicity ($Z = 0.2Z_⊙$) and a Chabrier IMF. Nebular continuum and emission lines are included, assuming solar and subsolar metallicity, and allowing electron density and ionization parameter values of $n_e = 10, 10^2, 10^3 \text{ cm}^{-3}$ and $\log_{10}(\mathcal{U}) = -3, -2, -1$, respectively. For dust attenuation, we adopt D. Calzetti et al. (2000), while far-IR emission is modeled with B. T. Draine et al. (2014). AGN emission is incorporated using the SKIRTOR models (M. Stalevski et al. 2012, 2016), following initial parameters from G. Yang et al. (2023) but allowing for both Type 1 (unobscured) and Type 2 (obscured) AGN.

For SYNTHESIZER-AGN, the SED is modeled as a composite stellar population (P. G. Pérez-González et al. 2003, 2008) with AGN emission from the accretion disk and dusty torus (P. G. Pérez-González et al. 2024a). Stellar emission consists

of a young and an older star formation event, each described by a delayed exponential function with $\tau$ between 1 Myr and 1 Gyr, and ages up to the cosmic time at the Virgil's redshift. Each stellar population experiences independent attenuation following D. Calzetti et al. (2000), with $A_V$ ranging from 0 to 10 mag. The stellar component is modeled with G. Bruzual & S. Charlot (2003) templates, assuming a Chabrier IMF with stellar masses between 0.1 and 100 $M_⊙$, including nebular emission (P. G. Pérez-González et al. 2003). The AGN emission follows a QSO average spectrum (D. E. Vanden Berk et al. 2001; E. Glikman et al. 2006), while its dust component is modeled using self-consistent AGN torus templates from R. Siebenmorgen et al. (2015).

Finally, for PROSPECTOR, we used a modified version built on the standard implementation, which incorporates the Flexible Stellar Population Synthesis (FSPS; C. Conroy & J. E. Gunn 2010) model for the stellar component. We assumed a P. Kroupa (2001) initial mass function and a
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Illustrative example of the modeling and subtraction procedure applied to the foreground LAE at $z \approx 4.77$ located near Virgil. Each row corresponds to a different MIRI filter: F560W (top), F770W (middle), and F1000W (bottom). The three columns show, from left to right, the original data, the best-fit model of the LAE, and the residual after model subtraction (data-model). No modeling was applied to Virgil, which is visible in all panels and highlighted in the F560W stamp. Note that the displayed scales are optimized independently for each panel for visual clarity and are not uniform across filters.</figcaption></figure>

delayed-$\tau$ star formation history. Nebular line and continuum emission were included, as preconfigured in FSPS (N. Byler et al. 2017). We adopted the Calzetti attenuation curve with a flexible slope, following M. Kriek & C. Conroy (2013). For galaxy dust emission, given that Virgil has $z = 6.64$, we used the empirical IR SED model of Haro 11, a low-metallicity, starbursting dwarf galaxy believed to exhibit typical features of first-generation galaxies in the early Universe (J. Lyu et al. 2016; M. E. De Rossi et al. 2018).

This code employs a set of semiempirical AGN SED models optimized for AGN identification and characterization (J. Lyu et al. 2022, 2024). We adopted a model configuration similar to that in J. Lyu et al. (2024) for the SMILES+JADES AGN identification. The AGN component includes both the AGN-powered continuum from the UV to the far-IR and the narrow and broad emission lines from the UV to the NIR, derived from empirical observations. The continuum shape and line strengths of the AGN SED can be adjusted using a hybrid attenuation model, featuring an SMC-like curve for the typical UV-optical attenuation in Type 1 AGNs and an empirical attenuation law for IR obscuration. We refer the reader to J. Lyu et al. (2024) for further details.

## 5.2. Results from SED Fitting

We now discuss the results of the SED fitting. E. Iani et al. (2025) showed that most SED-fitting codes struggle to reproduce the MIRI fluxes unless an AGN component is included, with the exception of SYNTHESIZER-AGN. This is mainly because, by design, SYNTHESIZER-AGN always fits two independent stellar populations—regardless of the inclusion or not of an AGN—providing greater flexibility to account for complex SED shapes.

We first fit Virgil’s photometry without including an AGN component, adopting a model with two stellar populations. As shown in E. Iani et al. (2025), only SYNTHESIZER-AGN (used

here without including an AGN component) was originally capable of such a fit due to its design. Here, we introduced this flexibility in BAGPIPES as well.<fnref n="46" /> In both cases, the inferred total $M_*$ is significantly higher than in AGN-inclusive models (see the next paragraph), with $\log_{10}(M_*/M_\odot)$ ranging from 9.80 to 10.50. The older population consistently exhibits $A_V \approx 2-4$. However, the solutions from BAGPIPES and SYNTHESIZER-AGN differ: the latter yields two very young populations ($t_{age} \lesssim 50$ Myr) with distinct extinctions ($A_V \approx 3.60$ for the most obscured), consistent with E. Iani et al. (2025), while BAGPIPES retrieves two markedly distinct populations—one extremely young ($t_{age} < 3$ Myr), dust-free, and low-mass ($\log_{10}(M_*/M_\odot) \approx 7.24$), and one older (approaching the age of the Universe at $z \approx 6.64$), massive ($\log_{10}(M_*/M_\odot) \approx 10.48$), and heavily obscured ($A_V \approx 4$). The estimated $M_*$, in this case, would overlap with the maximum $M_*$ expected from the stellar mass-halo mass relation of P. Behroozi et al. (2020) at these redshifts.

When including an AGN component, both SYNTHESIZER-AGN and PROSPECTOR consistently yield $\log_{10}(M_*/M_\odot) \approx 9$, although BAGPIPES retrieves a significantly lower stellar mass by $\approx 1$ dex: $\log_{10}(M_*/M_\odot) = 7.92$. Similarly, CIGALE finds $\log_{10}(M_*/M_\odot) = 8.50$. Notably, BAGPIPES predicts a stellar age of $t_{age} = 400^{+260}_{-260}$ Myr, consistent with both CIGALE and PROSPECTOR within the uncertainties ($t_{age} \approx 400-600$ Myr). However, its mass-weighted age, from BAGPIPES, is only $4^{+2}_{-3}$ Myr, indicating that a substantial fraction of the stellar mass formed recently, which is consistent with a very recent burst of star formation.

Interestingly, both SYNTHESIZER-AGN and PROSPECTOR, when invoking an AGN, predict $A_V \approx 2-4$. In particular, SYNTHESIZER-AGN predicts $A_V \approx 0.74$ mag for the young population ($t_{age} \approx 10$ Myr), in line with our estimate from NIRSpec/PRISM data, and $A_V \approx 2.40$ for the older population ($t_{age} \approx 32$ Myr). Similarly, BAGPIPES and CIGALE yield $A_V$ values consistent with those inferred from Balmer lines tracing the young population. Across all fits, the metallicity remains consistent at 10%–20% $Z_\odot$, in agreement with our spectroscopic constraints.

Overall, all models yield a reduced $\chi^2$ between 6 and 12, with BAGPIPES achieving the best $\chi^2$ when including an AGN, while SYNTHESIZER-AGN favors the model of two stellar populations without an AGN. As discussed in E. Iani et al. (2025), fitting this object remains challenging, reflecting the well-known difficulty in modeling LRDs in general. Nevertheless, most SED-fitting codes used in this work—as it was for E. Iani et al. (2025)—favor the (embedded) AGN interpretation.

Interestingly, our modified version of PROSPECTOR tailored for obscured AGNs (J. Lyu et al. 2024) fails to reproduce the observed SED of this object when the AGN component is excluded. We show the best-fit results in Figure 11 along with other models from the literature (see the next section).

# 6. Hide and Seek: What is the Nature of Virgil?

## 6.1. Star-forming Properties

E. Iani et al. (2025) introduced Virgil as the first LRD with a clearly detected host galaxy and a peculiar rising SED between

<footnote n="46">⁴⁶ This flexibility was implemented in BAGPIPES only for the model without the AGN contribution.</footnote>
<figure><img src="image_14.png" /><figcaption>**Figure 11.** The best-fit results for Virgil are shown, with photometric data indicated by black squares and $3\sigma$ upper limits by black inverted triangles. The reduced chi-squared values are also reported and computed consistently across the different SED-fitting codes to account for their varying treatments of upper limits. The red arrow highlights the IR excess relative to the typical SFG template, which both NIRCam and NIRSpec data would bias the SED fitting toward. In contrast, the blue arrow indicates the UV excess compared to the average LRD template proposed by H. B. Akins et al. (2025). For comparison, we include the LRD model from H. B. Akins et al. (2025) and Virgil's SED, excluding the MIRI bands, as derived from BAGPIPES. Additionally, we show the average BL GP SED (green circles) from R. Lin et al. (2025), considering only sources with $M_*$, similar to the estimated value for Virgil. Finally, we compare it with the Blue-excess HotDOG (HSC J1202; red circles) from A. Noboriguchi et al. (2022). RUBIES-BLAGN-1 model from B. Wang et al. (2025) is shown in blue. The vertical dashed red line marks the transition between the NIRCam and MIRI domains. All models are normalized to F444W in the Virgil's reference system.</figcaption></figure>

4 and 10 µm. At that time, their analysis was limited to photometric data from HST and JWST (NIRCam and MIRI). Leveraging multiple SED-fitting codes, we identified two plausible scenarios for Virgil's nature: either a dusty starburst galaxy or an SFG hosting an obscured AGN. However, in the previous sections, we demonstrated that incorporating spectroscopic data from NIRSpec/PRISM strengthens the interpretation of Virgil as a typical SFG during the EoR—or, at the very least, that its AGN nature becomes less clear when the redshift evolution of diagnostic diagrams is taken into account. This is supported by its relatively low $M_*$, metal-poor nature, and extreme $EW_0$ for [O III]λλ4959, 5007. Interestingly, Virgil shares similarities with GPs and BBs, which are often considered the best low-redshift analogs of high-$z$ galaxies.

In general, when focusing on NIRSpec/PRISM data, we find that Virgil appears to be a typical system at $z \approx 6-7$ undergoing a bursty star formation episode. At the time of observation, it is either experiencing or fading out of a starburst, while also harboring an older stellar population that dominates its $M_*$, with formation dating back a few hundred million years, as inferred from BAGPIPES, CIGALE, and PROSPECTOR fits. Virgil exhibits extreme emission-line ratios, including a high R3 index ($\approx 0.8$), consistent with other galaxies at similar redshifts (e.g., A. J. Cameron et al. 2023; M. Tang et al. 2023). This is likely driven by intense star formation, with a surface SFR of $\Sigma_{SFR(H\alpha)} \approx 5.25 M_\odot \text{ yr}^{-1} \text{ kpc}^{-2}$ at a moderately low metallicity

($Z \approx 0.1-0.2 Z_\odot$). Whether it contains an AGN remains ambiguous if the emission-line ratios are interpreted in the context of its redshift and the evolution of these ratios, consistent with trends observed in other LRDs by P. Rinaldi et al. (2025a). If this interpretation holds, it may also apply to other claimed Type 1 and Type 2 AGN in the literature (e.g., V. Kokorev et al. 2023; J. Scholtz et al. 2025).

## 6.2. Virgil's AGN

However, multiple SED-fitting codes that incorporate wavelengths beyond 5 µm suggest that this galaxy hosts a deeply obscured AGN, a characteristic commonly observed in many LRDs. Notably, Virgil's infrared SED rises even more steeply than those of typical LRDs, making it one of the most extreme cases identified to date at $z \approx 6-7$. While Virgil may appear exceptional, it has lower-redshift analogs, such as the obscured AGNs identified in the HUDF by J. Lyu et al. (2024).

Interestingly, R. Lin et al. (2025) recently identified a subset of GP galaxies that exhibit the V-shaped SED commonly found in high-$z$ LRDs, along with broad emission lines, suggesting a possible evolutionary link between these two populations (see Figure 11).

Further exploring analogies with lower redshift sources provides key insights into the nature of Virgil. E. Iani et al. (2025) measured a compactness of $c_{F444W} \approx 1.5$ for Virgil, consistent with the LRD criterion ($c_{F444W} < 2$; V. Kokorev et al. 2024). However, this
threshold is somewhat arbitrary. If Virgil were at $z \approx 3.5$, for instance, the angular scale would be $\approx 7.3 \text{ kpc arcsec}^{-1}$ instead of $5.4 \text{ kpc arcsec}^{-1}$ at $z \approx 6.64$, and the $4.44 \ \mu\text{m}$ emission would be even more dominated by the obscured AGN, which would make Virgil appear even more compact.

Another intriguing analogy is with the very red sources at $z \approx 2$ identified in Spitzer observations, commonly known as Dust Obscured Galaxies (DOGs; e.g., D. W. Weedman et al. 2006; R. S. Bussmann et al. 2012). The most extreme DOGs exhibit infrared SEDs that rise as steeply as that of Virgil, such as SST24 J142648.9+332927 (R. S. Bussmann et al. 2012). DOGs are typically classified into two subgroups based on their distinct SEDs: power-law (PL) DOGs, indicative of AGN dominance, and “bump” DOGs, which show a flattening near rest-frame $1.6 \ \mu\text{m}$, suggesting a significant contribution from a stellar population (e.g., A. Dey et al. 2008). Therefore, “bump” DOGs are considered to correspond to galaxies in a star-forming phase (R. S. Bussmann et al. 2011), while PL DOGs are associated with galaxies in an AGN phase (e.g., R. S. Bussmann et al. 2009). A. Dey et al. (2008) showed that the fraction of PL DOGs among all DOGs increases with increasing MIR flux density, similar to the luminosity dependence of the AGN fraction in ultraluminous infrared galaxies (ULIRGs; D. B. Sanders & I. F. Mirabel 1996; S. Veilleux et al. 1999). Notably, some of these features are commonly observed in LRDs as well (P. G. Pérez-González et al. 2024a; C. C. Williams et al. 2024).

Hence, Virgil is likely experiencing a transition phase at the moment we observe it, allowing us to witness a starburst event with an embedded AGN. As shown in E. Iani et al. (2025), the MIRI part of Virgil’s SED cannot be solely attributed to star formation heated dust emission. Otherwise, the required dust temperature would approach the typical sublimation temperature for the dust ($T_{\text{subl}} \approx 1500-2000$ K, depending on dust composition; M. J. Temple et al. 2021), a scenario that would be difficult to reconcile with “normal” star formation behavior (e.g., L. Sommovigo et al. 2020).

However, the DOG family is more diverse. For instance, HotDOGs are among the most luminous objects in the Universe, with bolometric luminosities often exceeding $10^{13-14} L_{\odot}$ (see an example of a HotDOG, HSC J1202, in Figure 11). Discovered with the Wide-field Infrared Survey Explorer (WISE; E. L. Wright 2010), they exhibit extremely red MIR SEDs. Their emission is dominated by hyperluminous, highly obscured AGNs that power their intense mid-to-far-infrared radiation (e.g., P. R. M. Eisenhardt et al. 2012). Interestingly, X-ray studies indicate that their obscuration levels are exceptionally high, often reaching or exceeding the Compton-thick limit ($N_H > 1.5 \times 10^{24} \text{ cm}^{-2}$) (e.g., C. Ricci et al. 2017). Due to this heavy obscuration, the host galaxy could dominate the UV and optical wavelengths, rather than direct AGN emission.

While most HotDOGs exhibit extreme dust obscuration that suppresses UV emission, a subset displays significant excess UV/optical light compared to the average population—these are known as Blue-excess HotDOGs (BHDs; e.g., R. J. Assef et al. 2016; A. Noboriguchi et al. 2022). R. J. Assef et al. (2016) proposed several plausible scenarios to explain this excess, including the possibility that a small fraction of AGN light leaks into our line of sight via scattering by dust, gas, or both—a scenario recently studied also by M. Stepney et al. (2024) in the case of ULASJ2315+0143, a source around

Cosmic Noon. In this scenario, the intrinsic emission from the accretion disk and broad-line region is largely absorbed, with only $\lesssim 1\%$ scattered into view on larger spatial scales, where dust grains efficiently reflect shorter-wavelength radiation (B. T. Draine 2003).

A similar contrast between the UV and optical colors has been already observed in LRDs (A. Noboriguchi et al. 2023). While most LRDs exhibit extremely red optical-to-NIR colors (0.4–1 $\ \mu\text{m}$), they show a blue excess in the rest-frame UV (0.2–0.4 $\ \mu\text{m}$), reminiscent of Blue-Excess HotDOGs found at $z \approx 2-3$. However, LRDs represent a less extreme manifestation of this phenomenon compared to Blue-Excess HotDOGs, as they are found at higher redshifts and exhibit lower $M_*$ and Eddington ratios ($\lambda_{\text{Edd}} \approx 0.1-1$), suggesting they may be earlier-stage counterparts.

Nonetheless, Virgil stands out as an even more extreme case among LRDs, exhibiting both an exceptionally steep SED between F444W and F1500W and an enhanced UV excess. Indeed, Virgil’s UV excess is more pronounced compared to the average LRD template proposed by H. B. Akins et al. (2025), suggesting a greater contribution from either scattered AGN light or, more likely, a strong star-forming component. At the same time, its SED rises more steeply than seen in most LRDs at wavelengths beyond $0.6 \ \mu\text{m}$ (rest frame), resembling the most obscured DOGs. This dual nature—both an extreme red excess in the MIRI regime and an unusually strong blue excess in the UV—potentially places Virgil as a high-$z$ analog of Blue-excess HotDOGs, possibly with a strong contribution at the UV regime from its host.

A. Noboriguchi et al. (2023) suggested that the connection between LRDs and Blue-excess HotDOGs reflects an evolutionary sequence in which gas-rich mergers trigger dusty starbursts, giving rise to heavily obscured DOGs. Interestingly, when modeling Virgil with FORCEPHO, it appears to be composed of two distinct knots in the NIRCam wavelength range, suggesting a possible merger-driven scenario. As AGN feedback disperses the surrounding dust, some of these systems develop a blue excess, akin to what is observed in Blue-excess HotDOGs. In this context, LRDs may represent an even earlier stage in this process, where lower metallicities and dust-to-gas ratios at high redshifts lead to a higher fraction of such objects exhibiting blue excess. Yet, Virgil appears even more extreme, with a UV excess likely arising from a combination of AGN activity and star formation, though primarily driven by the host galaxy, while its obscured SMBH shapes an extreme near-IR SED (rest frame).

As a final analogy with sources at lower redshifts, in Figure 11, we compared Virgil to another LRD at a much lower redshift ($z_{\text{spec}} = 3.1$), RUBIES-BLAGN-1, studied by B. Wang et al. (2025). RUBIES-BLAGN-1 closely follows Virgil’s SED at wavelengths beyond $4.4 \ \mu\text{m}$. While both objects exhibit a similar trend at IR wavelengths, they diverge significantly in the UV, further emphasizing the striking UV excess observed in Virgil compared to the average LRD population. This divergence highlights once again the peculiar nature of Virgil.

## 7. Summary and Conclusions

In this paper, we extend the analysis of Virgil, first presented in E. Iani et al. (2025), using recent data collected by the MIDIS and PAHSPECS teams for MIRI, and NIRSpec/
PRISM data by the OASIS team. Virgil is an LAE at $z_{spec} = 6.6379 \pm 0.003$, detected with VLT/MUSE (R. Bacon et al. 2023). It shares key photometric properties with the LRD population, including its compactness in F444W ($c_{F444W} \approx 1.5$; E. Iani et al. 2025), although it does not strictly meet all their color selection criteria (V. Kokorev et al. 2024; P. Rinaldi et al. 2025a).

Our main findings can be summarized as follows:

1. Virgil exhibits low to moderate dust attenuation based on its Balmer decrement, assuming the SMC reddening law from K. D. Gordon et al. (2003). Its metallicity aligns with galaxies of similar $M_*$ at $z \gtrsim 6$ (Figure 2). By comparing the UV- and Hα-based SFRs, we find that Virgil may be transitioning into—or fading out of—a bursty phase. Its $f_{esc, LyC}$ and $\xi_{ion}$ indicate a limited role in Cosmic Reionization, further supported by its moderate $EW_0(H\alpha)$ compared to other emitters at similar redshifts (e.g., P. Rinaldi et al. 2023, 2024; C. Simmonds et al. 2024b). Overall, Virgil’s spectral properties align with the average galaxy population during the EoR (e.g., A. J. Cameron et al. 2023). A full summary of its spectral properties is provided in Table 2.

2. The combination of NIRSpec/PRISM data with extensive photometric coverage from HST and JWST, including MIRI detections up to 15 µm, offers a unique opportunity to further investigate Virgil’s nature. Without an a priori assumption that it could host an AGN—as previously suggested by E. Iani et al. (2025) and based solely on photometry—we analyze its ISM properties using the R23 versus O32 and Ne3O2Hd versus Ne3O2 diagrams (Figure 3). By considering CLOUDY models for SFGs (following the approach outlined in A. Calabrò et al. 2024a), we find that Virgil’s ISM does not significantly deviate from the typical high-z SFGs during the EoR. Interestingly, its properties resemble well-studied GPs and BBs, which are among the best local analogs of high-z galaxies. Notably, despite its moderate $EW_0(H\alpha)$, Virgil exhibits a huge $EW_0$ for [O III] $\lambda\lambda4959, 5007$ (1514 ± 20 Å), comparable to the Extreme Line Emitters studied by K. Boyett et al. (2024) at similar redshifts.

3. To assess the AGN nature suggested by E. Iani et al. (2025), we search for broad Balmer lines (a feature already found in other claimed LRDs), modeling $H\alpha$ with and without a broad component. However, the NIRSpec/PRISM data do not provide a firm conclusion (Figure 4). We note that the source’s position within the slit (Figure 1, right panel)—with the red knot lying near the edge or possibly in the shutter gap—may limit our sensitivity to broad-line emission. Nonetheless, assuming the broadening is real, the inferred $M_*$, $L_{bol}$, and Eddington ratio align with expectations for LRDs (Figure 5) and with values already reported in E. Iani et al. (2025) by leveraging photometric data only. We examine emission-line diagnostics and find that while Virgil lies above the SFG-AGN separation line in the “OHNO” diagram (B. E. Backhaus et al. 2022), it falls below it when accounting for redshift evolution. Similarly, O3Hg-based diagnostics (G. Mazzolari et al. 2024) place it among AGNs at high $z$ (e.g., K. Nakajima et al. 2023; J. Scholtz et al. 2025), but when accounting for redshift evolution, its classification becomes

ambiguous (Figure 7). Overall, while these high-z diagnostics would classify Virgil as an AGN, they do not provide a definitive classification when the redshift evolution is taken into account, instead pointing to a likely mixed nature.

4. The newly acquired MIRI data from MIDIS and PAHSPECS, including a clear 15 µm detection of Virgil, enabled improved photometry and SED fitting. Using FORCEPHO, we identified two distinct knots, reinforcing evidence of complex UV morphology in LRDs (P. Rinaldi et al. 2025a), with one component accounting for the MIRI excess. Virgil appears even redder in F444W–F1500W color ($\approx 2.84$) than typical LRDs and, while it does not meet all standard LRD color criteria (V. Kokorev et al. 2024), it satisfies the V-shaped SED selection from D. D. Kocevski et al. (2025), suggesting that MIRI bands may be crucial for identifying extreme cases (Figure 8). To contextualize its spectral properties, we applied multiple SED-fitting codes (Figure 11), extending the analysis of E. Iani et al. (2025) by incorporating two independent stellar populations in BAGPIPES and utilizing an AGN-optimized version of PROSPECTOR (J. Lyu et al. 2024). Interestingly, the tailored version of PROSPECTOR fails to reproduce the observed SED when an AGN component is not included. While modeling this source remains challenging, our results show that the fits consistently require the presence of a dust-obscured AGN, consistent with the findings of E. Iani et al. (2025).

5. Combining photometric and spectroscopic data, we find Virgil to be among the most extreme LRDs identified to date, with a steeply rising SED beyond 4.44 µm (observed frame) and a pronounced UV excess. Its SED closely resembles Cosmic Noon objects observed with Spitzer, such as SST24 J142648.9+33292, a well-studied DOG, and shares similarities with HotDOG, whose IR properties align with those of LRDs. A. Noboriguchi et al. (2023) linked JWST/EROs to DOGs, particularly Blue-Excess HotDOG, which exhibit extreme IR colors and a flat UV continuum. The UV excess in Blue-Excess HotDOG has been attributed to AGN light ($<1\%$; R. J. Assef et al. 2016, 2020) scattered by dust or gas, a scenario that could also contribute to Virgil’s properties. A similar finding was also reported for ULASJ2315+014, a Cosmic Noon source, by M. Stepney et al. (2024). However, Virgil’s SED and emission-line features suggest a mix of AGN and star formation activity, with the UV emission likely dominated by the host galaxy. Finally, we compare Virgil to RUBIES-BLAGN-1 ($z_{spec} = 3.1$) from B. Wang et al. (2025), noting that their SEDs diverge at shorter wavelengths, with Virgil exhibiting a stronger UV excess.

This study underscores the challenges of studying such systems at high redshift and highlights the necessity of a multiwavelength approach to identify and characterize extreme sources during Cosmic Reionization. Deep MIRI imaging at $z > 6$ is crucial for uncovering highly dust-obscured AGNs that might otherwise be misclassified as normal SFGs, especially in studies limited to wavelengths below 5 µm. In such cases, both photometric and spectroscopic information can lead to ambiguous classifications, particularly when
redshift evolution is taken into account. This reinforces the need to expand MIRI coverage to systematically identify more of these extreme objects and refine our understanding of early black hole growth.

## Acknowledgments

The authors are deeply grateful to Antonello Calabrò for valuable insights on CLOUDY and PYCLOUDY, and for publicly sharing their SFG and AGN models, which were used as a reference to verify the consistency of our photoionization models. The authors also thank Adam Carnall for insightful input on BAGPIPES and for assistance with the implementation of the two-population model adopted in this work. Finally, they also thank Camilla Pacifici, Vasily Kokorev, and Cristian Vignali for their insightful discussions.

This work is based on observations made with the NASA/ESA/CSA JWST. The data were obtained from the Mikulski Archive for Space Telescopes (MAST) 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. These observations are associated with JWST programs GTO #1180, GO #1210, GTO#1283, GO #1963, GO #1895, GO# 3215, and GO#6511.

The authors acknowledge the FRESCO, JEMS, and #3215 teams led by co-PIs P. Oesch, C. C. Williams, M. Maseda, D. Eisenstein, and R. Maiolino for developing their observing program with a zero-exclusive-access period. Processing for the JADES NIRCam data release was performed on the lux cluster at the University of California, Santa Cruz, funded by NSF MRI grant AST 1828315. Also based on observations made with the NASA/ESA Hubble Space Telescope obtained from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 526555. The data presented in this article were obtained from MAST at the Space Telescope Science Institute. The specific observations analyzed can be accessed via doi: 10.17909/1rq3-8048 P. Oesch & D. Magee (2023), C. Williams et al. (2023), G. Illingworth (2015), and M. Rieke et al. (2023).

A.J.B. acknowledges funding from the “FirstGalaxies” Advanced Grant from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement No. 789056).

P.G.P.-G. acknowledges support from grant PID2022-139567NB-I00 funded by the Spanish Ministerio de Ciencia e Innovación MCIN/AEI/10.13039/501100011033, FEDER, UE.

B.E.R. acknowledges support from the NIRCam Science Team contract to the University of Arizona, NAS5-02015, and JWST Program 3215.

S.T. acknowledges support by the Royal Society Research Grant G125142.

The research of C.C.W. is supported by NOIRLab, which is managed by the Association of Universities for Research in Astronomy (AURA) under a cooperative agreement with the National Science Foundation.

J.W. gratefully acknowledges support from the Cosmic Dawn Center through the DAWN Fellowship. The Cosmic Dawn Center (DAWN) is funded by the Danish National Research Foundation under grant No. 140.

Y.Z., Z.J., and P.L. gratefully acknowledge the JWST/NIRCam contract to the University of Arizona NAS5-02015.

The work of G.H.R. and P.L. was also supported by grant 80NSSC18K0555, from the NASA Goddard Space Flight Center to the University of Arizona.

H.Ü. acknowledges funding by the European Union (ERC APEX, 101164796). Views and opinions expressed are however those of the authors only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them.

G.C.J. acknowledges support by the Science and Technology Facilities Council (STFC), ERC Advanced grant 695671 “QUENCH.”

A.C.G. acknowledges support by JWST contract B0215/JWST-GO-02926.

G.O. acknowledges support from the Swedish National Space Agency (SNSA).

H.I. acknowledges support from JSPS KAKENHI grant No. JP21H01129.

M.A. gratefully acknowledges support from ANID Basal Project FB210003 and ANID MILENIO NCN2024_112.

T.D.S. acknowledges the research project was supported by the Hellenic Foundation for Research and Innovation (HFRI) under the “2nd Call for HFRI Research Projects to Support Faculty Members and Researchers” (project No.: 03382).

R.M. acknowledges support by the Science and Technology Facilities Council (STFC), by the ERC through Advanced grant 695671 “QUENCH,” and by the UKRI Frontier Research grant RISEandFALL. R.M. also acknowledges funding from a research professorship from the Royal Society.

I.S. acknowledges funding from the Atracción de Talento grant No. 2022-T1/TIC-20472 of the Comunidad de Madrid, Spain, and the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant No. 101117541, DistantDust).

K.I.C. acknowledges funding from the Dutch Research Council (NWO) through the award of the Vici grant VI. C.212.036.

**Facilities:** HST, JWST.

**Software:** ASTROPY (Astropy Collaboration et al. 2022), BAGPIPES (A. C. Carnall et al. 2019), MSAEXP (G. Brammer 2023) NUMPY (C. R. Harris et al. 2020), PANDAS (The pandas development team 2024) PHOTUTILS (L. Bradley et al. 2016), TOPCAT (M. Taylor 2022).

## ORCID iDs

Pierluigi Rinaldi https://orcid.org/0000-0002-5104-8245

Pablo G. Pérez-González https://orcid.org/0000-0003-4528-5639

George H. Rieke https://orcid.org/0000-0003-2303-6519

Jianwei Lyu https://orcid.org/0000-0002-6221-1829

Francesco D'Eugenio https://orcid.org/0000-0003-2388-8172

Zihao Wu https://orcid.org/0000-0002-8876-5248

Stefano Carniani https://orcid.org/0000-0002-6719-380X

Tobias J. Looser https://orcid.org/0000-0002-3642-2446

Irene Shivaei https://orcid.org/0000-0003-4702-7561

Leindert A. Boogaard https://orcid.org/0000-0002-3952-8588

Tanio Diaz-Santos https://orcid.org/0000-0003-0699-6083

Luis Colina https://orcid.org/0000-0002-9090-4227
Göran Östlin https://orcid.org/0000-0002-3005-1349

Stacey Alberts https://orcid.org/0000-0002-8909-8782

Javier Álvarez-Márquez https://orcid.org/0000-0002-7093-1877

Marianna Annuziatella https://orcid.org/0000-0002-8053-8040

Manuel Aravena https://orcid.org/0000-0002-6290-3198

Rachana Bhatawdekar https://orcid.org/0000-0003-0883-2226

Andrew J. Bunker https://orcid.org/0000-0002-8651-9879

Karina I. Caputi https://orcid.org/0000-0001-8183-1460

Stéphane Charlot https://orcid.org/0000-0003-3458-2275

Alejandro Crespo Gómez https://orcid.org/0000-0003-2119-277X

Mirko Curti https://orcid.org/0000-0002-2678-2560

Andreas Eckart https://orcid.org/0000-0000-0000-0000

Steven Gillman https://orcid.org/0000-0001-9885-4589

Kevin Hainline https://orcid.org/0000-0003-4565-8239

Nimisha Kumari https://orcid.org/0000-0002-5320-2568

Jens Hjorth https://orcid.org/0000-0002-4571-2306

Edoardo Iani https://orcid.org/0000-0001-8386-3546

Hanae Inami https://orcid.org/0000-0003-4268-0393

Zhiyuan Ji https://orcid.org/0000-0001-7673-2257

Benjamin D. Johnson https://orcid.org/0000-0002-9280-7594

Gareth C. Jones https://orcid.org/0000-0002-0267-9024

Álvaro Labiano https://orcid.org/0000-0002-0690-8824

Roberto Maiolino https://orcid.org/0000-0002-4985-3819

Jens Melinder https://orcid.org/0000-0003-0470-8754

Thibaud Moutard https://orcid.org/0000-0002-3305-9901

Florian Peissker https://orcid.org/0000-0002-9850-2708

Marcia Rieke https://orcid.org/0000-0002-7893-6170

Brant Robertson https://orcid.org/0000-0002-4271-0364

Jan Scholtz https://orcid.org/0000-0001-6010-6809

Sandro Tacchella https://orcid.org/0000-0002-8224-4505

Paul P. van der Werf https://orcid.org/0000-0001-5434-5942

Fabian Walter https://orcid.org/0000-0003-4793-7880

Christina C. Williams https://orcid.org/0000-0003-2919-7495

Chris Willott https://orcid.org/0000-0002-4201-7367

Joris Witstok https://orcid.org/0000-0002-7595-121X

Hannah Übler https://orcid.org/0000-0003-4891-0794

Yongda Zhu https://orcid.org/0000-0003-3307-7525

## References

Abazajian, K. N., Adelman-McCarthy, J. K., Agüeros, M. A., et al. 2009, ApJS, 182, 543

Akhlaghi, M., & Ichikawa, T. 2017, in ASP, Astronomical Data Analysis Software and Systems XXV, 512 (ASP), 213

Akins, H. B., Casey, C. M., Lambrides, E., et al. 2025, ApJ, 991, 37

Alberts, S., Lyu, J., Shivaei, I., et al. 2024, ApJ, 976, 224

Amorín, R., Pérez-Montero, E., Contini, T., et al. 2015, A&A, 578, A105

Anders, P., & Fritze-v. Alvensleben, U. 2003, A&A, 401, 1063

Armah, M., Dors, O. L., Aydar, C. P., et al. 2021, MNRAS, 508, 371

Armus, L., Charmandaris, V., & Soifer, B. T. 2020, NatAs, 4, 467

Arrabal Haro, P., Dickinson, M., Finkelstein, S. L., et al. 2023, Natur, 622, 707

Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, ARA&A, 47, 481

Assef, R. J., Brightman, M., Walton, D. J., et al. 2020, ApJ, 897, 112

Assef, R. J., Walton, D. J., Brightman, M., et al. 2016, ApJ, 819, 111

Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167

Atek, H., Furtak, L. J., Oesch, P., et al. 2022, MNRAS, 511, 4464

Backhaus, B. E., Bridge, J. S., Trump, J. R., et al. 2023, ApJ, 943, 37

Backhaus, B. E., Trump, J. R., Cleri, N. J., et al. 2022, ApJ, 926, 161

Backhaus, B. E., Trump, J. R., Pirzkal, N., et al. 2024, ApJ, 962, 195

Backhaus, B. E., Cleri, N. J., Trump, J. R., et al. 2025, arXiv:250203519B

Bacon, R., Brinchmann, J., Conseil, S., et al. 2023, A&A, 670, A4

Bagley, M. B., Pirzkal, N., Finkelstein, S. L., et al. 2024, ApJ, 965, L6

Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5

Barro, G., Pérez-González, P. G., Kocevski, D. D., et al. 2024a, ApJ, 963, 128

Barro, G., Perez-Gonzalez, P. G., Kocevski, D. D., et al. 2024b, arXiv:2412.01887

Behroozi, P., Conroy, C., Wechsler, R. H., et al. 2020, MNRAS, 499, 5702

Berg, D. A., Skillman, E. D., Marble, A. R., et al. 2012, ApJ, 754, 98

Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393

Bian, F., Kewley, L. J., & Dopita, M. A. 2018, ApJ, 859, 175

Binette, L. 1985, A&A, 143, 334

Binette, L., Dopita, M. A., & Tuohy, I. R. 1985, ApJ, 297, 476

Binette, L., Zovaro, H. R. M., Villar Martín, M., et al. 2024, A&A, 684, A53

Blain, A. W., Smail, I., Ivison, R. J., Kneib, J. P., & Frayer, D. T. 2002, PhR, 369, 111

Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103

Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2011, ApJ, 737, 90

Boyett, K., Bunker, A. J., Curtis-Lake, E., et al. 2024, MNRAS, 535, 1796

Bradač, M. 2020, NatAs, 4, 478

Bradley, L., Sipocz, B., Robitaille, T., et al., 2016 Photutils: Photometry tools, Astrophysics Source Code Library, ascl:1609.011

Bradley, L., Sipőcz, B., Robitaille, T., et al., 2022 astropy/photutils, v1.4.0, Zenodo

Brammer, G., 2023 msaexp: NIRSpec analyis tools, v0.6.17, Zenodo

Brinchmann, J. 2023, MNRAS, 525, 2087

Brunker, S. W., Salzer, J. J., Janowiecki, S., Finn, R. A., & Helou, G. 2020, ApJ, 898, 68

Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000

Bunker, A. J. & NIRSPEC Instrument Science TeamJAESs Collaboration 2020, in IAU, Symp. Uncovering Early Galaxy Evolution in the ALMA and JWST Era Place, 352 (Cambridge Univ. Press), 342

Bunker, A. J., Cameron, A. J., Curtis-Lake, E., et al. 2024, A&A, 690, A288

Bunker, A. J., Saxena, A., Cameron, A. J., et al. 2023, A&A, 677, A88

Bussmann, R. S., Dey, A., Armus, L., et al. 2012, ApJ, 744, 150

Bussmann, R. S., Dey, A., Borys, C., et al. 2009, ApJ, 705, 184

Bussmann, R. S., Dey, A., Lotz, J., et al. 2011, ApJ, 733, 21

Byler, N., Dalcanton, J. J., Conroy, C., & Johnson, B. D. 2017, ApJ, 840, 44

Calabrò, A., Castellano, M., Zavala, J. A., et al. 2024a, ApJ, 975, 245

Calabrò, A., Pentericci, L., Feltre, A., et al. 2023, A&A, 679, A80

Calabrò, A., Pentericci, L., Santini, P., et al. 2024b, A&A, 690, A290

Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682

Cameron, A. J., Saxena, A., Bunker, A. J., et al. 2023, A&A, 677, A115

Caputi, K. I., Caminha, G. B., Fujimoto, S., et al. 2021, ApJ, 908, 146

Caputi, K. I., Deshmukh, S., Ashby, M. L. N., et al. 2017, ApJ, 849, 45

Cardamone, C. N., Schawinski, K., Sarzi, M., et al. 2009, MNRAS, 399, 1191

Carnall, A. C., McLure, R. J., Dunlop, J. S., et al. 2019, MNRAS, 490, 417

Carnall, A. C., McLure, R. J., Dunlop, J. S., et al. 2023, Natur, 619, 716

Carniani, S., Hainline, K., D'Eugenio, F., et al. 2024a, Natur, 633, 318

Carniani, S., Hainline, K., D'Eugenio, F., et al. 2024b, Natur, 633, 318

Chisholm, J., Saldana-Lopez, A., Flury, S., et al. 2022, MNRAS, 517, 5104

Cleri, N. J., Olivier, G. M., Hutchison, T. A., et al. 2023, ApJ, 953, 10

Colina, L., Crespo Gómez, A., Álvarez Márquez, J., et al. 2023, A&A, 673, L6
Conroy, C., & Gunn, J. E., 2010 FSPS: Flexible Stellar Population Synthesis, Astrophysics Source Code Library, ascl:1010.043

Curti, M., Cresci, G., Mannucci, F., et al. 2017, MNRAS, 465, 1384

Curti, M., D'Eugenio, F., Carniani, S., et al. 2023, MNRAS, 518, 425

Curti, M., Maiolino, R., Curtis-Lake, E., et al. 2024, A&A, 684, A75

Curti, M., Mannucci, F., Cresci, G., & Maiolino, R. 2020, MNRAS, 491, 944

Curti, M., Carniani, S., Cameron, A., et al. 2023, NatAs, 7, 622

de Graaff, A., Rix, H.-W., Carniani, S., et al. 2024, A&A, 684, A87

De Rossi, M. E., Rieke, G. H., Shivaei, I., Bromm, V., & Lyu, J. 2018, ApJ, 869, 4

D'Eugenio, F., Cameron, A. J., Scholtz, J., et al. 2025, ApJS, 277, 4

Dey, A., Soifer, B. T., Desai, V., et al. 2008, ApJ, 677, 943

Dickinson, M., Giavalisco, M. & GOODS Team 2003, in Proc. of the ESO/ USM Workshop "The Mass of Galaxies at Low and High Redshift", ed. R. Bender & A. Renzini (ESO), 324

Dors, O. L., Maiolino, R., Cardaci, M. V., et al. 2020, MNRAS, 496, 3209

Drape, B. T. 2003, ApJ, 598, 1017

Drape, B. T., Aniano, G., Krause, O., et al. 2014, ApJ, 780, 172

Eisenhardt, P. R. M., Wu, J., Tsai, C.-W., et al. 2012, ApJ, 755, 173

Eisenstein, D. J., Johnson, B. D., Robertson, B., et al. 2023a, arXiv:2310.12340

Eisenstein, D. J., Willott, C., Alberts, S., et al. 2023b, arXiv:2306.02465
Eldridge, J. J., Stanway, E. R., Xiao, L., et al. 2017, PASA, 34, e058

Elston, R., Rieke, G. H., & Rieke, M. J. 1988, ApJ, 331, L77

Elston, R., Rieke, M. J., & Rieke, G. H. 1989, ApJ, 341, 80

Faisst, A. L., Capak, P. L., Emami, N., Tacchella, S., & Larson, K. L. 2019, ApJ, 884, 133

Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10

Feltre, A., Charlot, S., & Gutkin, J. 2016, MNRAS, 456, 3354

Ferland, G. J., Chatzikos, M., Guzmán, F., et al. 2017, RMXAA, 53, 385

Ferland, G. J., Porter, R. L., van Hoof, P. A. M., et al. 2013, RMXAA, 49, 137

Ferruit, P., Jakobsen, P., Giardino, G., et al. 2022, A&A, 661, A81

Feuillet, L. M., Meléndez, M., Kraemer, S., et al. 2024, ApJ, 962, 104

Furtak, L. J., Zitrin, A., Plat, A., et al. 2023, ApJ, 952, 124

Gardner, J. P., Mather, J. C., Abbott, R., et al. 2023, PASP, 135, 068001

Gehrz, R. D., Roellig, T. L., Werner, M. W., et al. 2007, RSCI, 78, 011302

Giavalisco, M., Ferguson, H. C., Koekemoer, A. M., et al. 2004, ApJ, 600, L93

Gillman, S., Pye, J. P., Alonso-Herrero, A., et al. 2025, arXiv:2501.11491

Glikman, E., Djorgovski, S. G., Stern, D., Bogosavljević, M., & Mahabal, A. 2007, ApJ, 663, L73

Glikman, E., Helfand, D. J., & White, R. L. 2006, ApJ, 640, 579

Gordon, K. D., Clayton, G. C., Misselt, K. A., Landolt, A. U., & Wolff, M. J. 2003, ApJ, 594, 279

Greene, J. E., Labbe, I., Goulding, A. D., et al. 2024, ApJ, 964, 39

Greene, J. E., Seth, A., Kim, M., et al. 2016, ApJ, 826, L32

Groves, B. A., Dopita, M. A., & Sutherland, R. S. 2004, ApJS, 153, 75

Hainline, K. N., Johnson, B. D., Robertson, B., et al. 2024, ApJ, 964, 71

Hainline, K. N., Maiolino, R., Juodzbalis, P. A., et al. 2025, ApJ, 979, 138

Harikane, Y., Zhang, Y., Nakajima, K., et al. 2023, ApJ, 959, 39

Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Natur, 585, 357

Heintz, K. E., Brammer, G. B., Watson, D., et al. 2025, A&A, 693, A60

Helou, G., & Kessler, M. F. 1995, in ASP Conf. Ser. 73, From Gas to Stars to Dust, ed. M. R. Haas, J. A. Davidson, & E. F. Erickson (ASP), 659

Henry, A., Scarlata, C., Martin, C. L., & Erb, D. 2015, ApJ, 809, 19

Hirschmann, M., Charlot, S., Feltre, A., et al. 2023, MNRAS, 526, 3610

Hopkins, P. F., Strauss, M. A., Hall, P. B., et al. 2004, AJ, 128, 1112

Huang, K.-H., Bradač, M., Lemaux, B. C., et al. 2016, AJ, 817, 11

Iani, E., Caputi, K. I., Rinaldi, P., et al. 2024, ApJ, 963, 97

Iani, E., Rinaldi, P., Caputi, K. I., et al. 2025, ApJ, 989, 160

Ilingworth, G. 2015, Hubble Legacy Fields ("HLF"), STScI/MAST,

Isobe, Y., Ouchi, M., Nakajima, K., et al. 2023, ApJ, 956, 139

Izotov, Y. I., Schaerer, D., Thuan, T. X., et al. 2016, MNRAS, 461, 3683

Izotov, Y. I., Stasińska, G., Meynet, G., Guseva, N. G., & Thuan, T. X. 2006, A&A, 448, 955

Izotov, Y. I., Thuan, T. X., Guseva, N. G., & Liss, S. E. 2018, MNRAS, 473, 1956

Jakobsen, P., Ferruit, P., Alves de Oliveira, C., et al. 2022, A&A, 661, A80

Jaskot, A. E., & Oey, M. S. 2013, ApJ, 766, 91

Johnson, B. D., Leja, J., Conroy, C., & Speagle, J. S. 2021, ApJS, 254, 22

Juodžbalis, I., Maiolino, R., Baker, W. M., et al. 2024, Natur, 636, 594

Kauffmann, G., Heckman, T. M., White, S. D. M., et al. 2003, MNRAS, 341, 33

Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531

Kewley, L. J., Nicholls, D. C., Sutherland, R., et al. 2019, ApJ, 880, 16

Killi, M., Watson, D., Brammer, G., et al. 2024, A&A, 691, A52

Kocevski, D. D., Finkelstein, S. L., Barro, G., et al. 2025, ApJ, 986, 126

Kocevski, D. D., Onoue, M., Inayoshi, K., et al. 2023, ApJ, 954, L4

Kokorev, V., Caputi, K. I., Greene, J. E., et al. 2024, ApJ, 968, 38

Kokorev, V., Fujimoto, S., Labbe, I., et al. 2023, ApJL, 957, L7

Koyama, Y., Shimakawa, R., Yamamura, I., Kodama, T., & Hayashi, M. 2019, PASJ, 71, 8

Kriek, M., & Conroy, C. 2013, ApJ, 775, L16

Kron, R. G. 1980, ApJS, 43, 305

Kroupa, P. 2001, MNRAS, 322, 231

Kuruvanthodi, A., Schaerer, D., Marques-Chaves, R., et al. 2024, A&A, 691, A310

Labbé, I., van Dokkum, P., Nelson, E., et al. 2023, Natur, 616, 266

Lacy, M., & Sajina, A. 2020, NatAs, 4, 352

Lam, D., Bouwens, R. J., Coe, D., et al. 2019, arXiv:1903.08177

Langan, I., Ceverino, D., & Finlator, K. 2020, MNRAS, 494, 1988

Langeroodi, D., & Hjorth, J. 2024a, arXiv:2409.07455

Langeroodi, D., & Hjorth, J. 2024b, arXiv:2404.13045

Langeroodi, D., Hjorth, J., Chen, W., et al. 2023, ApJ, 957, 39

Langeroodi, D., Hjorth, J., Ferrara, A., & Gall, C. 2024, arXiv:2410.14671

Larson, R. L., Finkelstein, S. L., Kocevski, D. D., et al. 2023, ApJ, 953, L29

Leja, J., Carnall, A. C., Johnson, B. D., Conroy, C., & Speagle, J. S. 2019, ApJ, 876, 3

Li, M., Cai, Z., Bian, F., et al. 2023, ApJ, 955, L18

Libralato, M., Argyriou, I., Dicken, D., et al. 2024, PASP, 136, 034502

Liddle, A. R. 2007, MNRAS, 377, L74

Lin, R., Zheng, Z.-Y., Jiang, C., et al. 2025, ApJL, 980, L34

Lovell, C. C., Vijayan, A. P., Thomas, P. A., et al. 2021, MNRAS, 500, 2127

Luo, B., Brandt, W. N., Xue, Y. Q., et al. 2017, ApJS, 228, 2

Luridiana, V., Morisset, C., & Shaw, R. A. 2015, A&A, 573, A42

Lyu, J., Alberts, S., Rieke, G. H., & Rujopakarn, W. 2022, ApJ, 941, 191

Lyu, J., Alberts, S., Rieke, G. H., et al. 2024, ApJ, 966, 229

Lyu, J., Rieke, G. H., & Alberts, S. 2016, ApJ, 816, 85

Maiolino, R., Scholtz, J., Witstok, J., et al. 2024, Natur, 627, 59

Mascia, S., Pentericci, L., Calabrò, A., et al. 2023, A&A, 672, A155

Matsuoka, Y., Onoue, M., Kashikawa, N., et al. 2019, ApJ, 872, L2

Matthee, J., Feltre, A., Maseda, M., et al. 2022, A&A, 660, A10

Matthee, J., Naidu, R. P., Brammer, G., et al. 2024, ApJ, 963, 129

Matthee, J., Sobral, D., Best, P., et al. 2017, MNRAS, 465, 3637

Mazzolari, G., Übler, H., Maiolino, R., et al. 2024, A&A, 691, A345

McClymont, W., Tacchella, S., D'Eugenio, F., et al. 2025, MNRAS, 540, 190

Mingozzi, M., James, B. L., Arellano-Córdova, K. Z., et al. 2022, ApJ, 939, 110

Mobasher, B., Dickinson, M., Ferguson, H. C., et al. 2005, ApJ, 635, 832

Nagao, T., Murayama, T., & Taniguchi, Y. 2001, ApJ, 549, 155

Nakajima, K., & Maiolino, R. 2022, MNRAS, 513, 5134

Nakajima, K., Ouchi, M., Isobe, Y., et al. 2023, ApJS, 269, 33

Navarro-Carrera, R., Rinaldi, P., Caputi, K. I., et al. 2024, arXiv:2410.23249

Nelson, D., Pillepich, A., Springel, V., et al. 2019, MNRAS, 490, 3234

Netzer, H. 2009, MNRAS, 399, 1907

Neugebauer, G., Habing, H. J., van Duinen, R., et al. 1984, ApJ, 278, L1

Noboriguchi, A., Inoue, A. K., Nagao, T., Toba, Y., & Misawa, T. 2023, ApJ, 959, L14

Noboriguchi, A., Nagao, T., Toba, Y., et al. 2022, ApJ, 941, 195

Oesch, P., & Magee, D. 2023, The JWST FRESCO Survey, STScI/MAST,

Oesch, P. A., Bouwens, R. J., Illingworth, G. D., et al. 2014, ApJ, 786, 108

Oesch, P. A., Brammer, G., Naidu, R. P., et al. 2023, MNRAS, 525, 2864

Oesch, P. A., Brammer, G., van Dokkum, P. G., et al. 2016, ApJ, 819, 129

Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713

Osterbrock, D. E., & Ferland, G. J. 2006, Astrophysics of Gaseous Nebulae and Active Galactic Nuclei (Univ. Science Books)

Östlin, G., Pérez-González, P. G., Melinder, J., et al. 2025, A&A, 696, A570

Pacucci, F., Nguyen, B., Carniani, S., Maiolino, R., & Fan, X. 2023, ApJ, 957, L3

Papovich, C., Moustakas, L. A., Dickinson, M., et al. 2006, ApJ, 640, 92

Papovich, C., Simons, R. C., Estrada-Carpenter, V., et al. 2022, ApJ, 937, 22

Perna, M., Lanzuisi, G., Brusa, M., Cresci, G., & Mignoli, M. 2017, A&A, 606, A96

Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1

Pustilnik, S. A., Egorova, E. S., Kniazev, A. Y., et al. 2021, MNRAS, 507, 944

Pérez-González, P. G., Gil de Paz, A., Zamorano, J., et al. 2003, MNRAS, 338, 508

Pérez-González, P. G., Barro, G., Rieke, G. H., et al. 2024a, ApJ, 968, 4

Pérez-González, P. G., Costantin, L., Langeroodi, D., et al. 2023, ApJ, 951, L1

Pérez-González, P. G., Rieke, G. H., Villar, V., et al. 2008, ApJ, 675, 234

Pérez-González, P. G., Rinaldi, P., Caputi, K. I., et al. 2024b, ApJ, 969, L10

Ranalli, P., Comastri, A., Vignali, C., et al. 2013, A&A, 555, A42

Reddy, N. A., Shapley, A. E., Sanders, R. L., et al. 2018, ApJ, 869, 92

Reines, A. E., Greene, J. E., & Geha, M. 2013, ApJ, 775, 116

Reines, A. E., & Volonteri, M. 2015, ApJ, 813, 82

Rhoads, J. E., Wold, I. G. B., Harish, S., et al. 2017, ApJ, 942, L14

Ricci, C., Bauer, F. E., Treister, E., et al. 2017, MNRAS, 468, 1273

Rieke, G. H. 2006, The Last of the Great Observatories: Spitzer and the Era of Faster, Better, Cheaper at NASA (3rd edn.; Univ. Arizona Press)

Rieke, G. H., Alberts, S., Shivaei, I., et al. 2024, ApJ, 975, 83

Rieke, G. H., Young, E. T., Engelbracht, C. W., et al. 2004, ApJS, 154, 25

Rieke, M., Robertson, B., Tacchella, S., et al. 2023, Data from the JWST Advanced Deep Extragalactic Survey (JADES), STScI/MAST,

Rieke, M. J., Kelly, D. M., Misselt, K., et al. 2023a, PASP, 135, 028001

Rieke, M. J., Robertson, B., Tacchella, S., et al. 2023b, ApJS, 269, 16

Rinaldi, P., Bonaventura, N., Rieke, G. H., et al. 2025a, ApJ, 992, 71

Rinaldi, P., Caputi, K. I., Costantin, L., et al. 2023, ApJ, 952, 143

Rinaldi, P., Caputi, K. I., Iani, E., et al. 2024, ApJ, 969, 12

Rinaldi, P., Navarro-Carrera, R., Caputi, K. I., et al. 2025b, ApJ, 981, 161

Rinaldi, P., Pérez-González, P. G., Rieke, G. H., et al. 2022, ApJ, 930, 128

Roberts-Borsani, G., Treu, T., Shapley, A., et al. 2024, ApJ, 976, 193

Robertson, B. E., Furlanetto, S. R., Schneider, E., et al. 2013, ApJ, 768, 71

Sanders, D. B., & Mirabel, I. F. 1996, ARA&A, 34, 749
Sanders, R. L., Shapley, A. E., Topping, M. W., Reddy, N. A., & Brammer, G. B. 2024, ApJ, 962, 24

Sanders, R. L., Shapley, A. E., Jones, T., et al. 2021, ApJ, 914, 19

Sandles, L., D'Eugenio, F., Maiolino, R., et al. 2024, A&A, 691, A305

Scarlata, C., Hayes, M., Panagia, N., et al. 2024, arXiv:2404.09015

Scholtz, J., Maiolino, R., D'Eugenio, F., et al. 2025, A&A, 697, 175

Seyfert, C. K. 1943, ApJ, 97, 28

Shapley, A. E., Sanders, R. L., Topping, M. W., et al. 2025, ApJ, 980, 242

Shivaei, I., Reddy, N., Rieke, G., et al. 2020, ApJ, 899, 117

Siebenmorgen, R., Heymann, F., & Efstathiou, A. 2015, A&A, 583, A120

Simmonds, C., Tacchella, S., Hainline, K., et al. 2024a, MNRAS, 535, 2998

Simmonds, C., Tacchella, S., Hainline, K., et al. 2024b, MNRAS, 527, 6139

Simons, R. C., Papovich, C., Momcheva, I. G., et al. 2023, ApJS, 266, 13

Skillman, E. D., Salzer, J. J., Berg, D. A., et al. 2013, AJ, 146, 3

Sommovigo, L., Ferrara, A., Pallottini, A., et al. 2020, MNRAS, 497, 956

Stalevski, M., Fritz, J., Baes, M., Nakos, T., & Popović, L. Č. 2012, MNRAS, 420, 2756

Stalevski, M., Ricci, C., Ueda, Y., et al. 2016, MNRAS, 458, 2288

Steidel, C. C., Giavalisco, M., Dickinson, M., & Adelberger, K. L. 1996, AJ, 112, 352

Stepney, M., Banerji, M., Tang, S., et al. 2024, MNRAS, 533, 2948

Stern, D., Chary, R.-R., Eisenhardt, P. R. M., & Moustakas, L. A. 2006, AJ, 132, 1405

Stern, J., & Laor, A. 2012, MNRAS, 423, 600

Strait, V., Bradač, M., Coe, D., et al. 2020, ApJ, 888, 124

Sérsic, J. L. 1968, Atlas de Galaxies Australes (Observatorio Astronomico)

Tang, M., Stark, D. P., Chen, Z., et al. 2023, MNRAS, 526, 1657

Taylor, M. 2022, in ASP Conf. Ser. 532, Astronomical Data Analysis Software and Systems, ed. J. E. Ruiz, F. Pierfedereci, & P. Teuben (ASP), 3

Temple, M. J., Banerji, M., Hewett, P. C., Rankine, A. L., & Richards, G. T. 2021, MNRAS, 501, 3061

The pandas development team 2024, pandas-dev/pandas: Pandas, v2.2.3, Zenodo

Theios, R. L., Steidel, C. C., Strom, A. L., et al. 2019, ApJ, 871, 128

Topping, M. W., Sanders, R. L., Shapley, A. E., et al. 2025, MNRAS, 541, 1707

Topping, M. W., Stark, D. P., Senchyna, P., et al. 2024, MNRAS, 529, 3301

Tripodi, R., D'Eugenio, F., Maiolino, R., et al. 2024, A&A, 692, A184

Trouille, L., Barger, A., & Tremonti, C. 2011, AAS Meeting, 217, 326.03

Trump, J. R., Arrabal Haro, P., Simons, R. C., et al. 2023, ApJ, 945, 35

Übler, H., D'Eugenio, F., Perna, M., et al. 2024a, MNRAS, 533, 4287

Übler, H., Maiolino, R., Curtis-Lake, E., et al. 2023, A&A, 677, A145

Übler, H., Maiolino, R., Pérez-González, P. G., et al. 2024b, MNRAS, 531, 355

Ucci, G., Dayal, P., Hutter, A., et al. 2023, MNRAS, 518, 3557

Vanden Berk, D. E., Richards, G. T., Bauer, A., et al. 2001, AJ, 122, 549

Veilleux, S., Kim, D. C., & Sanders, D. B. 1999, ApJ, 522, 113

Wang, B., de Graaff, A., Davies, R. L., et al. 2025, ApJ, 984, 121

Wang, B., Fujimoto, S., Labbé, I., et al. 2023, ApJL, 957, L34

Weedman, D. W., Soifer, B. T., Hao, L., et al. 2006, ApJ, 651, 101

Whitaker, K. E., Ashas, M., Illingworth, G., et al. 2019, ApJS, 244, 16

Williams, C., Tacchella, S., Maseda, M., et al. 2023, Data from the JWST Extragalactic Medium-band Survey (JEMS), STScI/MAST,

Williams, C. C., Alberts, S., Ji, Z., et al. 2024, ApJ, 968, 34

Williams, C. C., Tacchella, S., Maseda, M. V., et al. 2023, ApJS, 268, 64

Wilson, G., Huang, J. S., Fazio, G. G., et al. 2007, ApJ, 660, L59

Wilson, G., Huang, J. S., Pérez-González, P. G., et al. 2004, ApJS, 154, 107

Witstok, J., Smit, R., Maiolino, R., et al. 2021, MNRAS, 508, 1686

Wright, E. L. 2010, AAS Meeting, 216, 104.01

Wright, G. S., Rieke, G. H., Glasse, A., et al. 2023, PASP, 135, 048003

Yan, H., Dickinson, M., Eisenhardt, P. R. M., et al. 2004, ApJ, 616, 63

Yan, H., Dickinson, M., Giavalisco, M., et al. 2006, ApJ, 651, 24

Yang, G., Caputi, K. I., Papovich, C., et al. 2023, ApJ, 950, L5

Yang, H., Malhotra, S., Gronke, M., et al. 2017b, ApJ, 844, 171

Yang, H., Malhotra, S., Rhoads, J. E., & Wang, J. 2017a, ApJ, 847, 38

York, D. G., Adelman, J., Anderson, J. E., Jr., et al. 2000, AJ, 120, 1579

Zavala, J. A., Castellano, M., Akins, H. B., et al. 2025, NatAs, 9, 155

Zeimann, G. R., Ciardullo, R., Gebhardt, H., et al. 2015, ApJ, 798, 29

Zhang, X. 2024, ApJ, 960, 108