Probing interstellar medium conditions at z 5.5-9.5 with ultra-deep JWST/NIRSpec spectroscopy

We present emission-line ratios from a sample of 27 Lyman-break galaxies from z ∼ 5 . 5 − 9 . 5 with − 17 . 0 < M 1500 < − 20 . 4, measured from ultra-deep JWST / NIRSpec multi-object spectroscopy from the JWST Advanced Deep Extragalactic Survey (JADES). We used a combination of 28h deep PRISM / CLEAR and 7h deep G 395 M / F 290 LP observations to measure, or place strong constraints on, ratios of widely studied rest-frame optical emission lines including H α , H β , [O ii ] λλ 3726, 3729, [Ne iii ] λ 3869, [O iii ] λ 4959, [O iii ] λ 5007, [O i ] λ 6300, [N ii ] λ 6583, and [S ii ] λλ 6716, 6731 in individual z > 5 . 5 spectra. We ﬁnd that the emission-line ratios exhibited by these z ∼ 5 . 5 − 9 . 5 galaxies occupy clearly distinct regions of line-ratio space compared to typical z ∼ 0 − 3 galaxies, instead being more consistent with extreme populations of lower-redshift galaxies. This is best illustrated by the [O iii ] / [O ii ] ratio, tracing interstellar medium (ISM) ionisation, in which we observe more than half of our sample to have [O iii ] / [O ii ] > 10. Our high signal-to-noise spectra reveal more than an order of magnitude of scatter in line ratios such as [O ii ] / H β and [O iii ] / [O ii ], indicating signiﬁcant diversity in the ISM conditions within the sample. We ﬁnd no convincing detections of [N ii ] λ 6583 in our sample, either in individual galaxies, or a stack of all G 395 M / F 290 LP spectra. The emission-line ratios observed in our sample are generally consistent with galaxies with extremely high ionisation parameters (log U ∼ − 1 . 5), and a range of metallicities spanning from ∼ 0 . 1 × Z (cid:12) to higher than ∼ 0 . 3 × Z (cid:12) , suggesting we are probing low-metallicity systems undergoing periods of rapid star formation, driving strong radiation ﬁelds. These results highlight the value of deep observations in constraining the properties of individual galaxies, and hence probing diversity within galaxy population.


Introduction
Emission lines can be some of the most prominent features of galaxy spectra.Regions of photoionised gas around young stars (Hii regions) routinely exhibit bright emission from hydrogen and helium recombination, as well as collisionally excited emission from numerous metal ions, especially O + , O 2+ , Ne 2+ , S + , N + (e.g.Peimbert et al. 2017).Ratios of emission line fluxes are sensitive to abundances of metals and their ionisation states, but also the physical conditions of the emitting nebulae, such as temperature and density, which are in turn linked to the nature of the ionising sources powering the emission (e.g.Osterbrock & Ferland 2006).Emission line ratios therefore have great diagnostic power for understanding properties of the interstellar medium (ISM) and stellar populations in galaxies (see Kewley et al. 2019 for a review).
However, observed emission lines in galaxies can arise from other phenomena, including active galacitc nuclei (AGN; Kewley et al. 2006;Feltre et al. 2016) or shock heated gas (Dopita & Sutherland 1996;Allen et al. 2008;Sutherland & Dopita 2017).Meanwhile, stellar populations other than the young O and B E-mail: alex.cameron@physics.ox.ac.uk stars typically associated with Hii regions, such as high-mass Xray binaries (HMXBs; Senchyna et al. 2020) or post AGB stars (Byler et al. 2019), can power emission in photoionised regions of gas.Each of these heating sources yields characteristically different emission spectra.An observed global emission spectrum of a star-forming galaxy will typically be a superposition of emission from many Hii regions with different conditions as well as possible contributions from these additional sources, making the interpretation of emission spectra challenging.
Despite this complexity, galaxies have been observed to follow remarkably tight trends in emission line-ratio space.The so-called BPT and VO87 diagrams, which relate the ratio of [O iii] λ5007/Hβ to ratios of either [N ii] λ6583/Hα, [S ii] λλ6716,6731/Hα, or [O i] λ6300/Hα (Baldwin et al. 1981;Veilleux & Osterbrock 1987), are perhaps the best known examples of such trends 1 .With large surveys like SDSS (York et al. 2000), star-forming galaxies at z ∼ 0 have been observed to form a tight sequence in the [N ii]-BPT diagram, with a general anti-correlation between the two ratios (e.g.Kauffmann et al. 2003;Kewley et al. 2006).This tight correlation is the result of the fact that the two parameters that most strongly dictate the emission spectra of galaxies, metallicity and ionisation parameter2 (U), are themselves generally correlated in galaxies (Dopita & Evans 1986;Carton et al. 2017).Position within this sequence is thus dominated by these two parameters, although variations in other properties including N/O abundance ratio, electron density and radiation hardness also contribute (e.g.Kewley et al. 2013a;Bian et al. 2016;Curti et al. 2022).These correlations are not unique to the BPT and VO87 diagrams, and indeed other combinations of emission lines, such as ([O ii]+[O iii])/Hβ vs.
Existing studies of earlier cosmic epochs show that, compared to the loci inhabited by local galaxies in line-ratio space, z ∼ 2 galaxies are generally offset to higher excitation and higher ionisation (e.g.Strom et al. 2017;Runco et al. 2021).This is perhaps best exhibited by the widely studied offset in the [N ii]-BPT diagram (Kewley et al. 2013b;Steidel et al. 2014Steidel et al. , 2016;;Runco et al. 2022), where z ∼ 2 galaxies are observed with typically higher [N ii]/Hα at a fixed [O iii]/Hβ.This is often explained as a combination of higher ionisation parameters, higher densities, and harder ionising spectra at high redshift (Brinchmann et al. 2008;Kewley et al. 2013a;Steidel et al. 2016;Sanders et al. 2016;Hirschmann et al. 2017).Many authors have suggested that super solar α/Fe abundance ratios at these earlier times play an important role by driving harder stellar ionisation fields (e.g.Strom et al. 2017;Topping et al. 2020b,a).However, the N/O abundance ratio is clearly another important parameter to consider, with some studies suggesting this offset can be explained as the result of higher N/O at fixed [O iii]/Hβ (Masters et al. 2016).
Despite the wealth of rest-frame optical studies highlighting the different ISM conditions in galaxies at z ∼ 2, emission line studies of galaxies at z 6 were, until recently, very limited.Broad band excesses had already suggested an abundance of galaxies at this epoch with very high equivalent width emission from [O iii]+Hβ (e.g.Endsley et al. 2021Endsley et al. , 2022;;Fujimoto et al. 2023).However constraints on their ISM conditions were limited to studies of far-infrared lines, such as [O iii] 88 µm and [C ii] 158 µm from ALMA, which suggested high ionisation parameters and enhanced α/Fe abundance ratios (e.g.Hashimoto et al. 2019;Carniani et al. 2020;Harikane et al. 2020;Fujimoto et al. 2022;Witstok et al. 2022).
The recent commissioning of JWST is transforming our ability to probe ISM physics in the early Universe.Already, emission line measurements from medium-depth observations from the initial months of JWST operations have observed z 6 galaxies to exhibit generally low metallicity, high ionisation, high excitation, and high temperature (Curti et al. 2023;Katz et al. 2023;Sun et al. 2022a;Rhoads et al. 2022;Tacchella et al. 2022;Trump et al. 2022;Langeroodi et al. 2022;Williams et al. 2022;Heintz et al. 2022;Matthee et al. 2022;Fujimoto et al. 2023;Mascia et al. 2023;Sanders et al. 2023;Nakajima et al. 2023;Tang et al. 2023).Furthermore, the observed emission spectra of some of these targets have been found to be consistent with the presence of AGN activity (Brinchmann 2022), while Katz et al. (2023) suggested that the presence of high mass X-ray binaries could help to explain the observed emission line ratios.
Although much insight has already been gained from these early observations, they have largely been based on relatively shallow spectroscopy.Here we leverage data taken from the 'Deep' spectroscopic tier of the JWST Advanced Deep Extragalactic Survey (JADES), the deepest spectroscopic observations yet taken with NIRSpec, to provide a more detailed look at the emission line ratios of galaxies at z ∼ 5.5 − 9.5.These MSA observations reached exposure times of up to 28 hours in the PRISM/CLEAR (R ∼ 100) and up to 7 hours in each of the 3 medium resolution gratings (R ∼ 1000), providing unprecedented new insights into the ISM properties of galaxies within the first Gyr of the Universe's history.The structure of this paper is as follows.In Section 2 we describe the observations, data reduction and measurement of emission line fluxes.In Section 3 we present our emission line ratios on various diagnostic diagrams and provide comparisons with literature samples.Section 4 then presents a discussion of the implications of these observations.We summarise in Section 5. Throughout this paper we adopt the Planck Collaboration et al. ( 2020) cosmology: Ω Λ = 0.6847, Ω m = 0.3153, H 0 = 67.36km s −1 Mpc −1 .All magnitues are quoted in the AB magnitude system (Oke & Gunn 1983).

