<figure><img src="image_1.png" /></figure>

# Metallicity of Galactic RR Lyrae from Optical and Infrared Light Curves. I. Period–Fourier–Metallicity Relations for Fundamental-mode RR Lyrae

Joseph P. Mullen¹, Massimo Marengo¹, Clara E. Martínez-Vázquez², Jillian R. Neeley³, Giuseppe Bono⁴,⁵, Massimo Dall’Ora⁶, Brian Chaboyer⁷, Frédéric Thévenin⁸, Vittorio F. Braga⁵,⁹, Juliana Crestani⁴,⁵,¹⁰, Michele Fabrizio⁵,⁹, Giuliana Fiorentino⁵, Christina K. Gilligan⁷, Matteo Monelli¹¹,¹², and Peter B. Stetson¹³

¹Department of Physics and Astronomy, Iowa State University, Ames, IA 50011, USA; jpmullen@iastate.edu

²Cerro Tololo Inter-American Observatory, NSF’s National Optical-Infrared Astronomy Research Laboratory, Casilla 603, La Serena, Chile

³Department of Physics, Florida Atlantic University, 777 Glades Road, Boca Raton, FL 33431, USA

⁴Dipartimento di Fisica, Università di Roma Tor Vergata, via della Ricerca Scientifica 1, I-00133 Roma, Italy

⁵INAF—Osservatorio Astronomico di Roma, via Frascati 33, I-00078 Monte Porzio Catone, Italy

⁶INAF—Osservatorio Astronomico di Capodimonte, Salita Moiariello 16, I-80131 Napoli, Italy

⁷Department of Physics and Astronomy, Dartmouth College, 6127 Wilder Laboratory, Hanover, NH 03755, USA

⁸Université de Nice Sophia-antipolis, CNRS, Observatoire de la Côte d’Azur, Laboratoire Lagrange, BP 4229, F-06304 Nice, France

⁹Space Science Data Center—ASI, via del Politecnico snc, I-00133 Roma, Italy

¹⁰Departamento de Astronomia, Universidade Federal do Rio Grande do Sul, Av. Bento Gonçalves 6500, Porto Alegre 91501-970, Brazil

¹¹IAC- Instituto de Astrofísica de Canarias Calle Vía Lactea s/n, E-38205 La Laguna, Tenerife, Spain

¹²Departamento de Astrofísica, Universidad de La Laguna, E-38206 La Laguna, Tenerife, Spain

¹³Herzberg Astronomy and Astrophysics, National Research Council, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada

Received 2020 October 23; revised 2021 March 16; accepted 2021 March 16; published 2021 May 14

## Abstract

We present newly calibrated period–$\phi_{31}$–[Fe/H] relations for fundamental-mode RR Lyrae stars in the optical and, for the first time, mid-infrared. This work’s calibration data set provides the largest and most comprehensive span of parameter space to date, with homogeneous metallicities from $-3 \lesssim [\text{Fe/H}] \lesssim 0.4$ and accurate Fourier parameters derived from 1980 ASAS-SN (V band) and 1083 WISE (NEOWISE extension, W1 and W2 bands) RR Lyrae stars with well-sampled light curves. We compare our optical period–$\phi_{31}$–[Fe/H] relation with those available in the literature and demonstrate that our relation minimizes systematic trends in the lower and higher metallicity range. Moreover, a direct comparison shows that our optical photometric metallicities are consistent with both those from high-resolution spectroscopy and globular clusters, supporting the good performance of our relation. We found an intrinsic scatter in the photometric metallicities (0.41 dex in the V band and 0.50 dex in the infrared) by utilizing large calibration data sets covering a broad metallicity range. This scatter becomes smaller when optical and infrared bands are used together (0.37 dex). Overall, the relations derived in this work have many potential applications, including large-area photometric surveys with James Webb Space Telescope in the infrared and LSST in the optical.

Unified Astronomy Thesaurus concepts: RR Lyrae variable stars (1410); RRab variable stars (1413); Pulsating variable stars (1307); Globular star clusters (656); Metallicity (1031)

Supporting material: machine-readable table

## 1. Introduction

RR Lyrae stars (RRLs) are the most widely used tracers of old (age >10 Gyr; Walker 1989) stellar populations in the Milky Way and Local Group galaxies (see, e.g., Catelan & Smith 2015 for a review). They can also be used as standard candles thanks to a well-defined $M_V$ versus iron abundance relation (Sandage 1990; Caputo 1998). The recent calibrations of theoretical (Marconi et al. 2015; Neeley 2017) and observational (Dambis et al. 2013; Muraveva et al. 2018; Neeley et al. 2019), period–Wesenheit and period–luminosity relations (in the optical and infrared, respectively) have revealed the true potential of these stars as high-precision distance indicators. These same studies have highlighted the role of metallicity in determining the absolute brightness of these variables, hence requiring period–luminosity–metallicity (PLZ) and period–Wesenheit–metallicity (PWZ) relations for a reliable estimate of their distances. Accurate measurements of RRL metal abundances, however, are hard to come by.

Nearly all currently available catalogs (see, e.g., Dambis et al. 2013 for a recent large compilation) tend to list [Fe/H] abundances derived with heterogeneous methods and often calibrated with different scales. Metallicities derived from

high-resolution (HR) spectra ($R \gtrsim 20,000$) offer the highest level of precision (~0.1 dex; e.g., Clementini et al. 1995; Fernley & Barnes 1996; Lambert et al. 1996; For et al. 2011; Nemec et al. 2013; Govea et al. 2014; Pancino et al. 2015; Chadid et al. 2017; Magurno et al. 2018, 2019; Gilligan et al. 2021) but require large amounts of telescope and analysis time and therefore exist only for a small number of stars. Medium-resolution spectroscopic surveys ($R \sim 2000$), however, are available for larger samples of RRLs. For example, over 2300 sources found in the Sloan Extension for Galactic Understanding and Exploration (SEGUE) survey (Yanny et al. 2009) and another ~3000 targets from the Large Scale Area Multi-Object Spectroscopic Telescope (LAMOST) DR2 survey (Deng et al. 2012; Liu et al. 2014) can also be leveraged to derive iron abundances. This is traditionally achieved with the $\Delta S$ method (Preston 1959), which relies on ratios between the equivalent widths of Ca and H lines. Recent applications of this method are shown in both Liu et al. (2020), which is based on a fit of theoretical models and utilizes LAMOST data, and Crestani et al. (2021, hereafter C21), which is based on empirical calibrators and is shown to be valid for both SEGUE and LAMOST data. These applications demonstrate the
<table><caption>**Table 1**
Calibration Data Sets</caption><thead><tr><th></th><th>ASAS-SN</th><th>WISE</th><th>Joint Sample</th></tr></thead><tbody><tr><td>Bands</td><td>V</td><td>W1</td><td>W2</td></tr><tr><td>RRab stars</td><td>1980</td><td>1083</td><td>707</td></tr><tr><td>Period in days (range)</td><td>0.36–0.89</td><td>0.36–0.85</td><td>0.36–0.85</td></tr><tr><td>Period in days (mean value)</td><td>0.57</td><td>0.57</td><td>0.57</td></tr><tr><td>[Fe/H] (range)</td><td>-3.06-(+0.36)</td><td>-3.06-(+0.36)</td><td>-3.06-(+0.36)</td></tr><tr><td>[Fe/H] (mean value)</td><td>-1.47</td><td>-1.46</td><td>-1.44</td></tr><tr><td>Number of epochs (range)<sup>a</sup></td><td>69–892</td><td>153–879</td><td>153–879</td></tr><tr><td>Number of epochs (mean value)</td><td>270</td><td>231</td><td>233</td></tr><tr><td>Magnitude (range)</td><td>9.55–17.41</td><td>7.85–14.76</td><td>7.87–14.11</td></tr></tbody></table>

**Note.**

a The distribution of epochs is recorded prior to removing any spurious photometric measurement, as described in Section 3.

viability of the ΔS method in providing homogeneous [Fe/H] abundances with an uncertainty of 0.2–0.3 dex.

While the ΔS method allows a dramatic expansion to the sample of RRLs with known metallicity to thousands of sources, it still requires spectra. As such, it is not directly applicable to the much larger number of RRLs that will be discovered within the Milky Way and the other Local Group galaxies by upcoming large-area photometric surveys such as the Rubin observatory Legacy Survey of Space and Time (LSST; Ivezić et al. 2019) in the optical and surveys that will be executed for the Roman telescope (Akeson et al. 2019) at near-infrared wavelength. A reliable and precise method to derive metallicities from photometric time series is necessary to enable an expansion of such measurements to distances where taking spectra is not possible at all, or in high-extinction environments that can only be probed photometrically by mid-infrared telescopes such as the James Webb Space Telescope (JWST; Gardner et al. 2006).

