11institutetext: Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126
Pisa, Italy 22institutetext: INAF - Osservatorio Astrofisico di Arcetri, Largo
E. Fermi 5, I-50125 Firenze, Italy 33institutetext: European Southern
Observatory, Karl-Schwarzschild-Strasse 2, 85748, Garching, Germany
44institutetext: INAF - Osservatorio di Astrofisica e Scienza dello Spazio,
via Gobetti 93/3, I-40129, Bologna, Italy 55institutetext: Centro de
Astrobiología (CAB), CSIC-INTA, Cra. de Ajalvir Km. 4, 28850, Torrejón de
Ardoz, Madrid, Spain 66institutetext: Department of Physics, University of
Oxford, Denys Wilkinson Building, Keble Road, Oxford, OX1 3RH, UK
77institutetext: Cosmic Dawn Center (DAWN), Jagtvej 128, 2200, Copenhagen N,
Denmark 88institutetext: Niels Bohr Institute, University of Copenhagen,
Jagtvej 128, 2200, Copenhagen N, Denmark 99institutetext: Kavli Institute for
Cosmology, University of Cambridge, Madingley Road, Cambridge, CB3 0HA, UK
1010institutetext: Cavendish Laboratory - Astrophysics Group, University of
Cambridge, 19 JJ Thomson Avenue, Cambridge, CB3 0HE, UK 1111institutetext:
Department of Physics and Astronomy, University College London, Gower Street,
London WC1E 6BT, UK 1212institutetext: INAF - Osservatorio Astronomico di
Roma, via di Frascati 33, I-00078 Monte Porzio Catone, Italy
1313institutetext: Space Telescope Science Institute, 3700 San Martin Drive,
Baltimore, MD 21218, USA 1414institutetext: Dipartimento di Fisica e
Astronomia, Università degli Studi di Firenze, Largo E. Fermi 1, 50125,
Firenze, Italy

# Gas-phase metallicity gradients in galaxies at $z\sim 6-8$

G. Venturi<EMAIL_ADDRESS>metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$ S. Carniani , Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$ E. Parlanti Gas-phase metallicity gradients in galaxies at $z\sim
6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$ M. Kohandel
Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$ M. Curti Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$ A. Pallottini Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$ L. Vallini Gas-phase metallicity gradients in galaxies
at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$ S.
Arribas Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$ A. J. Bunker Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$ A. J. Cameron Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$ M. Castellano Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$ A. Ferrara Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-
phase metallicity gradients in galaxies at $z\sim 6-8$ A. Fontana Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$ S. Gallerani Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$ V. Gelli , Gas-phase metallicity gradients in galaxies
at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-
phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$ R. Maiolino ,, Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$ E. Ntormousi Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$ C. Pacifici Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$ L. Pentericci Gas-phase metallicity gradients in galaxies at $z\sim
6-8$Gas-phase metallicity gradients in galaxies at $z\sim 6-8$ S. Salvadori ,
Gas-phase metallicity gradients in galaxies at $z\sim 6-8$Gas-phase
metallicity gradients in galaxies at $z\sim 6-8$Gas-phase metallicity
gradients in galaxies at $z\sim 6-8$Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$ E. Vanzella Gas-phase metallicity gradients in
galaxies at $z\sim 6-8$Gas-phase metallicity gradients in galaxies at $z\sim
6-8$

(Received ; accepted )

The study of gas-phase metallicity and its spatial distribution at high
redshift is crucial to understand the processes that shaped the growth and
evolution of galaxies in the early Universe. Here we study the spatially
resolved metallicity in three systems at $z$ $\sim$ 6–8, namely A2744-YD4,
BDF-3299, and COSMOS24108, with JWST NIRSpec IFU low-resolution ($R$ $\sim$
100) spectroscopic observations. These are among the highest-$z$ sources in
which metallicity gradients have been probed so far. Each of these systems
hosts several spatial components in the process of merging within a few
kiloparsecs, identified from the rest-frame UV and optical stellar continuum
and ionised gas emission line maps. The sources have heterogeneous properties,
with stellar masses log($M_{*}/M_{\odot}$) $\sim$ 7.6–9.3, star formation
rates (SFRs) $\sim$ 1–15 $M_{\odot}$ yr-1, and gas-phase metallicities
12+log(O/H) $\sim$ 7.7–8.3, which exhibit a large scatter within each system.
Their properties are generally consistent with those of the highest-redshift
samples to date ($z$ $\sim$ 3–10), though the sources in A2744-YD4 and
COSMOS24108 are at the high end of the mass-metallicity relation (MZR) defined
by the $z$ $\sim$ 3–10 sources. Moreover, the targets in this work follow the
predicted slope of the MZR at $z$ $\sim$ 6–8 from most cosmological
simulations. The gas-phase metallicity gradients are consistent with being
flat in the main sources of each system. Flat metallicity gradients are
thought to arise from gas mixing processes on galaxy scales, such as mergers
or galactic outflows and SN winds driven by intense stellar feedback, which
wash out any gradient formed in the galaxy. The existence of flat gradients at
$z$ $\sim$ 6–8 sets also important constraints on cosmological simulations,
whose predictions on the cosmic evolution of metallicity gradients often
differ significantly, especially at high redshift.

###### Key Words.:

Galaxies: high-redshift – Galaxies: abundances – Galaxies: ISM – Galaxies:
evolution – Techniques: imaging spectroscopy – Techniques: high angular
resolution

##  1 Introduction

Figure 1: Simple sketch summarising the three main classes of radial gas-phase
metallicity gradients and possible physical origin. Negative (decreasing)
gradients can be interpreted as the result of the inside-out growth of a
galaxy, in which star formation and the following chemical enrichment start
earlier in the inner regions. Flat gradients can result from radial mixing
processes, such as galaxy mergers as well as SN winds redistributing the
metals in the ISM and galactic outflows of metal-enriched gas evacuated (or
re-accreted) in the outer regions. Positive (increasing) “inverted” gradients
may arise from accretion of pristine gas to the central regions of the galaxy;
this can also lead to flattened gradients.

Gas-phase metallicity, its scaling relations with other galactic properties,
and its spatial distribution are fundamental tools to study and understand the
evolution of galaxies (e.g. Maiolino & Mannucci, 2019). The interplay between
gas accretion, star formation, outflows, and mergers regulates the growth of
galaxies and the build up of their metal content, leading to the correlation
between gas-phase metallicity and stellar mass in galaxies known as mass-
metallicity relation (MZR; e.g. Davé et al. 2011; Pallottini et al. 2014;
Somerville & Davé 2015). The spatial distribution of metals in galaxies,
usually described through radial metallicity gradients, bears the imprint of
these underlying processes. Therefore, spatially resolved studies of
metallicity in the early Universe, where the formation and primordial growth
of galaxies were taking place, are key to better constrain the main mechanisms
that have contributed to shape them.

Negative (radially decreasing) gradients are usually interpreted as resulting
from an inside-out galaxy formation scenario, in which stars start forming
earlier in the inner parts of galaxies and thus have more time to chemically
enrich the inner regions than the outer ones (e.g. Samland et al., 1997;
Prantzos & Boissier, 2000; Davé et al., 2011; Pilkington et al., 2012; Gibson
et al., 2013; Hemler et al., 2021; Tissera et al., 2022). Flattened gradients
could arise from radial mixing of gas and redistribution of metals on galaxy
scales, induced by SN winds within the galaxy or galactic outflows of metal-
enriched material expelled from the galaxy and re-accreted in the outer
regions (so-called “galactic fountains”), driven by intense stellar feedback
(e.g. Gibson et al., 2013; Ma et al., 2017), as well as by galaxy merging and
interactions (e.g. Rupke et al., 2010b, a; Rich et al., 2012; Torres-Flores et
al., 2014). Finally, positive (radially increasing; also called “inverted”)
gradients may be produced by accretion of external pristine gas towards the
central regions of the galaxy (e.g. Ceverino et al., 2016); this could also
contribute to flatten the metallicity gradient, if not invert it. Particularly
strong metal-loaded galactic outflows could also contribute to produce a
positive gradient, by moving the metal-rich gas from the central starburst
regions to the outskirts (e.g. Tissera et al., 2022). A schematic cartoon that
summarises the above framework is shown in Fig. 1.

In the local Universe, most spiral galaxies show negative metallicity
gradients, with more metal-enriched gas in the inner than in the outer regions
of the galaxy (e.g. Magrini et al., 2009, 2017; Stanghellini & Haywood, 2010;
Luck et al., 2011 for the Milky Way; Zaritsky et al., 1994; Magrini et al.,
2010; Kewley et al., 2010; Bresolin, 2011; Sánchez et al., 2014; Berg et al.,
2015; Ho et al., 2015; Belfiore et al., 2017 for other local galaxies). By
using metallicity diagnostics that trace past metal enrichment, primarily
planetary nebulae (PNe; up to 5–10 Gyr ago, i.e. $z$ $\lesssim$ 2) some
studies have suggested that the Milky Way and nearby galaxies had flatter
gradients in the past which grew progressively steeper with cosmic time
(Stanghellini et al., 2014; Stanghellini & Haywood, 2018; Magrini et al.,
2016). A dependence of gradient slope with mass has been observed in the local
Universe, with more massive galaxies exhibiting progressively steeper
gradients, while low-mass galaxies having almost flat gradients (Belfiore et
al., 2017). This relation between mass and metallicity gradient may be the
result of an evolutionary sequence in mass, with more massive, more evolved
galaxies having steeper gradients and low-mass, less evolved ones, analogues
of high-$z$ galaxies, having flatter gradients (e.g. Maiolino & Mannucci,
2019).

Cosmological simulations of galaxy evolution make very different predictions
on metallicity gradients at high redshift (up to $z$ $\sim$ 3), which can
strongly vary depending on prescriptions on star formation and stellar
feedback models and on their relative contribution to gas mixing as compared
to mergers (e.g. Gibson et al., 2013; Ma et al., 2017; Hemler et al., 2021;
Tissera et al., 2022). For example, according to the FIRE simulations (Hopkins
et al., 2014), mergers and rapid variations of metallicity gradients in
high-$z$ galaxies, induced by starburst episodes which drive strong outflows
and disrupt the gas disc, would flatten any negative gradients previously
developed, especially in low-mass galaxies where feedback mechanisms are more
efficient (Ma et al., 2017). Instead, the TNG50 simulations (Pillepich et al.,
2019; Nelson et al., 2019) which feature less bursty feedback, predict
negative gradients at all redshifts, steeper (more negative) with increasing
$z$ (Hemler et al., 2021). The metallicity gradients resulting from the MUGS
(“normal” feedback; Stinson et al., 2010) and MaGICC (“enhanced” feedback;
Brook et al., 2012) simulation suites are compared in Gibson et al. (2013).
The “enhanced” feedback (including pre-SN early stellar feedback from massive
stars) re-distributes energy and re-cycled ISM material over large scales,
through re-accretion in the outer parts of the galaxy of gas expelled via
outflows. In this case, relatively flat and temporally invariant abundance
gradients are predicted, in contrast to the steeper negative gradients at
increasing $z$ emerging from the normal feedback scenario. Finally, the EAGLE
simulations (Schaye et al., 2015) find that the median gradient is zero (flat)
at all redshifts, but with the scatter around the median increasing with $z$
due to individual galaxies transitioning between steep (either negative or
positive) and flatter gradients (Tissera et al., 2022). The authors report
that this behaviour results from the higher frequency of major mergers at
higher $z$ leading to episodes of enhanced accretion of low-metallicity gas
which trigger intense star formation and ejection of metal-enriched gas.

