The Clouds in the "Middle Ages" I. Star formation history and chemical abundances in the LMC

Introduction

The formation and evolution of galaxies is driven by two processes, the progress of star formation from the cold gas and the accretion resulting from the interaction with other systems. The Large and Small Magellanic Clouds (LMC and SMC) are the most massive and brightest dwarf galaxies and probably approached the Milky Way only a few Gyr ago (e.g. Patel et al. [48], but see Vasiliev [63] for a second passage scenario). The Clouds have also interacted with each other multiple times (e.g. Bekki et al. [1], Diaz & Bekki [14], Hammer et al. [26], Cullinane et al. [10]) and due to their proximity (50–60 kpc; de Grijs et al. [11], de Grijs & Bono [11]) they represent a laboratory to study the consequences of dwarf-dwarf galaxies interactions in detail, and the interaction with a massive host. Due to their relatively low average metallicity (about one-half and one-fifth of the solar iron abundance for the LMC and SMC, respectively) they are also fundamental benchmarks for stellar properties and stellar evolution.

The most intense period of star formation in the LMC occurred 0.5–4 Gyr ago (e.g. Ruiz-Lara et al. [54], Mazzi et al. [41], Burhenne et al. [4]) and a steep increase in the number of star clusters confirmed a burst 2–3 Gyr ago (e.g. Pieres et al. 2016). During this time, distinct peaks of star formation also occurred in the SMC (Massana et al. [39]) and its star clusters (Dhanush et al. [13]), suggesting a link between enhanced star formation and tidal interactions, whereas at ages > 5 Gyr the LMC and SMC might have evolved independently. Most of the SMC mass was assembled prior to 3.5 Gyr ago (Rubele et al. [53]). Morphological maps of high spatial resolution (about 0.15 kpc) show a well-defined LMC bar with a clear central density embedded in a disc-like distribution, a major SMC inner over-density and several protuberances (El Youssoufi et al. [15]), depending on the stellar population. The periphery of the Clouds reveals also a wealth of substructures in kinematics and metallicity which were mostly produced by multiple LMC-SMC interactions (e.g., Belokurov & Erkal [2]). The formation of northern and southern LMC arms are likely associated to an interaction at about 1 Gyr ago whereas the western LMC truncation and eastern SMC extensions could be associated to interactions at about 400 and 150 Myr ago, respectively (e.g. Cullinane et al. [10]). Metallicity maps confirmed a metallicity gradient in the LMC (e.g. Grady et al. [24], Omkumar et al. [45], Povick et al. [50]) and various asymmetries, including a metal-rich ring-like feature possibly of intermediate age (Choi et al. [6]). The stellar motion within morphological and metallicity features depicts the disruption of the SMC as imprinted both on its outer and inner regions (e.g. Vijayasree et al. 2026). Asymptotic giant branch (AGB) stars, are typical intermediate-age tracers of the Clouds, as the majority of the red giant branch (RGB) and red clump stars in the Clouds are on average $>3$ Gyr old, see El Youssoufi et al. [15] for the median ages of the different stellar populations. AGB stars are low- to intermediate-mass stars experiencing shell He and H burning with an atmosphere enriched in carbon and/or oxygen products and a large envelope where dust may also form. Atomic and molecular features dominate the spectra of carbon- (or C-rich) and oxygen- (or O-rich) AGB stars and unmistakably allow to distinguish between the two types. The formation of carbon stars is directly linked to the number and efficiency of the third dredge up process which brings elements synthesised near the core to the surface. This process is more efficient at low metallicity (see e.g. Karakas et al. [33]) and in stars with masses between $\sim1.5$ and $3\,M_{\odot}$ (Pastorelli et al. [47]), which means that carbon stars trace a specific age interval (from about 300 Myr to 2 Gyr). On the contrary, O-rich AGB stars span a larger age range (up to 10 Gyr). Most AGB stars intrinsically pulsate with long periods (from a few tens to several hundred days) and large amplitudes (from a fraction to a few magnitudes depending on waveband). Their pulsations are powered by convective processes (see e.g. Freytag et al. [20]). Furthermore, they are notable producers of dust and contribute, due to their strong stellar winds, to both the chemical enrichment of galaxies and the integrated luminosity of intermediate-age stellar populations (e.g., [32, 40, 37]). These properties corroborate their importance as tracers of stellar evolution, stellar populations and distances (e.g. [47, 35, 31, 59]).