Early work by Jurcsik & Kovacs (1996, hereafter JK96) has demonstrated that at optical wavelengths the shape of the light curve of fundamental-mode RRLs (RRab) is related to their metallic abundance. In particular, they derived a linear relation connecting RRab's [Fe/H] abundance with period and low-order parameters in the Fourier decomposition of the star's light curve in the V band, with the phase parameter $\phi_{31} = \phi_3 - 3 \cdot \phi_1$ providing the most sensitive diagnostics. Further work by Nemec et al. (2013, hereafter N13), Martínez-Vázquez et al. (2016, hereafter MV16), Smolec (2005), Ngeow et al. (2016), and Iorio & Belokurov (2021, hereafter IB21) extended this analysis, respectively, to well-sampled RRab light curves obtained with the Kepler space telescope (Koch et al. 2010), to include stars in globular clusters (GCs), to the Optical Gravitational Lensing Experiment (OGLE; Udalski et al. 1992) I band, to the Palomar Transient Factory (PTF; Law et al. 2009) R band, and by using Gaia DR2 (Gaia Collaboration et al. 2018; Holl et al. 2018; Clementini et al. 2019) G-band light curves.

In this work, we take advantage of both the large sample of homogeneous HR [Fe/H] abundances in C21 and apply the ΔS calibration of C21 to the full medium-resolution LAMOST DR6 and SDSS-SEGUE data sets to build an extensive HR +ΔS metallicity catalog. We have cross-matched the variables in the HR+ΔS metallicity catalog with well-sampled photometric time series in the All-Sky Automated Survey for Supernovae (ASAS-SN; Shappee et al. 2014; Jayasinghe et al. 2018) and the Near-Earth Objects reactivation mission (NEOWISE; Mainzer et al. 2011) of the Wide-field Infrared

Survey Explorer (WISE; Wright et al. 2010). We have then derived novel period-$\phi_{31}$-[Fe/H] relations in the optical (V band) and, for the first time, mid-infrared (W1 and W2 bands). Our work shows that these relations can indeed be extended to the thermal infrared, where the light curves are mostly determined by the radius variation during the star's pulsation rather than the effective temperature changes that dominate in the optical wavelengths. As mentioned above, this will be crucial to allow the determination of reliable metallicities in upcoming space infrared surveys.

This paper is structured as follows. In Section 2, we describe in detail the data sets we adopt for our work: the HR+ΔS metallicity catalog utilizing the work of C21, the ASAS-SN and WISE time-series catalogs, and the light curves for a sample of Galactic GCs with known metallicity that we will use to validate our relations. In Section 3, we explain how our period-$\phi_{31}$-[Fe/H] relations are calibrated and validated. Our results are discussed in Section 4, where we assess the precision of the infrared and optical relations, compare our relations with previous ones found in the literature, and apply our method to measure the [Fe/H] abundance in the sample of Galactic GCs. Our conclusions are presented in Section 5.

## 2. Field and Globular Cluster Fundamental-mode RR Lyrae Data Sets

In this section, we describe the properties of the sample of RRLs with known spectroscopic [Fe/H] abundances that we have adopted to calibrate our period-$\phi_{31}$-metallicity relations. From this catalog, we have derived three calibration samples: one with V-band time series (ASAS-SN sample), one with photometric data available in the thermal infrared (WISE sample), and one with data in both wavelengths (joint sample). The properties of these calibration samples are described in the following sections and listed in Table 1.

### 2.1. Calibration Sample

As mentioned before, the sources selected for this analysis are chosen from both the HR metallicity catalog of C21 and the full medium-resolution LAMOST DR6 and SDSS-SEGUE data sets, from which the spectrum selection criteria and ΔS metallicity calibration of C21 have been applied. The resultant HR+ΔS metallicity catalog is composed of 8660 fundamental-mode field RRLs that have also been cross-matched with the Gaia EDR3 database (Gaia Collaboration et al. 2021), as well as a number of other publicly available data sets, using an
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Period distribution (left) and spectroscopic [Fe/H] distribution (right) of the different calibration data sets. The histogram labeled "joint sample" (hatched) corresponds to those stars in common between the ASAS-SN (optical V band, in black) and WISE (infrared W1, W2 bands, in gray) data sets.</figcaption></figure>

<figure><img src="image_3.png" /></figure>

algorithm specifically developed for sparse catalogs (Marrese et al. 2019).

The HR+ΔS metallicity catalog provides both a homogenized sample of RRab iron abundances gathered from various sources in literature (170 of which are derived from HR spectra) and similarly homogeneous new ΔS metallicity estimates, as the ΔS calibration of C21 is based in part on the metallicity of 111 of the aforementioned HR RRab stars. For a complete and detailed description of the HR metallicity catalog’s demographics, the ΔS calibration, and the spectrum selection criterion, we refer the reader to the C21 paper.

All metallicities in this catalog are based on the metallicity scale utilized by C21, which is based on HR spectroscopy, and the most updated iron line parameters, most of which have transition parameters derived in laboratory studies. This is the same scale utilized by For et al. (2011), Chadid et al. (2017), and Sneden et al. (2017) and will be the default scale used throughout this paper unless otherwise declared. Note that other literature field RRL HR works can be brought to the same metallicity scale with the addition of a simple offset. A full analysis of the offsets between various HR [Fe/H] scales is offered in C21. It is worth mentioning here that the often-used (Carretta et al. 2009, hereafter C09) [Fe/H] scale can be converted to this work’s scale with the addition of a small rigid shift of 0.08 dex.

In this work, we focus on the stars in the HR+ΔS metallicity catalog that have a match in the ASAS-SN and WISE surveys. Figure 1 shows the distribution of period and metallicity for the subset of stars with available optical (ASAS-SN) or infrared (WISE) good-quality light curves (i.e., passing the stringent photometric and Fourier decomposition criteria described in Section 3). The joint sample, also shown in the figure, comprises the smaller subset of RRLs with light curves available at both optical and infrared wavelengths. All samples cover the entire period range expected for RRab variables and are well representative of the metallicity of Galactic halo RRLs, with a mean [Fe/H] abundance of ≈−1.45 in both the ASAS-SN and WISE samples. Figure 1 also shows that all histograms retain both the low-metallicity [Fe/H] ≲−2.2 and the high-metallicity [Fe/H] ≳−0.7 tail present in the C21 catalog. For a discussion of the significance of populations in the Galactic halo, we refer to Fabrizio et al. (2019, hereafter F19). Here, we want to remark that the broad range in metallicity is an

important feature of our samples, as it ensures broad leverage for accurate calibration of the metallicity slope in our period–$\phi_{31}$–[Fe/H] relation.

## 2.2. Field RRL Optical Data

The optical (V-band) time-series data for this analysis have been extracted from the ASAS-SN survey. The first telescope of ASAS-SN came online in 2013, and telescopes have gradually been added for a total of 24 telescopes scattered across the world as of early 2020. At its current capacity, ASAS-SN can survey the entire sky every night, providing high-cadence V-band photometry, ideally suited to acquire long-term, densely populated light curves of RRLs.

Out of the 6079 variables in the HR+ΔS catalog with a match in ASAS-SN, we were able to extract good-quality ASAS-SN light curves for 1980 stars using the procedure described in Section 3. Of these, 967 are in the joint sample, also having good-quality light curves in the infrared. The left panel of Figure 2 shows the distribution of average V-band apparent magnitudes for both the entire optical data set (solid black) and the stars in the joint sample (hatched). The cutoff at $V \simeq 17.4$ is due to the limiting magnitude of the ASAS-SN survey. The joint sample is instead truncated at $V \simeq 16$ owing to the shallower photometric depth of the WISE survey. For a detailed description of the generation of the light curves from the ASAS-SN time series and calculation of average magnitudes for this data set, we refer to Section 3.2.

## 2.3. Field RRL Infrared Data

We obtained archival infrared photometric time series from the WISE mission W1 and W2 bands (3.4 and 4.6 µm, respectively). The primary WISE mission surveyed the entire sky in four infrared bands every 6 months from 2010 January to the end of the post-cryogenic phase in 2011 February. The spacecraft was then reactivated in 2013 September in just the W1 and W2 bands as the NEOWISE mission. Still having sensitivity similar to the full cryogenic phase, NEOWISE, as of early 2020, has yielded 12 additional full-sky survey epochs (having a minimum of 12 measurements per survey epoch). With a minimum of ~156 individual epochs (~12 from WISE and 144 from NEOWISE) for each position in the sky, collected over a baseline of almost 10 yr, the combined WISE
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Distribution of the average V-band apparent magnitude (left) and the average W1-band apparent magnitude (right). Each panel shows both the entire calibration sample (solid fill) and the joint sample (hatched fill, subset of stars having both ASAS-SN and WISE data).</figcaption></figure>

and NEOWISE photometry provides the most comprehensive full-sky catalog in the infrared. Hereafter in the paper, the combination of photometry from both the primary WISE and ongoing NEOWISE mission will be referred to for brevity as WISE.