Observations
The data presented in this paper were obtained via multiobject spectroscopic observations from JWST/NIRSpec using the micro-shutter assembly (MSA; Jakobsen et al. 2022;Ferruit et al. 2022).Observations were carried out as part of JADES in three visits between 21-25 Oct 2022 (Program ID: 1210; PI: N. Luetzgendorf) in the Great Observatories Origins Deep Survey South (GOODS-S) legacy field in a region overlapping with the Hubble Ultra Deep Field (HUDF; Koekemoer et al. 2013).Each visit consisted of 33,613 s integration in the PRISM/CLEAR low-resolution setting and 8,403 s integration in each of G140M/F070LP, G235M/F170LP, G395M/F290LP, and G395H/F290LP filter/grating settings.Across three visits, this totals 28 hours of integration in the PRISM, which provides continuous spectral coverage from 0.6 -5.3 µm at R ∼ 30 − 300, and ∼7 hours in each of the medium resolution gratings, which combine to provide R ∼ 700 − 1300 across the full spectral range of NIRSpec, plus 7 hours in the high-resolution grating which provides R ∼ 2700 from ∼2.8 -5.1 µm.In this paper we present only measurements from the PRISM/CLEAR and G395M/F290LP observations.Results obtained from other configurations will be presented in forthcoming work.
Observations within each visit were performed as a 3-shutter nod.The central pointing of each visit was dithered (by < 1 arcsec) such that common targets were observed in different shutters and different detector real-estate.Thus, each visit had a unique MSA configuration, although target allocation (performed with the eMPT3 ) was optimised for maximising target commonality between all three dither positions.A total of 253 unique targets were observed in the PRISM configuration with the three dithers featuring 145, 155, and 149 targets respectively.67 targets were triply exposed in the PRISM and reached full depth, while 62 targets featured in only two dithers and a further 124 received only one-third total depth.The hierarchical allocation of shutters means that the double-and triple-exposed targets are biased toward the highest priority (and hence highest redshift) targets.Indeed, of the 27 galaxies presented in this work, 18 received triple coverage, 6 received double, and 3 received only single coverage.Details of the target selection scheme are presented in Bunker et al. (2023).
The low dispersion associated with the PRISM mode means that all targets are observed with non-overlapping spectra.However, in the higher resolution modes, individual spectra are dispersed over a much larger extent of detector real-estate.To minimise the possibility of contaminating emission, in our grating observations we isolate our highest priority targets by closing the shutters of low-priority targets on the same row (i.e.targets that could cause overlapping spectra).Thus, for our grating spectra we observe only 198 unique targets (119, 121, and 111 in each dither).However, the high-priority z ∼ 5.5 − 9.5 galaxies presented in this work are almost unaffected by this measure with only two exceptions.One target (JADES-GS+53.11351-27.77284)that was triply exposed in PRISM was only singly exposed in the gratings, while one target (JADES-GS+53.17582-27.77446)was not observed in the medium resolution setting at all.The total integration times obtained on each target are summarised in Table 1.

Data reduction
These observations were processed by adopting algorithms developed by the ESA NIRSpec Science Operations Team (SOT) and the NIRSpec GTO Team.Details of the data-processing workflow will be presented in a forthcoming NIRSpec GTO collaboration paper.As this work uses the observations carried out with only two spectral configurations of the program, we describe the main steps for the PRISM/CLEAR and G395M/F290LP filter/grating settings.Upon retrieving the level-1a data from the MAST archive, we estimated the count rate per pixel by using the unsaturated groups in the ramp and removing jumps due to cosmic rays identified by estimating the slope of the individual ramps.During this first stage, we also performed the master bias and dark subtraction, corrected snowball artifacts, and flagged saturated pixels.We then performed the pixel-bypixel background subtraction by combining the three nod exposures of each pointing.We note that for some targets we exclude one of the 3-shutter nods in the background subtraction stage due to the serendipitous presence of a source contaminating the open background shutter.We then created two-dimensional (2D) cutouts of each 3-shutter slit and performed the flat-field, spectrograph optics, and disperser corrections.We then ran the ab- solute calibration and corrected the 2D spectra for path-losses.
Path-loss corrections are dependent on the morphology of the target and its relative position within the shutter.In this paper, we apply the path-loss corrections assuming a point-like source at the location of the target centroid.The high-redshift targets presented in this paper are generally smaller, or comparable in size, to the angular resolution of the telescope at the observed wavelength of emission lines used in this study (λ 3 µm).We estimate that for the sample presented here, even the ratio with the longest wavelength baseline ([O iii] λ5007 / [O ii] λλ3726, 3729) would be affected only at the 1 % (0.01 dex) level by adopting an extended morphology based on fitting to broad-band imaging.Since it is not a priori known how closely the spatial distribution of emission lines trace the starlight observed in broadband imaging, we decide to adopt the point-like assumption, but acknowledge that this effect introduces some systematic uncertainty to our measurements.
We rectified and interpolated the 2D continuum map onto a regular grid for the G395M/F290LP observation, and an irregular grid for the PRISM/CLEAR to avoid an oversampling of the line spread function at short wavelengths.Finally, the 1D spectra were extracted from the 2D map adopting a box-car aperture as large as the shutter size, centred on the relative position of the target in the shutter.For each target, we combined all 1D spectra and removed any bad pixels with sigma-clipping.Cutouts of five example representative 1D spectra are shown in Figure 1.
As a verification of our flux calibration, we compared the NIRSpec spectrophotometry to broadband photometry from NIRCam and found the differences to be on average less than 10 % even before considering the effects of slit-losses due to the MSA.Details of this comparison are presented in Appendix B.

Emission line ratio sample
Although 253 unique targets were observed in these observations, this study focuses only on a rest-UV selected sample of confirmed z > 5.5 galaxies with robustly detected emission lines.The criteria used to select targets for this program are described in detail in Bunker et al. (2023).We describe here only the aspects most relevant to this work.
The highest priorities in target assignment were reserved for galaxies at the highest redshifts, with a preference for the UV-bright galaxies.Candidate z > 5.5 galaxies were primarily drawn from a compilation of existing literature catalogs based on Hubble Space Telescope (HST) broad-band imaging (e.g Bunker et al. 2004;Bouwens et al. 2015Bouwens et al. , 2021;;Finkelstein et al. 2015;Harikane et al. 2016), with a number of the highest redshift galaxies presented here having been identified with HST (e.g.Bouwens et al. 2011Bouwens et al. , 2016;;Ellis et al. 2013;McLure et al. 2013;Oesch et al. 2013).In fact, only five galaxies presented in this work were newly identified from JWST/NIRCam imaging (see Bunker et al. 2023 for details).In total, 47 galaxies were assigned MSA shutters as z > 5.7 candidates out of 299 possible  2. Top: Redshift histogram of galaxies observed in this study.We divide the sample into two sub-samples 'z ∼6' and 'z ∼ 8', based simply on cutting at z = 7 since this is the redshift at which Hα is redshifted beyond the spectral coverage of NIRSpec.Thus, for JADES 'z ∼ 8' we only have coverage of lines from [O ii] λλ3726, 3729 to [O iii] λ5007, while for JADES 'z ∼ 6' we can in principle observe all the rest-frame optical lines out to [S ii] λλ6716, 6731.The orange hatched histogram indicates galaxies for which we can constrain the [S ii] λλ6716, 6731 flux.Green hatching is the same but for [N ii] λ6583.Bottom: Redshift vs. M 1500 for the full z ∼ 5.5 − 9.5 sample.We derived M 1500 directly from the PRISM spectra as described in Section 2.4.1.
such candidates within the MSA footprint.Of these, 28 passed redshift-dependent rest-UV cuts to place them in our highest priority classes PC = 1 and PC = 4 (PC = 2 and 3 contained only 5 targets between them, none of which featured here; see Bunker et al. 2023 for class descriptions).For this analysis, we consider only galaxies with robust spectroscopically confirmed redshifts in which Hβ is detected with S /N > 5. Of the 28 rest-UV-selected targets in PC = 1 and PC = 4, 18 met these criteria.We also include 6 (1) galaxies selected as z > 5.7 candidates that had failed the rest-UV selection cut, but were included in PC = 6.1 (PC = 6.2) as a result of meeting the magnitude criterion of m F160W < 29 (m F160W > 29) and turned out to match our Hβ S /N criterion for this paper.Finally, three galaxies selected in PC = 7.5 (4.5 < z < 5.7, m F160W < 29) had z spec > 5.5 and S /N Hβ > 5, and we included these in our sample as well.One galaxy exhibited a broad-component under Hα in R ∼ 1000 spectroscopy, suggestive of nuclear activity, and was discarded so as to avoid contamination from emission from AGN.Thus, our final sample included 27 star-forming galaxies, drawn from a largely rest-UV-selected sample.
This work focuses predominantly on the low-resolution PRISM/CLEAR spectra ('PRISM' hereafter), since these are significantly deeper than our medium-resolution grating spectra.The spectral resolution of the PRISM varies significantly with wavelength from R ∼ 30 around 1.2 µm up to R ∼ 300 at the long wavelength cutoff around 5 µm (e.g. Figure 5 in Jakobsen et al. 2022).By limiting ourselves to rest-frame optical emission lines at z 5.5, the measurements presented in this paper are made at observed wavelengths of λ obs 2.4 µm where the resolution is R 100, and improving steeply as a function of wavelength.Thus, in all cases from the PRISM, the [O iii] λλ4959, 5007 doublet is completely resolved from Hβ.The doublet itself is always at least partially resolved, and is completely resolved in our highest redshift targets (Figure 1).
[O ii] λλ3726, 3729 and [S ii] λλ6716, 6731 are never resolved into two components, but the latter can easily be fully resolved from the Hα+[N ii] complex over this redshift range.However, Hα and [N ii] λ6583 require resolutions of R ∼ 500 to be fully resolved and thus are significantly blended across the full redshift range.Constraints on this ratio are instead derived from our higher resolution (but shallower) G395M/F290LP grating spectra.
Figure 2 shows the redshift histogram of galaxies in this sample.We divide galaxies into two subsamples, based on the observability of Hα.The purple 'z ∼ 6' sample spans over 5.507 ≤ z ≤ 6.931 (z median = 5.943) and includes 21 galaxies for which Hα is still within the wavelength coverage of NIRSpec.There are an additional six galaxies in our 'z ∼ 8' sample, spanning over 7.05 ≤ z ≤ 9.45 (z median = 7.629), for which we study only emission lines from Since the G395M/F290LP spectra are dispersed across a larger area of detector real estate, they do not provide continuous spectral coverage over the nominal spectral range.For 2/21 galaxies in the z ∼ 6 sample, we do not have [N ii]/Hα constraints because these lines fall into the detector gap of NIRSpec.We additionally have one galaxy (JADES-GS+53.17582-27.77446)which was only observed in the PRISM mode (see Section 2.1).For one galaxy in the z ∼ 6 sample, the [S ii] doublet falls beyond the observed spectral range.

