QuanTAlib/lib/statistics/stderr/tests/Stderr.Validation.Tests.cs

using Tulip;
using Xunit;

namespace QuanTAlib.Tests;

/// <summary>
/// Stderr cross-validation against pure-C# reference implementation.
/// The reference exactly replicates the OLS formula in the pine script.
/// Also cross-validated against Tulip <c>stderr</c> (Standard Error of Linear Regression)
/// — exact formula match: sqrt(SSR / (n-2)).
/// </summary>
public class StderrValidationTests
{
    // ─────────────────────────────────────────────────────────────
    // Reference: brute-force OLS over an explicit window array
    // ─────────────────────────────────────────────────────────────
    private static double ReferenceStderr(double[] window)
    {
        int n = window.Length;
        if (n < 3)
        {
            return 0;
        }

        double sumX = 0, sumY = 0, sumXY = 0, sumX2 = 0;
        for (int i = 0; i < n; i++)
        {
            sumX += i;
            sumY += window[i];
            sumXY += i * window[i];
            sumX2 += (double)i * i;
        }

        double denom = n * sumX2 - sumX * sumX;
        if (denom == 0)
        {
            return 0;
        }

        double slope = (n * sumXY - sumX * sumY) / denom;
        double intercept = (sumY - slope * sumX) / n;

        double ssr = 0;
        for (int i = 0; i < n; i++)
        {
            double predicted = slope * i + intercept;
            double res = window[i] - predicted;
            ssr += res * res;
        }

        return Math.Sqrt(ssr / (n - 2.0));
    }

    [Fact]
    public void Stderr_KnownLinearData_IsZero()
    {
        // Perfect linear trend → residuals = 0 → Stderr = 0
        var se = new Stderr(5);
        for (int i = 0; i < 5; i++)
        {
            se.Update(new TValue(DateTime.UtcNow, i * 3.0 + 2.0));
        }
        Assert.Equal(0.0, se.Last.Value, precision: 8);
    }

    [Fact]
    public void Stderr_KnownData_Manual()
    {
        // y = {2, 4, 5}: reference computed in test B
        double expected = ReferenceStderr(new double[] { 2, 4, 5 });
        var se = new Stderr(3);
        se.Update(new TValue(DateTime.UtcNow, 2.0));
        se.Update(new TValue(DateTime.UtcNow, 4.0));
        se.Update(new TValue(DateTime.UtcNow, 5.0));
        Assert.Equal(expected, se.Last.Value, precision: 10);
    }

    [Fact]
    public void Stderr_Batch_Matches_Reference_GBM()
    {
        const int period = 14;
        var gbm = new GBM(seed: 12345);
        var closes = new List<double>();
        var series = new TSeries();
        for (int i = 0; i < 300; i++)
        {
            var bar = gbm.Next();
            closes.Add(bar.Close);
            series.Add(new TValue(bar.Time, bar.Close));
        }

        var result = Stderr.Batch(series, period);

        for (int i = period - 1; i < closes.Count; i++)
        {
            double[] window = closes.Skip(i - period + 1).Take(period).ToArray();
            double expected = ReferenceStderr(window);
            Assert.Equal(expected, result[i].Value, precision: 8);
        }
    }

    [Fact]
    public void Stderr_Streaming_Matches_Reference_GBM()
    {
        const int period = 20;
        var gbm = new GBM(seed: 54321);
        var closes = new List<double>();
        var se = new Stderr(period);

        for (int i = 0; i < 200; i++)
        {
            var bar = gbm.Next();
            closes.Add(bar.Close);
            se.Update(new TValue(bar.Time, bar.Close));

            if (i >= period - 1)
            {
                double[] window = closes.Skip(i - period + 1).Take(period).ToArray();
                double expected = ReferenceStderr(window);
                Assert.Equal(expected, se.Last.Value, precision: 8);
            }
        }
    }

    [Fact]
    public void Stderr_Span_Matches_Reference_GBM()
    {
        const int period = 10;
        var gbm = new GBM(seed: 999);
        var closes = new List<double>();
        for (int i = 0; i < 100; i++)
        {
            closes.Add(gbm.Next().Close);
        }

        var src = closes.ToArray();
        var dst = new double[src.Length];
        Stderr.Batch(src.AsSpan(), dst.AsSpan(), period);

        for (int i = period - 1; i < closes.Count; i++)
        {
            double[] window = closes.Skip(i - period + 1).Take(period).ToArray();
            double expected = ReferenceStderr(window);
            Assert.Equal(expected, dst[i], precision: 8);
        }
    }

    [Fact]
    public void Stderr_SlidingWindow_CorrectlyDropsOldest()
    {
        // Feed 6 values with period=4. Verify last two windows.
        const int period = 4;
        double[] data = { 1, 3, 2, 5, 4, 6 };
        var se = new Stderr(period);
        for (int i = 0; i < data.Length; i++)
        {
            se.Update(new TValue(DateTime.UtcNow, data[i]));
        }
        double expected = ReferenceStderr(new double[] { 2, 5, 4, 6 });
        Assert.Equal(expected, se.Last.Value, precision: 8);
    }

    [Fact]
    public void Stderr_AlwaysNonNegative()
    {
        const int period = 14;
        var gbm = new GBM(seed: 42);
        var se = new Stderr(period);
        for (int i = 0; i < 500; i++)
        {
            var bar = gbm.Next();
            se.Update(new TValue(bar.Time, bar.Close));
            Assert.True(se.Last.Value >= 0.0, $"Stderr < 0 at bar {i}: {se.Last.Value}");
        }
    }

    [Fact]
    public void Stderr_IsNonNegative_GBM()
    {
        // SE is always non-negative by definition (sqrt of a variance-like quantity)
        const int period = 14;
        var gbm = new GBM(seed: 1);
        var series = new TSeries();
        for (int i = 0; i < 200; i++)
        {
            var bar = gbm.Next();
            series.Add(new TValue(bar.Time, bar.Close));
        }

        var seResult = Stderr.Batch(series, period);

        for (int i = 0; i < series.Count; i++)
        {
            Assert.True(seResult[i].Value >= 0.0,
                $"Stderr < 0 at bar {i}: {seResult[i].Value}");
        }
    }

    // Note: Tulip `stderr` is NOT the standard error of linear regression.
    // Tulip formula: stddev(x, n) / sqrt(n) = standard error of the mean.
    // QuanTAlib Stderr: sqrt(SSR / (n-2)) = standard error of OLS regression.
    // These are different statistics — no cross-validation is possible.
    // QuanTAlib is validated against its own brute-force OLS reference above.
}