Observations at the cosmic noon epoch (i.e., $z$ $\sim$ 1–4) have found
heterogeneous results, with both negative, positive, and flat metallicity
gradients (e.g. Cresci et al., 2010; Yuan et al., 2011; Queyrel et al., 2012;
Swinbank et al., 2012; Stott et al., 2014; Troncoso et al., 2014; Jones et
al., 2013, 2015; Leethochawalit et al., 2016; Wuyts et al., 2016; Wang et al.,
2017, 2022; Carton et al., 2018; Förster Schreiber et al., 2018; Curti et al.,
2020a). Nevertheless, the majority of high-$z$ measurements is consistent with
little or no cosmic evolution of metallicity gradients, which are found to be
flat or only moderately negative or positive ($\lesssim$ $|\pm 0.1|$ dex
kpc-1; see e.g. the compilation in Curti et al. 2020a).

At higher redshift ($z$ $\sim$ 7), Vallini et al. (2024) also found flat
gradients by making use of a physically motivated Bayesian model to derive
metallicities from rest-frame far-infrared lines ([O iii] 88$\mu$m and [C ii]
158$\mu$m) from Atacama Large Millimeter/submillimeter Array (ALMA)
observations. However, as the authors point out as a caveat, there may be
systematics in comparing their results with gradients of O/H abundance
obtained with the more standard methods based on optical lines, since the
adopted model does not take into account the likely enhancement of O/C at low
metallicities (see e.g. Maiolino & Mannucci, 2019), which would result in the
model returning higher values of metallicity. Therefore, while their results
consistent with flat gradients at $z$ $\sim$ 7 are very significant, they
require independent confirmation.

The James Webb Space Telescope (JWST) has recently opened up the possibility
of measuring metallicity gradients in the first 1–2 Gyr of the Universe by
using rest-frame optical emission lines (Wang et al., 2022; Rodríguez Del Pino
et al., 2023; Arribas et al., 2023) which are extensively used and calibrated
at lower redshifts (e.g. Curti et al., 2017, 2020b).

In this work, we study the spatially resolved gas-phase metallicity and
investigate the shape of metallicity gradients at high redshift ($z$ $\gtrsim$
6), by tracing the warm ($T$ $\sim$ $10^{4}$ K) ionised gas emission.
Specifically, we observed three high-$z$ systems with the Near-InfraRed
Spectrograph (NIRSpec) on board JWST in its Integral Field Unit (IFU) mode
with the low-spectral resolution PRISM/CLEAR ($R$ $\sim$ 100) disperser-filter
combination. These are A2744-YD4, part of the proto-cluster A2744-z7p9OD ($z$
$\sim$ 7.88), BDF-3299 ($z$ $\sim$ 7.11), and COSMOS24108 ($z$ $\sim$ 6.36).
Basic information on the targets is given in Table 1. We trace gas-phase
metallicity by making use of strong-line calibrators relying on rest-frame
optical and near-UV emission lines. We also present the basic properties of
each source, specifically stellar mass, obtained from SED fitting, and SFR,
from H$\beta$ or H$\alpha$. The star formation history of these targets is
studied in detail in a separate paper (Kohandel et al., in prep.).

Throughout this work, the reported wavelengths are in vacuum and we adopt a
flat $\Lambda$CDM cosmology with $H_{0}$ $\simeq$ 67.7 km s-1 Mpc-1,
$\Omega_{\mathrm{M}}$ $\simeq$ 0.31, and $\Omega_{\mathrm{\Lambda}}$ $\simeq$
0.69 (Planck Collaboration et al., 2020).

##  2 Description of observations, data reduction, and analysis

The observations employed in this work were carried out between October and
December 2022 as part of the JWST GO program ID 1893 (PI: S. Carniani;
Carniani et al. 2021). The three sources were observed for 15.102 ks ($\sim$4
h) on-source each (22.059 ks with overheads; $\sim$6 h) with 1.678 ks of
background exposure time each (2.830 ks with overheads). The observations were
performed with the PRISM/CLEAR disperser-filter combination and an eight-point
dither “medium” cycling pattern (step $\sim$0.5′′). The JWST NIRSpec
$\sim$3${}^{\prime\prime}\times 3^{\prime\prime}$ IFU PRISM observations
simultaneously span the spectral range 0.6–5.2 $\mu$m with a spectral
resolution ranging between $R$ $\sim$ 30–330.

Table 1: Basic information on the observed targets.

Target name | RA [hh:mm:ss.ss] | Dec [dd:mm:ss.ss] |  $z$a | Scale [kpc arcsec-1]
---|---|---|---|---
A2744-YD4 | +00:14:24.90 | –30:22:56.10 | 7.879 | 3.512b
BDF-3299 | +22:28:12.31 | –35:10:00.39 | 7.114 | 5.289
COSMOS24108 | +10:00:47.33 | +02:28:43.14 | 6.361 | 5.647

  * a

From this work.

  * b

The reported scale is corrected for the lensing magnification factor of
$\sim$2 (Morishita et al., 2023; Bergamini et al., 2023).

###  2.1 Data reduction

We retrieved the raw data from the MAST archive and we ran the three stages of
the pipeline using version 1.11 with CRDS (calibration reference data system)
context “jwst_1094.pmap”. First, at stage 1 “calwebb_detector1”, the pipeline
applied the detector-level corrections (e.g. check for saturation, dark
exposure subtraction, flagging of bad pixels and cosmic-ray persistences) and
performed ramp fitting for individual exposures. We then calibrated the count
rate images by executing stage 2 “calwebb_spec2” of the pipeline, which
corrects for flat field and performs the wavelength calibration. The
background was subtracted from each exposure during this stage by using the
observations of the dedicated background for each target. We processed the
background targets up to stage 2 of the pipeline, then we applied the
background step for the science target. In the background step, the pipeline
subtract the background exposure from each target exposure in the detector
space. Finally, each calibrated exposure was combined in stage 3
“calwebb_spec3” by using “drizzle” weighting and a spaxel size of 0.05″ to
obtain the final cubes. During stage 3, we applied the outlier rejection step
built in the pipeline and then a sigma clipping to remove any residual
outliers in the final cube.

###  2.2 Data analysis

Table 2: List of emission line diagnostic ratios used in this work together with their compact notation. Name | Emission line ratio
---|---
R2 | [O ii] $\lambda\lambda$3727,30/H$\beta$
R3 | [O iii] $\lambda$5008/H$\beta$
O3O2 | [O iii] $\lambda$5008/[O ii] $\lambda\lambda$3727,30
$\hat{\mathrm{R}}$ | 0.47 R2 \+ 0.88 R3
Ne3O2 | [Ne iii] $\lambda$3870/[O ii] $\lambda\lambda$3727,30
RS32 | [O iii] $\lambda$5008/H$\beta$ \+ [S ii] $\lambda\lambda$6718,33/H$\alpha$
O3S2 | [O iii] $\lambda$5008/H$\beta$ / [S ii] $\lambda\lambda$6718,33/H$\alpha$

In this section we describe the analysis of the NIRSpec R100 IFU data. In
brief, we obtained the emission line fluxes from emission-line modelling of
the spectra, which were used to infer the gas-phase metallicity. This was done
on both an integrated basis (to get the integrated metallicity of each target)
and a spatially resolved one, in this case both spaxel-by-spaxel (to produce
maps) and in concentric radial annuli (to obtain radial gradients). From the
integrated spectra, we also obtained the stellar mass ($M_{*}$) of each source
from SED fitting and their star formation rate (SFR) from the emission-line
modelling. We provide more details in the following.

####  2.2.1 Emission line fitting

The main goal of this work is to measure radial gas-phase metallicity
gradients. To do so, we infer the oxygen abundance relative to hydrogen
(12+log(O/H)), a proxy of gas-phase metallicity, by making use of the optical
and near-UV strong-line diagnostic ratios reported in Table 2. We adopt the
new diagnostic ratio $\hat{\mathrm{R}}$ = 0.47 R2 \+ 0.88 R3, first introduced
in Laseter et al. (2024), in place of the more traditional R23 = ([O ii]
$\lambda\lambda$3727,30 + [O iii] $\lambda\lambda$4960,5008)/H$\beta$, since
the former is more suited for high-$z$ galaxies, while the latter is a
projection mostly driven by low-metallicity local analogues. We adopt the
best-fit polynomial calibrations from Curti et al. (2017, 2020b), slightly
revisited in Curti et al. (2023a, b) to better probe the low-O/H regime, to
infer the gas metallicity from the combination of the above ratios. We tested
other gas-phase metallicity calibrations (e.g. Nakajima et al., 2022; Sanders
et al., 2024). These gave 12+log(O/H) values similar to the Curti et al. ones
within 0.1 dex, consistent with the uncertainties, and virtually no difference
in the metallicity radial gradients.

The FWHM of the point spread function (PSF) of JWST depends on wavelength.
Specifically, the PSF FWHM of NIRSpec IFU ranges from around 0.09′′ at 1
$\mu$m up to around 0.16′′ at 4.5 $\mu$m, and the variation with wavelength is
stronger in the direction perpendicular to the IFU slices than along them
(D’Eugenio et al., 2023). Since the aim of this work is to obtain spatially
resolved metallicities from emission line ratios sometimes far in wavelength
(see Table 2), we applied a wavelength-dependent spatial smoothing to the data
cube prior to the line fitting with the aim of achieving the same spatial
resolution at all wavelengths. We adopted the PSF FWHM curves from D’Eugenio
et al. (2023), empirically calibrated by matching NIRSpec IFU observations
with NIRCam ones of the same target reported in their Eqs. 3 and 4, for the
along- and across-slice cases, respectively. The smoothing was done by
convolving the cube with a wavelength-dependent 2D Gaussian kernel in order to
obtain a smoothed cube having the same spatial resolution at each wavelength,
specifically the PSF FWHM reported in D’Eugenio et al. (2023) at the
wavelength of the highest-wavelength line of interest for the along-slice case
(which has the largest FWHM among the two cases). This corresponds to FWHMgoal
$\sim$ 0.16′′ (at [O iii]) for A2744-YD4, $\sim$ 0.14′′ (at [O iii]) for
BDF-3299, and $\sim$ 0.18′′ (at H$\alpha$-[S ii]) for COSMOS24108.
Specifically, the 2D Gaussian kernel to be used for the convolution, at a
given wavelength $\lambda$, was defined as follows:
FWHM${}^{2}_{\mathrm{kern;dir}}$($\lambda$) =
FWHM${}^{2}_{\mathrm{goal}}$($\lambda$) –
FWHM${}^{2}_{\mathrm{mod;dir}}$($\lambda$), where FWHMmod;dir is the model
FWHM from D’Eugenio et al. (2023) and “dir” indicates the direction, either
along or across the slices.

In general, we employ the original unsmoothed data cube for the flux maps of
continuum and emission lines and for integrated measurements from circular
apertures, while we adopt the smoothed data cube for the maps and radial
profiles involving line ratios and metallicity, which would otherwise be
affected by the wavelength-dependence of the PSF. In case of A2744-YD4, the
highest-$z$ target of our sample, we adopt the smoothed data cube for every
map, to obtain visually clearer maps as compared to the more noisy ones from
the unsmoothed cube.