Calculation of M 1500
The continuum is generally well detected in the rest-frame ultraviolet (UV) across our sample.In all except our faintest target, the signal-to-noise on the continuum reaches S /N > 3 for individual spectral pixels long-ward of the Lyman break, and in six spectra we reach in excess of S /N 10 per pixel in the rest-frame UV.For our faintest target, JADES-GS+53.16746-27.77201,we do observe faint continuum in the rest-UV, but the signal-to-noise per pixel is only S /N ≈ 1.8.We derive rest-frame UV absolute magnitudes (M 1500 ) for our sample directly from the spectra.We shift the observed spectra to the rest-frame using the spectroscopically confirmed redshift and integrate across a synthetic 200 Å wide top-hat filter centred on 1500 Å.This integrated flux is then converted into an absolute magnitude.Derived M 1500 values are reported in Table 1.Aside from our faintest target (JADES-GS+53.16746-27.77201with M 1500 = −17.01± 0.22), all other targets have magnitudes in the range −17.5 ≤ M 1500 ≤ −20.4.

PRISM/CLEAR emission line measurements
Emission line measurements and continuum modelling are performed simultaneously using the penalised pixel fitting algo-rithm, ppxf (Cappellari 2017(Cappellari , 2022)).ppxf models the continuum as a linear superposition of simple stellar-population (SSP) spectra, using non-negative weights and matching the spectral resolution of the observed spectrum.As input, we used the highresolution (R=10,000) SSP library combining MIST isochrones (Choi et al. 2016) and the C3K theoretical atmospheres (Conroy et al. 2019).The flux blue-ward of the Lyman break was manually set to 0. These templates are complemented by a 5 thdegree multiplicative Legendre polynomial, to take into account systematic differences between the SSPs and the data (e.g., dust, mismatch between the SSP models and high-redshift stellar populations, and residual flux calibration problems).The emission lines are modelled as pixel-integrated Gaussians, again matching the observed spectral resolution.To reduce the number of degrees of freedom, we divide the emission lines into four kinematic groups, constrained to have the same redshift and intrinsic broadening.These groups are: UV lines (blueward of 3000 Å in the rest-frame), the Balmer series of Hydrogen, non-Hydrogen optical lines (blueward of 9000 Å in the rest-frame), and near-infrared (NIR) lines.The stellar continuum has the same kinematics as the Balmer lines.Furthermore, we tie together doublets that have fixed ratios, and constrain variableratio doublets to their physical ranges.In particular, this study focuses on the following lines of interest: ] λλ6716, 6731.We note that, at the resolution of the PRISM/CLEAR mode, the two components of the [O ii] λλ3726, 3729 doublet are significantly blended and so are fit as a single Gaussian component.We visually inspect these fits to the 1D spectra to ensure they are reliable.For the faintest lines, we also individually inspect the 2D spectrum for each target and confirm that the emission line can be visually identified along the same trace to confirm the measured flux is not affected by uncorrected artifacts.Where a line is not detected, we derive an upper limit directly from the noise spectrum.The 1-σ upper limit is obtained by integrating the variance spectrum over an interval spanning 1.4×σ line either side of the expected centroid of the line (where σ line is the Gaussian line width), and multiplying this by the spectral extent of a pixel.
We note that, at the spectral resolution of our PRISM observations, [Ne iii] λ3869 is partially blended with the He i λ3889 line.From grating spectra, we expect the [Ne iii] line to be brighter than He i, typically by at least a factor of two, suggesting that performing a two-component fit can recover the [Ne iii] flux (e.g.Appendix C3 in Cameron et al. 2021).The fitting procedure outlined above simultaneously fits for both lines in this complex, and we confirm from visual inspection that these fits model the complex well in cases with a high signal-to-noise ratio.In this paper, the only ratio for which we consider the [Ne iii] line is the [Ne iii] λ3869 / [O ii] λλ3726, 3729 ratio.This ratio has a very short wavelength baseline and hence the uncertainties introduced by partial blending are somewhat offset by having almost no dependence on wavelength dependent corrections or calibrations.Nonetheless, we note that that the [Ne iii]-based measurements presented in this work should be treated with caution, and that any inference requiring high precision [Ne iii]-based measurements should be performed with medium-to high-resolution spectra.

Constraints on [N ii]/Hα from G395M/F290LP
As discussed in Section 2.3, the only line ratio for which we consider the data from our higher spectral resolution, but much shallower, G395M/F290LP spectra is [N ii] λ6583 / Hα.These lines are sufficiently close in wavelength that they remain at least partially blended in PRISM/CLEAR spectra across the entire redshift range considered here.However, in the R ∼ 1000 G395M/F290LP spectra, this complex is always completely resolved.
We performed a visual inspection simultaneously for the PRISM and G395M/F290LP 1D and 2D spectra of each target to identify the possible presence of [N ii].In fact, we find no convincing evidence for the detection of [N ii] λ6583 in any of our 18 galaxies, either from the resolved R ∼ 1000 spectra, or by way of the appearance of a 'red wing' in the Hα profile in the PRISM spectra.We derive upper-limits on the [N ii]/Hα ratio in the same way as described above, fitting the Hα line with pPXF, and integrating the variance spectrum across the expected spectral extent of [N ii] λ6583, based on the redshift and line-width measured for Hα.

Emission line ratio measurements
In this paper we focus on the following key diagnostic line ratios: (1) Several of these line ratios, such as N2, S 2 and R3, are calculated across short wavelength baselines and are thus largely insensitive to dust attenuation or wavelength-dependent calibration issues.However, to accommodate the analysis of the longer baseline ratios, we correct our emission line fluxes and upper limits for dust extinction according to the SMC dust law from Gordon et al. (2003).We calculate the E(B − V) from the decrements of Balmer lines.For the z ∼ 6 sample where Hα is available, we use the Hα/Hβ decrement, assuming an intrinsic ratio of Hα/Hβ = 2.86 (Osterbrock & Ferland 2006).For the z ∼ 8 sample, where possible, we use Hγ/Hβ assuming an intrinsic ratio of Hγ/Hβ = 0.468.Above z > 7, Hγ is sufficiently deblended from [O iii] λ4363 so as to avoid a biased measurement.In two cases, the S/N on Hγ is too low.For these galaxies, we adopt the attenuation obtained from SED fitting with BEA-GLE (Chevallard et al. in prep.), but find that only very small corrections (E(B − V) < 0.02) are required.Overall, the mean E(B − V) from this sample is 0.1.Our final dust-corrected emission line ratios are presented in Table 2.

