Chemical Abundances of Young Massive Clusters in NGC 1313^∗

The ISPy3 suite of Python routines relies on the information provided by the Hertzsprung–Russell diagram (HRD) of the star clusters under analysis. Three main approaches have been adopted for deriving the physical parameters from the HRD: (a) color–magnitude diagrams (CMD; Larsen et al. 2017), (b) CMD+isochrones (Larsen et al. 2012; Hernandez et al. 2017; Larsen et al. 2018), and (c) isochrones (Hernandez et al. 2018a, 2018b, 2019). Method (a) is optimally applied for nearby or local stellar clusters where resolved CMDs are available covering all evolutionary stages present. Method (b) is commonly used for extragalactic stellar clusters where the resolved CMDs typically include only the brightest stars in the cluster. The HRD is then complemented with theoretical isochrones covering fainter magnitudes. In the analysis presented here, we generate HRDs using solar-scaled stellar isochrones from PADOVA v.1.2.S (Bressan et al. 2012), adopting the isochrone-only approach, as in (c). The main driver for adopting this approach is the absence of CMDs for the YMCs in our sample. We note, however, that Hernandez et al. (2017) found that the differences between methods (b) and (c) have only a minor effect on the inferred abundances of YMC NGC 1313–379, of the order of ≲0.1 dex.

For each of the YMCs, we select the initial isochrones assuming a Large Magellanic Cloud (LMC)-like metallicity (Walsh & Roy 1997), [Z] = −0.3 dex, and the ages listed in Table 2 published in the Legacy Extragalactic UV Survey (LEGUS) cluster catalog (Adamo et al. 2017). From the corresponding isochrones we extract the stellar parameters, e.g., effective temperatures (T_eff), stellar masses (M), and surface gravities (log g), following a Salpeter (1955) initial mass function (IMF) and a lower mass limit of 0.4 M_⊙. We note that Hernandez et al. (2019) found that changing the choice of IMF (e.g., from Salpeter to Kroupa IMF) modifies the inferred metallicities on average by ≪0.1 dex. Lastly, similar to previous IL studies of YMCs we adopt microturbulent velocity values (v_t ) depending on the effective temperature of the stars: v_t = 2 km s⁻¹ for stars with T_eff < 6000 K, v_t = 4 km s⁻¹ for stars with 6000 < T_eff <22,000 K, and v_t = 8 km s⁻¹ for stars with T_eff > 22,000 K.

Using the parameters described above, we generate atmospheric models for each stellar type. Similar to the work by Hernandez et al. (2018a, 2018b), we make use of the local thermodynamic equilibrium (LTE) 1D plane-parallel atmosphere modeling code ATLAS9 (Kurucz 1970) for any stars with T_eff > 5000 K. For stars with T_eff < 5000 K we then use MARCS atmospheric models instead, as these have been primarily developed with a special focus on late-type stars (Gustafsson et al. 2008). These MARCS models are LTE 1D plane-parallel or spherical. One important distinction between the two modeling codes is that we generate the ATLAS9 atmospheres on the fly, whereas the MARCS models have been precomputed and downloaded from their official website. ⁵ At this point in the analysis we set the initial abundances of the stellar population. We note that in the work presented here we adopt the solar composition from Grevesse & Sauval (1998).

We also highlight that the entirety of our analysis is based on LTE modeling. We currently do not apply any non-LTE (NLTE) corrections as these are dependent on the individual stellar types, and our work focuses on the integrated light of the star clusters. Previous studies have looked into the effects of NLTE treatment in the analysis of IL spectroscopic observations of GCs and found that the overall metallicities can change by ∼0.05 dex (Young & Short 2019). A different study by Conroy et al. (2018), also focused on older stellar populations (>1 Gyr), and concluded that the assumption of LTE is adequate when studying Na and Ca. Lastly, Eitner et al. (2019) find that in many cases these NLTE corrections, specifically for Ba, Mg, and Mn, are significant (>0.1 dex) for metal-poor GCs. Most of these IL NLTE studies have so far focused on old populations of stars, where our focus here is on YMCs.

