RobPer (version 1.2.3)

tsgen: Artificial light curve generator

Description

This function generates light curves (special time series) with unequally sampled observation times, different periodicities both in sampling and observed values, with white and power law (red) noise in the observed values and possibly disturbed observations. See RobPer-package for more information about light curves and also Thieler, Fried and Rathjens (2016) for more details in general.

Usage

tsgen(ttype, ytype, pf, redpart, s.outlier.fraction = 0, interval, npoints,
 ncycles, ps, SNR, alpha = 1.5)

Arguments

ttype

character string: Specifying the sampling pattern. Possible choices are "equi" for equidistant sampling without gaps (unperiodic), "unif" for uniform non-equidistant unperiodic sampling, "sine" for sampling with periodic sine density, "trian" for sampling with periodic triangular density, both with period \(p_s\) (see Details and sampler).

ytype

character string: Specifying the shape of the periodic fluctuation with period \(p_f\). Possible choices are "const" for constantly being zero (so no periodicity), "sine" for a sine, "trian" for a periodic triangular function, "peak" for a peak function (see Details and signalgen for more details).

pf

positive number: Period \(p_f\) of the periodic fluctuation, argument of signalgen (see Details and signalgen).

redpart

numeric value in [0,1]: Proportion of the power law noise in noise components (see Details). The generated measurement accuracies \(s_i\) do not contain information about this noise component.

s.outlier.fraction

numeric value in [0,1]: Fraction of measurement accuracies to be replaced by outliers. A value of 0 means that no measurement accuracy is replaced by an outlier (for more details see disturber).

interval

logical: If TRUE, the observed values belonging to a random time interval of length 3\(p_s\) are replaced by atypical values (for more details see disturber).

npoints

integer value: Defines the sample size \(n\) for the generated light curve.

ncycles

integer value: number \(n_s\) of sampling cycles that is observed (see Details).

ps

positive number: Sampling period \(p_s\), influencing the sampling and how the light curve is disturbed (see Details and disturber).

SNR

positive number: Defines the relation between signal \(y_f\) and noise \(y_w+y_r\) (see Details).

alpha

numeric value: Power law index \(\alpha\) for the power law noise component \(y_r\) (see Details).

Value

tt

numeric vector: Generated observation times\(\ t_1,\ldots,t_n\).

y

numeric vector: Generated observation values\(\ y_1,\ldots,y_n\).

s

numeric vector: Generated measurement accuracies\(\ s_1,\ldots,s_n\).

Details

tsgen generates an artificial light curve consisting of observation times \(t_1,\ldots,t_n\), observation values \(y_1,\ldots,y_n\) and measurement accuracies \(s_1,\ldots,s_n\). It calls several subfunctions (see there for details):

sampler is used to sample observation times \(t_1,\ldots,t_n\) in the interval \([0,n_s*p_s]\) with a possibly periodic sampling of period \(p_s\).

signalgen generates periodically varying values \(y_{f;1},\ldots,y_{f;n}\) at time points \(t_1,\ldots,t_n\) with fluctuation period \(p_f\).

lc_noise samples measurement accuracies \(s_1,\ldots,s_n\) from a Gamma(3,10)-distribution and a white noise component \(y_{w;1},\ldots,y_{w;n}\) with from \(\mathcal N(0,s_i^2)\) distributions. A second noise component \(y_{r;1},\ldots,y_{r;n}\) does not depend on the \(s_i\). It is generated as red noise, i.e. following a power law with power law index \(\alpha\). For white noise choose \(\alpha=0\), for flicker noise (pink noise) \(\alpha=1\), for brown noise \(\alpha=2\). The power law noise is generated using TK95_uneq and TK95. The noise components are scaled so that the variance of the \(y_{r;i}\) has approximately the proportion redpart in the overall noise variance and that SNR is the ratio \(var(y_f)/var(y_w+y_r)\). The observed values are set to \(y_i= y_{f;i}+y_{w;i}+y_{r;i} \forall i\).

disturber disturbes the light curve replacing measurement accuracies \(s_i\) by outliers (if s.outlier.fraction>0) and observed values \(y_i\) by atypical values (if interval=TRUE). In case of s.outlier.fraction=0 and interval=FALSE, the function returns all values unchanged.

References

Thieler, A. M., Backes, M., Fried, R. and Rhode, W. (2013): Periodicity Detection in Irregularly Sampled Light Curves by Robust Regression and Outlier Detection. Statistical Analysis and Data Mining, 6 (1), 73-89

Thieler, A. M., Fried, R. and Rathjens, J. (2016): RobPer: An R Package to Calculate Periodograms for Light Curves Based on Robust Regression. Journal of Statistical Software, 69 (9), 1-36, <doi:10.18637/jss.v069.i09>

See Also

Applies sampler, signalgen, lc_noise, disturber, TK95, TK95_uneq.

Examples

# NOT RUN {
# Generate a light curve:
set.seed(22)
lightcurve<- tsgen(ttype="sine", ytype="peak" , pf=7, redpart=0.1, s.outlier.fraction=0, 
    interval=FALSE, npoints=200, ncycles=100, ps=5, SNR=3, alpha=0)

# Or do it step by step:
# First sampling observation times:
set.seed(22)
tt <- sampler(ttype="sine", npoints=200, ncycles=100, ps=5)

# show obviously irregular observation times as grey vertical bars on a red time line:
plot(tt, tt*0, type="n", axes=FALSE, xlab="Observation Times", ylab="")
    abline(v=tt, col="grey")
axis(1, pos=0, col="red", col.axis="red")

# Sampling period is 5, look at observation times modulo 10:
hist(tt%%5, xlab="Observation time modulo 5", 
    main="Sine Distribution for the phase (tt modulo 5)", freq=FALSE)
dsin <- function(tt) 0.2*(sin(2*pi*tt/5)+1)
curve(dsin, add=TRUE)

# Then generate periodic fluctuation
yf <- signalgen(tt, ytype="peak", pf=7)

plot(tt, yf, xlab="Observation Times", ylab="Periodic Fluctuation")
plot(tt%%7, yf, main="Phase Diagram (time modulo 7)", 
    xlab="Observation time modulo 7",  ylab="Periodic Fluctuation")

# Add noise and scale signal to the right SNR
temp <- lc_noise(tt,sig=yf, SNR=3, redpart=0.1, alpha=0)
y <- temp$y
s <- temp$s

# Plotting the light curve (vertical bars show measurement accuracies)
plot(tt, y, pch=16, cex=0.5, xlab="t", ylab="y", main="a Light Curve")
rect(tt, y+s, tt, y-s)

# The lightcurve has period 7:
plot(tt%%7, y, pch=16, cex=0.5, xlab="t", ylab="y", 
    main="Phase Diagram of a Light Curve")
rect(tt%%7, y+s, tt%%7, y-s)

# replace measurement accuracies by tiny outliers or include a peak
temp <- disturber(tt,y,s,ps=5, s.outlier.fraction=0, interval=FALSE)

# Phase diagram (observation times modulo 10)
plot(tt%%7, temp$y, pch=16, cex=0.5, xlab="t", ylab="y", 
    main="Phase Diagram of a Light Curve")
rect(tt%%7, temp$y+temp$s, tt%%7, temp$y-temp$s)

# The result is the same:
all(cbind(tt,temp$y,temp$s)==lightcurve)
# }