Stacked spectra
Although this paper focuses on presenting line ratio measurements made for individual objects, we also perform a stack of all the spectra in each of our two redshift bins to aid in assessing median quantities of the sample, and push the limits of emission line non-detections.In our z ∼ 6 bin, we perform this stacking for both the medium-resolution G395M/F290LP spectra (for the N2 ratio) and the low-resolution PRISM spectra (for all other line ratios).In the z ∼ 8 bin we stack only the PRISM spectra.
Each spectrum is de-redshifted to the rest-frame and then resampled onto a uniform wavelength grid with spectral pixels 3 Å wide in the rest-frame.For the z ∼ 6 stack, we renormalise each spectrum by the Hα flux, while for the z ∼ 8 stack we renormalise by Hβ.We then take the median across all spectra in each spectral bin to form our final composite spectrum.
In order to extract emission line fluxes, we first subtract off the continuum by masking the regions of bright emission lines and fitting a spline.We then fit each emission line in the continuum-subtracted spectra with single-component Gaussian profiles.We correct the measured emission line ratios for the effects of dust-reddening in the same way as for the individual measurements, described in Section 2.4.4,using the measured Hα/Hβ ratio to derive E(B − V) for the z ∼ 6 sample and the Hγ/Hβ ratio for the z ∼ 8 sample.
We repeat the stacking process and line measurements for each sample where each spectrum is instead normalised by the UV continuum flux (taken as the median flux between 1400 and 1600 Å in the rest-frame spectrum).Finally, we repeat these measurements removing the target with the highest average noise level.Although the ratios measurement from these alternative stacks do not vary so much as to change any of the key findings in this paper, the change is larger than the statistical uncertainties on the measurements.For the final dust-corrected emission line ratios reported in Table A.2, we adopt the values measured from the Balmer-emission-normalised stacks.However, rather than adopting the statistical uncertainty, we report the systematic uncertainty, which we take as the maximum difference between the values derived from any two of the stacking regimes described above.

Comparison samples
In this section, we compare the emission line ratios measured for our JADES z ∼ 6 and z ∼ 8 samples to various literature samples in several diagnostic diagrams.To form our baseline z ∼ 0 comparison sample, we use 0.03 < z < 0.1 galaxies from SDSS MPA-JHU catalogs4 (Aihara et al. 2011) to generate 2D PDFs.We consider only galaxies with spectra flagged as 'reliable' in the public release catalog, and select for only galaxies with prominent emission lines by making a signal-to-noise cut of S /N > 40 on Hα.For diagnostic diagrams that include the [Ne iii] λ3869 line, we additionally apply a cut of S /N λ3869 > 5.
In the N2-BPT, S 2-VO87, and O1-VO87 diagrams we plot all sources that pass these criteria (including any AGN), however for all other plots, we additionally select for only star-forming galaxies using the Kauffmann et al. (2003) criterion.
We consider z ∼ 2 galaxies from the MOSDEF survey public release emission line catalog5 (Kriek et al. 2015;Reddy et al. 2015).For this sample we impose signal-to-noise cuts of S /N > 5 on each emission line and plot only detections (i.e.we do not show any limits).We additionally select only galaxies with a sky line flag with < 0.2 for all lines considered, as recommended by the catalogue release.For lines that are flagged as being near the edge of the spectrum, we also select only those flagged as having reliable flux.For diagrams involving N2, S 2 or O1, the MOSDEF sample spans from 1.2 < z < 2.6 with median z = 2.1, while for diagrams involving only the bluer lines, it probes slightly higher redshifts: typically 2.0 z 3.6 with median of z ≈ 2.3.
A number of line ratio measurements have already been reported from NIRSpec at z 6.The JWST Early Release Observations (ERO; Pontoppidan et al. 2022) provided emission line measurements for three galaxies at z > 7.5, one of which was reported to have an extremely low metallicity of 12+log(O/H) = 6.99 ± 0.11 (Curti et al. 2023).We use these three galaxies as comparison points in oxygen-based diagrams, adopting the recently re-measured ratio values from Nakajima et al. (2023).We compare our individual galaxy measurements with individual 4.5 ≤ z ≤ 8.0 galaxies targeted in the GLASS survey, which leverages gravitational lensing (Mascia et al. 2023).Additionally, we compare to composite spectra presented in Sanders et al. (2023) and Tang et al. (2023) from the Cosmic Evolution Early Release Science (CEERS) survey.The z ∼ 5.6 and z ∼ 7.5 redshift bins in Sanders et al. (2023) are well matched in redshift to the sample presented here.Those authors report a median stellar mass of log(M * /M ) = 8.57+0.04 −0.13 for their z ∼ 5.6 composite, while they do not estimate masses in their z ∼ 7.5 bin.Tang et al. (2023) present one composite spectrum of 16 z 7 galaxies from CEERS with median redshift z = 7.7.Those authors report a median UV magnitude of M UV = −20.6.This is considerably brighter than our median magnitude of M 1500 = −18.71,indicating that the shallower CEERS spectroscopy is probing generally brighter galaxies that our JADES data.
We also consider a number of possible 'high-redshift analogs' identified from low-redshift samples by various authors.'Green pea' galaxies are a population of galaxies at z ∼ 0 with unusually high equivalent-width emission of [O iii] λ5007 (∼1000 Å) (Cardamone et al. 2009;Yang et al. 2017a).These galaxies are characterised by masses of M * ∼ 10 8.5−9.5 M and high star-formation rates (∼ 10M yr −1 ), analogous to galaxies at higher redshifts.We consider also 'blueberry' galaxies, which are a population of z ∼ 0 dwarfs with even more extreme properties than green peas, having masses of M * ∼ 10 6.5−7.5 M , very high ionisation ([O iii]/[O ii] 10) and very low metallicities (12+log(O/H) ≈ 7.1 − 7.8) (Yang et al. 2017b).We also consider samples of the most metal-poor galaxies known in the low-redshift Universe (Izotov et al. 2018(Izotov et al. , 2019)).
The z ∼ 0 comparison samples are corrected for dust assuming a Cardelli et al. (1989) law with R V = 3.1 where the A V is derived from the Hα/Hβ decrement assuming an intrinsic ratio of 2.86.Meanwhile, the z ∼ 2 MOSDEF sample is corrected in an identical manner to that used for our z ∼ 6 sample, described in Section 2.4.

