380 likes | 537 Views
Detecting Periodicities in RR Lyrae Stars. R. F. Stellingwerf Stellingwerf Consulting January, 24, 2011 RR Lyrae Stars, Metal-Poor Stars and the Galaxy Pasadena, CA. RR Lyrae Special Attributes. Highly non-sinusoidal light variations in many cases
E N D
Detecting Periodicities inRR Lyrae Stars R. F. Stellingwerf Stellingwerf Consulting January, 24, 2011 RR Lyrae Stars, Metal-Poor Stars and the Galaxy Pasadena, CA
RR Lyrae Special Attributes • Highly non-sinusoidal light variations in many cases • Changes in the period and amplitude due to evolution • Changes in amplitude due to Blazhko effect • Possible multiple modes with nonlinear modal interaction • Other variations caused by companions, eclipses, etc. • Gaps in data sets (yearly, nightly, historic) • Multiple overlapping data sets (different colors, observers)
Review of Methods Spectral: • Fourier (discrete, FFT, etc) • Periodogram (Vanicek/Lomb/Scargle) • PDM / PDM2 Period Search:
Discrete Fourier Method • Xkare coefficients of a Fourier decomposition of {xn} • Complex amplitude of transform produces power spectrum • Mathematically well understood • FFT available • Problems with unequally spaced data and data gaps • Limitations in sampling and frequency coverage
Periodogram / Least Squares Fit to Sin() • Find optimal fit to a sin() function at each frequency • Must adjust sin() amplitude, phase, average value and frequency
Least Squares Periodogram Formula • Sum of squares of deviations: • Pick A to minimize sum of squares… • Differentiate wrt A, set to 0, solve for A, get… • Note, min S usually implies max A. So, one approach is to varyw(and possibly f) looking for large values of A (or A2). • Alternate view: minimize
Periodogram Variations • Classic Form (Fourier form) • Lomb/Scargle Form (least squares form) • Possible questions about normalization, smoothing, and bandwidth effects
Phase Dispersion Minimization 8 sj Mean curve is determined by the data. Compute scatter in each bin. s …where s is sj summed over 10 bins • Minimize Q => least squares fit to mean curve through bin means • Very simple and efficient to compute • Additivity of the variances allows generalizations and variations
PDM Calculation 9 phase = frac( ti* fj ) bin = int( 10 * phase ) U s1 s2 s3 B s4 s5 s6 Time …where s is sj summed over 10 bins, all segments, all data sets • Use separate phases for time segments to get a quick trial estimate, use single phase across all times to refine the period
PDM Window Function 10 • Standard test case: 101 points cover 10 cycles of a pure sine wave, time 0->10, amplitude 1, f = 0->2, 200 frequency points Top: wide (5/2) bins Bottom: narrow (10/1) bins Frequency scans from 0 -> F will detect signals at all higher frequencies. (multi- periodic cases, or “subharmonics”).
PDM Treatment of Sigmas (error estimate) 11 • Deviation from mean is taken from edge of error bar • If sigma > distance from mean, point does not contribute to variance Bin points shown, radii are estimated sigmas. Deviations used in variance shown. More accurate points contribute more to result. Tests show this gives cleaner results than 1/sigma2 weighting.
PDM2 Options 12 • Classic PDM – step function fit 10 bins, double wide option for small data sets. • Linear mean curvefit. Only used for significant frequencies. • Spline mean curve fit. • Subharmonic Averaging (uses PDM window transform) • Beta Distribution and Non-Parametric significance tests.
Sin() Wave Test Case 13 • 10 cycles of sin() variation plus Gaussian noise • Lomb shows clear peak at f = 1.0 • PDM shows mins at f = 1.0, 1/2, 1/3, etc.
Sawtooth Test Case 14 • 10 cycles of saw() variation - No noise - Gap inserted in data • Lomb shows peak at f = 1.0, side lobes • PDM shows similar patterns at f = 1, ½, 1/3 • Main peak stronger and narrower
Sawtooth – Big Gap 15 • 10 cycles of saw() variation - No noise - Larger Gap inserted in data • Lomb shows peak at f = 1.0, many side lobes • PDM shows similar patterns at f = 1, ½, 1/3 Subharmonic averaging helps here. • This pattern is typical of a long time-scale data set
Pulse Test Case 16 • 3 Pulses at t = 3, 6, 9 • Lomb shows peak at f = 2/3, but not a clear winner. • PDM shows the f = 2/3 result more clearly. Note narrowness of main result.
Double Peak Test Case 17 • Alternating peaks, period = 2, f = 0.5 • Lomb shows major peak at f = 1, wrong frequency. • PDM shows f = 0.5 correctly, as well as a possible hit at f = 1.
Significance Testing 18 • PDM2 now includes the incomplete beta distribution (Schwarzenberg-Czerny, 1997, 1999) with multiple trial correction. This is the correct distribution for PDM in noise dominated cases. • The Lomb analysis follows an exponential distribution. • Note, however, that the “noise” present in many variable star signals is due to other modes, Blazhko variability, etc., and is not Gaussian, as assumed in any parametric statistical criteria. • Furthermore, the “bandwidth correction” depends on the data distribution in time, and may not be accurate. Used for Lomb also. • A more accurate estimate is now included using Fisher Randomization / Monte Carlo re: Nemec & Nemec (1985). This uses multiple permutations to derive a specific “noise” distribution unique to each data set and analysis parameters.
Significance Testing – Compare with Lomb 19 Test Case 1 : 51 points of Gaussian noise, random phases Lomb PDM
Significance Testing – Compare with Lomb 20 Test Case 2 : 511 points of Gaussian noise, random phases Conf Lev = 0.18 Monte Carlo Conf Lev = 0.19 Lomb PDM Conf Lev = 0.79 Monte Carlo Conf Lev = 0.94 Result better for PDM Much worse for Lomb
Significance Testing – PDM 21 Test Case 1: 51 points of Gaussian noise, random phases, wide bins 5/2 binning, 220 frequency points Correction is large Right curve = Beta distribution. Left curve = Beta with bandwidth correction. Arrow = data range. Noise distribution extends to Theta = 0.65 at the 0.01 significance level. Any values lower than this are not likely due to noise. Extreme value dist Right curve = Randomized Monte-Carlo derived distribution for Theta. Left curve = Distribution of Theta_min over 500 separate analyses. Results agree for both curves.
Significance Testing – PDM - 2 22 Test Case 1: 51 points of Gaussian noise, random phases, narrow bins 10/1 binning, 220 frequency points Beta changes slightly Right curve = Beta distribution. Left curve = Beta with bandwidth correction. Arrow = data range. Noise distribution extends to Theta = 0.62 at the 0.01 significance level. Right curve = Randomized Monte-Carlo derived distribution for Theta. Left curve = Distribution of Theta_min over 500 separate analyses. Results no longer agree. The narrower bins have caused a much larger noise level, extending to Theta = 0.50 at the 0.01 level. Actual
Significance Testing - Lomb 23 Test Case 1: 51 points of Gaussian noise, random phases, Lomb Correction is very large Left curve = Exp confidence distribution. Right curve = Exp with bandwidth correction. Noise distribution extends to Power = 7 at the 0.01 significance level. Arrow = data range. Power Max Dist Left curve = Randomized Monte-Carlo derived distribution for Power. Right curve = Distribution of Power_max over 500 separate analyses. Good agreement, Above noise at Pwr > 7
IC4499 - v04 24 • V Data, colored on sigma • PDM - First do a “rough cut” with nine data segments, wide frequency range. Best value is at f = 1.603. • Zero in on the main candidate with full data set. Get f = 1.60355. (published value = 1.60354) Data courtesy of G. Bono/ A. Kunder
Comparison of Beta Dist and Monte Carlo for IC4499 v04 (V data - 340 pts) Parametric Beta distribution and Corrected values 0.91 Monte Carlo distrib and Extreme value distrib 250 or 500 iterations 0.85
IC4499 - v04 26 Best period, data plotted versus phase. Colored on point number (note red points offset to right on descending branch)
IC4499 - v04 27 PDM2Variable Period Option • Redo the analysis with a range of beta, where P = P0 + b t • Agrees with O-C result (0.18 +/- 0.05, shown) • Get best alignment at about + 0.15 d/My
IC4499 - v04 28 Phase plot with period increasing at 0.20 d/My. Note much cleaner descending branch.
IC4499 – v48(Blazhko) 29 • V Data, colored on sigma • First do a “rough cut” with six data segments, wide frequency range. Best value is at f = 1.924. • Zero in on the main candidate with full data set. Get f = 1.92385. (O-C frequency = 1.92385) • Data courtesy of G. Bono/ A. Kunder
IC4499 – v48 (Blazhko) 30 Best period, data plotted versus phase Colored on point number . Note red offset to left on rising branch.
31 IC4499 – v48(Blazhko) PDM2 Variable Period Option • This gives a period change of -0.13 +/- 0.05 • Based on full light curve, so applicable to Blazhko variables.
32 IC4499 – v48(Blazhko) Phase plot with period decreasing at -0.13 d/My. Note much cleaner rising branch. Also note offsets of the max light phases.
IC4499 – v83(Blazhko) 33 • V Data, colored on sigma • First do a “rough cut” with six data segments, wide frequency range. Best value is at f = 1.5936. • Zero in on the main candidate with full data set. Get f = 1.593572. (O-C value = 1.593570) • Data courtesy of G. Bono/ A. Kunder
IC4499 – v83(Blazhko - ?) 34 Best period, data plotted versus phase Colored on point number . Doesn’t look like a Blazhko…
IC4499 – v83 (Blazhko) 35 PDM2Variable Period Option • This gives a period change of -0.20 +/- 0.10 • Based on full light curve, so applicable to Blazhko variables.
36 IC4499 – v83(Blazhko) Phase plot with period decreasing at -0.20 d/My. Cleaner rising branch. But Blazhko variations are still missing.
37 PDM2 / Lomb Application Available • Part of S_tran.exe application on www.stellingwerf.com • Typical PDM2 script shown here • #---read the data--- • read_lab data.csv // specify csv data file • for_i = 1 100000 • read_dat Xdat[i] Ydat[i] //[Sdat[i]] • if eof break • end_i • N = i-1 • #---PDM ANALYSIS--- • pdm_lpoints 3 // frequency fineness • pdm_f_range 0 4 // frequency range • pdm2 // do a PDM analysis • lomb // do a Lomb analysis • Source C code for PDM2 routine also available for download Data.csv