Variability of Red Supergiants in M31 from the Palomar Transient Factory
Monika D. Soraisam Lars Bildsten Maria R. Drout Evan B. Bauer Marat Gilfanov Thomas Kupfer Russ R. Laher Frank Masci Thomas A. Prince Shrinivas R. Kulkarni Thomas Matheson Abhijit Saha





National Optical Astronomy Observatory, Tucson, AZ 85719, USA

Department of Physics, University of California, Santa Barbara, CA 93106, USA
Kavli Institute for Theoretical Physics, University of California, Santa Barbara, CA 93106, USA

Hubble and Carnegie-Dunlap Fellow
The Observatories of the Carnegie Institution for Science, 813 Santa Barbara St., Pasadena, CA 91101, USA
0000-0002-4791-6724]Evan B. Bauer
Department of Physics, University of California, Santa Barbara, CA 93106, USA

Max Planck Institute for Astrophysics, Karl-Schwarzschild-Str. 1, 85748 Garching, Germany
Space Research Institute, Russian Academy of Sciences, Profsoyuznaya 84/32, 117997 Moscow, Russia

Division of Physics, Mathematics, and Astronomy, California Institute of Technology, Pasadena, CA 91125, USA

Infrared Processing and Analysis Center, California Institute of Technology, Pasadena, CA 91125, USA

Infrared Processing and Analysis Center, California Institute of Technology, Pasadena, CA 91125, USA

Division of Physics, Mathematics, and Astronomy, California Institute of Technology, Pasadena, CA 91125, USA

Division of Physics, Mathematics, and Astronomy, California Institute of Technology, Pasadena, CA 91125, USA

National Optical Astronomy Observatory, Tucson, AZ 85719, USA

National Optical Astronomy Observatory, Tucson, AZ 85719, USA
Monika Soraisam
msoraisam@noao.edu




Most massive stars end their lives as Red Supergiants (RSGs), a short-lived evolution phase when they are known to pulsate with varying amplitudes. The RSG period-luminosity (PL) relation has been measured in the Milky Way, the Magellanic Clouds and M33 for about 120 stars in total. Using over 1500 epochs of R band monitoring from the Palomar Transient Factory (PTF) survey over a five-year period, we study the variability of 255 spectroscopically cataloged RSGs in M31. We find that all RGSs brighter than M_K≈ -10 mag (log(L/L_⊙)>4.8) are variable at Δ m_R>0.05 mag. Our period analysis finds 63 with significant pulsation periods.  Using the periods found and the known values of M_K for these stars, we derive the RSG PL relation in M31 and show that it is consistent with those derived earlier in other galaxies of different metallicities.  We also detect, for the first time, a sequence of likely first-overtone pulsations. Comparison to stellar evolution models from  confirms the first overtone hypothesis and indicates that the variable stars in this sample have 12 M_⊙<M<24 M_⊙. As these RSGs are the immediate progenitors to Type II-P core-collapse supernovae (SNe), we also explore the implication of their variability in the initial-mass estimates for SN progenitors based on archival images of the progenitors. We find that this effect is small compared to the present measurement errors.








§ INTRODUCTION

Evolution of massive stars, those with zero age main sequence (ZAMS) spectral types O and B, populating the upper part of the Hertzsprung-Russell (HR) diagram, remains one of the most puzzling fields in stellar physics. Its details are influenced by many uncertain physical processes, which include, among others, mass loss, convection, and rotation (e.g., , , , , ; see  for a comparison of different models employing varied prescriptions for these physical processes). Furthermore, there are many unknowns surrounding the mapping of these stars to their final fates, with cases of observational results conflicting theory—e.g., what produces which types of core-collapse supernovae (SNe), and what results in black holes (; ).  


Despite the incomplete knowledge of these aspects of the evolution of massive stars, they remain from birth to death one of the key players moderating and tracing properties of galaxies, such as the star formation rates, energetics—both radiation and mechanical—and chemical enrichment (see  for a recent review). Red Supergiants (RSGs), which are core He-burning (or beyond) stars representing the red-most excursion in the HR diagram during evolution of massive stars in the ZAMS mass range ∼ 1030 M_⊙ (for Population I stars), are the largest and among the optically brightest, and thus easily detected stars in the Local Universe. As a result of direct detections in pre-explosion images, they are recognized as the progenitors of the most abundant type of core-collapse SNe, namely the Type II-P <cit.>. 

In addition, many RSGs are variable, the properties of which have long been interpreted as pulsations <cit.>. The driving mechanism for pulsation in RSGs is not yet fully understood but it is thought to be the κ mechanism in the hydrogen ionization zone, coupled with convection, whose feedback however remains unknown <cit.>. Despite this theoretical uncertainty, the existence of an observed period-luminosity (PL) relation for RSGs entails their potential use as extragalactic distance indicators aided by their high luminosities and their higher prevalence as compared to Classical Cepheids <cit.>.

Most early work on the PL relation focused on the infrared, particularly the K band, for determining the luminosity, as these objects suffer from strong absorption in the optical, for example due to presence of TiO in their atmospheres <cit.>, and the bolometric correction is small in the K band. From such a relation, <cit.> were able to obtain a distance to the galaxy M101 consistent with the Cepheid distance.  However, studies investigating the PL relation for RSGs have been just a countable few <cit.> in the almost four decades starting from the pioneering work of <cit.>.  