The data analysis consisted in modelling the rest-frame optical and near-UV
emission lines available in the observed spectral range. The emission lines
are spectrally unresolved in the PRISM spectra (the spectral resolution is
$\sigma_{\mathrm{res}}$ $\sim$ 1250 km s-1 at $\sim$3 $\mu$m and
$\sigma_{\mathrm{res}}$ $\sim$ 450 km s-1 at $\sim$5 $\mu$m), therefore a
single Gaussian function per line was used. We fitted the [O iii]-H$\beta$,
the [O ii]-[Ne iii]-H$\gamma$, and the H$\alpha$-[N ii]-[S ii] line complexes
separately. The [N ii] 6550,85 was included to allow for a better modelling of
the H$\alpha$ line profile, which was otherwise showing a marked residual in
its redward wing. This was done only when fitting integrated spectra from
circular apertures or concentric annuli for radial profiles, not in the
lower-S/N case of spaxel-by-spaxel fitting, for which no residual wing
indicative of [N ii] was detected above the noise. For each line complex, the
velocity was tied to be the same for all the lines, while we allowed the line
width to vary to match the PRISM spectral resolution, given that all lines are
spectrally unresolved as mentioned. We fixed the flux ratios [O iii]
$\lambda$5008/$\lambda$4960 and [N ii] $\lambda$6585/$\lambda$6550 to their
theoretical value of 3 (Storey & Zeippen, 2000). We included an underlying
first-order polynomial to model the continuum. We accounted for the Balmer
break of the continuum at around 3645 Å rest-frame, when needed, by employing
three (first-order) polynomials, one modelling the jump and the other two on
each side of it.

The line fluxes used to infer the metallicity or the SFR were first corrected
for dust extinction, when possible (see below), using the H$\gamma$/H$\beta$
ratio (for A2744-YD4 and BDF-3299) or H$\alpha$/H$\beta$ (for COSMOS24108),
only when both lines had S/N$>$3\. We adopted a Calzetti et al. (2000)
reddening curve, with an $R_{V}$ = 4.05, suitable for high-$z$ low-metallicity
(12+log(O/H) $<$ 8.5) star-forming galaxies at the wavelengths of interest
here of $\gtrsim$3000 Å (see e.g. Shivaei et al. 2020). We assumed the
theoretical extinction-free H$\gamma$/H$\beta$ and H$\alpha$/H$\beta$ ratios
of 0.466 and 2.87, respectively, valid for case-B recombination and an
electron temperature of $T_{\mathrm{e}}$ $\sim$ 104 K (Osterbrock & Ferland,
2006), typical of the warm ionised gas emitting rest-frame optical lines.

We selected the line ratio diagnostics for the metallicity estimation in an
adaptive way, based on the line S/N. Specifically, a certain ratio (e.g. [Ne
iii]/[O ii]) was employed only when all the lines involved in it had S/N$>$3\.
In this way, for each fitted spectrum we only selected the sub-sample of line
ratio diagnostics whose line fluxes exceeded the S/N threshold, instead of
using all the ratios (even those involving S/N$<$3 lines) or none of them
(when just one or a few of the lines were below the S/N threshold). Finally,
given a set of diagnostic ratios, the best-fit metallicity was obtained by
minimising the chi-squared defined simultaneously by the different observed
ratios and their relative calibration curves, weighted by the observed
uncertainties and the intrinsic dispersion of the calibration (added in
quadrature; see Curti et al. 2020b).

We first modelled the emission lines in integrated spectra extracted from
circular apertures centred on each sub-source present in the three studied
systems. From these, we obtained the integrated gas-phase metallicity (as
described above) and SFR of each sub-source. The SFR was obtained from the
extinction-corrected flux of H$\alpha$, using the relation from Kennicutt &
Evans (2012). In A2744-YD4 and BDF-3299, for which H$\alpha$ was not
available, the SFR was obtained from the extinction-corrected flux of
H$\beta$, converted to H$\alpha$ by adopting the theoretical
H$\alpha$/H$\beta$ ratio of 2.87. We employed circular apertures of 0.15″
radius for A2744-YD4 and COSMOS24108, and of 0.1″ for BDF-3299, due to the
spatial vicinity of the sources in this latter system.

The emission line modelling was then performed both on a spaxel-by-spaxel
basis (for both the unsmoothed and the smoothed data cubes) and in concentric
radial annuli centred on each source (for the smoothed data cube). The latter
was done to increase the S/N on the emission lines and get more robust
estimates of metallicity with the goal of obtaining metallicity gradients. For
this annular line modelling, we extracted integrated spectra at each radius by
collapsing the spaxels within concentric circular annuli. The annuli, having
radial width of 1 spaxel (0.05′′) each, were centred on each of the main
spatial components detected in each targeted system, with a variable maximum
aperture radius depending on the component extension (as described in the next
sections separately for each target).

The extinction-correction of emission line fluxes used to infer metallicity
was done only for the case of integrated spectra from circular apertures and
from concentric radial annuli. This was not possible for the spaxel-by-spaxel
case because, among the useful Balmer lines tracing extinction, H$\gamma$ is
detected in almost no spaxels in A2744-YD4 and BDF-3299. Therefore, the maps
of metallicity were obtained without accounting for possible extinction. For
COSMOS24108, also including H$\alpha$ in the spectral range, H$\beta$ is
detected in a large enough number of spaxels to attempt for an extinction
correction of the line fluxes for the metallicity. However, the resulting map
is very similar to the non-extinction-corrected one, only more noisy,
therefore we report the metallicity map obtained without the extinction
correction.

Table 3: Properties of the individual sources. (1) Source name. (2) Redshift from [O iii] $\lambda$5008\. (3) Stellar mass from SED fitting with bagpipes. (4-9) [O ii] $\lambda\lambda$3727,30, [Ne iii] $\lambda$3870, H$\beta$, [O iii] $\lambda$5008, H$\alpha$, and [S ii] $\lambda$6718,33 measured fluxes (in units of 10-20 erg s-1 cm-2), (10) SFR (in $M_{\odot}$ yr-1), and (10) gas-phase metallicity from emission line Gaussian fitting. For the sources in the A2744-YD4 system, the reported $M_{*}$ and SFR are corrected for the lensing magnification factor of $\sim$2 (Morishita et al., 2023; Bergamini et al., 2023); the line fluxes are instead the observed ones. All the properties are obtained from integrated spectra extracted from circular apertures (having radius of 0.15″ for A2744-YD4 and COSMOS24108 and of 0.1″ for BDF-3299) centred on each source. For the line fluxes and SFR we also report the $3\sigma$ upper limits in case of no detection. Source | $z$ | $\log(M_{*}/M_{\odot})$ | $F_{\rm[O\,\textsc{ii}]}$ | $F_{\rm[Ne\,\textsc{iii}]}$ | $F_{\rm H\beta}$ | $F_{\rm[O\,\textsc{iii}]}$ | $F_{\rm H\alpha}$ | $F_{\rm[S\,\textsc{ii}]}$ | SFR | 12+log(O/H)
---|---|---|---|---|---|---|---|---|---|---
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11)
A2744-YD4 |
YD4 | 7.879 | 8.69${}^{+0.18}_{-0.17}$ | 72$\pm$3 |  $<$21 | 36$\pm$3 | 180$\pm$4 | – | – | 2.2${}^{+0.2}_{-0.2}$ | 8.19${}_{-0.13}^{+0.09}$
YD6 | 7.880 | 8.40${}^{+0.27}_{-0.13}$ | 26$\pm$3 |  $<$20 |  $<$17 | 44$\pm$3 | – | – |  $<$1.0 | 8.27${}_{-0.20}^{+0.16}$
YD6-[O iii]-E | 7.881 | 7.83${}^{+0.46}_{-0.22}$ |  $<$24 |  $<$24 |  $<$14 | 57$\pm$3 | – | – |  $<$0.87 |
YD4-[O iii]-W | 7.877 | 7.82${}^{+0.41}_{-0.48}$ |  $<$14 |  $<$15 |  $<$19 | 46$\pm 3$ | – | – |  $<$1.1 |
ZD1 | 7.872 | 7.58${}^{+0.18}_{-0.37}$ |  $<$13 |  $<$13 |  $<$18 | 27$\pm$3 | – | – |  $<$1.1 |
YD1 | 7.883 | 8.58${}^{+0.28}_{-0.19}$ | 107$\pm$4 |  $<$24 | 47$\pm$3 | 346$\pm$3 | – | – | 2.85${}^{+0.18}_{-0.18}$ | 8.11${}_{-0.07}^{+0.05}$
YD1-E | 7.882 | 8.44${}^{+0.43}_{-0.23}$ | 54$\pm$3 |  $<$20 | 34$\pm$2 | 270$\pm$2 | – | – | 2.02${}^{+0.13}_{-0.13}$ | 7.98${}_{-0.09}^{+0.07}$
s1 | 7.878 | 8.01${}^{+0.28}_{-0.47}$ |  $<$24 |  $<$25 | 18$\pm$3 | 151$\pm$3 | – | – | 1.09${}^{+0.16}_{-0.16}$ | 7.72${}_{-0.20}^{+0.20}$
BDF-3299 |
a | 7.114 | 7.90${}^{+0.25}_{-0.13}$ | 24$\pm$2 |  $<$13 | 36$\pm$1 | 195$\pm$3 | – | – | 3.47${}^{+0.13}_{-0.13}$ | 7.68${}_{-0.09}^{+0.13}$
b | 7.114 | 8.21${}^{+0.30}_{-0.11}$ | 28$\pm$2 |  $<$11 | 44$\pm$2 | 285$\pm$4 | – | – | 4.2${}^{+0.2}_{-0.2}$ | 7.69${}_{-0.06}^{+0.10}$
c | 7.112 | 7.62${}^{+0.44}_{-0.14}$ |  $<$13 |  $<$12 | 15$\pm$2 | 93$\pm$4 | – | – | 1.5${}^{+0.2}_{-0.2}$ | 7.80${}_{-0.38}^{+0.15}$
COSMOS24108 |
a | 6.361 | 9.29${}^{+0.09}_{-0.08}$ | 179$\pm$5 | 48$\pm$6 | 53$\pm$5 | 358$\pm$5 | 219$\pm$4 | 34$\pm$4 | 14.4${}^{+1.3}_{-1.3}$ | 8.20${}_{-0.06}^{+0.05}$
b | 6.358 | 8.81${}^{+0.11}_{-0.17}$ | 150$\pm$6 | 47$\pm$6 | 72$\pm$4 | 470$\pm$4 | 245$\pm$3 |  $<$18 | 9.8${}^{+0.3}_{-0.3}$ | 8.11${}_{-0.07}^{+0.04}$
c | 6.359 | 8.89${}^{+0.10}_{-0.06}$ | 67$\pm$5 |  $<$36 | 28$\pm$4 | 92$\pm$4 | 68$\pm$4 |  $<$22 | 1.74${}^{+0.09}_{-0.09}$ | 8.33${}_{-0.10}^{+0.08}$
[O iii]-Ea | 6.362 | 8.53${}^{+0.30}_{-0.22}$ | 62$\pm$3 | 36$\pm$3 | 41$\pm$3 | 296$\pm$3 | 183$\pm$3 |  $<$21 | 14.5${}^{+1.3}_{-1.3}$ | 8.01${}_{-0.12}^{+0.08}$
[O iii]-Eb | 6.356 | 8.40${}^{+0.72}_{-0.32}$ | 44$\pm$4 |  $<$24 | 37$\pm$3 | 244$\pm$3 | 113$\pm$3 |  $<$22 | 3.43${}^{+0.11}_{-0.11}$ | 7.95${}_{-0.17}^{+0.14}$

####  2.2.2 SED fitting

With the goal of obtaining the stellar mass of each target, we fitted the
spatially integrated spectra extracted from the data cube by using the
spectral energy distribution (SED) fitting code bagpipes (Carnall et al.,
2018). We adopted the stellar population models by Bruzual & Charlot (2003)
and included the nebular emission with cloudy (Ferland et al., 2017) with the
ionisation parameter ($-3.0<\log U<0.0$) as a free parameter and adopting a
${\rm[C/O]=[C/O]}_{\odot}$. We assumed a Kroupa (2001) initial mass function
truncated at 0.01 and 100 M⊙ and a Calzetti et al. (2000) attenuation curve.
Finally, we used a non-parametric star-formation history model with continuity
priors (see Leja et al. 2019) and with four time bins: $0<t<10$ Myr, $10<t<50$
Myr, $50<t<100$ Myr, and $100<t<300$ Myr.