Among AGB pulsators, Mira variables offer a valuable means of determining ages (see e.g. Sanders et al. [55]). These stars experience thermal pulsations in the fundamental mode, with periods ranging from 80 to 1000 days and typical K-band amplitudes of at least 0.4 mag (Habing & Olofsson [25]). Miras are known to exhibit well-defined period–luminosity (PL) relationships (see e.g. [21, 17, 30, 29]), because stars of a fixed mass only enter the fundamental-mode pulsation stage within a restricted range of radii (Trabucchi et al. [60]). For analogous reasons, theory predicts that Miras should obey a period–age relation (see e.g. [18, 59]). Trabucchi & Mowlavi [59] combined nonlinear pulsation calculations with isochrone models. These models predict that the fundamental mode period decreases with increasing age and the authors give for O-rich and C-rich models the corresponding best-fit relations.

Observational support for such a relation has existed for decades, stemming from observed links between Mira periods and their kinematic properties, as well as between periods and population scale heights (Feast [16]). Since hotter kinematics generally correspond to older stellar populations, these correlations are understood as evidence of a period–age connection, further validated by studies of Mira variables in stellar clusters ([22, 67]). Empirical calibrations have been derived from data in the solar neighbourhood ([19, 5, 23]). More recently, Zhang & Sanders [67] fitted dynamical models to Gaia DR3 Mira variables, linking their radial and vertical velocity dispersions to their period, and then comparing these dispersions with age-dependent dispersions from main-sequence turn-off and subgiant stars. This allowed them to derive a period–age relation consistent with earlier work and with cluster member data. Based on this work, Sanders et al. [55] traced the SFH of Mira variables in the nuclear stellar disc (see also Schultheis et al. [56]).

Determining SFHs from long-period variables (LPVs) has the advantage that sources are relatively bright, but are less numerous, than for SFH studies based on fitting the colour–magnitude-diagram (CMD) resulting from deep imaging (e.g. [54, 41]). However, a time-domain survey is needed to obtain well-sampled light-curves for determining accurate periods – especially long ones – whereas for CMD-fitting single-epoch two-waveband observations do usually suffice. Previously, Rezaeikh et al. [52] derived the SFH from the K-band magnitude distribution of LPVs, which was converted to a mass-age relation using the stellar models from Marigo et al. [38]. In that case, sparsely sampled light-curves were used to identify the LPV nature of the sources. In this project, we will take advantage of period–age relations to trace the SFH in the LMC by applying a similar method as in Sanders et al. [55], since the calibrated period–age relations demonstrate that Mira variables span ages from 1–10 Gyr making them excellent probes of the formation and distribution of the stellar population. The two approaches also differ in their systematic sensitivities. The luminosity-function method of Rezaeikh et al. [52] requires an assumed initial mass function (IMF) to invert the observed luminosity function into a star-formation rate, and relies on a bolometric luminosity-age relation that is itself strongly metallicity dependent; it further requires deep, highly complete photometry, since the oldest, faintest AGB stars are the most readily lost to incompleteness. The period–age approach instead assigns each star an age directly from its pulsation period: IMF dependence enters only through the exponent weighting the period-dependent duration of the Mira phase (Sect. 2.2), which we show has a negligible effect on our results.

This is the first work from the Clouds in the "Middle Ages" project that aims to obtain line-sight-velocities, stellar parameters and elemental abundances of AGB stars in the Clouds for tracing kinematical and chemical patterns associated to their SFHs. The patterns of elemental abundances and their gradients, with respect to similar quantities derived from RGB stars and younger stars, will produce vital information to develop chemical evolution models and explain the nature of structural features across the galaxies, supported or contrasted by the dynamical behaviour of the stars derived from the line-of-sight velocities and the proper motions. Galactic chemical evolution tries to explain how the chemical elements have formed and their abundances have evolved in time and distributed in space. In particular, they aim at following the evolution of the gas and its chemical composition inside galaxies. Such an evolution depends on the SFH, stellar nucleosynthesis and possible gas flows.

The paper is organised as follows. Sect. 2 shows the calculation of stellar ages for a subset of AGB stars, Mira variables, and the resulting SFH. Section 3 presents the computation of Mg, Si and Ca abundances from APOGEE spectra. Finally, Sect. 4 provides a perspective on upcoming analyses of the Clouds.

Mira variables as indicators of stellar ages

Data

We use the Gaia DR3 Long-Period Variable (LPV) star catalog of Lebzelter et al. [34], which is the largest all-sky LPV catalog to date. Based on 34 months of Gaia data, this catalog contains about 1.7 million LPVs over the full sky for a period range between 35–1 000 days and amplitudes in the G band larger than 0.1 mag. Based on the Red Photometer (RP; 640–1050 nm and $R \sim50$) spectra, they performed a classification of C stars which are flagged in the Gaia LPV catalog.