RSGs have periods of typically several hundreds of days, and a semi-regular and complex variability; consequently, a sufficiently long baseline is needed for their studies. This has been largely made possible either through decades-long collections of the American Association of Variable Star Observers (e.g.,  for the Galactic RSGs) or recent surveys such as ASAS <cit.> and MACHO (; e.g.,  for the RSGs in the Magellanic Clouds). However, until now, there has not been a study of the variability of RSGs in the closest spiral galaxy to us, M31, even though a good number of them have been identified and cataloged (Sect. <ref>). This was due to a dearth of long-baseline time-domain surveys for M31. However, with recent surveys such as the (intermediate) Palomar Transient Factory (iPTF) survey <cit.>, the WeCAPP survey <cit.> of the bulge of M31 (2000–2003; ) and the PAndromeda project (2010–2012) of Pan-STARRS 1 <cit.>, such analyses are now possible. iPTF, in particular, possesses a baseline ≈ 1.52.5 times longer than the other surveys and probes sources down to m_R≈ 21 mag over the whole galaxy (cf. Sect. <ref>). In this work, we present the study of the variability properties of RSGs in M31 and a measurement of their PL relation using observations from iPTF.


Furthermore, while the derivation of physical parameters for SN progenitors based on direct detections or upper limits in archival images has become more complex in recent years, such analyses still currently neglect any possible variability of the progenitor star <cit.>. Our knowledge of the behavior of the star in the few years leading to the supernova is far from complete, with observational and theoretical work pointing to the possibility of pre-explosion eruptions, <cit.>, enhanced mass loss <cit.>, and growth of pulsations <cit.>. Even if pulsational variability similar to that observed in known RSG populations were still to be present during the archival imaging of the star, it could have a consequence in interpreting its initial mass. In particular, if pulsation amplitude is  luminosity-dependent, systematic effects could influence current observations, which point to a deficiency of higher-mass progenitors for SNe II-P (dubbed the “Red Supergiant Problem”; e.g., ; see Sect. <ref>). As an illustrative example, we also assess quantitatively the possible extent of such an effect in this paper. 

The paper is organized as follows. In Sect. <ref>, we introduce the sample of RSGs in M31. We describe the iPTF data in Sect. <ref> and show the RSG lightcurves constructed from these data in Sect. <ref> and the method to extract their periods in Sect. <ref>. We then derive the PL relation in Sect. <ref> and compare it with theoretical  models in Sect. <ref>, and touch upon the implication of the variability of RSGs for SN progenitors in Sect. <ref>. We end with our conclusions in Sect. <ref>. 













§ RED SUPERGIANTS IN M31



 §.§ The sample of RSGs and their physical parameters


For studies of resolved stellar populations of several galaxies in the Local Group, the Local Group Galaxies Survey (LGGS; ) has proved to be a milestone. The survey was conducted with the 4-meter telescopes at the Kitt Peak National Observatory (KPNO) and the Cerro Tololo Inter-American Observatory. In particular, 10 fields of M31 covering 2.2 deg^2, were imaged with the Mosaic CCD camera at the 4-meter Mayall telescope at KPNO, from which broadband photometric measurements in multiple optical filters (U,B,V,R,I) of over 350,000 stars in M31 have been published by <cit.>. They achieved a photometric precision of ≈ (12)% at 21 mag, thus providing an excellent catalog, particularly for studies of massive stars.   

<cit.> formulated a method for photometrically selecting RSGs, which serves to remove as far as possible the contaminating foreground red dwarfs. The author used the two-color diagram of B-V versus V-R, wherein the two classes of objects are found to separate out into two sequences. For these spectral types, the V-R color mostly traces the effective temperature (T_ eff) while the B-V color also traces the surface gravity. Lower surface gravity objects (i.e., RSGs) have significantly redder B-V values (by some tenths of a magnitude) due to increased importance of metal line blanketing, which is prominently expressed in the B filter. Using the LGGS catalog, <cit.> compiled a sample of 437 photometrically selected RSG candidates based on these two colors. From this photometric sample, <cit.> [hereafter ME16] then obtained spectra of 255 stars and measured their radial velocities, in order to confirm their membership of M31. This yielded a large, spectroscopically pure sample of RSGs, which forms the fiducial sample of RSGs for our study, coupled with the data from the iPTF survey (described below). 

For their sample of confirmed RSGs, ME16 derived physical parameters (e.g., luminosities, temperatures), which we utilize in the sections below (cf. Sect. <ref>). In order to derive luminosities for these stars, ME16 used available near-IR K band photometry, in order to take advantage of the small bolometric and extinction corrections as compared to the optical bands. K band photometric measurements were acquired from a combination of targeted observations (with FLAMINGOS; ) and the Two Micron All-Sky Survey (2MASS; ). The bolometric correction for K band was computed as a function of T_ eff, and values of T_ eff for the stars were determined by fitting the  atmospheric models to the observed spectra.  Finally, ME16 also obtained initial-mass estimates of the stars by comparing to the  stellar evolutionary tracks (see ME16 for more details), and spectral types following <cit.> using the strengths of TiO bands (late K and M type stars) and that of G band and Ca I λ 4226 (early and mid-K type stars).