##  3 A2744-YD4

Figure 2: Maps for A2744-YD4 from JWST NIRSpec IFU. All the maps are obtained
from the (wavelength-dependent) spatially smoothed data cube (details in Sect.
2.2). No extinction correction is applied. A cut of S/N$>3$ on the peak flux
of each line is applied. Median stellar continuum emission in the observed
spectral range 2–3 $\mu$m (0.23–0.34 $\mu$m rest-frame; left); [O iii]
integrated flux (centre); map of oxygen abundance, 12+log(O/H) (right).
Contours mark the continuum from first panel. NIRSpec PRISM/CLEAR spectrum
extracted from a circular aperture with radius of 0.15″ centred at the
location of the target YD4 (bottom).

A2744-YD4 is part of the proto-cluster A2744-z7p9OD at $z$ $\sim$ 7.883
located behind the strong lensing cluster Abell 2744. The proto-cluster
includes 22 sources at 7 ¡ $z$ ¡ 9 identified through combined deep HST and
Spitzer IRAC photometry (Laporte et al., 2014; Zheng et al., 2014; Atek et
al., 2015; Ishigaki et al., 2016) as part of the Hubble Frontier Fields
program (Lotz et al., 2017). Seven of them were recently spectroscopically
confirmed to be at $z$ = 7.88 through JWST NIRSpec MSA spectroscopy (Morishita
et al., 2023). This makes A2744-z7p9OD the most distant proto-cluster known so
far. The system is highly over-dense, with an excess of surface number density
from the field average $\delta$ = $(n-\overline{n})/\overline{n}$ $\sim$
130${}^{+66}_{-51}$ (Ishigaki et al., 2016).

The $\sim$3${}^{\prime\prime}\times 3^{\prime\prime}$ IFU FOV of our
observations includes five sources of the proto-cluster known as the
‘quintet’, namely YD1, YD4, YD6, YD7, and ZD1, as well as a new source
identified by Hashimoto et al. (2023) and named s1 (Fig. 2, top-left). We only
partially cover the East component of YD7 (named YD7-E), while its emission
toward to West lies outside the FOV. The lensing magnification factor of the
source in the FOV is $\mu$ $\simeq$ 2.0 (Morishita et al., 2023; Bergamini et
al., 2023). YD4, at the centre of our IFU observation, was reported to be at
$z$ = 8.38 based on Ly$\alpha$ and ALMA [C ii] 158 $\mu$m and [O iii] 88
$\mu$m line emission (Laporte et al., 2017, 2019; Carniani et al., 2020).
Recently, Morishita et al. (2023) measured $z$ = 7.88 for YD4 based on the
detection of high-S/N [O iii] 4959,5007 $\AA$ and H$\beta$ emission lines in
high-spectral resolution ($R$2700) NIRSpec MSA observations, thus ruling out
the previous redshift measurement.

The top-left panel of Fig. 2 shows the map of the continuum emission in the
range 2–3 $\mu$m observed wavelength ($\sim$0.23–0.34 $\mu$m rest-frame). In
addition to the sources named YD4, YD6, YD1, s1, YD7-E, and ZD1, the continuum
map also shows the presence of a source at the northernmost part of the FOV,
which is also present in the NIRCam image but is not labelled in Hashimoto et
al. (2023); we name this source s2. However, given that this is at the edge of
the FOV, where many artefacts are present in the NIRSpec IFU data, we consider
this source as tentative and report its label in parenthesis in the figure.
Moreover, in the map of the continuum obtained from the unsmoothed data cube
integrated over the spectral range 1.2–2 $\mu$m (Fig. 10, top-left), where the
spatial resolution is the highest, the eastern tail of YD1 appears as a
separate peaked spatial component. We label this extra source as YD1-E.

Fig. 2, bottom panel, illustrates the integrated spectrum (aperture radius =
0.15′′) associated with the main component of the A2744-YD4 system, that is,
YD4 itself. The rest-frame optical lines [O ii], H$\beta$, and [O iii] are
detected with high S/N ($\gtrsim$6 on the peak). The spectroscopic redshift
based on Lyman break and optical lines ($z$ $\sim$ 7.88) is consistent with
that found by Morishita et al. (2023) and Hashimoto et al. (2023) from the
NIRSpec $R$2700 data. We find a stellar mass of log($M_{*}/M_{\odot}$) $\sim$
8.7 from the SED fitting and a SFR $\sim$ 2 $M_{\odot}$ yr-1 from H$\beta$ for
YD4 (Table 3). These were corrected for the magnification factor of 2 due to
the lensing. The gas-phase metallicity, relying on line ratios, is instead not
affected by it. For the other spatial components in the system, we find
log($M_{*}/M_{\odot}$) $\sim$ 7.5–8.6 and SFR $\sim$ 1–3 $M_{\odot}$ yr-1.

The [O iii] flux map resulting from our spectral emission-line modelling is
displayed in Fig. 2, central panel (the maps of H$\beta$ and [O ii] are shown
in Fig. 10). From the [O iii] map, we identify two new emission regions which
are weak or absent in continuum, which we label as YD6-[OIII]-E and
YD4-[OIII]-W. Overall, the [O iii] maps obtained from the PRISM data presented
here are deeper and reveal fainter and more extended features than those from
high-resolution grating data presented in Hashimoto et al. (2023) (who do not
report the detection of any other emission line due to the lower S/N of their
data). We can count nine sources in total in the
$\sim$3${}^{\prime\prime}\times 3^{\prime\prime}$ NIRSpec IFU FOV, among those
detected in continuum and in line emission, namely YD4, YD6, YD6-[OIII]-E,
YD4-[OIII]-W, YD7-E, ZD1, YD1, YD1-E, and s1, and possibly a tenth source, s2.
All these sources lie at $z$ $\sim$ 7.88, either confirmed spectroscopically
(YD4, YD6, YD6-[OIII]-E, YD4-[OIII]-W, YD1, YD1-E, and s1; e.g. Morishita et
al. 2023, Hashimoto et al. 2023, and this work) or from photometry (ZD1),
except for YD7-E (and the tentative source s2) whose redshift is not assessed.
The spectra of all targets, but YD4, are reported in Figs. 8 and 9.

In the top-right panel of Fig. 2 we show the map of oxygen abundance, 12 +
$\log$ (O/H), a proxy for gas-phase metallicity. The metallicity is inferred
by making use of the strong-line diagnostics reported in Table 2 as described
in Sect. 2.2 (the maps of the [O iii] $\lambda$5008/H$\beta$ and [O iii]
$\lambda$5008/[O ii] $\lambda\lambda$3727,30 emission line ratios are shown in
Fig. 10). The ratios involving lines which are not close in wavelength (all
but [O iii] and H$\beta$) may be affected by extinction. Unfortunately, the
only robustly detected Balmer line at the spaxel level is H$\beta$, and only
in a handful of spaxels, therefore estimating the dust extinction spaxel-by-
spaxel is not possible. Based on this, the abundance map should not be taken
as a robust measurement; moreover, the map is very noisy. Nevertheless, the
map seems to suggest that, in the northern system, YD1 is embedded between two
sources at lower metallicity, s1 and YD1-E. No clear pattern of metallicity
shows up in the southern system (the YD4 one) from the map.

##  4 BDF-3299

Figure 3: Maps for BDF-3299 from JWST NIRSpec IFU. The flux maps are obtained
from the original data cube, while the line ratio and metallicity maps from
the spatially smoothed data cube (details in Sect. 2.2). Median stellar
continuum emission in the observed spectral range 1.1–2 $\mu$m
($\sim$0.14–0.25 $\mu$m rest-frame; top-left). The rest is as in Fig. 2. The
PRISM/CLEAR spectrum (bottom) is extracted from a circular aperture with
radius of 0.1″ centred at the location of the most luminous component in
continuum, BDF-3299-a.

BDF-3299 is a spectroscopically confirmed star-forming galaxy (SFRUV $\sim$ 6
$M_{\odot}$ yr-1, SFRdust $\lesssim$ 12 $M_{\odot}$/yr) at $z$ $\sim$ 7.109
(Castellano et al., 2010; Vanzella et al., 2011; Maiolino et al., 2015;
Carniani et al., 2017), located in an over-density of galaxies (Castellano et
al., 2016). We report the maps for BDF-3299 from the NIRSpec IFU observations
in Fig. 3.

The median continuum in the observed spectral range 1.1–2 $\mu$m
($\sim$0.14–0.25 $\mu$m rest-frame) shows that BDF-3299 is composed of three
spatial components. The two brighter ones are located at the NW (the strongest
of the two in the continuum) and at the centre of the system, respectively; a
third, fainter one resides in the E part. We label these sources as
BDF-3299-a, BDF-3299-b, and BDF-3299-c, respectively, in decreasing order of
continuum brightness, as labelled in Fig. 3, left panel. The other continuum
emitters in the FOV are lower-redshift sources.

Fig. 3, bottom panel, shows the integrated spectrum (aperture radius = 0.1′′)
associated with the brightest component in continuum emission, BDF-3299-a. The
rest-frame optical lines [O ii], H$\beta$, and [O iii] are detected with high
S/N ($\gtrsim$6 on the peak). From the lines, we estimate a redshift of 7.114,
slightly higher but roughly consistent with the previously reported one from
Ly$\alpha$ (7.109; Vanzella et al., 2011) and [C ii] (7.107; Carniani et al.,
2017). We obtain a stellar mass of log($M_{*}/M_{\odot}$) $\sim$ 7.9 from the
SED fitting and a SFR $\sim$ 3.5 $M_{\odot}$ yr-1 from H$\beta$ for BDF-3299-a
(Table 3). The spectra of the other spatial components in the system, namely
BDF-3299-b and c, are reported in Fig. 11. From these, we find
log($M_{*}/M_{\odot}$) $\sim$ 7.6–8.2 and SFR $\sim$ 1.5–4 $M_{\odot}$ yr-1.

The [O iii] emission (Fig. 3, central panel) peaks on BDF-3299-b rather than
on BDF-3299-a (the strongest in continuum emission). The line ratio maps are
quite noisy and only few spaxels are above the S/N threshold of 3 (Fig. 12).
Therefore, the same applies to the metallicity map, which generally shows low
values, of 12+log(O/H) $\lesssim$ 7.8, in the few spaxels where it can be
estimated (Fig. 3, right).

##  5 COSMOS-24108

Figure 4: Same as in Fig. 3 for COSMOS24108. The continuum map (left) is the
median in the 1–2 $\mu$m observed range ($\sim$0.14–0.27 $\mu$m rest-frame).
The PRISM/CLEAR spectrum (bottom) is extracted from a circular aperture with
radius of 0.15″ centred at the location of the most luminous component in
continuum, COSMOS24108-a.

COSMOS-24108 ($z$ $\sim$ 6.36; SFRUV $\sim$ 29 $M_{\odot}$/yr, SFRdust
$\lesssim$ 6.2 $M_{\odot}$/yr) shows two, or possibly three, spatial
components in H-band HST rest-frame UV imaging (Pentericci et al., 2016). In
Fig. 4, left panel, we show the 1–2 $\mu$m observed continuum ($\sim$0.14–0.27
$\mu$m rest-frame) from the NIRSpec IFU data. This shows two main spatial
components, and a third, fainter one in the SE part of the system, consistent
with those seen with HST in a similar spectral band. We label these three
spatial components as COSMOS24108-a, b, and c, in order of continuum
brightness. The continuum at redder wavelengths, between 3–3.5 $\mu$m
($\sim$0.39–0.46 $\mu$m rest-frame; Fig. 14), is dominated by the southernmost
of the two main components, and only an extended emission towards the
northernmost source is present at these wavelengths, instead of a more clearly
separate spatial component as at lower wavelengths. Therefore, the
northernmost source has bluer continuum than the southernmost one. In the SW
corner of the FOV, another galaxy at a lower redshift is also present.

