Skip to main content

quantwave_core/indicators/incremental/
bbands.rs

1//! Native Bollinger Bands — TA-Lib parity, incremental for every `matype`.
2//!
3//! Both paths are streaming and bounded-memory:
4//! * `MaType::Sma` — O(1) amortised per bar via [`RollingVariance`], a
5//!   shifted-data accumulator with a deterministic periodic exact refresh. It
6//!   replaced a `sum` / `sum_sq` pair whose `E[X²] - E[X]²` evaluation lost
7//!   ~1e-6 relative precision on slowly-varying inputs (quantwave-qkft).
8//! * every other `MaType` — the middle band comes from [`MaStream`], and the
9//!   deviation is a two-pass sum of squared deviations over a rolling window of
10//!   `timeperiod` inputs (O(timeperiod) per bar), exactly as the batch
11//!   `talib_rs::overlap::bbands` non-SMA path computes it.
12//!
13//! The non-SMA path used to push every bar onto an ever-growing `history` and
14//! re-run the *batch* `bbands` over the whole thing, i.e. O(n) per bar / O(n²)
15//! per series plus an unbounded memory leak in long-running streams.
16
17use crate::indicators::incremental::ma_stream::MaStream;
18use crate::indicators::incremental::rolling_variance::RollingVariance;
19use crate::indicators::ma_type::MaType;
20use crate::traits::Next;
21
22const NAN_TRIPLE: (f64, f64, f64) = (f64::NAN, f64::NAN, f64::NAN);
23
24/// Bollinger Bands — matches `talib_rs::overlap::bbands` for every `matype`.
25#[derive(Debug, Clone)]
26#[allow(non_camel_case_types)]
27pub struct BBANDS {
28    pub timeperiod: usize,
29    pub nbdevup: f64,
30    pub nbdevdn: f64,
31    pub matype: MaType,
32    /// Rolling window of the last `timeperiod` inputs (both paths). Carries the
33    /// stable mean/variance accumulator used by the SMA fast path; the non-SMA
34    /// path only reads back the window contents.
35    window: RollingVariance,
36    // --- non-SMA path ---
37    ma: Option<MaStream>,
38    bars_seen: usize,
39}
40
41impl BBANDS {
42    pub fn new(timeperiod: usize, nbdevup: f64, nbdevdn: f64, matype: MaType) -> Self {
43        Self {
44            timeperiod,
45            nbdevup,
46            nbdevdn,
47            matype,
48            window: RollingVariance::new(timeperiod),
49            ma: (matype != MaType::Sma).then(|| MaStream::new(timeperiod, matype)),
50            bars_seen: 0,
51        }
52    }
53
54    #[inline]
55    fn bands_from_window(&self) -> (f64, f64, f64) {
56        let ma_val = self.window.mean();
57        let stddev = self.window.stddev();
58        let upper = ma_val + self.nbdevup * stddev;
59        let lower = ma_val - self.nbdevdn * stddev;
60        (upper, ma_val, lower)
61    }
62
63    fn next_sma(&mut self, input: f64) -> (f64, f64, f64) {
64        let tp = self.timeperiod;
65        if tp == 0 {
66            return NAN_TRIPLE;
67        }
68
69        if !self.window.push(input) {
70            return NAN_TRIPLE;
71        }
72
73        self.bands_from_window()
74    }
75
76    fn next_non_sma(&mut self, input: f64) -> (f64, f64, f64) {
77        let tp = self.timeperiod;
78        if tp == 0 {
79            return NAN_TRIPLE;
80        }
81        let i = self.bars_seen;
82        self.bars_seen += 1;
83
84        self.window.push(input);
85        let ma_val = match self.ma {
86            Some(ref mut ma) => ma.next(input),
87            None => return NAN_TRIPLE,
88        };
89
90        if i + 1 < tp {
91            // Before `lookback = timeperiod - 1` the batch fills NaN everywhere.
92            return NAN_TRIPLE;
93        }
94        if ma_val.is_nan() {
95            // Past the BBANDS lookback but still inside the MA's own (longer)
96            // warmup — e.g. DEMA/TEMA/T3/MAMA. The batch leaves the zero-
97            // initialised band values in place here and writes NaN only into
98            // the middle band, so reproduce exactly that rather than inventing
99            // a different warmup convention. See the `continue` in the non-SMA
100            // loop of `talib_rs::overlap::bbands`.
101            return (0.0, f64::NAN, 0.0);
102        }
103
104        // Two-pass deviation about the MA value over the rolling window. Matches
105        // the batch `sum_sq_diff(window, ma_val)` and is the numerically stable
106        // form quantwave-qkft wants everywhere (no E[X²] - E[X]² cancellation).
107        let mut sum_sq_diff = 0.0;
108        for &v in self.window.iter() {
109            let d = v - ma_val;
110            sum_sq_diff += d * d;
111        }
112        let stddev = (sum_sq_diff / tp as f64).max(0.0).sqrt();
113        (
114            ma_val + self.nbdevup * stddev,
115            ma_val,
116            ma_val - self.nbdevdn * stddev,
117        )
118    }
119}
120
121impl Next<f64> for BBANDS {
122    type Output = (f64, f64, f64);
123
124    fn next(&mut self, input: f64) -> Self::Output {
125        if self.matype == MaType::Sma {
126            self.next_sma(input)
127        } else {
128            self.next_non_sma(input)
129        }
130    }
131}
132
133#[cfg(test)]
134mod tests {
135    use super::*;
136    use proptest::prelude::*;
137
138    fn batch_bbands(
139        input: &[f64],
140        period: usize,
141        nbdevup: f64,
142        nbdevdn: f64,
143        matype: MaType,
144    ) -> (Vec<f64>, Vec<f64>, Vec<f64>) {
145        talib_rs::overlap::bbands(input, period, nbdevup, nbdevdn, matype.into()).unwrap_or_else(
146            |_| {
147                (
148                    vec![f64::NAN; input.len()],
149                    vec![f64::NAN; input.len()],
150                    vec![f64::NAN; input.len()],
151                )
152            },
153        )
154    }
155
156    fn assert_parity(streaming: &[(f64, f64, f64)], batch: &(Vec<f64>, Vec<f64>, Vec<f64>)) {
157        let (b_upper, b_middle, b_lower) = batch;
158        for (i, &(s_upper, s_middle, s_lower)) in streaming.iter().enumerate() {
159            for (s, b, name) in [
160                (s_upper, b_upper[i], "upper"),
161                (s_middle, b_middle[i], "middle"),
162                (s_lower, b_lower[i], "lower"),
163            ] {
164                if s.is_nan() {
165                    assert!(b.is_nan(), "bar {i} {name}: streaming NaN, batch {b}");
166                } else {
167                    assert!(!b.is_nan(), "bar {i} {name}: streaming {s}, batch NaN");
168                    approx::assert_relative_eq!(s, b, epsilon = 1e-6);
169                }
170            }
171        }
172    }
173
174    proptest! {
175        #[test]
176        fn test_bbands_parity(input in prop::collection::vec(0.1..100.0, 1..100)) {
177            let period = 10;
178            let matype = MaType::Sma;
179            let mut bbands = BBANDS::new(period, 2.0, 2.0, matype);
180            let streaming: Vec<(f64, f64, f64)> =
181                input.iter().map(|&x| bbands.next(x)).collect();
182            assert_parity(&streaming, &batch_bbands(&input, period, 2.0, 2.0, matype));
183        }
184
185        /// Regression coverage for the non-SMA path: the old suite only ever
186        /// exercised `MaType::Sma`, which is why an O(n)-per-bar batch call hid
187        /// inside `next()` for every other matype. Lengths start past the
188        /// longest MA warmup (MAMA's 32, TEMA's 3*(p-1)) so the batch oracle
189        /// returns real values rather than an `InsufficientData` error.
190        #[test]
191        fn test_bbands_parity_all_matypes(
192            input in prop::collection::vec(0.1..100.0, 120..200),
193            idx in 0usize..9,
194        ) {
195            let matype = MaType::ALL[idx];
196            let period = 10;
197            let mut bbands = BBANDS::new(period, 2.0, 2.5, matype);
198            let streaming: Vec<(f64, f64, f64)> =
199                input.iter().map(|&x| bbands.next(x)).collect();
200            assert_parity(&streaming, &batch_bbands(&input, period, 2.0, 2.5, matype));
201        }
202    }
203
204    // ---- quantwave-qkft ----------------------------------------------------
205
206    /// BBANDS is the headline consumer of the stddev accumulator. On log-price
207    /// the SMA path's band half-width must track a compensated two-pass
208    /// reference to well inside the strict 1e-7 gate.
209    #[test]
210    fn sma_bands_match_two_pass_reference_on_log_price() {
211        use crate::test_utils::{log_price_series, reference_mean, reference_stddev};
212        let period = 20;
213        let nbdev = 2.0;
214        let data = log_price_series(600);
215        let mut bbands = BBANDS::new(period, nbdev, nbdev, MaType::Sma);
216        let mut worst_mid: f64 = 0.0;
217        let mut worst_half: f64 = 0.0;
218        for (i, &x) in data.iter().enumerate() {
219            let (upper, middle, lower) = bbands.next(x);
220            if middle.is_nan() {
221                continue;
222            }
223            let window = &data[i + 1 - period..=i];
224            let ref_mid = reference_mean(window);
225            let ref_half = nbdev * reference_stddev(window);
226            assert!(ref_half > 0.0);
227            worst_mid = worst_mid.max(((middle - ref_mid) / ref_mid).abs());
228            // Compare the half-width, not the band level: `upper` is
229            // dominated by the ~8.0 middle band and would hide the error.
230            worst_half = worst_half.max((((upper - lower) / 2.0 - ref_half) / ref_half).abs());
231        }
232        assert!(
233            worst_half < 1e-7,
234            "BBANDS half-width vs reference on log-price: {worst_half:e} (gate 1e-7)"
235        );
236        // Noise floor here is the *test's* extraction, not the accumulator:
237        // recovering a ~1e-4 half-width by differencing two bands sitting at
238        // level ~8 costs ~eps * 8 / 1e-4 ≈ 1e-11 relative. The accumulator
239        // itself is at ~1e-15 (see rolling_variance's own tests).
240        assert!(worst_half < 1e-10, "expected ~1e-11, got {worst_half:e}");
241        assert!(worst_mid < 1e-12, "middle band drift {worst_mid:e}");
242    }
243
244    /// Bit-identical batch vs chunked streaming, for every `matype`.
245    #[test]
246    fn batch_streaming_bitwise_parity_all_matypes() {
247        use crate::test_utils::log_price_series;
248        let data = log_price_series(400);
249        for &matype in MaType::ALL.iter() {
250            let mut batch = BBANDS::new(10, 2.0, 2.5, matype);
251            let batch_out: Vec<_> = data.iter().map(|&x| batch.next(x)).collect();
252            for chunk in [1usize, 7, 64] {
253                let mut streaming = BBANDS::new(10, 2.0, 2.5, matype);
254                let mut out = Vec::with_capacity(data.len());
255                for part in data.chunks(chunk) {
256                    for &x in part {
257                        out.push(streaming.next(x));
258                    }
259                }
260                for (i, (b, s)) in batch_out.iter().zip(out.iter()).enumerate() {
261                    assert_eq!(
262                        (b.0.to_bits(), b.1.to_bits(), b.2.to_bits()),
263                        (s.0.to_bits(), s.1.to_bits(), s.2.to_bits()),
264                        "{matype:?} bar {i} chunk {chunk}: {b:?} vs {s:?}"
265                    );
266                }
267            }
268        }
269    }
270
271    /// The non-SMA path must not retain more than `timeperiod` samples.
272    #[test]
273    fn non_sma_window_is_bounded() {
274        let mut bbands = BBANDS::new(10, 2.0, 2.0, MaType::Wma);
275        for i in 0..10_000 {
276            bbands.next(50.0 + (i as f64 * 0.01).sin());
277        }
278        assert_eq!(bbands.window.len(), 10);
279    }
280}