<cit.> contend that the T_ eff values obtained from fitting the TiO bands in optical spectra—the method employed by ME16—underestimate the true values by several hundred Kelvin, since these molecular lines form high up in the atmosphere with lower temperature. In contrast, <cit.> find T_ eff≈4150±150 K for all spectral types of RSGs in the Magellanic Clouds that they have analyzed, when fitting the continuum of the full spectral energy distribution. However, <cit.> argue that this uniformity in T_ eff is inconsistent with the observed variation of RSG spectral types with metallicity (for example, comparing the Milky way, LMC and SMC RSGs; ), and the fact that the Hayashi limit moves to higher T_ eff values with decrease in metallicity <cit.>. 
Thus, the uncertainty in the temperature scale of RSGs is as yet unresolved.   
Indeed, the extent of the atmospheres of these supergiants and its consequences for the occurrence of convection, stratification in the temperature, etc., undoubtedly add to the complexities in modeling them, which still remains a challenge. For this paper, we mostly use the K-band magnitudes of these sources from ME16, and turn to their luminosities for comparing to stellar evolutionary models (Sect. <ref>) and in discussing the SN progenitor mass estimates (Sect. <ref>). 





 §.§ The iPTF data

The iPTF survey, carried out with the 1.2-meter Samuel Oschin telescope at Palomar Observatory, continued imaging M31 (including the outskirts of this galaxy) from its predecessor, the PTF survey <cit.>, with the same wide-field detector covering more than 7  deg^2 with 7 active CCDs, and a typical spatial sampling of 2” full-width at half maximum (FWHM). The baseline of the data for M31 analyzed here extends from 2012 May until 2017 February. The imaging was done in two passbands, R and g, with the largest fraction (>80%) in R at a limiting magnitude reaching ≈ 21 mag and an average cadence of 1 day (typically taking 2 observations per night).  

We thus use only the R band PTF/iPTF (hereafter simply PTF) data set of about 10000 images. This set comprises over 1500 observation epochs of M31 covering 1.8×2.4  deg^2 with 6 CCDs (as shown in Fig. <ref>), over the approximately 5 years of baseline. Difference imaging was performed on all these images via the difference imaging pipeline of the PTF collaboration, the IPAC/iPTF Discovery Engine (PTFIDE; ), and this pipeline also carried out detections of varying sources on the difference images, which we term here as “raw detections”.




 §.§ Optical lightcurves of RSGs and their variability


The PTFIDE raw detections catalog, in principle, provides us with a means to determine whether a source in M31 is variable or not. In fact, we find from the results of cross-matching (with a search radius of 2”, consistent with the typical FWHM) that all RSGs lying within the PTF footprint (253) are flagged as variable by the PTFIDE pipeline. However, as already established by <cit.> and <cit.>, these “raw detections” are largely dominated by artifacts of image differencing, for example, dipoles from imperfect PSF-matching, edges of masked image regions due to saturated stars, CCD defects, etc. 

For our science goal, we do not want to compromise completeness by resorting to thresholding on output parameters of image differencing. Rather, we implement forced photometry to construct the RSG lightcurves, using the difference images themselves, and analyze their variability. Difference imaging significantly alleviates the problem of crowding in the M31 fields for measuring their fluxes. Subtraction artifacts in individual difference images will largely be reflected as rogue points in the lightcurves, which we effectively deal with (e.g., by masking those points) without discarding the source. 

In the same manner as <cit.>, we perform aperture photometry at the positions of the M31 RSGs and apply a curve-of-growth correction to the measured fluxes. The latter is a different factor for each star applied to correct for the flux missed due to the limited size of the aperture (taken to be the FWHM). The subtracted fluxes of the sources are then added using the template/reference image PSF-fit photometry catalog and calibrated using the relevant zero-point in the science images (see  for details). We drop 9 RSGs with bad photometry in a large section of their lightcurves, for example due to their location close to masked parts of images or in regions close to the bulge where the quality of image differencing is poor. Example calibrated lightcurves are shown in Fig. <ref>[The PTF lightcurves of all the RSGs studied here are available from the corresponding author.]. 