Fig. 4, bottom panel, shows the integrated spectrum (aperture radius = 0.15′′)
associated with the brightest component in continuum emission, COSMOS24108-a.
The rest-frame optical lines [O ii], H$\beta$, [O iii], H$\alpha$, and [S ii]
are detected with high S/N ($\gtrsim$4 on the peak). We obtain a stellar mass
of log($M_{*}/M_{\odot}$) $\sim$ 9.3 from SED fitting and a SFR $\sim$ 15
$M_{\odot}$ yr-1 from H$\alpha$ for COSMOS24108-a (Table 3). The integrated
spectra of the other spatial components in the system are shown in Fig. 13.
For these, we find log($M_{*}/M_{\odot}$) $\sim$ 8.4–8.9 and SFR $\sim$ 2–15
$M_{\odot}$ yr-1.

The [O iii] ionised gas line emission (Fig. 4, central panel; H$\beta$, [O
ii], and H$\alpha$ maps are reported in Fig. 14) is much more extended than
the continuum, revealing two additional bright clumps to the SE of the system,
one of the two also tentatively detected in continuum (see contours). We label
these as COSMOS24108-[OIII]-Ea and [OIII]-Eb. The H$\alpha$ map (Fig. 14,
bottom-left) shows a bridge of gas (visible also in [O iii], though weaker
relative the rest of the emission) connecting the main system to the closer of
these two clumps, COSMOS24108-[OIII]-Ea. By comparing the ionised gas line
emission with the 1–2 $\mu$m continuum, we see that the former peaks in
between the two main continuum components, with a preferential extension
towards the northernmost, weaker continuum component (COSMOS24108-b) rather
than to the southernmost, brighter one (COSMOS24108-a).

Fig. 4, right panel, reports the map of gas metallicity (the [O iii]
$\lambda$5008/H$\beta$ and [O iii] $\lambda$5008/[O ii]
$\lambda\lambda$3727,30 line ratio maps are shown in Fig. 14). The metallicity
exhibits differences of up to $\sim$0.5–0.6 dex among the different sources in
the system. The northernmost of the two main continuum sources, COSMOS24108-b,
appears to have lower metallicity (12 + $\log$ (O/H) $\sim$ 8.0–8.1) than the
southernmost one, COSMOS24108-a ($\sim$ 8.2), while the SE minor continuum
component, COSMOS24108-c, has larger metallicity ($\sim$ 8.4). Towards the two
ionised gas clumps to the SE, COSMOS24108-[OIII]-Ea and COSMOS24108-[OIII]-Eb,
the metallicity is the lowest, with values of $\sim$ 7.8–7.9.

Figure 5:  Left: Star formation rate versus stellar mass (left) and gas-phase
metallicity versus stellar mass (right) for the single spatial components in
the systems from this work (red, green, and blue circles). In grey circles we
show for comparison the compilation of $z$ $\sim$ 3–10 galaxies from JADES
(Curti et al. 2023b; including GNz-11, Bunker et al. 2023), CEERS (Nakajima et
al., 2023, as re-computed for consistency in Curti et al. 2023b), and EROs
(Curti et al., 2023a; Laseter et al., 2024). The main sequences of star
formation (SFMS) at z $\sim$ 0 and z $\sim$ 8 from Popesso et al. (2023) are
displayed for reference. Right: best-fit mass-metallicity relations (MZRs)
from Curti et al. (2023b) for these $z$ $\sim$ 3–10 targets, Li et al. (2023)
at $z$ $\sim$ 3 (GLASS), Sanders et al. (2021) at $z$ $\sim$ 3.3 and 2.3
(MOSDEF), and Curti et al. (2020b) at $z$ $\sim$ 0.08 (SDSS). The high-$z$
MZRs predicted by EAGLE ($z$ $\sim$ 8; Schaye et al., 2015), FIRE ($z$ $\sim$
6; Ma et al., 2016), IllustrisTNG ($z$ $\sim$ 8; Torrey et al., 2019),
FirstLight ($z$ $\sim$ 6 and $z$ $\sim$ 8; Langan et al., 2020), and SERRA
($z$ $\sim$ 8; Pallottini et al., 2022) cosmological simulations are also
displayed.

##  6 Integrated mass-metallicity relation

We obtained the properties of each sub-source belonging to each target, by
extracting integrated spectra from circular apertures of radius of 0.15′′ for
A2744-YD4 and COSMOS24108 and of 0.1′′ for BDF-3299 centred on each sub-source
(Table 3). We find emission line ratios of log([O iii]/H$\beta$) $\sim$
0.4–0.9 for A2744-YD4, log([O iii]/H$\beta$) $\sim$ 0.7–0.8 for BDF-3299, and
log([O iii]/H$\beta$) $\sim$ 0.5–0.9 and log([S ii]/H$\alpha$) $\lesssim$–0.8
for COSMOS24108. These are consistent with the range of values found for $z$
$\sim$ 6–8 star-forming galaxies in the JADES survey (Cameron et al., 2023),
of log([O iii]/H$\beta$) $\sim$ 0.5–0.8 and log([S ii]/H$\alpha$) $\lesssim$
–0.8 (when only considering detections and not upper limits), in a range of
stellar masses (log($M_{*}/M_{\odot}$) $\sim$ 6.5–9.0) and SFRs ($\sim$0.1–30)
comprising those of the sources in this work. In some cases, like
COSMOS24108-a with log([O iii]/H$\beta$) $\sim$ 0.8 and log([S ii]/H$\alpha$)
$\sim$ –0.8, the ratios are at the high end of the values reported in (Cameron
et al., 2023).

We inferred the gas-phase metallicity and the SFR from emission-line fitting
and the stellar mass from SED fitting, as described in Sect. 2.2. In Fig. 5 we
display the metallicity versus stellar mass diagram for the spatial components
identified in the images (Table 3). We compare these with the values reported
by Curti et al. (2023b) for the $z$ $\sim$ 3–10 galaxies from the JADES
(including GNz-11; Bunker et al. 2023), CEERS (Nakajima et al., 2023, re-
computed by Curti et al. 2023b for consistency), and EROs (Curti et al.,
2023a; Laseter et al., 2024) samples (grey circles). The best-fit relation for
this $z$ $\sim$ 3–10 compilation obtained by Curti et al. (2023b) is also
shown. We also report for reference the best-fit mass-metallicity relations
(MZRs) for galaxies at $z$ $\sim$ 0.08 from SDSS (Curti et al., 2020b), $z$
$\sim$ 2–3 from MOSDEF (Sanders et al., 2021), and $z$ $\sim$ 3 from GLASS (Li
et al., 2023).

Our targets are generally compatible with the values found for the $z$ $\sim$
3–10 galaxies, especially the sources in BDF-3299 which sit on the MZR defined
by this high-$z$ collection. We note that most of the sources in A2744-YD4 and
COSMOS24108 ($z$ $\sim$ 8 and 6, respectively) are at the high end of the
JADES, CEERS, and EROs points, where the $z$ $\sim$ 2–3 MZRs from MOSDEF and
GLASS lie. This suggests that the A2744-YD4 and COSMOS24108 systems may
comprise more evolved sources as compared to the majority of sources at the
same redshifts and be instead more similar to galaxies at cosmic noon.

We further display the predictions for the high-$z$ MZR from a number of
cosmological simulations, specifically EAGLE ($z$ $\sim$ 8; Schaye et al.,
2015), FIRE ($z$ $\sim$ 6; Ma et al., 2016), IllustrisTNG ($z$ $\sim$ 8;
Torrey et al., 2019), FirstLight ($z$ $\sim$ 6 and $z$ $\sim$ 8; Langan et
al., 2020), and SERRA ($z$ $\sim$ 8; Pallottini et al., 2022)111The MZR from
the SERRA simulation (Pallottini et al., 2022) was originally presented in
Curti et al. (2023a) and was computed using 202 simulated galaxies at $z=7.7$.
In the present work, the MZR from SERRA is recomputed from an extended sample
of 245 objects Pallottini & Ferrara (2023). This causes some differences at
the low-mass end of the reported MZR with respect to that in Curti et al.
(2023a).. We note that the sources analysed in this work, both as a whole and
within the individual A2744-YD4 and COSMOS24108 systems, are better aligned
with simulations predicting steeper MZR slopes than the median of the
observations at $z$ $\sim$ 3–10.

##  7 Gas-phase metallicity gradients

Figure 6: Metallicity gradients for the main sources in A2744-YD4 (top two rows), BDF-3299 (mid-bottom), and COSMOS24108 (bottom) from radial annuli having radial width of 1 spaxel (0.05″), from the wavelength-dependent spatially smoothed cube (details in Sect. 2.2). Metallicity gradients (right) and [O iii] maps (left; also from smoothed cube), where the concentric annuli are reported on the source. The vertical grey dashed line marks the spatial resolution ($\sigma_{\mathrm{res.}}$). The fluxes of emission lines far in wavelength are corrected for extinction based on H$\gamma$/H$\beta$ in A2744-YD4 and BDF-3299, while in COSMOS24108 H$\alpha$/H$\beta$ is used for that and also [S ii] and H$\alpha$ are used to infer the metallicity. A cut of S/N$>3$ on the peak flux of each line is applied. Table 4: Metallicity gradients inferred for the main sources in this work. Source name |  $\nabla_{\mathrm{r}}\,\log(Z)$ [dex kpc-1]
---|---
A2744-YD4 | 0.14${}_{-0.16}^{+0.16}$
A2744-YD1 | –0.05${}_{-0.13}^{+0.12}$
BDF-3299 | –0.03${}_{-0.12}^{+0.12}$
COSMOS24108 | 0.03${}_{-0.04}^{+0.04}$

To more robustly constrain the spatially resolved metallicity, we extracted
the line fluxes from spatially integrated spectra from concentric radial
annuli, as described in Sect. 2.2. This allowed us to increase the S/N of the
emission lines at large radii and determine the dust extinction from the faint
Balmer lines. We thus determined the gas-phase metallicity as a function of
the distance from the centre of each source, defined as explained in the
following.

The resulting metallicity gradients are reported in Fig. 6. In the system
A2744-YD4, we determined the metallicity gradients only for the two brightest
sources YD4 and YD1 (centred around their continuum peak) since the
metallicity measurements for the other sources suffered from large
uncertainties due to low S/N of emission lines. For BDF-3299, since we could
have only extracted two radial points for each spatial component (BDF-3299-a,
b, and c) due to the compactness of this system, we obtained the metallicity
gradient profile for the whole BDF-3299 system, centred on the brightest
source in line emission, BDF-3299-b. For the COSMOS24108 system, given that
the line emission peaks between the two main continuum components,
COSMOS24108-a and COSMOS24108-b, and not on either of them (see Fig. 4), we
opted to extract a single set of concentric annular apertures centred on the
line emission peak, rather than two sets centred on each of the continuum
sources (since we considered this latter approach not very sensible in this
case). We also report the metallicity gradient from two additional smaller
apertures centred on the emission line clumps to the SE of the main system,
COSMOS24108-Ea and COSMOS24108-Eb.

For all the sources, A2744-YD4, YD1, BDF-3299, and COSMOS24108, the
metallicity gradients span a range between $\nabla_{r}\log(Z)$ $\sim$
$-0.05\leavevmode\nobreak\ {\rm dex\leavevmode\nobreak\ kpc^{-1}}$ and
$\nabla_{r}\log(Z)$ $\sim$ $0.15\leavevmode\nobreak\ {\rm
dex\leavevmode\nobreak\ kpc^{-1}}$ (Table 4). We note that, within the
uncertainties, all metallicity gradients in our sample are flat.

