From 9a86ee8152213db78cb51ef9a54d83b34501b81c Mon Sep 17 00:00:00 2001 From: "joseph.lizier" Date: Wed, 1 Aug 2012 07:26:18 +0000 Subject: [PATCH] Added classes for joint entropy and multivariate mutual information calculation with gaussian assumption. Restructured MeasurementDistribution to have child classes EmpiricalMeasurementDistribution and AnalyticMeasurementDistribution - this resulted in many classes being altered. Added ChiSquare distribution methods to MathsUtils. Added more covariance methods to MatrixUtils. --- .../ActiveInfoStorageCalculator.java | 6 +- .../continuous/ChannelCalculatorCommon.java | 6 +- ...nfoCalculatorMultiVariateWithDiscrete.java | 6 +- ...culatorMultiVariateWithDiscreteSource.java | 6 +- ...nfoCalculatorMultiVariateWithDiscrete.java | 6 +- .../ActiveInfoStorageCalculatorKernel.java | 6 +- ...utualInfoCalculatorMultiVariateKernel.java | 8 +- ...culatorMultiVariateWithDiscreteKernel.java | 8 +- .../TransferEntropyCalculatorKernel.java | 8 +- .../TransferEntropyCalculatorKernelPlain.java | 6 +- ...EntropyCalculatorKernelPlainIterators.java | 6 +- ...ansferEntropyCalculatorKernelSeparate.java | 8 +- ...erEntropyCalculatorMultiVariateKernel.java | 8 +- ...tualInfoCalculatorMultiVariateKraskov.java | 8 +- ...ulatorMultiVariateWithDiscreteKraskov.java | 8 +- .../kraskov/MultiInfoCalculatorKraskov.java | 8 +- ...tualInfoCalculatorMultiVariateKraskov.java | 8 +- ...oCalculatorMultiVariateKraskovByMulti.java | 8 +- ...ulatorMultiVariateWithDiscreteKraskov.java | 8 +- ...ansferEntropyCalculatorKraskovByMulti.java | 12 +- ...yCalculatorMultiVariateLinearGaussian.java | 177 ++++++++++ ...oCalculatorMultiVariateLinearGaussian.java | 307 ++++++++++++++++++ ...latorMultiVariateWithDiscreteSymbolic.java | 6 +- ...latorMultiVariateWithDiscreteSymbolic.java | 6 +- .../TransferEntropyCalculatorSymbolic.java | 6 +- .../ApparentTransferEntropyCalculator.java | 6 +- .../measures/discrete/ChannelCalculator.java | 4 +- .../CompleteTransferEntropyCalculator.java | 6 +- .../discrete/MutualInformationCalculator.java | 8 +- .../InterregionalChannelMeasure.java | 6 +- ...ntDistributionPermutationsOverSubsets.java | 4 +- .../EmpiricalMeasurementDistribution.java | 88 +++++ .../source/infodynamics/utils/MathsUtils.java | 117 +++++++ .../infodynamics/utils/MatrixUtils.java | 106 +++++- .../utils/MeasurementDistribution.java | 129 +++----- 35 files changed, 938 insertions(+), 186 deletions(-) create mode 100755 java/source/infodynamics/measures/continuous/lineargaussian/EntropyCalculatorMultiVariateLinearGaussian.java create mode 100755 java/source/infodynamics/measures/continuous/lineargaussian/MutualInfoCalculatorMultiVariateLinearGaussian.java create mode 100755 java/source/infodynamics/utils/EmpiricalMeasurementDistribution.java diff --git a/java/source/infodynamics/measures/continuous/ActiveInfoStorageCalculator.java b/java/source/infodynamics/measures/continuous/ActiveInfoStorageCalculator.java index 611cc14..9f492aa 100755 --- a/java/source/infodynamics/measures/continuous/ActiveInfoStorageCalculator.java +++ b/java/source/infodynamics/measures/continuous/ActiveInfoStorageCalculator.java @@ -1,6 +1,6 @@ package infodynamics.measures.continuous; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; public interface ActiveInfoStorageCalculator { @@ -77,9 +77,9 @@ public interface ActiveInfoStorageCalculator { public double[] computeLocalOfPreviousObservations() throws Exception; - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception; public void setDebug(boolean debug); diff --git a/java/source/infodynamics/measures/continuous/ChannelCalculatorCommon.java b/java/source/infodynamics/measures/continuous/ChannelCalculatorCommon.java index 3843981..48de4b1 100755 --- a/java/source/infodynamics/measures/continuous/ChannelCalculatorCommon.java +++ b/java/source/infodynamics/measures/continuous/ChannelCalculatorCommon.java @@ -1,6 +1,6 @@ package infodynamics.measures.continuous; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * An abstract interface for calculators computing measures from a source to a destination. @@ -51,7 +51,7 @@ public abstract interface ChannelCalculatorCommon { * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; /** *

As per {@link computeSignificance(int) computeSignificance()} but supplies @@ -64,7 +64,7 @@ public abstract interface ChannelCalculatorCommon { * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception; public void setDebug(boolean debug); diff --git a/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscrete.java b/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscrete.java index 70f2bfa..febe6fe 100755 --- a/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscrete.java +++ b/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscrete.java @@ -1,6 +1,6 @@ package infodynamics.measures.continuous; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * A conditional mutual information calculator between a joint set of continuous variables, @@ -31,9 +31,9 @@ public interface ConditionalMutualInfoCalculatorMultiVariateWithDiscrete { public double[] computeLocalUsingPreviousObservations(double[][] contStates, int[] discreteStates, int[] conditionedStates) throws Exception; - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; public void setDebug(boolean debug); diff --git a/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSource.java b/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSource.java index 8aff536..fe415d3 100755 --- a/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSource.java +++ b/java/source/infodynamics/measures/continuous/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSource.java @@ -1,6 +1,6 @@ package infodynamics.measures.continuous; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * A conditional mutual information calculator between a joint set of continuous variables, @@ -32,9 +32,9 @@ public interface ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSource { public double[] computeLocalUsingPreviousObservations(double[][] contStates, int[] discreteStates, double[][] conditionedStates) throws Exception; - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; public void setDebug(boolean debug); diff --git a/java/source/infodynamics/measures/continuous/MutualInfoCalculatorMultiVariateWithDiscrete.java b/java/source/infodynamics/measures/continuous/MutualInfoCalculatorMultiVariateWithDiscrete.java index cfa2a4a..b372713 100755 --- a/java/source/infodynamics/measures/continuous/MutualInfoCalculatorMultiVariateWithDiscrete.java +++ b/java/source/infodynamics/measures/continuous/MutualInfoCalculatorMultiVariateWithDiscrete.java @@ -1,6 +1,6 @@ package infodynamics.measures.continuous; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; public interface MutualInfoCalculatorMultiVariateWithDiscrete { @@ -15,9 +15,9 @@ public interface MutualInfoCalculatorMultiVariateWithDiscrete { public double[] computeLocalUsingPreviousObservations(double[][] contStates, int[] discreteStates) throws Exception; - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception; - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception; public void setDebug(boolean debug); diff --git a/java/source/infodynamics/measures/continuous/kernel/ActiveInfoStorageCalculatorKernel.java b/java/source/infodynamics/measures/continuous/kernel/ActiveInfoStorageCalculatorKernel.java index 04a12a3..a746b36 100755 --- a/java/source/infodynamics/measures/continuous/kernel/ActiveInfoStorageCalculatorKernel.java +++ b/java/source/infodynamics/measures/continuous/kernel/ActiveInfoStorageCalculatorKernel.java @@ -2,7 +2,7 @@ package infodynamics.measures.continuous.kernel; import infodynamics.measures.continuous.ActiveInfoStorageCalculator; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * @@ -207,7 +207,7 @@ public class ActiveInfoStorageCalculatorKernel * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { return miKernel.computeSignificance(numPermutationsToCheck); } @@ -222,7 +222,7 @@ public class ActiveInfoStorageCalculatorKernel * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { return miKernel.computeSignificance(newOrderings); diff --git a/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateKernel.java b/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateKernel.java index 9117bae..016de01 100755 --- a/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateKernel.java +++ b/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateKernel.java @@ -2,7 +2,7 @@ package infodynamics.measures.continuous.kernel; import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; public class MutualInfoCalculatorMultiVariateKernel implements @@ -187,7 +187,7 @@ public class MutualInfoCalculatorMultiVariateKernel implements * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(observations1.length, numPermutationsToCheck); @@ -206,7 +206,7 @@ public class MutualInfoCalculatorMultiVariateKernel implements * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { computeAverageLocalOfObservations(); @@ -217,7 +217,7 @@ public class MutualInfoCalculatorMultiVariateKernel implements double[][] originalData2 = observations2; double[][] data2; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateWithDiscreteKernel.java b/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateWithDiscreteKernel.java index 5e487c3..bbbca01 100755 --- a/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateWithDiscreteKernel.java +++ b/java/source/infodynamics/measures/continuous/kernel/MutualInfoCalculatorMultiVariateWithDiscreteKernel.java @@ -4,7 +4,7 @@ import java.util.Arrays; import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariateWithDiscrete; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; public class MutualInfoCalculatorMultiVariateWithDiscreteKernel implements @@ -219,7 +219,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKernel implements * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(contObservations.length, numPermutationsToCheck); @@ -238,7 +238,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKernel implements * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { computeAverageLocalOfObservations(); @@ -247,7 +247,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKernel implements double actualMI = lastAverage; int[] originalDiscrete = discObservations; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernel.java b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernel.java index 4cc00de..555d02a 100755 --- a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernel.java +++ b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernel.java @@ -5,7 +5,7 @@ import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.measures.continuous.kernel.TransferEntropyKernelCounts; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Iterator; @@ -379,7 +379,7 @@ public class TransferEntropyCalculatorKernel * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); @@ -397,7 +397,7 @@ public class TransferEntropyCalculatorKernel * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; @@ -408,7 +408,7 @@ public class TransferEntropyCalculatorKernel double[] oldSourceValues = sourceValues; int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source in the destPastSourceVectors // and destNextPastSourceVectors vectors diff --git a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlain.java b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlain.java index 5c4dc21..9fc957b 100755 --- a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlain.java +++ b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlain.java @@ -4,7 +4,7 @@ import infodynamics.measures.continuous.TransferEntropyCalculator; import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import java.util.Iterator; @@ -420,12 +420,12 @@ public class TransferEntropyCalculatorKernelPlain return true; } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { throw new RuntimeException("Not implemented in this calculator"); } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { throw new RuntimeException("Not implemented in this calculator"); } diff --git a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlainIterators.java b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlainIterators.java index 24eb9ef..e6957aa 100755 --- a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlainIterators.java +++ b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelPlainIterators.java @@ -3,7 +3,7 @@ package infodynamics.measures.continuous.kernel; import infodynamics.measures.continuous.TransferEntropyCalculator; import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.utils.MathsUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import java.util.Iterator; @@ -426,12 +426,12 @@ public class TransferEntropyCalculatorKernelPlainIterators return true; } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { throw new RuntimeException("Not implemented in this calculator"); } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { throw new RuntimeException("Not implemented in this calculator"); } diff --git a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelSeparate.java b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelSeparate.java index 45d93bc..640549e 100755 --- a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelSeparate.java +++ b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorKernelSeparate.java @@ -5,7 +5,7 @@ import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.measures.continuous.kernel.KernelEstimatorMultiVariate; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Iterator; @@ -531,7 +531,7 @@ public class TransferEntropyCalculatorKernelSeparate * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); @@ -549,7 +549,7 @@ public class TransferEntropyCalculatorKernelSeparate * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; @@ -560,7 +560,7 @@ public class TransferEntropyCalculatorKernelSeparate double[] originalSourceValuesInJoint = MatrixUtils.selectColumn(destPastSourceVectors, k); int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Check that the length of the reorderings is OK if (newOrderings[p].length != totalObservations) { diff --git a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorMultiVariateKernel.java b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorMultiVariateKernel.java index a6bac7f..edabb4b 100755 --- a/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorMultiVariateKernel.java +++ b/java/source/infodynamics/measures/continuous/kernel/TransferEntropyCalculatorMultiVariateKernel.java @@ -5,7 +5,7 @@ import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.measures.continuous.kernel.TransferEntropyKernelCounts; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Iterator; @@ -578,7 +578,7 @@ public class TransferEntropyCalculatorMultiVariateKernel * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); @@ -596,7 +596,7 @@ public class TransferEntropyCalculatorMultiVariateKernel * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; @@ -607,7 +607,7 @@ public class TransferEntropyCalculatorMultiVariateKernel double[][] oldSourceValues = sourceVectors; int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source in the destPastSourceVectors // and destNextPastSourceVectors vectors diff --git a/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateKraskov.java b/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateKraskov.java index 27e7795..95b7131 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateKraskov.java +++ b/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateKraskov.java @@ -2,7 +2,7 @@ package infodynamics.measures.continuous.kraskov; import infodynamics.utils.EuclideanUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; /** @@ -174,7 +174,7 @@ public abstract class ConditionalMutualInfoCalculatorMultiVariateKraskov { * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(data1.length, numPermutationsToCheck); @@ -193,7 +193,7 @@ public abstract class ConditionalMutualInfoCalculatorMultiVariateKraskov { * @param newOrderings the specific new orderings to use * @return the proportion of conditional MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!condMiComputed) { computeAverageLocalOfObservations(); @@ -201,7 +201,7 @@ public abstract class ConditionalMutualInfoCalculatorMultiVariateKraskov { // Store the real observations and their MI: double actualMI = condMi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov.java b/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov.java index f2a66cf..e26c32f 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov.java +++ b/java/source/infodynamics/measures/continuous/kraskov/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov.java @@ -4,7 +4,7 @@ import infodynamics.measures.continuous.ConditionalMutualInfoCalculatorMultiVari import infodynamics.utils.EuclideanUtils; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; /** @@ -487,7 +487,7 @@ public class ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov impl * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(continuousDataX.length, numPermutationsToCheck); @@ -506,7 +506,7 @@ public class ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov impl * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { @@ -515,7 +515,7 @@ public class ConditionalMutualInfoCalculatorMultiVariateWithDiscreteKraskov impl // Store the real observations and their MI: double actualMI = condMi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/MultiInfoCalculatorKraskov.java b/java/source/infodynamics/measures/continuous/kraskov/MultiInfoCalculatorKraskov.java index 6bf406f..d1ec679 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/MultiInfoCalculatorKraskov.java +++ b/java/source/infodynamics/measures/continuous/kraskov/MultiInfoCalculatorKraskov.java @@ -2,7 +2,7 @@ package infodynamics.measures.continuous.kraskov; import infodynamics.measures.continuous.MultiInfoCalculator; import infodynamics.utils.EuclideanUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Vector; @@ -184,7 +184,7 @@ public abstract class MultiInfoCalculatorKraskov implements * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][][] newOrderings = new int[numPermutationsToCheck][][]; @@ -209,7 +209,7 @@ public abstract class MultiInfoCalculatorKraskov implements * third index is the reordered variable number for that position. * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { computeAverageLocalOfObservations(); @@ -217,7 +217,7 @@ public abstract class MultiInfoCalculatorKraskov implements // Store the real observations and their MI: double actualMI = mi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskov.java b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskov.java index 12f39f5..5a92fa2 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskov.java +++ b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskov.java @@ -3,7 +3,7 @@ package infodynamics.measures.continuous.kraskov; import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; import infodynamics.utils.EuclideanUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; /** @@ -175,7 +175,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskov implements * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(data1.length, numPermutationsToCheck); @@ -194,7 +194,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskov implements * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { computeAverageLocalOfObservations(); @@ -202,7 +202,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskov implements // Store the real observations and their MI: double actualMI = mi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskovByMulti.java b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskovByMulti.java index 7c6dc48..f885c40 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskovByMulti.java +++ b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateKraskovByMulti.java @@ -2,7 +2,7 @@ package infodynamics.measures.continuous.kraskov; import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Hashtable; @@ -262,7 +262,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskovByMulti implements * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(data1.length, numPermutationsToCheck); @@ -281,7 +281,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskovByMulti implements * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { computeAverageLocalOfObservations(); @@ -289,7 +289,7 @@ public abstract class MutualInfoCalculatorMultiVariateKraskovByMulti implements // Store the real observations and their MI: double actualMI = mi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateWithDiscreteKraskov.java b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateWithDiscreteKraskov.java index d61f191..7f5a298 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateWithDiscreteKraskov.java +++ b/java/source/infodynamics/measures/continuous/kraskov/MutualInfoCalculatorMultiVariateWithDiscreteKraskov.java @@ -4,7 +4,7 @@ import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariateWithDisc import infodynamics.utils.EuclideanUtils; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; /** @@ -361,7 +361,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKraskov implements Mutu * @param numPermutationsToCheck * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public synchronized MeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { + public synchronized EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(continuousData.length, numPermutationsToCheck); @@ -380,7 +380,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKraskov implements Mutu * @param newOrderings the specific new orderings to use * @return the proportion of MI scores from the distribution which have higher or equal MIs to ours. */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; if (!miComputed) { @@ -389,7 +389,7 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteKraskov implements Mutu // Store the real observations and their MI: double actualMI = mi; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); int countWhereMiIsMoreSignificantThanOriginal = 0; for (int i = 0; i < numPermutationsToCheck; i++) { diff --git a/java/source/infodynamics/measures/continuous/kraskov/TransferEntropyCalculatorKraskovByMulti.java b/java/source/infodynamics/measures/continuous/kraskov/TransferEntropyCalculatorKraskovByMulti.java index 5f73b68..5ab36c2 100755 --- a/java/source/infodynamics/measures/continuous/kraskov/TransferEntropyCalculatorKraskovByMulti.java +++ b/java/source/infodynamics/measures/continuous/kraskov/TransferEntropyCalculatorKraskovByMulti.java @@ -4,7 +4,7 @@ import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; import infodynamics.measures.continuous.TransferEntropyCalculator; import infodynamics.measures.continuous.TransferEntropyCommon; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Hashtable; @@ -351,7 +351,7 @@ public class TransferEntropyCalculatorKraskovByMulti * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { // Generate the re-ordered indices: RandomGenerator rg = new RandomGenerator(); @@ -373,7 +373,7 @@ public class TransferEntropyCalculatorKraskovByMulti * @return * @throws Exception */ - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; @@ -381,7 +381,7 @@ public class TransferEntropyCalculatorKraskovByMulti double actualTE = computeAverageLocalOfObservations(); int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source in the destPastSourceVectors // and destNextPastSourceVectors vectors @@ -423,7 +423,7 @@ public class TransferEntropyCalculatorKraskovByMulti * @return * @throws Exception */ - public MeasurementDistribution computeSignificanceExplicitlyReordering( + public EmpiricalMeasurementDistribution computeSignificanceExplicitlyReordering( int[][] newOrderings) throws Exception { int numPermutationsToCheck = newOrderings.length; @@ -437,7 +437,7 @@ public class TransferEntropyCalculatorKraskovByMulti } int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source in the destPastSourceVectors // and destNextPastSourceVectors vectors diff --git a/java/source/infodynamics/measures/continuous/lineargaussian/EntropyCalculatorMultiVariateLinearGaussian.java b/java/source/infodynamics/measures/continuous/lineargaussian/EntropyCalculatorMultiVariateLinearGaussian.java new file mode 100755 index 0000000..aad5566 --- /dev/null +++ b/java/source/infodynamics/measures/continuous/lineargaussian/EntropyCalculatorMultiVariateLinearGaussian.java @@ -0,0 +1,177 @@ +package infodynamics.measures.continuous.lineargaussian; + +import infodynamics.measures.continuous.EntropyCalculatorMultiVariate; +import infodynamics.utils.MatrixUtils; + +/** + *

Computes the differential entropy of a given multivariate set of observations, + * assuming that the probability distribution function for these observations is + * a multivariate Gaussian distribution.

+ * + *

+ * Usage: + *

    + *
  1. Construct
    2. + *
  2. initialise()
    3. + *
  3. setObservations(), or setCovariance().
    4. + *
  4. computeAverageLocalOfObservations() to return the average differential + * entropy based on either the set variance or the variance of + * the supplied observations, or computeLocalUsingPrevious
    5. + *
+ *

+ * + * @see Differential entropy for Gaussian random variables defined at + * {@link http://mathworld.wolfram.com/DifferentialEntropy.html} + * @author Joseph Lizier joseph.lizier_at_gmail.com + * + */ +public class EntropyCalculatorMultiVariateLinearGaussian implements EntropyCalculatorMultiVariate { + + /** + * Covariance matrix of the most recently supplied observations + */ + protected double[][] covariance; + + /** + * The set of observations, retained in case the user wants to retrieve the local + * entropy values of these + */ + protected double[][] observations; + + /** + * Number of dimenions for our multivariate data + */ + protected int dimensions; + + protected double lastAverage; + + protected boolean debug; + + /** + * Constructor + */ + public EntropyCalculatorMultiVariateLinearGaussian() { + // Nothing to do + } + + /** + * Initialise the calculator ready for reuse + */ + public void initialise(int dimensions) { + covariance = null; + observations = null; + this.dimensions = dimensions; + } + + /** + * Provide the multivariate observations from which to compute the entropy + * + * @param observations the observations to compute the entropy from. + * First index is time, second index is variable number. + */ + public void setObservations(double[][] observations) { + this.observations = observations; + covariance = MatrixUtils.covarianceMatrix(observations); + // Check that the observations was of the correct number of dimensions: + // (done afterwards since the covariance matrix computation checks that + // all rows had the right number of columns + if (covariance.length != dimensions) { + throw new RuntimeException("Supplied observations does not match initialised number of dimensions"); + } + } + + /** + * Set the covariance of the distribution for which we will compute the + * entropy. + * + * @param covariance covariance matrix + */ + public void setCovariance(double[][] covariance) throws Exception { + observations = null; + // Make sure the supplied covariance matrix is square: + int rows = covariance.length; + if (rows != dimensions) { + throw new Exception("Supplied covariance matrix does not match initialised number of dimensions"); + } + for (int r = 0; r < rows; r++) { + if (covariance[r].length != rows) { + throw new Exception("Cannot compute the determinant of a non-square matrix"); + } + } + + this.covariance = covariance; + } + + /** + *

The joint entropy for a multivariate Gaussian-distribution of dimension n + * with covariance matrix C is 0.5*\log_e{(2*pi*e)^n*|det(C)|}, + * where det() is the matrix determinant of C.

+ * + *

Here we compute the joint entropy assuming that the recorded estimation of the + * covariance is correct (i.e. we will not make a bias correction for limited + * observations here).

+ * + * @return the joint entropy of the previously provided observations or from the + * supplied covariance matrix. + */ + public double computeAverageLocalOfObservations() { + try { + lastAverage = 0.5 * (dimensions* (1 + Math.log(2.0*Math.PI)) + + Math.log(Math.abs(MatrixUtils.determinant(covariance)))); + return lastAverage; + } catch (Exception e) { + // Should not happen, since we check the validity of the supplied + // matrix beforehand; so we'll throw an Error in this case + throw new Error(e); + } + } + + public void setDebug(boolean debug) { + this.debug = debug; + } + + /** + *

Set the given property to the given value.

+ * + *

There are currently no properties to set for this calculator

+ * + * @param propertyName name of the property + * @param propertyValue value of the property. + * @throws Exception + */ + public void setProperty(String propertyName, String propertyValue) + throws Exception { + // No properties to set here + } + + /** + * @return the lastAverage + */ + public double getLastAverage() { + return lastAverage; + } + + public double[] computeLocalUsingPreviousObservations(double[][] states) + throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + + public double[] computeLocalOfPreviousObservations() { + + // TODO Implement this function + if (true) + throw new RuntimeException("Not implemented yet"); + + if (observations == null) { + throw new RuntimeException("Cannot compute local values since no observations were supplied"); + } + double[] localEntropy = new double[observations.length]; + for (int t=0; t < observations.length; t++) { + // Compute the probability for the given observation, based on + // the assumption of a multivariate Gaussian PDF: + + } + return null; + } +} diff --git a/java/source/infodynamics/measures/continuous/lineargaussian/MutualInfoCalculatorMultiVariateLinearGaussian.java b/java/source/infodynamics/measures/continuous/lineargaussian/MutualInfoCalculatorMultiVariateLinearGaussian.java new file mode 100755 index 0000000..5c4dd70 --- /dev/null +++ b/java/source/infodynamics/measures/continuous/lineargaussian/MutualInfoCalculatorMultiVariateLinearGaussian.java @@ -0,0 +1,307 @@ +package infodynamics.measures.continuous.lineargaussian; + +import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; +import infodynamics.utils.ChiSquareMeasurementDistribution; +import infodynamics.utils.MathsUtils; +import infodynamics.utils.MatrixUtils; +import infodynamics.utils.EmpiricalMeasurementDistribution; + +/** + *

Computes the differential mutual information of two given multivariate sets of + * observations, + * assuming that the probability distribution function for these observations is + * a multivariate Gaussian distribution.

+ * + *

+ * Usage: + *

    + *
  1. Construct {@link #MutualInfoCalculatorMultiVariateLinearGaussian()}
    2. + *
  2. {@link #initialise(int, int)}
    3. + *
  3. Provide the observations to the calculator using: + * {@link #setObservations(double[][], double[][])}, or + * {@link #setCovariance(double[][])}, or + * a sequence of: + * {@link #startAddObservations()}, + * multiple calls to either {@link #addObservations(double[][], double[][])} + * or @{link {@link #addObservations(double[][], double[][], int, int)}, and then + * @{link {@link #finaliseAddObservations()}.
    4. + *
  4. Compute the required information-theoretic results, primarily: + * @{link #computeAverageLocalOfObservations()} to return the average differential + * entropy based on either the set variance or the variance of + * the supplied observations; or other calls to compute + * local values or statistical significance.
    5. + *
+ *

+ * + * @see Differential entropy for Gaussian random variables defined at + * {@link http://mathworld.wolfram.com/DifferentialEntropy.html} + * @author Joseph Lizier joseph.lizier_at_gmail.com + * + */ +public class MutualInfoCalculatorMultiVariateLinearGaussian implements + MutualInfoCalculatorMultiVariate { + + /** + * Covariance matrix of the most recently supplied observations + */ + protected double[][] covariance; + + /** + * The set of source observations, retained in case the user wants to retrieve the local + * entropy values of these. + * They're held in the order in which they were supplied in the + * {@link addObservations(double[][], double[][])} functions. + */ + protected double[][] sourceObservations; + + /** + * The set of destination observations, retained in case the user wants to retrieve the local + * entropy values of these. + * They're held in the order in which they were supplied in the + * {@link addObservations(double[][], double[][])} functions. + */ + protected double[][] destObservations; + + /** + * Number of dimenions for each of our multivariate data sets + */ + protected int dimensionsDest; + protected int dimensionsSource; + + protected double lastAverage; + + protected boolean debug; + + public MutualInfoCalculatorMultiVariateLinearGaussian() { + // Nothing to do + } + + /** + * Clear any previously supplied probability distributions and prepare + * the calculator to be used again. + * + * @param sourceDimensions number of joing variables in the source + * @param destDimensions number of joing variables in the destination + */ + public void initialise(int sourceDimensions, int destDimensions) { + covariance = null; + sourceObservations = null; + destObservations = null; + dimensionsSource = sourceDimensions; + dimensionsDest = destDimensions; + } + + /** + * Provide the complete set of observations to use to compute the + * mutual information. + * One cannot use the {@link addObservations(double[][], double[][])} + * style methods after this without calling initialise again first. + * + */ + public void setObservations(double[][] source, double[][] destination) + throws Exception { + sourceObservations = source; + destObservations = destination; + covariance = MatrixUtils.covarianceMatrix(source, destination); + // Check that the observations was of the correct number of dimensions: + // (done afterwards since the covariance matrix computation checks that + // all rows had the right number of columns + if (covariance.length != dimensionsSource + dimensionsDest) { + throw new RuntimeException("Supplied observations do not match initialised number of dimensions"); + } + } + + public void addObservations(double[][] source, double[][] destination) + throws Exception { + // TODO implement these addObservations style functions. + // This will not be hard to implement - see the implementation + // for TE in TransferEntropyCommon. It might be useful + // to have a MutualInfoCommon which pulls the same functionality + // together for the MI calculators anyway. + throw new RuntimeException("Not implemented yet"); + } + + public void addObservations(double[][] source, double[][] destination, + int startTime, int numTimeSteps) throws Exception { + throw new RuntimeException("Not implemented yet"); + } + + public void setObservations(double[][] source, double[][] destination, + boolean[] sourceValid, boolean[] destValid) throws Exception { + throw new RuntimeException("Not implemented yet"); + } + + public void setObservations(double[][] source, double[][] destination, + boolean[][] sourceValid, boolean[][] destValid) throws Exception { + throw new RuntimeException("Not implemented yet"); + } + + public void startAddObservations() { + throw new RuntimeException("Not implemented yet"); + } + + public void finaliseAddObservations() { + throw new RuntimeException("Not implemented yet"); + } + + /** + * Set the covariance of the distribution for which we will compute the + * mutual information. + * + * @param covariance covariance matrix of the source and destination + * variables, considered together (variable indices start with the source + * and continue into the destination). + */ + public void setCovariance(double[][] covariance) throws Exception { + sourceObservations = null; + destObservations = null; + // Make sure the supplied covariance matrix is square: + int rows = covariance.length; + if (rows != dimensionsSource + dimensionsDest) { + throw new Exception("Supplied covariance matrix does not match initialised number of dimensions"); + } + for (int r = 0; r < rows; r++) { + if (covariance[r].length != rows) { + throw new Exception("Cannot compute the determinant of a non-square matrix"); + } + } + + this.covariance = covariance; + } + + /** + *

The joint entropy for a multivariate Gaussian-distribution of dimension n + * with covariance matrix C is 0.5*\log_e{(2*pi*e)^n*|det(C)|}, + * where det() is the matrix determinant of C.

+ * + *

Here we compute the mutual information from the joint entropies + * of the source variables (H_s), destination variables (H_d), and all variables + * taken together (H_sd), giving MI = H_s + H_d - H_sd. + * We assume that the recorded estimation of the + * covariance is correct (i.e. we will not make a bias correction for limited + * observations here).

+ * + * @return the mutual information of the previously provided observations or from the + * supplied covariance matrix. + */ + public double computeAverageLocalOfObservations() throws Exception { + try { + int[] sourceIndicesInCovariance = MatrixUtils.range(0, dimensionsSource - 1); + int[] destIndicesInCovariance = MatrixUtils.range(dimensionsSource, + dimensionsSource + dimensionsDest - 1); + double[][] sourceCovariance = + MatrixUtils.selectRowsAndColumns(covariance, + sourceIndicesInCovariance, sourceIndicesInCovariance); + double[][] destCovariance = + MatrixUtils.selectRowsAndColumns(covariance, + destIndicesInCovariance, destIndicesInCovariance); + double sourceEntropy = 0.5 * + Math.log(Math.abs(MatrixUtils.determinant(sourceCovariance))); + double destEntropy = 0.5 * + Math.log(Math.abs(MatrixUtils.determinant(destCovariance))); + double jointEntropy = 0.5 * + Math.log(Math.abs(MatrixUtils.determinant(covariance))); + lastAverage = sourceEntropy + destEntropy - jointEntropy; + return lastAverage; + } catch (Exception e) { + // Should not happen, since we check the validity of the supplied + // matrix beforehand; so we'll throw an Error in this case + throw new Error(e); + } + } + + public double[] computeLocalOfPreviousObservations() throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + + /** + *

Compute the statistical significance of the mutual information + * result analytically, without creating a distribution + * under the null hypothesis by bootstrapping.

+ * + *

Brillinger (see reference below) shows that under the null hypothesis + * of no source-destination relationship, the MI for two + * Gaussian distributions follows a chi-square distribution with + * degrees of freedom equal to the product of the number of variables + * in each joint variable.

+ * + * @return MeasurementDistribution object with only the + * {@link EmpiricalMeasurementDistribution#actualValue} and + * {@link EmpiricalMeasurementDistribution#pValue} fields filled out. + * This object contains the proportion of MI scores from the distribution + * which have higher or equal MIs to ours. + * + * @see Brillinger, "Some data analyses using mutual information", + * {@link http://www.stat.berkeley.edu/~brill/Papers/MIBJPS.pdf} + * @see Cheng et al., "Data Information in Contingency Tables: A + * Fallacy of Hierarchical Loglinear Models", + * {@link http://www.jds-online.com/file_download/112/JDS-369.pdf} + * @see Barnett and Bossomaier, "Transfer Entropy as a Log-likelihood Ratio" + * {@link http://arxiv.org/abs/1205.6339} + */ + public ChiSquareMeasurementDistribution computeSignificance() { + // TODO Check that the null distribution actually follows chi with + // these degrees of freedom + return new ChiSquareMeasurementDistribution(lastAverage, + dimensionsSource * dimensionsDest); + } + + public EmpiricalMeasurementDistribution computeSignificance( + int numPermutationsToCheck) throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) + throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + + /** + *

Set the given property to the given value.

+ * + *

There are currently no properties to set for this calculator

+ * + * @param propertyName name of the property + * @param propertyValue value of the property. + * @throws Exception + */ + public void setProperty(String propertyName, String propertyValue) + throws Exception { + // No properties to set here + } + + public void setDebug(boolean debug) { + this.debug = debug; + } + + /** + * @return the previously computed average mutual information + */ + public double getLastAverage() { + return lastAverage; + } + + /** + * @return the number of previously supplied observations for which + * the mutual information will be / was computed. + */ + public int getNumObservations() { + return destObservations.length; + } + + public double computeAverageLocalOfObservations(int[] newOrdering) + throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + + public double[] computeLocalUsingPreviousObservations(double[][] states1, + double[][] states2) throws Exception { + // TODO Implement me + throw new RuntimeException("Not implemented yet"); + } + +} diff --git a/java/source/infodynamics/measures/continuous/symbolic/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java b/java/source/infodynamics/measures/continuous/symbolic/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java index 40c3512..87d1b4f 100755 --- a/java/source/infodynamics/measures/continuous/symbolic/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java +++ b/java/source/infodynamics/measures/continuous/symbolic/ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java @@ -5,7 +5,7 @@ import infodynamics.measures.discrete.ConditionalMutualInformationCalculator; import infodynamics.utils.FirstIndexComparatorDouble; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; public class ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic implements @@ -131,12 +131,12 @@ public class ConditionalMutualInfoCalculatorMultiVariateWithDiscreteSymbolic imp throw new Exception("Local method not implemented yet"); } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { return condMiCalc.computeSignificance(numPermutationsToCheck); } - public MeasurementDistribution computeSignificance(int[][] newOrderings) + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { return condMiCalc.computeSignificance(newOrderings); } diff --git a/java/source/infodynamics/measures/continuous/symbolic/MutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java b/java/source/infodynamics/measures/continuous/symbolic/MutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java index 7d7d849..082e749 100755 --- a/java/source/infodynamics/measures/continuous/symbolic/MutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java +++ b/java/source/infodynamics/measures/continuous/symbolic/MutualInfoCalculatorMultiVariateWithDiscreteSymbolic.java @@ -5,7 +5,7 @@ import infodynamics.measures.discrete.MutualInformationCalculator; import infodynamics.utils.FirstIndexComparatorDouble; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; public class MutualInfoCalculatorMultiVariateWithDiscreteSymbolic implements @@ -139,12 +139,12 @@ public class MutualInfoCalculatorMultiVariateWithDiscreteSymbolic implements throw new Exception("Local method not implemented yet"); } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { return miCalc.computeSignificance(numPermutationsToCheck); } - public MeasurementDistribution computeSignificance(int[][] newOrderings) + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { return miCalc.computeSignificance(newOrderings); } diff --git a/java/source/infodynamics/measures/continuous/symbolic/TransferEntropyCalculatorSymbolic.java b/java/source/infodynamics/measures/continuous/symbolic/TransferEntropyCalculatorSymbolic.java index 87aadae..2137c60 100755 --- a/java/source/infodynamics/measures/continuous/symbolic/TransferEntropyCalculatorSymbolic.java +++ b/java/source/infodynamics/measures/continuous/symbolic/TransferEntropyCalculatorSymbolic.java @@ -6,7 +6,7 @@ import infodynamics.measures.discrete.ApparentTransferEntropyCalculator; import infodynamics.utils.FirstIndexComparatorDouble; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; import java.util.Iterator; @@ -300,7 +300,7 @@ public class TransferEntropyCalculatorSymbolic return locals; } - public MeasurementDistribution computeSignificance( + public EmpiricalMeasurementDistribution computeSignificance( int numPermutationsToCheck) throws Exception { return teCalc.computeSignificance(numPermutationsToCheck); } @@ -310,7 +310,7 @@ public class TransferEntropyCalculatorSymbolic * which is not strictly what this method is meant to do. * */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) throws Exception { System.out.println("TESymbolic.computeSignificance(): Not using the new orderings supplied"); return teCalc.computeSignificance(newOrderings.length); diff --git a/java/source/infodynamics/measures/discrete/ApparentTransferEntropyCalculator.java b/java/source/infodynamics/measures/discrete/ApparentTransferEntropyCalculator.java index 10618a6..71423f0 100755 --- a/java/source/infodynamics/measures/discrete/ApparentTransferEntropyCalculator.java +++ b/java/source/infodynamics/measures/discrete/ApparentTransferEntropyCalculator.java @@ -1,7 +1,7 @@ package infodynamics.measures.discrete; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; /** @@ -653,7 +653,7 @@ public class ApparentTransferEntropyCalculator extends ContextOfPastMeasureCalcu * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) { + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) { double actualTE = computeAverageLocalOfObservations(); // Reconstruct the source values (not necessarily in order) @@ -698,7 +698,7 @@ public class ApparentTransferEntropyCalculator extends ContextOfPastMeasureCalcu ate2.pastCount = pastCount; ate2.nextPastCount = nextPastCount; int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source int[] newSourceData = MatrixUtils.extractSelectedTimePoints(sourceValues, newOrderings[p]); diff --git a/java/source/infodynamics/measures/discrete/ChannelCalculator.java b/java/source/infodynamics/measures/discrete/ChannelCalculator.java index 016a040..2c9cac3 100755 --- a/java/source/infodynamics/measures/discrete/ChannelCalculator.java +++ b/java/source/infodynamics/measures/discrete/ChannelCalculator.java @@ -1,6 +1,6 @@ package infodynamics.measures.discrete; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * An interface for calculators computing measures from a source to a destination. @@ -47,5 +47,5 @@ public interface ChannelCalculator { * @param numPermutationsToCheck * @return */ - public MeasurementDistribution computeSignificance(int numPermutationsToCheck); + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck); } diff --git a/java/source/infodynamics/measures/discrete/CompleteTransferEntropyCalculator.java b/java/source/infodynamics/measures/discrete/CompleteTransferEntropyCalculator.java index a8ddf00..33d9055 100755 --- a/java/source/infodynamics/measures/discrete/CompleteTransferEntropyCalculator.java +++ b/java/source/infodynamics/measures/discrete/CompleteTransferEntropyCalculator.java @@ -2,7 +2,7 @@ package infodynamics.measures.discrete; import infodynamics.utils.MathsUtils; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; @@ -342,7 +342,7 @@ public class CompleteTransferEntropyCalculator extends InfoMeasureCalculator { * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) { + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) { double actualTE = computeAverageLocalOfObservations(); // Reconstruct the source values (not necessarily in order) @@ -400,7 +400,7 @@ public class CompleteTransferEntropyCalculator extends InfoMeasureCalculator { cte.pastOthersCount = pastOthersCount; cte.destPastOthersCount = destPastOthersCount; int countWhereTeIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the source int[] newSourceData = MatrixUtils.extractSelectedTimePoints(sourceValues, newOrderings[p]); diff --git a/java/source/infodynamics/measures/discrete/MutualInformationCalculator.java b/java/source/infodynamics/measures/discrete/MutualInformationCalculator.java index 7439d46..3946e48 100755 --- a/java/source/infodynamics/measures/discrete/MutualInformationCalculator.java +++ b/java/source/infodynamics/measures/discrete/MutualInformationCalculator.java @@ -1,7 +1,7 @@ package infodynamics.measures.discrete; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.RandomGenerator; @@ -185,7 +185,7 @@ public class MutualInformationCalculator extends InfoMeasureCalculator * @param numPermutationsToCheck number of new orderings of the source values to compare against * @return */ - public MeasurementDistribution computeSignificance(int numPermutationsToCheck) { + public EmpiricalMeasurementDistribution computeSignificance(int numPermutationsToCheck) { RandomGenerator rg = new RandomGenerator(); int[][] newOrderings = rg.generateDistinctRandomPerturbations(observations, numPermutationsToCheck); return computeSignificance(newOrderings); @@ -197,7 +197,7 @@ public class MutualInformationCalculator extends InfoMeasureCalculator * @param newOrderings the reorderings to use * @return */ - public MeasurementDistribution computeSignificance(int[][] newOrderings) { + public EmpiricalMeasurementDistribution computeSignificance(int[][] newOrderings) { double actualMI = computeAverageLocalOfObservations(); int numPermutationsToCheck = newOrderings.length; @@ -222,7 +222,7 @@ public class MutualInformationCalculator extends InfoMeasureCalculator mi2.iCount = iCount; mi2.jCount = jCount; int countWhereMIIsMoreSignificantThanOriginal = 0; - MeasurementDistribution measDistribution = new MeasurementDistribution(numPermutationsToCheck); + EmpiricalMeasurementDistribution measDistribution = new EmpiricalMeasurementDistribution(numPermutationsToCheck); for (int p = 0; p < numPermutationsToCheck; p++) { // Generate a new re-ordered data set for the i variable int[] newDataI = MatrixUtils.extractSelectedTimePoints(iValues, newOrderings[p]); diff --git a/java/source/infodynamics/networkinference/interregional/InterregionalChannelMeasure.java b/java/source/infodynamics/networkinference/interregional/InterregionalChannelMeasure.java index a749509..ab91bef 100755 --- a/java/source/infodynamics/networkinference/interregional/InterregionalChannelMeasure.java +++ b/java/source/infodynamics/networkinference/interregional/InterregionalChannelMeasure.java @@ -5,7 +5,7 @@ import infodynamics.measures.continuous.MutualInfoCalculatorMultiVariate; import infodynamics.measures.continuous.TransferEntropyCalculatorMultiVariate; import infodynamics.utils.ArrayFileReader; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; import infodynamics.utils.ParsedProperties; import infodynamics.utils.RandomGenerator; @@ -455,7 +455,7 @@ public abstract class InterregionalChannelMeasure { // Compute the measure for set s measureForEachSet[s] = channelCalc.computeAverageLocalOfObservations(); - MeasurementDistribution measDist; + EmpiricalMeasurementDistribution measDist; if (allValid || !validityForIndividualElements) { // Ask the TE calculator to work out the significance for us measDist = channelCalc.computeSignificance(reorderings); @@ -463,7 +463,7 @@ public abstract class InterregionalChannelMeasure { // We need to explicitly reorder including the individual validities. // Can't ask the calculator to do it, as it will only use the source-dest // pairings it originally used, which won't match across all subsets - measDist = new MeasurementDistribution(reorderingsForSignificance); + measDist = new EmpiricalMeasurementDistribution(reorderingsForSignificance); measDist.actualValue = measureForEachSet[s]; int countWhereReorderedIsMoreSignificantThanOriginal = 0; for (int p = 0; p < reorderingsForSignificance; p++) { diff --git a/java/source/infodynamics/networkinference/interregional/MeasurementDistributionPermutationsOverSubsets.java b/java/source/infodynamics/networkinference/interregional/MeasurementDistributionPermutationsOverSubsets.java index 2a3e796..04def6e 100755 --- a/java/source/infodynamics/networkinference/interregional/MeasurementDistributionPermutationsOverSubsets.java +++ b/java/source/infodynamics/networkinference/interregional/MeasurementDistributionPermutationsOverSubsets.java @@ -4,7 +4,7 @@ package infodynamics.networkinference.interregional; import infodynamics.utils.MatrixUtils; -import infodynamics.utils.MeasurementDistribution; +import infodynamics.utils.EmpiricalMeasurementDistribution; /** * Extends MeasurementDistribution for computations over @@ -17,7 +17,7 @@ import infodynamics.utils.MeasurementDistribution; * @author Joseph Lizier * */ -public class MeasurementDistributionPermutationsOverSubsets extends MeasurementDistribution { +public class MeasurementDistributionPermutationsOverSubsets extends EmpiricalMeasurementDistribution { // The true measurements for each subset s double[] actualValues; diff --git a/java/source/infodynamics/utils/EmpiricalMeasurementDistribution.java b/java/source/infodynamics/utils/EmpiricalMeasurementDistribution.java new file mode 100755 index 0000000..592af72 --- /dev/null +++ b/java/source/infodynamics/utils/EmpiricalMeasurementDistribution.java @@ -0,0 +1,88 @@ +package infodynamics.utils; + +/** + * + * Structure to hold a distribution of info-theoretic measurements, + * and a significance value for how an original measurement compared + * with these. + * + * @author Joseph Lizier + * + */ +public class EmpiricalMeasurementDistribution extends MeasurementDistribution { + + /** + * Distribution of surrogate measurement values + */ + public double[] distribution; + /** + * Whether the mean of the surrogate measurement distribution has + * been computed + */ + protected boolean computedMean = false; + /** + * Computed mean of the surrogate measurement distribution + */ + protected double meanOfDist; + /** + * Computed mean of the surrogate measurement distribution + */ + protected double stdOfDist; + + public EmpiricalMeasurementDistribution(int size) { + super(); // Creating the super class with mean and pValue 0 + // These value will be filled out by the caller later. + distribution = new double[size]; + } + + public EmpiricalMeasurementDistribution(double[] distribution, double actualValue) { + super(actualValue, 0); // Using pValue = 0 temporarily ... + this.distribution = distribution; + int countWhereActualIsNotGreater = 0; + for (int i = 0; i < distribution.length; i++) { + if (distribution[i] >= actualValue) { + countWhereActualIsNotGreater++; + } + } + pValue = (double) countWhereActualIsNotGreater / (double) distribution.length; + } + + // TODO Compute the significance under the assumption of a Gaussian distribution + /* + public double computeGaussianSignificance() { + // Need to conpute the significance based on the assumption of + // an underlying Gaussian distribution. + // Use the t distribution for analysis, since we have a finite + // number of samples to comptue the mean and std from. + return 0; + } + */ + + public double getTSscore() { + if (! computedMean) { + meanOfDist = MatrixUtils.mean(distribution); + stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); + computedMean = true; + } + double t = (actualValue - meanOfDist) / stdOfDist; + return t; + } + + public double getMeanOfDistribution() { + if (! computedMean) { + meanOfDist = MatrixUtils.mean(distribution); + stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); + computedMean = true; + } + return meanOfDist; + } + + public double getStdOfDistribution() { + if (! computedMean) { + meanOfDist = MatrixUtils.mean(distribution); + stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); + computedMean = true; + } + return stdOfDist; + } +} diff --git a/java/source/infodynamics/utils/MathsUtils.java b/java/source/infodynamics/utils/MathsUtils.java index ca1c957..7f67385 100755 --- a/java/source/infodynamics/utils/MathsUtils.java +++ b/java/source/infodynamics/utils/MathsUtils.java @@ -281,6 +281,98 @@ public class MathsUtils { return result; } + /** + *

Return the value of the cummulative distribution function of the + * chi-square distribution, evaluated at x, for k degrees of freedom.

+ * + *

Note that this relies on our approximation of the error function, + * which is the limiting part of the accuracy. Testing against + * values produced by octave indicates this is accurate to 5-6 + * decimal places.

+ * + * @param x value at which to evaluate the CDF + * @param k degrees of freedom (must have k>0) + * @return chi squared CDF evaluated at x given k degrees of freedom + * @see {@link http://en.wikipedia.org/wiki/Chi-squared_distribution} + */ + public static double chiSquareCdf(double x, int k) { + if (k <= 0) { + throw new IllegalArgumentException("k (" + k + ") must be > 0"); + } + return lowerIncompleteGammaFunctionOfArgsOn2(k,x) / + gammaOfArgOn2Plus1(k-2); // denominator is Gamma(k/2) + } + + /** + * Return the value of the lower Incomplete Gamma function, + * given arguments s/2 and x/2. + * We assume postive integer parameter s (s could be complex in general, + * with positive real part, but we restrict it to real and integer for + * this method). We make the evaluation using a recurrence relation, + * which terminates at s/2 = 1 or 1/2 (i.e. s = 2 or 1) + * + * @param s for parameter s/2 to lower incomplete gamma + * @param x for value x/2 to lower incomplete gamma + * @return value of lower gamma incomplete function + * @see {@link http://en.wikipedia.org/wiki/Incomplete_Gamma_function} + */ + public static double lowerIncompleteGammaFunctionOfArgsOn2(int s, double x) { + if (s <= 0) { + throw new IllegalArgumentException("s must be > 0"); + } + if (s == 2) { + // Terminating condition: evaluate lower gamma(1, x/2): + return 1 - Math.exp(-x/2.0); + } else if (s == 1) { + // Terminating condition: evaluate lower gamma(1/2, x/2): + return Math.sqrt(Math.PI) * erf(Math.sqrt(x/2.0)); + } else { + // Else evaluate recurrence relation: + return (s/2.0-1.0)*lowerIncompleteGammaFunctionOfArgsOn2(s-2,x) - + Math.pow(x/2.0, s/2.0 - 1.0) * Math.exp(-x/2.0); + } + } + + /** + * Return the value of the error function at a given x. + * We approximate the error function using elementary functions + * as described at the link below (quoting Abramowitz and Stegun). + * This approximation is quoted to + * have maximum error 1.5e-7 (and indeed this appears to be the + * case in comparison to values produced by octave). + * + * @param x value at which to evaluate the error function + * @return erf(x) + * @see {@link http://en.wikipedia.org/wiki/Error_function#Approximation_with_elementary_functions} + * @see Abramowitz, Milton; Stegun, Irene A., eds. (1972), + * "Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables", + * New York: Dover Publications, ISBN 978-0-486-61272-0 + */ + public static double erf(double x) { + // Constants: + double p = 0.3275911; + double[] a = {0.254829592, -0.284496736, 1.421413741, + -1.453152027, 1.061405429}; + boolean negArg = (x < 0); + + if (negArg) { + // The rest of the method requires x >= 0, but since erf(x) + // is an odd function, we just reflect x. + x = -x; + } + + double t = 1.0 / (1 + p * x); + double multiplier = 0.0; + double tToPower = t; + for (int i = 0; i < 5; i++) { + multiplier += a[i] * tToPower; + tToPower *= t; + } + double retVal = 1.0 - multiplier * Math.exp(-x*x); + // Remember that erf(x) was an odd function: + return negArg ? - retVal: retVal; + } + /** * Return the number of possible combinations of p from n (i.e. n choose p) * @@ -359,6 +451,12 @@ public class MathsUtils { return upToSetNum; } + /** + * Perform some testing: + * + * @param args + * @throws Exception + */ public static void main(String args[]) throws Exception { /* System.out.println(numOfSets(158,4)); @@ -368,6 +466,7 @@ public class MathsUtils { // int[][] sets = generateAllSets(6,4); // MatrixUtils.printMatrix(System.out, sets); + /* System.out.printf("digamma() digammaOld()\n"); for (int n = 0; n < 100; n++) { System.out.printf("%d %.3f %.3f\n", n, MathsUtils.digamma(n), MathsUtils.digammaByDefinition(n)); @@ -375,5 +474,23 @@ public class MathsUtils { for (int n = 0; n < 101; n++) { System.out.printf("%d %.3f %.3f\n", n, MathsUtils.digamma(n), MathsUtils.digammaByDefinition(n)); } + */ + + /* + System.out.println("erf(" + 1 + ")= " + MathsUtils.erf(1)); + System.out.println("erf(" + 2 + ")= " + MathsUtils.erf(2)); + System.out.println("erf(" + 0 + ")= " + MathsUtils.erf(0)); + for (int n=0; n<100; n++) { + System.out.println("erf(" + n*0.1 + ")= " + MathsUtils.erf(n*0.1)); + }*/ + + int degFree = 10; + System.out.println("chi2cdf(1," + degFree +")= " + MathsUtils.chiSquareCdf(1, degFree)); + System.out.println("chi2cdf(2," + degFree +")= " + MathsUtils.chiSquareCdf(2, degFree)); + System.out.println("chi2cdf(3," + degFree +")= " + MathsUtils.chiSquareCdf(3, degFree)); + for (int n=0; n<100; n++) { + System.out.println("chi2cdf(" + n*0.1 + "," + degFree +")= " + MathsUtils.chiSquareCdf(n*0.1, degFree)); + } + } } diff --git a/java/source/infodynamics/utils/MatrixUtils.java b/java/source/infodynamics/utils/MatrixUtils.java index efed5d0..bff17fb 100755 --- a/java/source/infodynamics/utils/MatrixUtils.java +++ b/java/source/infodynamics/utils/MatrixUtils.java @@ -13,6 +13,21 @@ import java.util.Vector; */ public class MatrixUtils { + /** + * Return an array with values enumerated through the given range + * + * @param startValue first value for the array + * @param endValue last value for the array + * @return + */ + public static int[] range(int startValue, int endValue) { + int[] array = new int[endValue - startValue + 1]; + for (int i = 0; i < endValue - startValue + 1; i++) { + array[i] = startValue + i; + } + return array; + } + public static double sum(double[] input) { double total = 0; for (int i = 0; i < input.length; i++) { @@ -2249,6 +2264,32 @@ public class MatrixUtils { return c / (double) data.length; } + /** + *

Returns the covariance between two columns of data in + * two multivariate arrays.

+ *

See - Mathworld + *

+ * + * @param data1 first multivariate array of data; first index is time, second is + * variable number + * @param data2 second multivariate array of data; first index is time, second is + * variable number + * @param col1 variable number 1 to compute the covariance to + * @param col2 variable number 2 to compute the covariance to + * @param mean1 mean of variable 1 + * @param mean2 mean of variable 2 + * @return the covariance + */ + public static double covarianceTwoColumns( + double[][] data1, double[][] data2, int col1, int col2, + double mean1, double mean2) { + double c = 0; + for (int t = 0; t < data1.length; t++) { + c += (data1[t][col1] - mean1)*(data2[t][col2]-mean2); + } + return c / (double) data1.length; + } + /** * Compute the covariance matrix between all column pairs (variables) in the * multivariate data set @@ -2274,9 +2315,72 @@ public class MatrixUtils { covariances[c][r] = covariances[r][c]; } } - return null; + return covariances; } + /** + * Compute the covariance matrix between all column pairs (variables) in the + * multivariate data set, which consists of two separate + * multivariate vectors. + * + * @param data1 multivariate array of data; first index is time, second is + * variable number + * @param data2 a second multivariate array of data, which can be though + * of as extensions of rows of the first. + * @return covariance matrix, where the columns of dat1 are numbered + * first, and the columns of data2 after that. + */ + public static double[][] covarianceMatrix( + double[][] data1, double[][] data2) { + int numVariables1 = data1[0].length; + int numVariables2 = data2[0].length; + int numVariables = numVariables1 + numVariables2; + double[][] covariances = new double[numVariables][numVariables]; + // Compute means of each variable once up front to save time + double[] means1 = new double[numVariables1]; + double[] means2 = new double[numVariables2]; + for (int r = 0; r < numVariables1; r++) { + means1[r] = mean(data1, r); + } + for (int r = 0; r < numVariables2; r++) { + means2[r] = mean(data2, r); + } + // Now compute the covariances: + for (int r = 0; r < numVariables1; r++) { + // Compute the covariances internal to data1: + for (int c = r; c < numVariables1; c++) { + // Compute the covariance between variable r and c: + covariances[r][c] = covarianceTwoColumns(data1, r, c, + means1[r], means1[c]); + // And of course this is symmetric between c and r: + covariances[c][r] = covariances[r][c]; + } + // Compute the covariances between data1 and data2 + for (int c = 0; c < numVariables2; c++) { + // Compute the covariance between variable r and c: + covariances[r][numVariables1 + c] = + covarianceTwoColumns(data1, data2, + r, c, means1[r], means2[c]); + // And of course this is symmetric between c and r: + covariances[numVariables1 + c][r] = + covariances[r][numVariables1 + c]; + } + } + // Now compute the covariances internal to data2: + for (int r = 0; r < numVariables2; r++) { + for (int c = r; c < numVariables2; c++) { + // Compute the covariance between variable r and c: + covariances[numVariables1 + r][numVariables1 + c] = + covarianceTwoColumns(data2, r, c, + means2[r], means2[c]); + // And of course this is symmetric between c and r: + covariances[numVariables1 + c][numVariables1 + r] = + covariances[numVariables1 + r][numVariables1 + c]; + } + } + return covariances; + } + /** *

Returns the correlation between the two arrays of data.

*

The arrays are asssumed to have the same lengths

diff --git a/java/source/infodynamics/utils/MeasurementDistribution.java b/java/source/infodynamics/utils/MeasurementDistribution.java index 7c48e0e..b5360a6 100755 --- a/java/source/infodynamics/utils/MeasurementDistribution.java +++ b/java/source/infodynamics/utils/MeasurementDistribution.java @@ -1,85 +1,44 @@ -package infodynamics.utils; - -/** - * - * Structure to hold a distribution of info-theoretic measurements, - * and a significance value for how an original measurement compared - * with these. - * - * @author Joseph Lizier - * - */ -public class MeasurementDistribution { - - /** - * Distribution of surrogate measurement values - */ - public double[] distribution; - /** - * Actual observed value of the measurement - */ - public double actualValue; - /** - * Probability that surrogate measurement is greater than - * the observed value - */ - public double pValue; - - protected boolean computedMean = false; - protected double meanOfDist; - protected double stdOfDist; - - public MeasurementDistribution(int size) { - distribution = new double[size]; - } - - public MeasurementDistribution(double[] distribution, double actualValue) { - this.actualValue = actualValue; - this.distribution = distribution; - int countWhereActualIsNotGreater = 0; - for (int i = 0; i < distribution.length; i++) { - if (distribution[i] >= actualValue) { - countWhereActualIsNotGreater++; - } - } - pValue = (double) countWhereActualIsNotGreater / (double) distribution.length; - } - - /* - public double computeGaussianSignificance() { - // Need to conpute the significance based on the assumption of - // an underlying Gaussian distribution. - // Use the t distribution for analysis, since we have a finite - // number of samples to comptue the mean and std from. - return 0; - } - */ - - public double getTSscore() { - if (! computedMean) { - meanOfDist = MatrixUtils.mean(distribution); - stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); - computedMean = true; - } - double t = (actualValue - meanOfDist) / stdOfDist; - return t; - } - - public double getMeanOfDistribution() { - if (! computedMean) { - meanOfDist = MatrixUtils.mean(distribution); - stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); - computedMean = true; - } - return meanOfDist; - } - - public double getStdOfDistribution() { - if (! computedMean) { - meanOfDist = MatrixUtils.mean(distribution); - stdOfDist = MatrixUtils.stdDev(distribution, meanOfDist); - computedMean = true; - } - return stdOfDist; - } -} +package infodynamics.utils; + +/** + *

Structure to hold a distribution of info-theoretic measurements, + * and a significance value for how an original measurement compared + * with these.

+ * + *

While in theory this class could be directly used, it is it's + * children which are intended to be used.

+ * + * @author Joseph Lizier + * + */ +public class MeasurementDistribution { + /** + * Actual observed value of the measurement + */ + public double actualValue; + /** + * Probability that surrogate measurement is greater than + * the observed value. + * (Small pValue means that the observed value is highly significant). + */ + public double pValue; + + /** + * Allow empty constructor for internal use only when actualValue and + * pValue will be supplied later + */ + protected MeasurementDistribution() { + } + + /** + * Construct with supplied actual value and p-value for it. + * + * @param actualValue actual observed value + * @param pValue p-value that the surrogate measurement is larger + * than the observed value + */ + public MeasurementDistribution(double actualValue, double pValue) { + this.actualValue = actualValue; + this.pValue = pValue; + } +}