| Title: | Lowpass Filtering |
|---|---|
| Description: | Creates lowpass filters which are commonly used in ion channel recordings. It supports generation of random numbers that are filtered, i.e. follow a model for ion channel recordings, see <doi:10.1109/TNB.2018.2845126>. Furthermore, time continuous convolutions of piecewise constant signals with the kernel of lowpass filters can be computed. |
| Authors: | Pein Florian [aut, cre], Thomas Hotz [ctb], Inder Tecuapetla-Gómez [ctb] |
| Maintainer: | Pein Florian <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0-2 |
| Built: | 2024-11-11 03:43:48 UTC |
| Source: | https://github.com/cran/lowpassFilter |
Creates lowpass filters and offers further functionalities around them. Lowpass filters are commonly used in ion channel recordings.
The main function of this package is lowpassFilter which creates lowpass filters, currently only Bessel filters are supported. randomGeneration and randomGenerationMA allow to generate random numbers that are filtered, i.e. follow a model for ion channel recordings, see (Pein et al., 2018, 2020). getConvolution, getConvolutionJump, and getConvolutionPeak allow to compute the convolution of a signal with the kernel of a lowpass filter.
Pein, F., Bartsch, A., Steinem, C., and Munk, A. (2020) Heterogeneous idealization of ion channel recordings - Open channel noise. Submitted.
Pein, F., Tecuapetla-Gómez, I., Schütte, O., Steinem, C., Munk, A. (2018) Fully-automatic multiresolution idealization for filtered ion channel recordings: flickering event detection. IEEE Trans. Nanobioscience, 17(3):300-320.
Pein, F. (2017) Heterogeneous Multiscale Change-Point Inference and its Application to Ion Channel Recordings. PhD thesis, Georg-August-Universität Göttingen. http://hdl.handle.net/11858/00-1735-0000-002E-E34A-7.
Hotz, T., Schütte, O., Sieling, H., Polupanow, T., Diederichsen, U., Steinem, C., and Munk, A. (2013) Idealizing ion channel recordings by a jump segmentation multiresolution filter. IEEE Trans. Nanobioscience, 12(4):376-386.
lowpassFilter, randomGeneration, randomGenerationMA, getConvolution, getConvolutionJump, getConvolutionPeak
# creates a lowpass filter filter <- lowpassFilter(type = "bessel", param = list(pole = 4, cutoff = 0.1), sr = 1e4) time <- 1:4000 / filter$sr # creates a piecewise constant signal with a single peak stepfun <- getSignalPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40) # computes the convolution of the signal with the kernel of the lowpass filter signal <- getConvolutionPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40, filter = filter) # generates random numbers that are filtered data <- randomGenerationMA(n = 4000, filter = filter, signal = signal, noise = 1.4) # generated data plot(time, data, pch = 16) # zoom into the single peak plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) # use of data randomGeneration instead data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = 1.4) # similar result plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2)# creates a lowpass filter filter <- lowpassFilter(type = "bessel", param = list(pole = 4, cutoff = 0.1), sr = 1e4) time <- 1:4000 / filter$sr # creates a piecewise constant signal with a single peak stepfun <- getSignalPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40) # computes the convolution of the signal with the kernel of the lowpass filter signal <- getConvolutionPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40, filter = filter) # generates random numbers that are filtered data <- randomGenerationMA(n = 4000, filter = filter, signal = signal, noise = 1.4) # generated data plot(time, data, pch = 16) # zoom into the single peak plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) # use of data randomGeneration instead data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = 1.4) # similar result plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2)
For developers only; computes a time discrete convolution.
.convolve(val, kern).convolve(val, kern)
val |
a numeric vector giving the values |
kern |
a numeric vector giving the time discrete kernel |
A numeric vector giving the convolution.
For developers only; computes the deconvolution of a single jump or an isolated peak assuming that the observations are lowpass filtered. More details are given in (Pein et al., 2018).
.deconvolveJump(grid, observations, time, leftValue, rightValue, typeFilter, inputFilter, covariances) .deconvolvePeak(gridLeft, gridRight, observations, time, leftValue, rightValue, typeFilter, inputFilter, covariances, tolerance).deconvolveJump(grid, observations, time, leftValue, rightValue, typeFilter, inputFilter, covariances) .deconvolvePeak(gridLeft, gridRight, observations, time, leftValue, rightValue, typeFilter, inputFilter, covariances, tolerance)
grid, gridLeft, gridRight
|
numeric vectors giving the potential time points of the single jump, of the left and right jump points of the peak, respectively |
observations |
a numeric vector giving the observed data |
time |
a numeric vector of length |
leftValue, rightValue
|
single numerics giving the value (conductance level) before and after the jump / peak, respectively |
typeFilter, inputFilter
|
a description of the assumed lowpass filter, usually computed by |
covariances |
a numeric vector giving the (regularized) covariances of the observations |
tolerance |
a single numeric giving a tolerance for the decision whether the left jump point is smaller than the right jump point |
For .deconvolveJump a single numeric giving the jump point. For .deconvolvePeak a list containing the entries left, right and value giving the left and right jump point and the value of the peak, respectively.
Pein, F., Tecuapetla-Gómez, I., Schütte, O., Steinem, C., Munk, A. (2018) Fully-automatic multiresolution idealization for filtered ion channel recordings: flickering event detection. IEEE Trans. Nanobioscience, 17(3):300-320.
Creates piecewise constant signals with a single jump / peak. Computes the convolution of piecewise constant signals with the kernel of a lowpass filter.
getConvolution(t, stepfun, filter, truncated = TRUE) getSignalJump(t, cp, leftValue, rightValue) getConvolutionJump(t, cp, leftValue, rightValue, filter, truncated = TRUE) getSignalPeak(t, cp1, cp2, value, leftValue, rightValue) getConvolutionPeak(t, cp1, cp2, value, leftValue, rightValue, filter, truncated = TRUE)getConvolution(t, stepfun, filter, truncated = TRUE) getSignalJump(t, cp, leftValue, rightValue) getConvolutionJump(t, cp, leftValue, rightValue, filter, truncated = TRUE) getSignalPeak(t, cp1, cp2, value, leftValue, rightValue) getConvolutionPeak(t, cp1, cp2, value, leftValue, rightValue, filter, truncated = TRUE)
t |
a numeric vector giving the time points at which the signal / convolution should be computed |
stepfun |
specification of the piecewise constant signal, i.e. a |
cp, cp1, cp2
|
a single numeric giving the location of the single, first and second jump point, respectively |
value, leftValue, rightValue
|
a single numeric giving the function value at, before and after the peak / jump, respectively |
filter |
an object of class |
truncated |
a single logical (not NA) indicating whether the signal should be convolved with the truncated or the untruncated filter kernel |
a numeric of length length(t) giving the signal / convolution at time points t
Pein, F., Bartsch, A., Steinem, C., and Munk, A. (2020) Heterogeneous idealization of ion channel recordings - Open channel noise. Submitted.
Pein, F., Tecuapetla-Gómez, I., Schütte, O., Steinem, C., Munk, A. (2018) Fully-automatic multiresolution idealization for filtered ion channel recordings: flickering event detection. IEEE Trans. Nanobioscience, 17(3):300-320.
Pein, F. (2017) Heterogeneous Multiscale Change-Point Inference and its Application to Ion Channel Recordings. PhD thesis, Georg-August-Universität Göttingen. http://hdl.handle.net/11858/00-1735-0000-002E-E34A-7.
# creating and plotting a signal with a single jump at 0 from 0 to 1 time <- seq(-2, 13, 0.01) signal <- getSignalJump(time, 0, 0, 1) plot(time, signal, type = "l") # setting up the filter filter <- lowpassFilter(param = list(pole = 4, cutoff = 0.1)) # convolution with the truncated filter convolution <- getConvolutionJump(time, 0, 0, 1, filter) lines(time, convolution, col = "red") # without truncating the filter, looks almost equal convolution <- getConvolutionJump(time, 0, 0, 1, filter, truncated = FALSE) lines(time, convolution, col = "blue") # creating and plotting a signal with a single peak with jumps # at 0 and at 3 from 0 to 1 to 0 time <- seq(-2, 16, 0.01) signal <- getSignalPeak(time, 0, 3, 1, 0, 0) plot(time, signal, type = "l") # convolution with the truncated filter convolution <- getConvolutionPeak(time, 0, 3, 1, 0, 0, filter) lines(time, convolution, col = "red") # without truncating the filter, looks almost equal convolution <- getConvolutionPeak(time, 0, 3, 1, 0, 0, filter, truncated = FALSE) lines(time, convolution, col = "blue") # doing the same with getConvolution # signal can also be an object of class stepblock instead, # e.g. constructed by stepR::stepblock signal <- data.frame(value = c(0, 1, 0), leftEnd = c(-2, 0, 3), rightEnd = c(0, 3, 16)) convolution <- getConvolution(time, signal, filter) lines(time, convolution, col = "red") convolution <- getConvolution(time, signal, filter, truncated = FALSE) lines(time, convolution, col = "blue") # more complicated signal time <- seq(-2, 21, 0.01) signal <- data.frame(value = c(0, 10, 0, 50, 0), leftEnd = c(-2, 0, 3, 6, 8), rightEnd = c(0, 3, 6, 8, 21)) convolution <- getConvolution(time, signal, filter) plot(time, convolution, col = "red", type = "l") convolution <- getConvolution(time, signal, filter, truncated = FALSE) lines(time, convolution, col = "blue")# creating and plotting a signal with a single jump at 0 from 0 to 1 time <- seq(-2, 13, 0.01) signal <- getSignalJump(time, 0, 0, 1) plot(time, signal, type = "l") # setting up the filter filter <- lowpassFilter(param = list(pole = 4, cutoff = 0.1)) # convolution with the truncated filter convolution <- getConvolutionJump(time, 0, 0, 1, filter) lines(time, convolution, col = "red") # without truncating the filter, looks almost equal convolution <- getConvolutionJump(time, 0, 0, 1, filter, truncated = FALSE) lines(time, convolution, col = "blue") # creating and plotting a signal with a single peak with jumps # at 0 and at 3 from 0 to 1 to 0 time <- seq(-2, 16, 0.01) signal <- getSignalPeak(time, 0, 3, 1, 0, 0) plot(time, signal, type = "l") # convolution with the truncated filter convolution <- getConvolutionPeak(time, 0, 3, 1, 0, 0, filter) lines(time, convolution, col = "red") # without truncating the filter, looks almost equal convolution <- getConvolutionPeak(time, 0, 3, 1, 0, 0, filter, truncated = FALSE) lines(time, convolution, col = "blue") # doing the same with getConvolution # signal can also be an object of class stepblock instead, # e.g. constructed by stepR::stepblock signal <- data.frame(value = c(0, 1, 0), leftEnd = c(-2, 0, 3), rightEnd = c(0, 3, 16)) convolution <- getConvolution(time, signal, filter) lines(time, convolution, col = "red") convolution <- getConvolution(time, signal, filter, truncated = FALSE) lines(time, convolution, col = "blue") # more complicated signal time <- seq(-2, 21, 0.01) signal <- data.frame(value = c(0, 10, 0, 50, 0), leftEnd = c(-2, 0, 3, 6, 8), rightEnd = c(0, 3, 6, 8, 21)) convolution <- getConvolution(time, signal, filter) plot(time, convolution, col = "red", type = "l") convolution <- getConvolution(time, signal, filter, truncated = FALSE) lines(time, convolution, col = "blue")
Creates a lowpass filter.
lowpassFilter(type = c("bessel"), param, sr = 1, len = NULL, shift = 0.5) ## S3 method for class 'lowpassFilter' print(x, ...)lowpassFilter(type = c("bessel"), param, sr = 1, len = NULL, shift = 0.5) ## S3 method for class 'lowpassFilter' print(x, ...)
type |
a string specifying the type of the filter, currently only Bessel filters are supported |
param |
a |
sr |
a single numeric giving the sampling rate |
len |
a single integer giving the filter length of the truncated and digitised filter, see Value for more details. By default ( |
shift |
a single numeric between |
x |
the object |
... |
for generic methods only |
An object of class lowpassFilter, i.e. a list that contains
"type", "param", "sr", "len"
the corresponding arguments
"kernfun"the kernel function of the filter, obtained as the Laplace transform of the corresponding transfer function
"stepfun"the step-response of the filter, i.e. the antiderivative of the filter kernel
"acfun"the autocorrelation function, i.e. the convolution of the filter kernel with itself
"acAntiderivative"the antiderivative of the autocorrelation function
"truncatedKernfun"the kernel function of the at len / sr truncated filter, i.e. kernfun truncated and rescaled such that the new kernel still integrates to 1
"truncatedStepfun"the step-response of the at len / sr truncated filter, i.e. the antiderivative of the kernel of the truncated filter
"truncatedAcfun"the autocorrelation function of the at len / sr truncated filter, i.e. the convolution of the kernel of the truncated filter with itself
"truncatedAcAntiderivative"the antiderivative of the autocorrelation function of the at len / sr truncated filter
"kern"the digitised filter kernel normalised to one, i.e. kernfun((0:len + shift) / sr) / sum(kernfun((0:len + shift) / sr))
"step"the digitised step-response of the filter, i.e. stepfun((0:len + shift) / sr)
"acf"the discrete autocorrelation, i.e. acfun(0:len / sr)
"jump"the last index of the left half of the filter, i.e. min(which(ret$step >= 0.5)) - 1L, it indicates how much a jump is shifted in time by a convolution of the signal with the digitised kernel of the lowpassfilter; if all values are below 0.5, len is returned with a warning
"number"for developers; an integer indicating the type of the filter
"list"for developers; a list containing precomputed quantities to recreate the filter in C++
Pein, F., Bartsch, A., Steinem, C., and Munk, A. (2020) Heterogeneous idealization of ion channel recordings - Open channel noise. Submitted.
Pein, F., Tecuapetla-Gómez, I., Schütte, O., Steinem, C., Munk, A. (2018) Fully-automatic multiresolution idealization for filtered ion channel recordings: flickering event detection. IEEE Trans. Nanobioscience, 17(3):300-320.
Pein, F. (2017) Heterogeneous Multiscale Change-Point Inference and its Application to Ion Channel Recordings. PhD thesis, Georg-August-Universität Göttingen. http://hdl.handle.net/11858/00-1735-0000-002E-E34A-7.
Hotz, T., Schütte, O., Sieling, H., Polupanow, T., Diederichsen, U., Steinem, C., and Munk, A. (2013) Idealizing ion channel recordings by a jump segmentation multiresolution filter. IEEE Trans. Nanobioscience, 12(4):376-386.
filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), sr = 1e4) # filter kernel, truncated version plot(filter$kernfun, xlim = c(0, 20 / filter$sr)) t <- seq(0, 20 / filter$sr, 0.01 / filter$sr) # truncated version looks very similar lines(t, filter$truncatedKernfun(t), col = "red") # filter$len (== 11) is chosen automatically # this ensures that filter$acf < 1e-3 for this lag and at all larger lags plot(filter$acfun, xlim = c(0, 20 / filter$sr), ylim = c(-0.003, 0.003)) abline(h = 0.001, lty = "22") abline(h = -0.001, lty = "22") abline(v = (filter$len - 1L) / filter$sr, col = "grey") abline(v = filter$len / filter$sr, col = "red") # filter with sr == 1 filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4)) # filter kernel and its truncated version plot(filter$kernfun, xlim = c(0, 20 / filter$sr)) t <- seq(0, 20 / filter$sr, 0.01 / filter$sr) # truncated version looks very similar lines(t, filter$truncatedKernfun(t), col = "red") # digitised filter points((0:filter$len + 0.5) / filter$sr, filter$kern, col = "red", pch = 16) # without a shift filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), shift = 0) # filter$kern starts with zero points(0:filter$len / filter$sr, filter$kern, col = "blue", pch = 16) # much shorter filter filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), len = 4L) points((0:filter$len + 0.5) / filter$sr, filter$kern, col = "darkgreen", pch = 16)filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), sr = 1e4) # filter kernel, truncated version plot(filter$kernfun, xlim = c(0, 20 / filter$sr)) t <- seq(0, 20 / filter$sr, 0.01 / filter$sr) # truncated version looks very similar lines(t, filter$truncatedKernfun(t), col = "red") # filter$len (== 11) is chosen automatically # this ensures that filter$acf < 1e-3 for this lag and at all larger lags plot(filter$acfun, xlim = c(0, 20 / filter$sr), ylim = c(-0.003, 0.003)) abline(h = 0.001, lty = "22") abline(h = -0.001, lty = "22") abline(v = (filter$len - 1L) / filter$sr, col = "grey") abline(v = filter$len / filter$sr, col = "red") # filter with sr == 1 filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4)) # filter kernel and its truncated version plot(filter$kernfun, xlim = c(0, 20 / filter$sr)) t <- seq(0, 20 / filter$sr, 0.01 / filter$sr) # truncated version looks very similar lines(t, filter$truncatedKernfun(t), col = "red") # digitised filter points((0:filter$len + 0.5) / filter$sr, filter$kern, col = "red", pch = 16) # without a shift filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), shift = 0) # filter$kern starts with zero points(0:filter$len / filter$sr, filter$kern, col = "blue", pch = 16) # much shorter filter filter <- lowpassFilter(type = "bessel", param = list(pole = 4L, cutoff = 1e3 / 1e4), len = 4L) points((0:filter$len + 0.5) / filter$sr, filter$kern, col = "darkgreen", pch = 16)
Generate random numbers that are filtered. Both, signal and noise, are convolved with the given lowpass filter, see details. Can be used to generate synthetic data resembling ion channel recordings, please see (Pein et al., 2018, 2020) for the exact models.
randomGeneration(n, filter, signal = 0, noise = 1, oversampling = 100L, seed = n, startTime = 0, truncated = TRUE) randomGenerationMA(n, filter, signal = 0, noise = 1, seed = n, startTime = 0, truncated = TRUE)randomGeneration(n, filter, signal = 0, noise = 1, oversampling = 100L, seed = n, startTime = 0, truncated = TRUE) randomGenerationMA(n, filter, signal = 0, noise = 1, seed = n, startTime = 0, truncated = TRUE)
n |
a single positive integer giving the number of observations that should be generated |
filter |
an object of class |
signal |
either a numeric of length 1 or of length |
noise |
for |
oversampling |
a single positive integer giving the factor by which the errors should be oversampled, see Details |
seed |
will be passed to |
startTime |
a single finite numeric giving the time at which sampling should start |
truncated |
a single logical (not NA) indicating whether the signal should be convolved with the truncated or the untruncated filter kernel |
As discussed in (Pein et al., 2018) and (Pein et al., 2020), in ion channel recordings the recorded data points can be modelled as equidistant sampled at rate filter$sr from the convolution of a piecewise constant signal perturbed by Gaussian white noise scaled by the noise level with the kernel of an analogue lowpass filter. The noise level is either constant (homogeneous noise, see (Pein et al., 2018)) or itself varying (heterogeneous noise, see (Pein et al., 2020)). randomGeneration and randomGenerationMA generate synthetic data from such models. randomGeneration allows homogeneous and heterogeneous noise, while randomGenerationMA only allows homogeneous noise, i.e. noise has to be a single numeric giving the constant noise level. The resulting observations represent the conductance at time points startTime + 1:n / filter$sr.
The generated observations are the sum of a convolved signal evaluated at those time points plus centred Gaussian errors that are correlated (coloured noise), because of the filtering, and scaled by the noise level. The convolved signal evaluated at those time points can either by specified in signal directly or signal can specify a piecewise constant signal that will be convolved with the filter using getConvolution and evaluated at those time points. randomGenerationMA computes a moving average process with the desired autocorrelation to generate random errors. randomGeneration oversamples the error, i.e. generates errors at time points startTime + (seq(1 - filter$len + 1 / oversampling, n, 1 / oversampling) - 1 / 2 / oversampling) / filter$sr, which will then be convolved with the filter. For this function noise can either give the noise levels at those oversampled time points or specify a piecewise constant function that will be automatically evaluated at those time points.
a numeric vector of length n giving the generated random numbers
Pein, F., Bartsch, A., Steinem, C., and Munk, A. (2020) Heterogeneous idealization of ion channel recordings - Open channel noise. Submitted.
Pein, F., Tecuapetla-Gómez, I., Schütte, O., Steinem, C., Munk, A. (2018) Fully-automatic multiresolution idealization for filtered ion channel recordings: flickering event detection. IEEE Trans. Nanobioscience, 17(3):300-320.
Pein, F. (2017) Heterogeneous Multiscale Change-Point Inference and its Application to Ion Channel Recordings. PhD thesis, Georg-August-Universität Göttingen. http://hdl.handle.net/11858/00-1735-0000-002E-E34A-7.
filter <- lowpassFilter(type = "bessel", param = list(pole = 4, cutoff = 0.1), sr = 1e4) time <- 1:4000 / filter$sr stepfun <- getSignalPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40) signal <- getConvolutionPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40, filter = filter) data <- randomGenerationMA(n = 4000, filter = filter, signal = signal, noise = 1.4) # generated data plot(time, data, pch = 16) # zoom into the single peak plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) # use of randomGeneration instead data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = 1.4) # similar result plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) ## heterogeneous noise # manual creation of an object of class 'stepblock' # instead the function stepblock in the package stepR can be used noise <- data.frame(leftEnd = c(0, 0.2, 0.2 + 3 / filter$sr), rightEnd = c(0.2, 0.2 + 3 / filter$sr, 0.4), value = c(1, 30, 1)) attr(noise, "x0") <- 0 class(noise) <- c("stepblock", class(noise)) data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = noise) plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2)filter <- lowpassFilter(type = "bessel", param = list(pole = 4, cutoff = 0.1), sr = 1e4) time <- 1:4000 / filter$sr stepfun <- getSignalPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40) signal <- getConvolutionPeak(time, cp1 = 0.2, cp2 = 0.2 + 3 / filter$sr, value = 20, leftValue = 40, rightValue = 40, filter = filter) data <- randomGenerationMA(n = 4000, filter = filter, signal = signal, noise = 1.4) # generated data plot(time, data, pch = 16) # zoom into the single peak plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) # use of randomGeneration instead data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = 1.4) # similar result plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2) ## heterogeneous noise # manual creation of an object of class 'stepblock' # instead the function stepblock in the package stepR can be used noise <- data.frame(leftEnd = c(0, 0.2, 0.2 + 3 / filter$sr), rightEnd = c(0.2, 0.2 + 3 / filter$sr, 0.4), value = c(1, 30, 1)) attr(noise, "x0") <- 0 class(noise) <- c("stepblock", class(noise)) data <- randomGeneration(n = 4000, filter = filter, signal = signal, noise = noise) plot(time, data, pch = 16, xlim = c(0.199, 0.202), ylim = c(19, 45)) lines(time, stepfun, col = "blue", type = "s", lwd = 2) lines(time, signal, col = "red", lwd = 2)