We have selected our LMC sample as in Nidever et al. [42], using equatorial coordinates and a radial velocity selection ($160 < v_{\mathrm{rad}} < 350\ \mathrm{km\,s^{-1}}$) to ensure membership to the LMC. Our sample consists of 9588 LPVs which we cross-matched with the 2MASS catalogue. We include C-rich and O-rich stars using the Gaia DR3 classification flag is_cstar.

LPVs occupy different locations in the PL relation based on their pulsation mechanism. Wood [66] demonstrated that the different PL sequences correspond to different radial pulsation modes. Miras were unambiguously identified as fundamental mode pulsators on sequence C, while semiregular variables on sequences A and B were found to pulsate in the higher overtones. Subsequent refinements introduced additional sequences (A′, C′) and suggested that sequences such as A and B are composed of multiple, closely-spaced sub-sequences (e.g. Ita et al. [30]). Soszyński et al. [57] found a more precise characterisation of the PL relation based on the Wesenheit index ($W_{\mathrm{K}}$) which is a reddening independent quantity and therefore the period-$W_{\mathrm{K}}$ relations can be used as distance indicators even for heavily reddened Mira variables.

Figure 1 shows the period-$W_{\mathrm{K}}$ diagram of our LMC sample revealing the different pulsational sequences, and in particular the fundamental mode sequences (indicated by the green line for O-rich and the orange line for C-rich LPVs, Soszyński et al. [57] – with their uncertainties). We have selected our sample of fundamental pulsating LPVs within $\pm0.5$ mag around these relations (indicated as the filled areas in Figure 1). Note that these sequences also include fundamental-mode pulsating semiregular variables (SRVs). This sample consists of 2047 stars with a mean period of $\sim290$ days and a mean peak-to-peak amplitude in G of $\sim0.6$ mag.

As we will use the period-age relations for Mira variables to retrieve the SFH (see Sect. 2.2), we cross-matched this sample against the OGLE-III LMC catalogue of LPVs (Soszyński et al. [58]), which classifies sequence-C pulsators into genuine Miras and SRVs directly from well-sampled I-band light curves, using the standard peak-to-peak amplitude boundary $\Delta I > 0.8$ mag. We derived an empirical relation: $\Delta G = 0.149 + 0.959 \times\Delta I$ with a residual scatter of 0.18 mag. Applying this relation to OGLE’s own $\Delta I > 0.8$ mag Mira boundary gives an equivalent boundary of $\Delta G \gtrsim0.91$ mag for our sample (see Sect. A for further details). We find that only 251 stars ($\sim12\%$) meet this Mira criterion.

Star formation history

We make use of the period-age relation for Mira variables to construct the SFH of the LMC following the approach of Sanders et al. [55]. We employ two distinct period-age calibrations. The first is the polynomial relation of Trabucchi & Mowlavi [59], based on theoretical stellar evolution models:

where $\tau$ is the stellar age in years and $P$ is the pulsation period in days. The second is the empirical sigmoidal calibration of Zhang & Sanders [67]:

where $\tau$ is in Gyr. This relation was calibrated exclusively on O-rich, kinematically-selected Miras in the Milky Way disc; the Zhang & Sanders [67] sample explicitly excludes SRVs via an amplitude and period selection as well as C-rich stars. We therefore evaluate this calibration only on the subset of our sample satisfying Zhang & Sanders [67]’s own published O-rich amplitude and period selection criteria, whereas the Trabucchi & Mowlavi [59] relation, whose calibrating sample included both Miras and SRVs as well as O- and C-rich stars, is applied to the full chemistry-inclusive sample (Sect. 2.1).

It is worth noting that, independently of the amplitude selection, the two calibrations also differ in the population and host galaxy they were calibrated on. Trabucchi & Mowlavi [59] calibrate their relation on $\sim95$ fundamental-mode pulsators in star clusters spanning the Galaxy, the LMC (19 clusters), and the SMC (NGC 419), with ages taken from published cluster isochrone fits; the authors explicitly state the relation is "environment-dependent, and not necessarily universal". Zhang & Sanders [67], by contrast, calibrate their relation from the kinematics of $\sim8300$ field Mira candidates in the Milky Way disc, using a dynamical model to convert velocity dispersion in period bins into ages; their own spatial selection explicitly excludes the LMC’s sky position and distance range from that sample. Their relation was therefore calibrated on zero LMC stars, using Milky Way disc dynamics, and applying it to the LMC assumes the period-age relation is portable across host galaxies and environments — an assumption independent of, and larger than, the amplitude selection.

To properly propagate observational and systematic uncertainties, we implement a Monte Carlo (MC) simulation with $N_{\mathrm{MC}} = 100$ iterations. Each iteration applies three independent perturbations to the data.

First, period measurement uncertainties are propagated through Gaussian perturbations:

