ungrading thresholding and fixing RCFCaster errors for initial values (#385)

* ungrading thresholding and fixing RCFCaster errors for initial values

* optimizations

* tweaks and comments

* changes

* cleanup and optimizations

* fixes based on comments

* fixes and changes

* off by one correction

Author: sudiptoguha

Java/core/src/main/java/com/amazon/randomcutforest/state/Version.java

@@ -20,4 +20,5 @@ public class Version {
     public static final String V2_1 = "2.1";
     public static final String V3_0 = "3.0";
     public static final String V3_5 = "3.5";
+    public static final String V3_7 = "3.7";
 }

Java/examples/src/main/java/com/amazon/randomcutforest/examples/parkservices/LowNoisePeriodic.java

@@ -44,7 +44,7 @@ public String description() {
     public void run() throws Exception {
         // Create and populate a random cut forest
 
-        int shingleSize = 4;
+        int shingleSize = 8;
         int numberOfTrees = 50;
         int sampleSize = 256;
         int dataSize = 100000;
@@ -53,17 +53,21 @@ public void run() throws Exception {
         double[] reference = new double[] { 1.0f, 3.0f, 5.0f, 7.0f, 9.0f, 11.0f, 9.5f, 8.5f, 7.5f, 6.5f, 6.0f, 6.5f,
                 7.0f, 7.5f, 9.5f, 11.0f, 12.5f, 10.5f, 8.5f, 7.0f, 5.0f, 3.0f, 2.0f, 1.0f };
 
-        // change this to control the percent deviation
-        // NOTE that if the noise is smaller than 0.003 times the actual value then
-        // it would be difficult to detect the anomalies unless the slope is 0
-
-        double noise = 2.0;
+        // the noise should leave suffient gap between the consecutive levels
+        double noise = 0.25;
+        // the noise will be amplified by something within [factorRange, 2*factorRange]
+        // increase should lead to increased precision--recall; likewise decrease must
+        // also
+        // lead to decreased precision recall; if the factor = 1, then the anomalies are
+        // information theoretically almost non-existent
+        double anomalyFactor = 10;
 
         double slope = 0.2 * sampleSize
                 * (Arrays.stream(reference).max().getAsDouble() - Arrays.stream(reference).min().getAsDouble()) / 50000;
 
-        // to analyse without linear shift
-        // slope = 0;
+        // to analyse without linear shift; comment out the line below and change the
+        // slope above as desired
+        slope = 0;
 
         double anomalyRate = 0.005;
         long seed = new Random().nextLong();
@@ -74,12 +78,14 @@ public void run() throws Exception {
         int correct = 0;
         int late = 0;
 
-        // change the transformation below to experiment
+        // change the transformation below to experiment;
+        // if slope != 0 then NONE will have poor result
+        // both of the difference operations also introduce many errors
         TransformMethod method = TransformMethod.NORMALIZE;
 
         int dimensions = shingleSize;
-        ThresholdedRandomCutForest forest = ThresholdedRandomCutForest.builder().compact(true).dimensions(dimensions)
-                .randomSeed(0).numberOfTrees(numberOfTrees).shingleSize(shingleSize).sampleSize(sampleSize)
+        ThresholdedRandomCutForest forest = ThresholdedRandomCutForest.builder().dimensions(dimensions).randomSeed(0)
+                .numberOfTrees(numberOfTrees).shingleSize(shingleSize).sampleSize(sampleSize)
                 .internalShinglingEnabled(true).anomalyRate(0.01).forestMode(ForestMode.STANDARD).startNormalization(32)
                 .transformMethod(method).outputAfter(32).initialAcceptFraction(0.125)
                 // for 1D data weights should not alter results significantly (if in reasonable
@@ -99,7 +105,7 @@ public void run() throws Exception {
         // missed current value 3.0 (say X), intended 1.0 (equiv., X - noise), because
         // the shift up in the actual was not 2*noise
 
-        // forest.setIgnoreNearExpectedFromAbove( new double [] {2 * noise});
+        // forest.setIgnoreNearExpectedFromAbove( new double [] {2*noise});
 
         // or to suppress all anomalies that are shifted up from predicted
         // for any sequence; using Double.MAX_VALUE may cause overflow
@@ -110,7 +116,7 @@ public void run() throws Exception {
         // the shift down in the actual was not 2*noise, in effect we suppress all
         // anomalies
 
-        // forest.setIgnoreNearExpectedFromBelow(new double [] {2*noise});
+        // forest.setIgnoreNearExpectedFromBelow(new double [] {noise*2});
 
         // the following suppresses all anomalies that shifted down compared to
         // predicted
@@ -126,12 +132,15 @@ public void run() throws Exception {
             boolean anomaly = false;
 
             double intendedValue = reference[(count + 4) % reference.length] + slope * count;
-            // extremely periodic signal
+            // extremely periodic signal -- note that there is no periodicity detection
             value[0] = intendedValue;
             if (rng.nextDouble() < anomalyRate && count > initialSegment) {
-                value[0] += (rng.nextDouble() < 0.5) ? -noise : noise;
+                double anomalyValue = noise * anomalyFactor * (1 + rng.nextDouble());
+                value[0] += (rng.nextDouble() < 0.5) ? -anomalyValue : anomalyValue;
                 anomaly = true;
                 ++numAnomalies;
+            } else {
+                value[0] += (2 * rng.nextDouble() - 1) * noise;
             }
 
             AnomalyDescriptor result = forest.process(new double[] { value[0] }, 0);

Java/examples/src/main/java/com/amazon/randomcutforest/examples/parkservices/ThresholdedMultiDimensionalExample.java

@@ -47,15 +47,15 @@ public String description() {
     public void run() throws Exception {
         // Create and populate a random cut forest
 
-        int shingleSize = 4;
+        int shingleSize = 8;
         int numberOfTrees = 50;
         int sampleSize = 256;
         Precision precision = Precision.FLOAT_32;
         int dataSize = 4 * sampleSize;
 
         // change this to try different number of attributes,
         // this parameter is not expected to be larger than 5 for this example
-        int baseDimensions = 2;
+        int baseDimensions = 3;
 
         int dimensions = baseDimensions * shingleSize;
         ThresholdedRandomCutForest forest = ThresholdedRandomCutForest.builder().compact(true).dimensions(dimensions)
@@ -66,19 +66,33 @@ public void run() throws Exception {
 
         long seed = new Random().nextLong();
         System.out.println("seed = " + seed);
+
+        // basic amplitude of the waves -- the parameter will be randomly scaled up
+        // betwee 0-20 percent
+        double amplitude = 100.0;
+
+        // the amplitude of random noise it will be +ve/-ve uniformly at random
+        double noise = 5.0;
+
+        // the following controls the ratio of anomaly magnitude to noise
+        // notice amplitude/noise would determine signal-to-noise ratio
+        double anomalyFactor = 5;
+
+        // the following determines if a random linear trend should be added
+        boolean useSlope = false;
+
         // change the last argument seed for a different run
         MultiDimDataWithKey dataWithKeys = ShingledMultiDimDataWithKeys.getMultiDimData(dataSize + shingleSize - 1, 50,
-                100, 5, seed, baseDimensions);
+                amplitude, noise, seed, baseDimensions, anomalyFactor, useSlope);
         int keyCounter = 0;
         int count = 0;
         for (double[] point : dataWithKeys.data) {
 
             AnomalyDescriptor result = forest.process(point, 0L);
 
-            if (keyCounter < dataWithKeys.changeIndices.length
-                    && count + shingleSize - 1 == dataWithKeys.changeIndices[keyCounter]) {
-                System.out.println("timestamp " + (count + shingleSize - 1) + " CHANGE "
-                        + Arrays.toString(dataWithKeys.changes[keyCounter]));
+            if (keyCounter < dataWithKeys.changeIndices.length && count == dataWithKeys.changeIndices[keyCounter]) {
+                System.out.println(
+                        "timestamp " + (count) + " CHANGE " + Arrays.toString(dataWithKeys.changes[keyCounter]));
                 ++keyCounter;
             }
 

Java/parkservices/src/main/java/com/amazon/randomcutforest/parkservices/ErrorHandler.java

@@ -148,7 +148,7 @@ public ErrorHandler(int errorHorizon, int forecastHorizon, int sequenceIndex, do
                 pastForecasts[i] = new RangeVector(values, upper, lower);
                 System.arraycopy(actualsFlattened, i * inputLength, actuals[i], 0, inputLength);
             }
-            calibrate();
+            calibrate(null);
         }
     }
 
@@ -171,7 +171,7 @@ public void update(ForecastDescriptor descriptor, Calibration calibrationMethod)
             actuals[errorIndex][i] = (float) input[i];
         }
         ++sequenceIndex;
-        calibrate();
+        calibrate(descriptor.deviations);
         if (calibrationMethod != Calibration.NONE) {
             if (calibrationMethod == Calibration.SIMPLE) {
                 adjust(descriptor.timedForecast.rangeVector, errorDistribution);
@@ -260,8 +260,8 @@ public RangeVector computeErrorPercentile(double percentile, int newHorizon,
                         double fracRank = percentile * len;
                         Arrays.sort(copy);
                         values[pos] = interpolatedMedian(copy);
-                        lower[pos] = interpolatedLowerRank(copy, fracRank);
-                        upper[pos] = interpolatedUpperRank(copy, len, fracRank);
+                        lower[pos] = interpolatedLowerRank(copy, fracRank, 0);
+                        upper[pos] = interpolatedUpperRank(copy, len, fracRank, 0);
                     }
                 }
             }
@@ -286,8 +286,11 @@ protected double[] getErrorVector(int len, int leadtime, int inputCoordinate, in
      * this method computes a lot of different quantities, some of which would be
      * useful in the future.In particular it splits the RMSE into positive and
      * negative contribution which is informative about directionality of error.
+     *
+     *
+     * @param errorDeviations the weighted standard deviations seen so far
      */
-    protected void calibrate() {
+    protected void calibrate(double[] errorDeviations) {
         int inputLength = actuals[0].length;
         int arrayLength = pastForecasts.length;
         int errorIndex = (sequenceIndex - 1 + arrayLength) % arrayLength;
@@ -327,16 +330,19 @@ protected void calibrate() {
                     errorRMSE.high[pos] = (positiveCount > 0) ? Math.sqrt(positiveSqSum / positiveCount) : 0;
                     errorRMSE.low[pos] = (positiveCount < len) ? -Math.sqrt(negativeSqSum / (len - positiveCount)) : 0;
                     Arrays.sort(medianError, 0, len);
+                    // medianError array is now sorted
                     errorDistribution.values[pos] = interpolatedMedian(medianError);
-                    errorDistribution.upper[pos] = interpolatedUpperRank(medianError, len, len * percentile);
-                    errorDistribution.lower[pos] = interpolatedLowerRank(medianError, len * percentile);
+                    double deviation = (errorDeviations == null) ? 0 : errorDeviations[j];
+                    errorDistribution.upper[pos] = interpolatedUpperRank(medianError, len, len * percentile, deviation);
+                    errorDistribution.lower[pos] = interpolatedLowerRank(medianError, len * percentile, deviation);
                     intervalPrecision[pos] = intervalPrecision[pos] / len;
                 } else {
                     errorMean[pos] = 0;
                     errorRMSE.high[pos] = errorRMSE.low[pos] = 0;
                     errorDistribution.values[pos] = 0;
-                    errorDistribution.upper[pos] = Float.MAX_VALUE;
-                    errorDistribution.lower[pos] = -Float.MAX_VALUE;
+                    double deviation = (errorDeviations == null) ? 0 : errorDeviations[j];
+                    errorDistribution.upper[pos] = (float) (1.3 * deviation);
+                    errorDistribution.lower[pos] = -(float) (1.3 * deviation);
                     adders.upper[pos] = adders.lower[pos] = adders.values[pos] = 0;
                     intervalPrecision[pos] = 0;
                 }
@@ -354,28 +360,30 @@ float interpolatedMedian(double[] array) {
         }
     }
 
-    float interpolatedLowerRank(double[] array, double fracRank) {
+    float interpolatedLowerRank(double[] ascendingArray, double fracRank, double deviation) {
         if (fracRank < 1) {
-            return -Float.MAX_VALUE;
+            return (float) (-1.3 * deviation * (1 - fracRank) + fracRank * ascendingArray[0]);
         }
         int rank = (int) Math.floor(fracRank);
         if (!RCFCaster.USE_INTERPOLATION_IN_DISTRIBUTION) {
             // turn off interpolation
             fracRank = rank;
         }
-        return (float) (array[rank - 1] + (fracRank - rank) * (array[rank] - array[rank - 1]));
+        return (float) (ascendingArray[rank - 1]
+                + (fracRank - rank) * (ascendingArray[rank] - ascendingArray[rank - 1]));
     }
 
-    float interpolatedUpperRank(double[] array, int len, double fracRank) {
+    float interpolatedUpperRank(double[] ascendingArray, int len, double fracRank, double deviation) {
         if (fracRank < 1) {
-            return Float.MAX_VALUE;
+            return (float) (1.3 * deviation * (1 - fracRank) + fracRank * ascendingArray[len - 1]);
         }
         int rank = (int) Math.floor(fracRank);
         if (!RCFCaster.USE_INTERPOLATION_IN_DISTRIBUTION) {
             // turn off interpolation
             fracRank = rank;
         }
-        return (float) (array[len - rank] + (fracRank - rank) * (array[len - rank - 1] - array[len - rank]));
+        return (float) (ascendingArray[len - rank]
+                + (fracRank - rank) * (ascendingArray[len - rank - 1] - ascendingArray[len - rank]));
     }
 
     void adjust(RangeVector rangeVector, RangeVector other) {