diff --git a/java/unittests/infodynamics/measures/continuous/MultiInfoAbstractTester.java b/java/unittests/infodynamics/measures/continuous/MultiInfoAbstractTester.java new file mode 100755 index 0000000..e05486c --- /dev/null +++ b/java/unittests/infodynamics/measures/continuous/MultiInfoAbstractTester.java @@ -0,0 +1,106 @@ +/* + * Java Information Dynamics Toolkit (JIDT) + * Copyright (C) 2012, Joseph T. Lizier + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +package infodynamics.measures.continuous; + +import infodynamics.utils.EmpiricalMeasurementDistribution; +import infodynamics.utils.MatrixUtils; +import infodynamics.utils.RandomGenerator; +import junit.framework.TestCase; + +public abstract class MultiInfoAbstractTester extends TestCase { + + protected double lastResult = 0.0; + + /** + * Confirm that the local values average correctly back to the average value + * + * @param miCalc a pre-constructed MultiInfoCalculator object + * @param dimensions number of dimensions for the data to use + * @param timeSteps number of time steps for the random data + */ + public void testLocalsAverageCorrectly(MultiInfoCalculator miCalc, + int dimensions, int timeSteps) + throws Exception { + + miCalc.initialise(dimensions); + + // generate some random data + RandomGenerator rg = new RandomGenerator(); + double[][] data = rg.generateNormalData(timeSteps, dimensions, + 0, 1); + + miCalc.setObservations(data); + + //miCalc.setDebug(true); + double mi = miCalc.computeAverageLocalOfObservations(); + lastResult = mi; + //miCalc.setDebug(false); + double[] miLocal = miCalc.computeLocalOfPreviousObservations(); + + System.out.printf("Average was %.5f\n", mi); + + assertEquals(mi, MatrixUtils.mean(miLocal), 0.00001); + } + + /** + * Confirm that significance testing doesn't alter the average that + * would be returned. + * + * @param miCalc a pre-constructed MultiInfoCalculator object + * @param dimensions number of dimensions for the data to use + * @param timeSteps number of time steps for the random data + * @throws Exception + */ + public void testComputeSignificanceDoesntAlterAverage(MultiInfoCalculator miCalc, + int dimensions, int timeSteps) throws Exception { + + miCalc.initialise(dimensions); + + // generate some random data + RandomGenerator rg = new RandomGenerator(); + double[][] data = rg.generateNormalData(timeSteps, dimensions, + 0, 1); + + miCalc.setObservations(data); + + //miCalc.setDebug(true); + double mi = miCalc.computeAverageLocalOfObservations(); + //miCalc.setDebug(false); + //double[] miLocal = miCalc.computeLocalOfPreviousObservations(); + + System.out.printf("Average was %.5f\n", mi); + + // Now look at statistical significance tests + EmpiricalMeasurementDistribution measDist = + miCalc.computeSignificance(2); + // Make sure that (the first) surrogate TE does not + // match the actual TE (it could possibly match but with + // an incredibly low probability) + assertFalse(mi == measDist.distribution[0]); + + // And compute the average value again to check that it's consistent: + for (int i = 0; i < 10; i++) { + double lastAverage = miCalc.getLastAverage(); + assertEquals(mi, lastAverage); + double averageCheck1 = miCalc.computeAverageLocalOfObservations(); + assertEquals(mi, averageCheck1); + } + } + +} diff --git a/java/unittests/infodynamics/measures/continuous/kraskov/MultiInfoTester.java b/java/unittests/infodynamics/measures/continuous/kraskov/MultiInfoTester.java new file mode 100755 index 0000000..ffa8593 --- /dev/null +++ b/java/unittests/infodynamics/measures/continuous/kraskov/MultiInfoTester.java @@ -0,0 +1,344 @@ +/* + * Java Information Dynamics Toolkit (JIDT) + * Copyright (C) 2012, Joseph T. Lizier + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +package infodynamics.measures.continuous.kraskov; + +import infodynamics.measures.continuous.MultiInfoAbstractTester; +import infodynamics.utils.ArrayFileReader; +import infodynamics.utils.MatrixUtils; + +public class MultiInfoTester extends MultiInfoAbstractTester { + + protected String NUM_THREADS_TO_USE_DEFAULT = MultiInfoCalculatorKraskov.USE_ALL_THREADS; + protected String NUM_THREADS_TO_USE = NUM_THREADS_TO_USE_DEFAULT; + + /** + * Utility function to create a calculator for the given algorithm number + * + * @param algNumber + * @return + */ + public MultiInfoCalculatorKraskov getNewCalc(int algNumber) { + MultiInfoCalculatorKraskov miCalc = null; + if (algNumber == 1) { + miCalc = new MultiInfoCalculatorKraskov1(); + } else if (algNumber == 2) { + miCalc = new MultiInfoCalculatorKraskov2(); + } + return miCalc; + } + + /** + * Confirm that the local values average correctly back to the average value + * + * + */ + public void checkLocalsAverageCorrectly(int algNumber, String numThreads) throws Exception { + + MultiInfoCalculatorKraskov miCalc = getNewCalc(algNumber); + + String kraskov_K = "4"; + + miCalc.setProperty( + MultiInfoCalculatorKraskov.PROP_K, + kraskov_K); + miCalc.setProperty( + MultiInfoCalculatorKraskov.PROP_NUM_THREADS, + numThreads); + + super.testLocalsAverageCorrectly(miCalc, 2, 10000); + } + public void testLocalsAverageCorrectly() throws Exception { + checkLocalsAverageCorrectly(1, NUM_THREADS_TO_USE); + checkLocalsAverageCorrectly(2, NUM_THREADS_TO_USE); + } + + /** + * Confirm that significance testing doesn't alter the average that + * would be returned. + * + * @throws Exception + */ + public void checkComputeSignificanceDoesntAlterAverage(int algNumber) throws Exception { + + MultiInfoCalculatorKraskov miCalc = getNewCalc(algNumber); + + String kraskov_K = "4"; + + miCalc.setProperty( + MutualInfoCalculatorMultiVariateKraskov.PROP_K, + kraskov_K); + miCalc.setProperty( + MutualInfoCalculatorMultiVariateKraskov.PROP_NUM_THREADS, + NUM_THREADS_TO_USE); + + super.testComputeSignificanceDoesntAlterAverage(miCalc, 2, 100); + } + public void testComputeSignificanceDoesntAlterAverage() throws Exception { + checkComputeSignificanceDoesntAlterAverage(1); + checkComputeSignificanceDoesntAlterAverage(2); + } + + /** + * Utility function to run Kraskov MI for data with known results + * + * @param var1 + * @param kNNs array of Kraskov k nearest neighbours parameter to check + * @param expectedResults array of expected results for each k + */ + protected void checkMIForGivenData(double[][] data, + int[] kNNs, double[] expectedResults) throws Exception { + + // The Kraskov MILCA toolkit MIhigherdim executable + // uses algorithm 2 by default (this is what it means by rectangular): + MultiInfoCalculatorKraskov miCalc = getNewCalc(2); + + for (int kIndex = 0; kIndex < kNNs.length; kIndex++) { + int k = kNNs[kIndex]; + miCalc.setProperty( + MutualInfoCalculatorMultiVariateKraskov.PROP_K, + Integer.toString(k)); + miCalc.setProperty( + MutualInfoCalculatorMultiVariateKraskov.PROP_NUM_THREADS, + NUM_THREADS_TO_USE); + // No longer need to set this property as it's set by default: + //miCalc.setProperty(MutualInfoCalculatorMultiVariateKraskov.PROP_NORM_TYPE, + // EuclideanUtils.NORM_MAX_NORM_STRING); + miCalc.initialise(data[0].length); + miCalc.setObservations(data); + miCalc.setDebug(true); + double mi = miCalc.computeAverageLocalOfObservations(); + miCalc.setDebug(false); + + System.out.printf("k=%d: Average Multi-info %.8f (expected %.8f)\n", + k, mi, expectedResults[kIndex]); + // Dropping required accuracy by one order of magnitude, due + // to faster but slightly less accurate digamma estimator change + // 0.0000001 is fine for all but last test, so dropping again + // to 0.000001 + assertEquals(expectedResults[kIndex], mi, 0.000001); + } + } + + /** + * Test the computed Multi info for 2 variables (i.e. should be a regular MI!) + * against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIxnyn 1 1 3000 0 + * + * @throws Exception if file not found + * + */ + public void testUnivariateMIforRandomVariablesFromFile() throws Exception { + + // Test set 1: + + ArrayFileReader afr = new ArrayFileReader("demos/data/2randomCols-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 3, 4, 5, 6, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA = {-0.05294175, -0.03944338, -0.02190217, + 0.00120807, -0.00924771, -0.00316402, -0.00778205, -0.00565778}; + + System.out.println("Kraskov comparison 1 - univariate random data 1"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1}), + kNNs, expectedFromMILCA); + + //------------------ + // Test set 2: + + // We'll just take the first two columns from this data set + afr = new ArrayFileReader("demos/data/4randomCols-1.txt"); + data = afr.getDouble2DMatrix(); + + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {-0.04614525, -0.00861460, -0.00164540, + -0.01130354, -0.01339670, -0.00964035, -0.00237072, -0.00096891}; + + System.out.println("Kraskov comparison 2 - univariate random data 2"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1}), + kNNs, expectedFromMILCA_2); + + } + + /** + * Test the computed multivariate multi-info against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIhigherdim 4 1 1 3000 0 + * + * @throws Exception if file not found + * + */ + public void testMultivariateMIforRandomVariablesFromFile() throws Exception { + + // Test set 3: + + // We'll just take the first two columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/4randomCols-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 3, 4, 5, 6, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {0.03229833, -0.01146200, -0.00691358, + 0.00002149, -0.01056322, -0.01482730, -0.01223885, -0.01461794}; + + System.out.println("Kraskov comparison 3 - multivariate random data 1"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1, 2, 3}), + kNNs, expectedFromMILCA_2); + + } + + /** + * Test the computed multivariate multi-info against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIhigherdim 4 1 1 3000 0 + * + * @throws Exception if file not found + * + */ + public void testMultivariateMIVariousNumThreads() throws Exception { + + // Test set 3: + + // We'll just take the first two columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/4randomCols-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {3, 4}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {-0.00691358, + 0.00002149}; + + System.out.println("Kraskov comparison 3a - single threaded"); + NUM_THREADS_TO_USE = "1"; + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1, 2, 3}), + kNNs, expectedFromMILCA_2); + System.out.println("Kraskov comparison 3b - dual threaded"); + NUM_THREADS_TO_USE = "2"; + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1, 2, 3}), + kNNs, expectedFromMILCA_2); + NUM_THREADS_TO_USE = NUM_THREADS_TO_USE_DEFAULT; + } + + /** + * Test the computed multivariate MI against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIhigherdim 4 1 1 3000 0 + * + * @throws Exception if file not found + * + */ + public void testMultivariateMIforDependentVariablesFromFile() throws Exception { + + // Test set 6: + + // We'll just take the first two columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/4ColsPairedDirectDependence-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 3, 4, 5, 6, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {8.44056282, 7.69813699, 7.26909347, + 6.97095249, 6.73728113, 6.53105867, 5.96391264, 5.51627278}; + + System.out.println("Kraskov comparison 6 - multivariate dependent data 1"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1, 2, 3}), + kNNs, expectedFromMILCA_2); + + } + + /** + * Test the computed multivariate MI against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIhigherdim 4 1 1 3000 0 + * + * @throws Exception if file not found + * + */ + public void testMultivariateMIforNoisyDependentVariablesFromFile() throws Exception { + + // Test set 7: + + // We'll just take the first two columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/4ColsPairedNoisyDependence-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 3, 4, 5, 6, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {0.31900665, 0.37304998, 0.37213228, + 0.37982388, 0.37304217, 0.36802502, 0.36353436, 0.35095074}; + + System.out.println("Kraskov comparison 7 - multivariate dependent data 1"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0, 1, 2, 3}), + kNNs, expectedFromMILCA_2); + + } + + /** + * Test the computed multivariate MI against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIhigherdim 10 1 1 10000 0 + * + * @throws Exception if file not found + * + */ + public void testMultivariateMIforRandomGaussianVariablesFromFile() throws Exception { + + // Test set 8: + + // We'll take the columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/10ColsRandomGaussian-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 4, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {0.00932984, 0.00662195, 0.01697033, + 0.00397984, 0.00212609}; + + System.out.println("Kraskov comparison 8 - multivariate uncorrelated Gaussian data 1"); + checkMIForGivenData(MatrixUtils.selectColumns(data, + new int[] {0, 1, 2, 3, 4, 5, 6, 7, 8, 9}), + kNNs, expectedFromMILCA_2); + + } +} diff --git a/java/unittests/infodynamics/measures/continuous/kraskov/MutualInfoMultiVariateTester.java b/java/unittests/infodynamics/measures/continuous/kraskov/MutualInfoMultiVariateTester.java index 746ceaa..a2a5100 100755 --- a/java/unittests/infodynamics/measures/continuous/kraskov/MutualInfoMultiVariateTester.java +++ b/java/unittests/infodynamics/measures/continuous/kraskov/MutualInfoMultiVariateTester.java @@ -46,8 +46,6 @@ public class MutualInfoMultiVariateTester /** * Confirm that the local values average correctly back to the average value * - * TODO Add a test with say 10000 time steps, after we introduce fast nearest - * neighbour searching. * */ public void checkLocalsAverageCorrectly(int algNumber, String numThreads) throws Exception { @@ -63,7 +61,7 @@ public class MutualInfoMultiVariateTester MutualInfoCalculatorMultiVariateKraskov.PROP_NUM_THREADS, numThreads); - super.testLocalsAverageCorrectly(miCalc, 2, 100); + super.testLocalsAverageCorrectly(miCalc, 2, 10000); } public void testLocalsAverageCorrectly() throws Exception { checkLocalsAverageCorrectly(1, NUM_THREADS_TO_USE); @@ -297,8 +295,8 @@ public class MutualInfoMultiVariateTester * where the file has the first 30 rows repeated. * * Kraskov et al recommend that a small amount of noise should be - * added to the data to avoid issues with repeated scores; we have not - * implemented this yet. + * added to the data to avoid issues with repeated scores; this + * can be done in our toolkit by setting the relevant property * * @throws Exception if file not found * @@ -430,4 +428,36 @@ public class MutualInfoMultiVariateTester kNNs, expectedFromMILCA_2); } + + /** + * Test the computed multivariate MI against that calculated by Kraskov's own MILCA + * tool on the same data. + * + * To run Kraskov's tool (http://www.klab.caltech.edu/~kraskov/MILCA/) for this + * data, run: + * ./MIxnyn 1 1 10000 0 + * + * @throws Exception if file not found + * + */ + public void testMIforRandomGaussianVariablesFromLargeFile() throws Exception { + + // Test set 9: + + // We'll take the columns from this data set + ArrayFileReader afr = new ArrayFileReader("demos/data/10ColsRandomGaussian-1.txt"); + double[][] data = afr.getDouble2DMatrix(); + + // Use various Kraskov k nearest neighbours parameter + int[] kNNs = {1, 2, 4, 10, 15}; + // Expected values from Kraskov's MILCA toolkit: + double[] expectedFromMILCA_2 = {0.01542004, 0.01137151, 0.00210945, + 0.00159921, 0.00031277}; + + System.out.println("Kraskov comparison 9 - uncorrelated Gaussian data 1 - large file"); + checkMIForGivenData(MatrixUtils.selectColumns(data, new int[] {0}), + MatrixUtils.selectColumns(data, new int[] {1}), + kNNs, expectedFromMILCA_2); + + } } diff --git a/java/unittests/infodynamics/utils/KdTreeTest.java b/java/unittests/infodynamics/utils/KdTreeTest.java index 9665bf8..4f6c0e1 100755 --- a/java/unittests/infodynamics/utils/KdTreeTest.java +++ b/java/unittests/infodynamics/utils/KdTreeTest.java @@ -1,3 +1,21 @@ +/* + * Java Information Dynamics Toolkit (JIDT) + * Copyright (C) 2012, Joseph T. Lizier + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + package infodynamics.utils; import java.util.Calendar; @@ -5,7 +23,6 @@ import java.util.PriorityQueue; import java.util.Vector; import infodynamics.utils.KdTree.KdTreeNode; -import infodynamics.utils.KdTree.NeighbourNodeData; import junit.framework.TestCase; public class KdTreeTest extends TestCase { @@ -512,7 +529,7 @@ public class KdTreeTest extends TestCase { ((double) (endTimeValidate - startTime)/1000.0)); } - public void testCountNeighboursWithinRSeparateArraysWithDuplicates() { + public void testCountNeighboursWithinRSeparateArraysWithDuplicates() throws Exception { int variables = 3; int dimensionsPerVariable = 3; int samples = 2000; diff --git a/java/unittests/infodynamics/utils/UnivariateNearestNeighbourTest.java b/java/unittests/infodynamics/utils/UnivariateNearestNeighbourTest.java new file mode 100755 index 0000000..2b88089 --- /dev/null +++ b/java/unittests/infodynamics/utils/UnivariateNearestNeighbourTest.java @@ -0,0 +1,206 @@ +/* + * Java Information Dynamics Toolkit (JIDT) + * Copyright (C) 2012, Joseph T. Lizier + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +package infodynamics.utils; + +import java.util.Calendar; +import java.util.PriorityQueue; + +import junit.framework.TestCase; + +public class UnivariateNearestNeighbourTest extends TestCase { + + RandomGenerator rg = new RandomGenerator(); + + public void testSmallConstruction() throws Exception { + // Testing with an example from + // http://en.wikipedia.org/wiki/K-d_tree + + double[] data = { 2, 5, 9, 4, 8, 7 }; + + UnivariateNearestNeighbourSearcher searcher = + new UnivariateNearestNeighbourSearcher(data); + validateAllPointsInSearcher(searcher, data); + } + + public void testLargerConstruction() throws Exception { + double[] data = rg.generateNormalData(1000, 0, 1); + UnivariateNearestNeighbourSearcher searcher = + new UnivariateNearestNeighbourSearcher(data); + validateAllPointsInSearcher(searcher, data); + } + + public void validateAllPointsInSearcher(UnivariateNearestNeighbourSearcher searcher, + double[] data) { + // Check that the data are sorted in increasing order: + for (int t = 1; t < data.length; t++) { + assert(data[searcher.sortedArrayIndices[t]] >= + data[searcher.sortedArrayIndices[t-1]]); + } + // Check that the backreferences all work fine: + for (int t = 0; t < data.length; t++) { + assert(searcher.sortedArrayIndices[searcher.indicesInSortedArray[t]] == t); + } + } + + public void testNearestNeighbourSearch() throws Exception { + double[] data = rg.generateNormalData(10000, 0, 1); + UnivariateNearestNeighbourSearcher searcher = + new UnivariateNearestNeighbourSearcher(data); + int[] nearestNeighbourIndices = new int[data.length]; + long startTime = Calendar.getInstance().getTimeInMillis(); + for (int t = 0; t < data.length; t++) { + NeighbourNodeData nodeData = searcher.findNearestNeighbour(t); + nearestNeighbourIndices[t] = nodeData.sampleIndex; + } + long endTimeNNs = Calendar.getInstance().getTimeInMillis(); + System.out.printf("Found all nearest neighbours for %d points in: %.3f sec\n", + data.length, ((double) (endTimeNNs - startTime)/1000.0)); + + // Now do brute force: + int[] bruteForceNeighbourIndices = new int[data.length]; + for (int t = 0; t < data.length; t++) { + int neighbour = 0; + double minDist = Double.POSITIVE_INFINITY; + for (int t2 = 0; t2 < data.length; t2++) { + if (t2 == t) { + continue; + } + double norm = Math.abs(data[t] - data[t2]); + if (norm < minDist) { + minDist = norm; + neighbour = t2; + } + } + bruteForceNeighbourIndices[t] = neighbour; + } + long endTimeBruteForce = Calendar.getInstance().getTimeInMillis(); + System.out.printf("Found all nearest neighbours for %d points by brute force in: %.3f sec\n", + data.length, ((double) (endTimeBruteForce - endTimeNNs)/1000.0)); + + // Now validate the nearest neighbours: + for (int t = 0; t < data.length; t++) { + assertEquals(nearestNeighbourIndices[t], bruteForceNeighbourIndices[t]); + } + } + + public void testRangeFinderStrictWithin() throws Exception { + checkRangeFinder(true); + } + + public void testRangeFinderWithinOrEqual() throws Exception { + checkRangeFinder(false); + } + + public void checkRangeFinder(boolean strict) throws Exception { + int numTimeStepsInitial = 10000; + int duplicateSteps = 100; + int numTimeSteps = numTimeStepsInitial+duplicateSteps; + double[] dataRaw = rg.generateNormalData(numTimeStepsInitial, 0, 1); + double[] data = new double[numTimeSteps]; + // Now duplicate some of the time steps as a test: + System.arraycopy(dataRaw, 0, data, 0, numTimeStepsInitial); + System.arraycopy(dataRaw, 0, data, numTimeStepsInitial, duplicateSteps); + + UnivariateNearestNeighbourSearcher searcher = + new UnivariateNearestNeighbourSearcher(data); + int[] counts = new int[data.length]; + double r = 0.2; + long startTime = Calendar.getInstance().getTimeInMillis(); + for (int t = 0; t < data.length; t++) { + counts[t] = searcher.countPointsWithinR(t, r, !strict); + } + long endTimeNNs = Calendar.getInstance().getTimeInMillis(); + System.out.printf("Found all neighbours within %.3f (mean = %.3f) in: %.3f sec\n", + r, MatrixUtils.mean(counts), + ((double) (endTimeNNs - startTime)/1000.0)); + + // Now do brute force: + int[] bruteForceCounts = new int[data.length]; + for (int t = 0; t < data.length; t++) { + int count = 0; + for (int t2 = 0; t2 < data.length; t2++) { + if (t2 == t) { + continue; + } + double norm = Math.abs(data[t] - data[t2]); + if ((strict && (norm < r) ) || + (!strict && (norm <= r))) { + count++; + } + } + bruteForceCounts[t] = count; + } + long endTimeBruteForce = Calendar.getInstance().getTimeInMillis(); + System.out.printf("Found all neighbours within %.3f (mean = %.3f) by brute force in: %.3f sec\n", + r, MatrixUtils.mean(bruteForceCounts), + ((double) (endTimeBruteForce - endTimeNNs)/1000.0)); + + // Now validate the nearest neighbour counts: + for (int t = 0; t < data.length; t++) { + assertEquals(bruteForceCounts[t], counts[t]); + } + } + + public void testFindKNearestNeighbours() throws Exception { + for (int K = 1; K < 5; K++) { + double[] data = rg.generateNormalData(1000, 0, 1); + + long startTime = Calendar.getInstance().getTimeInMillis(); + UnivariateNearestNeighbourSearcher searcher = + new UnivariateNearestNeighbourSearcher(data); + long endTimeTree = Calendar.getInstance().getTimeInMillis(); + System.out.printf("Searcher of %d points for %d NNs constructed in: %.3f sec\n", + data.length, K, ((double) (endTimeTree - startTime)/1000.0)); + + startTime = Calendar.getInstance().getTimeInMillis(); + for (int t = 0; t < data.length; t++) { + PriorityQueue nnPQ = + searcher.findKNearestNeighbours(K, t); + assertTrue(nnPQ.size() == K); + // Now find the K nearest neighbours with a naive all-pairs comparison + double[][] distancesAndIndices = new double[data.length][2]; + for (int t2 = 0; t2 < data.length; t2++) { + if (t2 != t) { + distancesAndIndices[t2][0] = searcher.norm(data[t], data[t2]); + } else { + distancesAndIndices[t2][0] = Double.POSITIVE_INFINITY; + } + distancesAndIndices[t2][1] = t2; + } + int[] timeStepsOfKthMins = + MatrixUtils.kMinIndices(distancesAndIndices, 0, K); + for (int i = 0; i < K; i++) { + // Check that the ith nearest neighbour matches for each method. + // Note that these two method provide a different sorting order + NeighbourNodeData nnData = nnPQ.poll(); + if (timeStepsOfKthMins[K - 1 - i] != nnData.sampleIndex) { + // We have an error: + System.out.printf("Erroneous match between indices %d (expected) " + + " and %d\n", timeStepsOfKthMins[K - 1 - i], nnData.sampleIndex); + } + assertEquals(timeStepsOfKthMins[K - 1 - i], nnData.sampleIndex); + } + } + long endTimeValidate = Calendar.getInstance().getTimeInMillis(); + System.out.printf("All %d nearest neighbours found in: %.3f sec\n", + K, ((double) (endTimeValidate - startTime)/1000.0)); + } + } + +}