As it has been extensively described in our previous studies, our analysis requires the creation of individual synthetic spectra for each stellar type in the population under study. We continue to use the line lists by Castelli & Hubrig (2004). For those atmospheres generated with ALTAS9 we use the SYNTHE software (Kurucz & Furenlid 1979; Kurucz & Avrett 1981), and TURBOSPECTRUM (Plez 2012) for the MARCS atmospheric models. We then coadd the individual spectra to create a single IL model spectrum.

The IL model produced by the spectral synthesis codes is created at resolutions of R ∼ 500,000. To match the resolution of the X-shooter data, R ∼ 5000–9000, we fit for the best Gaussian dispersion, σ_sm, and convolve it with our IL synthetic spectrum. This smoothing parameter, σ_sm, accounts for the instrumental resolution, σ_inst, and the line-of-sight velocity dispersion of the star cluster, σ_1D.

Once the resolution of the synthetic spectrum and that of the observations are similar, the ISPy3 tool matches the continuum of the X-shooter spectrum with that of the model using a cubic spline or polynomial, depending on the size of the wavelength window being analyzed. Our software allows us to define specific weights for a given pixel in the X-shooter spectrum, anything between 0.0 and 1.0. To avoid contamination from emission lines or other problematic features in the IL observations we have set the weights of these wavelength regions to 0.0, to exclude them from our analysis.

In previous publications we have assessed the various uncertainties and systematics introduced by the different choices made in the modeling. As mentioned in Section 3.1, Hernandez et al. (2017) showed that the differences in the inferred abundances applying method (b) and (c) are very minor, of the order of ≲0.1 dex. Additionally, studying the integrated light of a large sample of stellar clusters Hernandez et al. (2018a) found overall differences <0.1 dex when comparing abundances inferred using solar-scaled (PARSEC) and α-enhanced isochrones (Dartmouth, Dotter et al. 2007). Regarding the uncertainties introduced by the different IMF selection, Hernandez et al. (2019) calculate abundances adopting a Kroupa IMF (Kroupa 2001) and a Salpeter IMF (Salpeter 1955) and found that the uncertainties introduced by this selection is on average of the order of ≪ 0.1 dex. Overall, it has been found that the systematic uncertainties should be of the order of those quoted in Table 5 or smaller.

We note that three of the YMCs analyzed here displayed strong emission lines, particularly in Hα. The three different clusters, NGC 1313-463, NGC 1313-503, and NGC 1313-505, exhibit slightly different emission profiles (Figure 3). The Hα profile observed in YMC NGC 1313-503 is comparable to those identified in Hernandez et al. (2017) for YMCs NGC 1313–379 and NGC 1705-1. The nature of the broad profiles, similar to those observed in these three clusters, is briefly discussed by Melnick et al. (1985), where they proposed that one possibility for the origin of this emission in Hα, other than ionized gas, is the presence of Be stars. These objects are identified as B-type nonsupergiant stars with rotational velocities of several hundred km s⁻¹, which display strong and broad Hα emission features (Marlborough 1982; Townsend et al. 2004). Additionally, studies have found an enhancement of Be stars in stellar clusters with ages <100 Myr (e.g., Mathew et al. 2008). The broad Hα emission observed in NGC 1313-503, reaching velocities of ∼300 km s⁻¹, along with its estimated age of 20 Myr, hints at the presence of Be stars.

Zoom In Zoom Out Reset image size

The Hα profile observed in NGC 1313-463, a cluster with an estimated age of 40 Myr, appears to be slightly more complex. From visual inspection we believe the emission could be originating from both a population of Be stars (creating a slightly broad component, reaching velocities of ≲200 km s⁻¹), and from ionized gas surrounding the cluster (much narrower component). Lastly, the emission features observed in NGC 1313-505 appear to be originating primarily from gas, although a population of Be stars cannot be fully discarded. We highlight that any wavelength regions contaminated by emission features, irrespective of their origin, are masked in our abundance analysis.

5.1.1.  α-elements

