mirror of https://github.com/jlizier/jidt
Unit tests for Dynamic correlation exclusion, or Theiler window, in all Kraskov estimators, plus direct validation against values from using Theiler window in TRENTOOL. Validates issue 35
This commit is contained in:
joseph.lizier
parent
733da8265d
commit
33ea4f6136
|
|
@ -124,6 +124,92 @@ public class TransferEntropyTester
|
|||
assertEquals(99, teCalc.getNumObservations());
|
||||
}
|
||||
|
||||
/**
|
||||
* Test the computed univariate TE
|
||||
* against that calculated by Wibral et al.'s TRENTOOL
|
||||
* on the same data, adding dynamic correlation exclusion
|
||||
*
|
||||
* To run TRENTOOL (http://www.trentool.de/) for this
|
||||
* data, run its TEvalues.m matlab script on the multivariate source
|
||||
* and dest data sets as:
|
||||
* TEvalues(source, dest, 1, 1, 1, kraskovK, dynCorrExcl)
|
||||
* with these values ensuring source-dest lag 1, history k=1,
|
||||
* embedding lag 1, and dynamic correlation exclusion window dynCorrExcl
|
||||
*
|
||||
* @throws Exception if file not found
|
||||
*
|
||||
*/
|
||||
public void testUnivariateTEforCoupledVariablesFromFileDynCorrExcl() throws Exception {
|
||||
// Test set 1:
|
||||
|
||||
ArrayFileReader afr = new ArrayFileReader("demos/data/2coupledRandomCols-1.txt");
|
||||
double[][] data = afr.getDouble2DMatrix();
|
||||
double[] col0 = MatrixUtils.selectColumn(data, 0);
|
||||
double[] col1 = MatrixUtils.selectColumn(data, 1);
|
||||
// Need to normalise these ourselves rather than letting the calculator do it -
|
||||
// this ensures the extra values in the time series (e.g. last value in source)
|
||||
// are taken into account, in line with TRENTOOL
|
||||
col0 = MatrixUtils.normaliseIntoNewArray(col0);
|
||||
col1 = MatrixUtils.normaliseIntoNewArray(col1);
|
||||
|
||||
// Use various Kraskov k nearest neighbours parameter
|
||||
int kNNs = 4;
|
||||
// Expected values from TRENTOOL for correlation exclusion window 10:
|
||||
int exclWindow = 10;
|
||||
double expectedFromTRENTOOL0to1 = 0.2930714;
|
||||
double expectedFromTRENTOOL1to0 = -0.0387031;
|
||||
|
||||
System.out.println("Kraskov TE comparison 1b to TRENTOOL - univariate coupled data with dynamic correlation exclusion");
|
||||
TransferEntropyCalculatorKraskov teCalc =
|
||||
new TransferEntropyCalculatorKraskov();
|
||||
teCalc.setProperty(ConditionalMutualInfoCalculatorMultiVariateKraskov.PROP_K, Integer.toString(kNNs));
|
||||
// We already normalised above, and this will do a different
|
||||
// normalisation without taking the extra values in to account if we did it
|
||||
teCalc.setProperty(ConditionalMutualInfoCalculatorMultiVariateKraskov.PROP_NORMALISE, "false");
|
||||
// Set dynamic correlation exclusion window:
|
||||
teCalc.setProperty(ConditionalMutualInfoCalculatorMultiVariateKraskov.PROP_DYN_CORR_EXCL_TIME,
|
||||
Integer.toString(exclWindow));
|
||||
|
||||
teCalc.initialise(1);
|
||||
teCalc.setObservations(col0, col1);
|
||||
double result = teCalc.computeAverageLocalOfObservations();
|
||||
System.out.printf("From 2coupledRandomCols 0->1, Theiler window 10, expecting %.6f, got %.6f\n",
|
||||
expectedFromTRENTOOL0to1, result);
|
||||
assertEquals(expectedFromTRENTOOL0to1, result, 0.000001);
|
||||
|
||||
teCalc.initialise(1);
|
||||
teCalc.setObservations(col1, col0);
|
||||
result = teCalc.computeAverageLocalOfObservations();
|
||||
assertEquals(expectedFromTRENTOOL1to0, result, 0.000001);
|
||||
System.out.printf("From 2coupledRandomCols 1->0, Theiler window 10, expecting %.6f, got %.6f\n",
|
||||
expectedFromTRENTOOL1to0, result);
|
||||
|
||||
assertEquals(99, teCalc.getNumObservations());
|
||||
|
||||
// Change dynamic correlation exclusion window:
|
||||
exclWindow = 20;
|
||||
expectedFromTRENTOOL0to1 = 0.2995997;
|
||||
expectedFromTRENTOOL1to0 = -0.0381608;
|
||||
teCalc.setProperty(ConditionalMutualInfoCalculatorMultiVariateKraskov.PROP_DYN_CORR_EXCL_TIME,
|
||||
Integer.toString(exclWindow));
|
||||
|
||||
teCalc.initialise(1);
|
||||
teCalc.setObservations(col0, col1);
|
||||
result = teCalc.computeAverageLocalOfObservations();
|
||||
System.out.printf("From 2coupledRandomCols 0->1, Theiler window 20, expecting %.6f, got %.6f\n",
|
||||
expectedFromTRENTOOL0to1, result);
|
||||
assertEquals(expectedFromTRENTOOL0to1, result, 0.000001);
|
||||
|
||||
teCalc.initialise(1);
|
||||
teCalc.setObservations(col1, col0);
|
||||
result = teCalc.computeAverageLocalOfObservations();
|
||||
assertEquals(expectedFromTRENTOOL1to0, result, 0.000001);
|
||||
System.out.printf("From 2coupledRandomCols 1->0, Theiler window 20, expecting %.6f, got %.6f\n",
|
||||
expectedFromTRENTOOL1to0, result);
|
||||
|
||||
assertEquals(99, teCalc.getNumObservations());
|
||||
}
|
||||
|
||||
/**
|
||||
* Test the computed univariate TE
|
||||
* against that calculated by Wibral et al.'s TRENTOOL
|
||||
|
|
|
|||
|
|
@ -335,6 +335,11 @@ public class KdTreeTest extends TestCase {
|
|||
PriorityQueue<NeighbourNodeData> nnPQ =
|
||||
kdTree.findKNearestNeighbours(K, t);
|
||||
assertTrue(nnPQ.size() == K);
|
||||
// Also check the result is still correct if we search with
|
||||
// zero size exclusion window:
|
||||
PriorityQueue<NeighbourNodeData> nnPQ_zeroExclWindow =
|
||||
kdTree.findKNearestNeighbours(K, t, 0);
|
||||
assertTrue(nnPQ_zeroExclWindow.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++) {
|
||||
|
|
@ -345,6 +350,58 @@ public class KdTreeTest extends TestCase {
|
|||
}
|
||||
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();
|
||||
NeighbourNodeData nnData_zeroWindow = nnPQ_zeroExclWindow.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);
|
||||
assertEquals(timeStepsOfKthMins[K - 1 - i], nnData_zeroWindow.sampleIndex);
|
||||
}
|
||||
}
|
||||
long endTimeValidate = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All %d nearest neighbours found in: %.3f sec\n",
|
||||
K, ((double) (endTimeValidate - startTime)/1000.0));
|
||||
}
|
||||
}
|
||||
|
||||
public void testFindKNearestNeighboursWithExclusionWindow() throws Exception {
|
||||
int dimension = 4;
|
||||
int numSamples = 1000;
|
||||
int exclusionWindow = 100;
|
||||
|
||||
for (int K = 1; K < 5; K++) {
|
||||
double[][] data = rg.generateNormalData(numSamples, dimension, 0, 1);
|
||||
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
KdTree kdTree = new KdTree(data);
|
||||
long endTimeTree = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("Tree of %d points for %d NNs constructed in: %.3f sec\n",
|
||||
data.length, K, ((double) (endTimeTree - startTime)/1000.0));
|
||||
|
||||
EuclideanUtils normCalculator = new EuclideanUtils(EuclideanUtils.NORM_MAX_NORM);
|
||||
startTime = Calendar.getInstance().getTimeInMillis();
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
PriorityQueue<NeighbourNodeData> nnPQ =
|
||||
kdTree.findKNearestNeighbours(K, t, exclusionWindow);
|
||||
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 (Math.abs(t2 - t) > exclusionWindow) {
|
||||
distancesAndIndices[t2][0] = normCalculator.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++) {
|
||||
|
|
@ -637,6 +694,199 @@ public class KdTreeTest extends TestCase {
|
|||
((double) (endTimeValidate - endTimeCount)/1000.0));
|
||||
}
|
||||
|
||||
public void testCountNeighboursWithinRExclusionWindow() {
|
||||
int dimensionsPerVariable = 3;
|
||||
int samples = 1000;
|
||||
double[][] data = rg.generateNormalData(samples, dimensionsPerVariable, 0, 1);
|
||||
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
KdTree kdTree = new KdTree(data);
|
||||
long endTimeTree = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("Tree of %d points constructed in: %.3f sec\n",
|
||||
samples, ((double) (endTimeTree - startTime)/1000.0));
|
||||
|
||||
double[] rs = {0.2, 0.6, 0.8, 1.0};
|
||||
for (int i = 0; i < rs.length; i++) {
|
||||
verifyCountNeighboursWithinRExclusionWindow(EuclideanUtils.NORM_MAX_NORM,
|
||||
data, samples, kdTree, rs[i], true);
|
||||
verifyCountNeighboursWithinRExclusionWindow(EuclideanUtils.NORM_MAX_NORM,
|
||||
data, samples, kdTree, rs[i], false);
|
||||
// and check using the array methods as well
|
||||
verifyCountNeighboursWithinRExclusionWindowArrayArgs(EuclideanUtils.NORM_MAX_NORM,
|
||||
data, samples, kdTree, rs[i], true);
|
||||
verifyCountNeighboursWithinRExclusionWindowArrayArgs(EuclideanUtils.NORM_MAX_NORM,
|
||||
data, samples, kdTree, rs[i], false);
|
||||
}
|
||||
}
|
||||
|
||||
private void verifyCountNeighboursWithinRExclusionWindow(int normType,
|
||||
double[][] data, int samples, KdTree kdTree, double r,
|
||||
boolean allowEqualToR) {
|
||||
|
||||
int exclWindow = 100;
|
||||
|
||||
// Set up the given norm type:
|
||||
EuclideanUtils normCalculator = new EuclideanUtils(normType);
|
||||
kdTree.setNormType(normType);
|
||||
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
int totalCount = 0;
|
||||
int[] counts = new int[samples];
|
||||
for (int t = 0; t < samples; t++) {
|
||||
int count =
|
||||
kdTree.countPointsWithinR(t, r, exclWindow, allowEqualToR);
|
||||
assertTrue(count >= 0);
|
||||
counts[t] = count;
|
||||
totalCount += count;
|
||||
}
|
||||
long endTimeCount = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All neighbours within %.3f (average of %.3f) found in: %.3f sec\n",
|
||||
r, (double) totalCount / (double) samples,
|
||||
((double) (endTimeCount - startTime)/1000.0));
|
||||
// Now find the neighbour count with a naive all-pairs comparison
|
||||
for (int t = 0; t < samples; t++) {
|
||||
int naiveCount = 0;
|
||||
for (int t2 = 0; t2 < samples; t2++) {
|
||||
double norm = Double.POSITIVE_INFINITY;
|
||||
if (Math.abs(t2 - t) > exclWindow) {
|
||||
norm = normCalculator.norm(
|
||||
data[t], data[t2]);
|
||||
}
|
||||
if ((!allowEqualToR && (norm < r)) ||
|
||||
( allowEqualToR && (norm <= r))) {
|
||||
naiveCount++;
|
||||
}
|
||||
}
|
||||
if (naiveCount != counts[t]) {
|
||||
// Add some debug prints:
|
||||
System.out.printf("All pairs : count %d\n", naiveCount);
|
||||
System.out.printf("kdTree search: count %d\n", counts[t]);
|
||||
Collection<NeighbourNodeData> pointsWithinR =
|
||||
kdTree.findPointsWithinR(t, r, allowEqualToR);
|
||||
for (NeighbourNodeData nnData : pointsWithinR) {
|
||||
System.out.printf("kd: t=%d: time diff %d, norm %.3f\n",
|
||||
t, Math.abs(nnData.sampleIndex - t),
|
||||
normCalculator.norm(
|
||||
data[t], data[nnData.sampleIndex]));
|
||||
}
|
||||
// and compare to all pairs:
|
||||
for (int t2 = 0; t2 < samples; t2++) {
|
||||
double norm = Double.POSITIVE_INFINITY;
|
||||
if (Math.abs(t2 - t) > exclWindow) {
|
||||
norm = normCalculator.norm(
|
||||
data[t], data[t2]);
|
||||
}
|
||||
if ((!allowEqualToR && (norm < r)) ||
|
||||
( allowEqualToR && (norm <= r))) {
|
||||
System.out.printf("naive: t=%d: time diff %d, norm %.3f\n",
|
||||
t, Math.abs(t2 - t), norm);
|
||||
}
|
||||
}
|
||||
}
|
||||
assertEquals(naiveCount, counts[t]);
|
||||
}
|
||||
long endTimeValidate = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All neighbours within %.3f (average of %.3f) validated in: %.3f sec\n",
|
||||
r, (double) totalCount / (double) samples,
|
||||
((double) (endTimeValidate - endTimeCount)/1000.0));
|
||||
}
|
||||
|
||||
private void verifyCountNeighboursWithinRExclusionWindowArrayArgs(int normType,
|
||||
double[][] data, int samples, KdTree kdTree, double r,
|
||||
boolean allowEqualToR) {
|
||||
|
||||
int exclWindow = 100;
|
||||
|
||||
// Set up the given norm type:
|
||||
EuclideanUtils normCalculator = new EuclideanUtils(normType);
|
||||
kdTree.setNormType(normType);
|
||||
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
int totalCount = 0;
|
||||
int[] counts = new int[samples];
|
||||
boolean[] withinR = new boolean[samples];
|
||||
int[] indicesWithinR = new int[samples];
|
||||
for (int t = 0; t < samples; t++) {
|
||||
// int count = kdTree.countPointsWithinR(t, r, allowEqualToR);
|
||||
kdTree.findPointsWithinR(t, r, exclWindow, allowEqualToR,
|
||||
withinR, indicesWithinR);
|
||||
int count = 0;
|
||||
// Run through the list of points returned as within r and check them
|
||||
for (int t2 = 0; indicesWithinR[t2] != -1; t2++) {
|
||||
count++;
|
||||
assertTrue(withinR[indicesWithinR[t2]]);
|
||||
// Reset this marker
|
||||
withinR[indicesWithinR[t2]] = false;
|
||||
// Check that it's really within r
|
||||
double maxNorm = normCalculator.norm(
|
||||
data[t], data[indicesWithinR[t2]]);
|
||||
assertTrue((!allowEqualToR && (maxNorm < r)) ||
|
||||
( allowEqualToR && (maxNorm <= r)));
|
||||
}
|
||||
// Check that there were no other points marked as within r:
|
||||
int count2 = 0;
|
||||
for (int t2 = 0; t2 < samples; t2++) {
|
||||
if (withinR[t2]) {
|
||||
count2++;
|
||||
}
|
||||
}
|
||||
assertEquals(0, count2); // We set all the ones that should have been true to false already
|
||||
// Now reset the boolean array
|
||||
counts[t] = count;
|
||||
totalCount += count;
|
||||
}
|
||||
long endTimeCount = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All neighbours within %.3f (average of %.3f) found in: %.3f sec\n",
|
||||
r, (double) totalCount / (double) samples,
|
||||
((double) (endTimeCount - startTime)/1000.0));
|
||||
// Now find the neighbour count with a naive all-pairs comparison
|
||||
for (int t = 0; t < samples; t++) {
|
||||
int naiveCount = 0;
|
||||
for (int t2 = 0; t2 < samples; t2++) {
|
||||
double norm = Double.POSITIVE_INFINITY;
|
||||
if (Math.abs(t2 - t) > exclWindow) {
|
||||
norm = normCalculator.norm(
|
||||
data[t], data[t2]);
|
||||
}
|
||||
if ((!allowEqualToR && (norm < r)) ||
|
||||
( allowEqualToR && (norm <= r))) {
|
||||
naiveCount++;
|
||||
}
|
||||
}
|
||||
if (naiveCount != counts[t]) {
|
||||
// Add some debug prints:
|
||||
System.out.printf("All pairs : count %d\n", naiveCount);
|
||||
System.out.printf("kdTree search: count %d\n", counts[t]);
|
||||
Collection<NeighbourNodeData> pointsWithinR =
|
||||
kdTree.findPointsWithinR(t, r, allowEqualToR);
|
||||
for (NeighbourNodeData nnData : pointsWithinR) {
|
||||
System.out.printf("kd: t=%d: time diff %d, norm %.3f\n",
|
||||
t, Math.abs(nnData.sampleIndex - t),
|
||||
normCalculator.norm(
|
||||
data[t], data[nnData.sampleIndex]));
|
||||
}
|
||||
// and compare to all pairs:
|
||||
for (int t2 = 0; t2 < samples; t2++) {
|
||||
double norm = Double.POSITIVE_INFINITY;
|
||||
if (Math.abs(t2 - t) > exclWindow) {
|
||||
norm = normCalculator.norm(
|
||||
data[t], data[t2]);
|
||||
}
|
||||
if ((!allowEqualToR && (norm < r)) ||
|
||||
( allowEqualToR && (norm <= r))) {
|
||||
System.out.printf("naive: t=%d: time diff %d, norm %.3f\n",
|
||||
t, Math.abs(t2 - t), norm);
|
||||
}
|
||||
}
|
||||
}
|
||||
assertEquals(naiveCount, counts[t]);
|
||||
}
|
||||
long endTimeValidate = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All neighbours within %.3f (average of %.3f) validated in: %.3f sec\n",
|
||||
r, (double) totalCount / (double) samples,
|
||||
((double) (endTimeValidate - endTimeCount)/1000.0));
|
||||
}
|
||||
|
||||
public void testCountNeighboursWithinRSeparateArraysWithDuplicates() throws Exception {
|
||||
int variables = 3;
|
||||
int dimensionsPerVariable = 3;
|
||||
|
|
|
|||
|
|
@ -100,14 +100,24 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
}
|
||||
|
||||
public void testRangeFinderStrictWithin() throws Exception {
|
||||
checkRangeFinder(true);
|
||||
checkRangeFinder(true, 0);
|
||||
}
|
||||
|
||||
public void testRangeFinderWithinOrEqual() throws Exception {
|
||||
checkRangeFinder(false);
|
||||
checkRangeFinder(false, 0);
|
||||
}
|
||||
|
||||
public void checkRangeFinder(boolean strict) throws Exception {
|
||||
public void testRangeFinderDynCorrExclusion() throws Exception {
|
||||
checkRangeFinder(true, 100);
|
||||
checkRangeFinder(false, 100);
|
||||
}
|
||||
|
||||
public void testRangeFinderArrayCallsWithDynCorrExclusion() throws Exception {
|
||||
checkRangeFinderArrayCalls(true, 100);
|
||||
checkRangeFinderArrayCalls(false, 100);
|
||||
}
|
||||
|
||||
public void checkRangeFinder(boolean strict, int exclWindow) throws Exception {
|
||||
int numTimeStepsInitial = 10000;
|
||||
int duplicateSteps = 100;
|
||||
int numTimeSteps = numTimeStepsInitial+duplicateSteps;
|
||||
|
|
@ -123,7 +133,7 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
double r = 0.2;
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
counts[t] = searcher.countPointsWithinR(t, r, !strict);
|
||||
counts[t] = searcher.countPointsWithinR(t, r, exclWindow, !strict);
|
||||
}
|
||||
long endTimeNNs = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("Found all neighbours within %.3f (mean = %.3f) in: %.3f sec\n",
|
||||
|
|
@ -135,7 +145,7 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
for (int t = 0; t < data.length; t++) {
|
||||
int count = 0;
|
||||
for (int t2 = 0; t2 < data.length; t2++) {
|
||||
if (t2 == t) {
|
||||
if (Math.abs(t2 - t) <= exclWindow) {
|
||||
continue;
|
||||
}
|
||||
double norm = Math.abs(data[t] - data[t2]);
|
||||
|
|
@ -157,6 +167,81 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
}
|
||||
}
|
||||
|
||||
public void checkRangeFinderArrayCalls(boolean strict, int exclWindow) 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];
|
||||
boolean[] withinR = new boolean[data.length];
|
||||
int[] indicesWithinR = new int[data.length];
|
||||
double r = 0.2;
|
||||
long startTime = Calendar.getInstance().getTimeInMillis();
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
searcher.findPointsWithinR(t, r, exclWindow, !strict,
|
||||
withinR, indicesWithinR);
|
||||
int count = 0;
|
||||
// Run through the list of points returned as within r and check them
|
||||
for (int t2 = 0; indicesWithinR[t2] != -1; t2++) {
|
||||
count++;
|
||||
assertTrue(withinR[indicesWithinR[t2]]);
|
||||
// Reset this marker
|
||||
withinR[indicesWithinR[t2]] = false;
|
||||
// Check that it's really within r
|
||||
double maxNorm = Math.abs(data[t] - data[indicesWithinR[t2]]);
|
||||
assertTrue((strict && (maxNorm < r)) ||
|
||||
( !strict && (maxNorm <= r)));
|
||||
}
|
||||
// Check that there were no other points marked as within r:
|
||||
int count2 = 0;
|
||||
for (int t2 = 0; t2 < data.length; t2++) {
|
||||
if (withinR[t2]) {
|
||||
count2++;
|
||||
}
|
||||
}
|
||||
assertEquals(0, count2); // We set all the ones that should have been true to false already
|
||||
// Now reset the boolean array
|
||||
counts[t] = count;
|
||||
}
|
||||
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 (Math.abs(t2 - t) <= exclWindow) {
|
||||
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);
|
||||
|
|
@ -172,12 +257,88 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
@SuppressWarnings("unchecked")
|
||||
PriorityQueue<NeighbourNodeData>[] pqs =
|
||||
(PriorityQueue<NeighbourNodeData>[]) new PriorityQueue[data.length];
|
||||
@SuppressWarnings("unchecked")
|
||||
PriorityQueue<NeighbourNodeData>[] pqsNoWindow =
|
||||
(PriorityQueue<NeighbourNodeData>[]) new PriorityQueue[data.length];
|
||||
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
PriorityQueue<NeighbourNodeData> nnPQ =
|
||||
searcher.findKNearestNeighbours(K, t);
|
||||
assertTrue(nnPQ.size() == K);
|
||||
pqs[t] = nnPQ;
|
||||
// Also check the result is still correct if we search with
|
||||
// zero size exclusion window:
|
||||
PriorityQueue<NeighbourNodeData> nnPQ_zeroWindow =
|
||||
searcher.findKNearestNeighbours(K, t, 0);
|
||||
assertTrue(nnPQ_zeroWindow.size() == K);
|
||||
pqsNoWindow[t] = nnPQ_zeroWindow;
|
||||
}
|
||||
long nnEndTime = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All %d nearest neighbours found in: %.3f sec\n",
|
||||
K, ((double) (nnEndTime - startTime)/1000.0));
|
||||
|
||||
// Now find the K nearest neighbours with a naive all-pairs comparison
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
double[][] distancesAndIndices = new double[data.length][2];
|
||||
for (int t2 = 0; t2 < data.length; t2++) {
|
||||
if (t2 != t) {
|
||||
distancesAndIndices[t2][0] =
|
||||
UnivariateNearestNeighbourSearcher.norm(data[t], data[t2],
|
||||
searcher.normTypeToUse);
|
||||
} else {
|
||||
distancesAndIndices[t2][0] = Double.POSITIVE_INFINITY;
|
||||
}
|
||||
distancesAndIndices[t2][1] = t2;
|
||||
}
|
||||
int[] timeStepsOfKthMins =
|
||||
MatrixUtils.kMinIndices(distancesAndIndices, 0, K);
|
||||
// Compare to what our NN technique picked out:
|
||||
PriorityQueue<NeighbourNodeData> nnPQ = pqs[t];
|
||||
PriorityQueue<NeighbourNodeData> nnPQNoWindow = pqsNoWindow[t];
|
||||
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);
|
||||
NeighbourNodeData nnDataNoWindow = nnPQNoWindow.poll();
|
||||
assertEquals(nnData.sampleIndex, nnDataNoWindow.sampleIndex);
|
||||
}
|
||||
}
|
||||
long endTimeValidate = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All %d nearest neighbours validated in: %.3f sec\n",
|
||||
K, ((double) (endTimeValidate - nnEndTime)/1000.0));
|
||||
}
|
||||
}
|
||||
|
||||
public void testFindKNearestNeighboursWithExclusionWindow() throws Exception {
|
||||
int dataLength = 1000;
|
||||
int exclWindow = 100;
|
||||
|
||||
for (int K = 1; K < 5; K++) {
|
||||
double[] data = rg.generateNormalData(dataLength, 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();
|
||||
@SuppressWarnings("unchecked")
|
||||
PriorityQueue<NeighbourNodeData>[] pqs =
|
||||
(PriorityQueue<NeighbourNodeData>[]) new PriorityQueue[data.length];
|
||||
|
||||
for (int t = 0; t < data.length; t++) {
|
||||
PriorityQueue<NeighbourNodeData> nnPQ =
|
||||
searcher.findKNearestNeighbours(K, t, exclWindow);
|
||||
assertTrue(nnPQ.size() == K);
|
||||
pqs[t] = nnPQ;
|
||||
}
|
||||
long nnEndTime = Calendar.getInstance().getTimeInMillis();
|
||||
System.out.printf("All %d nearest neighbours found in: %.3f sec\n",
|
||||
|
|
@ -187,7 +348,7 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
for (int t = 0; t < data.length; t++) {
|
||||
double[][] distancesAndIndices = new double[data.length][2];
|
||||
for (int t2 = 0; t2 < data.length; t2++) {
|
||||
if (t2 != t) {
|
||||
if (Math.abs(t2 - t) > exclWindow) {
|
||||
distancesAndIndices[t2][0] =
|
||||
UnivariateNearestNeighbourSearcher.norm(data[t], data[t2],
|
||||
searcher.normTypeToUse);
|
||||
|
|
@ -217,5 +378,4 @@ public class UnivariateNearestNeighbourTest extends TestCase {
|
|||
K, ((double) (endTimeValidate - nnEndTime)/1000.0));
|
||||
}
|
||||
}
|
||||
|
||||
}
|
||||
|
|
|
|||
Loading…
Reference in New Issue