//! EWMA (Exponentially Weighted Moving Average) Calculator
//!
//! Implements adaptive threshold calculation using exponentially weighted moving averages.
//! EWMA provides smooth tracking of trends while being responsive to recent changes.
//!
//! # Formula
//!
//! `EWMA_t` = α * `value_t` + (1 - α) * EWMA_{t-1}
//!
//! where α = 2 / (span + 1)
//!
//! # Usage
//!
//! ```rust
//! use ml::features::ewma::EWMACalculator;
//!
//! let mut calculator = EWMACalculator::new(100);
//!
//! // Update with new values
//! let ewma1 = calculator.update(100.0);
//! let ewma2 = calculator.update(105.0);
//! let ewma3 = calculator.update(102.0);
//!
//! // Get current EWMA
//! if let Some(current) = calculator.current() {
//!     println!("Current EWMA: {}", current);
//! }
//! ```
//!
//! # Span Selection
//!
//! - **Small span (10-20)**: High responsiveness, tracks recent changes closely
//! - **Medium span (50-100)**: Balanced smoothing and responsiveness
//! - **Large span (200+)**: Heavy smoothing, slower to respond to changes
//!
//! Common spans:
//! - **12**: Very responsive for short-term trends
//! - **26**: Medium-term trends (common in MACD)
//! - **50**: Balanced for most use cases
//! - **100**: Longer-term smoothing
//! - **200**: Very smooth, for long-term trends

use serde::{Deserialize, Serialize};

/// EWMA calculator for adaptive thresholds
///
/// Tracks exponentially weighted moving average of a time series.
/// The span parameter controls how much weight is given to recent vs historical values.
#[derive(Debug, Clone, Serialize, Deserialize)]
pub struct EWMACalculator {
    /// Number of periods for EWMA calculation (e.g., 100)
    span: usize,

    /// Smoothing factor: α = 2 / (span + 1)
    alpha: f64,

    /// Current EWMA value
    ewma: Option<f64>,
}

impl EWMACalculator {
    /// Create a new EWMA calculator with specified span
    ///
    /// # Arguments
    ///
    /// * `span` - Number of periods for EWMA calculation (e.g., 100)
    ///
    /// # Formula
    ///
    /// α = 2 / (span + 1)
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::EWMACalculator;
    ///
    /// // Create calculator with 100-period span
    /// let calculator = EWMACalculator::new(100);
    /// ```
    pub fn new(span: usize) -> Self {
        let alpha = 2.0 / (span as f64 + 1.0);
        Self {
            span,
            alpha,
            ewma: None,
        }
    }

    /// Update EWMA with a new value
    ///
    /// # Arguments
    ///
    /// * `value` - New value to incorporate into EWMA
    ///
    /// # Returns
    ///
    /// Updated EWMA value
    ///
    /// # Formula
    ///
    /// First update: EWMA = value (initialization)
    /// Subsequent updates: EWMA = α * value + (1 - α) * `EWMA_prev`
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::EWMACalculator;
    ///
    /// let mut calculator = EWMACalculator::new(100);
    /// let ewma1 = calculator.update(100.0);  // Initialize
    /// let ewma2 = calculator.update(105.0);  // Update
    /// ```
    pub fn update(&mut self, value: f64) -> f64 {
        self.ewma = Some(match self.ewma {
            Some(prev) => self.alpha * value + (1.0 - self.alpha) * prev,
            None => value, // Initialize on first value
        });
        // ewma is guaranteed to be Some after the match above sets it
        self.ewma.unwrap_or(value)
    }

    /// Get current EWMA value
    ///
    /// # Returns
    ///
    /// Current EWMA if initialized, None otherwise
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::EWMACalculator;
    ///
    /// let mut calculator = EWMACalculator::new(100);
    /// assert!(calculator.current().is_none());  // Not initialized
    ///
    /// calculator.update(100.0);
    /// assert!(calculator.current().is_some());  // Initialized
    /// ```
    pub const fn current(&self) -> Option<f64> {
        self.ewma
    }

    /// Get the span parameter
    ///
    /// # Returns
    ///
    /// Number of periods for EWMA calculation
    pub const fn span(&self) -> usize {
        self.span
    }

    /// Get the alpha (smoothing factor)
    ///
    /// # Returns
    ///
    /// Smoothing factor: α = 2 / (span + 1)
    pub const fn alpha(&self) -> f64 {
        self.alpha
    }

    /// Reset EWMA to uninitialized state
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::EWMACalculator;
    ///
    /// let mut calculator = EWMACalculator::new(100);
    /// calculator.update(100.0);
    /// assert!(calculator.current().is_some());
    ///
    /// calculator.reset();
    /// assert!(calculator.current().is_none());
    /// ```
    pub const fn reset(&mut self) {
        self.ewma = None;
    }

    /// Check if EWMA is initialized
    ///
    /// # Returns
    ///
    /// true if at least one value has been added, false otherwise
    pub const fn is_initialized(&self) -> bool {
        self.ewma.is_some()
    }
}

/// Adaptive threshold using EWMA
///
/// Provides dynamic threshold calculation based on EWMA and standard deviation.
/// Useful for detecting anomalies or regime changes in time series data.
#[derive(Debug, Clone)]
pub struct AdaptiveThreshold {
    /// EWMA calculator for mean tracking
    ewma: EWMACalculator,

    /// EWMA calculator for variance tracking
    variance_ewma: EWMACalculator,