Known as α elements, Ca, Ti, and Mg, are primarily produced by high-mass stars and ejected to the interstellar medium through core-collapse supernovae (Woosley & Weaver 1995). Due to their nature, α-elements are typically used to gain insight into the period of time when core-collapse supernovae shaped the chemical evolution of the host galaxy. Since Fe-peak elements are also produced in massive stars, the result is a constant [α/Fe] ratio during those epochs dominated by core-collapse supernovae. Typically, the [α/Fe] ratios begin to decline when Type Ia supernovae (SNIa) start to drive the chemical enrichment, as these release large amounts of Fe-peak elements (Timmes et al. 2003). Overall, these are the general trends we observed in the MW. Additionally, Galactic abundances of α elements appear to correlate with each other. Such a trend can be clearly observed in Figure 4, where we show in blue stars the abundances of field stars in the MW as measured by Reddy et al. (2003, 2006). We also show in Figure 4 as orange stars the observed abundances in the LMC, a high-mass dwarf galaxy, published by Van der Swaelmen et al. (2013). We note that the LMC and NGC 1313 share a few properties. Both galaxies are classified as spiral galaxies, the LMC as SBm and NGC 1313 as SB(s) (Sandage & Tammann 1981), with comparable metallicities (Russell & Dopita 1992; Walsh & Roy 1997).

Zoom In Zoom Out Reset image size

The α-element abundances presented here, are the first-ever measured in stellar populations in NGC 1313. The abundances for all the clusters in our sample, including those from Hernandez et al. (2017), are plotted as orange squares in Figure 4 as a function of [Fe/H]. The [Ca/Fe] ratios for all of the clusters in NGC 1313 tend to decrease with increasing [Fe/H], following a trend comparable to that observed in the LMC, and with slightly more depleted Ca abundances than those from the MW field stars.

The [Ti/Fe] values observed in the YMC sample show a much higher scatter than that of Ca; however, a similarly high scatter is seen in the field stars of the LMC. This larger scatter in the LMC [Ti/Fe] ratios is primarily introduced by the stars in the inner disk, compared to the abundances in the LMC bar (Van der Swaelmen et al. 2013). It is also clear from Figure 4 that overall, the Ti abundances are lower in NGC 1313, compared to those in the LMC for a given metallicity.

We observe a comparable scatter in the [Mg/Fe] ratios in NGC 1313, to that in the [Ti/Fe] abundances. However, we also highlight that the uncertainties in the Mg abundances are much higher than those from Ti, particularly for NGC 1313-505. For this cluster we infer [Mg/Fe] = +0.70 ± 0.7 dex, with abundances compatible to those in the field in the MW and LMC, within the uncertainties.

5.1.2. Fe-peak Elements

In Figure 5 we show the NGC 1313 YMC abundances for the Fe-peak elements Cr and Ni, and compared them with those observed in stars in the MW (Reddy et al. 2003, 2006) and the LMC (Van der Swaelmen et al. 2013). Both of these Fe-peak elements are primarily produced by Type Ia SNe. However, in spite of sharing a common origin, they exhibit slightly different patterns.

Zoom In Zoom Out Reset image size

Cr abundances, both in the LMC and the MW, display a much higher scatter than that observed in Ni. The stellar abundances in these two environments are typically flat around [Cr/Fe] = 0.0 dex. The [Cr/Fe] ratios measured for our sample of YMCs in NGC 1313 appear to follow a comparable trend to those observed in the MW and LMC. The only cluster that appears to slightly deviate from this trend is the most metal-poor YMC NGC 1313–379, studied in Hernandez et al. (2017), displaying a subtle Cr enhancement.

Ni also exhibits a flat distribution in both the MW and the LMC. The trend observed in the LMC, however, is consistently subsolar. Overall, the Ni abundances measured in NGC 1313 follow the flat trends of the MW and LMC, again with the exception of the most metal-poor YMCs, as well as the most metal-rich target.

