Exploring the Gas-phase Metallicity Gradients of Star-forming Galaxies at Cosmic Noon

We explore the relationships between the [O/H] gas-phase metallicity radial gradients and multiple galaxy properties for 238 star-forming galaxies at 0.6 < z < 2.6 selected from the CANDELS Lyα Emission at Reionization survey with stellar mass 8.5<logM*/M⊙<10.5 . The gradients cover the range from −0.11 to 0.22 dex kpc−1, with the median value close to 0. We reconstruct the nonparametric star formation histories (SFHs) of the galaxies with spectral energy distribution modeling using Prospector with more than 40 photometric bands from the Hubble Space Telescope, Spitzer, and ground-based facilities. In general, we find weak or no correlations between the metallicity gradients and most galaxy properties, including the mass-weighted age, recent star formation rate, dust attenuation, and morphology as quantified by both parametric and nonparametric diagnostics. We find a significant but moderate correlation between the gradients and the “evolutionary time,” a temporal metric that characterizes the evolutionary status of a galaxy, with flatter gradients observed in more evolved galaxies. Also, there is evidence that galaxies with multiple star formation episodes in their SFHs tend to develop more negative gas-phase metallicity gradients (higher [O/H] at the center). We conclude that gas kinematics, e.g., radial inflows and outflows, is likely an important process in setting the gas-phase metallicity gradients, in addition to the evolution of the SFH radial profile. Since the gradients are largely independent of the galaxies’ physical properties and only weakly dependent on their SFH, it would appear that the timescale of the gas kinematics is significantly shorter than the evolution of star formation.


