Description Usage Arguments Value References Examples
Calculates the Lomb Scargle Periodogram using the algorithm described by K. Hocke and N. Kämpfer (2009).
1 2 |
t |
Numeric vector of timepoints. |
y |
Numeric vector of values corresponding to t. |
ofac |
Oversampling factor. Recommend >=2. Higher improves granularity, but may cause artefacts. |
mc_sim |
monte carlo simulation object returned by function "monte_carlo_lsp". If not NULL will return p-values. If NULL will use spd for thresholding; using the "thresh" parameter |
zero_factor |
Integer, zero padding factor. Pads fourier series to increase resolution of approximated function. Number of zeros is equal to zero_factor*length of (unpadded) fourier series. Value of 0 means no padding. |
thresh |
Threshhold for p-value or Spectral Power Density. Filter fourier series based on peak strength p value if mc_sim not NULL, otherise on spd. This will affect the reconstructed series. |
verbose |
Verbosity level. 0=silent, 1=printed results. |
A list with the following elements:
t - original t.
y - original y.
t2 - approximated t.
y2 - reconstructed y. plot against t2.
frequency - frequency range used.
spectral_power_density - Plot against Frequency for the periodogram. unitless.
phase_radian - phase relative to 0 in radians. range 0:2*pi.
phase_default - phase relative to 0 in the same unit as t. This is the distance from 0 to the nearest peak after t[1].
peak_info - data.frame with index, frequency, spd, phase_default, and (optional) pvalue of each detected peak.
peak_idx - index of highest peak
peak_spd - highest spd
peak_period - period corresponding to peak_idx
peak_phase - phase_default corresponding to peak_idx
Hocke, K., and N. Kämpfer. "Gap filling and noise reduction of unevenly sampled data by means of the Lomb-Scargle periodogram." Atmospheric Chemistry and Physics 9.12 (2009): 4197-4206.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 | library(LSP)
#set up testing parameters
ofac = 8;
amp=3;
sd = 1.0;
n = 100; #initial series length
n_out = 50; #sereis length after random sampling
period = 45;
freq = 1/period;
offset = 25;
phase = pi/2;
p_thresh = 0.05;
N_mc_runs = 1000;
zero_factor = 2;
N_cores = 1;
#create testing data
x = seq(0,n-1);
x = x + offset;
y = amp*cos(2*pi*freq*x - phase);
#random sampling to simulate missing data/irregular timepoints
idx = seq(1,n);
idx = sample(idx,n_out);
idx = sort(idx)
#add some noise
y = rnorm(n,0,sd) + y;
#final data
y = y[idx];
x = x[idx];
#monte carlo
mc_sim = monte_carlo_lsp(x,N_mc_runs,ofac,N_cores);
#LSP
# result = lsp(x,y,ofac=ofac,mc_sim = mc_sim,zero_factor = zero_factor, thresh = p_thresh,verbose=1);
result = lsp(x,y,ofac=ofac,mc_sim = mc_sim,zero_factor = zero_factor, thresh = 0.05,verbose=1);
#plot input data vs reconstructed signal
plot(result$t,result$y,col="red")
lines(result$t2,result$y2,col="green")
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.