Over the last couple of decades, studies investigating the chemical variations as a function of galactocentric distances in spiral galaxies have provided essential constraints on the evolution of these systems. The majority of these studies, greatly relying on oxygen abundances of H ii regions, hint at an inside-out growth of galactic disks (Pilyugin et al. 2014; Bresolin & Kennicutt 2015; Ho et al. 2015; Lian et al. 2018). In this section we present our results on the stellar abundance gradients for Fe, and the α elements measured in this work.

In Figure 6 we show in different panels the abundances as a function galactocentric distance for the individual elements, Fe, Ca, Ti, and Mg. We compute linear regressions for each element and show the inferred slopes with a solid line. Given that the theorems underlying least square regressions assume asymptotic normality (meaning large-N samples), a condition that does not apply to our limited sample, we estimate our abundance gradients and their uncertainties applying a linear regression with nonparametric bootstrapping, where we apply the bootstrapping on the residuals, not the fitted parameters themselves. Briefly summarized, we find the optimal linear regression on our original data set, extract the residuals from the best fit, generate new equal-size samples using the residuals, and fit linear regressions on the new samples. This resampling method allows for a more accurate estimate of the uncertainties in our inferred parameters.

Zoom In Zoom Out Reset image size

We measure an Fe gradient of −0.124 ± 0.034 dex kpc⁻¹.We find a slightly shallower gradient for Ca, −0.091 ±0.026 dex kpc⁻¹. Finally, the gradients for Ti and Mg show a much larger scatter, with values of −0.117 ± 0.044 dex kpc⁻¹ and −0.069 ± 0.089 dex kpc⁻¹, respectively. We estimate a weighted average α-abundance gradient of −0.096 ± 0.009 dex kpc⁻¹, which is only slightly flatter than that observed for Fe. It is important to highlight this subtle difference as it is expected that these elements, Fe-peak and α, with different nucleosynthetic origins, imprint a different signature in the observed gradients. For example, in the past such a behavior has been predicted and observed in the MW (Chiappini et al. 2001; Palla et al. 2020). Overall, negative abundance gradients are obtained when the star formation has been more efficient in the inner regions of the galaxy, compared to those in the outer regions of the disk, which based on our study, appears to be the case in NGC 1313.

Lastly, we perform a final test as the abundance gradients in Figure 6 appear to be primarily driven by the individual point at R/R₂₅ ∼ 0.9. We remove this individual point and estimate new gradients with the remaining points at R/R₂₅ < 0.8. Overall we find that the abundance gradients for most elements remain within the uncertainties quoted in Figure 6, with the exception of Mg. The gradient for this particular element becomes slightly positive with much larger uncertainties, strongly highlighting the effects of the large scatter on the inferred Mg gradient.

Past studies have investigated the distribution of metals in NGC 1313. Pagel et al. (1980) analyzed the observations of six H ii regions and found a uniform gas-phase abundance distribution, with no clear gradient. Similarly, through the spectroscopic analysis of a much larger sample of 33 H ii regions, Walsh & Roy (1997) reported a flat gas-phase abundance distribution across the disk. These results made NGC 1313 the most massive barred galaxy to lack a radial abundance gradient in the ionized-gas component. We highlight that these abundance studies relied on strong-line calibrations, which are based on the flux ratios of some of the strongest forbidden lines. It is well known that different strong-line calibrations return different abundance gradients with systematic offsets in the measured abundances (Kewley & Ellison 2008; Maiolino & Mannucci 2019), particularly in high-metallicity environments (Bresolin et al. 2016). Gas-phase metallicities estimated using the T_e -based method are considered to be more reliable than those inferred from strong-line calibrations (Maiolino & Mannucci 2019).

In Figure 7 we compare our stellar weighted average α-abundance gradient to the much shallower gradient inferred from the H ii regions studied in Walsh & Roy (1997; hereinafter WR97), −0.016 ± 0.008 dex kpc⁻¹. The nebular abundances were estimated using the strong-line [O iii]/[N ii] calibration by Pettini & Pagel (2004) along with the measured fluxes by WR97. From this figure, it is clear that the overall trends imprinted in the stellar and gas-phase abundances are remarkably different. A puzzling effect as both of these components are expected to describe the present-day metallictiy distribution of NGC 1313. Studies such as that by Bresolin et al. (2016) and Hernandez et al. (2021) have previously shown that strong-line indicators produce nebular abundance gradients that are significantly shallower than those inferred from young stellar populations.