To examine the variability of these RSGs over the full baseline, we compute the root-mean-square (RMS) deviations from the mean for all of the extracted lightcurves, and compare them to those of static stars. We define a static star as a star in the catalog of the reference image used in the image subtraction, having no detection from the PTFIDE pipeline in all of the input images used in the subtraction (also termed science images by ). We extract the lightcurves for these static stars in the same manner as those of the RSGs. Since the surface brightness of M31 is not uniform over the face of the whole galaxy, there is a gradient in the local background noise of a resulting difference image, which deteriorates severely particularly toward the bulge.  The RSGs, on the other hand, are largely distributed in the disk of the galaxy. To obtain a sample of static stars with a similar background noise distribution as the population of RSGs, we select the static stars that are located in proximity (within about 1') of any of the RSGs.



In Fig. <ref> we plot the distribution of the RMS deviations found for these static stars as a function of their average R band magnitudes over the baseline of the survey (blue points). RMS values increase toward fainter magnitudes, and the distribution we find for stars in M31 agrees well with that for M-dwarf stars in PTF derived by <cit.>. Also shown in Fig. <ref> are the M31 RSGs (red circles) and we see that the variabilities (as determined by the RMS deviations) of some of the RSGs, around the mean <m_R> range of 18.0–19.0, are comparable to those of static stars, i.e.,  their variability is consistent with noise. Filled circles indicate RSGs with measured periods (Sect. <ref>).  

We visually separate the RSGs with significant variability from the RSGs whose RMS deviations are consistent with noise. To this end, we draw a straight line in the plane of <m_R>Δ m_R, with the slope similar to the linear fit through the static stars and try different intercepts to most efficiently select significantly variable RGSs. This demarcation is shown with the orange line in Fig. <ref>. This method is not very stringent, but we verified that the main results of this paper remain the same with small shifts in this line. Adopting this demarcation, we find that 167 RSGs, constituting ≈ 70% of our sample, have measurable variability from PTF. These are examined in the following section to determine what periodicity, if any, is present in their multi-year light curves. 





In Fig. <ref>, we plot the RMS deviations for our sample of M31 RSGs versus their M_K magnitudes (left panel) and luminosities (right panel), whenever available from ME16. RSGs with significant variability are plotted in red, while those with PTF upper limits on variability (those below the orange line in Fig. <ref>) are plotted in black. 

As can be seen, at low luminosities there is an overlap between the RSGs with detected variability and those below our detection threshold, despite the fact that they fully separate in the <m_R>Δ m_R plane. This can be attributed to the varying degree of extinction in the R band suffered by individual RSGs, which subsequently modulates their signal-to-noise ratio (S/N). However, despite this overlap, it is notable that all of the “non-varying” RSGs are characterized by low intrinsic luminosity. All RSGs in our sample brighter than M_K≈-10 (and log L/L_⊙≈4.8) show variability and many of them at the level Δ m_R>0.1. We assume this variability of RSGs to be associated with pulsations <cit.>.
Further, our result points to variability/pulsations occurring already for  log L/L_⊙≳4.6. Theoretical calculations of <cit.> could only find pulsations above log L/L_⊙>4.95 (e.g., their Fig. 1). This is, however, due to the fact that the implicit scheme adopted in stellar evolution calculations introduces substantial damping and thus, they could only place a lower limit on the existence of pulsations in RSGs[Private communication with M. Cantiello.]. Our result on the other hand, shows that stars with mass lower than the 1617 M_⊙ limit found by <cit.>, can also pulsate.
In addition, we find the hint of a positive correlation between variability amplitude and brightness for the RSGs with detected variability. In Sect. <ref>, we will explore the influence that this observed variability can have on progenitor mass estimates of core collapse supernovae, in particular the trend of larger amplitude variations for more luminous RSGs.
 











 §.§ Period determination


The lightcurves of RSGs are known to be semi-regular, and some of them are completely irregular (see Fig. <ref>). The periodicity of unevenly-sampled regular time-series/lightcurves has been efficiently handled by conventional standard algorithms, for example, Lomb-Scargle periodogram ( LS; ), phase dispersion minimization (PDM; ), analysis of variance <cit.>, and the hybrid algorithm of <cit.>.  However, for the RSGs, these methods tend to fall short, as these objects are characterized by one of the most complex lightcurve morphologies. This includes a strong red-noise component speculated to arise from the convective cells on the surface and multi-periodicity. For example, <cit.> found two distinct periods for a significant fraction of their Galactic RSG sample over a baseline extending to more than 5 decades: one on the order of a few hundred days, which is typically associated with RSG pulsations, and the other greater than 1000 days, similar to the long secondary periods in Miras with nature unknown (see also ). Compared to this, our baseline will not be sensitive to pick up these longer (>1000 days) periods. 





Bound by the limitations of the popular standard methods for our application, we turn to a relatively recent approach based on Gaussian Process (GP) modeling (a non-parametric statistical model that assumes the distribution over a set of random function values as a multi-variate Gaussian, typically with mean 0; see ). In principle, this Bayesian technique provides a statistically robust and powerful way to model any phenomenology, particularly ones for which prior knowledge may not be available.   
For a GP, determining the covariance matrix is the crux of the modeling, which in turn is determined by the power spectrum of the signal. In almost all of the existing off-the-shelf GP-based methods used by the stellar community, a functional form is assumed for generating the covariance matrix elements either assuming stationarity for the lightcurve (i.e., its statistical properties are assumed constant) or a single frequency (for periodic signals, e.g., ). Given the complexities of the RSG lightcurves as mentioned above, these algorithms are not accurate for their analysis. 

For our case, the GP modeling is performed by using the so-called critical filter algorithm of signal reconstruction as implemented by <cit.> using the  package of <cit.> aimed toward cosmological signals, and applied to period analysis by Oppermann et al. (2018, to be submitted).
Simply put, the critical filter provides a methodology to reconstruct the underlying signal, whose power spectrum is unknown, by simultaneously also reconstructing this spectrum from the observed data. The computation is iterative. We start with an initial guess of the power spectrum that defines the prior distribution, and reconstruct the signal thus obtaining a posterior distribution conditioned on the observed data and the given power spectrum, which is again a Gaussian. The latter distribution gives an improved estimate of the power spectrum, which is used to update the prior distribution, and then this is used to reconstruct the signal and so on. The iteration is stopped when the new estimate of the power spectrum converges.  Since <cit.> have illustrated the presence of a 1/f component in the power spectra of these RSGs, this provides a natural choice as our initial guess, but we also impose the condition that there has been at least five iterations before convergence. We obtain a similar reconstruction even when starting with a flat spectrum.  

Examples of the resulting power spectra of the lightcurves are shown in the right panels of Fig. <ref>, with the red-noise component evident. To extract the peaks in a spectrum, we fit a power law to the background red-noise and use that to set a S/N threshold of 7 above the background. We obtain the slopes for these RSGs in the range -1.51 to -0.86. We then select the peaks above this threshold with values greater than 10^-3  day^-1 (constrained by our baseline). For cases with more than one peak, we choose the one with the maximum power (i.e., the dominant mode) and for a few sources with blended peaks, we take their excess power-weighted average. We adopt these frequencies to define the periods of the lightcurves used in deriving the PL-relation (Sect. <ref>). 

The methodology we have used reconstructs the lightcurve along the complete time axis, fully taking into account the observational sampling. We have verified by shuffling the magnitude measurements while keeping the same irregular time sampling and calculating the power spectrum of the resulting lightcurve, that no artificial periods corresponding to the observational pattern emerge as power peaks. 

A detailed analysis of how the results from other period-finding algorithms fare against that of the GP modeling is beyond the scope of this paper (it will be presented by Oppermann et al., 2018; to be submitted). Nevertheless, for a simple comparison, we use available  implementations of common period-finding algorithms, specifically LS, supersmoother[<h>ttps://www.astroml.org/gatspy/periodic/supersmoother.html] and PDM, to estimate the periods of some example lightcurves of M31 RSGs. We consider the RSGs with IDs 227 and 189 shown in Fig. <ref>, with quite complicated shapes. The different algorithms give a range of periods: 139 days (LS), 1113 days (supersmoother),  472 days (PDM),  505 days (this work) for RSG_227;  263 days (LS), 1329 days (supersmoother), 786 days (PDM), 331 days (this work) for RSG_189. As can be seen from the power spectra of these two RSGs in Fig. <ref>, the results from the other methods, though associated with power peaks in most cases, do not pick up the dominant peaks barring PDM for RSG_227.


We successfully obtain the pulsation periods (less than ≈ 1000 days) for 63 RSGs (the remaining 104 do not have a significant peak), and these are shown in Table <ref>. 
We mark these RSGs in Fig. <ref> by the filled red circles, and as is evident, some bright, highly variable sources do not have detectable periods, an example of which is shown in the last panel of Fig. <ref>. It is quite likely that with an extended baseline (for example with the Zwicky Transient Facility, ), a period >1000 days could be found for such cases (see also Sect. <ref>). Another example of a bright, but less variable, non-periodic source is also shown in Fig. <ref>. 
A detailed investigation of why some of the sources in our sample show detectable periods and others, though variable, show none is beyond the scope of this paper. As described above, the mechanism of pulsation in RSGs and its interaction with convection is not well understood. 













§ PERIOD-LUMINOSITY RELATION OF RSGS IN M31


Figure <ref> shows the periods of the RSGs obtained in the previous section and their absolute K-band magnitudes from ME16 (Sect. <ref>). Of the 63 RSGs with identified periods, 51 appear as a coherent band between log(P/ days)≈2.4 at M_K≈-9.5 and log(P/ days)≈3.0 at M_K≈-11.5. However, 12 RSGs follow a separate sequence (copper circles, leftmost portion of the plot) with 10 of them appearing at similar M_K<-9.5 as the bulk of the RSGs in the adjacent band, but with an offset of log(P/ days)≈0.3. The nature of these stars will be examined in Sect. <ref>.

The remaining two stars of the 12 have periods <100 days (IDs 52 and 166, cf. Table <ref>). Short periods of around 100–150 days have been found for RSGs <cit.>, but no RSGs in the literature have pulsation periods as short as these stars, prompting us to examine them in detail. Both stars possess extremely regular lightcurves as compared to the other RSGs (see Fig. <ref>), relatively high T_ eff values (4300 K on the temperature scale of ME16), and relatively low luminosities (M_K≥-9.5 mag). These physical properties overlap with those observed for ultra-long period Cepheids <cit.>. Indeed, star RSG 166 was classified as a Cepheid by <cit.>. As the RSG sample of ME16 was simply selected to have temperature cooler than around 4300 K and be a member of M31 (as opposed to strictly stars on the RSG branch), it is not wholly unexpected that a small number of stars on the luminous and cool edge of the Cepheid instability strip would be identified. 


Ignoring the 12 RSGs from above, we obtain a simple linear fit through the rest as 

    M_K=(-3.38±0.27)×log P+(-1.32±0.75),

shown by the blue line in Fig. <ref>. 
The RSG PM_K relation we have derived here is the first for M31. We see a dispersion around this relation of 0.29 mag. A part of it is possibly due to the errors in the adopted K band magnitudes, which could reach a tenth of a magnitude (for example, the 2MASS values given by  and used in ME16), and the fact that these K band magnitudes are from single (nightly) observation epochs, while the amplitude in this band could be as high as 0.25 mag <cit.> (thus, a contribution of <0.18 mag to the dispersion). Also, circumstellar reddening may contribute to this dispersion; ME16, however, applied an average correction for M31 RSGs obtained by <cit.> to the M_K values of these stars that should reduce this effect. But a possible part of the cause could also be that the dispersion is intrinsic for these objects—an apparent effect of additional parameters (cf. Sect. <ref>).





Given the broad trend of increasing variability with brightness (Fig. <ref>) discussed in Sect. <ref>, and the positive correlation of period with brightness (Fig. <ref>), the same broad trend can be expected between the periods and the RMS deviations. This is indicated in Fig. <ref> by the gradient in the face-color of the points. In light of such a trend, it is feasible that a period longer than the range probed in this study is present in the highly variable sources without detectable periods (Fig. <ref>, Sect. <ref>). Furthermore, the theoretical result of increasing pulsation amplitude with initial masses or luminosities <cit.> appears consistent with our observational result. 


Figure <ref> shows the comparison of our results for the PM_K relation with those of others for different galaxies. Our slope and intercept are consistent with those of <cit.> for the Galactic RSG sample with 13 sources and a dispersion of 0.46 mag. Furthermore, they are in good agreement with those of the LMC, SMC, and M33 derived by <cit.> using 47, 21 and 40 sources, respectively, with dispersion value of ≈ 0.3 mag (they claim their relation for the Galactic sample is less accurate due to uncertainties in distances)[<cit.> adopted distance moduli of 18.50 (LMC), 18.91 (SMC) and 24.93 (M33).]. Thus, the data from several galaxies characterized by different metallicities (from sub-solar by a factor of ≈ 4 in the Small Magellanic Cloud to super-solar by a factor of ≈ 2 in M31, e.g., ), are all consistent with a universal period-luminosity relation for RSGs, although the uncertainties still allow a wide range of parameters. Such universality may not be unexpected given current ideas on what drives the pulsation in RSGs (κ mechanism in hydrogen, Sect. <ref>), and this potentially could have important implications for the life of very low-metallicity/Population III stars <cit.>.




§  MODELS: PHYSICAL PARAMETERS OF RSGS

We investigate the pulsation mode represented by the periods we have found for the RSGs with theoretical models constructed using Modules for Experiments in Stellar Astrophysics (; ), a standard tool for stellar evolutionary calculations. One of the latest updates presented by <cit.> includes fully incorporating existing methodology for asteroseismology computations ( software of ), and thus facilitates comparison of our observational results with theory. 


To produce theoretical predictions for RSG pulsation periods, we use version 9793 of the  software. We consider stellar models with ZAMS masses of 1224 M_⊙, which are evolved from ZAMS to depletion of Carbon in the core. Our  settings follow the recent work of <cit.>, including a calibrated mixing length parameter α_ MLT = 2.7 to match inferred RSG temperatures in M31 at Z = 0.04. We also include step overshoot with parameter f_ ov = 0.15 and mass loss according to the “Dutch” wind prescription <cit.>. More details on  settings can be found in <cit.>. Our models do not include any rotation.


After the models leave the main sequence, we include adiabatic  calculations in each step for T_ eff<4500 K. These calculations find all radial modes with periods in the range 10–10000 days, more than adequate to identify the fundamental mode and first overtone in these RSG models. The left panel of Fig. <ref> shows the resulting period tracks for several different mass models, along with points representing 20 equally spaced time intervals along each of the tracks for the fundamental mode periods. When comparing to the observed RSG data in the right panel of the same figure, we show the same tracks only after the onset of core helium burning, where the models spend more than 95% of their RSG lifetime.

Our modeling makes no attempt to discern which modes should be excited. The work of <cit.> suggests that RSG pulsations are excited by the κ mechanism in the hydrogen
ionization zone. <cit.> also find that pulsation periods roughly scale as P ∝ L/M in agreement with the analytic predictions of <cit.>. Our  models roughly agree with this scaling.


The right panel of Fig. <ref> shows that the bulk of the pulsating RSGs in our sample are consistent with the fundamental mode pulsation of the theoretical models. Furthermore, the sequence of 10 stars offset to shorter periods for a given M_K, is broadly consistent with the first overtone mode. 





The Hertzsprung-Russell (HR) diagram for cool supergiants is shown in Fig. <ref>. Colored lines represent the same  models used to produce Fig. <ref>, and the RSGs with measured pulsation periods are shown in grey (other M31 RSGs are plotted as blue circles for reference). The observed RSGs lie in proximity to  model tracks for broadly the same initial masses, whether period or T_ eff is used as parameter (right panel of Fig. <ref> and Fig. <ref>).  
The 10 stars with measured periods offset from the bulk of the sample in the period-luminosity plane that are not Cepheids are shown by the star symbols in the HR diagram. These stars with shorter periods overlap with the bulk of the RSG population with measured periods in their temperatures and luminosities. This strengthens our association of the variability observed in these stars with the first overtone mode of stars with the same properties as the bulk population.


As is already known, there is no single-valued mass-luminosity relation for the red supergiants. Consequently, tracks of different masses at the same luminosity can smear out the period for the RSGs, and may be a reason for the scatter seen in Fig. <ref>. Thus, the mass of the star, and more importantly its evolutionary stage, naturally appear as additional parameters contributing to some of the dispersion in the P-M_K relation. As a result, the precision of any extragalactic distance estimate based on this relation can never be better than the limit imposed by the intrinsic dispersion.  




§ EFFECT OF INTRINSIC VARIABILITY ON SN II-P PROGENITOR MASS ESTIMATES


RSGs as progenitors of Type II-P/L SNe <cit.> have been established based, in part, on the direct detection of some progenitor stars in archival images of nearby galaxies hosting these SNe <cit.>. A recent comprehensive review of such observations spanning a 15 year period, and their implications for the progenitor masses of hydrogen-rich core-collapse supernovae, was presented by <cit.> [hereafter SS15]. The general approach to determining progenitor masses from these pre-explosion observations is to convert the photometric measurements to bolometric luminosities and then compare to the initial mass-final luminosity relation from stellar models. One of the important results from such an analysis is the maximum-mass for these progenitors is claimed to be at a much lower value (≈ 18 M_⊙)  (; SS15) than the high-mass limit of RSGs (2530 M_⊙) predicted by stellar evolution theory. 

However, various systematic errors, often resulting from the uncertain final stages in RSG evolution, can influence this result. In particular, when converting observed fluxes in a limited (often single) photometric band to bolometric luminosities, treatments of extinction and bolometric corrections are critical. For example, <cit.> showed that the spectral type of RSGs evolves to later type as they approach core collapse, resulting in a larger bolometric correction than typically assumed for the imaged progenitor stars. They thus pointed out that failing to account for the spectral evolution of RSGs can introduce a systematic error in the adopted bolometric correction.  

Most of the archival detections of supernova progenitors are from single-epoch observations, 
and we have demonstrated in Sect. <ref> (see Fig. <ref>) that while lower luminosity RSGs have observed variability of ≲ (510)%, a large fraction of the RSGs with luminosities above log(L/L_⊙)≈4.8 show a significant variability. Thus, if our observed sample of M31 RSGs is representative, intrinsic variability of RSGs will preferentially affect the initial mass estimates of high luminosity (and therefore high initial mass) progenitor stars. It is therefore pertinent to examine how big such an effect can be. 


Recently <cit.> (see also ) examined the variability of four SN II-P progenitor stars, from approximately four to eight years prior to the supernovae. They derived upper limits on their variability at (5–10)% of their R band luminosities, corresponding to Δ m_R<0.050.10. Three of these SN progenitors had sufficient observations to derive bolometric luminosities, all of which were in the range log(L/L_⊙)=4.5-4.9 (). 
Thus, as is evident from Fig. <ref> (right panel), the upper limits on their variability are also consistent with our results for the RSGs in M31. 





In this section, we make a quantitative assessment of whether intrinsic variability of RSGs can have any significant effect on the SN progenitor mass estimates.  

In converting a photometric measurement in some filter (m_λ) to the bolometric luminosity (ℒ=log (L/L_⊙)), a correction (dm_λ to obtain the mean value <m_λ>) for possible variability of the progenitor can be applied as

    2.5·<ℒ>=M_⊙-(m_λ+ dm_λ-A_λ+ BC_λ),

where the term inside the brackets in the right-hand side is the bolometric magnitude of the star, A_λ and BC_λ are the extinction and bolometric corrections, respectively, for the star in the given filter, and  M_⊙ is the bolometric magnitude of the sun. In absence of color information to constrain the BC_λ of the star, typically the value for an M0 supergiant is assumed with its 1σ error as 0.3 mag corresponding to the standard deviation of the bolometric corrections spanning the spectral types from late K to late M (see , SS15). 


The detections and upper limits of the progenitor stars have been acquired mostly in the I and R bands. Since we do not have information on the intrinsic variability of our RSG sample in the I band, we assume a similar level of  variability as in the neighboring R band (Fig. <ref>). Then, given a random-phase measurement m_λ from the pre-explosion image of the star, and hence, a random-phase luminosity ℒ computed from it, we are interested in the mean luminosity <ℒ> for the star. To this end, we use Bayes' Theorem to derive the posterior probability distribution for the mean luminosity <ℒ>, given the observed random-phase luminosity ℒ, according to

    P(<ℒ>|ℒ)=P(ℒ|<ℒ>)· P(<ℒ>)/P(ℒ).


The luminosity distribution of our RSG sample is taken as the prior P(<ℒ>). We obtain a similar result when using a flat prior. To construct the likelihood P(ℒ|<ℒ>), we use our full sample of M31 RSGs with luminosity available from ME16, i.e., all the RSGs in Fig. <ref> both below and above the PTF variability sensitivity limit. We bin the stars in log L  with bin-width ≈ 0.1 (Fig. <ref>), and compute the distribution of dm=m-<m> magnitudes using the lightcurves of all the stars in each bin. Our construction of relative fluxes in the PTF R band assumed these relative fluxes are similar in the other filters where most progenitor images were taken. 
We scale the resulting distribution in the corresponding bin containing <ℒ> by a factor 2.5 and shift it by <ℒ> (cf. Eq. <ref>) to obtain the likelihood. The denominator in Eq. <ref> is just a normalization constant. With all these ingredients, we compute the posterior distribution P(<ℒ>|ℒ) and then obtain the mean and standard deviation of this distribution. 


The result is shown in Fig. <ref> by the red points. The diagonal dotted line in the figure shows the 1:1 relation between the random-phase ℒ and the mean <ℒ>, and thus represents the typical assumption made when computing supernova progenitor masses from pre-explosion images. 


From this figure, we see that <ℒ> vs. ℒ for our RSG sample is consistent with the 1:1 relation. 
The vertical bars of the red points represent standard deviations of the <ℒ> distributions, and as can be seen, they become larger with increasing luminosity, reflecting the fact that the observed variability in our sample is larger for higher luminosity stars (Fig. <ref>). We also plot the measured luminosities along with the measurement errors of the SN II-P progenitor stars given by SS15 with blue points based on the 1:1 relation. As is evident, the measurement errors (from the photometry, extinction and bolometric corrections, cf. Eq. <ref>) of the random-phase luminosities of the progenitors are much larger than the standard deviations of the <ℒ> distributions.
The overall shift in the masses of the progenitors accounting for the intrinsic variability (a mere 12 M_⊙) is not significant, given the rather large uncertainty on the applied corrections (Δ A_λ and Δ BC_λ) and photometric measurements (Δ m_λ). The uncertainty due to variability (dm), however, is an inevitable physical effect that cannot be overcome by taking more precise measurements and thus will eventually become the dominant effect.

The analysis presented above rests on the assumption that our full sample of RSGs in M31 is representative of the variability that RSGs exhibit in the final years before core-collapse. If these properties change as a RSG enter the final nuclear burning stages, our results and conclusions could be influenced. Further, we have assumed that the variability of RSGs in the PTF R band is similar in the filters (largely I) that the SN II-P progenitors were imaged. If RSGs were markedly more variable in the I band, then our results could be affected. Nevertheless, there is no obvious argument for variability amplitude to be larger in I than in R.

While intrinsic variability of RSGs with respect to luminosity estimates of SN II-P progenitors is shown here to be not significant, it could, however, play an important role in other physical aspects of RSGs of consequence to their initial mass estimates, in particular a possible enhancement of the mass loss <cit.>. A similar effect is known to take place in Miras <cit.> and in Cepheids <cit.>, wherein variability due to radial pulsations and the associated lifting and lowering of extended atmospheric layers is regarded to be an important driver for mass loss.  
Mutli-epoch monitoring, in multiple filters, for variability of potential SN II-P progenitors, like the efforts of <cit.>, <cit.>, <cit.> will continue to be beneficial and have implications beyond the photometric uncertainty of the progenitor stars.




§ CONCLUSIONS

The Palomar Transient Factory's nearly 2000 days of monitoring of M31 has allowed for the first measurements of the prevalence of RSG variability in this important nearby galaxy. In this first study, we confined our efforts to the 255 bright RSGs that were cataloged and typed by the LGGS group (ME16). As was known for our galaxy, the
Magellanic Clouds and M33, these evolved massive stars can be highly variable, and often periodic. Based on the RMS deviations of the PTF R band lightcurves, 167 of the 255 RSGs exhibit variability above the noise limit of PTF data (Fig. <ref>). For these sources, we see a broad positive correlation between the variability amplitudes, represented by the RMS values, and their brightness.

We found that the lightcurves of these variables are characterized by semi-regular (and irregular for some cases) long-period pulsations over a strong red noise <cit.>.  Using a robust, probabilistically motivated methodology for reconstruction of their power spectra, we determined significant pulsational periods for 63 of the variable RSGs.  Using the absolute magnitudes measured in the K band by ME16, we derived the PM_K relation for the RSGs in M31 and showed that it is consistent with that of the Galaxy, the Magellanic Clouds, and the M33 RSG samples, which represent a wide range of metallicities. This points toward a universality of this relation, at least for the Local Universe. It will be interesting to explore this relation with a large sample of RSGs from a larger number of distant resolved galaxies with the upcoming large telescopes like the Large Synoptic Survey Telescope (LSST; ) and the European Extremely Large Telescope <cit.>. 

Comparison of the measured Period Luminosity relation with the theoretical one based on the  models shows that the pulsations are consistent with the fundamental radial mode. At the same time, a clear sequence of outliers is seen in the distribution with our observed sample, which is found to be consistent with the first overtone mode based on the models—another novel feat of this work.

Finally, we explored whether the measured RSG variability could be large enough to confound the inference of stellar mass from the observed supernova progenitor luminosity. We found that this effect is small compared to current observational uncertainties. 
 
Future insights on RSG variability in M31 with the existing PTF data could be gained from large-scale variability studies of bright stars not necessarily in the ME16 catalog. In addition, with the recent commissioning of the ZTF instrument, continued monitoring of all RSGs will be possible and hopefully reveal longer period pulsations than possible with the original data set. Moreover, variability observations of different RSG populations spanning near and far host galaxies will be an asset for conducting many important studies of these stars. LSST in the near future, with its superb photometry and expected cadence of three days, has the potential to generate such a product from monitoring a select group of galaxies.






Acknowledgments: This work was supported by the National Science Foundation through grants PHY 11-25915, PHY 17-148958 and ACI 16-63688 and funded in part by the Gordon and Betty Moore Foundation through Grant GBMF5076.  
MS is grateful to Niels Oppermann for the detailed discussions on statistical approaches and methods of proper signal reconstruction. MS also thanks Mathew Graham and Ashish Mahabal for helpful discussions on the period finding algorithms, and Matteo Cantiello for the helpful suggestions and comments on the manuscript. EB thanks Sang-Hyun Chun for providing  inlists necessary to reproduce their RSG models, and Matteo Cantiello for discussion of RSG evolution in . PyFITS is a product of the Space Telescope Science Institute, which is operated by AURA for NASA. Support for this work was provided to MRD by NASA through Hubble Fellowship grant NSG-HF2-51373 awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., for NASA, under contract NAS5-26555. MRD acknowledges support from the Dunlap Institute at the University of Toronto. 




 <cit.>,  <cit.>,  <cit.>, ,  <cit.>,  <cit.>,  <cit.>,  <cit.>,  <cit.>, <cit.>.

yahapj