    /// Number of standard deviations for threshold
    num_std: f64,
}

impl AdaptiveThreshold {
    /// Create a new adaptive threshold calculator
    ///
    /// # Arguments
    ///
    /// * `span` - Number of periods for EWMA calculation
    /// * `num_std` - Number of standard deviations for threshold (e.g., 2.0 for 95% confidence)
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::AdaptiveThreshold;
    ///
    /// // Create threshold with 100-period span and 2 standard deviations
    /// let threshold = AdaptiveThreshold::new(100, 2.0);
    /// ```
    pub fn new(span: usize, num_std: f64) -> Self {
        Self {
            ewma: EWMACalculator::new(span),
            variance_ewma: EWMACalculator::new(span),
            num_std,
        }
    }

    /// Update threshold with new value and return (`lower_bound`, `upper_bound`)
    ///
    /// # Arguments
    ///
    /// * `value` - New value to incorporate
    ///
    /// # Returns
    ///
    /// Tuple of (`lower_bound`, `upper_bound`) for adaptive threshold.
    /// Returns (value, value) if not initialized.
    ///
    /// # Examples
    ///
    /// ```rust
    /// use ml::features::ewma::AdaptiveThreshold;
    ///
    /// let mut threshold = AdaptiveThreshold::new(100, 2.0);
    /// let (lower, upper) = threshold.update(100.0);
    /// ```
    pub fn update(&mut self, value: f64) -> (f64, f64) {
        let mean = self.ewma.update(value);

        // Update variance EWMA with squared deviation
        let deviation = value - mean;
        let variance = self.variance_ewma.update(deviation * deviation);
        let std_dev = variance.sqrt();

        // Calculate bounds
        let lower_bound = mean - self.num_std * std_dev;
        let upper_bound = mean + self.num_std * std_dev;

        (lower_bound, upper_bound)
    }

    /// Get current mean (EWMA)
    pub const fn mean(&self) -> Option<f64> {
        self.ewma.current()
    }

    /// Get current standard deviation
    pub fn std_dev(&self) -> Option<f64> {
        self.variance_ewma.current().map(|v| v.sqrt())
    }

    /// Reset threshold to uninitialized state
    pub const fn reset(&mut self) {
        self.ewma.reset();
        self.variance_ewma.reset();
    }
}

#[cfg(test)]
mod tests {
    use super::*;
    use approx::assert_relative_eq;

    #[test]
    fn test_ewma_initialization() {
        let calculator = EWMACalculator::new(100);
        assert_eq!(calculator.span(), 100);
        assert_relative_eq!(calculator.alpha(), 2.0 / 101.0, epsilon = 1e-10);
        assert!(!calculator.is_initialized());
        assert!(calculator.current().is_none());
    }

    #[test]
    fn test_ewma_first_value() {
        let mut calculator = EWMACalculator::new(100);
        let first_value = 100.0;
        let result = calculator.update(first_value);

        assert_relative_eq!(result, first_value, epsilon = 1e-10);
        assert!(calculator.is_initialized());
        assert_relative_eq!(calculator.current().unwrap(), first_value, epsilon = 1e-10);
    }

    #[test]
    fn test_ewma_formula() {
        let span = 10;
        let alpha = 2.0 / 11.0;
        let mut calculator = EWMACalculator::new(span);

        let v1 = 100.0;
        let ewma1 = calculator.update(v1);
        assert_relative_eq!(ewma1, v1, epsilon = 1e-10);

        let v2 = 110.0;
        let expected_ewma2 = alpha * v2 + (1.0 - alpha) * ewma1;
        let ewma2 = calculator.update(v2);
        assert_relative_eq!(ewma2, expected_ewma2, epsilon = 1e-10);
    }

    #[test]
    fn test_ewma_reset() {
        let mut calculator = EWMACalculator::new(100);
        calculator.update(100.0);
        assert!(calculator.is_initialized());

        calculator.reset();
        assert!(!calculator.is_initialized());
        assert!(calculator.current().is_none());
    }

    #[test]
    fn test_adaptive_threshold_basic() {
        let mut threshold = AdaptiveThreshold::new(100, 2.0);

        // Initialize with first value
        let (lower1, upper1) = threshold.update(100.0);
        assert_eq!(lower1, upper1); // First value has zero variance

        // Add more values
        threshold.update(105.0);
        threshold.update(102.0);
        let (lower, upper) = threshold.update(98.0);

        // Bounds should be different and mean should be between them
        assert!(lower < upper);
        if let Some(mean) = threshold.mean() {
            assert!(lower < mean && mean < upper);
        }
    }

    #[test]
    fn test_adaptive_threshold_volatility() {
        let mut threshold = AdaptiveThreshold::new(50, 2.0);

        // Low volatility period
        for _ in 0..20 {
            threshold.update(100.0);
        }
        let (low_vol_lower, low_vol_upper) = threshold.update(100.0);
        let low_vol_range = low_vol_upper - low_vol_lower;

        // High volatility period
        let high_vol_values = vec![100.0, 120.0, 90.0, 130.0, 80.0];
        for value in high_vol_values {
            threshold.update(value);
        }
        let (high_vol_lower, high_vol_upper) = threshold.update(100.0);
        let high_vol_range = high_vol_upper - high_vol_lower;

        // High volatility should produce wider bounds
        assert!(high_vol_range > low_vol_range);
    }
}