Out of 6264 sources in the HR+ΔS with a WISE catalog match, we were able to generate 1083 good-quality light curves in at least one of the two adopted WISE bands (1083 and 707 stars in the W1 and W2 bands, respectively). Figure 2 (right panel) shows the WISE distribution of the average W1 apparent magnitude, for both the entire infrared data set (solid gray) and those stars in the joint sample (hatched). Although the WISE sensitivity limit is ~16 mag in the two bands of interest, the figure shows a sharp drop at W1 ≃ 14.7; fainter RRLs, while detected in the WISE and NEOWISE catalogs, tend to have noisy light curves that are then rejected by the quality control procedures described in Section 3. The histogram for the W2 band is similar, although the cutoff in the magnitude distribution happens at W2 ≃ 14.1 owing to the lower sensitivity in this band, resulting in noisier light curves. The magnitude range and other characteristics in the joint sample subset of WISE bands are nearly the same as the complete infrared data set owing to the large overlap with available ASAS-SN data.

## 2.4. Globular Cluster Data Set

Lastly, we select a data set separate from those used in deriving our period-φ₃₁-[Fe/H] relations to use as an independent check of this work (see Section 4.2). This data set is composed of eight GCs homogeneously spread between [Fe/H] = −1.1 and −2.3 dex with a sizable number of RRab stars. These clusters are listed in Table 2. The GC data set comes mainly from the homogeneous photometric database of P. B. Stetson<fnref n="14" /> (hereafter PBS), except for NGC 3201 data, which comes from Piersimoni et al. (2002). The data in the PBS database were collected from ground-based telescopes using archival data from 1984 to the present. The main telescopes and cameras that contributed the most to the V-band data used in this work include the following: AAO LCOGT 1 m (CCD), CTIO 0.9 m

<table><caption>**Table 2**
Globular Clusters</caption><thead><tr><th>Clusters</th><th>[Fe/H]<sub>C09</sub></th><th>RRab Stars<sup>a</sup></th><th>Epochs</th></tr></thead><tbody><tr><td>NGC 7078 (M15)</td><td>-2.33</td><td>64 (64)</td><td>223</td></tr><tr><td>NGC 4590 (M68)</td><td>-2.27</td><td>14 (13)</td><td>41</td></tr><tr><td>NGC 4833</td><td>-1.89</td><td>11 (11)</td><td>72</td></tr><tr><td>NGC 5286</td><td>-1.70</td><td>30 (25)</td><td>111</td></tr><tr><td>NGC 3201</td><td>-1.51</td><td>72 (50)</td><td>80</td></tr><tr><td>NGC 5272 (M3)</td><td>-1.50</td><td>177 (175)</td><td>167</td></tr><tr><td>NGC 5904 (M5)</td><td>-1.33</td><td>90 (67)</td><td>87</td></tr><tr><td>NGC 6362</td><td>-1.07</td><td>18 (18)</td><td>80</td></tr></tbody></table>

**Note.**

(Tek2K), CTIO 1 m (Y4KCam), CTIO 1.3 m (ANDICAM), CTIO LCOGT 1 m (CCD), Cerro Pachón SOAR 4.1 m (SOI), La Silla NTT 3.6 m (EMMI TK2048EB), La Silla ESO/MPI 2.2 m (WFI), LCO Warsaw 1.3 m (8k-MOSAIC), Maunakea CFHT 3.6 m (CFH12K), ORM La Palma JKT 1 m (EEV7), ORM La Palma INT 2.5 m (WFC), SAAO LCOGT 1 m (CCD). For further information about the summary of observing runs for each cluster, the bands observed, and the imagers and telescopes used, we direct the reader to the PBS database (previously mentioned).

The spectroscopic metallicities of these GCs (second column in Table 2) are listed in the scale of C09, and the third column of Table 2 shows the total number of fundamental-mode RRLs in each GC according to Clement's catalog<fnref n="15" /> (Clement et al. 2001). In parentheses, we list the number of RRab actually available in the PBS photometry. In addition, the average number of epochs per V-band light curve is shown in the fourth column.

## 3. Calibration of Period–Fourier–Metallicity Relation

Before folding the ASAS-SN and WISE time series into phased light curves, we ensured that they had at least 30 available epochs, had a signal-to-noise ratio > 5, and no quality

<footnote n="14">¹⁴ https://www.canfar.net/storage/list/STETSON/homogeneous/Latest_photometry_for_targets_with_at_least_BVI. For additional information on the PBS database, readers may also contact Peter B. Stetson directly at peter.stetson@nrc-cnrc.gc.ca.</footnote>

<footnote>a The number of RRab listed in each GC is according to Clement's catalog (Clement et al. 2001), where the actual number of RRab available with the PBS photometry is in parentheses.</footnote>

<footnote n="15">¹⁵ http://www.astro.utoronto.ca/~cclement/read.html</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Multiple-band light curves for a typical short- and long-period star. The V (left), W1 (middle), and W2 bands (right) are shown for the star DM Cyg (period = 0.419866 days) in the top row and NSVS 13688631 (period = 0.737646 days) in the bottom row. The Fourier fit (black solid line) to the GLOESS light curve (see text) is plotted on top of the phased data (red). The average magnitude with its associated error is shown as horizontal black dotted lines. Points automatically rejected by the GLOESS fitting procedure have been removed.</figcaption></figure>

or contamination flags in their original archive. Most time series passed this check.

Figure 3 shows the V, W1, and W2 light curve typical for two of our stars with good-quality photometry and optimal phase sampling. While the V band has the characteristic sawtooth shape of RRab variables, with a sharp minimum and maximum connected by a quick rise, the infrared curves are more symmetric and sinusoidal, with a broad maximum, still retaining a sharp minimum. As mentioned above, this is a consequence of optical light curves being more sensitive to changes in the star's effective temperature, while the variations in radius are the determinant factor for the infrared light curves. We take advantage of these wavelength-dependent properties by deriving separate calibrations of the period–$\phi_{31}$–[Fe/H] in the optical and infrared, allowing us to probe different aspects of RRLs' stellar atmospheric physics.

This section describes the procedure we followed to fit and validate our relations: (a) we first refine the nominal period of each variable found in the literature by taking advantage of the large time coverage of our ASAS-SN and WISE time series, (b) we then smooth the phased light curves to prepare them for efficient Fourier analysis and perform quality control on the light curves to remove noisy data, (c) we calculate the Fourier decomposition and identify the parameters that are best correlated with metallicity, and (d) we finally fit the period–$\phi_{31}$–[Fe/H] relations.

### 3.1. Period Determination

The large temporal baseline of the ASAS-SN and WISE surveys makes our time series very sensitive to even small errors in the nominal period of the stars. Even inaccuracies of $10^{-6}$ days, over the 10 yr span of WISE and its extended mission, result in a ~0.7% phase shift when the photometric time series of a 0.55-day period of our typical RRab star is phased, readily detectable in our data sets. To avoid this issue, we have rederived the periods of all our stars on the basis of their V- and W1-band<fnref n="16" /> photometry, using the Lomb–Scargle method (Lomb 1976; Scargle 1982) within a search window from 0.2 to 1.5 days.

Typical rederived periods differ from the ones in the literature by $10^{-5}$ days or less. Light curves with period discrepancy larger than $10^{-5}$ days between the V and W1 data sets, or with the nominal period from the literature, were inspected manually for further quality checks and possible rejection. Approximately 5% of derived periods fell in this category and required manual inspection. The inability to find a reliable period was also used as a rejection criterion for an individual light curve. Note that this quality assurance process causes the removal of some high-amplitude Blazhko stars (Blažko 1907), due to their naturally larger dispersion in the light curve causing a high “false-alarm”

<footnote n="16">¹⁶ We did not use the W2 band for period determination owing to its larger photometric error, often leading to less accurate periods.</footnote>
probability in the period determination returned by the Lomb–Scargle routine. This removal is intended, as the photometric modulation of Blazhko stars can affect the Fourier decomposition of their light curves, complicating the dependence of [Fe/H] from period and $\phi_{31}$ (see, e.g., Skarka 2014).

## 3.2. GLOESS Light-curve Smoothing and Quality Control

Rather than directly calculate the Fourier decomposition of the phased light curve, we elected to smooth it first to better deal with the uneven sampling of the observing epochs, and primarily to allow the removal of data points with excessive noise. The smoothing was performed using a Gaussian locally weighted regression smoothing algorithm (GLOESS; Persson et al. 2004). This method places the phased light curve on an interpolated grid, to which a second-degree polynomial is locally fit to the photometric data. Fitting weights depend on both the photometric error and the Gaussian distance of the phased photometry from the interpolation point. The procedure was repeated multiple times in order to apply an iterative sigma clipping rejection scheme, designed for noisy data point removal. Only time series with a minimum of 30 valid photometric data points after sigma clipping were retained for processing, with the rest excluded from further consideration. A full description of our implementation of the GLOESS method is available in Neeley et al. (2015).

Figure 3 shows example light curves of the phased data (red points) for both a typical short- and long-period RRL star (DM Cyg and NSVS 13688631, respectively). Data points automatically rejected by the GLOESS fitting procedure are removed. From the GLOESS light curve, we have measured the mean magnitude (calculated as the mean of the smoothed light curve in flux units, converted back to magnitude), the amplitude, and the epoch of maximum (separately in each band). The mean magnitude uncertainty was calculated as the sum in quadrature of the photometric uncertainty and the uncertainty in the fit (see Neeley et al. 2015 for details). The error associated with the amplitude was instead defined as the standard deviation of the data points’ residual with respect to the smoothed light curve. Stars with an amplitude equal to less than three times the amplitude error were found to correspond to excessively noisy light curves and were excluded from further analysis. Note that this process also excluded large-modulation Blazhko stars, which would appear as noisy light curves with large uncertainty in amplitude.