INTRODUCTION
Galaxies form and evolve in a gas ecosystem, with metals in the gas phase continuously produced, transported and consumed throughout their evolution.Star formation, supernovae and neutron star mergers (e.g, Thielemann et al. 2017) contribute to the local metallic-ity enrichment of the interstellar medium (ISM).Meanwhile, the stellar remnants and the evolution of low-mass stars lead to the consumption of metals (e.g, Kobayashi et al. 2006;Sukhbold et al. 2016).The distribution and evolution of gas-phase metallicity gradients provide important insights into the star formation history (SFH) of galaxies as well as how they evolve over time.
If a galaxy were a closed box without any gas exchange, and thus metal exchange with the external environment, one would expect that the gas-phase metallicity gradient would tend to follow the radial distribu-tion of past generations of stars as they formed.However, metals can be lost through galactic winds and outflows (e.g, Rupke 2018;McQuinn et al. 2019;Wang et al. 2022a), accreted from the surrounding intergalactic medium (IGM) and/or the CGM (e.g, Cen & Chisari 2011;Muratov et al. 2017;Scholz-Díaz et al. 2021;Sánchez Almeida et al. 2022), and re-distributed either through internal dynamical mechanisms or through interactions with other galaxies such as mergers (e.g, Boco et al. 2021;Mapelli et al. 2021).All these factors can effectively change the distribution and evolution of metals in galaxies, and thus add to the complexity of interpreting the shape of the gas-phase metallicity gradients and their relationship with the galaxies' stellar mass and star-formation rate (SFR) gradients.
Gas-phase metallicity gradients have been extensively studied through observations in the local universe since the 1970s (Searle 1971;Shaver et al. 1983).In the local universe, early works showed that the metallicity decreases approximately exponentially with the galactocentric radius (namely showing a negative gradient) in disk-like galaxies (e.g.Zaritsky 1992;Henry & Howard 1995).With the development of contemporary instruments and large integral-field spectroscopic surveyse.g.CALIFA (Sánchez et al. 2012) and MaNGA (Bundy et al. 2015) -the measurements of gas-phase metallicity gradients have been dramatically expanded.Local galaxies are consistently reported to have negative radial gradients in gas-phase metallicity, with typical values ranging from 0 to -0.1 dex kpc −1 (e.g, Kewley et al. 2019;Kreckel et al. 2019;Sánchez 2020).Flat or positive gradients are also found occasionally, typically for low-mass and late-type galaxies where the estimation of age and metallicity is highly degenerated, or for galaxies in dense environments (e.g, Lian et al. 2019;Sánchez 2020).By measuring the gas fraction and the local escape velocity, Barrera-Ballesteros et al. (2018) were able to reproduce the observed radial metallicity gradients for a large sample of local disk galaxies.This result indicated that the gas-phase metallicity is a consequence of local star formation history, and momentum-driven outflows help shape the internal metallicity of star-forming galaxies.
Measuring metallicity gradients in high-z galaxies is challenging because galaxies were smaller in the past and the available sensitivity and angular resolution are still limited compared to observations in the local universe.To quantify the systematics encountered by observing poorly-resolved galaxies, Acharyya et al. (2020) used synthetic observations of numerical galaxy simulations to demonstrate that the observed metallicity gradient is artificially flattened at low spatial resolu-tion and low signal-to-noise ratio.For this reason, the sizes of the samples at high redshift with robust measurements of metallicity gradients have traditionally remained small.Thanks to the sensitive and highresolution near-infrared spectrographs onboard the Hubble Space Telescope and, more recently, with JWST, as well as many large ground-based facilities, the observations of gas-phase metallicity gradients in high-z galaxies have been expanding in the last decade (Wuyts et al. 2016;Wang et al. 2017Wang et al. , 2019;;Curti et al. 2020;Wang et al. 2020;Gillman et al. 2021;Simons et al. 2021;Wang et al. 2022b).Besides the metallicity gradient itself, these studies have also investigated the relationship between the metallicity gradient and galaxy properties, including the stellar mass, size, and morphology.
Gas-phase metallicity gradients have been reported to correlate with the stellar mass especially in the local universe, but how exactly the two properties are correlated with each other remains inconclusive.Some suggested that high-mass (typically log(M/M ⊙ ) ≳ 10.5) galaxies tend to have slightly more negative metallicity gradients compared with low-mass galaxies (e.g.Belfiore et al. 2017;Carton et al. 2018;Poetrodjojo et al. 2018).In contrast, Mingozzi et al. (2020) reported a slight steepening in metallicity gradient with stellar mass at log(M/M ⊙ ) = 9∼10.25,and then a flattening towards log(M/M ⊙ ) = 11.Sharda et al. (2021a) and Gillman et al. (2022) reached similar conclusion at z∼0 and z∼1.5 respectively.An analysis of the Illustris TNG highresolution cosmological simulations also suggested that at z<2, galaxies with larger stellar mass generally have flatter (i.e. less negative) gas-phase metallicity gradients (Hemler et al. 2021).Simons et al. (2021) reported that massive galaxies have flatter gradients and lower mass galaxies have more positive gradients.This dependence may indicate the redistribution of metals in the early stages of galaxy formation, such as powerful feedback effects, radial flows, and stochastic accretion from minor mergers in the early universe (Maiolino & Mannucci 2019).
The dependence of gas-phase metallicity gradients on galaxy sizes is also disputed.Some works reported that larger galaxies tend to have steeper (more negative) metallicity gradients (Carton et al. 2018;Boardman et al. 2020), while some others reported no dependence on size of gas metallicity gradients for their samples (e.g.Sánchez-Menguiano et al. 2018;Simons et al. 2021).
Because the structure of galaxies depends on the stellar mass, it is natural to test the connection between morphology and metallicity gradients.Such a test has been carried out in the local universe.For example, Sánchez-Menguiano et al. (2016) reported that earlier type spirals show flatter gradients than the later type ones.How this dependence is affected by the redshift evolution of morphology, however, remains poorly understood.
Gas-phase metallicity gradients also provide crucial constraints on theoretical models of galaxy formation and evolution.For instance, the observed negative radial gradients in both the gas-phase and stellar metallicities at low redshift has been interpreted as strong evidence for the scenario of inside-out galaxy formation.In this model, the center of the galaxy forms first, undergoes more generations of star formation, and thus accumulates more metals compared with the outskirts.Although the negative trend is qualitatively easy to understand, quantitatively reproducing the observed gradients as well as their dependence on galaxy properties requires detailed modeling that includes assumptions on the dependence of star formation activity on the local gas distribution (e.g. the star formation law or Kennicutt-Schmidt Law, Kennicutt 1998) as well as the SFR radial dependence (e.g, Belfiore et al. 2019;Kang et al. 2021;Sharda et al. 2021b).
Recent theoretical efforts to understand and reproduce the observed metallicity gradients in galaxies have relied mostly on semi-analytic models (e.g., Forbes et al. 2019;Yates et al. 2021) and cosmological hydrodynamical simulations (e.g., Tissera et al. 2019;Hemler et al. 2021;Bellardini et al. 2021).Many works have succeeded in reproducing the present-day metallicity gradient in our Milky Way (e.g.Hayden et al. 2014;Minchev et al. 2018;Queiroz et al. 2021), but analyzing the metallicity in high-z galaxies with differing evolutionary histories is far more complicated.The results can be affected by which processes are being considered or emphasized, but a general conclusion of these works is that there is only a weak steepening trend in metallicity gradients as galaxies evolve towards the present day from z ≤ 4 (Curti et al. 2020).Sharda et al. (2021b) presented a theory of the evolution of gas-phase metallicity gradients in galaxies, derived from first principles, which is based on the unified galactic disc model given by Krumholz et al. 2018.This theory is able to explain the dependence of metallicity gradients on the stellar mass and kinematics of galaxies and reproduces the observations in the local Universe.
Recently, Simons et al. (2021) measured gas-phase metallicity gradients for a large sample of over 200 star forming galaxies at 0.6 < z < 2.6.The data were obtained from HST /WFC3 G102 slitless grism spectra as part of the Cycle 23 CANDELS Lyman α Emission at Reionization survey (CLEAR, Simons et al. 2023).The Figure 1.The histograms of the gas-phase metallicity gradients taken from Simons et al. 2021 for the whole sample as well as for three stellar mass bins of equal sample size (8.5<log M * /M⊙ <9.5 for the low mass bin, 9.5< log M * /M⊙ <9.9 for the medium mass bin, and 9.9< log M * /M⊙ <10.5 for the high mass bin).There is evidence that galaxies with higher stellar mass tend to grow slightly more negative metallicity gradients.
authors find a broad dispersion of values of the gradients, ranging from -0.11 to 0.22 dex kpc −1 with essentially no apparent trend or correlation (Pearson correlation coefficient |ρ| < 0.2) with galaxy properties such as stellar mass, SFR and morphology.
In this paper, we follow-up on the work by Simons et al. (2021) and expand the investigation to explore the dependence of the gradients with the star-formation history of the galaxies, which we have derived in nonparametric form using the software package P rospector (Johnson et al. 2021).The idea behind this investigation is that the presence or absence of such relationships will help provide a better physical understanding on the formation and assembly history of high-z galaxies.
The rest of the paper is organized as follows.Section 2 describes our data and sample selection.Section 3 shows the spectral energy distribution fitting with P rospector.Section 4, 5 and 6 investigate the relationship between the gas-phase metallicity gradients and multiple galaxy properties, including the SFH, colors and the morphology.Section 7 discusses the implications of this study in terms of galaxy evolution, and Section 8 provides a brief summary of our main conclusions.
Throughout this paper, we adopt a flat ΛCDM cosmology model with H 0 = 70 km s −1 M pc −1 and Ω m0 = 0.3.
In this paper we use the measures of gas-phase metallicity gradients of high-z star-forming galaxies in the GOODS-S and GOODS-N fields made by Simons et al. (2021) and expand the analysis of their correlation with the star formation history of the galaxies.A full description of the sample selection, observations, data reduction and analysis, as well as the creation of the line emission maps and measure of the metallicity gradients can be founds in Section 2 of Simons et al. (2021).The sample and the data of this work are exactly the same as those in that paper.To illustrate the data, in Appendix A we provide the emission line maps of one galaxy as an example, color coded by pixel SNR, together with extracted 1D spectra from selected extraction boxes.Simons et al. ( 2021) provide a complete description of the line fitting algorithms and of the azimuthal averaging procedures adopted to measure the metallicity gradients.
The data were obtained from deep near-infrared Hubble Space Telescope slitless spectroscopy as a part of the CLEAR survey (see Simons et al. 2023).The footprint of this survey covers 12 pointings in the CANDELS GOODS-N and GOODS-S fields (Grogin et al. 2011;Koekemoer et al. 2011), providing deep HST/WFC3 G102 grism spectroscopy and companion HST/WFC3 F105W direct imaging.All publicly available HST/WFC3 G102 and G141 grism observations that overlap with the CLEAR footprint are also included to extend the spectral coverage and maximize the depth.The redshift fitting and grism spectral analysis are conducted with the software Grizli (Brammer 2019), and the emission line fluxes are derived from Flexible Stellar Population Synthesis (FSPS) templates (Conroy et al. 2009).Emission line maps are then created by drizzling the contamination-and continuum-subtracted 2D spectral beams.
The sample of galaxies were selected based on high S/N line ratio criteria (≥ 5σ integrated detection in at least 2 of the strong lines in [OIII] λλ4958, 5007, [OII] λ3727, and Hβ) and for having a radial extent that allowed spatially resolved emission line maps.The metallicity gradients were derived from line ratio diagnostics of these maps.Galaxies with detected X-ray emission consistent with a Active Galactic Nuclei (AGN) have been removed from the sample.The total sample contains 238 star-forming galaxies at 0.6 < z < 2.6 across a mass range of 8.5 < log M * /M ⊙ < 10.5.Around 96% of the sample are consistent with a flat or positive metallicity gradient within 3σ confidence while only 4% show a negative gradient.The number distribution of the derived metallicity gradients is shown in Figure 1.We show the distri-bution both for the total sample, and for three equally divided sub-samples according to the stellar mass.

