Skip to main content

quantwave_core/indicators/incremental/
statistics_ta.rs

1//! Native O(1) TA-Lib statistics: STDDEV, VAR, LINEARREG family, BETA, CORREL, TSF.
2
3use crate::traits::Next;
4use crate::utils::RingBuffer;
5
6#[inline]
7fn linreg_coeffs(ys: &[f64]) -> (f64, f64) {
8    let n = ys.len() as f64;
9    let mut sum_x = 0.0;
10    let mut sum_x2 = 0.0;
11    let mut sum_y = 0.0;
12    let mut sum_xy = 0.0;
13    for (i, &y) in ys.iter().enumerate() {
14        let x = i as f64;
15        sum_x += x;
16        sum_x2 += x * x;
17        sum_y += y;
18        sum_xy += x * y;
19    }
20    let denom = n.mul_add(sum_x2, -sum_x * sum_x);
21    if denom == 0.0 {
22        return (f64::NAN, f64::NAN);
23    }
24    let slope = n.mul_add(sum_xy, -sum_x * sum_y) / denom;
25    let intercept = (sum_y - slope * sum_x) / n;
26    (slope, intercept)
27}
28
29#[derive(Debug, Clone)]
30struct RollingWindow {
31    period: usize,
32    buf: RingBuffer<f64>,
33    sum: f64,
34    sum_sq: f64,
35}
36
37impl RollingWindow {
38    fn new(period: usize) -> Self {
39        Self {
40            period,
41            buf: RingBuffer::with_capacity(period),
42            sum: 0.0,
43            sum_sq: 0.0,
44        }
45    }
46
47    fn push(&mut self, v: f64) -> bool {
48        if self.buf.len() >= self.period {
49            if let Some(old) = self.buf.pop_front() {
50                self.sum -= old;
51                self.sum_sq -= old * old;
52            }
53        }
54        self.buf.push_back(v);
55        self.sum += v;
56        self.sum_sq += v * v;
57        self.buf.len() >= self.period
58    }
59
60    fn values(&self) -> Vec<f64> {
61        self.buf.iter().cloned().collect()
62    }
63}
64
65#[derive(Debug, Clone)]
66#[allow(non_camel_case_types)]
67pub struct TaSTDDEV {
68    pub timeperiod: usize,
69    pub nbdev: f64,
70    inner: RollingWindow,
71}
72
73impl TaSTDDEV {
74    pub fn new(timeperiod: usize, nbdev: f64) -> Self {
75        Self {
76            timeperiod,
77            nbdev,
78            inner: RollingWindow::new(timeperiod),
79        }
80    }
81}
82
83impl Next<f64> for TaSTDDEV {
84    type Output = f64;
85
86    fn next(&mut self, input: f64) -> Self::Output {
87        if !self.inner.push(input) {
88            return f64::NAN;
89        }
90        let n = self.timeperiod as f64;
91        let mean = self.inner.sum / n;
92        let var = (self.inner.sum_sq / n - mean * mean).max(0.0);
93        var.sqrt() * self.nbdev
94    }
95}
96
97#[derive(Debug, Clone)]
98#[allow(non_camel_case_types)]
99pub struct TaVAR {
100    pub timeperiod: usize,
101    pub nbdev: f64,
102    inner: RollingWindow,
103}
104
105impl TaVAR {
106    pub fn new(timeperiod: usize, nbdev: f64) -> Self {
107        Self {
108            timeperiod,
109            nbdev,
110            inner: RollingWindow::new(timeperiod),
111        }
112    }
113}
114
115impl Next<f64> for TaVAR {
116    type Output = f64;
117
118    fn next(&mut self, input: f64) -> Self::Output {
119        if !self.inner.push(input) {
120            return f64::NAN;
121        }
122        let n = self.timeperiod as f64;
123        let mean = self.inner.sum / n;
124        (self.inner.sum_sq / n - mean * mean).max(0.0)
125    }
126}
127
128#[derive(Debug, Clone)]
129struct LinRegCore {
130    period: usize,
131    inner: RollingWindow,
132}
133
134impl LinRegCore {
135    fn new(period: usize) -> Self {
136        Self {
137            period,
138            inner: RollingWindow::new(period),
139        }
140    }
141
142    fn push(&mut self, v: f64) -> Option<(f64, f64)> {
143        if !self.inner.push(v) {
144            return None;
145        }
146        Some(linreg_coeffs(&self.inner.values()))
147    }
148}
149
150#[derive(Debug, Clone)]
151#[allow(non_camel_case_types)]
152pub struct TaLINEARREG {
153    pub timeperiod: usize,
154    core: LinRegCore,
155}
156
157impl TaLINEARREG {
158    pub fn new(timeperiod: usize) -> Self {
159        Self {
160            timeperiod,
161            core: LinRegCore::new(timeperiod),
162        }
163    }
164}
165
166impl Next<f64> for TaLINEARREG {
167    type Output = f64;
168
169    fn next(&mut self, input: f64) -> Self::Output {
170        let Some((slope, intercept)) = self.core.push(input) else {
171            return f64::NAN;
172        };
173        let x = (self.timeperiod - 1) as f64;
174        intercept + slope * x
175    }
176}
177
178#[derive(Debug, Clone)]
179#[allow(non_camel_case_types)]
180pub struct TaLINEARREG_SLOPE {
181    pub timeperiod: usize,
182    core: LinRegCore,
183}
184
185impl TaLINEARREG_SLOPE {
186    pub fn new(timeperiod: usize) -> Self {
187        Self {
188            timeperiod,
189            core: LinRegCore::new(timeperiod),
190        }
191    }
192}
193
194impl Next<f64> for TaLINEARREG_SLOPE {
195    type Output = f64;
196
197    fn next(&mut self, input: f64) -> Self::Output {
198        match self.core.push(input) {
199            Some((slope, _)) => slope,
200            None => f64::NAN,
201        }
202    }
203}
204
205#[derive(Debug, Clone)]
206#[allow(non_camel_case_types)]
207pub struct TaLINEARREG_INTERCEPT {
208    pub timeperiod: usize,
209    core: LinRegCore,
210}
211
212impl TaLINEARREG_INTERCEPT {
213    pub fn new(timeperiod: usize) -> Self {
214        Self {
215            timeperiod,
216            core: LinRegCore::new(timeperiod),
217        }
218    }
219}
220
221impl Next<f64> for TaLINEARREG_INTERCEPT {
222    type Output = f64;
223
224    fn next(&mut self, input: f64) -> Self::Output {
225        match self.core.push(input) {
226            Some((_, intercept)) => intercept,
227            None => f64::NAN,
228        }
229    }
230}
231
232#[derive(Debug, Clone)]
233#[allow(non_camel_case_types)]
234pub struct TaLINEARREG_ANGLE {
235    pub timeperiod: usize,
236    core: LinRegCore,
237}
238
239impl TaLINEARREG_ANGLE {
240    pub fn new(timeperiod: usize) -> Self {
241        Self {
242            timeperiod,
243            core: LinRegCore::new(timeperiod),
244        }
245    }
246}
247
248impl Next<f64> for TaLINEARREG_ANGLE {
249    type Output = f64;
250
251    fn next(&mut self, input: f64) -> Self::Output {
252        match self.core.push(input) {
253            Some((slope, _)) => slope.atan().to_degrees(),
254            None => f64::NAN,
255        }
256    }
257}
258
259#[derive(Debug, Clone)]
260#[allow(non_camel_case_types)]
261pub struct TaTSF {
262    pub timeperiod: usize,
263    core: LinRegCore,
264}
265
266impl TaTSF {
267    pub fn new(timeperiod: usize) -> Self {
268        Self {
269            timeperiod,
270            core: LinRegCore::new(timeperiod),
271        }
272    }
273}
274
275impl Next<f64> for TaTSF {
276    type Output = f64;
277
278    fn next(&mut self, input: f64) -> Self::Output {
279        let Some((slope, intercept)) = self.core.push(input) else {
280            return f64::NAN;
281        };
282        intercept + slope * self.timeperiod as f64
283    }
284}
285
286#[derive(Debug, Clone)]
287struct DualWindow {
288    period: usize,
289    xs: RingBuffer<f64>,
290    ys: RingBuffer<f64>,
291}
292
293impl DualWindow {
294    fn new(period: usize) -> Self {
295        Self {
296            period,
297            xs: RingBuffer::with_capacity(period),
298            ys: RingBuffer::with_capacity(period),
299        }
300    }
301
302    fn push(&mut self, x: f64, y: f64) -> bool {
303        if self.xs.len() >= self.period {
304            let _ = self.xs.pop_front();
305            let _ = self.ys.pop_front();
306        }
307        self.xs.push_back(x);
308        self.ys.push_back(y);
309        self.xs.len() >= self.period
310    }
311
312    fn xs_values(&self) -> Vec<f64> {
313        self.xs.iter().cloned().collect()
314    }
315
316    fn ys_values(&self) -> Vec<f64> {
317        self.ys.iter().cloned().collect()
318    }
319}
320
321#[inline]
322fn correl_beta(xs: &[f64], ys: &[f64]) -> (f64, f64) {
323    let n = xs.len() as f64;
324    let mut sum_x = 0.0;
325    let mut sum_y = 0.0;
326    let mut sum_x2 = 0.0;
327    let mut sum_y2 = 0.0;
328    let mut sum_xy = 0.0;
329    for i in 0..xs.len() {
330        let x = xs[i];
331        let y = ys[i];
332        sum_x += x;
333        sum_y += y;
334        sum_x2 += x * x;
335        sum_y2 += y * y;
336        sum_xy += x * y;
337    }
338    let cov = n.mul_add(sum_xy, -sum_x * sum_y);
339    let var_x = n.mul_add(sum_x2, -sum_x * sum_x);
340    let var_y = n.mul_add(sum_y2, -sum_y * sum_y);
341    if var_x == 0.0 || var_y == 0.0 {
342        return (f64::NAN, f64::NAN);
343    }
344    let correl = cov / (var_x.sqrt() * var_y.sqrt());
345    // Linear regression Y = alpha + beta * X  =>  beta = cov(X,Y) / var(X).
346    let beta = cov / var_x;
347    (correl, beta)
348}
349
350#[derive(Debug, Clone)]
351#[allow(non_camel_case_types)]
352pub struct TaCORREL {
353    pub timeperiod: usize,
354    inner: DualWindow,
355}
356
357impl TaCORREL {
358    pub fn new(timeperiod: usize) -> Self {
359        Self {
360            timeperiod,
361            inner: DualWindow::new(timeperiod),
362        }
363    }
364}
365
366impl Next<(f64, f64)> for TaCORREL {
367    type Output = f64;
368
369    fn next(&mut self, (x, y): (f64, f64)) -> Self::Output {
370        if !self.inner.push(x, y) {
371            return f64::NAN;
372        }
373        correl_beta(&self.inner.xs_values(), &self.inner.ys_values()).0
374    }
375}
376
377/// Beta on percentage returns of two series (TA-Lib `statistic::beta`).
378#[derive(Debug, Clone)]
379#[allow(non_camel_case_types)]
380pub struct TaBETA {
381    pub timeperiod: usize,
382    prev0: Option<f64>,
383    prev1: Option<f64>,
384    rx: RingBuffer<f64>,
385    ry: RingBuffer<f64>,
386    sx: f64,
387    sy: f64,
388    sxx: f64,
389    sxy: f64,
390}
391
392impl TaBETA {
393    pub fn new(timeperiod: usize) -> Self {
394        Self {
395            timeperiod,
396            prev0: None,
397            prev1: None,
398            rx: RingBuffer::with_capacity(timeperiod),
399            ry: RingBuffer::with_capacity(timeperiod),
400            sx: 0.0,
401            sy: 0.0,
402            sxx: 0.0,
403            sxy: 0.0,
404        }
405    }
406
407    fn push_return(&mut self, x: f64, y: f64) -> Option<f64> {
408        let p = self.timeperiod;
409        if self.rx.len() >= p {
410            if let (Some(ox), Some(oy)) = (self.rx.pop_front(), self.ry.pop_front()) {
411                self.sx -= ox;
412                self.sy -= oy;
413                self.sxx -= ox * ox;
414                self.sxy -= ox * oy;
415            }
416        }
417        self.rx.push_back(x);
418        self.ry.push_back(y);
419        self.sx += x;
420        self.sy += y;
421        self.sxx += x * x;
422        self.sxy += x * y;
423        if self.rx.len() < p {
424            return None;
425        }
426        let n = p as f64;
427        let denom = n * self.sxx - self.sx * self.sx;
428        Some(if denom > 0.0 {
429            (n * self.sxy - self.sx * self.sy) / denom
430        } else {
431            0.0
432        })
433    }
434}
435
436impl Next<(f64, f64)> for TaBETA {
437    type Output = f64;
438
439    fn next(&mut self, (v0, v1): (f64, f64)) -> Self::Output {
440        let (Some(p0), Some(p1)) = (self.prev0, self.prev1) else {
441            self.prev0 = Some(v0);
442            self.prev1 = Some(v1);
443            return f64::NAN;
444        };
445        self.prev0 = Some(v0);
446        self.prev1 = Some(v1);
447        if p0 == 0.0 || p1 == 0.0 {
448            return f64::NAN;
449        }
450        let rx = (v0 - p0) / p0;
451        let ry = (v1 - p1) / p1;
452        self.push_return(rx, ry).unwrap_or(f64::NAN)
453    }
454}
455
456#[cfg(test)]
457mod tests {
458    use super::*;
459    use proptest::prelude::*;
460
461    proptest! {
462        #[test]
463        fn test_ta_beta_parity(
464            x in prop::collection::vec(1.0..100.0, 15..100),
465            y in prop::collection::vec(1.0..100.0, 15..100)
466        ) {
467            let len = x.len().min(y.len());
468            let period = 14;
469            let mut ind = TaBETA::new(period);
470            let streaming: Vec<f64> = (0..len).map(|i| ind.next((x[i], y[i]))).collect();
471            let batch = talib_rs::statistic::beta(&x[..len], &y[..len], period)
472                .unwrap_or_else(|_| vec![f64::NAN; len]);
473            for (s, b) in streaming.iter().zip(batch.iter()) {
474                if s.is_nan() { assert!(b.is_nan()); }
475                else if !b.is_nan() { approx::assert_relative_eq!(s, b, epsilon = 1e-5); }
476            }
477        }
478        #[test]
479        fn test_ta_stddev_parity(input in prop::collection::vec(0.1..100.0, 1..100)) {
480            let period = 10;
481            let nbdev = 1.0;
482            let mut s = TaSTDDEV::new(period, nbdev);
483            let streaming: Vec<f64> = input.iter().map(|&x| s.next(x)).collect();
484            let batch = talib_rs::statistic::stddev(&input, period, nbdev)
485                .unwrap_or_else(|_| vec![f64::NAN; input.len()]);
486            for (a, b) in streaming.iter().zip(batch.iter()) {
487                if a.is_nan() { assert!(b.is_nan()); }
488                else { approx::assert_relative_eq!(a, b, epsilon = 1e-6); }
489            }
490        }
491
492        #[test]
493        fn test_ta_linearreg_parity(input in prop::collection::vec(0.1..100.0, 1..100)) {
494            let period = 10;
495            let mut s = TaLINEARREG::new(period);
496            let streaming: Vec<f64> = input.iter().map(|&x| s.next(x)).collect();
497            let batch = talib_rs::statistic::linearreg(&input, period)
498                .unwrap_or_else(|_| vec![f64::NAN; input.len()]);
499            for (a, b) in streaming.iter().zip(batch.iter()) {
500                if a.is_nan() { assert!(b.is_nan()); }
501                else { approx::assert_relative_eq!(a, b, epsilon = 1e-6); }
502            }
503        }
504
505        #[test]
506        fn test_ta_correl_parity(
507            x in prop::collection::vec(0.1..100.0, 1..100),
508            y in prop::collection::vec(0.1..100.0, 1..100)
509        ) {
510            let len = x.len().min(y.len());
511            let period = 10;
512            let mut s = TaCORREL::new(period);
513            let streaming: Vec<f64> = (0..len).map(|i| s.next((x[i], y[i]))).collect();
514            let batch = talib_rs::statistic::correl(&x[..len], &y[..len], period)
515                .unwrap_or_else(|_| vec![f64::NAN; len]);
516            for (a, b) in streaming.iter().zip(batch.iter()) {
517                if a.is_nan() { assert!(b.is_nan()); }
518                else { approx::assert_relative_eq!(a, b, epsilon = 1e-6); }
519            }
520        }
521    }
522}