Zoom In Zoom Out Reset image size

We expand our nebular and stellar abundance gradient comparison to include the direct abundances for six H ii regions studied in Hadfield & Crowther (2007). In order to perform a meaningful comparison between the direct oxygen abundances and those measured from our YMC sample, we must correct the oxygen gas-phase metallicities for depletion onto interstellar dust grains. Overall, studies estimate depletion factors for oxygen that are between −0.08 and −0.2 dex (Jenkins 2009; Mesa-Delgado et al. 2009; Peimbert & Peimbert 2010; Simón-Díaz & Stasińska 2011). Similar to the work of Bresolin et al. (2016), for simplicity we adopt a depletion correction factor of −0.1 dex. In Figure 7 we show with square symbols the depletion-corrected [O/H] abundances by Hadfield & Crowther (2007). We note that due to the limited size of their H ii-region sample, Hadfield & Crowther (2007) report no evidence for a radial dependence in their inferred abundances. Figure 7, however, hints at a possibly steep gradient for the ionized gas primarily driven by the single abundance measurement at R ∼ 3 kpc. This figure clearly shows that the direct oxygen abundances appear to be much better aligned to the stellar abundance gradient, than to the strong-line abundance trend by WR97. To further confirm the presence of such a steep gradient in the gas-phase metallicities in NGC 1313, the sample of H ii regions with available [O ii] λ3727, [O iii] λ4363, and [O iii] λ5007 nebular emission features (to obtain direct abundances) needs to be expanded to cover galactocentric distances close to the nucleus, as well as in the outskirts of the galaxy. Given the limitations of the strong-line calibrations, it is critical that abundance studies comparing the gas-phase and stellar populations adopt whenever possible T_e -based nebular calibrations, otherwise, any information on their chemical evolution or mixing timescales might be misconstrued.

Age–metallicity relations are important for understanding the chemical evolution of individual galaxies. These types of relations are able to shed some light on correlations between the age of the stellar clusters and possible interactions or close encounters with other systems. Livanou et al. (2013) studied the age–metallicity relation in the LMC by analyzing a sample of small open clusters. They report a metallicity gradient with higher abundances for the younger clusters. Such a trend is of course expected in the normal evolution of stellar populations in galaxies. Livanou et al. (2013) also identify a considerable increase in metallicity around the ages of ∼600 Myr, which they propose could be the result of the most recent encounter between the LMC and the SMC.

Our sample of YMCs in NGC 1313 covers ages between 20 and 300 Myr. In Figure 8 we display the age–metallicity trend for NGC 1313 with orange triangles, and compare it to the relation observed in the LMC (shown as blue squares) for clusters with comparable ages. Overall, there is a clear trend where clusters in NGC 1313 with ages ≥100 Myr began to show a steady decline in their metallicities. One YMC appears to deviate from this subtle trend, NGC 1313–379, with an age of 55 Myr and a metallicity of [Fe/H] = −0.84 ± 0.07 dex. Interestingly, through the study of resolved stellar populations Silva-Villa & Larsen (2012) find that the star formation histories suggest a burst in the southwestern region, right at the location of NGC 1313–379. This in contrast to the rest of the galaxy, the bar and arms, which exhibits a constant star formation history. Their findings support a scenario where the southwest region of NGC 1313 is experiencing an interaction with a relatively small companion satellite. Previous to the work by Silva-Villa & Larsen (2012), several other studies had found evidence of satellite interactions in this same location (Blackman 1981; Peters et al. 1994; Larsen et al. 2007). The age–metallicity relation observed as part of this work (Figure 8) supports a scenario of constant star formation history throughout the galaxy, with a possible recent burst in star formation close to the region where NGC 1313–379 is located.

Zoom In Zoom Out Reset image size