Photometry and morphology data
We have modeled the spectral energy distributions (SEDs) of the galaxies in our sample and reconstructed their SFH in non-parametric form with the package P rospector using the extended CANDELS (Cosmic Assembly Near-infrared Deep Extragalactic Legacy Survey, Grogin et al. 2011) photometric catalogs in the GOODS fields (Giavalisco et al. 2004).These include over 40 photometric bands from both ground-based and space-borne instruments, ranging from the restframe ultraviolet (UV) portion of the spectrum to the mid-infrared (MIR) one.For the GOODS-N field, we use the photometric catalog from Barro et al. (2019).This includes broadband data from the observed UV (U band from KPNO and LBC), optical (HST/ACS F435W, F606W, F775W, F814W, and F850LP), NIR (HST/WFC3 F105W, F125W, F140W, and F160W; Subaru/MOIRCS Ks; CFHT/Megacam K; NIR and MIR Spitzer/IRAC 3.6, 4.5, 5.8, and 8.0 µm), plus 25 consecutive medium bands from the SHARDS survey at λ = 500-950 nm (Barro et al. 2019).This yields a total of 42 photometric bands used in the P rospector SED modeling.
For the GOODS-S field we adopt the ASTRODEEP-GS43 photometric catalog made by Merlin et al. (2021), which is built and improved upon the previously released CANDELS catalog.To add to the CANDELS photometry in 7 bands (CTIO U, Hubble Space Telescope WFC3 and ISAAC-K), they also obtained measures in the other 36 bands (VIMOS, HST ACS, HAWK-I Ks, Spitzer IRAC, and 23 from Subaru SuprimeCAM and Magellan-Baade Fourstar), all the aperture-matched with the template-fitting package TPHOT (Merlin et al. 2021).This yields a total of 43 wavelength bands (25 wide and 18 medium filters).
Finally, Huertas-Company et al. (2015) presented a catalog of visual-like H-band morphologies of CAN-DELS galaxies estimated from Convolutional Neural Networks.We use their morphology classifications along with the best-fit H-band Sérsic indexes provided by the CANDELS catalog (van der Wel et al. 2012) to explore the dependence of metallicity gradients on the morphological properties of galaxies in Section 6.

SED FITTING WITH PROSPECTOR
With high-quality, multi-band photometry data, we use P rospector (Johnson et al. 2021) to perform SED modeling and reconstruct the SFHs of selected galaxies.P rospector is a Python-based package to conduct prin-cipled Bayesian inference of stellar-population properties by fitting photometric and/or spectroscopic data to spectral population synthesis models.The code enables one to obtain the galaxies' SFH in a non-parametric form, as opposed to a pre-assumed analytical functional form, which has been shown to minimize systematics in the measures of stellar mass, SFR and stellar age (Iyer et al. 2019;Leja et al. 2019aLeja et al. , 2021;;Tacchella et al. 2021;Dong et al. 2022).The fitting is based on the FSPS code (Conroy et al. 2009), which uses the MIST stellar isochrone library (Choi et al. 2016) and the MILES spectral library (Falcón-Barroso et al. 2011).Also, the code self-consistently models nebular emission (line plus continuum) using a photo-ionization model described in Byler et al. 2017.
Specifically, we adopt the non-parametric SFH model which consists with a certain number of bins in lookback time from the time of observations where the past SFR is estimated, using a prior to control the smoothness of the transition of the SFR from one bin to the next.Following other authors (e.g.Johnson et al. 2021;Tacchella et al. 2022a;Ji & Giavalisco 2022a), who have extensively tested the procedure, we also adopted the continuity prior described in Leja et al. 2019a in our analysis to emphasize smoothness in the star formation rate (SFR).The code directly fits for ∆ log(SF R) between adjacent time bins and explicitly weights against sharp transitions (see Leja et al. 2019a for details).We have also experimented with the Dirichelet prior, however, finding that the choice of the prior only minimally affects our results and conclusions.
We adopt 9 lookback time bins in this study, with the SFR being constant within each bin.The time bins are adjusted based on the age of the universe at the given redshift for each galaxy.The first lookback time bin is fixed at 0 < t < 30 Myr to capture the recent episodes of star formations, and the last bin is set to be 85%−100% of the age of universe at any given redshift.All the intervening bins are evenly spaced in log(t) scale.To test the possible dependence of the SED fitting results on the time binning procedure, we reran the fitting for all the GOODS-N galaxies in our sample with 7 look-back time bins and compared with the 9-bin results (see Appendix B for details).After comparing the recovered SFH plots and galaxy physical properties in the two cases, we concluded that the fitting results with 7 bins and 9 bins are fully consistent with each other.This analysis is also in agreement with previous work reporting that the fitting results are generally insensitive to the number of time bins when there are more than 5 bins (Leja et al. 2019a).Following other groups (e.g, Tacchella et al. 2022a;Ji & Giavalisco 2022a) for the rest of the paper we will be pre-senting our results obtained using the 9-bin SFH.The detailed parameter settings are described as follows.

Other adopted assumptions
In this work we adopt the Kroupa's initial mass function (Kroupa 2001) and the Calzetti's dust attenuation law (Calzetti et al. 2000).The diffuse dust V-band optical depth (the 'dust2' parameter in FSPS) is fitted with an uniform prior between 0 and 4. We also fix the redshift at the best-fit grism redshift given by Simons et al. (2021).The stellar metallicity (the 'logzsol' parameter in FSPS) is left as a free parameter, but with a strong Gaussian prior centered at the average observed gas-phase metallicity, i.e. obtained by averaging the observed gradient within 2R e , and width equal to the corresponding uncertainty.
The sampling is performed with the dynamical nested sampling code Dynesty (Speagle 2020), which allows us to control the effective time resolution and to sample preferentially near the given posterior.Though costly in terms of time, the Dynesty sampling method has been proved to be more accurate than the traditional MCMC sampling method, especially when a nonparametric SFH is assumed (Leja et al. 2019a).

Results
We do not directly test for the robustness of P rospector's SFH reconstruction here, because such tests have been done in detail by several different teams of investigators, including the builders of the code themselves (e.g, Leja et al. 2017;Johnson et al. 2021;Tacchella et al. 2022b;Ji & Giavalisco 2022a and reference therein).The general result from these testing procedures, which rely on synthetic galaxies from cosmological hydrodynamics simulations, such as IllustrisTNG (Pillepich et al. 2018), is that P rospector generally yields accurate and relatively precise SFH.
In this paper, we have tested the stability of P rospector against variations in the quality and number of the input photometric bands to take into account the fact that not all galaxies have the same photometric coverage and that SNRs in the various bands vary from galaxy to galaxy.First, we test for the effect of missing bands by rerunning the fitting with less photometric bands.In comparison with the ∼40 bands adopted before, we use a set of ∼20 photometric bands provided by the 3DHST survey (Momcheva et al. 2016), including data from the HHDFN U, B, V, R bands; HST/ACS F435W, F606W, F775W, and F850LP; HST/WFC3 F125W, F140W, and F160W; MODS J, H, Ks; and Spitzer/IRAC 3.6, 4.5, 5.8, and 8.0 µm bands.Although the results are highly consistent for most galaxies in our  sample, a few of them show different features in the bestfit SFH (See Appendix B).Therefore, to secure highquality SFH reconstructions we have visually inspected the photometry of our sample to only select galaxies with high-quality, densely sampled photometry.Specifically, we have removed galaxies that have more than 5 missing bands and/or have large photometric error bars, i.e.SNR≈ 5 or less, in ground-based observations.This vetting procedure leaves us with a sample of 119 starforming galaxies, about half of the original sample size.
After running P rospector on this final sample, we have recorded the best-fit 9-bin SFH, the stellar mass and star formation rate in the most recent bin (rSF R), i.e. the mass and SFR at the time of observation, as well as the mass-weighted age (t mw ) and the dust obscuration.Figure 2 shows the best-fit SED for one of the galaxies as a demonstration, while Figure 3 shows the individual plots of the SFH for each of the 119 sample galaxies.While all the galaxies are actively forming stars, there is a big variety in the individual shape of the SFH.To identify galaxies with star-formation rejuvenation, we visually divide the sample into two subsets based on the number of major episodes of star formation in their SFH.Among the total selected sample, we visually identify 39 galaxies with prominent multiple star formation episodes in their SFH, i.e. an overall nonmonotonic function, and 80 with only one main star formation episode, i.e. an overall monotonic function.Figure 4 demonstrates the 'typical' SFH curves for galaxies with single and multiple star formation episode.In the following analysis, we will test for possible correlations using both the total sample and each of these two subsets.