## 3.3. Fourier Decomposition

Fourier decomposition is a widely adopted tool to quantify the shape of a light curve, as the lower-order terms in the expansion are usually sufficient to fully characterize its shape (see, e.g., Simon & Lee 1981). As justified above, for this work we elected to perform Fourier expansion of the evenly spaced smoothed light curves, rather than the individual data points, in the form

$$m(\Phi) = A_0 + \sum_{i=1}^{n} A_i \sin[2\pi i(\Phi + \Phi_0) + \phi_i], \quad (1)$$

where $m(\Phi)$ is the observed magnitude for either the ASAS-SN or WISE bands, $A_0$ is the mean magnitude, $n$ is the order of the expansion, $\Phi$ is the phase from the GLOESS light curve varying from 0 to 1, $\Phi_0$ is the phase that corresponds to the time of maximum light $T_0$, and $A_i$ and $\phi_i$ are the $i$th-order Fourier amplitude and phase coefficients, respectively.

We determined the Fourier coefficients with a weighted least-squares fit of Equation (1) to the smoothed GLOESS light curve and its locally calculated error. The locally calculated error was defined as the local photometric scatter around the smoothed light curve, estimated as the weighted standard deviation of the residual with the data convolved with the GLOESS smoothing kernel. We found that a fifth-order ($n=5$) Fourier expansion was sufficient to reproduce the shape of the light curve in each band. Examples are shown in Figure 3 for two typical RRLs. The Fourier decomposition fit (black solid line) is plotted on top of the actual photometric data in red.

Simon & Lee (1981) first showed that certain combinations of Fourier coefficients were directly related to some physical parameters of pulsating stars. These coefficients are typically defined either as linear combinations of Fourier phases,

$$\phi_{ij} = j \cdot \phi_i - i \cdot \phi_j, \quad (2)$$

where $\phi_{ij}$ is cyclic in nature and ranges from 0 to $2\pi$, or as ratios of the Fourier amplitudes,

$$R_{ij} = \frac{A_i}{A_j}. \quad (3)$$

In Appendix A, we discuss the correlation between various combinations of Fourier parameters among themselves, as well as with period, and how they help discriminate sources with different metallicities. Our analysis confirms the conclusions of early studies such as JK96, suggesting that a relation between period and $\phi_{31}$ is a good indicator of metallicity in the V band. We found this to be true also in the infrared. Uncertainties in the calculated $\phi_{31}$ values are several orders of magnitude less than our best uncertainties in metallicity and are therefore deemed negligible throughout the rest of this work.

Figure 4 shows that the W1 and W2 bands produce indistinguishable values of $\phi_{31}$ owing to the light curves in these bands being nearly identical. This is demonstrated quantitatively by Equation (4), which shows that the best-fit slope between $\phi_{31(W1)}$ and $\phi_{31(W2)}$ is within the errors close to unity, with a dispersion of 0.30:

$$\phi_{31(W2)} = (1.011 \pm 0.023) \cdot \phi_{31(W1)} + (0.066 \pm 0.043). \quad (4)$$

We take advantage of this strong correlation by averaging, whenever possible, the $\phi_{31}$ parameters calculated for the two bands, falling back on $\phi_{31(W1)}$ or $\phi_{31(W2)}$ when only one is available. This essentially doubles the signal of the WISE light curves by combining the data of independent measurements taken at different wavelengths.

## 3.4. Period–Fourier–[Fe/H] Fitting

The metallicity of each star in either our ASAS-SN or WISE data sets can be effectively represented by a plane in the period, $\phi_{31}$, and [Fe/H] space, as is demonstrated in Appendix B by performing principal component analysis (PCA) on each data set. To determine the orientation of this plane, we adopted the orthogonal distance regression (ODR) routine part of the SciPy package,<fnref n="17" /> which utilizes a modified trust-region Levenberg-Marquardt-type algorithm (Boggs & Rogers 1990) to estimate the best fitting-parameters. We chose ODR because it can be used similarly to PCA, with both techniques minimizing the perpendicular distance to the fit with no differentiation between dependent and independent variables. This allows ODR to

<footnote n="17">¹⁷ https://docs.scipy.org/doc/scipy/reference/odr.html</footnote>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Comparison between the W1 and W2 $\phi_{31}$ parameter. The dashed orange line shows the rms between the two different $\phi_{31}$ values, or the typical scatter around the ideal one-to-one relation. The data points shown correspond to the subset of the WISE data set with both W1 and W2 $\phi_{31}$ measurements available.</figcaption></figure>

produce an optimal fit despite the correlation between variables we found in Appendix B (especially the well-known correlation between period and metallicity; see, e.g., F19). These correlations, in the case of ordinary least-squares or similar fitting algorithms, would result in a trend of the fitted metallicity residuals with respect to the other variables.

We performed the ODR fit on both our ASAS-SN and WISE calibration sets using the equation

$$[\text{Fe/H}] = a + b \cdot (P - P_0) + c \cdot (\phi_{31} - \phi_{31_0}), \quad (5)$$

where $P_0$ and $\phi_{31_0}$ are pivot offsets necessary to add an extra element of robustness in the fitting procedure and reduce the fitting parameter uncertainties. We chose the pivot offsets to be near the mean period and $\phi_{31}$ value of each data set, equal to $P_0 = 0.58$ days for both data sets, and $\phi_{31_0} = 5.25$ and $1.90$ rad for ASAS-SN and WISE, respectively. To ensure an accurate fit, two rounds of fitting were performed with an intermediate $4\sigma$ clipping between the fitted and calibration [Fe/H] (removing $\lesssim 1\%$ of the stars in each calibration data set).

Our best-fit period-$\phi_{31}$-[Fe/H] relation based on the ASAS-SN V-band light curves is

$$[\text{Fe/H}] = (-1.22 \pm 0.01) + (-7.60 \pm 0.24) \cdot (P - 0.58) \\ + (1.42 \pm 0.05) \cdot (\phi_{31} - 5.25); \quad \text{rms} = 0.41, \qquad (6)$$

where, due to the $2\pi$ periodic ambiguity in the $\phi_{31}$ coefficient, some $\phi_{31}$ values required adding $2\pi$ to their phase in order to lie closer to the mean $\phi_{31}$ value, as suggested by JK96.

The period-$\phi_{31}$-[Fe/H] relation, based on the average $\phi_{31}$ parameters of the WISE light curves (averaged when possible between the W1 and W2 bands, as described in Section 3.3), is instead

$$[\text{Fe/H}] = (-1.47 \pm 0.02) + (-8.33 \pm 0.34) \cdot (P - 0.58) \\ + (0.92 \pm 0.05) \cdot (\phi_{31} - 1.90); \quad \text{rms} = 0.50. \qquad (7)$$

Parameter uncertainties have been checked with bootstrap resampling and are consistent with those via ODR. The errors from bootstrap resampling are those listed in the above relations.

The rms values of the two relations (0.41 and 0.50 dex for the ASAS-SN and WISE samples, respectively) are similar, showing that indeed accurate photometric metallicities can be obtained from infrared light curves. The rms values are also comparable to the dispersion that we found with a nonparametric regression scheme, based on the k-NN method, of the same data (0.33 and 0.40 dex, respectively; see Appendix C). This shows that our ODR fits provide an accurate description of the dependence of [Fe/H] on period and $\phi_{31}$, with uncertainty only slightly larger than the data's intrinsic scatter. Table 3 shows the derived photometric properties for both the V-band and mid-IR (WISE) calibration data sets, including the period, $\phi_{31}$ value, and photometric metallicity in each band. Following the light-curve fitting, quality control process, and the final plane fit described above, we were left with 1980 variables with good-quality ASAS-SN light curves and 1083 variables with good WISE (in at least one of the two W1 and W2 bands) light curves.

Figures 5 and 6 (for the ASAS-SN and WISE samples, respectively) demonstrate in graphical form the ability of our fits to provide photometric metallicities in agreement with the spectroscopic values in the HR+ΔS calibration sample. The top left panel shows the distribution of the sources in the $\phi_{31}$ versus period plane, superimposed with the relations in Equations (6) and (7) calculated for fixed values of [Fe/H]. Both data and fit lines are color-binned by metallicity. Note the good agreement between Note the good agreement between the distribution of the spectroscopic metallicities of individual sources with the locus of the photometric relations corresponding to the same metallicity bin. The histogram in the bottom left panel shows how the distribution of the photometric [Fe/H] from our fits closely reproduces the metallicity distribution, from spectroscopy, of the calibration sample. The two panels on the right instead show the distribution of the residuals between the photometric and spectroscopic metallicities: the top panel confirms that our choice of fitting the data with an ODR method indeed avoids residual trends, while the bottom panel shows a close-to-Gaussian distribution of the residuals.