where $\sigma_{P_i}$ is the measured period uncertainty for star $i$. This accounts for observational errors in the period determination.

Second, completeness corrections account for the sample’s detection efficiency. We assume a typical 90% completeness based on the recovery rate of Gaia LPVs with respect to OGLE and ASAS-SN data [34], modeled as:

clipped to the range [0.70, 1.00] to avoid unphysical values.

Third, bootstrap resampling accounts for statistical uncertainty arising from the finite sample size. For each MC iteration, we generate a new dataset of equal size to the original by randomly drawing stellar indices with replacement from the parent sample. This procedure naturally propagates Poisson noise through the entire SFH reconstruction pipeline—including period-age conversion, weighting, and binning—without making assumptions about the underlying age distribution.

To convert the number density in age into a star formation rate, we follow Sanders et al. [55], who showed that the number of Mira variables per unit age is proportional to $\tau^{0.48}$ times the star formation rate (SFR):

The $\tau^{0.48}$ factor arises from the combination of the initial mass function (IMF) and the mass-dependent duration of the Mira phase, which we assume to be constant with age. Because both the IMF slope and the mass-dependent Mira-phase duration are metallicity dependent, and this exponent was derived by Sanders et al. [55] for the Milky Way’s nuclear stellar disc – a different environment, star-formation regime, and metallicity from the LMC – we tested its impact by recomputing the SFH for $n \in[0.30, 0.60]$. The resulting star-formation-rate-weighted mean age changes by less than 10% for the Trabucchi & Mowlavi [59] calibration and by less than 1% for the Zhang & Sanders [67] calibration over this range, and none of our qualitative conclusions depend on the specific value adopted (Figure 2).

For each MC iteration, the total weight for each star is:

where $C$ is the survey completeness factor, and $s(\log P)$ is a spline-based systematic correction accounting for uncertainties in the period-age relation. The spline correction is modeled as a Gaussian process with covariance:

with parameters $\sigma_{\mathrm{spline}=0.03}$, $\sigma_{\mathrm{corr}}=0.005$, and correlation length $\ell=2$ in knot space. This captures both independent and correlated systematic uncertainties.

The total SFH in each age bin $j$ is then computed as:

where $\mathbb{I}(\tau_i \in\mathrm{bin}_j)$ is the indicator function that equals unity if the age $\tau_i$ falls within the age interval defining the $j$-th bin, and zero otherwise. The denominator ensures normalization where

represents the relative SFR per Gyr.

Uncertainties are propagated via the MC simulations, with final uncertainties given by the 16th and 84th percentiles of the resulting distributions:

where $\mathrm{SFR}_j^{(k)}$ is the SFH in bin $j$ computed from the $k$-th MC iteration.

Results

Figure 3 and 4 present the reconstructed SFH of the LMC, derived respectively from the full fundamental-mode-pulsator sample and from the amplitude-selected Mira-only sample (Sect. 2.1), comparing the two period-age calibrations with independent SFHs. Ruiz-Lara et al. [54] used photometric data from the optical Survey of the Magellanic Stellar History (SMASH; Nidever et al. [43]) reaching sources down the main sequence turnoff and adopting a CMD-fitting technique to derive the global SFH in the LMC. They found a distinct population in the spiral arms with an age of about 2 Gyr. Mazzi et al. [41] used the photometric near-infrared data from the VISTA Survey of the Magellanic Clouds system (VMC; Cioni et al. [8]) and compared the observed CMDs with synthetic ones derived from stellar population analysis with TRILEGAL¹. They found that the period of most intense star formation occurred between 0.5 and 4 Gyr, with an enhanced SFR at ∼ 1 Gyr.

A striking result of our analysis is the significant disagreement between the two period-age calibrations, with the Trabucchi & Mowlavi [59] and Zhang & Sanders [67] relations yielding substantially different SFH shapes, amplitudes, and peak ages. We find especially significant differences at the oldest ages where Zhang & Sanders [67] finds a peak in the SFR in the age-range of 8–10 Gyr.

Neither calibration reproduces the external CMD-fitting SFHs of Mazzi et al. [41] and Ruiz-Lara et al. [54] equally well across all ages: the Trabucchi & Mowlavi [59] relation more closely matches their intermediate-age (0.5–4 Gyr) SFR level and shape, while the Zhang & Sanders [67] relation is more consistent with the substantial old population implied by the short, ancient formation episode found by Harris & Zaritsky [27]. We therefore do not adopt either relation as a single fiducial calibration; instead we present results under both

throughout, and treat their difference as our primary estimate of the systematic uncertainty associated with the choice of period-age calibration.