SFH VS. METALLICITY GRADIENTS
As shown in Figure 5, there is no obvious correlation between the gas-phase metallicity gradient and t mw or rSF R. To reduce high-frequency noise and identify more clearly possible trends, we have binned the points by dividing the X axis into 10 equally-sized bins and  calculated the median Y values in each bin.The binned data points as well as the best linear fit are shown in black, while the grey shaded region shows the 1σ fitting error.To quantify the significance of the possible correlation, we calculate the Spearman's rank correlation coefficient for each case.Except that t mw shows a very weak correlation with the metallicity gradients for galaxies with multiple SF episodes (Figure 5, left panel), the sense being that for these galaxies the older means more positive gradient, no correlation is found in any of the other cases (with |Spearman's rank correlation coefficients| < 0.2).In terms of rSF R, there is no difference between galaxies with single or multiple star formation episodes.

Characterizing the SFH
It is customary to use some characteristic redshift as a parameter to describe in a synthetic form the shape of a galaxy's SFH, for example the redshifts when a certain percent of the stellar mass has formed.Here we use an additional set of parameters, proposed by Ji & Giavalisco (2022b), designed to provide quantitative information on the overall shape of the SFH.The four parameters used here are listed below, where t obs refers to the age of the Universe at the observed redshift.
• τ tot : the total time length of star-formation (0 < t < t obs ) • τ 1 : the time length of the late period of star formation (t mw < t < t obs ), referred to as τ SF in Ji & Giavalisco 2022b • τ 2 : the time length of the early period of star formation (0 < t < t mw ), referred to as τ Q in Ji & Giavalisco 2022b • Skewness: the third standardized moment of stellar ages ( 8 We refer the readers to Section 2.4 of Ji & Giavalisco (2022b) for detailed descriptions of these parameters.The ratio between τ 1 and τ 2 , which strongly correlated with Skewness, serves as a good metric for quantifying the asymmetry of SFHs.We note that galaxies that fall at the extreme ends of the Skewness range usually have multiple SF episodes, i.e. a more complex SFHs than one characterized by an increasing phase, a peak more or less extended, and a decrease.All these parameters enable us to better identify the overall shape of SFHs as well as capture some feature in the star formation history.
Finally, following Giavalisco et al. (2023, in prep.), in an effort to further characterize the SFH in physical terms, we also use a new time parameter, defined to describe the evolutionary state of galaxies at the time of observation, named the evolutionary time, t E : where sSF R is the specific star formation rate, t obs the time at observation, and t i is the time when star formation started in the galaxy, which in practice can be assumed coincident with the Big Bang in this work.The quantity inside the integral is, essentially, the mass doubling time of the galaxy at the time t in unit of the Hubble Time (e.g.see Tacchella et al. 2019).Thus, t E , the integral of this quantity over the life time of the galaxy at the time of observation, is the average amount of time the galaxy spent doubling its stellar mass.Approximately, galaxies with the same t E are in the same evolutionary phase.Therefore, the relation between t E and the average of galaxy properties at t obs provides a description on the evolution of such properties with time.

Relationship between metallicity gradients and the SFH
We now investigate the relationship between the gasphase metallicity gradients and the five SFH parameters defined above.
As shown in Figure 6, with all the Spearman's rank correlation coefficients less than 0.2, we find no clear correlation between any of the four parameters that describe the overall shape of the SFH and the metallicity gradients.There is a hint that galaxies with single or multiple star formation episodes might occupy different regions in the diagrams of τ 1 and τ 2 versus the gradients.However, when we look at the Skewness (analogous to the ratio τ 1 /τ 2 ), the distribution of points shows no dependence on the number of SF episodes.Similarly, the relationship between τ tot and the metallicity gradients is also insensitive to the number of SF episodes.The results are summarized in Table 1, the machine-readable version will be available online.Only part of the table is shown here, the full data will be available online alongside publication.
To have additional insight into the possible effects of multiple episodes of star formation during the SFH of galaxies on the gas-phase metallicity gradients, in Figure 7 we plot the fraction of galaxies whose SFH shows evidence of distinct and repeated episodes of star formation versus the gradients.Such a fraction is defined as the ratio of the number of galaxies with multiple SF episodes to the total number of galaxies, and to do so we adopt the same selection criteria of 'multiple SF episodes' as described in Ji & Giavalisco (2022c).To test for the effect of uncertainties in such a selection, we also develop a more inclusive sample of multi-SF galaxies by adding in all 'possible' cases, whenever the error bars in the reconstructed SFH make it difficult to count the number of star-formation episodes (see Fig. 11 of Ji & Giavalisco 2022c).As shown in the figure, for both the original 'safe' selection and the 'more inclusive' one a highlysignificant correlation of the multi-SF fraction as the metallicity gradient increases is observed, with a Spearman's coefficient of -0.80 and a p-value of 0.10.Therefore, galaxies with more negative gas-phase metallicity gradients are more likely to have experienced multiple SF episodes in their SFH, and this effect seems robust with the details of the selection of SFH with multiple SF episodes.
We then investigate the relationship between the gasphase metallicity gradients and the evolutionary time t E , which is directly calculated from the SFH by replacing the integral with the sum over the 9 bins of look-back time.Figure 8 shows the effective radius R e and the metallicity gradients as a function of t E in Gyr.The size of these galaxies increases with the evolutionary time with a Spearman's coefficient of 0.58 and a pvalue of 0.03, inferring that star-forming galaxies grow larger along their evolution process.This finding is consistent with previous work about the size evolution of star-forming galaxies (e.g.Allen et al. 2017;Jiménez-Andrade et al. 2019).The metallicity gradient is also significantly correlated with t E , with a Spearman's coefficient of -0.60 and a p-value of 0.01.The gradient goes from positive values at small t E to ≈ 0 at large t E .The meaning of this correlation is that less evolved galaxies tend to have more positive metallicity gradient.In time, these gradients become progressively smaller and smaller, and possibly negative, as the galaxies evolve.
To gain further insight into the physical meaning of these correlations, in the left panel of Figure 9, we plot the gas-phase metallicity gradient versus the stellar mass and color-code the points for having one major SF episode (blue) or multiple episodes (orange) in their SFH.Although with large scatter, there is a significant correlation (ρ = −0.32 with a p-value < 0.01) such that smaller galaxies have positive gradients and more massive galaxies mildly negative ones, with galaxies with log M * /M ⊙ ≈ 10 having flat gradients.Given the big uncertainties in the measurement of ρ, we calculate the bootstrap distribution of ρ by re-sampling the data with corresponding error bars 100,000 times.As shown in the right panel of Figure 9, the 1 − σ bootstrap confidence interval is derived as [-0.57, -0.28], which confirms that the negative correlation between metallicity gradient and stellar mass is robust.Galaxies with multiple SF episodes are mostly concentrated toward the high-mass end of the distribution, while single-peak SFHs populate the low-mass side.The fact the galaxies with multiple SF peaks are on average more massive than those with one peak only has been reported by Ji & Giavalisco (2022c).The novelty here is that the former on average shows positive gas-phase metallicity gradients, while the latter shows flat or mildly negative gradients.Thus, the picture that seems to be emerging from these analyses is that, on average, galaxies develop positive gradients early on in their evolution.As evolution proceeds, these positive gradients are progressively reduced, becoming flat or mildly negative for highly evolved galaxies.This process seems to be more effective in galaxies with multiple bursts in their star formation history.Taking all together these elements suggests that gas accretion and radial transport of gas reduce the gradients by progressively enriching the inner regions and diluting the chemical content of the outer ones with fresh gas.

MEASUREMENTS OF COLORS AND ATTENUATION
The rest-frame colors as well as the radial color gradients are closely related to the metallicity of galaxies.Recent studies at moderate and high redshifts (z ≈ 2) have shown that stellar age, dust attenuation and metallicity can all affect the color gradients of star-forming galaxies (e.g, Liu et al. 2017;Miller et al. 2022).In an effort to disentangle the degeneracy among these effects, we first test the relationship between the metallicity gradients and the rest-frame UV and optical colors.The best-fit SEDs are convolved with the Johnson/Bessel transmission filters to get the flux densities of the rest-frame U, V, J bands.As shown in Figure 10, both the U − V and V − J colors show no correlation with the metallicity gradients.We continue to derive the color gradients for galaxies in our sample.We measured the azimuth-averaged radial surface brightness profiles in the PSF-matched HST/ACS V band (F606W) and z band (F850LP) images for each galaxy using the photutils package (Bradley et al. 2020), corresponding approximately to the U V −U color in the rest frame.We then performed a linear fit to characterize the color profiles between these two bands in semi-log scale (e.g., Wu et al. 2005;Liu et al. 2016).We adopted an inner radius limit of 0.3 kpc to reduce the PSF effect, and an outer radius cut at 2R e .Then we expressed the color gradient as d(F (V ) − F (z))/dlog(R) (0.3 kpc < R < 2R e ).
The derived color gradients of our galaxies are generally flat, with some scattering on both positive and negative sides (Figure 11).We calculated the Spearman's rank coefficient of the correlation between the color gradients and metallicity gradients, finding no significant correlation between these two quantities.This suggests that the gas-phase metallicity gradient does not drive color gradients in these galaxies.
We have also considered the possibility that differential dust attenuation across the profile of galaxies could also give rise to color gradients and in turn, provide systematic bias in the measures of metallicity gradients.As discussed by Simons et al. (2021), the measures are obtained from spectral indexes based on the ratio of multiple emission lines located inside a well defined range of wavelength.Using Monte Carlo simulations applied on the line maps, they showed that a spatially varying extinction, if any, could only change the metallicity gradient measurements by 5-30%, depending on the diagnostic and the intrinsic metallicity gradient.Further-    extinction gradients of the galaxies in our sample should also be, on average, flat.
To test the robustness of the measures of A V , we have used both the outputs from our P rospector SED fittings and the measures from the 3D-HST catalogs (Skelton et al. 2014), which have been obtained with a different procedure that uses the FAST package (Kriek et al. 2009), assuming discrete grid values of A V .The results are compared in Figure 12.As seen in the figure, the two results are qualitatively consistent with each other, possibly showing a tentative, very weak anti-correlation.The Spearman's coefficient for the Prospector measure is ρ ∼ −0.1 while that for the 3D-HST measures is ρ ∼ −0.3 with a p-value < 0.01, suggesting a more significant, although very week anti-correlation.In any case, both plots suggest that more dusty galaxies are more likely to develop steeper (more negative) metallicity gradients, in qualitative agreement with the general idea that more evolved galaxies develop more negative gradients.
Finally, we conducted two additional tests of the robustness of the results above.In the first test, with the results from the SED modeling in hand, we looked for any possible dependence of the gradients on the dust attenuation (A V ) on the distribution of metallicity gradients, which could be evidence that the gradients themselves are the results of the dust-metallicity degeneracy in the SED fitting procedure.We split the sample into two equal-size bins in A V measurements (low attenuation bin vs. high attenuation bin) and studied the distribution of the gradients in each bin.A K-S test shows that the distributions of metallicity gradients in the two bins are statistically equivalent.This is not surprising, given the overall rather low values of A V in our sample.In the second test, we studied the dependence of the gradients with the size of the galaxies, looking for the possibility that stronger gradients could be found in larger galaxies simply because of the larger dynamic range of radii.We similarly split the sample into two equal-size bins of R e and studied the distributions of gradients in both bins, finding no significant difference between them.In conclusion, the metallicity gradient measurements discussed here are not significantly biased by the gradients in the dust extinction or size measurements.

GALAXY MORPHOLOGY VS. METALLICITY GRADIENTS
Another galaxy property that may be related to metallicity gradients is the morphology.We characterize the morphology of our sample of 238 galaxies using a number of diagnostics.We adopt the Sérsic parameters n and R e provided by van der Wel et al. (2012).As shown in Figure 13, the relationship between gas-phase metallicity gradient and the Sérsic index n is flat.
In the top right panel of Figure 8 in Simons et al. (2021), they show the correlation between gas-phase metallicity gradient and the circularized effective radius normalized by the population average at the mass and redshift of the galaxy, reporting null correlations between these two measurements.As shown in Figure 14, we revisit this test and find consistent results.The 68% confidence interval of the bootstrap distribution gives ρ ∈ [−0.23, −0.18].Though very extended galaxies usually show more negative metallicity gradients, there is no correlation within R e ∼ 3 kpc.
Subsequently, using the H-band Sérsic index, stellar mass and effective radius, we calculate Σ 1 , the projected stellar mass surface densities within 1 kpc radius (Σ 1 ).As plotted in Figure 15, there is no strong correlation between the central mass densities and the metallicity gradients, which is consistent with the lack of trend probed by Sérsic indexes.The bootstrap distribution is narrowly centered with a 68% confidence interval of [-0.114, -0.109].
Given that Σ 1 strongly depends on the stellar mass, with more massive galaxies having larger central density, we also compare the metallicity gradients to the fraction of stellar mass within the central 1 kpc, M * , whose dependence on stellar mass is much weaker than that of Σ 1 (Ji et al. 2022).The relationship between the metallicity gradients and is plotted in Figure 16, which similarly shows no strong correlation.
With a significant correlation with stellar mass and tentative correlation with size, we detect null correlations with different measurements of concentration.This can be attributed to the large scatter in R e and Σ1 at a fixed stellar mass for star-forming galaxies.Promising candidates for the origin of this scatter are halo concentration and/or halo formation time (Chen et al. 2020).
We have also used the visual classification of H-band morphologies given by Huertas-Company et al. 2015, which provides the probability of being spheroid/disklike/irregular for each galaxy.Here we combine the disk and irregular probabilities into one single category for easier and more reliable comparison with spheroid probabilities.From Figure 17, the dependence of metallicity gradients on the visual morphology is rather weak.
The plane defined by the Gini coefficient G and the second moment M 20 , (Gini, M 20 ), provides a widelyused non-parametric classification of galaxy morphology, both regarding morphology as well as indication on the probability that they are mergers observed during various stages of the merging event (see, e.g., Lotz et al. Figure 12.The relation between the gas-phase metallicity gradients and the V band attenuation AV .Marks and labels are the same as Figure 8.In the left panel, AV is derived from SED fittings using P rospector; in the right panel, AV is directly taken from the 3D-HST data release (only reported to one significant digit).
Figure 13.The binned scatter plot between the gas-phase metallicity gradients and the Sérsic indexes.Marks and labels are the same as Figure 8. 2004,2008).Figure 18 shows the distribution of our galaxies on the (Gini, M 20 ) plane.The lines mark the boundaries between mergers, disk-like galaxies and ellipticals/S0s, based on nearby galaxies and adjusted to z∼1 relative to the calibration of the coefficient made at z ≈ 0 (Lotz et al. 2004(Lotz et al. , 2008)).It has been found, however, that there is no significant change in this diagnostic at higher redshifts (e.g.Kartaltepe et al. 2022).As seen in the plot, most of the galaxies in our sample have disk-like morphologies and relatively few of them are located in the region of the plane where candidate mergers are found.The figure also shows that there is a large dispersion of metallicity gradients for both mergers and non-mergers (disk galaxies), which indicates that merging activities are not the dominant mechanism for the redistribution of gas and metals and the flattening of the metallicity gradients.2021), the gas-phase metallicity gradients in star-forming galaxies of approximately the Milky Way mass at 0.6 < z < 2.6 are essentially flat, with a broad scatter observed across the population (see Figure 1).This scatter accounts for galaxies with both negative and positive gradients, with the former being a minority of cases.The flat and mildly positive metallicity gradients of star-forming galaxies at the cosmic noon are in stark contrast with the measurements in the local universe, which in general show negative gradients (e.g, Kreckel et al. 2019;Sánchez 2020).
Overall, we find that the gas-phase metallicity gradients show rather weak correlations with the properties of galaxies, which is in agreement with the analysis by Simons et al. 2021.A summary of all the Spearman's rank correlation coefficients and p-values is given in Table 2. From new measures of the stellar mass (total mass formed) that we have obtained from SED modeling with P rospector, we find a noisy but significant anti-correlation between the metallicity gradient and the stellar mass.Namely, more massive star-forming galaxies have less positive (flatter and negative) metallicity gradients, although the scatter is large at lower stellar mass.The dependence of metallicity gradients on the stellar mass has been widely discussed (e.g., Barrera-Ballesteros et al. 2016;Hwang et al. 2019).Low-mass galaxies and active star-forming galaxies tend to have flat gradients, likely due to more efficient stellar feedback (Ma et al. 2017).Barrera-Ballesteros et al. 2016 argue that the resolved mass-metallicity relation can reproduce the gas-phase metallicity gradients for all galaxies except the least massive ones.The stellar-mass range of galaxies in our sample is relatively narrow (8.5 < log M * /M ⊙ < 10.5) and falls in the low-mass end of the present-day mass function, which may explain the   general flatness and lack of correlations in the gradients, given the relatively small dynamic range.
Given that star-forming galaxies increase their stellar mass as they evolve, another way to interpret this anti-correlation is that as the stellar mass of galaxies grows with time, the gas-phase metallicity gradients decrease.Additional evidence in support of this trend is that galaxies with more negative gas-phase metallicity gradients are more likely to have experienced multiple star forming episodes in their SFH, which also requires longer evolutionary time.
Further independent evidence comes from the weak anti-correlation between the metallicity gradient and the size of galaxy as quantified by R e , which increases with time in star-forming galaxies as they evolve (e.g.Allen et al. 2017;Jiménez-Andrade et al. 2019).Note that the correlation with R e is most likely diluted by the large scatter in R e at fixed stellar mass for star-forming galaxies at high redshifts (Chen et al. 2020).The growth of R e can be directly observed in the left panel of Figure 8, where R e is plotted against the evolutionary time t E , providing a metric of the amount of the evolution a given galaxy has undergone through at the epoch of observation (z obs ) of observation (Giavalisco et al. in prep; see also the 'mass doubling number' in Tacchella et al. 2022a).
Lastly, additional evidence that the evolutionary trends of metallicity gradients comes from the anticorrelation between the metallicity gradient and the evolutionary time t E itself, which is shown in the right panel of Figure 8.This anti-correlation directly reflects the 'flattening' or the decrease of the gas-phase metallicity gradients with time.In conclusion, the implications from all these anti-correlations are in good agreement with each other, indicating that star-forming galaxies develop flatter or more negative metallicity gradients as they evolve.
Except for the correlation illustrated above, we find no additional statistically-significant correlation between the gas-phase metallicity gradients and other stellar population properties of the galaxies, such as the recent SFR, mass-weighted age, dust attenuation, rest-frame UV and optical colors.This finding aligns with earlier results reported in Simons et al. (2021), where no strong correlations are found between the metallicity gradients and SFR, sizes, star-formation surface density, and starformation per potential energy.Clearly, the generally negative gradients observed in the local universe must be the result of evolutionary process that took place after the epochs of observation reported here, namely 0.6 < z < 2.6.The lack of correlation with rest-frame colors and color gradients described here also suggests that the stellar age or stellar activities in the broad picture can not fully explain the observed metallicity gradients.Dust attenuation plays a marginally important role to shape the metallicity gradients, where more dusty galaxies grow or retain more negative metallicity gradients, but the correlation is quite noisy, especially at the low dust attenuation end.
As mentioned in Section 1, due to the local enrichment from star-formation, a negative radial gradient in the gas-phase metallicity would be a natural consequence of inside-out galaxy formation scenarios.However, gas exchanges continuously take place in galaxies on rapid timescales (e.g.Somerville et al. 2015), including outflows of enriched gas and inflows of gas at low metallicity accreted from outside the virial radius.These gaseous exchanges evidently contribute to shaping the evolution of the metallicity gradients.The magnitude of such effects, of course, depends on the extent of radial gaseous flows.The phenomenology presented in this paper suggest that star-forming galaxies evolve from an early evolutionary phase characterized by flat or mildly positive gradients toward more evolved phases during which the metallicity gradients are essentially flat or negative.This dependence is not necessarily in contrast with the inside-out galaxy formation scenario, it just evidences the effects of large-scale redistribution of metals in the early stages of galaxy evolution.
The broad dispersion of metallicity gradients at high redshifts, in contrast with observations in the local universe, adds to the complexity of understanding how the metallicity gradients are affected by the assembly history of galaxies.Given the generally higher and sustained star formation rates, the big diversity of gasphase metallicity gradients observed at high-z would be naturally explained by the negative gradients being rapidly destroyed before fully established, so that the well defined negative gradients are rare at among the general population of star-forming galaxies.
Conversely, the flat and mildly positive metallicity gradients may trace accretion or radial inflow of nearpristine gas.Accretion of metal-depleted gas and efficient transport of this gas toward the center of the galaxies, which is intimately related to metal transports, can also efficiently flatten or even reverse gas-phase metallicity gradients, mainly through winds, accretion, stellar Figure 17.The binned scatter plot between the gas-phase metallicity gradients and the morphologically classification probabilities.The left panel refers to the probability of being disk-like or irregular galaxies, the right panel refers to the probability of being spheroid galaxies.Marks and labels are the same as Figure 8.
Figure 18.The H band Gini value as a function of M20 color-coded by the gas-phase metallicity gradients.The black lines show the boundaries between mergers (dashed), disks and elliptical galaxies (dash-dotted) given by Lotz et al. 2008.feedback, and merging events (Ma et al. 2017;Grand et al. 2019;Grossi et al. 2020).Note that our analysis is based on a sample of active star-forming galaxies which mostly include disks.Recently, Bílek et al. (2022) report that early-type galaxies with lower rotational support tend to have more negative metallicity gradients.
As the galaxies evolve, however, our analysis suggests that the rate of change of metallicity due to starformation is higher in the centers of galaxies than in their outskirts (see the toy model in Simons et al. 2021), and thus negative metallicity gradients would naturally develop.Losing metal-rich gas preferentially from the outskirt, where the escape velocity is lower, and accreting high-angular-momentum gas would also help decrease the gas-phase metallicity gradients.Genzel et al. (2023) find significant evidence of non-tangential flows in star-forming galaxies at z∼2, where the gravitational torques drive gas rapidly inward and result in the formation of central disks and large bulges.Similarly, Weldon et al. ( 2023) report the discovery of cool gas inflows towards three star-forming galaxies at z∼2.3, their analysis suggests that the inflows are due to recycling metalenriched gas from previous ejections.
We now attempt to develop some physical explanations for the lack of correlations with some of stellar population properties discussed above.First of all, the measured metallicity gradients in our sample don't show big deviation from flat gradients.The overall range is [-0.10,0.16] dex kpc −1 , with more than 60% galaxies have flat gradients (within [-0.05, 0.05] dex kpc −1 ).Considering the uncertainties of measurements, the small diversity of gradient values disfavors strong correlations with any galaxy properties.Some works claim that the total stellar mass plus the local stellar mass surface density are sufficient to qualitatively reproduce the observed trends of metallicity gradients (e.g.Boardman et al. 2021).However, whether we can use similar relations to understand the gas-phase metallicity gradients in high-z galaxies is questionable.At early times, it is suggested by the flat metallicity gradients that metals are more frequently and efficiently redistributed due to the prevalence of gas transport mechanisms such as gas accretion, outflows and mergers.In parallel, the elevated levels of gas turbulence detected in the interstellar gas of galaxies at z > 0.5 (Weiner et al. 2006;Förster Schreiber et al. 2006, 2009;Kassin et al. 2007Kassin et al. , 2012;;Wisnioski et al. 2015;Simons et al. 2016Simons et al. , 2017;;Übler et al. 2019;Price et al. 2020) will promote elevated mixing of the gas-phase metals in galaxies.
The kinematic properties of high-z galaxies can vary on a relatively fast time scale together with their starformation activities and feedback, and then significantly change the gas-phase metallicity gradient.Note that Genzel et al. (2023) reports evidence of radial gaseous flows based on the deviations of the gas kinematics from differential rotation.
Thus, the general lack of strong correlations between metallicity gradients and a number of properties related to the star-formation activity presented here suggests that the observed metallicity gradients could reflect more the instantaneous star formation distribution and rapid gaseous flows rather than the long-term growth history of high-z galaxies.This conclusion also finds support by numerical simulations (e.g., Ma et al. 2017).