# 4. Discussion

## 4.1. Optical versus Infrared Relations

The analysis in Section 3.4 shows that optical and infrared period-$\phi_{31}$-[Fe/H] relations provide photometric metallicities of comparable accuracy. However, we still need to test whether Equations (6) and (7) provide consistent values of [Fe/H] for individual stars. Figure 7 shows that this is indeed the case: we find an excellent agreement between the photometric metallicities derived at the two wavelength ranges for the joint sample on a per-star basis. The dispersion between the two data sets, of 0.44 dex, is comparable to the rms of the individual fits, as well as to the dispersion of the residuals from the k-NN method in Appendix C.

This tight correlation between optical and infrared photometric metallicities, with no apparent trends, led us to consider whether averaging the two photometric metallicities would yield a value significantly closer to the spectroscopic metallicity from the HR+ΔS sample. Furthermore, besides the obvious advantage of increasing the statistics of the photometric measurements, by combining optical and infrared data we can probe the effects that different aspects of stellar atmospheres have on the RRL light curves acquired in these two separate wavelength ranges (the optical light curves are
<figure><img src="image_7.png" /><figcaption>**Figure 5.** ASAS-SN V-band period-$\phi_{31}$-[Fe/H] fit. Top left: period vs. $\phi_{31}$ plane with stars binned and color-coded based on their spectroscopic metallicity. Solid lines represent the best-fit relation calculated for the average metallicity in each bin, from metal-rich in the bottom left (red) to metal-poor in the top right (blue). Bottom left: comparison of the spectroscopic calibration [Fe/H] with the best-fit photometric metallicity. Top right: spectroscopic [Fe/H] plotted as a function of the residuals of the photometric and spectroscopic metallicity. The dashed vertical lines represent the rms error between the photometric and spectroscopic metallicity. Bottom right: histogram of the difference between the photometric and spectroscopic [Fe/H]. Again, the vertical lines show the rms dispersion of the data.</figcaption></figure>

<table><caption>**Table 3**

Derived Photometric Properties</caption><thead><tr><th>Gaia ID<br>(DR3)</th><th>Period<sup>a</sup><br>(days)</th><th>φ<sub>31</sub> (V)<br>(rad)</th><th>[Fe/H]<sub>V</sub><br>(dex)</th><th>φ<sub>31</sub> (W)<br>(rad)</th><th>[Fe/H]<sub>W</sub><br>(dex)</th></tr></thead><tbody><tr><td>507222753405440</td><td>0.6088826</td><td></td><td></td><td>1.65842</td><td>-1.93</td></tr><tr><td>4235220006525184</td><td>0.4664141</td><td></td><td></td><td>1.18099</td><td>-1.18</td></tr><tr><td>5355313117668352</td><td>0.4945746</td><td>4.90759</td><td>-1.06</td><td></td><td></td></tr><tr><td>14233869512030080</td><td>0.6670034</td><td>4.96473</td><td>-2.29</td><td>1.62441</td><td>-2.44</td></tr><tr><td>15489408711727488</td><td>0.6511669</td><td>5.19695</td><td>-1.84</td><td>1.80280</td><td>-2.15</td></tr><tr><td>15891245851805568</td><td>0.5738989</td><td>5.04774</td><td>-1.46</td><td>1.51379</td><td>-1.77</td></tr><tr><td>18268974106572416</td><td>0.5541744</td><td>4.75948</td><td>-1.72</td><td>1.26200</td><td>-1.84</td></tr><tr><td>19606736160298112</td><td>0.5465583</td><td>5.01736</td><td>-1.30</td><td></td><td></td></tr><tr><td>20161096179157248</td><td>0.6042955</td><td></td><td></td><td>2.31260</td><td>-1.29</td></tr><tr><td>20357148550791168</td><td>0.4853687</td><td>4.31030</td><td>-1.84</td><td></td><td></td></tr></tbody></table>

**Note.**

<footnote><sup>a</sup> When both V-band and mid-IR (WISE) data are present, the period included was calculated from ASAS-SN (V-band) data, as the period is usually more accurate owing to the higher amplitude and steeper light curve.</footnote>

<footnote>(This table is available in its entirety in machine-readable form.)</footnote>
<figure><img src="image_8.png" /><figcaption>Figure 6. Same as Figure 5, but for the WISE light curves.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Comparison of the photometric [Fe/H] derived from ASAS-SN and WISE samples. The left panel directly compares the photometric metallicities derived for each star in the joint sample, while the right panel shows a histogram of the differences between these photometric metallicities. The dashed lines show the rms between the two data sets [Fe/H].</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Figure 8.** Left: comparison of the distribution of the spectroscopic calibration metallicity (gray histogram) and the fitted photometric metallicity averaged from both the ASAS-SN and WISE P − φ₃₁ − [Fe/H] fit (black histogram). Right: histogram of the difference between the average photometric and spectroscopic metallicities. The vertical orange lines represent the rms dispersion between the two metallicity values.</figcaption></figure>

<figure><img src="image_11.png" /></figure>

dominated by temperature, while infrared emission follows
more closely radius variations).

Based on these arguments, we average together the
photometric metallicities derived from the optical and infrared
relationships for each star in the joint sample. Figure 8 shows
that the average photometric metallicity (black histogram in left
panel) closely mirrors the spectroscopic metallicity values from
the HR+ΔS calibration sample (gray histogram) for the entire
metallicity range. The right panel shows that the residuals
between the average photometric [Fe/H] and the spectroscopic
metallicities fall fairly symmetrically around zero.

The rms dispersion between the two sets of values is ~0.37 dex —smaller than the individual error in the optical (~0.41 dex) or infrared (~0.50 dex) sample alone. This rms value is near the dispersion we found with the k-NN method in Appendix C (0.33 and 0.40 dex for the ASAS-SN and WISE data sets, respectively), which shows that combining these data sets allows one to approach the local scatter in the individual data sets.

## 4.2. Comparison with Globular Cluster Metallicity

In order to test the V-band relation obtained in Section 3.3 on
an independent sample, we selected a list of eight Galactic GCs
homogeneously spread between [Fe/H] = −1.1 and −2.3 dex.
The selected GCs are described in Section 2.4 and listed in
Table 2. For each cluster, we chose the RRab stars with the best-
sampled light curves in the V band, avoiding those RRab stars that
suffer the Blazhko effect because of the modulation of the
amplitude and the shape of their light curves. A Fourier
decomposition was performed on each light curve to obtain their
φ₃₁ parameter. The period-φ₃₁-[Fe/H] relation (Equation (6))
was then applied to estimate the metallicity of each cluster star.

Figure 9 shows the spectroscopic [Fe/H] versus the mean
photometric [Fe/H] values calculated for each GC with our
relation. For consistency with the field RRLs described in
Section 2.1, the spectroscopic metallicities of these GCs (second
column in Table 2) from C09 were converted into the scale of this
paper with the addition of an offset of 0.08 dex, as noted in C21
and described in Section 2.1. The error bars correspond to the

<figure><img src="image_12.png" /><figcaption>**Figure 9.** [Fe/H]$_{spec}$ vs. [Fe/H]$_{phot}$ obtained by using Equation (6) on a sample of GCs. The solid orange line is the one-to-one relation, while the dashed orange lines show the standard deviation (±0.09 dex). Error bars correspond to the uncertainties in spectroscopic metallicity and the statistical error in the photometric metallicity.</figcaption></figure>

spectroscopic and photometric metallicity uncertainties. The former
comes from C09, while the latter was assessed as the standard error
of the mean for star photometric metallicities in each cluster. The
good performance obtained using the period-φ₃₁-[Fe/H] relation
derived in this work is clearly noticeable in Figure 9. There are no
signals of possible systematic effects, and the predicted [Fe/H]
using the period-φ₃₁-[Fe/H] relation obtained in this work is
within ± 0.09 dex of the spectroscopic [Fe/H]. In fact, this
accuracy is similar to the metallicity uncertainties measured from
the HR spectroscopy on individual stars.

4.3. *Comparison with Other Relations*

In this section, we compare our optical period-φ₃₁-[Fe/H]
relation (Equation (6)) with previous relations found in the
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Comparison of the [Fe/H] derived with our V-band $P - \phi_{31}$-[Fe/H] relation with that derived using the relationships of JK96 (top left), N13 (top right), MV16 (bottom left), and IB21 (bottom right), applied to the period and $\phi_{31}$ values derived in the ASAS-SN calibration data set. The range of calibration metallicities used in each work to derive their respective relations is noted in each panel with the two red dashed-dotted lines. The blue dashed lines around the 1:1 secant (solid line) represent the rms dispersion between this work's and the literature [Fe/H] abundance for the stars within the red dashed-dotted calibration range. Note that in the top right panel, outside the calibration region of N13, the quadratic trend deviates significantly from this work for low metallicities (beyond the bounds of this panel).</figcaption></figure>