##  8 Discussion

We found flat gradients within the uncertainties for the sources in A2744-YD4,
BDF-3299, and COSMOS24108 analysed in this work. As introduced in Sect. 1 and
sketched in Fig. 1, negative (radially decreasing) gradients are expected as a
consequence of the inside-out galaxy formation scenario, while inverted
(radially increasing) gradients may arise as a consequence of pristine gas
accretion towards the central regions. Possible explanations for flat radial
gradients are radial gas mixing processes, occurring as a consequence of
mergers, SN-driven winds and/or large-scale gas circulation (“galactic
fountains”, i.e. metal-loaded outflows expelled from the galaxy and re-
accreted in the external regions), which redistribute the gas in the galaxy
and wash out any pre-existing radial gradient of metallicity.

All the targets in this work show multiple sources within a few kpc and
disturbed morphologies. This indicates that they are experiencing, or have
experienced, interactions and galaxy merging processes, which are expected to
be more frequent in the early Universe according to cosmological simulations
(see e.g. Kohandel et al., 2020; Pallottini & Ferrara, 2023). Specifically,
galaxies at $z$ $\sim$ 6 with stellar masses $\gtrsim$ 5$\times$108
$M_{\odot}$ living in dense environments, similar to those studied in this
work (see Table 3), are expected to have already experienced multiple merger
events (Gelli et al., 2020). On the other hand, while in the low-spectral
resolution ($R$ $\sim$ 100) NIRSpec PRISM data analysed in this work the
emission lines are spectrally unresolved, [O iii] and H$\beta$ do not show
evidence for asymmetric wings or any complex line profile indicative of
outflows in the higher resolution ($R$ $\sim$ 2700) NIRSpec IFU or MSA data
from Hashimoto et al. (2023) and Morishita et al. (2023), respectively, for
the case of YD4, YD1, and s1 in the A2744-YD4 system (YD4 having both IFU and
MSA $R$ $\sim$ 2700 data). Therefore, this indicates that strong outflows are
not occurring at present in this system. However, in principle we cannot
exclude that the redistribution of metals across the galaxy may also be the
result of past intense SN winds mixing the ISM or re-accretion in the outer
regions of metal-enriched gas expelled by past galactic outflows, whose trace
is absent in present-day spectra. All in all, mergers and/or possibly either
past galactic outflows or SN wind mixing seem the most likely mechanism
driving the observed flat metallicity gradients.

We stress that, even if we did our best in extracting the metallicity gradient
for each separate spatial component, complications arise because in some cases
these sub-sources are very close to each other and may contaminate each
other’s metallicity (as appears to be the case of YD1 and YD1-E; Fig. 2) or
there is confusion as to which source the detected gas belongs (as in the case
of COSMOS24108-a and COSMOS24108-b, where gas is located between the two
sources). Moreover, the sources are not settled in radially symmetric
metallicity distributions. All these aspects can contribute to yielding a flat
gradient. Nevertheless, the spatial vicinity of the sources and the
displacement of gas and stars are the result of the ongoing merging processes.
Therefore, ultimately, mergers appear to be the most likely cause for the
flatness of the metallicity gradients.

###  8.1 Metallicity gradients across cosmic time

Figure 7: Gas-phase metallicity gradient (in dex kpc-1) as a function of
redshift (bottom $x$ axis) and cosmic time (upper $x$ axis). The targets
analysed in this work are marked with red (A2744-YD4 and YD1), green
(BDF-3299), and blue (COSMOS24108). The other symbols are a compilation of
gas-phase metallicity gradients from literature. Specifically, at $z$ $\sim$
0, we report the measurements from Rupke et al. (2010c, full and empty squares
for isolated and merging galaxies, respectively) and Belfiore et al. (2017,
median). The evolution of metallicity gradients with cosmic time for the Milky
Way, M33, and M81 from H II regions ($z$ $\sim$ 0) and PNe (up to a few Gyr
ago), from Magrini et al. (2009, 2010), Stanghellini & Haywood (2010, 2018),
and Stanghellini et al. (2014) is shown. We report the gradients estimated at
various $z$ from Yuan et al. (2011), Swinbank et al. (2012), Queyrel et al.
(2012, median), Stott et al. (2014, median), Troncoso et al. (2014), Jones et
al. (2013, 2015), Leethochawalit et al. (2016), Wuyts et al. (2016), Carton et
al. (2018), Wang et al. (2017, 2019, 2020, 2022), Förster Schreiber et al.
(2018), Curti et al. (2020a), Gillman et al. (2021), Simons et al. (2021),
Vallini et al. (2024), and Arribas et al. (2023). Hatched regions are reported
in case of large samples. The predictions for the evolution of metallicity
gradients from cosmological simulations are also shown, specifically for MUGS
(enhanced feedback) and MaGICC (normal feedback) from Gibson et al. (2013),
FIRE (Ma et al., 2017), TNG50 (Hemler et al., 2021), and EAGLE (Tissera et
al., 2022).

In Fig. 7, we report the metallicity gradients from a compilation of
observational studies, from $z$ = 0 up to $z$ $\sim$ 8 (this study), namely
the highest redshift probed so far. The flat gradients within the
uncertainties that we find at $z$ $\sim$ 6–8 are compatible with the other
currently available measured gradients at these redshifts, namely Vallini et
al. (2024) from FIR lines with ALMA and Arribas et al. (2023) from rest-frame
optical lines with JWST/NIRSpec. Specifically, the results from Vallini et al.
(2024) showcase the potential for synergies between spatially resolved
observations with JWST and ALMA, though the spatial resolutions may be quite
different ($\sigma$ $\sim$ 0.2–0.4 kpc in our case versus $\sim$ 1.5 kpc in
the case of Vallini et al. 2024).

We also report the predictions for the evolution of metallicity gradients with
redshift from cosmological simulations. Some of them predict steeper negative
gradients with increasing redshift (MUGS, with normal feedback, Gibson et al.
2013; TNG50, Hemler et al. 2021), while others predict median flatter
gradients at higher $z$ though with increasing scatter (MAGICC, with enhanced
feedback, Gibson et al. 2013; FIRE, Ma et al. 2017; EAGLE, Tissera et al.
2022).

Unfortunately, predictions for the cosmic evolution of metallicity gradients
from simulations are limited to $z$ $\lesssim$ 3 so far. Nevertheless, if we
assume that the extrapolation at higher $z$ follows the same trend which is in
place up to $z$ $\sim$ 3, we can compare the predictions with our measurements
at $z$ $\sim$ 6–8. Specifically, a rough linear extrapolation of the trends
from TNG50 (Hemler et al., 2021) and MUGS (Gibson et al., 2013) would give a
metallicity gradient at $z$ $\sim$ 8 of $\nabla_{\mathrm{r}}\,\log(Z)$ $\sim$
–0.2 dex kpc-1 and $\lesssim$–0.5 dex kpc-1, respectively, which are both
incompatible with our estimates considering the uncertainties (Fig. 7).

Our results, showing flat gradients within the uncertainties, are more
compatible with expectations from those cosmological simulations predicting
flatter gradients at high $z$ (MaGICC, Gibson et al. 2013; FIRE, Ma et al.
2017; EAGLE, Tissera et al. 2022). This conclusion is based on the
extrapolation of current model predictions, therefore we stress the need for
updated metallicity gradient predictions from simulations at $z$ $\gtrsim$ 3
in order to match the redshifts reached by current observations.

We point out that the merger fraction increases with redshift (e.g.
theoretically Fakhouri et al. 2010; Rodriguez-Gomez et al. 2015; O’Leary et
al. 2021; observationally Duncan et al. 2019), therefore galaxies at high
redshift often have disturbed morphologies due to interactions and do not have
radially symmetric gas and chemical distributions (see e.g. Figs. 2, 3, and 4;
also e.g. observations in de Graaff et al. 2023; Arribas et al. 2023;
Rodríguez Del Pino et al. 2023 and simulations in Pallottini et al. 2017,
2019). Specifically, the metallicity maps in our work show in some cases non-
azimuthally symmetric structures which are averaged out when producing radial
gradients (as it is e.g. the case of YD1, having the highest values to the W
and the lowest to the E of its spatial emission peak). Therefore, radial
gradients may not be the optimal means to describe the distribution of metals
in high-$z$ galaxies and other alternative methods to quantify it (such as the
metallicity scatter in the galaxy, as suggested by Maiolino & Mannucci 2019)
should be sought and included in predictions from simulations.

In summary, the flat radial gradients and the asymmetric spatial distribution
of metallicity observed in some cases seem to support the merging scenario,
though ISM mixing due to past SN winds or re-accretion in the outer regions of
metal-enriched gas from past galactic outflows (galactic fountains) are also a
viable mechanism to explain the flatness of the gradients. Spatially
asymmetric low-metallicity accretion from the CGM and IGM could also
contribute to produce non-azimuthally symmetric metallicities as suggested by
Arribas et al. (2023).

##  9 Conclusions

In this work, we have presented new JWST NIRSpec IFU observations at low
spectral resolution ($R$ $\sim$ 100; PRISM/CLEAR) of three high-$z$ systems,
namely A2744-YD4 in the proto-cluster A2744-z7p9OD ($z$ $\sim$ 7.88), BDF-3299
($z$ $\sim$ 7.11), and COSMOS24108 ($z$ $\sim$ 6.36). At these redshifts, the
NIRSpec PRISM spectra (spanning $\sim$0.6–5.2 $\mu$m) cover the rest-frame UV
and optical spectral ranges and include the main optical emission lines from
warm ($T$ $\sim$ $10^{4}$ K) ionised gas, [O ii], H$\gamma$, H$\beta$, [O
iii], and, for COSMOS24108, also H$\alpha$ and [S ii]. We modelled the ionised
gas lines in the spectra with Gaussian functions and mapped their emission.
The main goal of this work is to study the spatially resolved gas-phase
metallicity and investigate the shape of metallicity gradients at high
redshift ($z$ $\gtrsim$ 6).

The targets have very disturbed morphologies in both the stellar continuum and
ionised gas line emission. We identify several spatial components concentrated
within a few kiloparsecs in all the three systems. We obtained the main
integrated properties of each spatial component in each target, specifically
stellar mass from SED fitting, SFR from H$\alpha$ or H$\beta$, and gas-phase
metallicity by means of flux ratios of rest-frame optical emission lines which
we modelled with Gaussian functions. We found log($M_{*}/M_{\odot}$) $\sim$
7.6–9.3, SFRs $\sim$ 1–15 $M_{\odot}$ yr-1, and gas-phase metallicities
12+log(O/H) $\sim$ 7.7–8.3 by extracting integrated spectra from circular
apertures centred on each spatial component.

We compared the stellar masses and gas-phase metallicities of the targets with
the mass-metallicity relations (MZRs) inferred observationally at different
redshifts and with those predicted by cosmological simulations. In general,
the sources in the systems studied in this work are consistent with the
distribution of the highest-$z$ galaxies to date ($z$ $\sim$ 3–10) from JADES,
CEERS, and EROs samples in the mass-metallicity plane. Nevertheless, most of
the sources in A2744-YD4 and COSMOS24108 lie at the upper end of this $z$
$\sim$ 3–10 distribution, in terms of metallicity at a given stellar mass, and
are closer to the best-fit MZRs measured at $z$ $\sim$ 2–3 (e.g. MOSDEF and
GLASS surveys) than to the MZR at $z$ $\sim$ 3–10. Relative to the MZRs from
cosmological simulations, the sources in the three systems studied in this
work are in good agreement with the slopes predicted at $z$ $\sim$ 6–8 by most
simulations.

