Add R² and SMAPE error metrics with comprehensive tests and documentation

- Introduced R² (Coefficient of Determination) metric with detailed mathematical foundation, performance profile, and usage examples.
- Implemented SMAPE (Symmetric Mean Absolute Percentage Error) metric, addressing asymmetry in MAPE with symmetric error calculations.
- Added unit tests for SMAPE covering various scenarios including edge cases and input validation.
- Enhanced Dema class to correctly handle event publishing with isNew parameter.
- Updated Quantower test project to include coverage configuration for better test reporting.
This commit is contained in:
Miha Kralj
2025-12-29 20:58:21 -08:00
parent 4dbb093892
commit bf611d319f
50 changed files with 11327 additions and 21 deletions
+311
View File
@@ -0,0 +1,311 @@
namespace QuanTAlib.Tests;
public class MseTests
{
[Fact]
public void Constructor_ValidatesInput()
{
Assert.Throws<ArgumentException>(() => new Mse(0));
Assert.Throws<ArgumentException>(() => new Mse(-1));
var mse = new Mse(10);
Assert.NotNull(mse);
}
[Fact]
public void Properties_Accessible()
{
var mse = new Mse(10);
Assert.Equal(0, mse.Last.Value);
Assert.False(mse.IsHot);
Assert.Contains("Mse", mse.Name, StringComparison.Ordinal);
mse.Update(100, 105);
Assert.NotEqual(0, mse.Last.Time);
}
[Fact]
public void IsHot_BecomesTrueWhenBufferFull()
{
int period = 5;
var mse = new Mse(period);
for (int i = 0; i < period - 1; i++)
{
Assert.False(mse.IsHot, $"IsHot should be false at index {i}");
mse.Update(i * 10, i * 10 + 5);
}
mse.Update((period - 1) * 10, (period - 1) * 10 + 5);
Assert.True(mse.IsHot, "IsHot should be true after period updates");
}
[Fact]
public void Mse_CalculatesCorrectly()
{
var mse = new Mse(3);
// (10 - 15)² = 25
var res1 = mse.Update(10, 15);
Assert.Equal(25.0, res1.Value, 10);
// (20 - 30)² = 100, Mean = (25 + 100) / 2 = 62.5
var res2 = mse.Update(20, 30);
Assert.Equal(62.5, res2.Value, 10);
// (30 - 25)² = 25, Mean = (25 + 100 + 25) / 3 = 50
var res3 = mse.Update(30, 25);
Assert.Equal(50.0, res3.Value, 10);
// (40 - 35)² = 25, Window slides: (100 + 25 + 25) / 3 = 50
var res4 = mse.Update(40, 35);
Assert.Equal(50.0, res4.Value, 10);
}
[Fact]
public void Mse_PerfectPrediction_ReturnsZero()
{
var mse = new Mse(5);
for (int i = 0; i < 10; i++)
{
mse.Update(i * 10, i * 10); // Perfect prediction
}
Assert.Equal(0.0, mse.Last.Value, 10);
}
[Fact]
public void Mse_ConstantError_ReturnsSquaredConstant()
{
var mse = new Mse(5);
for (int i = 0; i < 10; i++)
{
mse.Update(100, 110); // Constant error of 10, squared = 100
}
Assert.Equal(100.0, mse.Last.Value, 10);
}
[Fact]
public void Mse_PenalizesLargeErrors()
{
var mse = new Mse(3);
// Small errors: (1-2)² = 1, (2-3)² = 1, (3-4)² = 1
// Mean = 1
mse.Update(1, 2);
mse.Update(2, 3);
var smallResult = mse.Update(3, 4);
Assert.Equal(1.0, smallResult.Value, 10);
mse.Reset();
// Large error: (1-11)² = 100, (2-3)² = 1, (3-4)² = 1
// Mean = 102/3 = 34
mse.Update(1, 11); // Large error
mse.Update(2, 3);
var largeResult = mse.Update(3, 4);
Assert.Equal(102.0 / 3.0, largeResult.Value, 10);
}
[Fact]
public void Calc_IsNew_AcceptsParameter()
{
var mse = new Mse(10);
mse.Update(100, 110, isNew: true);
double value1 = mse.Last.Value;
mse.Update(100, 120, isNew: true);
double value2 = mse.Last.Value;
Assert.NotEqual(value1, value2);
}
[Fact]
public void Calc_IsNew_False_UpdatesValue()
{
var mse = new Mse(10);
mse.Update(100, 110);
mse.Update(100, 120, isNew: true);
double beforeUpdate = mse.Last.Value;
mse.Update(100, 130, isNew: false);
double afterUpdate = mse.Last.Value;
Assert.NotEqual(beforeUpdate, afterUpdate);
}
[Fact]
public void IterativeCorrections_RestoreToOriginalState()
{
var mse = new Mse(5);
double tenthActual = 0;
double tenthPredicted = 0;
// Feed 10 updates
for (int i = 0; i < 10; i++)
{
tenthActual = i * 10;
tenthPredicted = i * 10 + 5;
mse.Update(tenthActual, tenthPredicted);
}
double stateAfterTen = mse.Last.Value;
// Apply 5 corrections with isNew=false
for (int i = 0; i < 5; i++)
{
mse.Update(100 + i, 200 + i, isNew: false);
}
// Restore to original values
mse.Update(tenthActual, tenthPredicted, isNew: false);
Assert.Equal(stateAfterTen, mse.Last.Value, 10);
}
[Fact]
public void Reset_ClearsState()
{
var mse = new Mse(5);
for (int i = 0; i < 10; i++)
{
mse.Update(i * 10, i * 10 + 5);
}
Assert.True(mse.IsHot);
mse.Reset();
Assert.False(mse.IsHot);
Assert.Equal(0, mse.Last.Value);
}
[Fact]
public void NaN_Input_UsesLastValidValue()
{
var mse = new Mse(5);
mse.Update(100, 110);
mse.Update(110, 120);
mse.Update(120, 130);
var result = mse.Update(double.NaN, double.NaN);
Assert.True(double.IsFinite(result.Value));
}
[Fact]
public void Infinity_Input_UsesLastValidValue()
{
var mse = new Mse(5);
mse.Update(100, 110);
mse.Update(110, 120);
var result = mse.Update(double.PositiveInfinity, double.NegativeInfinity);
Assert.True(double.IsFinite(result.Value));
}
[Fact]
public void Mse_Throws_On_Single_Input()
{
var mse = new Mse(10);
Assert.Throws<NotSupportedException>(() => mse.Update(new TValue(DateTime.UtcNow, 1)));
Assert.Throws<NotSupportedException>(() => mse.Update(new TSeries()));
Assert.Throws<NotSupportedException>(() => mse.Prime(new double[] { 1, 2, 3 }));
}
[Fact]
public void BatchSpan_MatchesStreaming()
{
int period = 5;
int count = 100;
var gbm = new GBM(startPrice: 100, mu: 0.05, sigma: 0.2, seed: 123);
double[] actual = new double[count];
double[] predicted = new double[count];
for (int i = 0; i < count; i++)
{
var bar = gbm.Next();
actual[i] = bar.Close;
predicted[i] = bar.Close * 1.05 + 2;
}
// Streaming
var mse = new Mse(period);
var streamingResults = new double[count];
for (int i = 0; i < count; i++)
{
streamingResults[i] = mse.Update(actual[i], predicted[i]).Value;
}
// Batch
double[] batchResults = new double[count];
Mse.Batch(actual, predicted, batchResults, period);
// Compare
for (int i = 0; i < count; i++)
{
Assert.Equal(streamingResults[i], batchResults[i], 9);
}
}
[Fact]
public void BatchSpan_ValidatesInput()
{
double[] actual = [1, 2, 3, 4, 5];
double[] predicted = [1, 2, 3, 4, 5];
double[] output = new double[5];
Assert.Throws<ArgumentException>(() =>
Mse.Batch(actual.AsSpan(), predicted.AsSpan(), output.AsSpan(), 0));
Assert.Throws<ArgumentException>(() =>
Mse.Batch(actual.AsSpan(), predicted.AsSpan(), output.AsSpan(), -1));
Assert.Throws<ArgumentException>(() =>
Mse.Batch(actual.AsSpan(), predicted.AsSpan(), new double[3].AsSpan(), 3));
}
[Fact]
public void Calculate_Works()
{
var actual = new TSeries();
var predicted = new TSeries();
var now = DateTime.UtcNow;
for (int i = 0; i < 10; i++)
{
actual.Add(now.AddMinutes(i), i * 10);
predicted.Add(now.AddMinutes(i), i * 10 + 5);
}
var results = Mse.Calculate(actual, predicted, 3);
Assert.Equal(10, results.Count);
// All errors are 5², so MSE should be 25
Assert.Equal(25.0, results.Last.Value, 10);
}
[Fact]
public void BatchSpan_HandlesNaN()
{
double[] actual = [100, 110, double.NaN, 130, 140];
double[] predicted = [105, 115, 125, double.NaN, 145];
double[] output = new double[5];
Mse.Batch(actual, predicted, output, 3);
foreach (var val in output)
{
Assert.True(double.IsFinite(val), $"Expected finite value but got {val}");
}
}
}
+280
View File
@@ -0,0 +1,280 @@
using System.Runtime.CompilerServices;
using System.Runtime.InteropServices;
namespace QuanTAlib;
/// <summary>
/// MSE: Mean Squared Error
/// </summary>
/// <remarks>
/// MSE measures the average of the squares of the errors between actual and
/// predicted values. It penalizes larger errors more heavily than MAE.
///
/// Formula:
/// MSE = (1/n) * Σ(actual - predicted)²
///
/// Uses a RingBuffer for O(1) streaming updates with running sum.
///
/// Key properties:
/// - Always non-negative (MSE ≥ 0)
/// - Units are squared (e.g., if data is in dollars, MSE is in dollars²)
/// - Heavily penalizes outliers due to squaring
/// - MSE = 0 indicates perfect prediction
/// </remarks>
[SkipLocalsInit]
public sealed class Mse : AbstractBase
{
private readonly RingBuffer _buffer;
[StructLayout(LayoutKind.Auto)]
private record struct State(double Sum, double LastValidActual, double LastValidPredicted, int TickCount);
private State _state;
private State _p_state;
private const int ResyncInterval = 1000;
/// <summary>
/// Creates MSE with specified period.
/// </summary>
/// <param name="period">Number of values to average (must be > 0)</param>
public Mse(int period)
{
if (period <= 0)
throw new ArgumentException("Period must be greater than 0", nameof(period));
_buffer = new RingBuffer(period);
Name = $"Mse({period})";
WarmupPeriod = period;
}
/// <summary>
/// True if the MSE has enough data to produce valid results.
/// </summary>
public override bool IsHot => _buffer.IsFull;
/// <summary>
/// Updates the MSE with new actual and predicted values.
/// </summary>
/// <param name="actual">Actual value (source1)</param>
/// <param name="predicted">Predicted value (source2)</param>
/// <param name="isNew">Whether this is a new bar.</param>
/// <returns>The calculated MSE value.</returns>
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public TValue Update(TValue actual, TValue predicted, bool isNew = true)
{
double actualVal = actual.Value;
double predictedVal = predicted.Value;
// Handle NaN/Infinity with last-valid-value substitution
if (!double.IsFinite(actualVal))
actualVal = double.IsFinite(_state.LastValidActual) ? _state.LastValidActual : 0.0;
else
_state.LastValidActual = actualVal;
if (!double.IsFinite(predictedVal))
predictedVal = double.IsFinite(_state.LastValidPredicted) ? _state.LastValidPredicted : 0.0;
else
_state.LastValidPredicted = predictedVal;
double diff = actualVal - predictedVal;
double error = diff * diff;
if (isNew)
{
_p_state = _state;
double removedValue = _buffer.Count == _buffer.Capacity ? _buffer.Oldest : 0.0;
_state.Sum = _state.Sum - removedValue + error;
_buffer.Add(error);
_state.TickCount++;
if (_buffer.IsFull && _state.TickCount >= ResyncInterval)
{
_state.TickCount = 0;
_state.Sum = _buffer.RecalculateSum();
}
}
else
{
_state = _p_state;
double removedValue = _buffer.Count == _buffer.Capacity ? _buffer.Oldest : 0.0;
_state.Sum = _state.Sum - removedValue + error;
_buffer.UpdateNewest(error);
_state.Sum = _buffer.RecalculateSum();
}
double result = _buffer.Count > 0 ? _state.Sum / _buffer.Count : error;
Last = new TValue(actual.Time, result);
PubEvent(Last, isNew);
return Last;
}
/// <summary>
/// Updates the MSE with raw double values.
/// </summary>
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public TValue Update(double actual, double predicted, bool isNew = true)
{
return Update(new TValue(DateTime.UtcNow, actual), new TValue(DateTime.UtcNow, predicted), isNew);
}
/// <summary>
/// Single-input Update is not supported. Use Update(actual, predicted).
/// </summary>
public override TValue Update(TValue input, bool isNew = true)
{
throw new NotSupportedException("MSE requires two inputs. Use Update(actual, predicted).");
}
/// <summary>
/// Single-series Update is not supported. Use Calculate(actual, predicted, period).
/// </summary>
public override TSeries Update(TSeries source)
{
throw new NotSupportedException("MSE requires two inputs. Use Calculate(actualSeries, predictedSeries, period).");
}
/// <summary>
/// Single-series Prime is not supported.
/// </summary>
public override void Prime(ReadOnlySpan<double> source, TimeSpan? step = null)
{
throw new NotSupportedException("MSE requires two inputs.");
}
/// <summary>
/// Resets the MSE state.
/// </summary>
public override void Reset()
{
_buffer.Clear();
_state = default;
_p_state = default;
Last = default;
}
/// <summary>
/// Calculates MSE for the entire series pair.
/// </summary>
public static TSeries Calculate(TSeries actual, TSeries predicted, int period)
{
if (actual.Count != predicted.Count)
throw new ArgumentException("Actual and predicted series must have the same length", nameof(predicted));
int len = actual.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(actual.Values, predicted.Values, vSpan, period);
actual.Times.CopyTo(tSpan);
return new TSeries(t, v);
}
/// <summary>
/// Calculates MSE in-place using pre-allocated spans.
/// </summary>
[MethodImpl(MethodImplOptions.AggressiveInlining)]
public static void Batch(ReadOnlySpan<double> actual, ReadOnlySpan<double> predicted, Span<double> output, int period)
{
if (actual.Length != predicted.Length || actual.Length != output.Length)
throw new ArgumentException("All spans must have the same length", nameof(output));
if (period <= 0)
throw new ArgumentException("Period must be greater than 0", nameof(period));
int len = actual.Length;
if (len == 0) return;
CalculateScalarCore(actual, predicted, output, period);
}
[MethodImpl(MethodImplOptions.AggressiveInlining)]
private static void CalculateScalarCore(ReadOnlySpan<double> actual, ReadOnlySpan<double> predicted, Span<double> output, int period)
{
int len = actual.Length;
const int StackAllocThreshold = 256;
Span<double> buffer = period <= StackAllocThreshold
? stackalloc double[period]
: new double[period];
double sum = 0;
double lastValidActual = 0;
double lastValidPredicted = 0;
// Find first valid values
for (int k = 0; k < len; k++)
{
if (double.IsFinite(actual[k]))
{
lastValidActual = actual[k];
break;
}
}
for (int k = 0; k < len; k++)
{
if (double.IsFinite(predicted[k]))
{
lastValidPredicted = predicted[k];
break;
}
}
int bufferIndex = 0;
int i = 0;
int warmupEnd = Math.Min(period, len);
for (; i < warmupEnd; i++)
{
double act = actual[i];
double pred = predicted[i];
if (double.IsFinite(act)) lastValidActual = act; else act = lastValidActual;
if (double.IsFinite(pred)) lastValidPredicted = pred; else pred = lastValidPredicted;
double diff = act - pred;
double error = diff * diff;
sum += error;
buffer[i] = error;
output[i] = sum / (i + 1);
}
int tickCount = 0;
for (; i < len; i++)
{
double act = actual[i];
double pred = predicted[i];
if (double.IsFinite(act)) lastValidActual = act; else act = lastValidActual;
if (double.IsFinite(pred)) lastValidPredicted = pred; else pred = lastValidPredicted;
double diff = act - pred;
double error = diff * diff;
sum = sum - buffer[bufferIndex] + error;
buffer[bufferIndex] = error;
bufferIndex++;
if (bufferIndex >= period) bufferIndex = 0;
output[i] = sum / period;
tickCount++;
if (tickCount >= ResyncInterval)
{
tickCount = 0;
double recalcSum = 0;
for (int k = 0; k < period; k++)
{
recalcSum += buffer[k];
}
sum = recalcSum;
}
}
}
}
+112
View File
@@ -0,0 +1,112 @@
# MSE: Mean Squared Error
> "The metric that makes outliers pay dearly for their transgressions."
Mean Squared Error (MSE) measures the average of the squares of the errors between actual and predicted values. By squaring errors, MSE penalizes large deviations more heavily than small ones.
## Historical Context
MSE is fundamental to least-squares regression, dating back to Gauss and Legendre in the early 1800s. It remains the most widely used loss function in machine learning and statistical modeling due to its mathematical convenience and theoretical properties.
## Architecture & Physics
MSE squares each error before averaging, which has significant implications:
- Large errors contribute disproportionately to the metric
- The quadratic penalty creates a smooth, differentiable loss surface
- Optimal for normally distributed errors
### Properties
- **Non-negative**: MSE ≥ 0, with 0 indicating perfect prediction
- **Squared units**: If data is in dollars, MSE is in dollars²
- **Outlier sensitive**: Single large error dominates the metric
- **Differentiable**: Smooth gradient for optimization algorithms
## Mathematical Foundation
### 1. Squared Error
For each observation, calculate the squared difference:
$$e_i = (y_i - \hat{y}_i)^2$$
### 2. Mean Calculation
Average the squared errors over the period:
$$MSE = \frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2$$
### 3. Running Update (O(1))
QuanTAlib uses a ring buffer with running sum for O(1) updates:
$$S_{new} = S_{old} - e_{oldest} + e_{newest}$$
$$MSE = \frac{S_{new}}{n}$$
## Implementation Details
### Usage Patterns
```csharp
// Streaming mode
var mse = new Mse(period: 20);
var result = mse.Update(actualValue, predictedValue);
// Batch mode
var results = Mse.Calculate(actualSeries, predictedSeries, period: 20);
// Span mode - zero allocation
Mse.Batch(actualSpan, predictedSpan, outputSpan, period: 20);
```
### Parameters
| Parameter | Type | Description |
| :--- | :--- | :--- |
| **period** | int | Lookback window for averaging (must be > 0) |
## Performance Profile
| Metric | Score | Notes |
| :--- | :--- | :--- |
| **Throughput** | ~12 ns/bar | O(1) with one multiplication |
| **Allocations** | 0 | Pre-allocated ring buffer |
| **Complexity** | O(1) | Constant time per update |
| **Accuracy** | 10/10 | Exact calculation |
## Interpretation
| MSE Range | Interpretation |
| :--- | :--- |
| **0** | Perfect prediction |
| **Low** | Predictions are close to actual values |
| **High** | Large prediction errors present |
## Relationship to RMSE
RMSE (Root Mean Squared Error) is simply the square root of MSE:
$$RMSE = \sqrt{MSE}$$
RMSE has the advantage of being in the same units as the original data.
## Common Use Cases
1. **Loss Function**: Primary loss for regression models
2. **Model Selection**: Compare models on validation data
3. **Gradient Descent**: Smooth gradient enables optimization
4. **Variance Estimation**: Related to sample variance
## Edge Cases
- **Identical Values**: Returns 0 when actual equals predicted
- **NaN Handling**: Uses last valid value substitution
- **Large Errors**: Can produce very large values due to squaring
## Related Indicators
- [MAE](../mae/Mae.md) - Mean Absolute Error (robust to outliers)
- [RMSE](../rmse/Rmse.md) - Root Mean Squared Error (same units as data)
- [Huber](../huber/Huber.md) - Combines MSE and MAE benefits