From 644ea702de0f5e76d86ee11aa8ee0606ef1309e6 Mon Sep 17 00:00:00 2001 From: Marcus Cuda Date: Sun, 3 May 2015 16:33:38 +0300 Subject: [PATCH] added support for infinity --- src/Numerics/Statistics/MovingStatistics.cs | 218 ++++++++++++------ .../StatisticsTests/MovingStatisticsTests.cs | 155 ++++++++++++- 2 files changed, 300 insertions(+), 73 deletions(-) diff --git a/src/Numerics/Statistics/MovingStatistics.cs b/src/Numerics/Statistics/MovingStatistics.cs index 71c0ee44..44cf1662 100644 --- a/src/Numerics/Statistics/MovingStatistics.cs +++ b/src/Numerics/Statistics/MovingStatistics.cs @@ -46,11 +46,11 @@ namespace MathNet.Numerics.Statistics long _totalCountOffset; int _lastIndex; int _lastNaNTimeToLive; + int _lastPosInfTimeToLive; + int _lastNegInfTimeToLive; double _m1; double _m2; - double _m3; - double _m4; double _max = double.NegativeInfinity; double _min = double.PositiveInfinity; @@ -64,7 +64,8 @@ namespace MathNet.Numerics.Statistics _oldValues = new double[_windowSize]; } - public MovingStatistics(int windowSize, IEnumerable values) : this(windowSize) + public MovingStatistics(int windowSize, IEnumerable values) + : this(windowSize) { PushRange(values); } @@ -88,7 +89,20 @@ namespace MathNet.Numerics.Statistics /// public double Minimum { - get { return _count > 0 && _lastNaNTimeToLive == 0 ? _min : double.NaN; } + get + { + if (_lastNaNTimeToLive > 0) + { + return double.NaN; + } + + if (_lastNegInfTimeToLive > 0) + { + return double.NegativeInfinity; + } + + return (_count > 0 || _lastPosInfTimeToLive > 0) ? _min : double.NaN; + } } /// @@ -97,7 +111,20 @@ namespace MathNet.Numerics.Statistics /// public double Maximum { - get { return _count > 0 && _lastNaNTimeToLive == 0 ? _max : double.NaN; } + get + { + if (_lastNaNTimeToLive > 0) + { + return double.NaN; + } + + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } + + return (_count > 0 || _lastNegInfTimeToLive > 0) ? _max : double.NaN; + } } /// @@ -106,7 +133,25 @@ namespace MathNet.Numerics.Statistics /// public double Mean { - get { return _count > 0 && _lastNaNTimeToLive == 0 ? _m1 : double.NaN; } + get + { + if (_lastNaNTimeToLive > 0 || (_lastPosInfTimeToLive > 0 && _lastNegInfTimeToLive > 0)) + { + return double.NaN; + } + + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } + + if (_lastNegInfTimeToLive > 0) + { + return double.NegativeInfinity; + } + + return _count == 0 ? double.NaN : _m1; + } } /// @@ -116,7 +161,20 @@ namespace MathNet.Numerics.Statistics /// public double Variance { - get { return _count < 2 || _lastNaNTimeToLive > 0 ? double.NaN : _m2/(_count - 1); } + get + { + if (_lastNaNTimeToLive > 0 || _lastNegInfTimeToLive > 0 || (_lastPosInfTimeToLive > 0 && _lastNegInfTimeToLive > 0)) + { + return double.NaN; + } + + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } + + return _count < 2 ? double.NaN : _m2 / (_count - 1); + } } /// @@ -126,7 +184,20 @@ namespace MathNet.Numerics.Statistics /// public double PopulationVariance { - get { return _count < 2 || _lastNaNTimeToLive > 0 ? double.NaN : _m2/_count; } + get + { + if (_lastNaNTimeToLive > 0 || _lastNegInfTimeToLive > 0 || (_lastPosInfTimeToLive > 0 && _lastNegInfTimeToLive > 0)) + { + return double.NaN; + } + + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } + + return _count < 2 ? double.NaN : _m2 / _count; + } } /// @@ -136,7 +207,20 @@ namespace MathNet.Numerics.Statistics /// public double StandardDeviation { - get { return _count < 2 || _lastNaNTimeToLive > 0 ? double.NaN : Math.Sqrt(_m2/(_count - 1)); } + get + { + if (_lastNaNTimeToLive > 0 || _lastNegInfTimeToLive > 0 || (_lastPosInfTimeToLive > 0 && _lastNegInfTimeToLive > 0)) + { + return double.NaN; + } + + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } + + return _count < 2 ? double.NaN : Math.Sqrt(_m2 / (_count - 1)); + } } /// @@ -146,72 +230,48 @@ namespace MathNet.Numerics.Statistics /// public double PopulationStandardDeviation { - get { return _count < 2 || _lastNaNTimeToLive > 0 ? double.NaN : Math.Sqrt(_m2/_count); } - } - -/* /// - /// Estimates the unbiased population skewness from the provided samples. - /// Uses a normalizer (Bessel's correction; type 2). - /// Returns NaN if data has less than three entries or if any entry is NaN. - /// - public double Skewness - { - get { return Count < 3 ? double.NaN : (Count*_m3*Math.Sqrt(_m2/(Count - 1))/(_m2*_m2*(Count - 2)))*(Count - 1); } - } + get + { + if (_lastNaNTimeToLive > 0 || _lastNegInfTimeToLive > 0 || (_lastPosInfTimeToLive > 0 && _lastNegInfTimeToLive > 0)) + { + return double.NaN; + } - /// - /// Evaluates the population skewness from the full population. - /// Does not use a normalizer and would thus be biased if applied to a subset (type 1). - /// Returns NaN if data has less than two entries or if any entry is NaN. - /// - public double PopulationSkewness - { - get { return Count < 2 ? double.NaN : Math.Sqrt(Count)*_m3/Math.Pow(_m2, 1.5); } - } + if (_lastPosInfTimeToLive > 0) + { + return double.PositiveInfinity; + } - /// - /// Estimates the unbiased population kurtosis from the provided samples. - /// Uses a normalizer (Bessel's correction; type 2). - /// Returns NaN if data has less than four entries or if any entry is NaN. - /// - public double Kurtosis - { - get { return Count < 4 ? double.NaN : ((double) Count*Count - 1)/((Count - 2)*(Count - 3))*(Count*_m4/(_m2*_m2) - 3 + 6.0/(Count + 1)); } + return _count < 2 ? double.NaN : Math.Sqrt(_m2 / _count); + } } - /// - /// Evaluates the population kurtosis from the full population. - /// Does not use a normalizer and would thus be biased if applied to a subset (type 1). - /// Returns NaN if data has less than three entries or if any entry is NaN. - /// - public double PopulationKurtosis - { - get { return Count < 3 ? double.NaN : Count*_m4/(_m2*_m2) - 3.0; } - }*/ - /// /// Update the running statistics by adding another observed sample (in-place). /// public void Push(double value) { + DecrementTimeToLive(); + if (double.IsNaN(value)) { - _totalCountOffset += _count + 1; - _count = 0; _lastNaNTimeToLive = _windowSize; + Reset(double.PositiveInfinity, double.NegativeInfinity); + return; + } - _m1 = 0.0; - _m2 = 0.0; - _m3 = 0.0; - _m4 = 0.0; - _max = double.NegativeInfinity; - _min = double.PositiveInfinity; + if (double.IsPositiveInfinity(value)) + { + _lastPosInfTimeToLive = _windowSize; + Reset(_min, double.NegativeInfinity); return; } - if (_lastNaNTimeToLive > 0) + if (double.IsNegativeInfinity(value)) { - _lastNaNTimeToLive--; + _lastNegInfTimeToLive = _windowSize; + Reset(double.PositiveInfinity, _max); + return; } if (_count < _windowSize) @@ -219,13 +279,10 @@ namespace MathNet.Numerics.Statistics _oldValues[_count] = value; _count++; var d = value - _m1; - var s = d/_count; -// var s2 = s * s; - var t = d*s*(_count - 1); + var s = d / _count; + var t = d * s * (_count - 1); _m1 += s; -// _m4 += t * s2 * (Count * Count - 3 * Count + 3) + 6 * s2 * _m2 - 4 * s * _m3; -// _m3 += t * s * (Count - 2) - 3 * s * _m2; _m2 += t; if (value < _min) @@ -242,18 +299,14 @@ namespace MathNet.Numerics.Statistics { var oldValue = _oldValues[_lastIndex]; var d = value - oldValue; - var s = d/_count; -// var s2 = s * s; + var s = d / _count; var oldM1 = _m1; _m1 += s; var x = (value - _m1 + oldValue - oldM1); - var t = d*x; + var t = d * x; _m2 += t; -// _m4 += t * s2 * (Count * Count - 3 * Count + 3) + 6 * s2 * _m2 - 4 * s * _m3; -// _m3 += t * (x /(Count-1) - 3 * s * _m2; - _oldValues[_lastIndex] = value; _lastIndex++; if (_lastIndex == WindowSize) @@ -275,5 +328,32 @@ namespace MathNet.Numerics.Statistics Push(value); } } + + private void DecrementTimeToLive() + { + if (_lastNaNTimeToLive > 0) + { + _lastNaNTimeToLive--; + } + + if (_lastPosInfTimeToLive > 0) + { + _lastPosInfTimeToLive--; + } + + if (_lastNegInfTimeToLive > 0) + { + _lastNegInfTimeToLive--; + } + } + + private void Reset(double min, double max) + { + _totalCountOffset += _count + 1; + _count = 0; + _m1 = 0; + _max = max; + _min = min; + } } } diff --git a/src/UnitTests/StatisticsTests/MovingStatisticsTests.cs b/src/UnitTests/StatisticsTests/MovingStatisticsTests.cs index 3753d706..867ab4a6 100644 --- a/src/UnitTests/StatisticsTests/MovingStatisticsTests.cs +++ b/src/UnitTests/StatisticsTests/MovingStatisticsTests.cs @@ -28,6 +28,7 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.Distributions; using MathNet.Numerics.Random; using MathNet.Numerics.Statistics; @@ -38,6 +39,153 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests [TestFixture, Category("Statistics")] public class MovingStatisticsTests { + [Test] + public void PositiveInfinityTest() + { + var ms = new MovingStatistics(3); + + ms.Push(1.0); + ms.Push(2.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.EqualTo(double.PositiveInfinity)); + + ms.Push(double.PositiveInfinity); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(1.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(double.PositiveInfinity); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(2.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(3.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(4.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.EqualTo(double.PositiveInfinity)); + } + + [Test] + public void NegativeInfinityTest() + { + var ms = new MovingStatistics(3); + + ms.Push(1.0); + ms.Push(2.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.EqualTo(double.NegativeInfinity)); + + ms.Push(double.NegativeInfinity); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(1.0); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(double.NegativeInfinity); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(2.0); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(3.0); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(4.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.NaN); + } + + [Test] + public void MixedInfinityTest() + { + var ms = new MovingStatistics(3); + + ms.Push(1.0); + ms.Push(2.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.EqualTo(double.PositiveInfinity)); + + ms.Push(double.NegativeInfinity); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(1.0); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(double.PositiveInfinity); + Assert.That(ms.Minimum, Is.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.NaN); + Assert.That(ms.StandardDeviation, Is.NaN); + + ms.Push(2.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(3.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.EqualTo(double.PositiveInfinity)); + + ms.Push(4.0); + Assert.That(ms.Minimum, Is.Not.EqualTo(double.NegativeInfinity)); + Assert.That(ms.Maximum, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.Mean, Is.Not.EqualTo(double.PositiveInfinity)); + Assert.That(ms.StandardDeviation, Is.Not.EqualTo(double.PositiveInfinity)); + } + [Test] public void StabilityTest() { @@ -53,10 +201,9 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests Assert.AreEqual(33.33, ms.Mean, 1e-11); Assert.AreEqual(308.58025, ms.Variance, 1e-10); - - //AssertHelpers.AlmostEqualRelative(stats0.Mean, ms.Mean, 14); - //AssertHelpers.AlmostEqualRelative(stats0.Variance, ms.Variance, 14); - //AssertHelpers.AlmostEqualRelative(stats0.StandardDeviation, ms.StandardDeviation, 14); + Assert.AreEqual(17.5664524022354, ms.StandardDeviation, 1e-11); + Assert.AreEqual(246.8642, ms.PopulationVariance, 1e-10); + Assert.AreEqual(15.7119126779651, ms.PopulationStandardDeviation, 1e-10); } [Test]