We inferred the gas-phase metallicity radial gradients by extracting
integrated spectra from concentric radial annuli centred on the main sources
of each system for which the S/N allowed for a robust metallicity estimate.
The gas-phase metallicity gradients are flat within the uncertainties. Flat
gradients can be associated with processes which mix the gas and the metals on
galaxy scales and therefore wash out any gradient which may have formed in the
galaxy, such as mergers, SN wind mixing, and re-accretion of metal-loaded
galactic outflows in the outer regions.

These are among the very few measurements of spatially resolved gas
metallicity at these high redshifts ($z$ $\sim$ 6–8). In particular, YD4 (at
$z$ $\sim$ 8) constitutes the highest-$z$ source in which a metallicity
gradient has been probed so far.

Cosmological simulations make very different predictions regarding the cosmic
evolution of metallicity gradients, though they are limited to $z$ $\lesssim$
2–3. Some of them predict steeper negative gradients with increasing redshift
(e.g. MUGS and TNG50), while others predict median flatter gradients at higher
$z$ (e.g. MaGICC, FIRE, and EAGLE). The results of this work, indicating
mainly flat gradients, are more compatible with the predictions from this
latter group of simulations, when considering their extrapolation at $z$
$\sim$ 6–8.

All in all, the results of this work in terms of the MZR and gas-phase
metallicity gradients at $z$ $\sim$ 6–8 provide important constraints to guide
future cosmological simulations. In particular, they urge for specific
predictions on the cosmic evolution of metallicity gradients and metallicity
maps out to redshift of 8, given that such predictions are currently limited
to $z$ $\lesssim$ 3.

Moreover, galaxies at high $z$ tend to have more irregular morphologies and
may not have azimuthally symmetric chemical distributions, as a result of
frequent mergers and asymmetric low-mass gas infall. Radial gradients, which
average out any azimuthal information, may then not be the optimal means to
quantify the distribution of metals in high-$z$ galaxies. Therefore, other
alternative quantitative tracers of the spatially resolved metallicity should
be considered as part of predictions from simulations, to be compared with
observations.

###### Acknowledgements.

