diff --git a/src/Numerics/Statistics/ArrayStatistics.cs b/src/Numerics/Statistics/ArrayStatistics.cs index 991b9ce7..d2ec678f 100644 --- a/src/Numerics/Statistics/ArrayStatistics.cs +++ b/src/Numerics/Statistics/ArrayStatistics.cs @@ -279,16 +279,17 @@ namespace MathNet.Numerics.Statistics /// /// Estimates the median value from the unsorted data array. - /// Approximately median-unbiased regardless of the sample distribution (R8). /// WARNING: Works inplace and can thus causes the data array to be reordered. /// /// Sample array, no sorting is assumed. Will be reordered. public static double MedianInplace(double[] data) { - return QuantileInplace(data, 0.5d); + var k = data.Length/2; + return data.Length.IsOdd() + ? SelectInplace(data, k) + : (SelectInplace(data, k - 1) + SelectInplace(data, k))/2.0; } - /// /// Estimates the p-Percentile value from the unsorted data array. /// If a non-integer Percentile is needed, use Quantile instead. diff --git a/src/Numerics/Statistics/SortedArrayStatistics.cs b/src/Numerics/Statistics/SortedArrayStatistics.cs index f00e1244..07ae7b26 100644 --- a/src/Numerics/Statistics/SortedArrayStatistics.cs +++ b/src/Numerics/Statistics/SortedArrayStatistics.cs @@ -81,7 +81,12 @@ namespace MathNet.Numerics.Statistics /// Sample array, must be sorted ascendingly. public static double Median(double[] data) { - return Quantile(data, 0.5d); + if (data.Length == 0) return double.NaN; + + var k = data.Length/2; + return data.Length.IsOdd() + ? data[k] + : (data[k - 1] + data[k])/2.0; } /// diff --git a/src/UnitTests/StatisticsTests/StatisticsTests.cs b/src/UnitTests/StatisticsTests/StatisticsTests.cs index 8d7d4492..690afc8d 100644 --- a/src/UnitTests/StatisticsTests/StatisticsTests.cs +++ b/src/UnitTests/StatisticsTests/StatisticsTests.cs @@ -912,6 +912,20 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests Assert.AreEqual(21.578697, a.Variance(), 1e-5); Assert.AreEqual(21.578231, a.PopulationVariance(), 1e-5); } + + [Test] + public void MedianIsRobustOnCloseInfinities() + { + Assert.That(Statistics.Median(new[] { 2.0, double.NegativeInfinity, double.PositiveInfinity }), Is.EqualTo(2.0)); + Assert.That(Statistics.Median(new[] { 2.0, double.NegativeInfinity, 3.0, double.PositiveInfinity }), Is.EqualTo(2.5)); + Assert.That(ArrayStatistics.MedianInplace(new[] { 2.0, double.NegativeInfinity, double.PositiveInfinity}), Is.EqualTo(2.0)); + Assert.That(ArrayStatistics.MedianInplace(new[] { double.NegativeInfinity, 2.0, double.PositiveInfinity }), Is.EqualTo(2.0)); + Assert.That(ArrayStatistics.MedianInplace(new[] { double.NegativeInfinity, double.PositiveInfinity, 2.0 }), Is.EqualTo(2.0)); + Assert.That(ArrayStatistics.MedianInplace(new[] { double.NegativeInfinity, 2.0, 3.0, double.PositiveInfinity }), Is.EqualTo(2.5)); + Assert.That(ArrayStatistics.MedianInplace(new[] { double.NegativeInfinity, 2.0, double.PositiveInfinity, 3.0, }), Is.EqualTo(2.5)); + Assert.That(SortedArrayStatistics.Median(new[] { double.NegativeInfinity, 2.0, double.PositiveInfinity }), Is.EqualTo(2.0)); + Assert.That(SortedArrayStatistics.Median(new[] { double.NegativeInfinity, 2.0, 3.0, double.PositiveInfinity }), Is.EqualTo(2.5)); + } } }