literature for similar wavelength ranges. In particular, we focus on the relations found in JK96, N13, MV16, and IB21. Note that a similar comparison for our mid-infrared relation (Equation (7)) is not possible, since we could not find any previous work studying the relation between metallicity and Fourier parameters at wavelengths longer than the *I* band (Smolec 2005). The results are shown in Figure 10, which plots the [Fe/H] abundances calculated with the periods and $\phi_{31}$ parameters from the ASAS-SN sample, using our relation versus the relations published in the literature.

The top left panel shows the comparison with the formula found in JK96 (their Equation (3)). Their relation was derived using a total of 81 field RRab stars with V-band photometry from heterogeneous observations at various sites. Due to the lack of phase coverage or excessive noise, a direct Fourier fit in JK96 was often not available, and individual polynomial fits or small parabolas were fit to light-curve segments as needed. Since JK96 adopted metallicities based on the high-dispersion spectroscopy scale of Jurcsik (1995), for consistency we first

converted the metallicities derived with their relation to the C09 scale, using the relation from Kapakos et al. (2011): $[Fe/H]_{C09} = 1.001[Fe/H]_{JK96} - 0.112$. A scale offset of 0.08 was then added to transform the C09 scale to that adopted by this work’s HR+ΔS calibration sample (see Section 2.1). We found agreement within the calibration range of JK96 (red horizontal lines), with an rms value of 0.17 dex, which is smaller than the nominal uncertainty in the calculated [Fe/H] abundances, and only a mildly apparent trend.

N13 introduced a quadratic period-$\phi_{31}$-[Fe/H] relation (their Equation (2)), calibrated using stars observed in the Kepler photometric band during the first 970 days of the Kepler mission. Since Kepler, during its nominal mission, surveyed the same region of the sky nearly continuously, each star resulted in having ~350,000 data points spread over ~2500 pulsation cycles, yielding the best-sampled light curves (and Fourier decomposition) currently available for any sample of RRLs. However, due to the fixed field of view and shallow depth of Kepler’s field, their calibration data set only had a total
of 26 RRab stars with accurate metallicity measurements, nine of which were Blazhko. In order to use the relation of N13, derived in the Kepler photometric system (Kp), with our data, we had to convert the V-band $\phi_{31}$ value to the Kepler system by adding the systematic offset derived by Nemec et al. (2011):

$$ \phi_{31}(V) = \phi_{31}(Kp) - (0.151 \pm 0.026). \quad (8) $$

Furthermore, we converted the [Fe/H] abundances provided by N13 from their adopted metallicity scale C09 to the scale in this work, using the previously noted scale offset of 0.08 dex. The results are presented in the top right panel of Figure 10 and show an excellent agreement for higher metallicities in the N13 calibration range ($-1.5 \text{ dex} \lesssim [\text{Fe/H}] \lesssim 0.03 \text{ dex}$). For lower metallicities, however, the two relations rapidly diverge ($\text{rms} \approx 0.30 \text{ dex}$ if calculated over their entire calibration range), caused by the higher-order term in N13, and possibly due to the scarcity of calibrators in their samples for low [Fe/H] (they only have one RRL with $[\text{Fe/H}] \lesssim -2.0 \text{ dex}$).

Next, we compare our relationship with the one derived by MV16 (end of Section 2 in their paper), calibrated using a sample of 381 RRLs in GCs binned by period, with the addition of eight field RRLs chosen to extend the metallicity range of the sample. The C09 metallicity scale used by MV16 has also been converted to that of this work to allow a comparison with our sample. The result is presented in the bottom right panel of Figure 10, showing a general agreement with our fit within the calibration range of MV16, albeit with a small slope with respect to our relation. Note that we set the lower limit of the calibration range plotted to $[\text{Fe/H}] \approx -0.99$ since there was a discontinuity in the data beyond which only had two variable stars.

Finally, IB21 introduced a G-band period-$\phi_{31}$-[Fe/H] relation (their Equation (3)), based on Gaia DR2 light curves. The relation was calibrated with a sample of 84 stars with known spectroscopic metallicity. To allow comparison with our metallicities, we first had to convert the V-band $\phi_{31}$ value of our data set to the G-band system, using the relation derived by Clementini et al. (2016):

$$ \phi_{31}(G) = (0.104 \pm 0.020) + (1.000 \pm 0.008)\phi_{31}(V); \\ \text{rms} = 0.055. \quad (9) $$

Furthermore, an additional $\pi$ offset had to be subtracted from $\phi_{31}$ to set the coefficient on the same scale as IB21. Metallicity abundances from IB21 were on the scale of Zinn & West (1984, ZW84) and were thus converted using the following equation from C09 with a subsequent scale conversion to that of this work: $[\text{Fe/H}]_{C09} = 1.105[\text{Fe/H}]_{ZW84} + 0.160$.

The results are shown in the bottom right panel of Figure 10 and exhibit a clear trend between the two relations for the entire range of metallicity. It should be noted that although the stars used for the calibration of the IB21 relation are in the range of $-2.53 \lesssim [\text{Fe/H}] \lesssim 0.33$, the fit lacks a significant number of calibrators at the low- and high-metallicity ends. In comparison to this work, their relation tends to overestimate the metallicity at the metal-poor end and underestimate the metallicity at the metal-rich end.

### 4.4. Comparison with High-resolution Spectroscopic Metallicities

To assess the reliability of the period-$\phi_{31}$-metallicity relations, both derived in this paper and from the literature, in providing [Fe/H] abundances, we compare their predictions

with metallicities measured from HR spectroscopy in the HR +ΔS sample. The results are shown in Figure 11, plotting the difference between spectroscopic and photometric [Fe/H] abundance for each of the stars in C21 with metallicity from an HR spectrum. The rms scatter of each relation, calculated with respect to zero residuals over the entire spectroscopic metallicity range, is indicated for each photometric relation.

The figure shows that all relations have similar residuals ($\text{rms} \approx 0.3-0.4 \text{ dex}$). The scatter in our relations (top row), in particular, is consistent with the expected uncertainties of our photometric metallicities. There is marginal evidence for a systematic shift downward in panels of Figure 11 that we attribute to the minor differences between the ΔS and HR metallicities as shown in Figure 9 of C21. All relations from the literature show a trend of the residuals with HR spectroscopic metallicity due to their tendency to overestimate [Fe/H] abundances at the low metallicity range and underestimate high metallicities, with the exception of N13 showing a significantly larger scatter ($\text{rms} \approx 0.57 \text{ dex}$). This apparent trend may result from an inadequate sampling of the metal-rich and metal-poor tails or possibly non-LTE effects in each of the calibration data set's metallicities that are not fully taken into account. More investigation is needed in order to clarify the mechanisms behind this trend. This trend is largest in the IB21, which underestimates the metallicity of all stars with $[\text{Fe/H}] \lesssim -2.0$ and overestimates the metallicity for $[\text{Fe/H}] \gtrsim 0.7 \text{ dex}$, which was already partially noted by IB21 within their calibration range.

It is worth noting that the rms values presented in Figure 11 are generally larger than those quoted inside these respective works of literature. The rms values quoted in other papers are a measure of their fit's dispersion with respect to their calibration data set and should not be construed as a measurement of accuracy between two relations. A comparison between works needs to be made on an identical validation data set (as was done in this section). We consider the larger dispersion in our relations as denoting a more realistic representation of the intrinsic scatter found in the period-$\phi_{31}$-[Fe/H] relation, which is not apparent in smaller calibration samples, usually covering a smaller metallicity range. Our calibration sample is $76 \times$ larger than N13 (26 RRab) and $24 \times$ larger than JK96 (81 RRab). We could decrease the dispersion of this work's relations by artificially cutting the tails of the metallicity distribution (very metal-poor/metal-rich). However, this approach would introduce severe systematic drifts in the metallicity estimates, i.e., the same problem affecting other current calibrations, which is what we have strived to resolve in this work.

## 5. Conclusions

In this work, we have calibrated a new relation to provide RRab photometric metallicities based on the $\phi_{31}$ Fourier parameter of optical and, for the first time, infrared light curves. Our relations are based on a set of homogeneous spectroscopic [Fe/H] abundances derived by using the ΔS method of C21 and validated with a sample of [Fe/H] metallicities from HR spectra. The photometric time series of our calibration stars were extracted from the ASAS-SN (V band) and WISE (NEOWISE extension, W1 and W2 bands) surveys, providing well-sampled light curves that allow for reliable Fourier expansions and accurate determination of the RRL pulsation periods (better than $10^{-6}$ days).
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Comparison of the difference between the HR spectroscopic [Fe/H] and that derived with the photometric relations discussed in Section 4.3. The black dashed lines around the zero residual solid line represent the rms dispersion between the respective work's [Fe/H] abundance and the spectroscopic metallicities from HR data. The range of calibration metallicities used in each work to derive their respective relations is noted in each panel with the two red dashed-dotted lines. Note that the mid-infrared comparison (top right panel) can only be shown for a smaller subset of HR stars for which the WISE light curve is available.</figcaption></figure>