SUMMARY
In this paper, we use the Bayesian SED fitting code P rospector to reconstruct the non-parametric SFHs of 238 star-forming galaxies, and explore possible relations between a number of galaxy properties and the gasphase metallicity gradients provided by the CLEAR survey (Simons et al. 2021).Our main conclusions are summarized below: • In our sample, there is a weak but significant anticorrelation between the metallicity gradient and the stellar mass (see Figure 9), suggesting that the gas-phase metallicity gradients decrease as the stellar mass of galaxies grows with time.
• There is an anti-correlation between the metallicity gradient and the evolutionary time (Figure 8), suggesting that less evolved galaxies tend to have more positive metallicity gradients, these gradients become progressively smaller and possibly negative as the galaxies evolve.
• Galaxies with smaller metallicity gradients are more likely to have experienced multiple star forming episodes in their SFH (Figure 7).
• We find no statistically-significant correlation between the gas-phase metallicity gradients and other stellar population properties of the galaxies (|Spearman's rank correlation coefficients| < 0.2), including the mass-weighted age, recent star formation rate, the overall SFH shape (e.g., total time width, asymmetry), dust attenuation, restframe UV/optical colors, and the morphology (see analysis in Section 4, 5, and 6).The lack of correlations could be partly due to the narrow stellar mass range, which also extends toward the lowmass end.
While the metallicity gradients of nearby galaxies can be understood by certain relations, the situation is much more complex in the high-z universe.At early time, gas kinematics can frequently change the distribution of metals in and around galaxies with a relatively short timescale, which means that the observed metallicity gradients may only reflect the instantaneous state of galaxies at the time of observation (see e.g., Acharyya et al. in prep).
These observed correlations or, better, lack of correlations, can help constrain theoretical models of galaxy evolution and be compared with numerical simulations.The GOODS-N and GOODS-S fields have been probed with several JWST surveys, including the JWST Advanced Deep Extragalactic Survey (JADES; Eisenstein et al. ( 2023)), the JWST Extragalactic Medium-band Survey (JEMS; Williams et al. (2023)), and the First Reionization Epoch Spectroscopic COmplete Survey (FRESCO; Oesch et al. (2023)).In the near future, by incorporating the JWST photometry, we will be able to revisit our analysis with larger data sample and more reliable SED modeling.Given high-quality spatial resolved metallicity measurements at high redshifts, robust comparison between the data and models can be realised.

ACKNOWLEDGEMENT
We thank our collaborators from the CLEAR team for many valuable discussions and contributions.This work has made use of the Rainbow Cosmological Surveys Database, which is operated by the Universidad Complutense de Madrid (UCM), partnered with the University of California Observatories at Santa Cruz (UCO/Lick,UCSC).
The parallel SED fittings are performed by the Unity cluster, which is a collaborative, multi-institutional high-performance computing cluster located at the Massachusetts Green High Performance Computing Center (MGHPCC).

A. SPECTRAL LINE DETECTION AND FITTING
The data used in this paper are presented and described in the paper by Simons et al. (2021), including data reduction, analysis, extraction of the spectra and measures of the metallicity gradients.Here we provide some examples of the emission line maps and extracted spectra to illustrate the data.In Figure 19, we show 2D maps of the optical lines of one galaxy from our sample, color coded by pixel SNR, as well as the 1D extracted spectra of the lines in selected extraction boxes.Simons et al. (2021) also detail the fitting procedures to measure the line properties, as well as the azimuthal averaging algorithm to derive the metallicity gradients.

B. ROBUSTNESS TEST OF PROSPECTOR FITTINGS
Different choices of the input photometry bands and the sampling of the lookback time could affect the SED fitting results, and in turn the reconstructed SFHs.To test the stability and robustness of our P rospector fittings, we compare the output SFHs with different numbers of input photometry bands and lookback time bins.
As described in Section 3, we adopt the ∼20 photometry bands provided by the 3DHST survey to test the SFH dependence on the input photometry bands.The results are highly consistent for most galaxies in our sample, while a few of them show different features in the best-fit SFH.Particularly, comparing with the original choice of ∼40 photometry bands, there is less variation in the reconstructed SFH among different galaxies when using ∼20 bands.As shown in Figure 20, there is no systematic trend in the SFH difference between the two.We conclude that the optimization could be slightly undermined by local minima when fewer photometry bands are available.Therefore, we decide on using all available photometry bands and filter those sources with scarce data points.For the test of time bins, aside from the original setting of 9 bins, we adopt the 7 time bins defined by Leja et al. 2019a and adjust them based on the age of the universe at the known redshift.To compare the 7-bin SFH of the test runs with the 9-bin SFH from our primary runs, we perform a linear interpolation of the 9-bin SFHs to make equivalent 7-bin results.The differences in the reconstructed SFHs are shown in the right panel of Figure 20.The results are highly consistent for different choices of time resolution and there is no systematic trend in the SFH differences.
In summary, the test runs indicate that the reconstruction of SFHs from P rospector is stable against the input photometry and the time resolution, supporting the robustness of our results.
Furthermore, to study the differences between SED fitting codes, we compare our fitting results with those reported by Simons et al. (2021).Their SED modeling is conducted by EAZY (Brammer et al. 2008), which is optimized for producing high-quality photometric redshifts over 0 < z < 4. The comparison of output stellar masses, star formation rates, as well as the dust attenuation A V is shown in Figure 21.P rospector tends to give higher stellar mass than EAZY at the high mass end, and provides relatively lower SFR, while the dust attenuation is uniformly scattered.The differences between multiple SED fitting codes have long been identified and discussed in literature (e.g., Leja et al. 2019b;Pacifici et al. 2023).Here we find no significative differences bias in the results of SED modeling from the two codes.

Figure 2 .
Figure 2.An example of the best-fit SED for source No.16607 in the GOODSN field.The wavelengths are in units of Å, the flux densities are in units of maggies (Jy/3631).The red circles show the observed data, and the gray lines show the transmission curves of the filters.A random model drawn from the sampling chain is plotted in blue, while the model at the maximum posterior probability (MAP) is plotted in green.

Figure 3 .
Figure 3.The stacked star formation history plots with 9 lookback time bins

Figure 4 .
Figure 4. Example plots for galaxies with single (left) or multiple (right) star formation episodes in their SFHs.

Figure 5 .
Figure5.The gas-phase metallicity gradients as a function of mass-weighted age (left) and the recent star formation rate (right).Blue represents galaxies with only one main SF episode in their SFH, orange represents galaxies with multiple SF episodes.The binned scatter plot with 1σ error-bars are plotted in black.The black solid line shows the best linear fit, and the grey shaded region shows the 1σ fitting error.The text label shows the Spearman's rank correlation coefficient for the total sample (black, ρtot), the single-SF sample (blue, ρ single ) and the multi-SF sample (orange, ρ multi ).No strong correlation is found in both cases.

Figure 6 .
Figure6.The relationship between the gas-phase metallicity gradients and four SFH characteristics.Marks and labels are the same as Figure5.According to the low Spearman's coefficients, no strong correlations are identified here.

Figure 7 .
Figure 7.The fraction of galaxies showing multiple starforming episodes in their reconstructed SFH at different gasphase metallicity gradients.The blue points show the safe selection of multi-SF galaxies, while the orange points show the more inclusive selection.An obvious anti-correlation is seen in each of the case.

Figure 8 .
Figure 8.The effective radius Re (left) and the gas-phase metallicity gradients (right) as a function of the evolution time tE.The binned scatter points are shown in black squares with 1σ error-bars.The black solid line shows the best linear fit, and the grey shaded region shows the 1σ fitting error.The unbinned data points are plotted in light blue.

Figure 9 .
Figure 9.The gas-phase metallicity gradients as a function of the stellar mass.The best linear fit as well as the 1σ error range of the binned scatter plot are shown in black.The number distribution of the correlation coefficient after bootstrapping is shown in the right panel.The red and black dashed lines show the center and the 68% confidence interval of the distribution.

Figure 10 .
Figure 10.The relation between the gas-phase metallicity gradients and the rest-frame U − V and V − J colors.

Figure 11 .
Figure 11.The relation between the gas-phase metallicity gradients and the radial color gradients.Marks and labels are the same as Figure 8.

Figure 14 .
Figure 14.The binned scatter plot between the gas-phase metallicity gradients and the effective radius Re.The number distribution of the correlation coefficient after bootstrapping is shown in the right panel.The red and black dashed lines show the center and the 68% confidence interval of the distribution.

Figure 15 .
Figure 15.The binned scatter plot between the gas-phase metallicity gradients and the central stellar mass surface densities.The number distribution of the correlation coefficient after bootstrapping is shown in the right panel.Marks and labels are the same as Figure 14.

Figure 16 .
Figure 16.The binned scatter plot between the gas-phase metallicity gradients and the normalized central stellar mass ( M 1kpc M * ).Marks and labels are the same as Figure 8.

Software:Figure 19 .
Figure 19.Examples of line maps and spatially-resolved 1D extracted spectra from the CLEAR observations.The data have been presented and discussed by Simons et al. (2021), who also provide a complete description of the data reduction, analysis and spectral extraction procedures.In this paper we have used exactly the same data set.The lower panel shows the 2D spectra and S/N maps for one of the sample galaxies.The magenta and orange boxes mark examples of the extraction regions for the center and outskirt respectively.The gradients have been obtained after azimuthal averaging.The extracted 1D spectra are plotted in the upper panel, where the spectral features are clearly identified.

Figure 20 .
Figure 20.The differences between the reconstructed SFHs in terms of different photometry band choices (left) and lookback time bin choices (right).The orange scattered points represent individual galaxies, the black squares indicate the median Y axis values for each time bin, and the shaded areas show the 68th and 95th percentiles.In the right panel, to directly compare the 9-bin results with the 7-bin results, the SFR values have been interpolated from adjacent bins.

Figure 21 .
Figure 21.Comparison of SED fitting outputs between EAZY (Simons et al. (2021)) and P rospector (this work).The dashed line marks the one-to-one relation.

Table 1 .
Simons et al. (2021)s and SFH characteristics derived from SED fittings Note-The IDs are matched to the CANDLES photometric catalog.The metallicity gradient measurements are taken fromSimons et al. (2021).

Table 2 .
The Spearman's rank correlation coefficients and p-values between the gas-phase metallicity gradients and different galaxy physical properties.