| Title: | High-Performance Wavelet Lifting Transforms |
|---|---|
| Description: | Performs Wavelet Lifting Transforms focusing on signal denoising and functional data analysis (FDA). Implements a hybrid architecture with a zero-allocation 'C++' core for high-performance processing. Features include unified offline (batch) denoising, causal (real-time) filtering using a ring buffer engine, and adaptive recursive thresholding. |
| Authors: | Moises da Silva [aut, cre] |
| Maintainer: | Moises da Silva <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.9.1 |
| Built: | 2026-06-03 08:34:22 UTC |
| Source: | https://github.com/mkyou/rlifting |
Comparison of execution time between rLifting's optimized causal mode and
a naive sliding-window implementation using wavethresh.
data(benchmark_causal)data(benchmark_causal)
A list containing:
Average time (seconds) for rLifting.
Time (seconds) for naive sliding window.
Ratio of Naive Time to rLifting Time.
Comparison of execution time and reconstruction error (MSE) between rLifting and other packages
(wavethresh, wavelets) using Haar wavelet.
data(benchmark_offline)data(benchmark_offline)
A data frame with the following columns:
Package name.
Execution time in seconds.
Mean Squared Error.
Estimates the optimal noise threshold based on current window statistics. Implements the recursive formula from Liu et al. (2014). Accelerated with 'C++'.
compute_adaptive_threshold(lwt_obj, alpha = 0.3, beta = 1.2)compute_adaptive_threshold(lwt_obj, alpha = 0.3, beta = 1.2)
lwt_obj |
Object returned by |
alpha |
Recursive adjustment parameter (Eq. 9). |
beta |
Initial threshold scale factor (Eq. 9). |
Object of class adaptive_thresholds (a list of thresholds).
Wrapper to create a lifting_scheme object from manual steps.
custom_wavelet(name, steps, norm = c(1, 1))custom_wavelet(name, steps, norm = c(1, 1))
name |
Identifier name for the wavelet. |
steps |
List of steps created via |
norm |
Normalization vector c(K, 1/K). |
An object of class lifting_scheme.
p1 = lift_step("predict", c(1), position = "center") u1 = lift_step("update", c(0.5), position = "center") w = custom_wavelet("HaarManual", list(p1, u1), c(1.41, 0.707))p1 = lift_step("predict", c(1), position = "center") u1 = lift_step("update", c(0.5), position = "center") w = custom_wavelet("HaarManual", list(p1, u1), c(1.41, 0.707))
Processes a complete signal simulating the sequential arrival of data.
Uses the specialized 'C++' class WaveletEngine to perform causal
filtering efficiently on a historical dataset.
denoise_signal_causal( signal, scheme, levels = 1, window_size = 256, alpha = 0.3, beta = 1.2, method = "semisoft", extension = "symmetric", update_freq = 1 )denoise_signal_causal( signal, scheme, levels = 1, window_size = 256, alpha = 0.3, beta = 1.2, method = "semisoft", extension = "symmetric", update_freq = 1 )
signal |
Complete vector of the noisy signal. |
scheme |
|
levels |
Decomposition levels. |
window_size |
Window size. |
alpha |
Threshold decay parameter (Eq 9). |
beta |
Threshold gain factor (Eq 9). |
method |
Thresholding method ("soft", "hard", "semisoft"). |
extension |
Boundary treatment ('symmetric', 'periodic'). |
update_freq |
Frequency of threshold updates. |
Filtered vector (same length as input).
Performs denoising on the entire signal at once using a non-causal approach. Uses global statistics for recursive threshold calculation (Eq. 9). This function is fully optimized in 'C++' (Zero-Allocation).
denoise_signal_offline( signal, scheme, alpha = 0.3, beta = 1.2, levels = 3, method = "semisoft", extension = "symmetric" )denoise_signal_offline( signal, scheme, alpha = 0.3, beta = 1.2, levels = 3, method = "semisoft", extension = "symmetric" )
signal |
Numeric vector containing the complete signal. |
scheme |
A |
alpha |
Recursive threshold parameter. |
beta |
Threshold scale factor. |
levels |
Number of decomposition levels. |
method |
Thresholding method ("hard", "soft", "semisoft"). |
extension |
Extension mode ("symmetric", "periodic", "zero"). |
Filtered numeric vector (same length as input).
Runs a battery of physical and mathematical tests on a wavelet.
diagnose_wavelet(wavelet_name, config, verbose = TRUE, plot = TRUE)diagnose_wavelet(wavelet_name, config, verbose = TRUE, plot = TRUE)
wavelet_name |
Name string or a |
config |
Configuration list (is_ortho, vm_degrees, max_taps). |
verbose |
Print results to console handling? (Defaults to TRUE). |
plot |
Boolean. Visualize basis functions during diagnosis? (Defaults to TRUE). |
An object of class wavelet_diagnosis (S3), which is a list containing
the results of each test. The object has a dedicated print method.
A synthetic dataset containing a Doppler signal contaminated with Gaussian noise. Used in the "General Usage" vignette.
data(doppler_example)data(doppler_example)
A data frame with 2048 rows and 3 columns:
Time index.
The pure Doppler signal.
The signal with added Gaussian noise (sd=0.5).
Reconstructs the original signal from wavelet coefficients. Optimized with 'C++' backend.
ilwt(lwt_obj, scheme = NULL)ilwt(lwt_obj, scheme = NULL)
lwt_obj |
Object of class 'lwt' returned by |
scheme |
(Optional) |
Numeric vector containing the reconstructed signal.
s = c(1, 2, 3, 4) sch = lifting_scheme("haar") fwd = lwt(s, sch) rec = ilwt(fwd) print(rec) # Should match ss = c(1, 2, 3, 4) sch = lifting_scheme("haar") fwd = lwt(s, sch) rec = ilwt(fwd) print(rec) # Should match s
Measurement of energy leakage into the "past" when processing an impulse signal. Used to demonstrate the zero-lookahead property of the causal mode.
data(leakage_results)data(leakage_results)
A data frame with:
Method description (e.g. "rLifting causal").
Sum of squared differences (leakage energy).
Helper function to create prediction (P) or update (U) steps, abstracting the complexity of index management.
lift_step( type = c("predict", "update"), coeffs, start_idx = NULL, position = "center" )lift_step( type = c("predict", "update"), coeffs, start_idx = NULL, position = "center" )
type |
Step type: "predict" (P) or "update" (U). |
coeffs |
Numeric vector containing the filter coefficients. |
start_idx |
(Optional) Manual start index. If provided, ignores the
|
position |
Automatic index adjustment
(used only if
|
A list formatted for the internal lifting engine.
Creates an S3 object containing the prediction (P) and update (U) steps required for the Lifting Transform.
lifting_scheme(wavelet = "haar", custom_steps = NULL, custom_norm = NULL)lifting_scheme(wavelet = "haar", custom_steps = NULL, custom_norm = NULL)
wavelet |
Wavelet name (string). Options: "haar", "db2", "cdf53", "cdf97", "dd4", "lazy". |
custom_steps |
List of custom steps (optional). If provided, ignores internal lookup. |
custom_norm |
Normalization vector (optional). |
An object of class lifting_scheme.
Performs the Forward Wavelet Transform using the Lifting Scheme. Optimized with 'C++' backend.
lwt(signal, scheme, levels = 1, extension = "symmetric")lwt(signal, scheme, levels = 1, extension = "symmetric")
signal |
Numeric vector containing the input signal. |
scheme |
A |
levels |
Integer. Number of decomposition levels. |
extension |
Boundary extension mode: "symmetric" (default), "periodic", or "zero". |
An object of class 'lwt'. It is a list containing 'coeffs' (list of details d1..dn and approximation an) and 'scheme' (the scheme object used).
data = c(1, 2, 3, 4, 5, 6, 7, 8) sch = lifting_scheme("haar") res = lwt(data, sch, levels = 2) print(res)data = c(1, 2, 3, 4, 5, 6, 7, 8) sch = lifting_scheme("haar") res = lwt(data, sch, levels = 2) print(res)
Generates a stateful function backed by a high-performance 'C++' Ring Buffer engine. It implements Sliding Window + Lifting Decomposition + Adaptive Thresholding in highly efficient time per sample.
new_wavelet_stream( scheme, window_size = 256, levels = 1, alpha = 0.3, beta = 1.2, method = "semisoft", extension = "symmetric", update_freq = 1 )new_wavelet_stream( scheme, window_size = 256, levels = 1, alpha = 0.3, beta = 1.2, method = "semisoft", extension = "symmetric", update_freq = 1 )
scheme |
A |
window_size |
Sliding window size (W). Must be > 8. |
levels |
Decomposition levels (default 1). |
alpha |
Threshold decay parameter (Eq 9). |
beta |
Threshold gain factor (Eq 9). |
method |
Shrinkage method: "hard", "soft", "semisoft". |
extension |
Boundary handling ('symmetric', 'periodic', 'zero'). |
update_freq |
How often to recompute threshold statistics (default 1). |
A closure function processor(new_sample) that accepts
a single numeric value and returns the filtered value immediately.
Plot method for Adaptive Thresholds
## S3 method for class 'adaptive_thresholds' plot(x, ...)## S3 method for class 'adaptive_thresholds' plot(x, ...)
x |
Object of class |
... |
Additional arguments. |
Invisibly returns NULL.
Plot method for Lifting Scheme
## S3 method for class 'lifting_scheme' plot(x, ...)## S3 method for class 'lifting_scheme' plot(x, ...)
x |
An object of class |
... |
Additional arguments passed to |
Invisibly returns NULL.
Plot method for LWT Decomposition
## S3 method for class 'lwt' plot(x, ...)## S3 method for class 'lwt' plot(x, ...)
x |
An object of class |
... |
Additional arguments. |
Invisibly returns NULL.
Print method for Adaptive Thresholds
## S3 method for class 'adaptive_thresholds' print(x, ...)## S3 method for class 'adaptive_thresholds' print(x, ...)
x |
Object of class |
... |
Additional arguments. |
Invisibly returns x.
Print method
## S3 method for class 'lifting_scheme' print(x, ...)## S3 method for class 'lifting_scheme' print(x, ...)
x |
object of class |
... |
additional arguments. |
Invisibly returns NULL. Called for side effects (printing).
Print method for LWT
## S3 method for class 'lwt' print(x, ...)## S3 method for class 'lwt' print(x, ...)
x |
An object of class lwt. |
... |
Additional arguments. |
Invisibly returns NULL. Called for side effects (printing).
Print method for Wavelet Diagnosis
## S3 method for class 'wavelet_diagnosis' print(x, ...)## S3 method for class 'wavelet_diagnosis' print(x, ...)
x |
Object of class |
... |
Additional arguments. |
Invisibly returns x.
Print method for Wavelet Stream Processor
## S3 method for class 'wavelet_stream' print(x, ...)## S3 method for class 'wavelet_stream' print(x, ...)
x |
Object of class |
... |
Additional arguments. |
Invisibly returns x.
A unified framework for Wavelet Transforms using the Lifting Scheme. It provides robust tools for offline signal analysis and functional data analysis (FDA), while also enabling high-performance causal processing for real-time applications via a specialized 'C++' core.
Maintainer: Moises da Silva [email protected]
Useful links:
General Thresholding Wrapper
threshold(x, lambda, method = "soft")threshold(x, lambda, method = "soft")
x |
Input vector. |
lambda |
Threshold value. |
method |
Method: "hard", "soft" or "semisoft". |
Numeric vector of the same length as x with thresholded coefficients.
Sets coefficients below the threshold to zero, keeping others unchanged. Known as the "keep or kill" policy.
threshold_hard(x, lambda)threshold_hard(x, lambda)
x |
Vector of coefficients (details). |
lambda |
Positive threshold value. |
Processed vector.
Implementation based on Liu et al. (2014).
Combines the stability of Soft Thresholding with the amplitude
precision of Hard Thresholding.
Function: sign(x) * sqrt(x^2 - lambda^2) for values above lambda.
threshold_semisoft(x, lambda)threshold_semisoft(x, lambda)
x |
Vector of coefficients. |
lambda |
Positive threshold value. |
Processed vector.
Sets coefficients below the threshold to zero and shrinks others towards zero. Reduces noise variance but introduces amplitude bias.
threshold_soft(x, lambda)threshold_soft(x, lambda)
x |
Vector of coefficients. |
lambda |
Positive threshold value. |
Processed vector.
Verifies if the impulse response is finite (FIR Filter).
validate_compact_support(scheme, max_width)validate_compact_support(scheme, max_width)
scheme |
Object of class |
max_width |
Maximum expected width (number of taps). |
List with status and number of active taps.
Verifies Parseval's Theorem. Only applicable for orthogonal wavelets.
validate_orthogonality(scheme, expected = TRUE, tol = 1e-09)validate_orthogonality(scheme, expected = TRUE, tol = 1e-09)
scheme |
Object of class |
expected |
Boolean. If TRUE, expects orthogonality. |
tol |
Tolerance (default 1e-9). |
List with status and energy ratio (Out/In).
Verifies wavelet invertibility against a battery of signals.
validate_perfect_reconstruction(scheme, tol = 1e-09)validate_perfect_reconstruction(scheme, tol = 1e-09)
scheme |
Object of class |
tol |
Numerical error tolerance (default 1e-9). |
List with global status and maximum error found.
Decimated wavelets are not translation invariant. This test quantifies the variation in detail energy when shifting the input signal by 1 sample.
validate_shift_sensitivity(scheme)validate_shift_sensitivity(scheme)
scheme |
Object of class |
List with status and percentage variation.
Verifies if the wavelet cancels polynomials of a specific degree.
validate_vanishing_moments(scheme, degree = 0, tol = 1e-09)validate_vanishing_moments(scheme, degree = 0, tol = 1e-09)
scheme |
Object of class |
degree |
Polynomial degree (0=Constant, 1=Ramp, 2=Parabola...). |
tol |
Residual energy tolerance (default 1e-9). |
List with status and residual energy.
Plots the waveform by iterating the reconstruction over several levels.
visualize_wavelet_basis(scheme, plot = TRUE, levels = 8)visualize_wavelet_basis(scheme, plot = TRUE, levels = 8)
scheme |
Object of class |
plot |
Boolean. |
levels |
Number of cascade levels. |
Invisibly returns NULL. Called for side effects (plotting).