Comparisons with other optical photometric metallicity relations available in the literature show that our relation provides reliable [Fe/H] abundances without noticeable trends over the broadest metallicity range (from [Fe/H] ≲ −3 to solar). We have also shown that there is an intrinsic scatter in the period−*φ*<sub>31</sub>−[Fe/H] plane that becomes apparent with large calibration data sets covering a broad metallicity range. Our V-band relation (Equation (6)) is consistent within ≃0.40 dex

with [Fe/H] abundances from HR spectroscopy over the entire metallicity range. This relation allows for a quick determination of reliable metallicities for large photometric samples of RRLs that will be observed in upcoming optical wide-area time-domain surveys. It offers a good compromise between efficiency and accuracy, with uncertainties only 2 or 3 times larger than metallicities from labor-intensive line-fitting HR spectroscopy. Our relation can also be applied to existing data
sets of well-sampled optical light curves obtained by Kepler and Gaia, by converting the $\phi_{31}$ parameters from their respective $K_p$ and G bands, using transformations available in the literature. Finally, we show that our V-band relation, applied to ensembles of RRLs in Galactic GCs, provides estimates of the clusters' metallicity with accuracy comparable to HR spectroscopy measurements (within $\pm 0.09$ dex).

In addition, we have obtained for the first time a mid-infrared period-$\phi_{31}$-[Fe/H] relation using WISE W1 and W2 bands (Equation (7)). Despite having a smaller number of RRLs to derive this relation, the rms obtained was similar to that we obtained in the optical. While our infrared relations are slightly less accurate (rms $\simeq 0.50$ dex), they can still be used to obtain statistically representative metallicities for large ensembles of RRLs. This will be crucial with the advent of sensitive telescopes in the mid-infrared (such as JWST), which will allow observations of extragalactic RRLs across the Local Group of galaxies, for which spectral observations will not be feasible. Further new large optical (WEAVE, 18 4MOST$^{19}$) and near-infrared (MOONS$^{20}$ at VLT) spectroscopic surveys will provide the unique opportunity to improve the current calibration of the period-Fourier-metallicity relations using both medium-resolution and HR spectra. However, by combining optical and infrared photometric metallicities, we have shown that it is possible to further improve the measurements' reliability.

Having a larger sample size of RRLs with known metallicity will enable a more in-depth study of the metallicity distributions of RRLs (tracers of older stellar populations) in our local neighborhood.

This publication makes use of data products from WISE, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory (JPL)/California Institute of Technology (Caltech), funded by the National Aeronautics and Space Administration (NASA), and from NEOWISE, which is a JPL/Caltech project funded by NASA's Planetary Science Division.

This publication also makes use of data products from the ASAS-SN project, which has their telescopes hosted by Las Cumbres Observatory. ASAS-SN is supported by the Gordon and Betty Moore Foundation through grant GBMF5490 and the NSF

by grants AST-151592 and AST-1908570. Development of ASAS-SN has been supported by Peking University, Mt. Cuba Astronomical Foundation, Ohio State University Center for Cosmology and AstroParticle Physics, the Chinese Academy of Sciences South America Center for Astronomy (CASSACA), the Villum Foundation, and George Skestos.

This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement.

M. Marengo and J.P.M. are supported by the National Science Foundation under grant No. AST-1714534. M. Monelli has been supported by the Spanish Ministry of Economy and Competitiveness (MINECO) under grant AYA2017-89076-P.

**Facilities:** WISE, ASAS-SN, Gaia.

**Software:** Astropy (Astropy Collaboration et al. 2013), SciPy (Virtanen et al. 2020), scikit-learn (Pedregosa et al. 2011).

## Appendix A: Fourier Parameters

To visualize whether the shape of an infrared RRL light curve retains information about the star's metallicity (as is the case in the optical), we have analyzed several low-order Fourier parameters, plotted one against the other, or versus the period. Figures 12 and 13 plot individual RRLs, color-coded on their [Fe/H], on the basis of $R_{21}$, $R_{31}$, $\phi_{21}$, $\phi_{31}$, and period for the V and W1 bands, respectively. In all panels where a Fourier parameter is plotted as a function of period, we can readily see a gradient in the [Fe/H] distribution of the stars.

The stronger separation is apparent in the $\phi_{21}$ and $\phi_{31}$ versus period plots, as already discussed in JK96 for V-band light curves. We found that the same is true in the WISE W1 and W2 bands. In accordance with previous literature, and based on the observation that the $\phi_{31}$ parameter provides a tighter relation with period and metallicity at both optical and infrared wavelengths, we have adopted $\phi_{31}$ as the parameter of choice for the analysis presented in Section 3.3.

<footnote n="18">18 https://ingconfluence.ing.iac.es:8444/confluence/display/WEAV/The+WEAVE+Project</footnote>

<footnote n="19">19 https://www.4most.eu/cms/facility/overview/</footnote>

<footnote n="20">20 https://www.eso.org/sci/facilities/develop/instruments/MOONS.html</footnote>
<figure><img src="image_15.png" /><figcaption>**Figure 12.** ASAS-SN V-band Fourier amplitude parameter ratios (top row) and linear phase parameter combinations (bottom row). RRLs are color-coded based on their spectroscopic calibration metallicities. In the middle and left panels, [Fe/H] generally goes from metal-rich in the left (red) to metal-poor in the right (blue).</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Figure 13. Same as Figure 12, but for our WISE sample.</figcaption></figure>

## Appendix B
### Principal Component Analysis

Previous works have successfully determined linear relations between period, $\phi_{31}$, and [Fe/H] abundance in RRab stars (except for N13, which preferred a quadratic form for the relation). In the case of uncorrelated variables, coefficients of such relations can be estimated with a simple least-squares fit of the available data. It is well known, however, that there exists a correlation between period and metallicity in RRab stars (see, e.g., Equation (14) in F19). To explore this correlation and test whether the distribution of metallicity can be effectively described with a plane in period and $\phi_{31}$ not just at optical wavelengths but also in the infrared, we have performed a PCA to these three dimensions of our data.

PCA (see, e.g., Shlens 2014) is a machine-learning algorithm often utilized for reducing the number of variables needed to describe a data set. Given a data set with *p* different variables, PCA finds a vector (principal component) in this *p*-dimensional space that can explain the most variance. Subsequent components are found orthogonal to the prior components in a direction that explains the next-highest amount of variance. The percent variance of each of these components quantifies their relative importance. In other words, PCA gives the ability to locate an *n*-dimensional hyperplane in the *p*-dimensional variable space ($n < p$, with the smallest possible *n*) that characterizes the largest amount of variance in the data.

The first step for PCA analysis is to renormalize each variable (period, $\phi_{31}$, and [Fe/H] abundance) in both samples (optical and infrared) so that they have zero mean and unitary

<table><caption>Table 4
Variance Attributed to Each Principal Component Axis</caption><thead><tr><td>Principal Component</td><td>ASAS-SN</td><td>WISE</td></tr></thead><tbody><tr><td>Axis 1</td><td>54.23%</td><td>56.51%</td></tr><tr><td>Axis 2</td><td>38.67%</td><td>35.40%</td></tr><tr><td>Axis 3</td><td>7.10%</td><td>8.09%</td></tr></tbody></table>

variance. This avoids common issues in PCA when different variables range in scale. The principal components were then found by utilizing the PCA decomposition subroutine from Python's scikit-learn package.<fnref n="21" />

Table 4 shows the PCA variance for each of the ASAS-SN and WISE data sets. It is clear that both samples can be effectively described with just the first two axes, with very little scatter associated with the third component, orthogonal to the plane described by axes 1 and 2. This shows that our data can be efficiently represented by a plane. For both data sets, the first axis is well aligned with the period of the stars, while the second axis is an almost equally weighted combination of $\phi_{31}$ and [Fe/H]. Based on this result, we feel confident that our data set can be represented by a linear fit of period, $\phi_{31}$, and [Fe/H]. We determine the exact coefficients using the ODR method described in Section 3.4. This method still provides rotated axes for the fit, but it allows us to do so while including uncertainties in the fitting variables and data.

<footnote n="21">²¹ https://scikit-learn.org/stable/modules/generated/sklearn.decomposition.PCA.html</footnote>
<figure><img src="image_17.png" /><figcaption>**Figure 14.** Difference in the [Fe/H] predicted by the *k*-NN method to the metallicity of the calibration sample for the ASAS-SN data set (left) and WISE data set (right). The rms error between the *k*-NN and calibration metallicity is listed in the upper right corners and shown by the vertical dashed lines.</figcaption></figure>

## Appendix C
### K-Nearest Neighbor Analysis

In this work, we chose to represent the relation between [Fe/H], period, and $\phi_{31}$ of RRab stars using a linear function of these three variables (see Section 3.4). An alternative approach consists of adopting nonparametric techniques that could provide a similar result without assumptions about the functional form of the relation to fit. One such technique is the *k*-Nearest Neighbor (k-NN) method (Cover & Hart 1967).