As a further check, we compared the SFH derived from our amplitude-selected Mira-only sample (Figure 4) with that obtained from the full fundamental-mode-pulsator sample (Figure 3), i.e. without applying the amplitude criterion used to isolate genuine Miras from other fundamental-mode variables. The two samples agree on the qualitative pattern of disagreement between calibrations: under Trabucchi & Mowlavi [59] both show their dominant contribution near 1 Gyr, and under Zhang & Sanders [67] both show enhanced star formation in the oldest age bins. The amplitude selection concentrates these features and the two calibrations show a better agreement for ages <2 Gyr.

We note that, for the amplitude-selected Mira sample, the old age bins contain fewer than five stars each under both calibrations (marked explicitly in Figure 4); while the old-age excess is qualitatively recovered independently in both the Mira-only and the full fundamental-mode samples, its precise amplitude in the Mira-only sample should be interpreted with this small-number caveat in mind. We nonetheless adopt the amplitude-selected Mira sample as our primary result, since the period-age calibrations used here (Trabucchi & Mowlavi [59]; Zhang & Sanders [67]) are most directly applicable to fundamental-mode Mira pulsators specifically.

The Mira-based SFH reveals three distinct phases of star formation activity, which we discuss in detail below. The oldest component of the LMC’s SFH contributes approximately 20% of the total SFR formed over the LMC’s history. The presence of a significant old population in the LMC is also consistent with the SFH derived from Harris & Zaritsky [27], where a significant portion of the LMC’s stellar mass was formed during a short period. Rezaeikh et al. [52], using also LPVs, found a 10 Gyr old single star formation epoch, while in the LMC bar they identified a second distinct star formation episode starting around 3 Gyr ago. Both events are captured in our study. We also examined whether our sample could be split spatially (e.g. bar versus outer disc), as suggested by the spatially-resolved SFHs of Ruiz-Lara et al. [54] and by this distinct bar episode of Rezaeikh et al. [52]; however, our classification-confirmed Mira sample is too small to support such a subdivision without individual age bins becoming dominated by only a handful of stars, and we defer this test to future work with a larger sample (Sect. 4). Interestingly, we see that for the intermediate-age population ($\sim$3–6 Gyr) the SFR is low with no obvious evidence of a star burst in line with Mazzi et al. [41] and Ruiz-Lara et al. [54]. One major feature of our reconstructed SFH is a recent starburst with a peak at $\sim$1–2 Gyr coinciding very nicely with the discovery of a young population at $\sim$ 2 Gyr in Mazzi et al. [41]. A comprehensive study based on data from the Hubble Space Telescope, whilst suggesting significant spatial variations, confirms a star formation activity at $\sim$2 Gyr (Cohen et al. 2024).

To summarise, based on the period-age relation of fundamental pulsating AGB stars, we can reconstruct the SFH of the LMC consistently with other studies using global CMD-fitting techniques. The two period-age calibrations show broad consistency in their overall shape but differ in the relative contributions of the ancient and intermediate-age populations, reflecting the dominant systematic uncertainty in this analysis. This discrepancy, which exceeds the statistical uncertainties by a factor of 2–3, shows that the choice of period-age calibration is the dominant systematic uncertainty in the Mira-based SFH reconstruction.

Chemical abundances of AGB stars in the LMC

Spectroscopic surveys of the Clouds based on AGB stars remain scarce. Olsen et al. [44] targeted $\sim$6 000 giants, including many AGB stars, revealing a counter-rotating population of metal-poor stars potentially accreted from the SMC. van der Marel et al. [61] exploited $\sim$1 000 spectra of C-rich stars to provide key constraints on the LMC’s kinematics, geometrical structure, total mass, and tidal radius.

The Apache Point Observatory Galactic Evolution Experiment (APOGEE; Majewski et al. [36]) produced high-resolution near-infrared spectra for 4 000 AGB stars across the Clouds (Nidever et al. [42]), but the abundance of multiple elements, their trends and interpretation using a chemical evolution model have been obtained only for RGB stars (Nidever et al. [42]; Hasselquist et al. [28]).

The APOGEE Stellar Parameter and Chemical Abundances Pipeline (ASPCAP) has been specifically designed to determine stellar parameters and chemical abundances from APOGEE observations. It uses two steps to analyze the data: First, the pipeline performs a global fit to the entire APOGEE spectrum to determine the key stellar atmospheric parameters. This involves simultaneously solving for effective temperature, surface gravity, microturbulence, overall metallicity ([M/H]), $\alpha$-elements ([α/M]), carbon ([C/M]), and nitrogen ([N/M]). With the stellar parameters fixed from the first stage, ASPCAP then derives the individual chemical abundances. It refocuses the fit on specific small wavelength regions ("windows") in the spectrum that are dominated by the spectral lines of a particular element. These windows are chosen because they contain spectral features that are dominated by transitions of the element of interest. The pipeline applies a mask to analyze only these specific portions of the spectrum for the abundance determination

