diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj
index 4cd23169..262f8e2a 100644
--- a/src/Numerics/Numerics.csproj
+++ b/src/Numerics/Numerics.csproj
@@ -111,6 +111,7 @@
+
diff --git a/src/Numerics/Statistics/Percentile.cs b/src/Numerics/Statistics/Percentile.cs
index 674c03d4..65f3bb8a 100644
--- a/src/Numerics/Statistics/Percentile.cs
+++ b/src/Numerics/Statistics/Percentile.cs
@@ -71,7 +71,7 @@ namespace MathNet.Numerics.Statistics
///
/// Holds the data.
///
- private readonly List _data;
+ private readonly double[] _data;
///
/// Gets or sets the method used to calculate the percentiles.
@@ -94,9 +94,9 @@ namespace MathNet.Numerics.Statistics
{
throw new ArgumentNullException("data");
}
-
- _data = new List(data);
- _data.Sort();
+
+ _data = data.ToArray();
+ Array.Sort(_data);
}
///
@@ -106,39 +106,19 @@ namespace MathNet.Numerics.Statistics
/// the requested percentile.
public double Compute(double percentile)
{
- if (percentile < 0 || percentile > 1 || _data.Count == 0)
- {
- return double.NaN;
- }
-
- if (percentile == 0.0 || _data.Count == 1)
- {
- return _data[0];
- }
-
- if (percentile == 1.0)
- {
- return _data[_data.Count - 1];
- }
-
- var result = double.NaN;
switch (Method)
{
case PercentileMethod.Nist:
- result = Nist(percentile);
- break;
+ return SortedArrayStatistics.QuantileCompatible(_data, percentile, QuantileCompatibility.Nist);
case PercentileMethod.Nearest:
- result = Nearest(percentile);
- break;
+ return SortedArrayStatistics.QuantileCompatible(_data, percentile, QuantileCompatibility.R3);
case PercentileMethod.Interpolation:
- result = Interpolation(percentile);
- break;
+ return SortedArrayStatistics.QuantileCompatible(_data, percentile, QuantileCompatibility.R5);
case PercentileMethod.Excel:
- result = Excel(percentile);
- break;
+ return SortedArrayStatistics.QuantileCompatible(_data, percentile, QuantileCompatibility.Excel);
+ default:
+ return SortedArrayStatistics.Quantile(_data, percentile);
}
-
- return result;
}
///
@@ -155,64 +135,5 @@ namespace MathNet.Numerics.Statistics
return percentiles.Select(Compute).ToList();
}
-
- ///
- /// Computes the percentile using the nearest value.
- ///
- /// The percentile.
- /// the percentile using the nearest value.
- private double Nearest(double percentile)
- {
- var n = (int)Math.Round((_data.Count * percentile) + 0.5, 0);
- return _data[n - 1];
- }
-
- ///
- /// Computes the percentile using Excel's method.
- ///
- /// The percentile.
- /// the percentile using Excel's method.
- private double Excel(double percentile)
- {
- var tmp = 1 + (percentile * (_data.Count - 1.0));
- var k = (int)tmp;
- var d = tmp - k;
-
- return _data[k - 1] + (d * (_data[k] - _data[k - 1]));
- }
-
- ///
- /// Computes the percentile using interpolation.
- ///
- /// The percentile.
- /// the percentile using the interpolation.
- private double Interpolation(double percentile)
- {
- var k = (int)(_data.Count * percentile);
- var pk = (k - 0.5) / _data.Count;
- if(k == 0)
- return _data[0];
-
- return _data[k - 1] + (_data.Count * (percentile - pk) * (_data[k] - _data[k - 1]));
- }
-
- ///
- /// Computes the percentile using NIST's method.
- ///
- /// The percentile.
- /// the percentile using NIST's method.
- private double Nist(double percentile)
- {
- var tmp = percentile * (_data.Count + 1.0);
- var k = (int)tmp;
- if(k == 0)
- return _data[0];
- if(k == _data.Count)
- return _data[k - 1];
-
- var d = tmp - k;
-
- return _data[k - 1] + (d * (_data[k] - _data[k - 1]));
- }
}
}
diff --git a/src/Numerics/Statistics/SortedArrayStatistics.cs b/src/Numerics/Statistics/SortedArrayStatistics.cs
new file mode 100644
index 00000000..00641e57
--- /dev/null
+++ b/src/Numerics/Statistics/SortedArrayStatistics.cs
@@ -0,0 +1,137 @@
+//
+// Math.NET Numerics, part of the Math.NET Project
+// http://numerics.mathdotnet.com
+// http://github.com/mathnet/mathnet-numerics
+// http://mathnetnumerics.codeplex.com
+//
+// Copyright (c) 2009-2013 Math.NET
+//
+// Permission is hereby granted, free of charge, to any person
+// obtaining a copy of this software and associated documentation
+// files (the "Software"), to deal in the Software without
+// restriction, including without limitation the rights to use,
+// copy, modify, merge, publish, distribute, sublicense, and/or sell
+// copies of the Software, and to permit persons to whom the
+// Software is furnished to do so, subject to the following
+// conditions:
+//
+// The above copyright notice and this permission notice shall be
+// included in all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
+// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
+// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
+// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
+// OTHER DEALINGS IN THE SOFTWARE.
+//
+
+using System;
+
+namespace MathNet.Numerics.Statistics
+{
+ public enum QuantileCompatibility
+ {
+ Default=0,
+ Nist,Nearest,Excel,
+ R1,R2,R3,R4,R5,R6,R7,R8,R9,
+ SAS1,SAS2,SAS3,SAS4,SAS5
+ }
+
+ public static class SortedArrayStatistics
+ {
+ const double Third = 1d / 3d;
+ const double Half = 1d / 2d;
+
+ ///
+ /// R-8, SciPy-(1/3,1/3):
+ /// Linear interpolation of the approximate medians for order statistics.
+ /// When tau < (2/3) / (N + 1/3), use x1. When tau >= (N - 1/3) / (N + 1/3), use xN.
+ ///
+ public static double Quantile(double[] data, double tau)
+ {
+ if (tau < 0d || tau > 1d || data == null || data.Length == 0) return double.NaN;
+ if (tau == 0d || data.Length == 1) return data[0];
+ if (tau == 1d) return data[data.Length - 1];
+
+ double h = (data.Length + Third)*tau + Third;
+ var hf = (int) h;
+ return data[hf - 1] + (h - hf)*(data[hf] - data[hf - 1]);
+ }
+
+ public static double QuantileCompatible(double[] data, double tau, QuantileCompatibility compatibility)
+ {
+ if (tau < 0d || tau > 1d || data == null || data.Length == 0) return double.NaN;
+ if (tau == 0d || data.Length == 1) return data[0];
+ if (tau == 1d) return data[data.Length - 1];
+
+ switch (compatibility)
+ {
+ case QuantileCompatibility.R1:
+ case QuantileCompatibility.SAS3:
+ {
+ double h = data.Length*tau + Half;
+ return data[(int) Math.Ceiling(h - Half) - 1];
+ }
+ case QuantileCompatibility.R2:
+ case QuantileCompatibility.SAS5:
+ {
+ double h = data.Length * tau + Half;
+ return (data[(int) Math.Ceiling(h - Half) - 1] + data[(int) (h + Half) - 1])*Half;
+ }
+ case QuantileCompatibility.R3:
+ case QuantileCompatibility.SAS2:
+ case QuantileCompatibility.Nearest:
+ {
+ double h = data.Length*tau;
+ return data[(int) Math.Round(h) - 1];
+ }
+ case QuantileCompatibility.R4:
+ case QuantileCompatibility.SAS1:
+ {
+ double h = data.Length*tau;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ case QuantileCompatibility.R5:
+ {
+ double h = data.Length*tau + Half;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ case QuantileCompatibility.R6:
+ case QuantileCompatibility.SAS4:
+ case QuantileCompatibility.Nist:
+ {
+ double h = (data.Length + 1)*tau;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ case QuantileCompatibility.R7:
+ case QuantileCompatibility.Excel:
+ {
+ double h = (data.Length - 1)*tau + 1d;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ case QuantileCompatibility.R8:
+ case QuantileCompatibility.Default:
+ {
+ double h = (data.Length + Third) * tau + Third;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ case QuantileCompatibility.R9:
+ {
+ double h = (data.Length + 1d/4d) * tau + 3d/8d;
+ var hf = (int)h;
+ return data[hf - 1] + (h - hf) * (data[hf] - data[hf - 1]);
+ }
+ default:
+ throw new NotSupportedException();
+ }
+ }
+ }
+}
diff --git a/src/Portable/Portable.csproj b/src/Portable/Portable.csproj
index c8bd617c..bf01491b 100644
--- a/src/Portable/Portable.csproj
+++ b/src/Portable/Portable.csproj
@@ -1059,6 +1059,9 @@
Statistics\Percentile.cs
+
+ Statistics\SortedArrayStatistics.cs
+
Statistics\Statistics.cs