While traditionally used in machine-learning classification problems (see, e.g., Marengo & Sanchez 2009 for an application in astronomy), the *k*-NN method also serves as a powerful nonparametric regression technique to estimate the value of a given variable based on the values of the closest neighbors, in a properly defined *n*-dimensional space. In this work, the location on the period-$\phi_{31}$ plane is used as a predictor for [Fe/H]. Following Ball et al. (2007), we calculate the Euclidean distance, in the period-$\phi_{31}$ plane, of each star from every other star in the sample. We then estimate the *k*-NN metallicity of each star as the [Fe/H] weighted average of its nearest *k* neighbors, where *k* is a suitable integer number. The weights used are proportional to the inverse distance in the period-$\phi_{31}$ plane from the test star so that nearer stars contribute more to the average than the farthest sources.

The number of nearest neighbors *k* is chosen in order to balance the need of averaging out the natural scatter in the neighbors' [Fe/H], while at the same time still preserving predictive power at the sparse edges of the samples' distribution in the period-$\phi_{31}$ plane. Since our goal is to capture the global trend over the entire range of values, and not only the highest-density region, we chose *k* = 5 as the optimum number of near neighbors for both our samples.

Figure 14 shows the difference between the metallicity derived with the *k*-NN method and the spectroscopic [Fe/H] for our optical and infrared samples. The rms [Fe/H] error of the residuals is 0.33 and 0.40 dex for the ASAS-SN and WISE samples, respectively. Since the *k*-NN method, in contrast to parametric fits, does not rely on any particular functional form, we can assume that the quoted rms errors are representative of the best possible overall [Fe/H] uncertainty (averaged over the entire range of metallicity) that a specific fit could achieve with the available data.

## ORCID iDs

Joseph P. Mullen https://orcid.org/0000-0002-1650-2764

Massimo Marengo https://orcid.org/0000-0001-9910-9230

Clara E. Martínez-Vázquez https://orcid.org/0000-0002-9144-7726

Jillian R. Neeley https://orcid.org/0000-0002-8894-836X

Giuseppe Bono https://orcid.org/0000-0002-4896-8841

Massimo Dall'Ora https://orcid.org/0000-0001-8209-0449

Brian Chaboyer https://orcid.org/0000-0003-3096-4161

Frédéric Thévenin https://orcid.org/0000-0002-5032-2476

Vittorio F. Braga https://orcid.org/0000-0001-7511-2830

Juliana Crestani https://orcid.org/0000-0001-8926-3496

Michele Fabrizio https://orcid.org/0000-0001-5829-111X

Giuliana Fiorentino https://orcid.org/0000-0003-0376-6928

Christina K. Gilligan https://orcid.org/0000-0003-4510-0964

Matteo Monelli https://orcid.org/0000-0001-5292-6380

Peter B. Stetson https://orcid.org/0000-0001-6074-6830

## References

Akeson, R., Armus, L., Bachelet, E., et al. 2019, arXiv:1902.05569

Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33

Ball, N. M., Brunner, R. J., Myers, A. D., et al. 2007, ApJ, 663, 774

Blažko, S. 1907, AN, 175, 325

Boggs, P., & Rogers, J. 1990, Contemp. Math., 112, 183

Caputo, F. 1998, A&ARv, 9, 33

Carretta, E., Bragaglia, A., Gratton, R., D'Orazi, V., & Lucatello, S. 2009, A&A, 508, 695, (C09)

Catelan, M., & Smith, H. A. 2015, Pulsating Stars (Weinheim: Wiley-VCH)

Chadid, M., Sneden, C., & Preston, G. W. 2017, ApJ, 835, 187

Clement, C. M., Muzzin, A., Dufton, Q., et al. 2001, AJ, 122, 2587

Clementini, G., Carretta, E., Gratton, R., et al. 1995, AJ, 110, 2319

Clementini, G., Ripepi, V., Leccia, S., et al. 2016, A&A, 595, A133

Clementini, G., Ripepi, V., Molinaro, R., et al. 2019, A&A, 622, A60

Cover, T., & Hart, P. 1967, ITIT, 13, 21

Crestani, J., Fabrizio, M., Braga, V. F., et al. 2021, ApJ, 908, 20, (C21)

Dambis, A. K., Berdnikov, L. N., Kniazev, A. Y., et al. 2013, MNRAS, 435, 3206

Deng, L.-C., Newberg, H. J., Liu, C., et al. 2012, RAA, 12, 735

Fabrizio, M., Bono, G., Braga, V. F., et al. 2019, ApJ, 882, 169, (F19)

Fernley, J., & Barnes, T. G. 1996, A&A, 312, 957

For, B.-Q., Sneden, C., & Preston, G. W. 2011, ApJS, 197, 29

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2018, A&A, 616, A1

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, 649, A1

Gardner, J. P., Mather, J. C., Clampin, M., et al. 2006, SSRv, 123, 485
Gilligan, C. K., Chaboyer, B., Marengo, M., et al. 2021, MNRAS, 503, 4719

Govea, J., Gomez, T., Preston, G. W., & Sneden, C. 2014, ApJ, 782, 59

Holl, B., Audard, M., Nienartowicz, K., et al. 2018, A&A, 618, A30

Iorio, G., & Belokurov, V. 2021, MNRAS, 502, 5686, (IB21)

Ivezić, Ž, Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111

Jayasinghe, T., Kochanek, C. S., Stanek, K. Z., et al. 2018, MNRAS, 477, 3145

Jurcsik, J. 1995, AcA, 45, 653

Jurcsik, J., & Kovacs, G. 1996, A&A, 312, 111, (JK96)

Kapakos, E., Hatzidimitriou, D., & Soszyński, I. 2011, MNRAS, 415, 1366

Koch, D. G., Borucki, W. J., Basri, G., et al. 2010, ApJL, 713, L79

Lambert, D. L., Heath, J. E., Lemke, M., & Drake, J. 1996, ApJS, 103, 183

Law, N. M., Kulkarni, S. R., Dekany, R. G., et al. 2009, PASP, 121, 1395

Liu, G. C., Huang, Y., Zhang, H. W., et al. 2020, ApJS, 247, 68

Liu, X. W., Yuan, H. B., Huo, Z. Y., et al. 2014, in IAU Symp. 298, Setting the Scene for Gaia and LAMOST, ed. S. Feltzing et al. (Cambridge: Cambridge Univ. Press), 310

Lomb, N. R. 1976, Ap&SS, 39, 447

Magurno, D., Sneden, C., Bono, G., et al. 2019, ApJ, 881, 104

Magurno, D., Sneden, C., Braga, V. F., et al. 2018, ApJ, 864, 57

Mainzer, A., Bauer, J., Grav, T., et al. 2011, ApJ, 731, 53

Marconi, M., Coppola, G., Bono, G., et al. 2015, ApJ, 808, 50

Marengo, M., & Sanchez, M. C. 2009, AJ, 138, 63

Marrese, P. M., Marinoni, S., Fabrizio, M., & Altavilla, G. 2019, A&A, 621, A144

Martínez-Vázquez, C. E., Monelli, M., Bono, G., et al. 2016, CoKon, 105, 53, (MV16)

Muraveva, T., Delgado, H. E., Clementini, G., Sarro, L. M., & Garofalo, A. 2018, MNRAS, 481, 1195

Neeley, J. R., Marengo, M., Bono, G., et al. 2015, ApJ, 808, 11

Neeley, J. R., Marengo, M., Bono, G., et al. 2017, ApJ, 841, 84

Neeley, J. R., Marengo, M., Freedman, W. L., et al. 2019, MNRAS, 490, 4254

Nemec, J. M., Cohen, J. G., Ripepi, V., et al. 2013, ApJ, 773, 181, (N13)

Nemec, J. M., Smolec, R., Benkő, J. M., et al. 2011, MNRAS, 417, 1022

Ngeow, C.-C., Yu, P.-C., Bellm, E., et al. 2016, ApJS, 227, 30

Pancino, E., Britavskiy, N., Romano, D., et al. 2015, MNRAS, 447, 2404

Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2011, J. Mach. Learn. Res., 12, 2825, http://jmlr.org/papers/v12/pedregosa11a.html

Persson, S. E., Madore, B. F., Krzemiński, W., et al. 2004, AJ, 128, 2239

Piersimoni, A. M., Bono, G., & Ripepi, V. 2002, AJ, 124, 1528

Preston, G. W. 1959, ApJ, 130, 507

Sandage, A. 1990, ApJ, 350, 603

Scargle, J. D. 1982, ApJ, 263, 835

Shappee, B. J., Prieto, J. L., Grupe, D., et al. 2014, ApJ, 788, 48

Shlens, J. 2014, arXiv:1404.1100

Simon, N. R., & Lee, A. S. 1981, ApJ, 248, 291

Skarka, M. 2014, A&A, 562, A90

Smolec, R. 2005, AcA, 55, 59

Sneden, C., Preston, G. W., Chadid, M., & Adamów, M. 2017, ApJ, 848, 68

Udalski, A., Szymanski, M., Kaluzny, J., Kubiak, M., & Mateo, M. 1992, AcA, 42, 253

Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, NatMe, 17, 261

Walker, A. R. 1989, PASP, 101, 570

Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868

Yanny, B., Rockosi, C., Newberg, H. J., et al. 2009, AJ, 137, 4377

Zinn, R., & West, M. J. 1984, ApJS, 55, 45