diff --git a/src/Numerics/Statistics/Correlation.cs b/src/Numerics/Statistics/Correlation.cs index 11566965..1faf845f 100644 --- a/src/Numerics/Statistics/Correlation.cs +++ b/src/Numerics/Statistics/Correlation.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// 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 @@ -28,19 +28,21 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; +using System.Collections.Generic; +using System.Linq; +using MathNet.Numerics.LinearAlgebra; +using MathNet.Numerics.Properties; + namespace MathNet.Numerics.Statistics { - using System; - using System.Collections.Generic; - using System.Linq; - /// /// A class with correlation measures between two datasets. /// public static class Correlation { /// - /// Computes the Pearson product-moment correlation coefficient. + /// Computes the Pearson Product-Moment Correlation coefficient. /// /// Sample data A. /// Sample data B. @@ -62,7 +64,7 @@ namespace MathNet.Numerics.Statistics { if (!ieB.MoveNext()) { - throw new ArgumentOutOfRangeException("dataB", "Datasets dataA and dataB need to have the same length. dataB is shorter."); + throw new ArgumentOutOfRangeException("dataB", Resources.ArgumentArraysSameLength); } double currentA = ieA.Current; double currentB = ieB.Current; @@ -82,7 +84,7 @@ namespace MathNet.Numerics.Statistics } if (ieB.MoveNext()) { - throw new ArgumentOutOfRangeException("dataA", "Datasets dataA and dataB need to have the same length. dataA is shorter."); + throw new ArgumentOutOfRangeException("dataA", Resources.ArgumentArraysSameLength); } } @@ -90,24 +92,78 @@ namespace MathNet.Numerics.Statistics } /// - /// Computes the Spearman Ranked Correlation Coefficient. + /// Computes the Pearson Product-Moment Correlation matrix. + /// + /// Array of sample data vectors. + /// The Pearson product-moment correlation matrix. + public static Matrix PearsonMatrix(params double[][] vectors) + { + var m = Matrix.Build.DenseIdentity(vectors.Length); + for (int i = 0; i < vectors.Length; i++) + for (int j = i + 1; j < vectors.Length; j++) + { + var c = Pearson(vectors[i], vectors[j]); + m.At(i, j, c); + m.At(j, i, c); + } + return m; + } + + /// + /// Computes the Pearson Product-Moment Correlation matrix. + /// + /// Enumerable of sample data vectors. + /// The Pearson product-moment correlation matrix. + public static Matrix PearsonMatrix(IEnumerable vectors) + { + return PearsonMatrix(vectors as double[][] ?? vectors.ToArray()); + } + + /// + /// Computes the Spearman Ranked Correlation coefficient. /// /// Sample data series A. /// Sample data series B. - /// The Spearman Ranked Correlation Coefficient. + /// The Spearman ranked correlation coefficient. public static double Spearman(IEnumerable dataA, IEnumerable dataB) { - return Pearson(RankedSeries(dataA.ToList()), RankedSeries(dataB.ToList())); + return Pearson(Rank(dataA), Rank(dataB)); + } + + /// + /// Computes the Spearman Ranked Correlation matrix. + /// + /// Array of sample data vectors. + /// The Spearman ranked correlation matrix. + public static Matrix SpearmanMatrix(params double[][] vectors) + { + return PearsonMatrix(vectors.Select(Rank).ToArray()); } - private static IEnumerable RankedSeries(ICollection series) + /// + /// Computes the Spearman Ranked Correlation matrix. + /// + /// Enumerable of sample data vectors. + /// The Spearman ranked correlation matrix. + public static Matrix SpearmanMatrix(IEnumerable vectors) { - if (series == null || series.Count == 0) - return Enumerable.Empty(); + return PearsonMatrix(vectors.Select(Rank).ToArray()); + } - var rankedSamples = series.Select((sample, index) => new { Sample = sample, RankIndex = index }).OrderBy(s => s.Sample).ToList(); + private static double[] Rank(IEnumerable series) + { + if (series == null) + { + return new double[0]; + } - var rankedArray = new double[series.Count]; + var rankedSamples = series.Select((sample, index) => new {Sample = sample, RankIndex = index}).OrderBy(s => s.Sample).ToArray(); + if (rankedSamples.Length == 0) + { + return new double[0]; + } + + var rankedArray = new double[rankedSamples.Length]; var previousSample = rankedSamples.Select((sampleIndex, index) => new { SampleIndex = sampleIndex, LoopIndex = index }).First(); foreach (var rankedSampleIndex in rankedSamples.Select((sampleIndex, index) => new { SampleIndex = sampleIndex, LoopIndex = index })) @@ -115,18 +171,24 @@ namespace MathNet.Numerics.Statistics var currentSample = rankedSampleIndex; if (Math.Abs(currentSample.SampleIndex.Sample - previousSample.SampleIndex.Sample) <= 0) + { continue; + } var rankedValue = (currentSample.LoopIndex + previousSample.LoopIndex - 1) / 2d + 1; foreach (var index in Enumerable.Range(previousSample.LoopIndex, currentSample.LoopIndex - previousSample.LoopIndex)) + { rankedArray[rankedSamples[index].RankIndex] = rankedValue; + } previousSample = currentSample; } - var finalValue = (rankedSamples.Count + previousSample.LoopIndex - 1) / 2d + 1; - foreach (var index in Enumerable.Range(previousSample.LoopIndex, rankedSamples.Count - previousSample.LoopIndex)) + var finalValue = (rankedSamples.Length + previousSample.LoopIndex - 1) / 2d + 1; + foreach (var index in Enumerable.Range(previousSample.LoopIndex, rankedSamples.Length - previousSample.LoopIndex)) + { rankedArray[rankedSamples[index].RankIndex] = finalValue; + } return rankedArray; }