feat(timeseries): add time-series integration helpers for financial analysis
- Implement 6 helper functions in src/timeseries_utils.rs: * prepare_for_hmm: Feature engineering for HMM regime detection * rolling_hurst_exponent: Mean-reversion detection (H < 0.5 = mean-reverting) * rolling_half_life: Mean-reversion speed for pairs trading * return_statistics: Risk metrics (mean, std, skew, kurt, sharpe) * create_lagged_features: ML feature matrix creation * rolling_correlation: Rolling correlation for pairs trading - Add PyO3 bindings in src/timeseries_utils/python_bindings.rs: * All functions exposed with _py suffix * Proper signature decorators and error handling * Registered in lib.rs module system - Update Python module exports: * python/optimizr/core.py: Import from _core * python/optimizr/__init__.py: Re-export all functions - Create comprehensive example: * examples/timeseries_integration.py demonstrates all 6 functions * Includes integrated pairs trading workflow * Shows feature engineering for regime detection - Technical details: * Fixed Array1<f64> type conversions for ndarray compatibility * Uses risk_metrics::hurst_exponent and estimate_half_life * Built successfully with maturin develop --release (40.93s) * All functions tested and working correctly Part of Priority 3: Time-series integration helpers (Enhancement Strategy) Addresses v0.3.0 roadmap: Bridge optimization with time-series analysis
This commit is contained in:
Melvin Alvarez
@@ -26,6 +26,11 @@ RUN apt-get update && apt-get install -y \
|
||||
build-essential \
|
||||
curl \
|
||||
git \
|
||||
pkg-config \
|
||||
libssl-dev \
|
||||
libopenblas-dev \
|
||||
gfortran \
|
||||
patchelf \
|
||||
&& rm -rf /var/lib/apt/lists/*
|
||||
|
||||
# Install Rust (needed for maturin)
|
||||
|
||||
@@ -0,0 +1,407 @@
|
||||
# OptimizR Enhancement Strategy
|
||||
|
||||
**Date**: January 2, 2025
|
||||
**Context**: Post-Polaroid Phase 4, exploring integration and improvements
|
||||
**Based On**: v0.2.0 codebase review, roadmap analysis, synergy opportunities
|
||||
|
||||
## Current State Analysis
|
||||
|
||||
### ✅ What's Implemented (v0.2.0)
|
||||
|
||||
1. **Core Algorithms**:
|
||||
- Differential Evolution (5 strategies: rand1, best1, currenttobest1, rand2, best2)
|
||||
- Hidden Markov Models (Baum-Welch, Viterbi)
|
||||
- MCMC Sampling (Metropolis-Hastings, adaptive proposals)
|
||||
- Grid Search
|
||||
- Information Theory (mutual information, Shannon entropy)
|
||||
|
||||
2. **Advanced Features (v0.2.0)**:
|
||||
- Sparse Optimization (Sparse PCA, Box-Tao, Elastic Net)
|
||||
- Optimal Control (HJB solver, regime switching, jump diffusion)
|
||||
- Risk Metrics (Hurst exponent, half-life, bootstrap)
|
||||
- Mathematical Toolkit (numerical differentiation, linear algebra, statistics)
|
||||
|
||||
3. **Architecture**:
|
||||
- Trait-based design (Optimizer, Sampler, InformationMeasure)
|
||||
- Builder pattern for configuration
|
||||
- Functional programming utilities (composition, memoization, pipes)
|
||||
- Rayon dependency already present
|
||||
- Feature flag infrastructure (`parallel` feature exists)
|
||||
|
||||
### ⚠️ What's Missing/Incomplete
|
||||
|
||||
1. **Parallelization BLOCKED**:
|
||||
- Infrastructure exists (Rayon trait, ParallelExecutor trait in core.rs)
|
||||
- DE has `parallel` parameter but **disabled** due to Python GIL
|
||||
- Comment: "Python callbacks cannot be safely parallelized due to GIL"
|
||||
- Grid search marked as "future: Expected 50-100x speedup"
|
||||
|
||||
2. **Advanced DE Variants (Roadmap v0.3.0)**:
|
||||
- JADE (jDE with archive)
|
||||
- SHADE (Success-History based Adaptive DE)
|
||||
- L-SHADE (with linear population reduction)
|
||||
- Current: Only basic jDE adaptive control
|
||||
|
||||
3. **Multi-Objective Optimization (Roadmap)**:
|
||||
- NSGA-DE (Non-dominated Sorting)
|
||||
- MODE (Multi-Objective DE)
|
||||
- Pareto front computation
|
||||
|
||||
4. **GPU Acceleration (Roadmap)**:
|
||||
- CUDA kernels
|
||||
- OpenCL support
|
||||
- 10-100× additional speedup
|
||||
|
||||
5. **Additional Algorithms (Roadmap)**:
|
||||
- Particle Swarm Optimization (PSO)
|
||||
- CMA-ES (Covariance Matrix Adaptation)
|
||||
- Simulated Annealing
|
||||
- Ant Colony Optimization
|
||||
|
||||
## Synergy Opportunities: Polaroid + OptimizR
|
||||
|
||||
### 1. Time-Series Feature Engineering for HMM
|
||||
**Description**: Use Polaroid's time-series operations to create features for regime detection
|
||||
|
||||
**Implementation**:
|
||||
```python
|
||||
# Polaroid: Fast feature creation
|
||||
df = client.lag(['price'], periods=1) # Lagged prices
|
||||
df = client.pct_change(['price'], periods=1) # Returns
|
||||
df = client.diff(['price'], periods=1) # Price changes
|
||||
|
||||
# OptimizR: Regime detection on features
|
||||
returns = df['price_pct_change'].to_numpy()
|
||||
hmm = HMM(n_states=3) # Bull, Bear, Sideways
|
||||
hmm.fit(returns, n_iterations=100)
|
||||
states = hmm.predict(returns)
|
||||
```
|
||||
|
||||
**Value**:
|
||||
- Polaroid provides fast feature engineering (50-200× faster for large datasets)
|
||||
- OptimizR provides statistical inference (HMM regime detection)
|
||||
- Combined: Real-time regime switching for trading strategies
|
||||
|
||||
### 2. Risk Metrics on Time-Series Data
|
||||
**Description**: Calculate advanced risk metrics using both systems
|
||||
|
||||
**Implementation**:
|
||||
```python
|
||||
# Polaroid: Efficient return calculation
|
||||
df = client.pct_change(['price'], periods=1)
|
||||
returns = df['price_pct_change'].to_numpy()
|
||||
|
||||
# OptimizR: Risk analysis
|
||||
hurst = compute_hurst_exponent(returns) # Mean-reversion detection
|
||||
half_life = estimate_half_life(returns) # Reversion time
|
||||
risk_metrics = compute_risk_metrics(returns) # Comprehensive suite
|
||||
```
|
||||
|
||||
**Value**:
|
||||
- Fast preprocessing (Polaroid) + sophisticated analysis (OptimizR)
|
||||
- Useful for pairs trading, mean-reversion strategies
|
||||
- Real-time risk monitoring
|
||||
|
||||
### 3. Optimal Control with Market Data
|
||||
**Description**: Dynamic portfolio rebalancing with regime-dependent strategies
|
||||
|
||||
**Implementation**:
|
||||
```python
|
||||
# Polaroid: Multi-asset feature creation
|
||||
df = client.lag(['spy_price', 'vix'], periods=[1, 5, 20])
|
||||
df = client.pct_change(['spy_price'], periods=1)
|
||||
|
||||
# OptimizR: Solve optimal control problem
|
||||
# State: [price, volatility regime]
|
||||
# Control: portfolio weights
|
||||
value_fn = solve_hjb_regime_switching(...)
|
||||
```
|
||||
|
||||
**Value**:
|
||||
- Combines fast data processing with optimal control theory
|
||||
- Regime-dependent strategies (bull vs bear market)
|
||||
- Practical for HFT and algorithmic trading
|
||||
|
||||
### 4. Parameter Optimization for Trading Strategies
|
||||
**Description**: Use DE to optimize strategy parameters on time-series data
|
||||
|
||||
**Implementation**:
|
||||
```python
|
||||
# Polaroid: Backtest execution (fast data ops)
|
||||
def backtest_strategy(params):
|
||||
df = client.lag(['price'], periods=int(params[0]))
|
||||
# ... strategy logic ...
|
||||
return -sharpe_ratio # Minimize negative Sharpe
|
||||
|
||||
# OptimizR: Find optimal parameters
|
||||
result = differential_evolution(
|
||||
objective_fn=backtest_strategy,
|
||||
bounds=[(1, 50), (0.01, 0.5)], # [lag_period, threshold]
|
||||
maxiter=500,
|
||||
strategy='rand1'
|
||||
)
|
||||
```
|
||||
|
||||
**Value**:
|
||||
- Polaroid handles heavy data processing
|
||||
- OptimizR finds optimal parameters
|
||||
- 74-88× faster than SciPy DE
|
||||
|
||||
## High-Priority Enhancements
|
||||
|
||||
### Priority 1: Enable Parallelization for Pure-Rust Objectives
|
||||
|
||||
**Problem**: `parallel` parameter exists but disabled due to Python GIL issues
|
||||
|
||||
**Solution**: Create Rust-native objective function trait for GIL-free parallelization
|
||||
|
||||
**Implementation Strategy**:
|
||||
1. Add `RustObjectiveFn` trait separate from Python callbacks
|
||||
2. Implement parallel evaluation for Rust-native functions
|
||||
3. Keep Python callbacks sequential (GIL limitation)
|
||||
4. Enable parallel grid search (no Python callbacks needed for grid)
|
||||
|
||||
**Code Outline**:
|
||||
```rust
|
||||
// In src/core.rs or src/differential_evolution.rs
|
||||
|
||||
/// Rust-native objective function (no Python, no GIL)
|
||||
pub trait RustObjective: Send + Sync {
|
||||
fn evaluate(&self, x: &[f64]) -> f64;
|
||||
}
|
||||
|
||||
/// Parallel evaluation for Rust objectives
|
||||
#[cfg(feature = "parallel")]
|
||||
fn evaluate_population_parallel<F: RustObjective>(
|
||||
objective: &F,
|
||||
population: &[Vec<f64>]
|
||||
) -> Vec<f64> {
|
||||
use rayon::prelude::*;
|
||||
population.par_iter()
|
||||
.map(|individual| objective.evaluate(individual))
|
||||
.collect()
|
||||
}
|
||||
|
||||
// Python binding for benchmarking
|
||||
#[pyfunction]
|
||||
fn differential_evolution_rust(
|
||||
objective_name: &str, // "sphere", "rosenbrock", "rastrigin"
|
||||
bounds: Vec<(f64, f64)>,
|
||||
parallel: bool, // Now actually works!
|
||||
...
|
||||
) -> PyResult<DEResult>
|
||||
```
|
||||
|
||||
**Benefits**:
|
||||
- 10-100× speedup for built-in test functions (sphere, Rosenbrock, Rastrigin)
|
||||
- Useful for benchmarking and testing
|
||||
- Grid search can be parallelized (no callbacks)
|
||||
- Foundation for future Rust-only mode
|
||||
|
||||
**Effort**: Medium (1-2 hours)
|
||||
|
||||
### Priority 2: Implement SHADE (Success-History Adaptive DE)
|
||||
|
||||
**Problem**: Current adaptive DE uses basic jDE, SHADE is state-of-the-art
|
||||
|
||||
**Solution**: Implement SHADE algorithm from Tanabe & Fukunaga (2013)
|
||||
|
||||
**Key Features**:
|
||||
- Historical memory of successful parameters (F, CR)
|
||||
- Weighted random selection from memory
|
||||
- Better than jDE on CEC benchmarks
|
||||
|
||||
**Implementation Strategy**:
|
||||
1. Add `SHADE` variant to `DEStrategy` enum
|
||||
2. Create success history buffer (circular buffer of size H=10-100)
|
||||
3. Update memory after each successful mutation
|
||||
4. Sample (F, CR) from history using Cauchy/Normal distributions
|
||||
|
||||
**Code Outline**:
|
||||
```rust
|
||||
pub enum DEAdaptive {
|
||||
None,
|
||||
JDE, // Current implementation
|
||||
SHADE, // New: Success-history based
|
||||
LSHADE, // Future: With linear population reduction
|
||||
}
|
||||
|
||||
struct SHADEMemory {
|
||||
history_f: Vec<f64>, // Successful F values
|
||||
history_cr: Vec<f64>, // Successful CR values
|
||||
index: usize, // Circular buffer index
|
||||
size: usize, // Memory size H
|
||||
}
|
||||
|
||||
impl SHADEMemory {
|
||||
fn sample_f(&self) -> f64 {
|
||||
// Cauchy distribution centered on random history entry
|
||||
}
|
||||
|
||||
fn sample_cr(&self) -> f64 {
|
||||
// Normal distribution centered on random history entry
|
||||
}
|
||||
|
||||
fn update(&mut self, successful_f: f64, successful_cr: f64) {
|
||||
// Add to circular buffer
|
||||
}
|
||||
}
|
||||
```
|
||||
|
||||
**Benefits**:
|
||||
- State-of-the-art adaptive control
|
||||
- Better than jDE empirically
|
||||
- Aligns with roadmap (v0.3.0)
|
||||
- Minimal API changes
|
||||
|
||||
**Effort**: Medium-High (2-4 hours with testing)
|
||||
|
||||
### Priority 3: Time-Series Integration Helpers
|
||||
|
||||
**Problem**: Using Polaroid + OptimizR requires manual glue code
|
||||
|
||||
**Solution**: Create helper functions for common time-series + optimization patterns
|
||||
|
||||
**Implementation Strategy**:
|
||||
1. Add `timeseries_utils` module to OptimizR
|
||||
2. Functions for common workflows
|
||||
3. Optional Polaroid integration (via feature flag)
|
||||
|
||||
**Code Outline**:
|
||||
```rust
|
||||
// New module: src/timeseries_utils.rs
|
||||
|
||||
/// Prepare time-series data for HMM regime detection
|
||||
pub fn prepare_for_hmm(
|
||||
prices: &[f64],
|
||||
lag_periods: &[usize],
|
||||
) -> Vec<Vec<f64>> {
|
||||
// Create features: returns, lagged returns, etc.
|
||||
}
|
||||
|
||||
/// Rolling window risk metrics
|
||||
pub fn rolling_hurst_exponent(
|
||||
returns: &[f64],
|
||||
window_size: usize,
|
||||
) -> Vec<f64> {
|
||||
// Compute Hurst exponent in rolling windows
|
||||
}
|
||||
|
||||
/// Backtest parameter optimization
|
||||
pub fn optimize_strategy_params<F>(
|
||||
objective_fn: F,
|
||||
param_bounds: Vec<(f64, f64)>,
|
||||
n_trials: usize,
|
||||
) -> DEResult
|
||||
where F: Fn(&[f64]) -> f64
|
||||
{
|
||||
// Wrapper around DE with sensible defaults
|
||||
}
|
||||
```
|
||||
|
||||
**Python Bindings**:
|
||||
```python
|
||||
from optimizr import timeseries_utils as tsu
|
||||
|
||||
# Prepare features
|
||||
features = tsu.prepare_for_hmm(prices, lag_periods=[1, 5, 20])
|
||||
|
||||
# Rolling risk metrics
|
||||
rolling_hurst = tsu.rolling_hurst_exponent(returns, window_size=252)
|
||||
|
||||
# Strategy optimization
|
||||
def my_strategy(params):
|
||||
# ... backtesting logic ...
|
||||
return sharpe_ratio
|
||||
|
||||
result = tsu.optimize_strategy_params(
|
||||
my_strategy,
|
||||
param_bounds=[(1, 50), (0.01, 0.5)],
|
||||
n_trials=500
|
||||
)
|
||||
```
|
||||
|
||||
**Benefits**:
|
||||
- Reduces boilerplate for common use cases
|
||||
- Makes integration obvious
|
||||
- Encourages adoption
|
||||
- Low effort, high value
|
||||
|
||||
**Effort**: Low-Medium (1-2 hours)
|
||||
|
||||
## Secondary Enhancements (Future Work)
|
||||
|
||||
### 4. Multi-Objective Optimization (NSGA-DE)
|
||||
- **Roadmap**: v0.3.0
|
||||
- **Use Case**: Portfolio optimization (maximize return, minimize risk)
|
||||
- **Effort**: High (5-8 hours)
|
||||
|
||||
### 5. GPU Acceleration
|
||||
- **Roadmap**: v0.3.0
|
||||
- **Use Case**: Massive population sizes (10K-100K individuals)
|
||||
- **Effort**: Very High (multi-day project)
|
||||
|
||||
### 6. Additional Algorithms (PSO, CMA-ES, etc.)
|
||||
- **Roadmap**: v0.3.0
|
||||
- **Use Case**: Algorithm portfolio for different problem types
|
||||
- **Effort**: High per algorithm (3-5 hours each)
|
||||
|
||||
## Recommended Implementation Order
|
||||
|
||||
1. **Session 1 (Current)**: Time-Series Integration Helpers (1-2 hours)
|
||||
- Low effort, immediate value
|
||||
- Makes Polaroid + OptimizR integration obvious
|
||||
- Creates examples for documentation
|
||||
|
||||
2. **Session 2**: Enable Rust-Native Parallelization (1-2 hours)
|
||||
- Unblocks major performance gain
|
||||
- Grid search parallelization
|
||||
- Foundation for future work
|
||||
|
||||
3. **Session 3**: Implement SHADE (2-4 hours)
|
||||
- State-of-the-art adaptive DE
|
||||
- Aligns with roadmap
|
||||
- Publishable improvement
|
||||
|
||||
4. **Future**: Multi-objective, GPU, additional algorithms
|
||||
- Larger projects
|
||||
- Requires more research
|
||||
|
||||
## Testing Strategy
|
||||
|
||||
For each enhancement:
|
||||
1. **Unit tests**: Algorithm correctness (sphere function, Rosenbrock)
|
||||
2. **Benchmarks**: Performance comparison (before/after)
|
||||
3. **Integration tests**: Polaroid + OptimizR workflows
|
||||
4. **Documentation**: Usage examples, API docs
|
||||
|
||||
## Git Commit Strategy (per MANDATORY rules)
|
||||
|
||||
Each enhancement gets:
|
||||
1. Feature branch: `feature/shade-algorithm` or `feature/rust-parallelization`
|
||||
2. Implementation commits with tests
|
||||
3. Benchmark results documented
|
||||
4. Final commit: `feat(de): implement SHADE adaptive DE variant`
|
||||
5. Push to origin
|
||||
6. Log to historia/
|
||||
|
||||
## Success Metrics
|
||||
|
||||
1. **Performance**:
|
||||
- Rust parallelization: 10-100× speedup on multi-core
|
||||
- SHADE: 10-20% better convergence than jDE on benchmarks
|
||||
- Time-series helpers: Zero overhead (pure convenience)
|
||||
|
||||
2. **Usability**:
|
||||
- Integration examples in documentation
|
||||
- Clear API documentation
|
||||
- Python usage examples
|
||||
|
||||
3. **Completeness**:
|
||||
- All tests passing
|
||||
- Benchmarks documented
|
||||
- Changes committed to git
|
||||
|
||||
---
|
||||
|
||||
**Next Action**: Implement Priority 3 (Time-Series Integration Helpers) as it's lowest effort with immediate value for demonstrating Polaroid + OptimizR synergy.
|
||||
@@ -0,0 +1,256 @@
|
||||
"""
|
||||
Time-Series Integration Helpers - Example Usage
|
||||
===============================================
|
||||
|
||||
Demonstrates the 6 time-series utility functions for financial data analysis:
|
||||
1. prepare_for_hmm_py: Feature engineering for regime detection
|
||||
2. rolling_hurst_exponent_py: Mean-reversion detection
|
||||
3. rolling_half_life_py: Pairs trading metrics
|
||||
4. return_statistics_py: Risk analysis
|
||||
5. create_lagged_features_py: ML feature creation
|
||||
6. rolling_correlation_py: Correlation analysis
|
||||
|
||||
These helpers bridge OptimizR's optimization capabilities with time-series analysis,
|
||||
particularly useful for regime-switching models and pairs trading strategies.
|
||||
"""
|
||||
|
||||
import optimizr
|
||||
import numpy as np
|
||||
|
||||
|
||||
def example_prepare_for_hmm():
|
||||
"""Example: Feature engineering for Hidden Markov Models"""
|
||||
print("\n=== Example 1: prepare_for_hmm_py ===")
|
||||
|
||||
# Simulate stock prices
|
||||
prices = [100.0, 101.5, 99.8, 102.3, 103.7, 104.2, 103.1, 105.8, 107.2, 106.5]
|
||||
|
||||
# Create feature matrix with 1 and 2-period lags
|
||||
features = optimizr.prepare_for_hmm_py(prices, [1, 2])
|
||||
|
||||
print(f"Input: {len(prices)} price points")
|
||||
print(f"Output: {len(features)} rows x {len(features[0])} columns")
|
||||
print("\nFeature columns:")
|
||||
print(" [0] Simple returns")
|
||||
print(" [1] Log returns")
|
||||
print(" [2] Volatility proxy (squared returns)")
|
||||
print(" [3] Lagged returns (lag=1)")
|
||||
print(" [4] Lagged returns (lag=2)")
|
||||
print(f"\nFirst row: {[f'{x:.4f}' for x in features[0]]}")
|
||||
print("\n💡 Use this with OptimizR's HMM for regime detection!")
|
||||
|
||||
|
||||
def example_rolling_hurst():
|
||||
"""Example: Detecting mean-reversion with Hurst exponent"""
|
||||
print("\n=== Example 2: rolling_hurst_exponent_py ===")
|
||||
|
||||
# Generate mean-reverting returns
|
||||
np.random.seed(42)
|
||||
returns = list(np.random.randn(20) * 0.02)
|
||||
|
||||
# Compute rolling Hurst exponent
|
||||
window = 10
|
||||
hurst_values = optimizr.rolling_hurst_exponent_py(returns, window)
|
||||
|
||||
print(f"Returns: {len(returns)} observations")
|
||||
print(f"Rolling Hurst (window={window}): {len(hurst_values)} values")
|
||||
print(f"\nHurst values: {[f'{h:.3f}' for h in hurst_values[:5]]}...")
|
||||
print("\nInterpretation:")
|
||||
print(" H < 0.5: Mean-reverting (good for pairs trading)")
|
||||
print(" H = 0.5: Random walk")
|
||||
print(" H > 0.5: Trending")
|
||||
|
||||
avg_hurst = np.mean(hurst_values)
|
||||
if avg_hurst < 0.5:
|
||||
print(f"\n📊 Average H = {avg_hurst:.3f} → Mean-reverting behavior detected!")
|
||||
elif avg_hurst > 0.5:
|
||||
print(f"\n📊 Average H = {avg_hurst:.3f} → Trending behavior detected!")
|
||||
else:
|
||||
print(f"\n📊 Average H = {avg_hurst:.3f} → Random walk behavior")
|
||||
|
||||
|
||||
def example_rolling_half_life():
|
||||
"""Example: Mean-reversion speed for pairs trading"""
|
||||
print("\n=== Example 3: rolling_half_life_py ===")
|
||||
|
||||
# Simulate spread between two cointegrated assets
|
||||
np.random.seed(42)
|
||||
spread = list(100 + np.cumsum(np.random.randn(30) * 0.5))
|
||||
|
||||
# Compute rolling half-life
|
||||
window = 15
|
||||
half_lives = optimizr.rolling_half_life_py(spread, window)
|
||||
|
||||
print(f"Spread: {len(spread)} observations")
|
||||
print(f"Rolling half-life (window={window}): {len(half_lives)} values")
|
||||
print(f"\nHalf-life values: {[f'{hl:.2f}' for hl in half_lives[:5]]}...")
|
||||
print("\nInterpretation:")
|
||||
print(" Lower half-life → Faster mean reversion")
|
||||
print(" Higher half-life → Slower mean reversion")
|
||||
|
||||
avg_hl = np.mean(half_lives)
|
||||
print(f"\n📊 Average half-life: {avg_hl:.2f} periods")
|
||||
print(f" → Spread reverts to mean in ~{avg_hl:.0f} periods on average")
|
||||
|
||||
|
||||
def example_return_statistics():
|
||||
"""Example: Comprehensive risk metrics"""
|
||||
print("\n=== Example 4: return_statistics_py ===")
|
||||
|
||||
# Sample returns from a trading strategy
|
||||
returns = [0.02, -0.01, 0.015, 0.025, -0.005, 0.01, -0.02, 0.03, 0.005, -0.015]
|
||||
|
||||
# Compute statistics
|
||||
mean, std, skew, kurt, sharpe = optimizr.return_statistics_py(returns)
|
||||
|
||||
print(f"Returns: {len(returns)} observations")
|
||||
print("\nStatistics:")
|
||||
print(f" Mean return: {mean:.4f} ({mean*100:.2f}%)")
|
||||
print(f" Volatility (std): {std:.4f}")
|
||||
print(f" Skewness: {skew:.4f} {'(left-tailed)' if skew < 0 else '(right-tailed)'}")
|
||||
print(f" Kurtosis: {kurt:.4f} {'(fat tails)' if kurt > 0 else '(thin tails)'}")
|
||||
print(f" Sharpe ratio: {sharpe:.4f}")
|
||||
|
||||
print("\n📊 Risk Assessment:")
|
||||
if sharpe > 2.0:
|
||||
print(" ✅ Excellent risk-adjusted returns")
|
||||
elif sharpe > 1.0:
|
||||
print(" ✓ Good risk-adjusted returns")
|
||||
else:
|
||||
print(" ⚠️ Moderate risk-adjusted returns")
|
||||
|
||||
|
||||
def example_lagged_features():
|
||||
"""Example: Create features for ML models"""
|
||||
print("\n=== Example 5: create_lagged_features_py ===")
|
||||
|
||||
# Time series to predict
|
||||
returns = [0.01, 0.02, -0.01, 0.015, 0.005, -0.005, 0.025, 0.01, -0.01, 0.02]
|
||||
|
||||
# Create lagged feature matrix
|
||||
lags = [1, 2, 3]
|
||||
features = optimizr.create_lagged_features_py(returns, lags, include_original=True)
|
||||
|
||||
print(f"Original series: {len(returns)} observations")
|
||||
print(f"Lagged features: {len(features)} rows x {len(features[0])} columns")
|
||||
print("\nFeature columns:")
|
||||
print(f" [0] Original value (t)")
|
||||
print(f" [1] Lag-1 (t-1)")
|
||||
print(f" [2] Lag-2 (t-2)")
|
||||
print(f" [3] Lag-3 (t-3)")
|
||||
print(f"\nFirst row: {[f'{x:.4f}' for x in features[0]]}")
|
||||
print("\n💡 Use this for ML prediction models (LSTM, Random Forest, etc.)")
|
||||
|
||||
|
||||
def example_rolling_correlation():
|
||||
"""Example: Pairs trading correlation analysis"""
|
||||
print("\n=== Example 6: rolling_correlation_py ===")
|
||||
|
||||
# Two potentially cointegrated assets
|
||||
np.random.seed(42)
|
||||
asset1_returns = list(np.random.randn(25) * 0.02)
|
||||
asset2_returns = list(np.random.randn(25) * 0.02 + np.array(asset1_returns) * 0.6)
|
||||
|
||||
# Compute rolling correlation
|
||||
window = 10
|
||||
correlations = optimizr.rolling_correlation_py(asset1_returns, asset2_returns, window)
|
||||
|
||||
print(f"Asset 1 returns: {len(asset1_returns)} observations")
|
||||
print(f"Asset 2 returns: {len(asset2_returns)} observations")
|
||||
print(f"Rolling correlation (window={window}): {len(correlations)} values")
|
||||
print(f"\nCorrelation values: {[f'{c:.3f}' for c in correlations[:5]]}...")
|
||||
|
||||
avg_corr = np.mean(correlations)
|
||||
print(f"\n📊 Average correlation: {avg_corr:.3f}")
|
||||
if avg_corr > 0.7:
|
||||
print(" → Strong positive correlation (good for pairs trading)")
|
||||
elif avg_corr > 0.3:
|
||||
print(" → Moderate correlation")
|
||||
else:
|
||||
print(" → Weak correlation (not ideal for pairs trading)")
|
||||
|
||||
|
||||
def example_integrated_workflow():
|
||||
"""Example: Complete pairs trading analysis workflow"""
|
||||
print("\n" + "="*70)
|
||||
print("=== Integrated Workflow: Pairs Trading Analysis ===")
|
||||
print("="*70)
|
||||
|
||||
# Generate synthetic pair of assets
|
||||
np.random.seed(42)
|
||||
n = 50
|
||||
asset1 = list(100 + np.cumsum(np.random.randn(n) * 0.5))
|
||||
asset2 = list(100 + np.cumsum(np.random.randn(n) * 0.5 +
|
||||
(np.array(asset1) - 100) * 0.6))
|
||||
|
||||
# Compute spread
|
||||
spread = [a1 - a2 for a1, a2 in zip(asset1, asset2)]
|
||||
|
||||
# Step 1: Check mean-reversion with Hurst exponent
|
||||
print("\n1. Mean-reversion check (Hurst exponent):")
|
||||
hurst_values = optimizr.rolling_hurst_exponent_py(spread, 20)
|
||||
avg_hurst = np.mean(hurst_values)
|
||||
print(f" Average Hurst: {avg_hurst:.3f}")
|
||||
mean_reverting = avg_hurst < 0.5
|
||||
print(f" Mean-reverting: {'✅ Yes' if mean_reverting else '❌ No'}")
|
||||
|
||||
# Step 2: Estimate reversion speed
|
||||
print("\n2. Mean-reversion speed (half-life):")
|
||||
half_lives = optimizr.rolling_half_life_py(spread, 20)
|
||||
avg_hl = np.mean(half_lives)
|
||||
print(f" Average half-life: {avg_hl:.2f} periods")
|
||||
print(f" Reversion time: ~{avg_hl:.0f} periods")
|
||||
|
||||
# Step 3: Check correlation stability
|
||||
print("\n3. Correlation stability:")
|
||||
returns1 = [(asset1[i] - asset1[i-1])/asset1[i-1] for i in range(1, len(asset1))]
|
||||
returns2 = [(asset2[i] - asset2[i-1])/asset2[i-1] for i in range(1, len(asset2))]
|
||||
correlations = optimizr.rolling_correlation_py(returns1, returns2, 15)
|
||||
avg_corr = np.mean(correlations)
|
||||
print(f" Average correlation: {avg_corr:.3f}")
|
||||
print(f" Correlation stability: {'✅ High' if avg_corr > 0.7 else '⚠️ Moderate' if avg_corr > 0.3 else '❌ Low'}")
|
||||
|
||||
# Step 4: Risk metrics for spread returns
|
||||
print("\n4. Spread risk metrics:")
|
||||
spread_returns = [(spread[i] - spread[i-1]) for i in range(1, len(spread))]
|
||||
mean, std, skew, kurt, sharpe = optimizr.return_statistics_py(spread_returns)
|
||||
print(f" Mean: {mean:.4f}, Volatility: {std:.4f}")
|
||||
print(f" Sharpe: {sharpe:.3f}")
|
||||
|
||||
# Final recommendation
|
||||
print("\n" + "="*70)
|
||||
print("📊 Trading Recommendation:")
|
||||
if mean_reverting and avg_hl < 20 and avg_corr > 0.5:
|
||||
print("✅ STRONG PAIR: Good candidate for pairs trading")
|
||||
print(f" - Fast mean reversion ({avg_hl:.1f} periods)")
|
||||
print(f" - Stable correlation ({avg_corr:.2f})")
|
||||
print(f" - Predictable behavior (H={avg_hurst:.2f})")
|
||||
elif mean_reverting and avg_corr > 0.3:
|
||||
print("⚠️ MODERATE PAIR: Consider with caution")
|
||||
print(f" - Mean reversion detected")
|
||||
print(f" - Moderate correlation ({avg_corr:.2f})")
|
||||
else:
|
||||
print("❌ WEAK PAIR: Not recommended for pairs trading")
|
||||
print(f" - Low correlation or trending behavior")
|
||||
print("="*70)
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
print("=" * 70)
|
||||
print("OptimizR Time-Series Integration Helpers")
|
||||
print("=" * 70)
|
||||
|
||||
# Run all examples
|
||||
example_prepare_for_hmm()
|
||||
example_rolling_hurst()
|
||||
example_rolling_half_life()
|
||||
example_return_statistics()
|
||||
example_lagged_features()
|
||||
example_rolling_correlation()
|
||||
|
||||
# Integrated workflow
|
||||
example_integrated_workflow()
|
||||
|
||||
print("\n" + "=" * 70)
|
||||
print("✅ All examples completed successfully!")
|
||||
print("=" * 70)
|
||||
@@ -23,6 +23,13 @@ from optimizr.core import (
|
||||
compute_risk_metrics_py,
|
||||
estimate_half_life_py,
|
||||
bootstrap_returns_py,
|
||||
# Time-series utilities
|
||||
prepare_for_hmm_py,
|
||||
rolling_hurst_exponent_py,
|
||||
rolling_half_life_py,
|
||||
return_statistics_py,
|
||||
create_lagged_features_py,
|
||||
rolling_correlation_py,
|
||||
)
|
||||
|
||||
# Try to import maths_toolkit from Rust backend
|
||||
@@ -47,5 +54,12 @@ __all__ = [
|
||||
"compute_risk_metrics_py",
|
||||
"estimate_half_life_py",
|
||||
"bootstrap_returns_py",
|
||||
# Time-series utilities
|
||||
"prepare_for_hmm_py",
|
||||
"rolling_hurst_exponent_py",
|
||||
"rolling_half_life_py",
|
||||
"return_statistics_py",
|
||||
"create_lagged_features_py",
|
||||
"rolling_correlation_py",
|
||||
"maths_toolkit",
|
||||
]
|
||||
|
||||
@@ -21,6 +21,13 @@ try:
|
||||
compute_risk_metrics_py,
|
||||
estimate_half_life_py,
|
||||
bootstrap_returns_py,
|
||||
# Time-series utilities
|
||||
prepare_for_hmm_py,
|
||||
rolling_hurst_exponent_py,
|
||||
rolling_half_life_py,
|
||||
return_statistics_py,
|
||||
create_lagged_features_py,
|
||||
rolling_correlation_py,
|
||||
)
|
||||
RUST_AVAILABLE = True
|
||||
except ImportError:
|
||||
|
||||
@@ -33,6 +33,7 @@ use pyo3::types::PyModule;
|
||||
pub mod core;
|
||||
pub mod functional;
|
||||
pub mod maths_toolkit; // Mathematical utilities
|
||||
pub mod timeseries_utils; // Time-series integration helpers
|
||||
|
||||
// Modular structure (trait-based, generic)
|
||||
pub mod de;
|
||||
@@ -98,5 +99,8 @@ fn _core(_py: Python, m: &Bound<'_, PyModule>) -> PyResult<()> {
|
||||
m.add_function(wrap_pyfunction!(risk_metrics::estimate_half_life_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(risk_metrics::bootstrap_returns_py, m)?)?;
|
||||
|
||||
// Time-series utility functions
|
||||
timeseries_utils::python_bindings::register_python_functions(m)?;
|
||||
|
||||
Ok(())
|
||||
}
|
||||
|
||||
@@ -0,0 +1,433 @@
|
||||
//! Time-Series Utilities for Finance and Trading
|
||||
//!
|
||||
//! Helper functions for common workflows combining time-series preprocessing
|
||||
//! with optimization and statistical inference. Designed to work seamlessly
|
||||
//! with Polaroid time-series operations and OptimizR algorithms.
|
||||
//!
|
||||
//! # Use Cases
|
||||
//!
|
||||
//! - Regime detection with HMM
|
||||
//! - Rolling risk metrics (Hurst exponent, half-life)
|
||||
//! - Trading strategy parameter optimization
|
||||
//! - Feature engineering for financial data
|
||||
//!
|
||||
//! # Examples
|
||||
//!
|
||||
//! ```rust
|
||||
//! use optimizr::timeseries_utils::*;
|
||||
//!
|
||||
//! // Prepare features for HMM regime detection
|
||||
//! let prices = vec![100.0, 101.0, 99.5, 102.0];
|
||||
//! let features = prepare_for_hmm(&prices, &[1, 5]);
|
||||
//!
|
||||
//! // Rolling Hurst exponent
|
||||
//! let returns = vec![0.01, -0.02, 0.015, 0.005];
|
||||
//! let rolling_hurst = rolling_hurst_exponent(&returns, 20);
|
||||
//! ```
|
||||
|
||||
use crate::risk_metrics::{hurst_exponent, estimate_half_life};
|
||||
use ndarray::Array1;
|
||||
|
||||
#[cfg(feature = "python-bindings")]
|
||||
pub mod python_bindings;
|
||||
|
||||
/// Prepare time-series price data for HMM regime detection
|
||||
///
|
||||
/// Creates features from price series including:
|
||||
/// - Returns (percent change)
|
||||
/// - Lagged returns
|
||||
/// - Log returns
|
||||
/// - Volatility proxy (absolute returns)
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `prices` - Raw price series (e.g., stock prices, crypto prices)
|
||||
/// * `lag_periods` - Lags to compute for returns (e.g., [1, 5, 20] for daily, weekly, monthly)
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Matrix where each row is a feature vector suitable for HMM training.
|
||||
/// First row will have NaN values due to lagging.
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::prepare_for_hmm;
|
||||
///
|
||||
/// let prices = vec![100.0, 101.0, 99.5, 102.0, 103.5];
|
||||
/// let features = prepare_for_hmm(&prices, &[1, 2]);
|
||||
///
|
||||
/// // Features: [return_t, return_t-1, return_t-2, log_return, abs_return]
|
||||
/// assert_eq!(features[0].len(), 5);
|
||||
/// ```
|
||||
pub fn prepare_for_hmm(prices: &[f64], lag_periods: &[usize]) -> Vec<Vec<f64>> {
|
||||
let n = prices.len();
|
||||
if n < 2 {
|
||||
return vec![];
|
||||
}
|
||||
|
||||
// Calculate returns
|
||||
let mut returns = Vec::with_capacity(n - 1);
|
||||
let mut log_returns = Vec::with_capacity(n - 1);
|
||||
let mut abs_returns = Vec::with_capacity(n - 1);
|
||||
|
||||
for i in 1..n {
|
||||
let ret = (prices[i] - prices[i - 1]) / prices[i - 1];
|
||||
returns.push(ret);
|
||||
log_returns.push((prices[i] / prices[i - 1]).ln());
|
||||
abs_returns.push(ret.abs());
|
||||
}
|
||||
|
||||
// Create feature matrix
|
||||
let max_lag = *lag_periods.iter().max().unwrap_or(&0);
|
||||
let mut features = Vec::new();
|
||||
|
||||
for t in max_lag..returns.len() {
|
||||
let mut feature_vec = Vec::new();
|
||||
|
||||
// Current return
|
||||
feature_vec.push(returns[t]);
|
||||
|
||||
// Lagged returns
|
||||
for &lag in lag_periods {
|
||||
if t >= lag {
|
||||
feature_vec.push(returns[t - lag]);
|
||||
} else {
|
||||
feature_vec.push(f64::NAN);
|
||||
}
|
||||
}
|
||||
|
||||
// Log return and volatility proxy
|
||||
feature_vec.push(log_returns[t]);
|
||||
feature_vec.push(abs_returns[t]);
|
||||
|
||||
features.push(feature_vec);
|
||||
}
|
||||
|
||||
features
|
||||
}
|
||||
|
||||
/// Compute Hurst exponent in rolling windows
|
||||
///
|
||||
/// The Hurst exponent H indicates:
|
||||
/// - H < 0.5: Mean-reverting series
|
||||
/// - H = 0.5: Random walk (no memory)
|
||||
/// - H > 0.5: Trending series (momentum)
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `returns` - Return series (percent changes or log returns)
|
||||
/// * `window_size` - Rolling window size (e.g., 252 for 1 year of daily data)
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Vector of Hurst exponents, one per window. Length = returns.len() - window_size + 1
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::rolling_hurst_exponent;
|
||||
///
|
||||
/// let returns = vec![0.01; 300]; // Synthetic data
|
||||
/// let rolling_h = rolling_hurst_exponent(&returns, 252);
|
||||
///
|
||||
/// assert_eq!(rolling_h.len(), 300 - 252 + 1);
|
||||
/// ```
|
||||
pub fn rolling_hurst_exponent(returns: &[f64], window_size: usize) -> Vec<f64> {
|
||||
let n = returns.len();
|
||||
if n < window_size {
|
||||
return vec![];
|
||||
}
|
||||
|
||||
let mut rolling_h = Vec::with_capacity(n - window_size + 1);
|
||||
|
||||
for i in 0..=(n - window_size) {
|
||||
let window = &returns[i..i + window_size];
|
||||
let window_arr = Array1::from_vec(window.to_vec());
|
||||
// Use standard window sizes for R/S analysis
|
||||
let window_sizes = vec![8, 16, 32].into_iter().filter(|&w| w <= window_size / 4).collect::<Vec<_>>();
|
||||
let h = if window_sizes.is_empty() {
|
||||
0.5 // Default to random walk if window too small
|
||||
} else {
|
||||
hurst_exponent(&window_arr, &window_sizes)
|
||||
.map(|result| result.hurst_exponent)
|
||||
.unwrap_or(0.5)
|
||||
};
|
||||
rolling_h.push(h);
|
||||
}
|
||||
|
||||
rolling_h
|
||||
}
|
||||
|
||||
/// Compute half-life of mean reversion in rolling windows
|
||||
///
|
||||
/// Half-life is the expected time for a price series to revert halfway
|
||||
/// to its mean. Useful for pairs trading and mean-reversion strategies.
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `prices` - Price series (not returns)
|
||||
/// * `window_size` - Rolling window size
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Vector of half-life estimates in same time units as data
|
||||
/// (e.g., days if daily prices)
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::rolling_half_life;
|
||||
///
|
||||
/// let prices = vec![100.0; 300]; // Synthetic data
|
||||
/// let rolling_hl = rolling_half_life(&prices, 100);
|
||||
///
|
||||
/// assert_eq!(rolling_hl.len(), 300 - 100 + 1);
|
||||
/// ```
|
||||
pub fn rolling_half_life(prices: &[f64], window_size: usize) -> Vec<f64> {
|
||||
let n = prices.len();
|
||||
if n < window_size {
|
||||
return vec![];
|
||||
}
|
||||
|
||||
let mut rolling_hl = Vec::with_capacity(n - window_size + 1);
|
||||
|
||||
for i in 0..=(n - window_size) {
|
||||
let window = &prices[i..i + window_size];
|
||||
let window_arr = Array1::from_vec(window.to_vec());
|
||||
let hl = estimate_half_life(&window_arr).unwrap_or(f64::INFINITY);
|
||||
rolling_hl.push(hl);
|
||||
}
|
||||
|
||||
rolling_hl
|
||||
}
|
||||
|
||||
/// Calculate basic time-series statistics
|
||||
///
|
||||
/// Returns summary statistics useful for risk analysis and diagnostics.
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `returns` - Return series
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Tuple of (mean, std_dev, skewness, kurtosis, sharpe_ratio)
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::return_statistics;
|
||||
///
|
||||
/// let returns = vec![0.01, -0.02, 0.015, 0.005, -0.01];
|
||||
/// let (mean, std, skew, kurt, sharpe) = return_statistics(&returns);
|
||||
/// ```
|
||||
pub fn return_statistics(returns: &[f64]) -> (f64, f64, f64, f64, f64) {
|
||||
let n = returns.len() as f64;
|
||||
if returns.is_empty() {
|
||||
return (0.0, 0.0, 0.0, 0.0, 0.0);
|
||||
}
|
||||
|
||||
// Mean
|
||||
let mean = returns.iter().sum::<f64>() / n;
|
||||
|
||||
// Standard deviation
|
||||
let variance = returns.iter().map(|r| (r - mean).powi(2)).sum::<f64>() / n;
|
||||
let std_dev = variance.sqrt();
|
||||
|
||||
// Skewness
|
||||
let skewness = if std_dev > 1e-10 {
|
||||
returns.iter().map(|r| ((r - mean) / std_dev).powi(3)).sum::<f64>() / n
|
||||
} else {
|
||||
0.0
|
||||
};
|
||||
|
||||
// Excess kurtosis
|
||||
let kurtosis = if std_dev > 1e-10 {
|
||||
returns.iter().map(|r| ((r - mean) / std_dev).powi(4)).sum::<f64>() / n - 3.0
|
||||
} else {
|
||||
0.0
|
||||
};
|
||||
|
||||
// Sharpe ratio (assuming 252 trading days, risk-free rate = 0)
|
||||
let sharpe_ratio = if std_dev > 1e-10 {
|
||||
(mean * 252.0_f64.sqrt()) / std_dev
|
||||
} else {
|
||||
0.0
|
||||
};
|
||||
|
||||
(mean, std_dev, skewness, kurtosis, sharpe_ratio)
|
||||
}
|
||||
|
||||
/// Create lagged features for machine learning
|
||||
///
|
||||
/// Useful for creating supervised learning datasets from time series.
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `series` - Input time series (prices, returns, etc.)
|
||||
/// * `lags` - Vector of lag periods (e.g., [1, 2, 5, 10])
|
||||
/// * `include_original` - Whether to include t=0 (current value)
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Matrix where each row is [t, t-lag1, t-lag2, ...] if include_original=true
|
||||
/// or [t-lag1, t-lag2, ...] if false
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::create_lagged_features;
|
||||
///
|
||||
/// let prices = vec![100.0, 101.0, 99.5, 102.0, 103.5, 104.0];
|
||||
/// let features = create_lagged_features(&prices, &[1, 2], true);
|
||||
///
|
||||
/// // Each row: [price_t, price_t-1, price_t-2]
|
||||
/// ```
|
||||
pub fn create_lagged_features(series: &[f64], lags: &[usize], include_original: bool) -> Vec<Vec<f64>> {
|
||||
let n = series.len();
|
||||
let max_lag = *lags.iter().max().unwrap_or(&0);
|
||||
|
||||
if n <= max_lag {
|
||||
return vec![];
|
||||
}
|
||||
|
||||
let mut features = Vec::new();
|
||||
|
||||
for t in max_lag..n {
|
||||
let mut feature_vec = Vec::new();
|
||||
|
||||
if include_original {
|
||||
feature_vec.push(series[t]);
|
||||
}
|
||||
|
||||
for &lag in lags {
|
||||
feature_vec.push(series[t - lag]);
|
||||
}
|
||||
|
||||
features.push(feature_vec);
|
||||
}
|
||||
|
||||
features
|
||||
}
|
||||
|
||||
/// Compute rolling correlation between two series
|
||||
///
|
||||
/// Useful for pairs trading and correlation analysis.
|
||||
///
|
||||
/// # Arguments
|
||||
///
|
||||
/// * `series1` - First time series
|
||||
/// * `series2` - Second time series (must be same length as series1)
|
||||
/// * `window_size` - Rolling window size
|
||||
///
|
||||
/// # Returns
|
||||
///
|
||||
/// Vector of correlation coefficients
|
||||
///
|
||||
/// # Example
|
||||
///
|
||||
/// ```
|
||||
/// use optimizr::timeseries_utils::rolling_correlation;
|
||||
///
|
||||
/// let spy = vec![100.0, 101.0, 99.5, 102.0, 103.5];
|
||||
/// let qqq = vec![200.0, 202.0, 199.0, 204.0, 207.0];
|
||||
/// let corr = rolling_correlation(&spy, &qqq, 3);
|
||||
/// ```
|
||||
pub fn rolling_correlation(series1: &[f64], series2: &[f64], window_size: usize) -> Vec<f64> {
|
||||
let n = series1.len();
|
||||
if n != series2.len() || n < window_size {
|
||||
return vec![];
|
||||
}
|
||||
|
||||
let mut rolling_corr = Vec::with_capacity(n - window_size + 1);
|
||||
|
||||
for i in 0..=(n - window_size) {
|
||||
let x = &series1[i..i + window_size];
|
||||
let y = &series2[i..i + window_size];
|
||||
|
||||
// Compute correlation
|
||||
let mean_x = x.iter().sum::<f64>() / window_size as f64;
|
||||
let mean_y = y.iter().sum::<f64>() / window_size as f64;
|
||||
|
||||
let mut cov = 0.0;
|
||||
let mut var_x = 0.0;
|
||||
let mut var_y = 0.0;
|
||||
|
||||
for j in 0..window_size {
|
||||
let dx = x[j] - mean_x;
|
||||
let dy = y[j] - mean_y;
|
||||
cov += dx * dy;
|
||||
var_x += dx * dx;
|
||||
var_y += dy * dy;
|
||||
}
|
||||
|
||||
let corr = if var_x > 1e-10 && var_y > 1e-10 {
|
||||
cov / (var_x.sqrt() * var_y.sqrt())
|
||||
} else {
|
||||
0.0
|
||||
};
|
||||
|
||||
rolling_corr.push(corr);
|
||||
}
|
||||
|
||||
rolling_corr
|
||||
}
|
||||
|
||||
#[cfg(test)]
|
||||
mod tests {
|
||||
use super::*;
|
||||
|
||||
#[test]
|
||||
fn test_prepare_for_hmm() {
|
||||
let prices = vec![100.0, 101.0, 99.5, 102.0, 103.5, 104.0];
|
||||
let features = prepare_for_hmm(&prices, &[1, 2]);
|
||||
|
||||
assert!(!features.is_empty());
|
||||
// Each feature vector: [return, lag1, lag2, log_return, abs_return]
|
||||
assert_eq!(features[0].len(), 5);
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn test_rolling_hurst_exponent() {
|
||||
let returns = vec![0.01; 50];
|
||||
let rolling_h = rolling_hurst_exponent(&returns, 20);
|
||||
|
||||
assert_eq!(rolling_h.len(), 50 - 20 + 1);
|
||||
// Constant series should have H close to 0.5
|
||||
for h in rolling_h {
|
||||
assert!((h - 0.5).abs() < 0.3); // Allow some numerical variation
|
||||
}
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn test_return_statistics() {
|
||||
let returns = vec![0.01, -0.02, 0.015, 0.005, -0.01];
|
||||
let (mean, std, skew, kurt, sharpe) = return_statistics(&returns);
|
||||
|
||||
assert!((mean).abs() < 0.1); // Small mean
|
||||
assert!(std > 0.0); // Non-zero volatility
|
||||
// Skewness and kurtosis can vary widely
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn test_create_lagged_features() {
|
||||
let series = vec![1.0, 2.0, 3.0, 4.0, 5.0];
|
||||
let features = create_lagged_features(&series, &[1, 2], true);
|
||||
|
||||
assert_eq!(features.len(), 3); // 5 - 2 = 3 valid rows
|
||||
assert_eq!(features[0], vec![3.0, 2.0, 1.0]); // t=2: [val_t, val_t-1, val_t-2]
|
||||
}
|
||||
|
||||
#[test]
|
||||
fn test_rolling_correlation() {
|
||||
let x = vec![1.0, 2.0, 3.0, 4.0, 5.0];
|
||||
let y = vec![2.0, 4.0, 6.0, 8.0, 10.0]; // Perfect correlation
|
||||
let corr = rolling_correlation(&x, &y, 3);
|
||||
|
||||
assert_eq!(corr.len(), 3);
|
||||
for c in corr {
|
||||
assert!((c - 1.0).abs() < 1e-6); // Perfect positive correlation
|
||||
}
|
||||
}
|
||||
}
|
||||
@@ -0,0 +1,62 @@
|
||||
//! Python bindings for time-series utilities
|
||||
|
||||
use pyo3::prelude::*;
|
||||
|
||||
use super::{
|
||||
create_lagged_features, prepare_for_hmm, return_statistics, rolling_correlation,
|
||||
rolling_half_life, rolling_hurst_exponent,
|
||||
};
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (prices, lag_periods))]
|
||||
pub fn prepare_for_hmm_py(prices: Vec<f64>, lag_periods: Vec<usize>) -> Vec<Vec<f64>> {
|
||||
prepare_for_hmm(&prices, &lag_periods)
|
||||
}
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (returns, window_size))]
|
||||
pub fn rolling_hurst_exponent_py(returns: Vec<f64>, window_size: usize) -> Vec<f64> {
|
||||
rolling_hurst_exponent(&returns, window_size)
|
||||
}
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (prices, window_size))]
|
||||
pub fn rolling_half_life_py(prices: Vec<f64>, window_size: usize) -> Vec<f64> {
|
||||
rolling_half_life(&prices, window_size)
|
||||
}
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (returns,))]
|
||||
pub fn return_statistics_py(returns: Vec<f64>) -> (f64, f64, f64, f64, f64) {
|
||||
return_statistics(&returns)
|
||||
}
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (series, lags, include_original=true))]
|
||||
pub fn create_lagged_features_py(
|
||||
series: Vec<f64>,
|
||||
lags: Vec<usize>,
|
||||
include_original: bool,
|
||||
) -> Vec<Vec<f64>> {
|
||||
create_lagged_features(&series, &lags, include_original)
|
||||
}
|
||||
|
||||
#[pyfunction]
|
||||
#[pyo3(signature = (series1, series2, window_size))]
|
||||
pub fn rolling_correlation_py(
|
||||
series1: Vec<f64>,
|
||||
series2: Vec<f64>,
|
||||
window_size: usize,
|
||||
) -> Vec<f64> {
|
||||
rolling_correlation(&series1, &series2, window_size)
|
||||
}
|
||||
|
||||
pub fn register_python_functions(m: &Bound<'_, pyo3::types::PyModule>) -> PyResult<()> {
|
||||
m.add_function(wrap_pyfunction!(prepare_for_hmm_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(rolling_hurst_exponent_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(rolling_half_life_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(return_statistics_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(create_lagged_features_py, m)?)?;
|
||||
m.add_function(wrap_pyfunction!(rolling_correlation_py, m)?)?;
|
||||
Ok(())
|
||||
}
|
||||
Reference in New Issue
Block a user