Excitation properties of z ∼ 6 galaxies: BPT and VO87 diagrams
The left panel of Figure 3 shows our z ∼ 6 JADES sample plotted onto the so-called N2-BPT diagram (Baldwin et al. 1981).As described in Section 2.4.3, we find no convincing detections of [N ii] λ6583 in any of our JADES z ∼ 6 sample (purple points).Figure 3 also shows comparison samples of z ∼ 0 galaxies from SDSS (background 2D PDF), extremely-metal poor galaxies at z ∼ 0 (Izotov et al. 2018(Izotov et al. , 2019, yellow pentagons), yellow pentagons), z ∼ 2 galaxies from MOSDEF (green circles), and the z ∼ 5.6 composite CEERS spectrum from Sanders et al. (2023).
Based on the R3 ratios alone, we can already see that these JADES z ∼ 6 galaxies span the upper extreme of what is observed for galaxies in SDSS and MOSDEF, with the 18 galaxies shown here having a median R3 of 0.74, and all lying above R3 > 0.53.Considering the N2 ratio, SDSS galaxies lying along the star-forming BPT sequence with 0.72 < R3 < 0.76, have a median value of N2 median = −1.4.None of our JADES galaxies have upper limits on the N2 ratio which are inconsistent with this value at the 3-σ level, while the majority of our sample have individual upper limits of only −0.9 (see also Table 2).Beyond our individual non-detections, we find that we are unable to robustly detect [N ii] λ6583 even in a stack of all 18 spectra, placing a 3-σ upper limit of N2 < −1.44.Previous studies have widely observed that z ∼ 2 galaxies are generally offset from the z ∼ 0 N2-BPT sequence, toward higher values of N2 at fixed R3 (Kewley et al. 2013b;Steidel et al. 2014;Strom et al. 2017;Runco et al. 2022; green MOSDEF points in Figure 3).Given the anticorrelation between N2 and R3 in this sequence, and the high R3 values measured in our sample, it is not possible to characterise from individual measurements what fraction of the JADES sample would lie significantly 'above' the z ∼ 0 sequence.However the upper limit imposed from our stacked spectrum suggests the median value at this redshift might be more aligned with the z ∼ 0 sequence than measurements from z ∼ 2 samples.Sun et al. (2022a) made emission line measurements at z ∼ 6 with NIRCam slitless spectroscopy, with two modest detections of [N ii] λ6583 showing much higher N2 values of N2 = −0.77and N2 = −0.69(see also Sun et al. 2022b).However, those two galaxies were more than an order of magnitude brighter (m F444W < 25) than the JADES galaxies probed here (m F444W;median = 28.1, for galaxies with NIRCam coverage), suggesting Sun et al. (2022a) are probing more massive systems.Sanders et al. (2023) present six individual galaxies at z > 5 with detections of [N ii] λ6583 from CEERS spectroscopy, also showing much higher N2 values (N2 −1.0).Taking a composite spectrum of all 38 galaxies in the CEERS sample from 5.0 < z < 6.5 they find N2 = −1.46+0.12 −0.09 , quoting a relatively  (2003).Left: Classical BPT or 'N2-BPT' -Purple triangles provide 3-σ upper limits on the locations of our JADES galaxies since we are unable to detect [N ii] λ6583 in any of the 18 galaxies for which we have G395M/F290LP spectral coverage of the Hα complex, although the R3 ratio is very well constrained in all cases from PRISM data.Even after stacking the grating spectra of all galaxies, we do not recover a robust detection of [N ii] λ6583.Centre: In the S 2-VO87 diagram (often referred to as 'S 2-BPT'), 3/20 of our JADES galaxies show detections of [S ii] in individual spectra.The tight 3-σ upper limits we place on our non-detections highlight that there must be more than an order of magnitude in scatter in the S 2 ratio within the sample.Right: We find one tentative detection of high median stellar mass of log (M * /M ) = 8.57+0.04 −0.13 .Our inability to detect N2 in any of our sample may indicate none of our galaxies have reached the same levels of metal-enrichment as some of the brighter individual galaxies observed in Sun et al. (2022b) or Sanders et al. (2023).However, comparison with lowredshift galaxies with extremely low metallicity (Izotov et al. 2018(Izotov et al. , 2019, yellow pentagons) , yellow pentagons) highlights that the 7 hour deep G395M/F290LP spectra presented here are well short of being sufficient to probe [N ii] emission from the most metal-poor systems in the early Universe.We see in the left panel of Figure 3 that, at z ∼ 0, these most metal poor systems often show N2 ratios well below N2 < −2, and have generally milder values of R3 (0.0 < R3 < 0.5).At z ∼ 0, the sequence of R3 in galaxies is known to be double-valued with metallicity, increasing from the lowest metallicities up to a turnover around 12 + log(O/H) ≈ 8, before decreasing toward higher metallicity (e.g.Curti et al. 2020).The fact that our JADES R3 ratios are notably higher than those from the Izotov et al. (2018) and Izotov et al. ( 2019) galaxies suggests we are not probing galaxies below this turnover and are more likely sampling galaxies only down to around ∼ 0.1 × Z , where the R3 sequence is relatively flat (Curti et al. 2020).Understanding N/O abundance ratios in metal-poor galaxies in the Universe will be critical to developing chemical evolution models of galaxies (Maiolino & Mannucci 2019;Hayden-Pawson et al. 2022).This study demonstrates that probing N/O abundances in the most metal poor galaxies at z 5.5 with emission lines will likely be very challenging even with very deep spectra.
The centre panel of Figure 3 shows our JADES sample plotted onto the S 2-VO87 diagram (Veilleux & Osterbrock 1987) 6 .Because the [S ii] λλ6716, 6731 doublet is more widely separated in wavelength from Hα, we can measure this ratio directly from the much deeper PRISM spectrum7 .Consequently, we observe three [S ii] detections at greater than 3-σ.These detections are broadly consistent with lying along the upper extension of the z ∼ 0 sequence and z ∼ 2 sequence.The deeper PRISM spectra also allow us to place much tighter upper limits on the S2 ratio for the 17/20 non-detections in our sample.In some cases the 3-σ upper limits are more than an order of magnitude lower in S 2 than our three detections, suggesting there is significant scatter in S 2 ratios within the z ∼ 6 population.The composite z ∼ 5.6 CEERS spectrum from Sanders et al. (2023) yields S 2 = −1.30+0.17 −0.02 , which is broadly in line with the three JADES galaxies with measurable [S ii].This suggests that stacks based on the shallower, but wider area, CEERS spectroscopy are biased toward galaxies that are more in line with the tip of the z ∼ 0 sequence, which are likely more evolved systems.
In the right panel of Figure 3 3.We additionally show comparison with two populations of z ∼ 0 extreme starbursts: green peas (Yang et al. 2017a; orange plusses) and blueberries (Yang et al. 2017b; blue crosses).The three well-studied z > 7.5 galaxies from the SMACS 0723 ERO observations are shown as lime green pentagons, adopting the ratio values presented by Nakajima et al. (2023).Unlike the BPT and VO87 diagrams, the majority of our sample have > 3σ detections of all lines in individual spectra, revealing a large diversity within the sample.Ratios measured from stacked JADES spectra are shown as the open magenta and red diamonds (see Section 2.5).
GS+53.15613-27.77584),indicating the presence of neutral gas in this galaxy.The measured O1 ratio for this galaxy places it above the theoretical maximum starburst excitation limit from Kewley et al. (2001), which could suggest the presence of shockheated gas in the ISM (Sutherland & Dopita 2017).Interestingly, this galaxy is not one of the three galaxies with a [S ii] detection from Figure 3.However, it does have the second lowest O32 ratio of galaxies in this sample (Table 2).Overall, the individual O1 constraints presented here are not strongly constraining relative to the Kewley et al. (2001) demarcation line.The measurement from our stacked spectrum suggests the median O1 in this sample could be consistent with the upper extension of the z ∼ 0 star-forming sequence.However, we showed above that there is significant diversity in S 2 ratios within our sample.If this diversity also holds for O1, it is clear that larger samples probing even deeper would be required to characterise the the sample properties of O1 at z ∼ 6.

Excitation properties from
In Figure 4 we instead consider the R2 = [O ii] λλ3726, 3729 / Hβ ratio.The longer wavelength baseline of this ratio has the drawback of being much more sensitive to the reddening correction applied compared to the N2-BPT, S 2-VO87, and O1-VO87 diagrams discussed in Section 3.2.However, the shorter rest-frame wavelength of these emission lines means that the R2 − R3 dia-gram can be used to explore the ISM conditions of galaxies up to z ∼ 9.5 with JWST/NIRSpec.[O ii] is typically intrinsically much brighter than [N ii], [S ii], or [O i].Moreover, like [S ii] and [O i], ratio measurements can be taken directly from the PRISM spectrum.As a result, not only does Figure 4 include the additional six galaxies from the z ∼ 8 sample, we also have a much higher detection rate, with 24/27 galaxies in our combined sample having detections of [O ii].
We see immediately in Figure 4 that there is more than an order of magnitude in scatter in the R2 ratio, suggesting there is significant diversity in the ISM conditions of galaxies within this sample.We find a median value of R2 = −0.23 with standard deviation of 0.38 among our 24 detections, significantly offset from the median value from MOSDEF of R2 = 0.38.Instead, we find that the range of R2 values observed for this sample are much more aligned with those observed in extreme z ∼ 0 dwarf starbursts such as green peas (Cardamone et al. 2009;Yang et al. 2017a; orange plusses) and blueberries (Yang et al. 2017b; blue crosses).
Both R2 and R3 are known to form a double-valued sequence in metallicity (e.g.Curti et al. 2020).In z ∼ 0 galaxies, R3 does not change strongly with metallicity between ∼ 0.1 − 0.5 × Z .On the other hand, R2 turns over at higher metallicities between ∼ 0.5 − 0.8 × Z , at a maximum value of R2 ≈ 0.5 (Curti et al. 2020).The combination of highly consistent R3 ratios but wildly varying R2 ratios indcates that our sample likely falls into the metallicity range of ∼ 0.1 − 0.5 × Z .We have five galaxies in our JADES sample with R2 > 0.0, overlapping with the region inhabited by typical z ∼ 2 MOSDEF galaxies.These may represent systems which are more evolved and more metal-rich.Below its metallicity turnover, R2 values in the z ∼ 0 sample from Curti et al. (2020) drop steeply with decreasing metallicity.The lowest R2 values measured in our sample lie below R2 < −0.75, below the range probed by the Curti et al. (2020) calibration.However we note that the low R2 portion of our sample exhibits a fair resemblence to blueberry galaxies from Yang et al. (2017b) which have been show to have very low-metallicties 0.1 × Z .However, given that R2 is strongly affected by the ionisation parameter, inferring metallicity values from this ratio is likely not robust.Nonetheless, there is clearly a large range of metallicity being probed in this sample.
We see no evidence for evolution of the R2 ratio between our z ∼ 8 and z ∼ 6 samples, albeit with only six z ∼ 8 galaxies.We note that all of our R2 > 0 galaxies are from the z ∼ 6 sample; however, the measurements from stacked spectra in the two samples actually yield a higher value for the z ∼ 8 sample in Figure 4 (open diamonds), likely affected by uncertainty due to small sample size.This suggests that there is no rapid time evolution in this diagram at this epoch and that the scatter is instead driven by sample diversity.