of each element. The “weights” for these windows, which can be assigned based on the reliability of the line or level of contamination from other features, are provided in the APOGEE data release.

Sample

We selected our sample of stars based on their equatorial coordinates and radial velocities in the same way it was done for the Gaia DR3 light curves (see Sect.2.1). In order to ensure to have high quality data, we excluded stars with STARFLAGS corresponding to BAD_PIXELS, PERSIST_HIGH, VERY_BRIGHT_NEIGHBOR, ASPCAP flags with NO_ASPCAP_RESULTS, STAR_BAD and SNR > 20. We refer the reader to the APOGEE DR17 documentation² for more information. Our final sample consists of 567 stars. We cross-matched our sample with the Gaia LPV catalogue from Lebzelter et al. [34] where we found 296 sources in common. We considered only those as AGB stars and classified the rest as RGB stars. As mentioned in Sect. 2.2, the typical completeness of Gaia is 90% for recovering O-rich AGB stars.

Figure 5 shows the Kiel diagram of our AGB/RGB star sample where the data reveal a clear separation between the two populations, with AGB stars generally occupying a region of lower effective temperatures (T_eff ≤ 4000K) compared to RGB stars. The mean metallicities of the RGB stars is [M/H] = −0.77±0.42 dex while it is [M/H] = −0.79±0.22 dex for the AGB stars showing no systematic difference between these two populations.

Chemical abundances

As in APOGEE we run the FERRE code with an improved version where we perform line-to-line individual abundances. We use the APOGEE input spectra together with their derived stellar parameters to feed in our pipeline. The algorithm begins by scanning through all wavelength positions in the mask file and identifying pixels where the target element’s weight meets the minimum threshold criterion. We have chosen the minimum weight to be 0.01 allowing to exploit weak lines. Secondly, the merging algorithm examines adjacent clusters and determines whether they should be merged based on the separation between their central wavelengths. Based on the instrumental profile of APOGEE, the threshold is set at 2.10 Å. Then we apply a quality filter where we calculate the line-to-line dispersion using a threshold of 0.2 dex. Lines with insufficient measurements or excessive scatter are rejected. Furthermore, measurements falling exactly on the boundaries of the synthetic grid (−0.5 dex and +0.75 dex, respectively) are also excluded. Since these edge values represent failed fits rather than true physical measurements, their removal is necessary not to affect the line-to-line dispersion computations. We specify in Tab. 1 the mean abundance for each line together with its standard deviation. Overall, we do not notice any significant differences in the mean abundance between RGB and AGB stars in our sample.

Figure 6 shows each individual Mg line plotted against [M/H]. The majority of the lines show a large intrinsic dispersion; hence those lines can not be trusted. Some lines, e.g. the Mg line at 15 419 Å or at 16 600 Å show just a large dispersion with no trend of [Mg/Fe] with metallicity. This could be due mainly to blends with strong molecules (e.g. CN, CO and OH) for these cool objects. Finally, for Mg only three lines (at 15 698 Å, 15 745 Å and 15 891 Å) fulfill our criteria and are used for our analysis (see Tab. 1). Figure 7 compares our Mg abundances, derived from selected lines, with the total Mg abundances from the optical sample of Van der Swaelmen et al. [62] and from APOGEE DR17 (lower left panel): the smaller dispersion in the Mg abundances at a given metallicity and the good agreement with the optical sample are clearly visible. Our line-by-line Mg abundances show a scatter about 25% smaller than that of the APOGEE DR17 pipeline values.

!Individual Mg lines as a function of [M/H] for RGB and AGB stars bbox=0.14,0.08,0.92,0.50

In addition, we do not notice any trends of the Mg abundances related to the effective temperature (Fig. D.1, left panels), indicating that AGB stars can be analysed in this way without introducing biases in the Mg abundance.

We refer to Fig. B.1 and Tab. 1 for details about the adopted Si lines. Compared with the total Si abundances from APOGEE DR17 (Figure 7, central lower panel), our abundances show a dispersion about 10% smaller than the APOGEE DR17 values. We find a trend of the Si abundances in the sense that cooler temperatures reveal lower Si abundances. Also, the temperature effect in APOGEE DR17 seems to be slightly larger than in our analysis (see Fig. D.1, middle panels).

