Wavelet Analysis for Trading: Multi-Scale Decomposition
Financial markets operate across multiple time scales simultaneously—high-frequency noise, daily fluctuations, weekly trends, monthly cycles, and longer-term regime shifts all coexist in price data. Wavelet analysis provides a powerful mathematical framework for decomposing price series into these constituent scales, revealing structure invisible to traditional time-domain or frequency-domain methods alone. Unlike Fourier analysis which assumes stationary sinusoids, wavelets are localized in both time and frequency—perfect for analyzing the non-stationary, multi-scale nature of financial markets.Understanding Wavelet Transforms
A wavelet is a wave-like oscillation that is localized in time. The wavelet transform decomposes a signal into wavelets scaled and shifted across time.Continuous Wavelet Transform (CWT)
The CWT of signal x(t) is defined as: ``
W(a,b) = (1/√a) ∫ x(t) ψ*((t-b)/a) dt
`
Where:
- a: scale parameter (inversely related to frequency)
- b: translation parameter (time location)
- ψ: mother wavelet function
- ψ*: complex conjugate of mother wavelet
The scale parameter a dilates (a > 1) or compresses (a < 1) the wavelet, corresponding to low or high frequency analysis.
Discrete Wavelet Transform (DWT)
DWT uses discrete scales and translations:
`
W(j,k) = (1/√2^j) Σ x[n] ψ((n-2^j k)/2^j)
`
Where j is the scale level and k is the translation index.
DWT decomposes signals into approximation (low-frequency) and detail (high-frequency) coefficients:
`
x[n] = Σ cA_j φ_j,k + Σ cD_j ψ_j,k
`
Key Takeaways
- Wavelets analyze signals across multiple time scales simultaneously
- Localized in both time and frequency (vs Fourier: frequency only)
- Ideal for non-stationary financial data with regime changes
- Enables multi-resolution analysis from intraday to long-term trends
- Can separate noise from signal at appropriate scales
- Reveals time-varying cyclical patterns invisible to other methods
Why Wavelets Matter for Trading
1. Multi-Scale Market Structure
Markets exhibit fractal-like behavior—patterns repeat across time scales. Wavelets naturally capture this:
`python
import pywt
import numpy as np
def multi_scale_decomposition(prices, wavelet='db4', levels=5):
"""
Decompose price series into multiple scales
"""
# Calculate returns
returns = np.log(prices / prices.shift(1)).dropna()
# Discrete wavelet decomposition
coeffs = pywt.wavedec(returns, wavelet, level=levels)
# Coefficients: [cA_n, cD_n, cD_n-1, ..., cD_1]
approximation = coeffs[0] # Trend
details = coeffs[1:] # Details at each scale
return approximation, details
Example usage
prices = data['close']
approx, details = multi_scale_decomposition(prices, levels=5)
Level 1: ~2-4 day cycles (highest frequency)
Level 2: ~4-8 day cycles
Level 3: ~8-16 day cycles
Level 4: ~16-32 day cycles
Level 5: ~32-64 day cycles
Approximation: trend (>64 days)
`
2. Noise Filtration
Separate signal from noise by reconstructing only significant scales:
`python
def denoise_signal(prices, wavelet='db4', level=3, threshold_factor=1.0):
"""
Denoise price series using wavelet thresholding
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Decompose
coeffs = pywt.wavedec(returns, wavelet, level=level)
# Threshold detail coefficients (soft thresholding)
sigma = np.median(np.abs(coeffs[-1])) / 0.6745 # Noise estimate
threshold = threshold_factor sigma np.sqrt(2 * np.log(len(returns)))
coeffs_thresh = coeffs.copy()
for i in range(1, len(coeffs)):
coeffs_thresh[i] = pywt.threshold(coeffs[i], threshold, mode='soft')
# Reconstruct
denoised = pywt.waverec(coeffs_thresh, wavelet)
# Trim to original length
return denoised[:len(returns)]
`
3. Time-Varying Cycles
Identify cycles that appear, disappear, and shift over time:
`python
import matplotlib.pyplot as plt
def plot_wavelet_scalogram(prices, wavelet='cmor1.5-1.0', scales=None):
"""
Create scalogram showing time-frequency energy distribution
"""
returns = np.log(prices / prices.shift(1)).dropna()
if scales is None:
# Auto-generate scales (2^1 to 2^8 ≈ 2 to 256 days)
scales = 2 ** np.arange(1, 9, 0.25)
# Continuous wavelet transform
coefficients, frequencies = pywt.cwt(returns, scales, wavelet)
power = np.abs(coefficients) ** 2
# Plot
plt.figure(figsize=(15, 8))
plt.imshow(power, extent=[0, len(returns), scales[-1], scales[0]],
cmap='jet', aspect='auto', vmax=np.percentile(power, 95))
plt.colorbar(label='Power')
plt.ylabel('Scale (approximate period in days)')
plt.xlabel('Time')
plt.title('Wavelet Scalogram - Time-Frequency Power')
return coefficients, frequencies, power
`
Practical Trading Applications
Trend-Following at Multiple Scales
Trade trends identified at different wavelet scales:
`python
class MultiScaleTrendStrategy:
def __init__(self, wavelet='db4', levels=[2, 3, 4]):
self.wavelet = wavelet
self.levels = levels
def generate_signals(self, prices):
"""
Generate trend signals at multiple scales
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Decompose
max_level = max(self.levels)
coeffs = pywt.wavedec(returns, self.wavelet, level=max_level)
signals = {}
# For each scale of interest
for level in self.levels:
# Reconstruct signal at this scale
# Zero out all other scales
coeffs_single = [np.zeros_like(c) for c in coeffs]
coeffs_single[-(level)] = coeffs[-(level)]
reconstructed = pywt.waverec(coeffs_single, self.wavelet)
reconstructed = reconstructed[:len(returns)]
# Signal: sign of wavelet coefficient
signal = np.sign(reconstructed)
period = 2 ** level
signals[f'scale_{period}d'] = signal
return pd.DataFrame(signals, index=returns.index)
def combined_signal(self, prices, weights=None):
"""
Combine signals across scales
"""
signals = self.generate_signals(prices)
if weights is None:
# Equal weight
weights = np.ones(len(self.levels)) / len(self.levels)
combined = np.zeros(len(signals))
for i, col in enumerate(signals.columns):
combined += weights[i] * signals[col].values
return pd.Series(combined, index=signals.index)
`
Mean Reversion at High Frequencies
Trade mean reversion in high-frequency components:
`python
def wavelet_mean_reversion(prices, wavelet='db4', detail_level=1, z_threshold=2.0):
"""
Mean reversion strategy using high-frequency wavelet details
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Decompose
coeffs = pywt.wavedec(returns, wavelet, level=5)
# Extract high-frequency detail
detail_coeffs = coeffs[-detail_level]
# Reconstruct only this detail
coeffs_single = [np.zeros_like(c) for c in coeffs]
coeffs_single[-detail_level] = detail_coeffs
hf_component = pywt.waverec(coeffs_single, wavelet)
hf_component = hf_component[:len(returns)]
# Z-score of high-frequency component
rolling_mean = pd.Series(hf_component).rolling(20).mean()
rolling_std = pd.Series(hf_component).rolling(20).std()
z_score = (pd.Series(hf_component) - rolling_mean) / rolling_std
# Mean reversion signals
signals = np.zeros(len(returns))
signals[z_score > z_threshold] = -1 # Short when too high
signals[z_score < -z_threshold] = 1 # Long when too low
return pd.Series(signals, index=returns.index)
`
Volatility Estimation
Estimate volatility at different time scales:
`python
def multi_scale_volatility(prices, wavelet='db4', levels=5):
"""
Compute volatility contribution from each scale
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Decompose
coeffs = pywt.wavedec(returns, wavelet, level=levels)
# Variance at each scale
total_var = np.var(returns)
scale_vars = {}
# Approximation variance
approx_var = np.var(coeffs[0])
scale_vars['trend'] = approx_var / total_var
# Detail variances
for i, detail in enumerate(coeffs[1:]):
scale = 2 ** (levels - i)
detail_var = np.var(detail)
scale_vars[f'{scale}d'] = detail_var / total_var
return scale_vars
Example
vol_decomp = multi_scale_volatility(prices)
Shows what % of volatility comes from each time scale
`
Regime Change Detection
Detect regime changes using wavelet variance:
`python
def detect_regime_change(prices, wavelet='db4', level=3, window=60):
"""
Detect regime changes via time-varying wavelet variance
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Rolling wavelet variance
variance_series = []
for i in range(window, len(returns)):
window_data = returns.iloc[i-window:i]
coeffs = pywt.wavedec(window_data, wavelet, level=level)
# Total wavelet variance
total_var = sum(np.var(c) for c in coeffs)
variance_series.append(total_var)
variance_series = pd.Series(variance_series, index=returns.index[window:])
# Regime change when variance shifts significantly
var_change = variance_series.pct_change()
regime_changes = np.abs(var_change) > 2 * var_change.std()
return regime_changes, variance_series
`
Advanced Techniques
Wavelet Packet Decomposition
Decompose both approximation and details for finer resolution:
`python
def wavelet_packet_analysis(prices, wavelet='db4', maxlevel=4):
"""
Full wavelet packet decomposition
"""
returns = np.log(prices / prices.shift(1)).dropna()
# Create wavelet packet tree
wp = pywt.WaveletPacket(data=returns, wavelet=wavelet, maxlevel=maxlevel)
# Extract all nodes at max level
nodes = wp.get_level(maxlevel, order='freq')
# Frequency band analysis
freq_bands = {}
for i, node in enumerate(nodes):
freq_bands[f'band_{i}'] = {
'coefficients': node.data,
'energy': np.sum(node.data ** 2),
'path': node.path
}
return freq_bands
Usage
bands = wavelet_packet_analysis(prices, maxlevel=4)
16 frequency bands at level 4
`
Complex Wavelets for Phase Analysis
Use complex wavelets to extract phase information:
`python
def wavelet_phase_analysis(prices, scales=None):
"""
Extract phase using complex Morlet wavelet
"""
returns = np.log(prices / prices.shift(1)).dropna()
if scales is None:
scales = np.arange(2, 128)
# Complex Morlet wavelet
coefficients, frequencies = pywt.cwt(
returns,
scales,
'cmor1.5-1.0' # Complex Morlet
)
# Extract magnitude and phase
magnitude = np.abs(coefficients)
phase = np.angle(coefficients)
# Instantaneous frequency (phase derivative)
inst_freq = np.diff(phase, axis=1) / (2 * np.pi)
return magnitude, phase, inst_freq
`
Multi-Dimensional Wavelet Analysis
Analyze correlations across assets using 2D wavelets:
`python
def cross_asset_wavelet_correlation(prices1, prices2, wavelet='db4', level=4):
"""
Wavelet coherence between two assets
"""
returns1 = np.log(prices1 / prices1.shift(1)).dropna()
returns2 = np.log(prices2 / prices2.shift(1)).dropna()
# Align series
common_idx = returns1.index.intersection(returns2.index)
returns1 = returns1.loc[common_idx]
returns2 = returns2.loc[common_idx]
# Decompose both
coeffs1 = pywt.wavedec(returns1, wavelet, level=level)
coeffs2 = pywt.wavedec(returns2, wavelet, level=level)
# Correlation at each scale
correlations = {}
for i in range(len(coeffs1)):
scale_name = 'approx' if i == 0 else f'detail_{level-i+1}'
corr = np.corrcoef(coeffs1[i], coeffs2[i])[0, 1]
correlations[scale_name] = corr
return correlations
`
Adaptive Wavelets
Choose optimal wavelet for specific data:
`python
def select_optimal_wavelet(prices, wavelet_families=['db', 'sym', 'coif']):
"""
Select wavelet with best entropy criterion
"""
from scipy.stats import entropy
returns = np.log(prices / prices.shift(1)).dropna()
best_wavelet = None
best_entropy = float('inf')
for family in wavelet_families:
wavelets = pywt.wavelist(family)
for wavelet in wavelets:
try:
coeffs = pywt.wavedec(returns, wavelet, level=4)
# Compute entropy of coefficients
all_coeffs = np.concatenate(coeffs)
ent = entropy(np.abs(all_coeffs))
if ent < best_entropy:
best_entropy = ent
best_wavelet = wavelet
except:
continue
return best_wavelet, best_entropy
`
Implementation Best Practices
1. Boundary Effects
Wavelets suffer from edge effects. Handle carefully:
`python
def handle_boundary_effects(signal, wavelet='db4', level=4, mode='periodic'):
"""
Mitigate boundary effects
modes: 'periodic', 'symmetric', 'reflect', 'zero'
"""
# Use appropriate boundary mode
coeffs = pywt.wavedec(signal, wavelet, level=level, mode=mode)
# Alternatively, pad signal
pad_length = 2 ** (level + 1)
padded = np.pad(signal, (pad_length, pad_length), mode='reflect')
coeffs_padded = pywt.wavedec(padded, wavelet, level=level)
reconstructed = pywt.waverec(coeffs_padded, wavelet)
# Remove padding
reconstructed = reconstructed[pad_length:-pad_length]
return reconstructed
`
2. Choosing Wavelet Type
Different wavelets for different applications:
`python
Daubechies (db4, db8): Good general purpose, smooth
Symlets (sym4, sym8): Nearly symmetric, good for financial data
Coiflets (coif1-coif5): Symmetric, good for trend extraction
Morlet: Complex wavelet for phase analysis
Mexican hat: Good for peak detection
def wavelet_characteristics(wavelet_name):
"""
Get properties of a wavelet
"""
wavelet = pywt.Wavelet(wavelet_name)
return {
'family': wavelet.family_name,
'orthogonal': wavelet.orthogonal,
'biorthogonal': wavelet.biorthogonal,
'symmetry': wavelet.symmetry,
'compact_support': wavelet.vanishing_moments_psi > 0,
'filter_length': len(wavelet.dec_lo)
}
`
3. Coefficient Thresholding
Apply thresholds to reduce noise:
`python
def threshold_wavelet_coeffs(coeffs, threshold_type='soft', threshold=None):
"""
Apply thresholding to wavelet coefficients
threshold_type: 'soft', 'hard', 'garrote'
"""
if threshold is None:
# Universal threshold
n = sum(len(c) for c in coeffs)
sigma = np.median(np.abs(coeffs[-1])) / 0.6745
threshold = sigma np.sqrt(2 np.log(n))
coeffs_thresh = coeffs.copy()
for i in range(1, len(coeffs)): # Skip approximation
if threshold_type == 'soft':
coeffs_thresh[i] = pywt.threshold(coeffs[i], threshold, mode='soft')
elif threshold_type == 'hard':
coeffs_thresh[i] = pywt.threshold(coeffs[i], threshold, mode='hard')
elif threshold_type == 'garrote':
coeffs_thresh[i] = pywt.threshold(coeffs[i], threshold, mode='garrote')
return coeffs_thresh
`
Real-World Case Study
Complete Multi-Scale Trading System
`python
class WaveletTradingSystem:
def __init__(self, wavelet='db4', trend_levels=[3, 4, 5],
mr_level=1, vol_window=60):
self.wavelet = wavelet
self.trend_levels = trend_levels
self.mr_level = mr_level
self.vol_window = vol_window
def decompose(self, prices):
"""Wavelet decomposition"""
returns = np.log(prices / prices.shift(1)).dropna()
max_level = max(self.trend_levels)
coeffs = pywt.wavedec(returns, self.wavelet, level=max_level)
return returns, coeffs
def trend_signals(self, coeffs):
"""Extract trend signals"""
signals = []
for level in self.trend_levels:
coeffs_single = [np.zeros_like(c) for c in coeffs]
coeffs_single[-(level)] = coeffs[-(level)]
reconstructed = pywt.waverec(coeffs_single, self.wavelet)
signal = np.sign(reconstructed)
signals.append(signal)
# Average across scales
combined = np.mean(signals, axis=0)
return combined
def mean_reversion_signals(self, coeffs, returns):
"""High-frequency mean reversion"""
# Extract HF detail
coeffs_hf = [np.zeros_like(c) for c in coeffs]
coeffs_hf[-self.mr_level] = coeffs[-self.mr_level]
hf = pywt.waverec(coeffs_hf, self.wavelet)
hf = hf[:len(returns)]
# Z-score
hf_series = pd.Series(hf, index=returns.index)
z = (hf_series - hf_series.rolling(20).mean()) / hf_series.rolling(20).std()
signals = np.zeros(len(returns))
signals[z > 2] = -1
signals[z < -2] = 1
return signals
def volatility_filter(self, coeffs):
"""Compute current volatility regime"""
# Total wavelet variance
total_var = sum(np.var(c) for c in coeffs)
# Rolling average variance (would need history)
# Simplified: return current variance
return np.sqrt(total_var)
def generate_signals(self, prices):
"""Combined signal generation"""
returns, coeffs = self.decompose(prices)
# Components
trend_sig = self.trend_signals(coeffs)
mr_sig = self.mean_reversion_signals(coeffs, returns)
vol = self.volatility_filter(coeffs)
# Combine: trend in low vol, mean reversion in high vol
vol_threshold = 0.02 # 2% daily vol
if vol < vol_threshold:
# Low vol: follow trends
combined = trend_sig[:len(returns)]
else:
# High vol: mean reversion
combined = mr_sig
return pd.Series(combined, index=returns.index)
def backtest(self, prices):
"""Backtest the system"""
signals = self.generate_signals(prices)
returns = prices.pct_change()
# Align
common_idx = signals.index.intersection(returns.index)
signals = signals.loc[common_idx]
returns = returns.loc[common_idx]
# Strategy returns
strategy_returns = signals.shift(1) * returns
return strategy_returns.dropna()
Usage
system = WaveletTradingSystem(wavelet='db4', trend_levels=[3, 4, 5])
strategy_returns = system.backtest(prices)
Performance
sharpe = strategy_returns.mean() / strategy_returns.std() * np.sqrt(252)
cumulative = (1 + strategy_returns).cumprod()
print(f"Sharpe: {sharpe:.2f}")
``