Ionisation-excitation diagrams: R23-O32 & R23-Ne3O2
Figure 5 shows the R23-O32 diagram with our JADES sample plotted against the same literature comparison samples shown in Figure 4. R23 is often taken as an indication of total excitation, as it encompasses emission lines from both singly and doubly ionised oxygen, while the O32 ratio is sensitive to ionisation parameter (e.g.Kewley et al. 2019).
It is in this diagram that the offset of z 6 from typical z ∼ 0 and z ∼ 2 galaxies is perhaps most striking: O32 ratios are clearly much higher in JADES z ∼ 5.5 − 9.5 galaxies compared to MOSDEF galaxies at z ∼ 2. Across our combined sample we find a median O32 value of 0.98 with a standard deviation of 0.36, with 18 of our galaxies have robust detections showing O32 > 0.75.This is much higher than MOSDEF galaxies plotted in Figure 5 which have a median O32 value 0.18 and are already tracing the upper end of the z ∼ 0 relation.This indicates that galaxies across our sample exhibit very high ionisation parameters, much larger than what are seen in typical galaxies at z ∼ 0 − 3.However, we can see in Figure 5 that our JADES sample largely overlaps with the region spanned by blueberries (Yang et al. 2017b) and green peas (Yang et al. 2017a), which are examples of z ∼ 0 galaxies with extreme star-formation activity.Comparable O32 ratios have also been observed in z ∼ 2 − 4 Lyman-α emitters (Nakajima et al. 2016;Tang et al. 2021) and z ∼ 1.3 − 2.4 extreme [O iii] emitters (Tang et al. 2019) therefore likely that the high O32 ratios observed in our sample are indicative of intense star formation activity.
We find that our z ∼ 6 sample is significantly offset from the z ∼ 5.6 stack presented in Sanders et al. (2023), with our composite O32 value found to be more than a factor of four larger.Indeed where individual detections from CEERS are reported in Sanders et al. (2023) they typically have O32 < 0.75.Mirroring what was observed in our R2 sequence from Figure 4, we find six galaxies with much lower O32, more in line with the composite spectrum from Sanders et al. (2023).Given the resemblance between the upper end of our O32 sample and dwarf starbursts such as blueberries, it is likely that this portion of our sample is probing higher star-formation rate and lower-metallicity galaxies (analogous to blueberries), while the handful of lower O32 galaxies presented here in alignment with the Sanders et al. (2023) z ∼ 5.6 composite could be probing more evolved systems.
On the other hand, the composite CEERS values at z ∼ 7.5 and z ∼ 7.7, from Sanders et al. (2023) and Tang et al. (2023) respectively, are in closer agreement with our composite spectra, albeit still offset toward slightly higher R23.This could suggest that the fraction of more evolved, metal-enriched systems is increasing rapidly from z ∼ 8 to z ∼ 6, despite the bulk of our UV-selected z ∼ 6 sample showing no significant evolution from our z ∼ 8 sample.Addressing this question in more detail requires larger samples of deep spectroscopy where properties of individual galaxies can be isolated.
By leveraging gravitational lensing, the measurements of 4.5 z 8 galaxies from Mascia et al. (2023) can more readily probe fainter and lower mass systems.These are similarly offset from the Sanders et al. (2023) z ∼ 5.6 composite and largely overlap with our JADES sample, spanning a large range in R23.Two galaxies are reported by Mascia et al. (2023) as having R23 0.6 despite having high O32 comparable with our JADES sample.This indicates that gravitational lensing may be important for study low metallicity galaxies (below the R23 turnover) and not simply ultra-deep observations like those presented here.
We note that, as with Figure 4, this diagram relies on ratios with long baselines in wavelength and that dust reddening can have a significant effect.The dust-correction (described in Section 2.4.4) is most significant on the O32 ratio here, with larger corrections serving to reduce O32 ratios.From the corrections we apply, all but two galaxies have corrections that shift O32 by less than 0.15.The offset between the Sanders et al. (2023) z ∼ 5.6 composite and our z ∼ 6 stack cannot be explained by the dust correction; 'correcting' our value down to O32 ≈ 0.5 would require E(B − V) 1, significantly larger than our measured value of E(B − V) = 0.06.Nonetheless, in Figure 6 we show our JADES z ∼ 6 and z ∼ 8 samples with Ne3O2 replacing O32 to provide the ionisation axis, which provides the bonus that Ne3O2 is essentially unaffected by dust extinction and the associated correc-tion uncertainties8 .Despite the lower detection rate, and generally much larger measurement uncertainties associated with switching from [O iii] λ5007 to the intrinsically much fainter [Ne iii] λ3869, we see the same story reflected here with the JADES sample exhibiting large scatter (median Ne3O2 = −0.11with standard deviation 0.29) and generally offset from SDSS and MOSDEF galaxies, instead being more aligned with blueberry dwarf starbursts.This again supports the conclusion that our JADES sample is probing a diverse sample of galaxies, many of which exhibit very high ionisation parameters analogous to those observed in local compact starbursts.Tang et al. (2023) presented Ne3O2 measurements from their z ∼ 7.7 composite CEERS spectrum.These are marginally lower than our Ne3O2 composite values.As in Figure 5, the larger offset is again the 0.1 dex difference in R23.This offset could reflect an increase in metallicity from our sample to the Tang et al. (2023) galaxies which are reported has having typically brighter UV magnitudes (median M UV = −20.6 from Tang et al. 2023).
Finally, we note that [Ne iii] λ3869 emission may become important for confirming spectroscopic redshifts at z 10 where [O iii] λλ 4959, 5007 is no longer observable within the wavelength coverage of NIRSpec.We find the [Ne iii] λ3869 emission to be brighter than the [O ii] doublet in at least 8/27 galaxies in our sample.Given that these cases should represent the more extreme examples (highest ionisation parameters and lowest metallicities), [Ne iii] λ3869 may be one of the most accessible emission lines in searches for the highest redshift galaxies.However, existing deep spectroscopic observations of such galaxies have already demonstrated that observing any emission lines at those redshifts can be challenging (Curtis-Lake et al. 2022;Robertson et al. 2022).

ISM conditions in the Epoch of Reionisation
Section 3 showed that our z ∼ 5.5 − 9.5 sample from JADES exhibits significant diversity in ISM conditions.The measured ratios highlight a clear difference between this high-redshift sample and typical galaxies at z ∼ 0 − 3, although parallels can be drawn instead with extreme starburst populations of low-redshift galaxies, such as blueberries and green peas.In this section we use photoionisation models to explore the parameter space spanned by our sample.
Figure 7 highlights that photoionisation models driven by SEDs of young stars can reproduce the parameter space populated by our JADES galaxies with only fairly standard assumptions.This indicates we do not need to invoke alternative heating sources (such as AGN or shocks) to explain our findings.We see from this simple comparison that the large scatter observed in our emission line ratio sample likely reflects a diversity in metallicity and ionisation parameter among these galaxies.Our most extreme galaxies (i.e.highest O32) are only reproduced by the models with very low-metallicity (0.07 − 0.15 × Z ) and very high ionisation parameters (−2.0 log U −1.0).These ionisation parameters are signifcantly higher than those typically observed in z ∼ 0 H ii regions (Dopita et al. 2000) or z ∼ 2 − 3 galaxies (Sanders et al. 2016).We note that we have not explored the effect of harder radiation fields, which would increase O32 at a fixed ionisation parameter (e.g.Sanders et al. 2016).Indeed, many studies at z ∼ 2−3 have argued that harder radiation fields, driven by super-solar α/Fe abundance ratios, are a dominant effect in driving the z ∼ 2 N2-BPT offset (Steidel et al. 2016;Strom et al. 2017;Topping et al. 2020b,a;Cullen et al. 2021).Given the long timescale of iron enrichment via Type Ia supernovae, one would of course expect that these z 5.5 galaxies would be α-enhanced, and models with non-solar abundance ratios would be more appropriate for parameter estimation.
We can see from these comparisons that our sample likely spans quite a range in metallicity.At lower O32 values, we see a turnover in the models where metallicity and ionisation parameter become highly degenerate.However, model metallicities of at least 0.45 × Z are required to approach our highest R2 and R23 values.The R2, R3, and R23 ratios are all known to exhibit double-valued relations (e.g.Curti et al. 2020) in metallicity.We tentatively observe that our highest O32 galaxies tend to have the the lowest R23, which reverses the trend seen in SDSS galaxies (grey 2D PDF), suggesting we are sampling the turnover of this relation.Some of our 'mildest' galaxies (i.e.highest R2, lowest O32) have R2 ≈ 0.4 and R23 ≈ 1.0 which are broadly consistent with the turnover observed for these ratios in z ∼ 0 galaxies from Curti et al. (2020).This suggests these galaxies are the most chemically evolved systems in our sample and may have metallicities as high as 0.3 − 0.5 × Z .Indeed, these are also the points that are in best agreement with the Sanders et al. (2023) z ∼ 5.6 composite spectrum.This suggests the shallower CEERS spectroscopy is probing on average more chemically evolved systems than those targeted in JADES.
Overall, it is clear from this simple model comparison that a range of ISM conditions is required to explain the diversity of emission line ratios observed in this sample.As noted in Section 3, our sample largely spans the same regions of emission-line-ratio-space inhabited by low-redshift extreme starbursts such as blueberries and green peas.Thus, the emerging picture is that the galaxies in our sample are likely lowmetallicity (∼ 0.1 − 0.3 × Z ) systems undergoing intense periods of star-formation driving strong radiation fields with high ionisation parameters (−2.0 log U −1.0) in the ISM.

