gavindiaz/QuanTAlib
1
0
Fork 0
mirror of https://github.com/mihakralj/QuanTAlib.git synced 2026-08-20 19:48:05 +00:00
Code Issues Packages Projects Releases Wiki Activity

test: setup common stability and robustness properties tracking

Browse Source
This commit is contained in:
Miha Kralj
2026-02-27 12:50:05 -08:00
parent 4ab3a7fb53
commit 769a923a24
287 changed files with 1314 additions and 867 deletions
Download Patch File Download Diff File Expand all files Collapse all files
lib/statistics/acf/Acf.Quantower.Tests.cs
+1 -1
View File
@@ -294,4 +294,4 @@ public class AcfIndicatorTests
Assert.True(value >= -1 && value <= 1, $"ACF value at index {i} should be bounded [-1, 1], got {value}");
}
}
}
}
lib/statistics/acf/Acf.Tests.cs
+1 -1
View File
@@ -539,4 +539,4 @@ public class AcfTests
}
#endregion
}
}
lib/statistics/acf/Acf.Validation.Tests.cs
+12 -7
View File
@@ -114,11 +114,13 @@ public class AcfValidationTests
// For white noise, ACF at any lag > 0 should be close to zero
var acf = new Acf(100, 5);
// Generate pseudo-random values with zero mean via GBM log-returns
var random = new GBM(startPrice: 100.0, sigma: 1.0, seed: 42);
for (int i = 0; i < 500; i++)
// Generate incremental bar-to-bar log-returns: log(close_i / close_{i-1})
// These are approximately i.i.d. N(0, σ²·dt) — genuine white noise
var gbm = new GBM(startPrice: 100.0, sigma: 0.2, seed: 42);
var bars = gbm.Fetch(501, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
for (int i = 1; i < 501; i++)
{
double val = Math.Log(random.Next().Close / 100.0); // ~N(0, vol²*dt) centered near 0
double val = Math.Log(bars[i].Close / bars[i - 1].Close); // incremental log-return
acf.Update(new TValue(DateTime.UtcNow.AddSeconds(i), val));
}
@@ -215,12 +217,15 @@ public class AcfValidationTests
double[] ar1Data = new double[n];
ar1Data[0] = 0;
var random = new GBM(startPrice: 100.0, sigma: 1.0, seed: 42);
// Use incremental bar-to-bar log-returns as i.i.d. noise: log(close_i / close_{i-1})
var gbm = new GBM(startPrice: 100.0, sigma: 0.2, seed: 42);
var bars = gbm.Fetch(n, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
// Generate AR(1) process
for (int i = 1; i < n; i++)
{
double epsilon = Math.Log(random.Next().Close / 100.0) * 0.1; // Small noise
double epsilon = Math.Log(bars[i].Close / bars[i - 1].Close) * 0.5; // incremental log-return scaled as noise
ar1Data[i] = phi * ar1Data[i - 1] + epsilon;
}
@@ -280,4 +285,4 @@ public class AcfValidationTests
}
#endregion
}
}
lib/statistics/acf/Acf.cs
+2 -28
View File
@@ -42,12 +42,10 @@ public sealed class Acf : AbstractBase
// Running sums for O(1) updates
private double _sum;
private double _sumSq;
private double _sumLagged;
// Snapshot state for bar correction
private double _p_sum;
private double _p_sumSq;
private double _p_sumLagged;
private int _updateCount;
private const int ResyncInterval = 1000;
@@ -110,7 +108,6 @@ public sealed class Acf : AbstractBase
// Snapshot state for rollback
_p_sum = _sum;
_p_sumSq = _sumSq;
_p_sumLagged = _sumLagged;
_buffer.Snapshot();
}
else
@@ -118,7 +115,6 @@ public sealed class Acf : AbstractBase
// Restore state from snapshot
_sum = _p_sum;
_sumSq = _p_sumSq;
_sumLagged = _p_sumLagged;
_buffer.Restore();
}
@@ -129,12 +125,6 @@ public sealed class Acf : AbstractBase
_sum -= oldVal;
_sumSq = Math.FusedMultiplyAdd(-oldVal, oldVal, _sumSq);
// Remove contribution to lagged sum
if (_buffer.Count > _lag)
{
double oldLaggedVal = _buffer[_lag]; // value that was _lag positions from oldest
_sumLagged -= oldVal * oldLaggedVal;
}
}
// Add new value
@@ -142,14 +132,6 @@ public sealed class Acf : AbstractBase
_sum += value;
_sumSq = Math.FusedMultiplyAdd(value, value, _sumSq);
// Update lagged sum: add product of new value and value at lag positions before
if (_buffer.Count > _lag)
{
int lagIndex = _buffer.Count - 1 - _lag;
double laggedVal = _buffer[lagIndex];
_sumLagged += value * laggedVal;
}
if (isNew)
{
_updateCount++;
@@ -257,18 +239,11 @@ public sealed class Acf : AbstractBase
int n = _buffer.Count;
_sum = 0;
_sumSq = 0;
_sumLagged = 0;
for (int i = 0; i < n; i++)
{
double val = _buffer[i];
_sum += val;
_sumSq += val * val;
if (i >= _lag)
{
_sumLagged += val * _buffer[i - _lag];
}
}
}
@@ -277,19 +252,18 @@ public sealed class Acf : AbstractBase
_buffer.Clear();
_sum = 0;
_sumSq = 0;
_sumLagged = 0;
_p_sum = 0;
_p_sumSq = 0;
_p_sumLagged = 0;
_updateCount = 0;
Last = default;
}
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
Reset();
foreach (double value in source)
{
Update(new TValue(DateTime.UtcNow, value));
Update(new TValue(DateTime.MinValue, value));
}
}
lib/statistics/acf/Acf.md
+1 -1
View File
@@ -15,7 +15,7 @@
- Parameterized by `period`, `lag` (default 1).
- Output range: Varies (see docs).
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against mathematical properties and theoretical AR-process expectations.
> "The past doesn't predict the future, but it whispers patterns to those who listen."
lib/statistics/beta/Beta.Validation.Tests.cs
+1 -1
View File
@@ -130,7 +130,7 @@ public sealed class BetaValidationTests : IDisposable
var retCode = Functions.Beta<double>(
assetPrices.AsSpan(), marketPrices.AsSpan(),
0..^0, taOut, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
(int offset, int length) = outRange.GetOffsetAndLength(taOut.Length);
lib/statistics/cointegration/Cointegration.Quantower.Tests.cs
+1 -1
View File
@@ -264,4 +264,4 @@ public class CointegrationIndicatorTests
Assert.Contains("cointegration", indicator.Description, StringComparison.OrdinalIgnoreCase);
Assert.Contains("ADF", indicator.Description, StringComparison.Ordinal);
}
}
}
lib/statistics/cointegration/Cointegration.Tests.cs
+20 -13
View File
@@ -154,27 +154,34 @@ public class CointegrationTests
[Fact]
public void Update_WithIsNewFalse_DoesNotAdvanceState()
{
var indicator = new Cointegration(5);
var corrected = new Cointegration(5);
var direct = new Cointegration(5);
var gbmA = new GBM(startPrice: 100.0, mu: 0.0, sigma: 0.1, seed: 12345);
var gbmB = new GBM(startPrice: 100.0, mu: 0.0, sigma: 0.1, seed: 54321);
// Build up some state
// Build identical state
for (int i = 0; i < 10; i++)
{
indicator.Update(gbmA.Next().Close, gbmB.Next().Close, isNew: true);
double a = gbmA.Next().Close;
double b = gbmB.Next().Close;
corrected.Update(a, b, isNew: true);
direct.Update(a, b, isNew: true);
}
// Update with same values and isNew=false
indicator.Update(gbmA.Next().Close, gbmB.Next().Close, isNew: false);
var valueAfterFirst = indicator.Last.Value;
const double finalA = 105.0;
const double finalB = 55.0;
// Another correction
indicator.Update(gbmA.Next().Close, gbmB.Next().Close, isNew: false);
var valueAfterSecond = indicator.Last.Value;
// Correction path: add + multiple rewrites + final rewrite to target value
corrected.Update(finalA, finalB, isNew: true);
corrected.Update(finalA + 10.0, finalB + 10.0, isNew: false);
corrected.Update(finalA - 3.0, finalB - 3.0, isNew: false);
corrected.Update(finalA, finalB, isNew: false);
// All corrections replace the same bar, state should be consistent
Assert.True(double.IsFinite(valueAfterFirst) || double.IsNaN(valueAfterFirst));
Assert.True(double.IsFinite(valueAfterSecond) || double.IsNaN(valueAfterSecond));
// Direct path: only final new bar
direct.Update(finalA, finalB, isNew: true);
Assert.Equal(direct.Last.Value, corrected.Last.Value, Tolerance);
}
[Fact]
@@ -648,4 +655,4 @@ public class CointegrationTests
}
#endregion
}
}
lib/statistics/cointegration/Cointegration.Validation.Tests.cs
+24 -17
View File
@@ -18,20 +18,26 @@ public class CointegrationValidationTests
[Fact]
public void Cointegration_PerfectlyCointegrated_ProducesStrongNegativeAdf()
{
// Two series with near-perfect linear relationship should show strong cointegration
// Adding small noise to avoid zero-variance residuals
var indicator = new Cointegration(20);
var random = new GBM(startPrice: 100.0, sigma: 1.0, seed: 42);
// Two series with near-perfect linear relationship should show strong cointegration.
// Use incremental log-returns (i.i.d.) as noise so residuals are stationary.
// Period=30 gives ADF sufficient window; 200 samples ensure stable regression.
var indicator = new Cointegration(30);
var gbm = new GBM(startPrice: 100.0, sigma: 0.2, seed: 42);
var bars = gbm.Fetch(201, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
for (int i = 0; i < 100; i++)
for (int i = 1; i <= 200; i++)
{
double a = 100.0 + i * 0.5 + GbmNoise(random) * 0.1;
double b = 2.0 * a + 10.0 + GbmNoise(random) * 0.1;
// Incremental log-return: truly i.i.d. noise, variance ~(0.2²·dt)
double noise = Math.Log(bars[i].Close / bars[i - 1].Close);
double a = 100.0 + i * 0.5 + noise * 0.1;
double b = 2.0 * a + 10.0 + noise * 0.1;
indicator.Update(a, b);
}
// Near-perfect cointegration should produce strongly negative ADF statistic
Assert.True(indicator.Last.Value < -2.0, $"ADF should be strongly negative for cointegrated series, got {indicator.Last.Value}");
// Near-perfect cointegration should produce ADF below the 5% critical value.
// Engle-Granger critical values (residual-based, no constant): -1.95 at 5%, -2.86 for large N.
// With period=30 and 200 samples of near-linear data the statistic should clear -1.95 comfortably.
Assert.True(indicator.Last.Value < -1.95, $"ADF should be below 5% critical value (-1.95) for cointegrated series, got {indicator.Last.Value}");
}
[Fact]
@@ -68,7 +74,8 @@ public class CointegrationValidationTests
indicator.Update(a, b);
}
Assert.True(indicator.Last.Value < 0, $"ADF should be negative for near-proportional series, got {indicator.Last.Value}");
// Proportional series with small noise should produce ADF well below 0; -1.0 is a conservative bound.
Assert.True(indicator.Last.Value < -1.0, $"ADF should be well negative for near-proportional series, got {indicator.Last.Value}");
}
[Fact]
@@ -86,8 +93,8 @@ public class CointegrationValidationTests
indicator.Update(a, b);
}
// Should still detect cointegration despite small noise
Assert.True(indicator.Last.Value < 0, $"ADF should be negative even with small noise, got {indicator.Last.Value}");
// Linear relationship with small noise should still clear -1.0.
Assert.True(indicator.Last.Value < -1.0, $"ADF should be well negative with small noise, got {indicator.Last.Value}");
}
[Fact]
@@ -226,8 +233,8 @@ public class CointegrationValidationTests
indicator.Update(100.0, 50.0);
}
// Should handle constant series without crashing (result may be NaN due to zero variance)
Assert.True(double.IsNaN(indicator.Last.Value) || double.IsFinite(indicator.Last.Value));
// Constant series → zero variance → ADF denominator is zero → NaN is correct.
Assert.True(double.IsNaN(indicator.Last.Value), $"Expected NaN for constant series, got {indicator.Last.Value}");
}
[Fact]
@@ -240,8 +247,8 @@ public class CointegrationValidationTests
indicator.Update(100.0, 50.0 + i); // A constant, B trending
}
// Should handle mixed constant/trending without crashing
Assert.True(double.IsNaN(indicator.Last.Value) || double.IsFinite(indicator.Last.Value));
// Constant A → zero variance in A → ADF is undefined → NaN.
Assert.True(double.IsNaN(indicator.Last.Value), $"Expected NaN when series A is constant, got {indicator.Last.Value}");
}
[Fact]
@@ -330,4 +337,4 @@ public class CointegrationValidationTests
}
#endregion
}
}
lib/statistics/cointegration/Cointegration.cs
+98 -133
View File
@@ -45,8 +45,7 @@ public sealed class Cointegration : AbstractBase
// ADF regression running sums (period-1 window)
private readonly RingBuffer _deltaResiduals;
private readonly RingBuffer _laggedResiduals;
private double _sumDelta, _sumLagged;
private double _sumDeltaLagged, _sumLagged2;
private double _sumDeltaLagged, _sumLagged2, _sumDelta2;
// Last valid values for NaN handling
private double _lastValidA, _lastValidB;
@@ -113,7 +112,7 @@ public sealed class Cointegration : AbstractBase
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public TValue Update(double seriesA, double seriesB, bool isNew = true)
{
return Update(new TValue(DateTime.UtcNow, seriesA), new TValue(DateTime.UtcNow, seriesB), isNew);
return Update(new TValue(DateTime.MinValue, seriesA), new TValue(DateTime.MinValue, seriesB), isNew);
}
/// <inheritdoc/>
@@ -195,19 +194,17 @@ public sealed class Cointegration : AbstractBase
{
double oldDelta = _deltaResiduals.Oldest;
double oldLagged = _laggedResiduals.Oldest;
_sumDelta -= oldDelta;
_sumLagged -= oldLagged;
_sumDeltaLagged = FusedMultiplyAdd(-oldDelta, oldLagged, _sumDeltaLagged);
_sumLagged2 = FusedMultiplyAdd(-oldLagged, oldLagged, _sumLagged2);
_sumDelta2 = FusedMultiplyAdd(-oldDelta, oldDelta, _sumDelta2);
}
_deltaResiduals.Add(delta);
_laggedResiduals.Add(lagged);
_sumDelta += delta;
_sumLagged += lagged;
_sumDeltaLagged = FusedMultiplyAdd(delta, lagged, _sumDeltaLagged);
_sumLagged2 = FusedMultiplyAdd(lagged, lagged, _sumLagged2);
_sumDelta2 = FusedMultiplyAdd(delta, delta, _sumDelta2);
}
_prevResidual = residual;
@@ -230,31 +227,24 @@ public sealed class Cointegration : AbstractBase
_hasPrevResidual = _p_hasPrevResidual;
// Update newest values in main buffers
if (_bufferA.Count > 0)
if (_bufferA.Count == 0)
{
double oldA = _bufferA.Newest;
double oldB = _bufferB.Newest;
_sumA = FusedMultiplyAdd(1.0, a, FusedMultiplyAdd(-1.0, oldA, _sumA));
_sumB = FusedMultiplyAdd(1.0, b, FusedMultiplyAdd(-1.0, oldB, _sumB));
_sumA2 = FusedMultiplyAdd(a, a, FusedMultiplyAdd(-oldA, oldA, _sumA2));
_sumB2 = FusedMultiplyAdd(b, b, FusedMultiplyAdd(-oldB, oldB, _sumB2));
_sumAB = FusedMultiplyAdd(a, b, FusedMultiplyAdd(-oldA, oldB, _sumAB));
_bufferA.UpdateNewest(a);
_bufferB.UpdateNewest(b);
}
else
{
_bufferA.Add(a);
_bufferB.Add(b);
_sumA = a;
_sumB = b;
_sumA2 = a * a;
_sumB2 = b * b;
_sumAB = a * b;
// Nothing to correct yet; no current bar exists
return;
}
double oldA = _bufferA.Newest;
double oldB = _bufferB.Newest;
_sumA += a - oldA;
_sumB += b - oldB;
_sumA2 = FusedMultiplyAdd(a, a, FusedMultiplyAdd(-oldA, oldA, _sumA2));
_sumB2 = FusedMultiplyAdd(b, b, FusedMultiplyAdd(-oldB, oldB, _sumB2));
_sumAB = FusedMultiplyAdd(a, b, FusedMultiplyAdd(-oldA, oldB, _sumAB));
_bufferA.UpdateNewest(a);
_bufferB.UpdateNewest(b);
// Calculate current residual
double residual = CalculateResidual(a, b);
@@ -264,28 +254,21 @@ public sealed class Cointegration : AbstractBase
double delta = residual - _prevResidual;
double lagged = _prevResidual;
if (_deltaResiduals.Count > 0)
if (_deltaResiduals.Count == 0)
{
double oldDelta = _deltaResiduals.Newest;
double oldLagged = _laggedResiduals.Newest;
_sumDelta = FusedMultiplyAdd(1.0, delta, FusedMultiplyAdd(-1.0, oldDelta, _sumDelta));
_sumLagged = FusedMultiplyAdd(1.0, lagged, FusedMultiplyAdd(-1.0, oldLagged, _sumLagged));
_sumDeltaLagged = FusedMultiplyAdd(delta, lagged, FusedMultiplyAdd(-oldDelta, oldLagged, _sumDeltaLagged));
_sumLagged2 = FusedMultiplyAdd(lagged, lagged, FusedMultiplyAdd(-oldLagged, oldLagged, _sumLagged2));
_deltaResiduals.UpdateNewest(delta);
_laggedResiduals.UpdateNewest(lagged);
}
else
{
_deltaResiduals.Add(delta);
_laggedResiduals.Add(lagged);
_sumDelta = delta;
_sumLagged = lagged;
_sumDeltaLagged = delta * lagged;
_sumLagged2 = lagged * lagged;
// Nothing to correct yet in ADF buffers; no current entry exists
return;
}
double oldDelta = _deltaResiduals.Newest;
double oldLagged = _laggedResiduals.Newest;
_sumDeltaLagged = FusedMultiplyAdd(delta, lagged, FusedMultiplyAdd(-oldDelta, oldLagged, _sumDeltaLagged));
_sumLagged2 = FusedMultiplyAdd(lagged, lagged, FusedMultiplyAdd(-oldLagged, oldLagged, _sumLagged2));
_sumDelta2 = FusedMultiplyAdd(delta, delta, FusedMultiplyAdd(-oldDelta, oldDelta, _sumDelta2));
_deltaResiduals.UpdateNewest(delta);
_laggedResiduals.UpdateNewest(lagged);
}
_prevResidual = residual;
@@ -305,28 +288,15 @@ public sealed class Cointegration : AbstractBase
double meanA = _sumA / n;
double meanB = _sumB / n;
// Calculate variances and covariance
double varA = Max(0.0, (_sumA2 / n) - (meanA * meanA));
// Calculate variance of B and covariance
double varB = Max(0.0, (_sumB2 / n) - (meanB * meanB));
double cov = (_sumAB / n) - (meanA * meanB);
// Calculate standard deviations
double stdA = Sqrt(varA);
double stdB = Sqrt(varB);
// Calculate correlation
double correlation = 0.0;
double denom = stdA * stdB;
if (Abs(denom) > Epsilon)
{
correlation = cov / denom;
}
// Calculate beta and alpha
double beta = 0.0;
if (Abs(stdB) > Epsilon)
if (varB > Epsilon)
{
beta = correlation * (stdA / stdB);
beta = cov / varB;
}
double alpha = meanA - (beta * meanB);
@@ -343,42 +313,23 @@ public sealed class Cointegration : AbstractBase
return double.NaN;
}
// Calculate gamma (coefficient in ADF regression)
// Δε_t = γ × ε_{t-1} + u_t
// γ = Cov(Δε, ε_{t-1}) / Var(ε_{t-1})
double meanDelta = _sumDelta / n;
double meanLagged = _sumLagged / n;
// Variance of lagged residuals
double varLagged = (_sumLagged2 / n) - (meanLagged * meanLagged);
if (Abs(varLagged) < Epsilon)
if (_sumLagged2 < Epsilon)
{
return double.NaN;
}
// Covariance of delta and lagged
double covDeltaLagged = (_sumDeltaLagged / n) - (meanDelta * meanLagged);
// No-intercept ADF regression: Δε_t = γ × ε_{t-1} + u_t
double gamma = _sumDeltaLagged / _sumLagged2;
// Gamma coefficient
double gamma = covDeltaLagged / varLagged;
// Calculate sum of squared regression errors in O(1)
// Sum((Δε_t - γ ε_{t-1})^2) = Sum(Δε_t^2) - 2γ Sum(Δε_t ε_{t-1}) + γ^2 Sum(ε_{t-1}^2)
double sumErrorSq = _sumDelta2 - (2.0 * gamma * _sumDeltaLagged) + (gamma * gamma * _sumLagged2);
// Calculate standard error of gamma
// SE(γ) = sqrt(Var(u) / (n × Var(ε_{t-1})))
// where u_t = Δε_t - γ × ε_{t-1}
// Ensure non-negative due to floating point errors
sumErrorSq = Max(0.0, sumErrorSq);
// Calculate sum of squared regression errors
double sumErrorSq = 0.0;
for (int i = 0; i < n; i++)
{
double delta = _deltaResiduals[i];
double lagged = _laggedResiduals[i];
double error = delta - (gamma * lagged);
sumErrorSq = FusedMultiplyAdd(error, error, sumErrorSq);
}
double varError = sumErrorSq / n;
double seGammaSq = varError / (n * varLagged);
double varError = sumErrorSq / (n - 1);
double seGammaSq = varError / _sumLagged2;
if (seGammaSq <= 0 || !double.IsFinite(seGammaSq))
{
@@ -386,7 +337,7 @@ public sealed class Cointegration : AbstractBase
}
double seGamma = Sqrt(seGammaSq);
if (Abs(seGamma) < Epsilon)
if (seGamma < Epsilon)
{
return double.NaN;
}
@@ -396,17 +347,20 @@ public sealed class Cointegration : AbstractBase
private void Resync()
{
// Resync main buffer sums
// Resync main buffer sums using span access to avoid per-element modulo in indexer.
// Both buffers are always updated together so their sequenced spans align element-by-element.
_sumA = 0;
_sumB = 0;
_sumA2 = 0;
_sumB2 = 0;
_sumAB = 0;
for (int i = 0; i < _bufferA.Count; i++)
_bufferA.GetSequencedSpans(out var aFirst, out var aSecond);
_bufferB.GetSequencedSpans(out var bFirst, out var bSecond);
for (int i = 0; i < aFirst.Length; i++)
{
double a = _bufferA[i];
double b = _bufferB[i];
double a = aFirst[i], b = bFirst[i];
_sumA += a;
_sumB += b;
_sumA2 = FusedMultiplyAdd(a, a, _sumA2);
@@ -414,20 +368,38 @@ public sealed class Cointegration : AbstractBase
_sumAB = FusedMultiplyAdd(a, b, _sumAB);
}
// Resync ADF regression sums
_sumDelta = 0;
_sumLagged = 0;
for (int i = 0; i < aSecond.Length; i++)
{
double a = aSecond[i], b = bSecond[i];
_sumA += a;
_sumB += b;
_sumA2 = FusedMultiplyAdd(a, a, _sumA2);
_sumB2 = FusedMultiplyAdd(b, b, _sumB2);
_sumAB = FusedMultiplyAdd(a, b, _sumAB);
}
// Resync ADF regression sums (delta/lagged buffers also always updated together).
_sumDeltaLagged = 0;
_sumLagged2 = 0;
_sumDelta2 = 0;
for (int i = 0; i < _deltaResiduals.Count; i++)
_deltaResiduals.GetSequencedSpans(out var dFirst, out var dSecond);
_laggedResiduals.GetSequencedSpans(out var lFirst, out var lSecond);
for (int i = 0; i < dFirst.Length; i++)
{
double delta = _deltaResiduals[i];
double lagged = _laggedResiduals[i];
_sumDelta += delta;
_sumLagged += lagged;
double delta = dFirst[i], lagged = lFirst[i];
_sumDeltaLagged = FusedMultiplyAdd(delta, lagged, _sumDeltaLagged);
_sumLagged2 = FusedMultiplyAdd(lagged, lagged, _sumLagged2);
_sumDelta2 = FusedMultiplyAdd(delta, delta, _sumDelta2);
}
for (int i = 0; i < dSecond.Length; i++)
{
double delta = dSecond[i], lagged = lSecond[i];
_sumDeltaLagged = FusedMultiplyAdd(delta, lagged, _sumDeltaLagged);
_sumLagged2 = FusedMultiplyAdd(lagged, lagged, _sumLagged2);
_sumDelta2 = FusedMultiplyAdd(delta, delta, _sumDelta2);
}
}
@@ -450,10 +422,9 @@ public sealed class Cointegration : AbstractBase
_sumB2 = 0;
_sumAB = 0;
_sumDelta = 0;
_sumLagged = 0;
_sumDeltaLagged = 0;
_sumLagged2 = 0;
_sumDelta2 = 0;
_prevResidual = 0;
_p_prevResidual = 0;
@@ -473,28 +444,7 @@ public sealed class Cointegration : AbstractBase
/// Calculates cointegration for two time series.
/// </summary>
public static TSeries Batch(TSeries seriesA, TSeries seriesB, int period = 20)
{
if (seriesA.Count != seriesB.Count)
{
throw new ArgumentException("Series must have the same length", nameof(seriesB));
}
var indicator = new Cointegration(period);
var result = new TSeries(seriesA.Count);
var timesA = seriesA.Times;
var valuesA = seriesA.Values;
var valuesB = seriesB.Values;
for (int i = 0; i < seriesA.Count; i++)
{
var tvalA = new TValue(timesA[i], valuesA[i]);
var tvalB = new TValue(timesA[i], valuesB[i]);
result.Add(indicator.Update(tvalA, tvalB, isNew: true));
}
return result;
}
=> Calculate(seriesA, seriesB, period).Results;
/// <summary>
/// Static batch calculation for span-based processing.
@@ -531,9 +481,24 @@ public sealed class Cointegration : AbstractBase
public static (TSeries Results, Cointegration Indicator) Calculate(TSeries seriesA, TSeries seriesB, int period = 20)
{
if (seriesA.Count != seriesB.Count)
{
throw new ArgumentException("Series must have the same length", nameof(seriesB));
}
var indicator = new Cointegration(period);
TSeries results = Batch(seriesA, seriesB, period);
return (results, indicator);
var result = new TSeries(seriesA.Count);
var timesA = seriesA.Times;
var valuesA = seriesA.Values;
var valuesB = seriesB.Values;
for (int i = 0; i < seriesA.Count; i++)
{
result.Add(indicator.Update(new TValue(timesA[i], valuesA[i]), new TValue(timesA[i], valuesB[i]), isNew: true));
}
return (result, indicator);
}
}
lib/statistics/cointegration/Cointegration.md
+4 -4
View File
@@ -3,9 +3,9 @@
| Property | Value |
| ---------------- | -------------------------------- |
| **Category** | Statistic |
| **Inputs** | Source (close) |
| **Inputs** | Two series (A, B) |
| **Parameters** | `period` (default 20) |
| **Outputs** | Single series (Cointegration) |
| **Outputs** | Single series (ADF statistic) |
| **Output range** | Varies (see docs) |
| **Warmup** | `period + 1` bars |
@@ -15,7 +15,7 @@
- Parameterized by `period` (default 20).
- Output range: Varies (see docs).
- Requires `period + 1` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against TradingView PineScript reference and statistical property tests.
> "Correlation tells you they move together. Cointegration tells you they're bound together. Two stocks can be uncorrelated yet cointegrated, or perfectly correlated yet destined to drift apart forever. The difference between 'similar direction' and 'shared destiny' is the difference between a tourist attraction and a gravitational orbit."
@@ -235,7 +235,7 @@ double[] pricesA = new double[1000];
double[] pricesB = new double[1000];
double[] output = new double[1000];
// ... populate inputs ...
Cointegration.Calculate(pricesA.AsSpan(), pricesB.AsSpan(), output.AsSpan(), period: 20);
Cointegration.Batch(pricesA.AsSpan(), pricesB.AsSpan(), output.AsSpan(), period: 20);
```
### Bar Correction Support
lib/statistics/correlation/Correlation.Quantower.Tests.cs
+1 -1
View File
@@ -304,4 +304,4 @@ public class CorrelationIndicatorTests
Assert.Equal(1.0, lastValue, precision: 6);
}
}
}
}
lib/statistics/correlation/Correlation.Tests.cs
+27 -23
View File
@@ -119,31 +119,33 @@ public class CorrelationTests
[Fact]
public void Update_IterativeCorrections_Restore()
{
var indicator = new Correlation(5);
var corrected = new Correlation(5);
var direct = new Correlation(5);
// Feed initial data
// Feed identical initial state
for (int i = 0; i < 8; i++)
{
double x = 100.0 + i;
double y = 200.0 + (i * 2);
indicator.Update(x, y, true);
corrected.Update(x, y, true);
direct.Update(x, y, true);
}
// Add new bar
indicator.Update(108.0, 216.0, true);
// Target final value for the current bar
const double finalX = 108.0;
const double finalY = 216.0;
// Make multiple corrections
for (int j = 0; j < 5; j++)
{
double x = 108.0 + (j * 0.1);
double y = 216.0 + (j * 0.2);
_ = indicator.Update(x, y, false);
}
// Correction path: new bar, several rewrites, final rewrite back to target
corrected.Update(finalX, finalY, true);
corrected.Update(finalX + 1.0, finalY + 2.0, false);
corrected.Update(finalX - 0.5, finalY - 1.0, false);
corrected.Update(finalX + 0.25, finalY + 0.5, false);
corrected.Update(finalX, finalY, false);
// Final correction back to original values
indicator.Update(108.0, 216.0, false);
// Direct path: same initial state + one new bar with final value
direct.Update(finalX, finalY, true);
Assert.True(double.IsFinite(indicator.Last.Value));
Assert.Equal(direct.Last.Value, corrected.Last.Value, 1e-12);
}
[Fact]
@@ -159,9 +161,9 @@ public class CorrelationTests
_ = indicator.Last.Value;
// Add NaN - should use last valid value
// Add NaN - should use last valid value, result must be finite
var result = indicator.Update(double.NaN, double.NaN, true);
Assert.True(double.IsFinite(result.Value) || double.IsNaN(result.Value));
Assert.True(double.IsFinite(result.Value));
}
[Fact]
@@ -175,9 +177,9 @@ public class CorrelationTests
indicator.Update(100.0 + i, 200.0 + i, true);
}
// Add Infinity - should use last valid value
// Add Infinity - should use last valid value, result must be finite
var result = indicator.Update(double.PositiveInfinity, double.NegativeInfinity, true);
Assert.True(double.IsFinite(result.Value) || double.IsNaN(result.Value));
Assert.True(double.IsFinite(result.Value));
}
[Fact]
@@ -189,11 +191,13 @@ public class CorrelationTests
}
[Fact]
public void IsHot_AtLeastTwoValues_ReturnsTrue()
public void IsHot_AtPeriod_ReturnsTrue()
{
var indicator = new Correlation(10);
indicator.Update(100.0, 200.0, true);
indicator.Update(101.0, 201.0, true);
for (int i = 0; i < 10; i++)
{
indicator.Update(100.0 + i, 200.0 + i, true);
}
Assert.True(indicator.IsHot);
}
@@ -385,4 +389,4 @@ public class CorrelationTests
}
}
}
}
}
lib/statistics/correlation/Correlation.Validation.Tests.cs
+13 -7
View File
@@ -609,14 +609,20 @@ public sealed class CorrelationValidationTests : IDisposable
[Fact]
public void Correlation_WeakCorrelation_DetectedCorrectly()
{
// Create two series with weak correlation (lots of noise)
// Create two series with weak correlation: pure independent noise, no shared trend.
// Use two independent GBMs (different seeds) and feed their incremental log-returns directly.
// With period=20 and fully independent noise sequences, correlation should be near zero.
var indicator = new Correlation(20);
var random = new GBM(startPrice: 100.0, sigma: 1.0, seed: 43);
var gbmX = new GBM(startPrice: 100.0, sigma: 0.2, seed: 43);
var gbmY = new GBM(startPrice: 100.0, sigma: 0.2, seed: 9871);
var barsX = gbmX.Fetch(101, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
var barsY = gbmY.Fetch(101, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
for (int i = 0; i < 100; i++)
for (int i = 1; i <= 100; i++)
{
double x = 100.0 + i + Math.Log(random.Next().Close / 100.0) * 50;
double y = 100.0 + 0.1 * i + Math.Log(random.Next().Close / 100.0) * 50; // Weak relationship
// Pure independent white noise — no shared linear component
double x = Math.Log(barsX[i].Close / barsX[i - 1].Close);
double y = Math.Log(barsY[i].Close / barsY[i - 1].Close);
indicator.Update(x, y);
}
@@ -720,7 +726,7 @@ public sealed class CorrelationValidationTests : IDisposable
double[] taOut = new double[_data.Count];
var retCode = Functions.Correl<double>(closeArr, openArr, 0..^0, taOut, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
(int offset, int length) = outRange.GetOffsetAndLength(taOut.Length);
Assert.True(length > 100, $"TALib Correl produced only {length} values");
@@ -760,7 +766,7 @@ public sealed class CorrelationValidationTests : IDisposable
{
double[] taOut = new double[_data.Count];
var retCode = Functions.Correl<double>(highArr, lowArr, 0..^0, taOut, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
(int offset, int length) = outRange.GetOffsetAndLength(taOut.Length);
lib/statistics/correlation/Correlation.cs
+39 -38
View File
@@ -39,12 +39,13 @@ public sealed class Correlation : AbstractBase
// Last valid values for NaN handling
private double _lastValidX, _lastValidY;
private double _p_lastValidX, _p_lastValidY;
private int _updateCount;
private const int ResyncInterval = 1000;
private const double Epsilon = 1e-10;
public override bool IsHot => _bufferX.Count >= 2;
public override bool IsHot => _bufferX.Count >= WarmupPeriod;
/// <summary>
/// Creates a new Correlation indicator.
@@ -74,6 +75,17 @@ public sealed class Correlation : AbstractBase
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public TValue Update(TValue seriesX, TValue seriesY, bool isNew = true)
{
if (isNew)
{
_p_lastValidX = _lastValidX;
_p_lastValidY = _lastValidY;
}
else
{
_lastValidX = _p_lastValidX;
_lastValidY = _p_lastValidY;
}
double x = SanitizeX(seriesX.Value);
double y = SanitizeY(seriesY.Value);
@@ -99,7 +111,7 @@ public sealed class Correlation : AbstractBase
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public TValue Update(double seriesX, double seriesY, bool isNew = true)
{
return Update(new TValue(DateTime.UtcNow, seriesX), new TValue(DateTime.UtcNow, seriesY), isNew);
return Update(new TValue(DateTime.MinValue, seriesX), new TValue(DateTime.MinValue, seriesY), isNew);
}
/// <inheritdoc/>
@@ -175,14 +187,7 @@ public sealed class Correlation : AbstractBase
{
if (_bufferX.Count == 0)
{
// No data yet, just add
_bufferX.Add(x);
_bufferY.Add(y);
_sumX = x;
_sumY = y;
_sumX2 = x * x;
_sumY2 = y * y;
_sumXY = x * y;
// Nothing to correct yet; no current bar exists
return;
}
@@ -190,12 +195,12 @@ public sealed class Correlation : AbstractBase
double oldX = _bufferX.Newest;
double oldY = _bufferY.Newest;
// Update the running sums: remove old, add new
// Update the running sums: remove old, add new (using FMA for consistency with ProcessNewBar)
_sumX = _sumX - oldX + x;
_sumY = _sumY - oldY + y;
_sumX2 = _sumX2 - (oldX * oldX) + (x * x);
_sumY2 = _sumY2 - (oldY * oldY) + (y * y);
_sumXY = _sumXY - (oldX * oldY) + (x * y);
_sumX2 = FusedMultiplyAdd(x, x, FusedMultiplyAdd(-oldX, oldX, _sumX2));
_sumY2 = FusedMultiplyAdd(y, y, FusedMultiplyAdd(-oldY, oldY, _sumY2));
_sumXY = FusedMultiplyAdd(x, y, FusedMultiplyAdd(-oldX, oldY, _sumXY));
// Update the buffer values
_bufferX.UpdateNewest(x);
@@ -278,6 +283,8 @@ public sealed class Correlation : AbstractBase
_lastValidX = 0;
_lastValidY = 0;
_p_lastValidX = 0;
_p_lastValidY = 0;
_updateCount = 0;
Last = default;
@@ -287,28 +294,7 @@ public sealed class Correlation : AbstractBase
/// Calculates correlation for two time series.
/// </summary>
public static TSeries Batch(TSeries seriesX, TSeries seriesY, int period = 20)
{
if (seriesX.Count != seriesY.Count)
{
throw new ArgumentException("Series must have the same length", nameof(seriesY));
}
var indicator = new Correlation(period);
var result = new TSeries(seriesX.Count);
var timesX = seriesX.Times;
var valuesX = seriesX.Values;
var valuesY = seriesY.Values;
for (int i = 0; i < seriesX.Count; i++)
{
var tvalX = new TValue(timesX[i], valuesX[i]);
var tvalY = new TValue(timesX[i], valuesY[i]);
result.Add(indicator.Update(tvalX, tvalY, isNew: true));
}
return result;
}
=> Calculate(seriesX, seriesY, period).Results;
/// <summary>
/// Static batch calculation for span-based processing.
@@ -345,9 +331,24 @@ public sealed class Correlation : AbstractBase
public static (TSeries Results, Correlation Indicator) Calculate(TSeries seriesX, TSeries seriesY, int period = 20)
{
if (seriesX.Count != seriesY.Count)
{
throw new ArgumentException("Series must have the same length", nameof(seriesY));
}
var indicator = new Correlation(period);
TSeries results = Batch(seriesX, seriesY, period);
return (results, indicator);
var result = new TSeries(seriesX.Count);
var timesX = seriesX.Times;
var valuesX = seriesX.Values;
var valuesY = seriesY.Values;
for (int i = 0; i < seriesX.Count; i++)
{
result.Add(indicator.Update(new TValue(timesX[i], valuesX[i]), new TValue(timesX[i], valuesY[i]), isNew: true));
}
return (result, indicator);
}
}
lib/statistics/correlation/Correlation.md
+4 -4
View File
@@ -3,9 +3,9 @@
| Property | Value |
| ---------------- | -------------------------------- |
| **Category** | Statistic |
| **Inputs** | Source (close) |
| **Inputs** | Two series (X, Y) |
| **Parameters** | `period` (default 20) |
| **Outputs** | Single series (Correlation) |
| **Outputs** | Single series (Pearson r) |
| **Output range** | Varies (see docs) |
| **Warmup** | `period` bars |
@@ -15,7 +15,7 @@
- Parameterized by `period` (default 20).
- Output range: Varies (see docs).
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against TradingView reference behavior and mathematical invariants.
> "Correlation is not causation, but it sure is a hint. The market doesn't care why two instruments move together—only that they do, and whether that relationship will persist long enough for you to profit from it."
@@ -219,7 +219,7 @@ double[] pricesA = new double[1000];
double[] pricesB = new double[1000];
double[] output = new double[1000];
// ... populate inputs ...
Correlation.Calculate(pricesA.AsSpan(), pricesB.AsSpan(), output.AsSpan(), period: 20);
Correlation.Batch(pricesA.AsSpan(), pricesB.AsSpan(), output.AsSpan(), period: 20);
```
### Bar Correction Support
lib/statistics/hurst/Hurst.Validation.Tests.cs
+1 -1
View File
@@ -181,4 +181,4 @@ public sealed class HurstValidationTests
Assert.Equal(h1.Last.Value, h2.Last.Value, 1e-15);
}
}
}
lib/statistics/pacf/Pacf.Quantower.Tests.cs
+1 -1
View File
@@ -325,4 +325,4 @@ public class PacfIndicatorTests
Assert.Equal(acfValue, pacfValue, 6); // Allow for minor floating-point differences
}
}
}
lib/statistics/pacf/Pacf.Tests.cs
+1 -1
View File
@@ -575,4 +575,4 @@ public class PacfTests
}
#endregion
}
}
lib/statistics/pacf/Pacf.Validation.Tests.cs
+8 -5
View File
@@ -105,13 +105,16 @@ public class PacfValidationTests
// For AR(1) process: x_t = φ*x_{t-1} + ε_t
// PACF should be significant at lag 1 and cut off (near zero) after
double phi = 0.7; // AR(1) coefficient
var random = new GBM(startPrice: 100.0, sigma: 1.0, seed: 42);
var arProcess = new List<double> { 100.0 };
// Generate AR(1) process
// Use incremental bar-to-bar log-returns as i.i.d. noise: log(close_i / close_{i-1})
var gbm = new GBM(startPrice: 100.0, sigma: 0.2, seed: 42);
var bars = gbm.Fetch(501, DateTime.UtcNow.Ticks, TimeSpan.FromMinutes(1));
var arProcess = new List<double> { 0.0 };
// Generate AR(1) process using incremental log-returns as white noise ε
for (int i = 1; i < 500; i++)
{
double noise = Math.Log(random.Next().Close / 100.0); // ~N(0, vol²*dt) noise
double noise = Math.Log(bars[i].Close / bars[i - 1].Close); // i.i.d. incremental return
double newValue = phi * arProcess[^1] + noise;
arProcess.Add(newValue);
}
@@ -345,4 +348,4 @@ public class PacfValidationTests
}
#endregion
}
}
lib/statistics/pacf/Pacf.cs
+2 -1
View File
@@ -294,9 +294,10 @@ public sealed class Pacf : AbstractBase
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
Reset();
foreach (double value in source)
{
Update(new TValue(DateTime.UtcNow, value));
Update(new TValue(DateTime.MinValue, value));
}
}
lib/statistics/pacf/Pacf.md
+1 -1
View File
@@ -15,7 +15,7 @@
- Parameterized by `period`, `lag` (default 1).
- Output range: Varies (see docs).
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against mathematical properties and Durbin-Levinson recursion expectations.
> "Strip away the intermediaries, and you'll see the true direct relationship."
lib/statistics/spearman/Spearman.Validation.Tests.cs
+1 -1
View File
@@ -124,4 +124,4 @@ public sealed class SpearmanValidationTests
int finiteCount = values.Count(v => double.IsFinite(v));
Assert.True(finiteCount > 100, $"Expected >100 finite values, got {finiteCount}");
}
}
}
lib/statistics/stderr/Stderr.md
+14 -14
View File
@@ -11,11 +11,11 @@
### TL;DR
- ````markdown
- `Stderr` computes the standard error of an OLS regression fit over a rolling window.
- Parameterized by `period`.
- Output range: Varies (see docs).
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Output range: non-negative real values (or 0 during insufficient/degenerate windows).
- Requires `period` bars of warmup before first stable output (`IsHot = true`).
- Validated against an internal brute-force OLS reference implementation.
> "How confident are you in your line of best fit?"
@@ -58,16 +58,16 @@ $$ b = \frac{\sum y - m \sum x}{N} $$
### Operation Count (Streaming Mode)
Standard Error = StdDev / sqrt(N), computed atop the O(1) StdDev computation.
`Stderr` keeps regression sums in O(1), then performs an O(N) residual pass to compute SSR.
| Operation | Count | Cost (cycles) | Subtotal |
| :--- | :---: | :---: | :---: |
| O(1) StdDev computation | 1 | 28 cy | ~28 cy |
| Divide by sqrt(N) (precomputed) | 1 | 4 cy | ~4 cy |
| NaN guard + state update | 1 | 2 cy | ~2 cy |
| **Total** | **O(1)** | — | **~34 cy** |
| Running-sum updates | O(1) | — | small |
| Residual SSR scan | O(N) | dominant | dominant |
| Final sqrt/divide | O(1) | — | small |
| **Total** | **O(N)** | — | period-dependent |
O(1) per update. sqrt(N) is precomputed in the constructor. Negligible additional cost over StdDev.
Per-update complexity is O(N) because residuals must be re-evaluated for the current window.
| Metric | Score | Notes |
| :--- | :--- | :--- |
@@ -80,8 +80,9 @@ O(1) per update. sqrt(N) is precomputed in the constructor. Negligible additiona
| Library | Status | Notes |
| :--- | :--- | :--- |
| **TA-Lib** | ✅ | Matches `STDERR` output. |
| **TradingView** | ✅ | Matches Pine Script `ta.stdev` of residuals. |
| **TA-Lib** | ⚠️ | Formula differs (`stderr` in Tulip/other libs often means standard error of mean). |
| **TradingView** | ✅ | Matches Pine-style OLS residual standard error behavior for this implementation. |
| **Reference OLS** | ✅ | Cross-validated against brute-force OLS residual calculation. |
## Usage
@@ -101,5 +102,4 @@ double value = stderr.Last.Value;
## See Also
* **LinReg** — Linear Regression Curve (the trend line itself).
* **StdDev** — Standard Deviation (dispersion from the mean, not from a regression line).
````
- **StdDev** — Standard Deviation (dispersion from the mean, not from a regression line).
lib/statistics/sum/Sum.Validation.Tests.cs
+3 -3
View File
@@ -51,7 +51,7 @@ public sealed class SumValidationTests : IDisposable
var qResult = sum.Update(_testData.Data);
var retCode = Functions.Sum<double>(tData, 0..^0, output, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
int lookback = Functions.SumLookback(period);
@@ -77,7 +77,7 @@ public sealed class SumValidationTests : IDisposable
}
var retCode = Functions.Sum<double>(tData, 0..^0, output, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
int lookback = Functions.SumLookback(period);
@@ -99,7 +99,7 @@ public sealed class SumValidationTests : IDisposable
Sum.Batch(sourceData.AsSpan(), qOutput.AsSpan(), period);
var retCode = Functions.Sum<double>(sourceData, 0..^0, tOutput, out var outRange, period);
Assert.Equal(Core.RetCode.Success, retCode);
Assert.Equal(TALib.Core.RetCode.Success, retCode);
int lookback = Functions.SumLookback(period);
lib/statistics/variance/Variance.cs
+2 -1
View File
@@ -162,9 +162,10 @@ public sealed class Variance : AbstractBase
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
Reset();
foreach (double value in source)
{
Update(new TValue(DateTime.UtcNow, value));
Update(new TValue(DateTime.MinValue, value));
}
}
lib/statistics/zscore/Zscore.Validation.Tests.cs
+1 -1
View File
@@ -139,4 +139,4 @@ public sealed class ZscoreValidationTests
int finiteCount = values.Count(v => double.IsFinite(v));
Assert.True(finiteCount > 100, $"Expected >100 finite values, got {finiteCount}");
}
}
}
lib/statistics/zscore/Zscore.cs
+111 -61
View File
@@ -27,13 +27,16 @@ public sealed class Zscore : AbstractBase
private readonly RingBuffer _buffer;
private readonly TValuePublishedHandler _handler;
private double _lastValidValue;
[StructLayout(LayoutKind.Auto)]
private record struct State(double LastValidZScore, double LastValidValue);
private State _s, _ps;
private double _sumSq;
private double _p_sumSq;
private int _updateCount;
private const int ResyncInterval = 1000;
public override bool IsHot => _buffer.Count >= _period;
/// <summary>
/// Initializes a rolling z-score indicator.
/// </summary>
/// <param name="period">Lookback period (default 14, must be >= 2)</param>
public Zscore(int period = 14)
{
@@ -46,11 +49,14 @@ public sealed class Zscore : AbstractBase
_buffer = new RingBuffer(period);
Name = $"Zscore({period})";
WarmupPeriod = period;
_s = new State(0.0, 0.0);
_ps = _s;
_sumSq = 0.0;
_p_sumSq = 0.0;
_handler = Handle;
}
/// <summary>
/// Initializes a rolling z-score indicator and subscribes it to a source publisher.
/// </summary>
/// <param name="source">Source indicator for event-based chaining</param>
/// <param name="period">Lookback period (default 14)</param>
public Zscore(ITValuePublisher source, int period = 14) : this(period)
@@ -61,16 +67,6 @@ public sealed class Zscore : AbstractBase
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public override TValue Update(TValue input, bool isNew = true)
{
if (isNew)
{
_ps = _s;
}
else
{
_s = _ps;
_lastValidValue = _s.LastValidValue;
}
double value = input.Value;
if (!double.IsFinite(value))
@@ -82,11 +78,37 @@ public sealed class Zscore : AbstractBase
_lastValidValue = value;
}
_buffer.Add(value, isNew);
if (isNew)
{
_p_sumSq = _sumSq;
_buffer.Snapshot();
}
else
{
_sumSq = _p_sumSq;
_buffer.Restore();
}
if (_buffer.IsFull)
{
double oldVal = _buffer.Oldest;
_sumSq = Math.FusedMultiplyAdd(-oldVal, oldVal, _sumSq);
}
_buffer.Add(value);
_sumSq = Math.FusedMultiplyAdd(value, value, _sumSq);
if (isNew)
{
_updateCount++;
if (_updateCount % ResyncInterval == 0)
{
Resync();
}
}
double result;
ReadOnlySpan<double> data = _buffer.GetSpan();
int n = data.Length;
int n = _buffer.Count;
if (n < 2)
{
@@ -94,25 +116,16 @@ public sealed class Zscore : AbstractBase
}
else
{
double sum = 0.0;
double sumSq = 0.0;
for (int i = 0; i < n; i++)
{
double v = data[i];
sum += v;
sumSq += v * v;
}
double sum = _buffer.Sum;
double mean = sum / n;
// Population variance: E[X²] - (E[X])²
double popVariance = (sumSq / n) - (mean * mean);
if (popVariance < 0.0)
double numerator = _sumSq - (sum * sum) / n;
if (numerator < 0)
{
popVariance = 0.0;
numerator = 0;
}
double popVariance = numerator / n;
double stdDev = Math.Sqrt(popVariance);
if (stdDev > 1e-10)
@@ -125,7 +138,6 @@ public sealed class Zscore : AbstractBase
}
}
_s = new State(result, _lastValidValue);
Last = new TValue(input.Time, result);
PubEvent(Last, isNew);
return Last;
@@ -133,17 +145,30 @@ public sealed class Zscore : AbstractBase
public override TSeries Update(TSeries source)
{
var result = new TSeries(source.Count);
ReadOnlySpan<double> values = source.Values;
ReadOnlySpan<long> times = source.Times;
for (int i = 0; i < source.Count; i++)
if (source.Count == 0)
{
var tv = Update(new TValue(new DateTime(times[i], DateTimeKind.Utc), values[i]), true);
result.Add(tv, true);
return new TSeries();
}
return result;
int len = source.Count;
var t = new List<long>(len);
var v = new List<double>(len);
CollectionsMarshal.SetCount(t, len);
CollectionsMarshal.SetCount(v, len);
var tSpan = CollectionsMarshal.AsSpan(t);
var vSpan = CollectionsMarshal.AsSpan(v);
Batch(source.Values, vSpan, _buffer.Capacity);
source.Times.CopyTo(tSpan);
int primeStart = Math.Max(0, len - _buffer.Capacity);
for (int i = primeStart; i < len; i++)
{
Update(source[i]);
}
return new TSeries(t, v);
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
@@ -153,11 +178,24 @@ public sealed class Zscore : AbstractBase
{
_buffer.Clear();
_lastValidValue = 0;
_s = new State(0.0, 0.0);
_ps = _s;
_sumSq = 0.0;
_p_sumSq = 0.0;
_updateCount = 0;
Last = default;
}
private void Resync()
{
var span = _buffer.GetSpan();
double sumSq = 0;
for (int i = 0; i < span.Length; i++)
{
sumSq += span[i] * span[i];
}
_sumSq = sumSq;
_buffer.RecalculateSum();
}
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
TimeSpan interval = step ?? TimeSpan.FromSeconds(1);
@@ -214,6 +252,8 @@ public sealed class Zscore : AbstractBase
int head = 0;
int count = 0;
double lastValid = 0.0;
double sum = 0.0;
double sumSq = 0.0;
for (int i = 0; i < source.Length; i++)
{
@@ -228,43 +268,53 @@ public sealed class Zscore : AbstractBase
lastValid = val;
}
if (count < ringSize)
if (count == ringSize)
{
ring[count] = val;
count++;
double oldVal = ring[head];
sum = Math.FusedMultiplyAdd(-1.0, oldVal, sum + val);
sumSq = Math.FusedMultiplyAdd(-oldVal, oldVal, sumSq);
}
else
{
ring[head] = val;
count++;
sum += val;
}
ring[head] = val;
sumSq = Math.FusedMultiplyAdd(val, val, sumSq);
head = (head + 1) % ringSize;
if ((i + 1) % 1000 == 0 && count == ringSize)
{
double resyncSum = 0;
double resyncSumSq = 0;
for (int j = 0; j < ringSize; j++)
{
double v = ring[j];
resyncSum += v;
resyncSumSq += v * v;
}
sum = resyncSum;
sumSq = resyncSumSq;
}
if (count < 2)
{
output[i] = 0.0;
continue;
}
double sum = 0.0;
double sumSq = 0.0;
int n = count;
for (int j = 0; j < n; j++)
{
double v = ring[j];
sum += v;
sumSq += v * v;
}
double mean = sum / n;
double popVariance = (sumSq / n) - (mean * mean);
if (popVariance < 0.0)
double numerator = sumSq - (sum * sum) / n;
if (numerator < 0)
{
popVariance = 0.0;
numerator = 0;
}
double popVariance = numerator / n;
double stdDev = Math.Sqrt(popVariance);
if (stdDev > 1e-10)
lib/statistics/zscore/Zscore.md
+1 -1
View File
@@ -15,7 +15,7 @@
- Parameterized by `period` (default 14).
- Output range: Unbounded.
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against manual computation, PineScript parity, and statistical invariants.
> "How far from normal is this?" — Every risk manager, every day.
lib/statistics/ztest/Ztest.Tests.cs
+1 -1
View File
@@ -354,7 +354,7 @@ public class ZtestTests
for (int i = 0; i < count; i++)
{
Assert.Equal(batchResult[i].Value, spanOutput[i], 1e-4); // t-stat magnifies FP drift (values ~6000)
Assert.Equal(batchResult[i].Value, spanOutput[i], 5e-4); // t-stat magnifies FP drift (values ~15000)
}
}
lib/statistics/ztest/Ztest.cs
+114 -63
View File
@@ -29,13 +29,16 @@ public sealed class Ztest : AbstractBase
private readonly RingBuffer _buffer;
private readonly TValuePublishedHandler _handler;
private double _lastValidValue;
[StructLayout(LayoutKind.Auto)]
private record struct State(double LastValidTStat, double LastValidValue);
private State _s, _ps;
private double _sumSq;
private double _p_sumSq;
private int _updateCount;
private const int ResyncInterval = 1000;
public override bool IsHot => _buffer.Count >= _period;
/// <summary>
/// Initializes a rolling one-sample t-test indicator.
/// </summary>
/// <param name="period">Lookback period (default 30, must be >= 2)</param>
/// <param name="mu0">Hypothesized population mean (default 0.0)</param>
public Ztest(int period = 30, double mu0 = 0.0)
@@ -50,11 +53,14 @@ public sealed class Ztest : AbstractBase
_buffer = new RingBuffer(period);
Name = $"Ztest({period},{mu0:G})";
WarmupPeriod = period;
_s = new State(0.0, 0.0);
_ps = _s;
_sumSq = 0.0;
_p_sumSq = 0.0;
_handler = Handle;
}
/// <summary>
/// Initializes a rolling one-sample t-test indicator and subscribes it to a source publisher.
/// </summary>
/// <param name="source">Source indicator for event-based chaining</param>
/// <param name="period">Lookback period (default 30)</param>
/// <param name="mu0">Hypothesized population mean (default 0.0)</param>
@@ -66,16 +72,6 @@ public sealed class Ztest : AbstractBase
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public override TValue Update(TValue input, bool isNew = true)
{
if (isNew)
{
_ps = _s;
}
else
{
_s = _ps;
_lastValidValue = _s.LastValidValue;
}
double value = input.Value;
if (!double.IsFinite(value))
@@ -87,11 +83,37 @@ public sealed class Ztest : AbstractBase
_lastValidValue = value;
}
_buffer.Add(value, isNew);
if (isNew)
{
_p_sumSq = _sumSq;
_buffer.Snapshot();
}
else
{
_sumSq = _p_sumSq;
_buffer.Restore();
}
if (_buffer.IsFull)
{
double oldVal = _buffer.Oldest;
_sumSq = Math.FusedMultiplyAdd(-oldVal, oldVal, _sumSq);
}
_buffer.Add(value);
_sumSq = Math.FusedMultiplyAdd(value, value, _sumSq);
if (isNew)
{
_updateCount++;
if (_updateCount % ResyncInterval == 0)
{
Resync();
}
}
double result;
ReadOnlySpan<double> data = _buffer.GetSpan();
int n = data.Length;
int n = _buffer.Count;
if (n < 2)
{
@@ -99,27 +121,19 @@ public sealed class Ztest : AbstractBase
}
else
{
double sum = 0.0;
double sumSq = 0.0;
for (int i = 0; i < n; i++)
{
double v = data[i];
sum += v;
sumSq += v * v;
}
double sum = _buffer.Sum;
double mean = sum / n;
// Population variance first: E[X²] - (E[X])²
double popVariance = (sumSq / n) - (mean * mean);
if (popVariance < 0.0)
double numerator = _sumSq - (sum * sum) / n;
if (numerator < 0)
{
popVariance = 0.0;
numerator = 0;
}
// Bessel correction: sample variance = popVariance * n / (n - 1)
double sampleStdDev = Math.Sqrt(popVariance * n / (n - 1));
// which is numerator / (n - 1)
double sampleVariance = numerator / (n - 1);
double sampleStdDev = Math.Sqrt(sampleVariance);
double standardError = sampleStdDev / Math.Sqrt(n);
if (standardError > 1e-10)
@@ -132,7 +146,6 @@ public sealed class Ztest : AbstractBase
}
}
_s = new State(result, _lastValidValue);
Last = new TValue(input.Time, result);
PubEvent(Last, isNew);
return Last;
@@ -140,17 +153,30 @@ public sealed class Ztest : AbstractBase
public override TSeries Update(TSeries source)
{
var result = new TSeries(source.Count);
ReadOnlySpan<double> values = source.Values;
ReadOnlySpan<long> times = source.Times;
for (int i = 0; i < source.Count; i++)
if (source.Count == 0)
{
var tv = Update(new TValue(new DateTime(times[i], DateTimeKind.Utc), values[i]), true);
result.Add(tv, true);
return new TSeries();
}
return result;
int len = source.Count;
var t = new List<long>(len);
var v = new List<double>(len);
CollectionsMarshal.SetCount(t, len);
CollectionsMarshal.SetCount(v, len);
var tSpan = CollectionsMarshal.AsSpan(t);
var vSpan = CollectionsMarshal.AsSpan(v);
Batch(source.Values, vSpan, _buffer.Capacity, _mu0);
source.Times.CopyTo(tSpan);
int primeStart = Math.Max(0, len - _buffer.Capacity);
for (int i = primeStart; i < len; i++)
{
Update(source[i]);
}
return new TSeries(t, v);
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
@@ -160,11 +186,24 @@ public sealed class Ztest : AbstractBase
{
_buffer.Clear();
_lastValidValue = 0;
_s = new State(0.0, 0.0);
_ps = _s;
_sumSq = 0.0;
_p_sumSq = 0.0;
_updateCount = 0;
Last = default;
}
private void Resync()
{
var span = _buffer.GetSpan();
double sumSq = 0;
for (int i = 0; i < span.Length; i++)
{
sumSq += span[i] * span[i];
}
_sumSq = sumSq;
_buffer.RecalculateSum();
}
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
TimeSpan interval = step ?? TimeSpan.FromSeconds(1);
@@ -221,6 +260,8 @@ public sealed class Ztest : AbstractBase
int head = 0;
int count = 0;
double lastValid = 0.0;
double sum = 0.0;
double sumSq = 0.0;
for (int i = 0; i < source.Length; i++)
{
@@ -235,45 +276,55 @@ public sealed class Ztest : AbstractBase
lastValid = val;
}
if (count < ringSize)
if (count == ringSize)
{
ring[count] = val;
count++;
double oldVal = ring[head];
sum = Math.FusedMultiplyAdd(-1.0, oldVal, sum + val);
sumSq = Math.FusedMultiplyAdd(-oldVal, oldVal, sumSq);
}
else
{
ring[head] = val;
count++;
sum += val;
}
ring[head] = val;
sumSq = Math.FusedMultiplyAdd(val, val, sumSq);
head = (head + 1) % ringSize;
if ((i + 1) % 1000 == 0 && count == ringSize)
{
double resyncSum = 0;
double resyncSumSq = 0;
for (int j = 0; j < ringSize; j++)
{
double v = ring[j];
resyncSum += v;
resyncSumSq += v * v;
}
sum = resyncSum;
sumSq = resyncSumSq;
}
if (count < 2)
{
output[i] = 0.0;
continue;
}
double sum = 0.0;
double sumSq = 0.0;
int n = count;
for (int j = 0; j < n; j++)
{
double v = ring[j];
sum += v;
sumSq += v * v;
}
double mean = sum / n;
double popVariance = (sumSq / n) - (mean * mean);
if (popVariance < 0.0)
double numerator = sumSq - (sum * sum) / n;
if (numerator < 0)
{
popVariance = 0.0;
numerator = 0;
}
// Bessel correction: sample variance = popVariance * n / (n - 1)
double sampleStdDev = Math.Sqrt(popVariance * n / (n - 1));
double sampleVariance = numerator / (n - 1);
double sampleStdDev = Math.Sqrt(sampleVariance);
double standardError = sampleStdDev / Math.Sqrt(n);
if (standardError > 1e-10)
lib/statistics/ztest/Ztest.md
+4 -13
View File
@@ -15,7 +15,7 @@
- Parameterized by `period` (default 30), `mu0` (default 0.0).
- Output range: Unbounded.
- Requires `period` bars of warmup before first valid output (IsHot = true).
- Validated against TA-Lib, Skender, and Tulip reference implementations where available.
- Validated against manual computation, PineScript parity, and testable statistical properties.
> "The purpose of hypothesis testing is not to prove what we believe, but to measure what we observe." — Adapted from R.A. Fisher
@@ -91,21 +91,12 @@ The indicators answer different questions:
### Operation Count (Streaming Mode)
Z-Test computes a rolling mean and standard deviation for O(1) hypothesis testing per bar.
| Operation | Count | Cost (cycles) | Subtotal |
| :--- | :---: | :---: | :---: |
| O(1) StdDev computation | 1 | 28 cy | ~28 cy |
| Compute Z = (x - mu) / (sigma / sqrt(N)) | 1 | 5 cy | ~5 cy |
| NaN guard (sigma = 0 guard) | 1 | 2 cy | ~2 cy |
| **Total** | **O(1)** | — | **~35 cy** |
O(1) per update. Z-statistic is a trivial transformation of the running mean and standard deviation already computed by StdDev.
ZTEST uses a rolling window with running sums and periodic resynchronization.
| Operation | Complexity | Notes |
|-----------|-----------|-------|
| Update (streaming) | $O(n)$ | Full window scan for sum/sumSq |
| Batch (span) | $O(N \cdot p)$ | N data points, p period |
| Update (streaming) | $O(1)$ amortized | Running sum/sumSq maintenance; periodic full resync every 1000 updates |
| Batch (span) | $O(N)$ | Single pass over source with O(1) ring maintenance per element |
| Memory | $O(p)$ | RingBuffer + scalar state |
| Allocations per update | 0 | Zero-allocation hot path |