We acknowledge support from European Union’s HE ERC Starting Grant No.
101040227 - WINGS (G.V., S.C.), INAF Minigrant “Reionization and fundamental
cosmology with high-redshift galaxies” (M.C.), “FirstGalaxies” Advanced Grant
from the European Research Council (ERC) under the European Union’s Horizon
2020 research and innovation program, Grant agreement No. 789056 (A.J.B.), ERC
Advanced Grant INTERSTELLAR H2020/740120 (A.F.), grant PID2021-127718NB-I00
funded by the Spanish Ministry of Science and Innovation/State Agency of
Research, MICIN/AEI/ 10.13039/501100011033 (S.A.). This work is based on
observations made with the NASA/ESA/CSA James Webb Space Telescope. The data
were obtained from the Mikulski Archive for Space Telescopes at the Space
Telescope Science Institute, which is operated by the Association of
Universities for Research in Astronomy, Inc., under NASA contract NAS 5-03127
for JWST. These observations are associated with program #1893. This research
has made use of NASA’s Astrophysics Data System Bibliographic Services. This
research has made use of the Python packages NumPy (Harris et al., 2020),
SciPy (Virtanen et al., 2020), IPython (Perez & Granger, 2007), Matplotlib
(Hunter, 2007), emcee (Foreman-Mackey et al., 2013), and Astropy
(http://www.astropy.org), a community-developed core Python package for
Astronomy (Astropy Collaboration et al., 2013, 2018). This research has made
use of the FITS files visualisation tool QFitsView
(https://www.mpe.mpg.de/~ott/QFitsView/).

## References

  * Arribas et al. (2023) Arribas, S., Perna, M., Rodríguez Del Pino, B., et al. 2023, arXiv e-prints, arXiv:2312.00899
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Atek et al. (2015) Atek, H., Richard, J., Kneib, J.-P., et al. 2015, ApJ, 800, 18
  * Belfiore et al. (2017) Belfiore, F., Maiolino, R., Tremonti, C., et al. 2017, MNRAS, 469, 151
  * Berg et al. (2015) Berg, D. A., Skillman, E. D., Croxall, K. V., et al. 2015, ApJ, 806, 16
  * Bergamini et al. (2023) Bergamini, P., Acebron, A., Grillo, C., et al. 2023, ApJ, 952, 84
  * Bresolin (2011) Bresolin, F. 2011, ApJ, 730, 129
  * Brook et al. (2012) Brook, C. B., Stinson, G., Gibson, B. K., Wadsley, J., & Quinn, T. 2012, MNRAS, 424, 1275
  * Bruzual & Charlot (2003) Bruzual, G. & Charlot, S. 2003, MNRAS, 344, 1000
  * Bunker et al. (2023) Bunker, A. J., Saxena, A., Cameron, A. J., et al. 2023, A&A, 677, A88
  * Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682
  * Cameron et al. (2023) Cameron, A. J., Saxena, A., Bunker, A. J., et al. 2023, A&A, 677, A115
  * Carnall et al. (2018) Carnall, A. C., McLure, R. J., Dunlop, J. S., & Davé, R. 2018, MNRAS, 480, 4379
  * Carniani et al. (2021) Carniani, S., Arribas, S., Bunker, A., et al. 2021, Galaxy assembly at z ¿ 6: unraveling the origin of the spatial offset between the UV and FIR emission, JWST Proposal. Cycle 1, ID. #1893
  * Carniani et al. (2020) Carniani, S., Ferrara, A., Maiolino, R., et al. 2020, MNRAS, 499, 5136
  * Carniani et al. (2017) Carniani, S., Maiolino, R., Pallottini, A., et al. 2017, A&A, 605, A42
  * Carton et al. (2018) Carton, D., Brinchmann, J., Contini, T., et al. 2018, MNRAS, 478, 4293
  * Castellano et al. (2016) Castellano, M., Dayal, P., Pentericci, L., et al. 2016, ApJ, 818, L3
  * Castellano et al. (2010) Castellano, M., Fontana, A., Paris, D., et al. 2010, A&A, 524, A28
  * Ceverino et al. (2016) Ceverino, D., Sánchez Almeida, J., Muñoz Tuñón, C., et al. 2016, MNRAS, 457, 2605
  * Cresci et al. (2010) Cresci, G., Mannucci, F., Maiolino, R., et al. 2010, Nature, 467, 811
  * Curti et al. (2017) Curti, M., Cresci, G., Mannucci, F., et al. 2017, MNRAS, 465, 1384
  * Curti et al. (2023a) Curti, M., D’Eugenio, F., Carniani, S., et al. 2023a, MNRAS, 518, 425
  * Curti et al. (2020a) Curti, M., Maiolino, R., Cirasuolo, M., et al. 2020a, MNRAS, 492, 821
  * Curti et al. (2023b) Curti, M., Maiolino, R., Curtis-Lake, E., et al. 2023b, arXiv e-prints, arXiv:2304.08516
  * Curti et al. (2020b) Curti, M., Mannucci, F., Cresci, G., & Maiolino, R. 2020b, MNRAS, 491, 944
  * Davé et al. (2011) Davé, R., Finlator, K., & Oppenheimer, B. D. 2011, MNRAS, 416, 1354
  * de Graaff et al. (2023) de Graaff, A., Rix, H.-W., Carniani, S., et al. 2023, arXiv e-prints, arXiv:2308.09742
  * D’Eugenio et al. (2023) D’Eugenio, F., Perez-Gonzalez, P., Maiolino, R., et al. 2023, arXiv e-prints, arXiv:2308.06317
  * Duncan et al. (2019) Duncan, K., Conselice, C. J., Mundy, C., et al. 2019, ApJ, 876, 110
  * Fakhouri et al. (2010) Fakhouri, O., Ma, C.-P., & Boylan-Kolchin, M. 2010, MNRAS, 406, 2267
  * Ferland et al. (2017) Ferland, G. J., Chatzikos, M., Guzmán, F., et al. 2017, Rev. Mexicana Astron. Astrofis., 53, 385
  * Foreman-Mackey et al. (2013) Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306
  * Förster Schreiber et al. (2018) Förster Schreiber, N. M., Renzini, A., Mancini, C., et al. 2018, ApJS, 238, 21
  * Gelli et al. (2020) Gelli, V., Salvadori, S., Pallottini, A., & Ferrara, A. 2020, MNRAS, 498, 4134
  * Gibson et al. (2013) Gibson, B. K., Pilkington, K., Brook, C. B., Stinson, G. S., & Bailin, J. 2013, A&A, 554, A47
  * Gillman et al. (2021) Gillman, S., Tiley, A. L., Swinbank, A. M., et al. 2021, MNRAS, 500, 4229
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357
  * Hashimoto et al. (2023) Hashimoto, T., Álvarez-Márquez, J., Fudamoto, Y., et al. 2023, ApJ, 955, L2
  * Hemler et al. (2021) Hemler, Z. S., Torrey, P., Qi, J., et al. 2021, MNRAS, 506, 3024
  * Ho et al. (2015) Ho, I. T., Kudritzki, R.-P., Kewley, L. J., et al. 2015, MNRAS, 448, 2030
  * Hopkins et al. (2014) Hopkins, P. F., Kereš, D., Oñorbe, J., et al. 2014, MNRAS, 445, 581
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90
  * Ishigaki et al. (2016) Ishigaki, M., Ouchi, M., & Harikane, Y. 2016, ApJ, 822, 5
  * Jones et al. (2013) Jones, T., Ellis, R. S., Richard, J., & Jullo, E. 2013, ApJ, 765, 48
  * Jones et al. (2015) Jones, T., Wang, X., Schmidt, K. B., et al. 2015, AJ, 149, 107
  * Kennicutt & Evans (2012) Kennicutt, R. C. & Evans, N. J. 2012, ARA&A, 50, 531
  * Kewley et al. (2010) Kewley, L. J., Rupke, D., Zahid, H. J., Geller, M. J., & Barton, E. J. 2010, ApJ, 721, L48
  * Kohandel et al. (2020) Kohandel, M., Pallottini, A., Ferrara, A., et al. 2020, MNRAS, 499, 1250
  * Kroupa (2001) Kroupa, P. 2001, MNRAS, 322, 231
  * Langan et al. (2020) Langan, I., Ceverino, D., & Finlator, K. 2020, MNRAS, 494, 1988
  * Laporte et al. (2017) Laporte, N., Ellis, R. S., Boone, F., et al. 2017, ApJ, 837, L21
  * Laporte et al. (2019) Laporte, N., Katz, H., Ellis, R. S., et al. 2019, MNRAS, 487, L81
  * Laporte et al. (2014) Laporte, N., Streblyanska, A., Clement, B., et al. 2014, A&A, 562, L8
  * Laseter et al. (2024) Laseter, I. H., Maseda, M. V., Curti, M., et al. 2024, A&A, 681, A70
  * Leethochawalit et al. (2016) Leethochawalit, N., Jones, T. A., Ellis, R. S., et al. 2016, ApJ, 820, 84
  * Leja et al. (2019) Leja, J., Carnall, A. C., Johnson, B. D., Conroy, C., & Speagle, J. S. 2019, ApJ, 876, 3
  * Li et al. (2023) Li, M., Cai, Z., Bian, F., et al. 2023, ApJ, 955, L18
  * Lotz et al. (2017) Lotz, J. M., Koekemoer, A., Coe, D., et al. 2017, ApJ, 837, 97
  * Luck et al. (2011) Luck, R. E., Andrievsky, S. M., Kovtyukh, V. V., Gieren, W., & Graczyk, D. 2011, AJ, 142, 51
  * Ma et al. (2016) Ma, X., Hopkins, P. F., Faucher-Giguère, C.-A., et al. 2016, MNRAS, 456, 2140
  * Ma et al. (2017) Ma, X., Hopkins, P. F., Feldmann, R., et al. 2017, MNRAS, 466, 4780
  * Magrini et al. (2016) Magrini, L., Coccato, L., Stanghellini, L., Casasola, V., & Galli, D. 2016, A&A, 588, A91
  * Magrini et al. (2017) Magrini, L., Randich, S., Kordopatis, G., et al. 2017, A&A, 603, A2
  * Magrini et al. (2009) Magrini, L., Sestito, P., Randich, S., & Galli, D. 2009, A&A, 494, 95
  * Magrini et al. (2010) Magrini, L., Stanghellini, L., Corbelli, E., Galli, D., & Villaver, E. 2010, A&A, 512, A63
  * Maiolino et al. (2015) Maiolino, R., Carniani, S., Fontana, A., et al. 2015, MNRAS, 452, 54
  * Maiolino & Mannucci (2019) Maiolino, R. & Mannucci, F. 2019, A&A Rev., 27, 3
  * Morishita et al. (2023) Morishita, T., Roberts-Borsani, G., Treu, T., et al. 2023, ApJ, 947, L24
  * Nakajima et al. (2023) Nakajima, K., Ouchi, M., Isobe, Y., et al. 2023, ApJS, 269, 33
  * Nakajima et al. (2022) Nakajima, K., Ouchi, M., Xu, Y., et al. 2022, ApJS, 262, 3
  * Nelson et al. (2019) Nelson, D., Pillepich, A., Springel, V., et al. 2019, MNRAS, 490, 3234
  * O’Leary et al. (2021) O’Leary, J. A., Moster, B. P., Naab, T., & Somerville, R. S. 2021, MNRAS, 501, 3215
  * Osterbrock & Ferland (2006) Osterbrock, D. E. & Ferland, G. J. 2006, Astrophysics of gaseous nebulae and active galactic nuclei
  * Pallottini & Ferrara (2023) Pallottini, A. & Ferrara, A. 2023, A&A, 677, L4
  * Pallottini et al. (2019) Pallottini, A., Ferrara, A., Decataldo, D., et al. 2019, MNRAS, 487, 1689
  * Pallottini et al. (2022) Pallottini, A., Ferrara, A., Gallerani, S., et al. 2022, MNRAS, 513, 5621
  * Pallottini et al. (2014) Pallottini, A., Ferrara, A., Gallerani, S., Salvadori, S., & D’Odorico, V. 2014, MNRAS, 440, 2498
  * Pallottini et al. (2017) Pallottini, A., Ferrara, A., Gallerani, S., et al. 2017, MNRAS, 465, 2540
  * Pentericci et al. (2016) Pentericci, L., Carniani, S., Castellano, M., et al. 2016, ApJ, 829, L11
  * Perez & Granger (2007) Perez, F. & Granger, B. E. 2007, Computing in Science and Engineering, 9, 21
  * Pilkington et al. (2012) Pilkington, K., Few, C. G., Gibson, B. K., et al. 2012, A&A, 540, A56
  * Pillepich et al. (2019) Pillepich, A., Nelson, D., Springel, V., et al. 2019, MNRAS, 490, 3196
  * Planck Collaboration et al. (2020) Planck Collaboration, Aghanim, N., Akrami, Y., et al. 2020, A&A, 641, A6
  * Popesso et al. (2023) Popesso, P., Concas, A., Cresci, G., et al. 2023, MNRAS, 519, 1526
  * Prantzos & Boissier (2000) Prantzos, N. & Boissier, S. 2000, MNRAS, 313, 338
  * Queyrel et al. (2012) Queyrel, J., Contini, T., Kissler-Patig, M., et al. 2012, A&A, 539, A93
  * Rich et al. (2012) Rich, J. A., Torrey, P., Kewley, L. J., Dopita, M. A., & Rupke, D. S. N. 2012, ApJ, 753, 5
  * Rodríguez Del Pino et al. (2023) Rodríguez Del Pino, B., Perna, M., Arribas, S., et al. 2023, arXiv e-prints, arXiv:2309.14431
  * Rodriguez-Gomez et al. (2015) Rodriguez-Gomez, V., Genel, S., Vogelsberger, M., et al. 2015, MNRAS, 449, 49
  * Rupke et al. (2010a) Rupke, D., Kewley, L., & Barnes, J. 2010a, in Astronomical Society of the Pacific Conference Series, Vol. 423, Galaxy Wars: Stellar Populations and Star Formation in Interacting Galaxies, ed. B. Smith, J. Higdon, S. Higdon, & N. Bastian, 355
  * Rupke et al. (2010b) Rupke, D. S. N., Kewley, L. J., & Barnes, J. E. 2010b, ApJ, 710, L156
  * Rupke et al. (2010c) Rupke, D. S. N., Kewley, L. J., & Chien, L. H. 2010c, ApJ, 723, 1255
  * Samland et al. (1997) Samland, M., Hensler, G., & Theis, C. 1997, ApJ, 476, 544
  * Sánchez et al. (2014) Sánchez, S. F., Rosales-Ortega, F. F., Iglesias-Páramo, J., et al. 2014, A&A, 563, A49
  * Sanders et al. (2021) Sanders, R. L., Shapley, A. E., Jones, T., et al. 2021, ApJ, 914, 19
  * Sanders et al. (2024) Sanders, R. L., Shapley, A. E., Topping, M. W., Reddy, N. A., & Brammer, G. B. 2024, ApJ, 962, 24
  * Schaye et al. (2015) Schaye, J., Crain, R. A., Bower, R. G., et al. 2015, MNRAS, 446, 521
  * Shivaei et al. (2020) Shivaei, I., Reddy, N., Rieke, G., et al. 2020, ApJ, 899, 117
  * Simons et al. (2021) Simons, R. C., Papovich, C., Momcheva, I., et al. 2021, ApJ, 923, 203
  * Somerville & Davé (2015) Somerville, R. S. & Davé, R. 2015, ARA&A, 53, 51
  * Stanghellini & Haywood (2010) Stanghellini, L. & Haywood, M. 2010, ApJ, 714, 1096
  * Stanghellini & Haywood (2018) Stanghellini, L. & Haywood, M. 2018, ApJ, 862, 45
  * Stanghellini et al. (2014) Stanghellini, L., Magrini, L., Casasola, V., & Villaver, E. 2014, A&A, 567, A88
  * Stinson et al. (2010) Stinson, G. S., Bailin, J., Couchman, H., et al. 2010, MNRAS, 408, 812
  * Storey & Zeippen (2000) Storey, P. J. & Zeippen, C. J. 2000, MNRAS, 312, 813
  * Stott et al. (2014) Stott, J. P., Sobral, D., Swinbank, A. M., et al. 2014, MNRAS, 443, 2695
  * Swinbank et al. (2012) Swinbank, A. M., Sobral, D., Smail, I., et al. 2012, MNRAS, 426, 935
  * Tissera et al. (2022) Tissera, P. B., Rosas-Guevara, Y., Sillero, E., et al. 2022, MNRAS, 511, 1667
  * Torres-Flores et al. (2014) Torres-Flores, S., Scarano, S., Mendes de Oliveira, C., et al. 2014, MNRAS, 438, 1894
  * Torrey et al. (2019) Torrey, P., Vogelsberger, M., Marinacci, F., et al. 2019, MNRAS, 484, 5587
  * Troncoso et al. (2014) Troncoso, P., Maiolino, R., Sommariva, V., et al. 2014, A&A, 563, A58
  * Vallini et al. (2024) Vallini, L., Witstok, J., Sommovigo, L., et al. 2024, MNRAS, 527, 10
  * Vanzella et al. (2011) Vanzella, E., Pentericci, L., Fontana, A., et al. 2011, ApJ, 730, L35
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261
  * Wang et al. (2022) Wang, X., Jones, T., Vulcani, B., et al. 2022, ApJ, 938, L16
  * Wang et al. (2020) Wang, X., Jones, T. A., Treu, T., et al. 2020, ApJ, 900, 183
  * Wang et al. (2019) Wang, X., Jones, T. A., Treu, T., et al. 2019, ApJ, 882, 94
  * Wang et al. (2017) Wang, X., Jones, T. A., Treu, T., et al. 2017, ApJ, 837, 89
  * Wuyts et al. (2016) Wuyts, E., Wisnioski, E., Fossati, M., et al. 2016, ApJ, 827, 74
  * Yuan et al. (2011) Yuan, T. T., Kewley, L. J., Swinbank, A. M., Richard, J., & Livermore, R. C. 2011, ApJ, 732, L14
  * Zaritsky et al. (1994) Zaritsky, D., Kennicutt, Robert C., J., & Huchra, J. P. 1994, ApJ, 420, 87
  * Zheng et al. (2014) Zheng, W., Shu, X., Moustakas, J., et al. 2014, ApJ, 795, 93

##  Appendix A Additional figures for A2744-YD4

Figure 8: NIRSpec PRISM/CLEAR spectra extracted from a circular aperture with
radius of 0.15″ centred at the location of the targets YD6, YD6-[O iii]-E,
YD4-[O iii]-W, and ZD1. Figure 9: NIRSpec PRISM/CLEAR spectra extracted from a
circular aperture with radius of 0.15″ centred at the location of the targets
YD1, YD1-E, and s1.

Figure 10: Additional maps for A2744-YD4. Highest-resolution continuum map of
A2744-YD4 from the original unsmoothed data cube in observed spectral range
1.2-2 $\mu$m, where the spatial resolution is the highest (top-left). All the
other maps are obtained from the (wavelength-dependent) spatially smoothed
data cube, as in Fig. 2 (details in Sect. 2.2). First row: H$\beta$ (centre)
and [O ii] (right) integrated fluxes. Second row: [O iii]/H$\beta$ (left) and
[O iii]/[O ii] (right) line flux ratios. Contours mark the 2–3 $\mu$m
(observed-frame) continuum.

##  Appendix B Additional figures for BDF-3299

Figure 11: NIRSpec PRISM/CLEAR spectra extracted from a circular aperture with
radius of 0.1″ centred at the location of the targets BDF-3299-b and
BDF-3299-c.

Figure 12: Additional maps for BDF-3299. As in Fig. 3, the flux maps are
obtained from the original data cube, while the line ratio maps from the
spatially smoothed data cube (details in Sect. 2.2). Top row: maps of H$\beta$
integrated flux (left), [O iii]/H$\beta$ (centre) and [O iii]/[O ii] (right)
line flux ratios.

##  Appendix C Additional figures for COSMOS24108

Figure 13: NIRSpec PRISM/CLEAR spectra extracted from a circular aperture with
radius of 0.15″ centred at the location of the targets COSMOS24108-b,
COSMOS24108-c, COSMOS24108-[OIII]-Ea, and COSMOS24108-[OIII]-Eb.

Figure 14: Additional maps for COSMOS24108. As in Fig. 4, the flux maps are
obtained from the original data cube, while the line ratio maps from the
spatially smoothed data cube (details in Sect. 2.2). Top row: median continuum
map in the observed spectral range 3–3.5 $\mu$m ($\sim$0.39–0.46 $\mu$m rest-
frame; left), with the contours of the continuum in the observed spectral
range 1–2 $\mu$m (same as in Fig. 4, left) superimposed; H$\beta$ (centre-
left), [O ii] (centre-right), and H$\alpha$ (right) flux maps. Middle row: [O
iii]/H$\beta$ (left) and [O iii]/[O ii] (centre-left) line ratio maps.

Figure 15: Metallicity gradients for the two minor sources in COSMOS24108, [O
iii]-Ea and [O iii]-Eb.