element	$\lambda$ (Å)	$\sigma$ (dex)	RGB	AGB
Mg	15698	0.187	0.114 ± 0.12	0.11 ± 0.09
Mg	15745	0.163	0.096 ± 0.09	0.039 ± 0.18
Mg	15891	0.137	0.101 ± 0.17	0.018 ± 0.09
Si	15889	0.173	0.121 ± 0.22	0.058 ± 0.16
Si	15893	0.141	0.225 ± 0.15	0.114 ± 0.13
Si	16247	0.167	0.027 ± 0.19	0.104 ± 0.12
Si	16385	0.191	0.106 ± 0.18	0.04 ± 0.14
Ca	16141	0.159	0.041 ± 0.15	0.043 ± 0.15
Ca	16155	0.137	0.082 ± 0.14	0.077 ± 0.15
Ca	16161	0.117	0.142 ± 0.10	0.140 ± 0.11

We were able to successfully exploit all three Ca lines within our observed spectra, since they provide good fits (see Fig. C.1 and Tab. 1 for more information). We see here similar patterns as for the other $\alpha$-elements: an improvement in the intrinsic dispersion of the Ca abundances compared to APOGEE DR17 (Figure 7, lower right panel), with a scatter about 10% smaller, and no signs of a temperature trend (Fig. D.1, right panels).

The individual-line Mg, Si, and Ca abundances for all 567 stars in our sample, together with the adopted [Mg/Fe], [Si/Fe], and [Ca/Fe] values used throughout this paper, are made available in their entirety as a machine-readable table at the CDS (Table 2). Table 2 is published in its entirety at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc.u-strasbg.fr/cgi-bin/qcat?/J/A+A/; a portion is shown here for guidance regarding its form and content.

Chemical abundances of AGB for tracing the intermediate-age populations

As AGB stars are large-amplitude pulsators and the stellar atmospheres are dominated by dynamical processes, we investigate here for the three $\alpha$-elements if we see any correlation with stellar parameters, pulsation period and amplitude of the star.

To study whether the chemical abundances of Mg, Si, and Ca depend on stellar pulsation period, we examined the relationship between these key $\alpha$-elements and the fundamental pulsation period derived from our data set. Figure 8 presents the abundance distributions as a function of period for Mg (top panel), Si (middle panel), and Ca (bottom panel). A running mean analysis reveals a relatively flat trend across the entire period range for all elements with no sign of a correlation between the period and the abundance. The mean abundances of Mg are 0.020 ± 0.103 dex, of Si 0.015 ± 0.093 dex, and of Ca 0.073 ± 0.088 dex showing similar trends within the uncertainties. We also do not see any correlation between the amplitude of the star and the abundances nor with the stellar parameter such as T$_{\rm eff}$, log g, and [M/H]. We therefore conclude that AGB stars (O-rich Miras) can be used to trace chemical abundances and therefore the SFH. However, a careful selection of lines suitable for the analysis is essential.

Nidever et al. [42] and Hasselquist et al. [28] produced chemical evolution models to fit the chemical abundances of RGB stars derived from APOGEE spectra and several elements. They found a SFH peak 2–3 Gyr ago in agreement with a previous chemical evolution study by Bekki & Tsujimoto [1] who only looked at Ba and Mg abundances, also for RGB stars. The latest chemical evolution models use as input the SFH from external data (e.g. from CMD fitting). For example, Vasini et al. [64] found that the SFH by Harris & Zaritsky [27] best reproduced the current SFR, supernova and nova rates, gas mass and metallicity of the LMC from their chemical evolution models. More recently, Palla et al. [46] fitted r- and s-process elements in the SMC adopting the SFH from the VMC survey [53]. In a future work we will compute specifically trimmed chemical evolution models including our derived SFH (Sect. 2.2) together with our improved analysis of chemical abundances of the intermediate-age population.

Conclusions and future perspectives

In this paper, we have presented a first study of the SFH and chemical abundances of the LMC, with a particular emphasis on the AGB stellar population—an intermediate-age tracer that has been largely overlooked in previous stellar population studies. This is the first work from the Clouds in the "Middle Ages" project that uses AGB stars to find and characterise kinematical and chemical patterns associated to the SFH of the Clouds. We present the SFH under both (Figure 3 and Figure 4) and treat their difference as our primary estimate of the systematic uncertainty.

ID	$T_{\mathrm{eff}}$ (K)	$\log g$	[M/H]	Type	Mg15698	[Mg/Fe]	Si16247	[Si/Fe]	Ca16141	[Ca/Fe]	...
2M05434720-6644121	4081.68	2.034	−0.523	RGB	0.516	0.041	0.231	0.249	0.097	0.067	...
2M04413697-6818530	3555.29	0.009	−0.873	AGB	0.282	0.043	−0.230	−0.167	0.020	0.111	...
2M04324400-6937420	3536.81	0.199	−0.961	AGB	0.163	0.080	−0.063	−0.038	0.023	0.088	...
...	...	...	...	...	...	...	...	...	...	...	...

