Skip to main content

quantwave_core/indicators/
correlation_cycle.rs

1use crate::indicators::metadata::{IndicatorMetadata, ParamDef};
2use crate::traits::Next;
3use crate::utils::RingBuffer as VecDeque;
4use std::f64::consts::PI;
5
6/// Correlation Cycle Indicator
7///
8/// Based on John Ehlers' "Correlation as a Cycle Indicator".
9/// Computes the correlation of price against a fixed-period Cosine and negative Sine wave
10/// to derive a phase angle. This angle is used to identify market cycles and trends.
11#[derive(Debug, Clone)]
12pub struct CorrelationCycle {
13    period: usize,
14    price_window: VecDeque<f64>,
15    cosine_wave: Vec<f64>,
16    sine_wave: Vec<f64>,
17    prev_angle: f64,
18    count: usize,
19}
20
21impl CorrelationCycle {
22    pub fn new(period: usize) -> Self {
23        let mut cosine_wave = Vec::with_capacity(period);
24        let mut sine_wave = Vec::with_capacity(period);
25        for n in 0..period {
26            let angle = 2.0 * PI * n as f64 / period as f64;
27            cosine_wave.push(angle.cos());
28            sine_wave.push(-angle.sin());
29        }
30
31        Self {
32            period,
33            price_window: VecDeque::with_capacity(period),
34            cosine_wave,
35            sine_wave,
36            prev_angle: 0.0,
37            count: 0,
38        }
39    }
40
41    fn pearson_correlation(n: usize, x: &VecDeque<f64>, y: &[f64]) -> f64 {
42        let mut sx = 0.0;
43        let mut sy = 0.0;
44        let mut sxx = 0.0;
45        let mut syy = 0.0;
46        let mut sxy = 0.0;
47
48        for i in 0..n {
49            let xi = x[i];
50            let yi = y[i];
51            sx += xi;
52            sy += yi;
53            sxx += xi * xi;
54            syy += yi * yi;
55            sxy += xi * yi;
56        }
57
58        let nf = n as f64;
59        let num = nf * sxy - sx * sy;
60        let den = ((nf * sxx - sx * sx) * (nf * syy - sy * sy)).sqrt();
61
62        if den > 0.0 { num / den } else { 0.0 }
63    }
64}
65
66impl Default for CorrelationCycle {
67    fn default() -> Self {
68        Self::new(20)
69    }
70}
71
72impl Next<f64> for CorrelationCycle {
73    type Output = (f64, f64, f64); // (Real, Imag, Angle)
74
75    fn next(&mut self, input: f64) -> Self::Output {
76        self.count += 1;
77        self.price_window.push_front(input);
78        if self.price_window.len() > self.period {
79            self.price_window.pop_back();
80        }
81
82        if self.price_window.len() < self.period {
83            return (0.0, 0.0, 0.0);
84        }
85
86        let real = Self::pearson_correlation(self.period, &self.price_window, &self.cosine_wave);
87        let imag = Self::pearson_correlation(self.period, &self.price_window, &self.sine_wave);
88
89        // Angle in degrees: 90 + atan(Real / Imag)
90        // Ehlers resolve ambiguity: if Imag > 0 then Angle = Angle - 180
91        let mut angle = if imag != 0.0 {
92            (real / imag).atan().to_degrees() + 90.0
93        } else {
94            90.0
95        };
96
97        if imag > 0.0 {
98            angle -= 180.0;
99        }
100
101        // Do not allow rate change of angle to go negative
102        // If Angle[1] - Angle < 270 and Angle < Angle[1] Then Angle = Angle[1];
103        if self.count > self.period + 1
104            && self.prev_angle - angle < 270.0
105            && angle < self.prev_angle
106        {
107            angle = self.prev_angle;
108        }
109
110        self.prev_angle = angle;
111
112        (real, imag, angle)
113    }
114}
115
116pub const CORRELATION_CYCLE_METADATA: IndicatorMetadata = IndicatorMetadata {
117    name: "CorrelationCycle",
118    description: "Determines cycle phase angle by correlating price with orthogonal sinusoids.",
119    usage: "Use to measure the dominant cycle period via autocorrelation in an amplitude-independent way. Prefer over DFT methods when price amplitude varies significantly across the measurement window.",
120    keywords: &["cycle", "dominant-cycle", "ehlers", "dsp", "spectral"],
121    ehlers_summary: "Ehlers introduces Correlation Cycle measurement in Cycle Analytics for Traders (2013) as an improvement on DFT. By normalizing autocorrelation coefficients to unity variance, the resulting periodogram is independent of price amplitude variations, producing more consistent cycle period estimates.",
122    params: &[ParamDef {
123        name: "period",
124        default: "20",
125        description: "Correlation wavelength",
126    }],
127    formula_source: "https://github.com/lavs9/quantwave/blob/main/references/Ehlers%20Papers/CORRELATION%20AS%20A%20CYCLE%20INDICATOR.pdf",
128    formula_latex: r#"
129\[
130R = \text{Corr}(Price, \cos(2\pi n/P)), I = \text{Corr}(Price, -\sin(2\pi n/P))
131\]
132\[
133\text{Angle} = 90 + \arctan(R/I) \text{ (with quadrant resolution)}
134\]
135"#,
136    gold_standard_file: "correlation_cycle.json",
137    category: "Ehlers DSP",
138};
139
140#[cfg(test)]
141mod tests {
142    use super::*;
143    use crate::traits::Next;
144    use proptest::prelude::*;
145
146    #[test]
147    fn test_correlation_cycle_basic() {
148        let mut cc = CorrelationCycle::new(20);
149        for i in 0..100 {
150            let (r, im, a) = cc.next((2.0 * PI * i as f64 / 20.0).sin());
151            assert!(!r.is_nan());
152            assert!(!im.is_nan());
153            assert!(!a.is_nan());
154        }
155    }
156
157    proptest! {
158        #[test]
159        fn test_correlation_cycle_parity(
160            inputs in prop::collection::vec(1.0..100.0, 100..150),
161        ) {
162            let period = 20;
163            let mut cc = CorrelationCycle::new(period);
164            let streaming_results: Vec<(f64, f64, f64)> = inputs.iter().map(|&x| cc.next(x)).collect();
165
166            // Batch implementation
167            let mut batch_results = Vec::with_capacity(inputs.len());
168            let mut cos_w = Vec::new();
169            let mut sin_w = Vec::new();
170            for n in 0..period {
171                let ang = 2.0 * PI * n as f64 / period as f64;
172                cos_w.push(ang.cos());
173                sin_w.push(-ang.sin());
174            }
175
176            let mut prev_a = 0.0;
177            for i in 0..inputs.len() {
178                if i < period - 1 {
179                    batch_results.push((0.0, 0.0, 0.0));
180                    continue;
181                }
182
183                let mut sx = 0.0;
184                let mut sy_c = 0.0;
185                let mut sy_s = 0.0;
186                let mut sxx = 0.0;
187                let mut syy_c = 0.0;
188                let mut syy_s = 0.0;
189                let mut sxy_c = 0.0;
190                let mut sxy_s = 0.0;
191
192                for j in 0..period {
193                    let xi = inputs[i - j];
194                    let yc = cos_w[j];
195                    let ys = sin_w[j];
196                    sx += xi;
197                    sy_c += yc;
198                    sy_s += ys;
199                    sxx += xi * xi;
200                    syy_c += yc * yc;
201                    syy_s += ys * ys;
202                    sxy_c += xi * yc;
203                    sxy_s += xi * ys;
204                }
205
206                let nf = period as f64;
207                let den_c = ((nf * sxx - sx * sx) * (nf * syy_c - sy_c * sy_c)).sqrt();
208                let real = if den_c > 0.0 { (nf * sxy_c - sx * sy_c) / den_c } else { 0.0 };
209
210                let den_s = ((nf * sxx - sx * sx) * (nf * syy_s - sy_s * sy_s)).sqrt();
211                let imag = if den_s > 0.0 { (nf * sxy_s - sx * sy_s) / den_s } else { 0.0 };
212
213                let mut angle = if imag != 0.0 {
214                    (real / imag).atan().to_degrees() + 90.0
215                } else {
216                    90.0
217                };
218                if imag > 0.0 { angle -= 180.0; }
219
220                if i > period {
221                    if prev_a - angle < 270.0 && angle < prev_a {
222                        angle = prev_a;
223                    }
224                }
225
226                prev_a = angle;
227                batch_results.push((real, imag, angle));
228            }
229
230            for (s, b) in streaming_results.iter().zip(batch_results.iter()) {
231                approx::assert_relative_eq!(s.0, b.0, epsilon = 1e-10);
232                approx::assert_relative_eq!(s.1, b.1, epsilon = 1e-10);
233                approx::assert_relative_eq!(s.2, b.2, epsilon = 1e-10);
234            }
235        }
236    }
237}