High O32: Implications for escape fraction and reionisation
The high O32 ratios presented in Figure 5 could have implications for interpreting the ionising photon escape fraction of these z ∼ 5.5 − 9.5 galaxies.When and how the Universe was reionised remains one of the most important open questions in high-redshift astrophysics.The efficiency with which ionising photons are produced within galaxies and the fraction of these  7. JADES z ∼ 5.5 − 9.5 galaxies compared to photoionisation models.Measured ratios are shown as purple points, as in Figures 3 -6 except that all galaxies are plotted in the same colour, regardless of redshift.SDSS galaxies are shown in the same manner as Figures 3 -6.Each line shows a set of models with constant metallicity with values of Z = [1,2,3,4,6,8,10,14,20,30, 40] ×10 −3 increasing from black to yellow.Assuming Z = 0.0134 (Asplund et al. 2009), these correspond to Z/Z = 0.07, 0.15, 0.22, 0.30, 0.45, 0.60, 0.75, 1.0, 1.5, 2.2, 3.0.Ionisation parameter varies along the grid line, increasing with marker size, with values of log U = −4.0,−3.5, −3.0, −2.5, −2.0, −1.5, −1.0, −0.5.Details of the photoionisation models are described in Section 4.1.
which escape ( f esc ) are two important quantities required to understand the contribution of galaxies to this process.Unfortunately, f esc can never be directly measured within the epoch of reionisation itself because any light escaping the galaxy blueward of the Lyman-limit will interact with intervening neutral hydrogen.Thus indirect probes of f esc are highly sought after.
High O32 has been proposed as being an indicator of higher f esc (e.g.Nakajima & Ouchi 2014;Faisst 2016).The reasoning is that the high O32 ratio selects for highly ionised systems which have a higher likelihood of having density-bounded channels through which ionising photons can escape.However, some studies have shown that O32 does not necessarily correlate well with f esc ; instead finding that results are very dependent on the observed line of sight, and that high O32 can be produced in systems without Lyman continuum (LyC) escape (Paalvast et al. 2018;Bassett et al. 2019;Katz et al. 2020).Nonetheless, samples of Lyman continuum (LyC) leaking galaxies at low-redshift generally show that the fraction of galaxies with high f esc increases toward higher O32, even if the correlation is not tight (e.g.Izotov et al. 2016;Flury et al. 2022).Indeed, the limited number of known LyC leakers at higher redshift typically have high O32 ratios (Vanzella et al. 2016(Vanzella et al. , 2020(Vanzella et al. , 2022;;Fletcher et al. 2019;Saxena et al. 2022).Extreme [O iii] emitters at z ∼ 2 − 3 have also been shown to exhibit high O32 ratios ([O iii]/[O ii] ∼ 10), while also preferentially selecting for high equivalent-width Lyman-α emission (Tang et al. 2019(Tang et al. , 2021)).This evidence therefore seems to suggest that high O32 is a necessary but not sufficient condition for high f esc .
We have shown in Section 3.4 that our JADES z ∼ 5.5 − 9.5 galaxies have characteristically high O32 ratios, with 21/27 of the sample having [O iii]/[O ii] > 5 and a median value of 9.4.Owing to the large scatter in O32-f esc relations, we do not at-tempt to derive constraints on the f esc of these galaxies.However, the high fraction of our sample that meets the criterion of high O32 suggests that the fraction of potential LyC leakers might be increasing to high redshift.

Conclusions
We have presented rest-frame optical emission line ratio measurements in 27 individual galaxies at 5.5 < z < 9.5 spanning over −17.0 < M 1500 < −20.4,observed as part of the 'Deep' tier of the JADES program, which represent the deepest integrations yet taken with JWST/NIRSpec.We find that galaxies in this sample occupy regions of line-ratio space that are offset from those inhabited by 'typical' galaxies at z ∼ 0 or z ∼ 2, although generally aligned with more extreme low-redshift populations such as 'blueberry' and 'green pea' dwarf starbursts (Yang et al. 2017a,b).
The high signal-to-noise achieved in individual spectra from our 28 hour deep PRISM/CLEAR integrations reveal significant intrinsic scatter in the observed line ratios, with our JADES z ∼ 5.5 − 9.5 galaxies spanning more than an order of magnitude in [O ii]/Hβ and [O iii]/[O ii] ratios.This highlights the diversity in ISM conditions of galaxies at this epoch.
Despite the depth of our observations, we find that low ionisation lines such as [N ii] λ6583, [S ii] λλ6716, 6731 and [O i] λ6300 are very challenging to detect in these systems.We find no convincing detections of [N ii] λ6583 in 7 hour deep G395M/F290LP spectra for any of the 18 galaxies in which we have adequate spectral coverage.Stacking these 18 spectra is still unable to deliver a nitrogen detection.However, the 3-σ upper limit suggests our sample median is not offset from the z ∼ 0 N2-BPT sequence contrary to what has been observed in z ∼ 2 samples (e.g.Steidel et al. 2014).In our PRISM spectra we find three modest detections of [S ii] λλ6716, 6731, while the tight upper limits we derive on other targets suggests more than an order of magnitude in intrinsic scatter in the S 2 ratio within the sample.The fainter [O i] λ6300 line yields only one tentative detection.
The bluer [O ii] λλ3726, 3729 and [Ne iii] λ3869 are much more readily detected throughout the sample.Measured R2, O32 and Ne3O2 ratios particularly highlight the diversity in this sample; we find median values and standard deviations of R2 = −0.23 ± 0.38, O32 = 0.98 ± 0.36, and Ne3O2 = −0.03± 0.27.We find that these ratios are generally consistent with very high ionisation parameters up to log U ∼ −1.0.The sample appears to span a significant range in metallicity of ∼ 0.1−0.3×Z, depending on metallicity calibrator assumptions.These ratios are found to be generally consistent with extreme low-redshift populations, particularly green pea and blueberry dwarf starburst samples.This is consistent with the suggestion that these z ∼ 5.5 − 9.5 galaxies are low-mass, low-metallicity galaxies undergoing periods of rapid star formation, driving strong radiation fields.
We compare our sample to composite spectra at z ∼ 5.6 and z ∼ 7.5 from CEERS (Sanders et al. 2023;Tang et al. 2023).We find significant offset between our z ∼ 6 sample median and the z ∼ 5.6 CEERS composite.We instead observe good agreement between this z ∼ 5.6 CEERS composite measurement and the low O32 end of our sample.Assuming these galaxies represent the more evolved, metal-enriched extent of our sample, it is clear that the shallower CEERS spectra are preferentially picking up on the more evolved portion of the z ∼ 6 population.We find no evidence for significant evolution between our z ∼ 6 and z ∼ 8 samples (albeit with a very small sample).Both our z ∼ 6 and z ∼ 8 sample medians are found to be in somewhat closer agreement with the CEERS z ∼ 7.5 composite spectra, which showed large offset from the CEERS z ∼ 5.6 composite spectrum.This could indicate that the fraction of moderately metal-enriched systems increases rapidly from z ∼ 8 to z ∼ 5.5, even if our sample median does not show evidence of significant redshift evolution.
Within its first six months of observations, JWST/NIRSpec has already shed considerable light on the ISM conditions of z 5.5 galaxies.The diversity of ISM conditions presented in this work highlights the need for observations that span the full range of galaxy properties.This will require a combination of deep observations reaching down to even fainter brightness limits than those presented here, as well as wider area surveys probing rarer, more evolved systems at early epochs.Characterising the ISM conditions in larger samples of galaxies across this broader parameter space will help chart galaxy assembly across the first ∼1 Gyr of cosmic time.F277W, 0 .116 in F356W 11 ).Thus, at the wavelengths considered in this work, the assumption of a point source profile should be reasonable for all but a few of our largest targets.The majority of the emission line ratios studied in this work have relatively short wavelength baselines, for which the precise slitloss corrections will make a negligible difference.The O32 and R23 ratios are the notable exceptions to this.However, even before slitloss corrections are applied, Figure B.2 demonstrates that, on average, there is only a ∼10 % effect across the full wavlength range.Given that [O ii] λλ3726,3729 and [O iii] λ5007 are only separated by ∼1.4 µm in the observed frame, the effect would be even smaller than this quoted 10 %.Furthermore, the point-source slitloss corrections we do apply should be adequate for most sources.Thus, the final flux calibration uncertainty on the O32 ratio is likely well below 10 % which is smaller than the quoted measurements uncertainties throughout this paper.

Appendix B.1: Comparison of grating flux calibration
This paper almost exclusively makes use of the PRISM/CLEAR mode.This is because, at redshifts beyond z > 5, the restframe optical emission lines fall beyond λ obs 2.5µm where the PRISM/CLEAR spectral resolution is R ∼ 200 − 300, meaning most of these lines are well resolved.The notable exception to this is the [N ii] λ6583 / Hα ratio, for which we require the R ∼ 1000 G395M/F290LP grating.Although reported ratios are only ever measured within the same spectral observation, we briefly consider the relative flux calibration of the PRISM/CLEAR and G395M/F290LP.Figure B.2 shows the ratio of Hα fluxes measured in the G395M/F290LP to those measured in the PRISM/CLEAR mode for galaxies in this sample which have spectral coverage in both settings.Note that the fluxes measured here are not those reported in the main text, where the continuum is modelled with stellar templates.Rather the continuum here is modelled with a simple spline fit.This is done to minimise any model-dependent effects in this comparison.The dashed magenta line shows the median value of 1.16, indicating that fluxes measured from the grating are systematically larger.This effect has been noted in Bunker et al. (2023).Note that in the PRISM/CLEAR, the reported Hα flux would be blended with [N ii], whereas these lines are resolved in the G395M/F290LP measurements.However, as pointed out in Section 3, [N ii] λ6583 is not clearly detected in any of these galaxies, and this would contribute at the < 10 % level.Furthermore, the effect of this blending would be to increase the PRISM/CLEAR flux, but it is the G395M/F290LP flux that is determined to be larger.As outlined above, the flux calibration in the PRISM/CLEAR mode is in good agreement with broadband photometry.Thus, we consider that the PRISM/CLEAR flux calibration is likely reliable, while the grating flux may be overestimated.However, given we only make use of the G395M/F290LP data to measure the ratio of two lines very close in wavelength within the same observation, the systematically high fluxes recorded should have no effect on any of the findings reported here.

Fig. 1 .
Fig. 1.Example spectra from five galaxies included in our sample.Left: Low-resolution PRISM/CLEAR spectra shifted to the rest-frame according to the observed redshift.Flux is shown in f λ but has been renormalised relative to the peak flux from the [O iii] λ5007 line.Solid horizontal lines show the zero flux of each spectrum.The shaded region around the 'zero' line indicates the 1-σ noise spectrum in the same renormalised units.We show only the subset of spectral coverage that includes the key emission lines used in this study.Vertical dotted lines show the expected centroids of the [O ii] λ3727, [Ne iii] λ3869, [Ne iii] λ3967, Hδ, Hγ, [O iii] λ4363, Hβ, [O iii] λ4959, [O iii] λ5007, [O i] λ6300, Hα, [N ii] λ6583, [S ii] λ6716, and [S ii] λ6731 lines respectively (moving left to right).Right: Zoom in on the Hα+[N ii] complex as observed in the G395M/F290LP grating spectra of these same galaxies.Observed spectra are shown in the same way as for the left panel (except renormalised to Hα, rather than [O iii]).Vertical dotted lines show the expected centroids of Hα and [N ii] λ6583.
Fig.2.Top: Redshift histogram of galaxies observed in this study.We divide the sample into two sub-samples 'z ∼6' and 'z ∼ 8', based simply on cutting at z = 7 since this is the redshift at which Hα is redshifted beyond the spectral coverage of NIRSpec.Thus, for JADES 'z ∼ 8' we only have coverage of lines from [O ii] λλ3726, 3729 to [O iii] λ5007, while for JADES 'z ∼ 6' we can in principle observe all the rest-frame optical lines out to [S ii] λλ6716, 6731.The orange hatched histogram indicates galaxies for which we can constrain the [S ii] λλ6716, 6731 flux.Green hatching is the same but for [N ii] λ6583.Bottom: Redshift vs. M 1500 for the full z ∼ 5.5 − 9.5 sample.We derived M 1500 directly from the PRISM spectra as described in Section 2.4.1.

Fig. 3 .
Fig.3.BPT and VO87 diagrams showing the comparison between JADES z ∼ 6 galaxies and samples across various redshifts.Individual measurements are shown as solid purple points, while values dervied from stacked JADES spectra are shown as the open magenta diamonds.The background grey 2D PDF shows z ∼ 0 galaxies from SDSS.Extremely metal-poor z ∼ 0 galaxies fromIzotov et al. (2018) andIzotov et al. (2019) are shown as yellow pentagons.MOSDEF z ∼ 2 galaxies are shown as green circles(Kriek et al. 2015).The z ∼ 5.6 composite CEERS spectra presented bySanders et al. (2023) are shown as the blue open circles.The solid navy lines in each panel show the theoretical maximum starbust lines fromKewley et al. (2001), while the the dotted line in the left panel shows the empirical demarcation derived byKauffmann et al. (2003).Left: Classical BPT or 'N2-BPT' -Purple triangles provide 3-σ upper limits on the locations of our JADES galaxies since we are unable to detect [N ii] λ6583 in any of the 18 galaxies for which we have G395M/F290LP spectral coverage of the Hα complex, although the R3 ratio is very well constrained in all cases from PRISM data.Even after stacking the grating spectra of all galaxies, we do not recover a robust detection of [N ii] λ6583.Centre: In the S 2-VO87 diagram (often referred to as 'S 2-BPT'), 3/20 of our JADES galaxies show detections of [S ii] in individual spectra.The tight 3-σ upper limits we place on our non-detections highlight that there must be more than an order of magnitude in scatter in the S 2 ratio within the sample.Right: We find one tentative detection of [O i] λ6300 in our O1-VO87 diagram.
Fig. 4. The R2 − R3 diagram showing our two JADES sub-samples at z ∼ 6 (purple points) and z ∼ 8 (maroon points).Comparison samples from SDSS and MOSDEF are shown as in Figure3.We additionally show comparison with two populations of z ∼ 0 extreme starbursts: green peas(Yang et al. 2017a; orange plusses) and blueberries(Yang et al. 2017b; blue crosses).The three well-studied z > 7.5 galaxies from the SMACS 0723 ERO observations are shown as lime green pentagons, adopting the ratio values presented byNakajima et al. (2023).Unlike the BPT and VO87 diagrams, the majority of our sample have > 3σ detections of all lines in individual spectra, revealing a large diversity within the sample.Ratios measured from stacked JADES spectra are shown as the open magenta and red diamonds (see Section 2.5).

Fig. 5 .
Fig. 5.The R23-O32 diagram, probing the ionisation and excitation of the ISM.The JADES sub-samples at z ∼ 6 (purple points) and z ∼ 8 (maroon points) lie at generally high excitation an high ionisation.Comparison samples, including SDSS, MOSDEF, green peas, blueberries, and SMACS 0723, are shown as in Figure 4. z ∼ 4.5 − 8 galaxies from GLASS are shown as grey diamond (Mascia et al. 2023).Blue and orange open circles show the z ∼ 5.6 and z ∼ 7.5 composite CEERS spectra from Sanders et al. (2023), while the z ∼ 7.7 composite CEERS spectrum from (Tang et al. 2023) is shown as the green open hexagon.Stacked JADES spectra are shown as the open magenta and red diamonds.

Fig. 6 .
Fig. 6.The R23-Ne3O2 diagram, probing the ionisation and excitation of the ISM.R23-Ne3O2 closely resembles the R23-O32 diagram with the additional advantage that the [Neiii] λ3869 / [O ii]λλ3726,329 ratio is less sensitive to the applied dust-correction, or any wavelength-dependent flux calibrations.JADES galaxies are plotted as in Figure5.Fewer literature samples present [Ne iii] ratios so we show only SDSS, blueberries, MOSDEF, and the z ∼ 7.7 composite CEERS spectrum from(Tang et al. 2023).These are each plotted in the same way as in previous figures.
Fig. B.2.Comparison of Hα flux measurements made from lowresolution PRISM/CLEAR and R ∼ 1000 G395M/F290LP measurements for galaxies in this sample.Fluxes reported in the R ∼ 1000 data are systematically higher by a median factor of 1.16.

Table 1 .
(Bunker et al. 2023tegration times on each of the targets presented in this sample.The six targets above the horizontal rule comprise the JADES 'z ∼ 8' sample and do not have coverage of any of Hα, [N ii] or [S ii] in either PRISM or G395M/F290LP data.The remainder are referred to as the 'z ∼ 6' sample.'NIRSpec_ID'refers to the integer ID used in the JADES data release(Bunker et al. 2023).
a Target does not have coverage of the Hα+[N ii] complex in G395M/F290LP observation.

Table 2 .
Table of dust-corrected emission line ratios measured for each of the 27 galaxies in our sample.Line ratio names are defined in Equation 1. Non-detections are quoted as 3-σ upper limits.a Target does not have coverage of the Hα+[N ii] complex in G395M/F290LP observation.
b Target was not observed in our G395M/F290LP MSA configuration.c Ne3O2 ratio is unconstrained as neither [O ii] nor [Ne iii] were robustly detected.A. J. Cameron et al.: ISM conditions from deep z 6 spectra (Rieke et al. 2023) radii from the JADES public catalog(Rieke et al. 2023)for the subset of our sample for which NIRCam imaging was available.'NIRCam ID' gives this cross-matched ID from those imaging catalogs.The horizontal rule marks the FWHM of the NIRCam PSF in the F356W filter (0 .116); i.e., targets below this rule have half-light radii larger than the NIRCam PSF.