We reconstructed the SFH of the LMC using a chemistry-inclusive sample of 2047 fundamental-mode pulsating AGB stars (879 O-rich, 1168 C-rich) identified from the Gaia DR3 LPV catalog, applying two independent period-age calibrations: the theoretical relation of Trabucchi & Mowlavi [59] and the empirical calibration of Zhang & Sanders [67]. Our analysis reveals significant systematic differences between these calibrations, particularly at ages $>3$ Gyr where the Trabucchi & Mowlavi [59] relation shows some star formation only 3 − 6 Gyr ago while the Zhang & Sanders [67] relation predicts a prominent event $>6$ Gyr ago. Neither calibration reproduces the independent CMD-based SFH determinations ([54, 41]) equally well across all ages, so we do not adopt either relation as a single fiducial calibration; instead

An old-age population, and the disagreement between calibrations more generally, is recovered consistently whether we use our amplitude-selected Mira-only sample or the full fundamental-mode-pulsator sample. We identify three distinct phases of star formation: an old population at $\sim7$–$10$ Gyr contributing to $\sim20\%$ of the total stellar mass formed, a quiescent intermediate-age phase at 3–6 Gyr with no evidence of major bursts, and a recent starburst at $\sim1$–$2$ Gyr. These findings demonstrate that Mira variables can serve as effective tracers of the LMC’s SFH when carefully calibrated, though the choice of period-age relation remains the dominant systematic uncertainty.

Our chemical abundance analysis of 567 RGB and AGB stars observed by APOGEE, using an improved line-by-line abundance pipeline, demonstrates that AGB stars can reliably trace chemical enrichment despite their pulsational variability and complex atmospheric dynamics. We identified high-quality spectral lines for Mg (15 698 Å, 15 745 Å, 15 891 Å), Si (15 889 Å, 15 893 Å, 16 247 Å, 16 385 Å), and Ca (16 141 Å, 16 155 Å, 16 161 Å) that yield significantly reduced intrinsic scatter compared to the APOGEE DR17 pipeline results. Crucially, we find no statistically significant correlations between the derived $\alpha$-element abundances ([Mg/Fe], [Si/Fe], [Ca/Fe]) and pulsation period, pulsation amplitude, or stellar parameters such as $T_{\mathrm{eff}}$, $\log g$, or [M/H]. The mean abundances for the AGB sample ([Mg/Fe] = 0.020 ± 0.103 dex, [Si/Fe] = 0.015 ± 0.093 dex, [Ca/Fe] = 0.073 ± 0.088 dex) are consistent with those of RGB stars and with the optical sample of Van der Swaelmen et al. [62], indicating that AGB stars preserve the chemical signatures of their progenitor populations without substantial bias introduced by the pulsation phenomenon. Looking ahead, the 4-metre Multi-Object Spectroscopic Telescope facility (4MOST; de Jong et al. [12]), mounted at the 4m-class telescope VISTA, will allow us to make a significant step forward in using AGB stars to probe the formation and evolution of the Clouds. In particular, the 1001 Magellanic fields survey (1001MC; Cioni et al. [7]), a survey of the 4MOST consortium, will provide high-resolution, multi-object spectroscopy for thousands of stars across the entire LMC and SMC. This will enable detailed abundance analysis of AGB stars across all spatial regions, including the disc, bar, and outer halo, allowing us to map the chemical signatures of past interaction events with unprecedented precision. Combined with precise Gaia DR4 proper motions, and an enlarged sample of LPVs, these data will also permit dynamical age estimates independent of period-age relations, offering a powerful cross-calibration of the Mira period-age relation and a robust constraint on the LMC-SMC-Milky Way interaction history. We expect the number of LPVs in Gaia DR4 to increase by a factor of three with a major improvement on the number of small amplitude variables beyond the DR3 lower limit of the amplitude in G ≥ 0.1 mag.

Acknowledgements

Acknowledgements. This work was supported by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – [Project number 550131653]. This work is in the honor of Prof. Paola Marigo, a leader in the study of AGB stars and esteemed colleague, who left us prematurely.

Amplitude classification

We show the selection of Mira variables from the amplitude in the $G$ band as translated from the amplitude defined in the OGLE-III study [58]. Figure A.1 shows the $G$-band amplitude vs. period while Figure A.2 vs. OGLE-$I$-band amplitude.

Si abundances

We show the line-to-line comparison of Si abundances.

Ca abundances

We show the line-to-line comparison of Ca abundances.

Temperature dependence of $\alpha$-elements

We show the temperature dependence of the elements [Mg/Fe], [Si/Fe] and [Ca/Fe] vs. [Fe/H] for the AGB and RGB stars in our sample.