Unit tests for Kraskov multi-info calculators using fast neighbour search and common multi-info class. Additional mutual info unit test added also on larger data set, and unit tests for the UnivariateNearestNeighbourSearcher.

This commit is contained in:
joseph.lizier 2014-10-30 04:24:42 +00:00
parent 2a56301b75
commit 9aa0e9ca0e
5 changed files with 710 additions and 7 deletions

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
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);
}
}
}

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
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 <dataFile> 1 1 3000 <kNearestNeighbours> 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 <dataFile> 4 1 1 3000 <kNearestNeighbours> 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 <dataFile> 4 1 1 3000 <kNearestNeighbours> 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 <dataFile> 4 1 1 3000 <kNearestNeighbours> 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 <dataFile> 4 1 1 3000 <kNearestNeighbours> 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 <dataFile> 10 1 1 10000 <kNearestNeighbours> 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);
}
}

View File

@ -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 <dataFile> 1 1 10000 <kNearestNeighbours> 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);
}
}

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
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;

View File

@ -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 <http://www.gnu.org/licenses/>.
*/
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<NeighbourNodeData> 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));
}
}
}