Spatial anomaly

src/examples/hotgym/HelloSPTP.cpp

@@ -106,8 +106,9 @@ Real64 BenchmarkHotgym::run(UInt EPOCHS, bool useSPlocal, bool useSPglobal, bool
     //Encode
     tEnc.start();
     x+=0.01f; //step size for fn(x)
-    enc.encode(sin(x), input); //model sin(x) function //TODO replace with CSV data
-//    cout << x << "\n" << sin(x) << "\n" << input << "\n\n";
+    const Real value = sin(x);
+    enc.encode(value, input); //model sin(x) function //TODO replace with CSV data
+//    cout << x << "\n" << value << "\n" << input << "\n\n";
     tEnc.stop();
 
     tRng.start();
@@ -117,13 +118,13 @@ Real64 BenchmarkHotgym::run(UInt EPOCHS, bool useSPlocal, bool useSPglobal, bool
     //SP (global x local)
     if(useSPlocal) {
     tSPloc.start();
-    spLocal.compute(input, true, outSPlocal);
+    spLocal.compute(input, true, outSPlocal, value /* optional for spatial anomaly*/);
     tSPloc.stop();
     }
 
     if(useSPglobal) {
     tSPglob.start();
-    spGlobal.compute(input, true, outSPglobal);
+    spGlobal.compute(input, true, outSPglobal, value /* optional for spatial anomaly */);
     tSPglob.stop();
     }
     outSP = outSPglobal; //toggle if local/global SP is used further down the chain (TM, Anomaly)
@@ -168,6 +169,7 @@ Real64 BenchmarkHotgym::run(UInt EPOCHS, bool useSPlocal, bool useSPglobal, bool
       // output values
       cout << "Epoch = " << e << endl;
       cout << "Anomaly = " << an << endl;
+      cout << "Anomaly (spatial) = " << spGlobal.anomaly << endl;
       cout << "Anomaly (avg) = " << avgAnom10.getCurrentAvg() << endl;
       cout << "Anomaly (Likelihood) = " << anLikely << endl;
       cout << "SP (g)= " << outSP << endl;
@@ -217,10 +219,12 @@ Real64 BenchmarkHotgym::run(UInt EPOCHS, bool useSPlocal, bool useSPglobal, bool
         if(useSPglobal) { NTA_CHECK(outSPglobal == goldSP) << "Deterministic output of SP (g) failed!\n" << outSP << "should be:\n" << goldSP; }
         if(useSPlocal) {  NTA_CHECK(outSPlocal == goldSPlocal) << "Deterministic output of SP (l) failed!\n" << outSPlocal << "should be:\n" << goldSPlocal; }
         if(useTM) {       NTA_CHECK(outTM == goldTM) << "Deterministic output of TM failed!\n" << outTM << "should be:\n" << goldTM; }
-        NTA_CHECK(static_cast<UInt>(an *10000.0f) == static_cast<UInt>(goldAn *10000.0f)) //compare to 4 decimal places
+        // anomalies
+	NTA_CHECK(static_cast<UInt>(an *10000.0f) == static_cast<UInt>(goldAn *10000.0f)) //compare to 4 decimal places
 		               << "Deterministic output of Anomaly failed! " << an << "should be: " << goldAn;
 	NTA_CHECK(static_cast<UInt>(avgAnom10.getCurrentAvg() * 10000.0f) == static_cast<UInt>(goldAnAvg * 10000.0f)) 
 		<< "Deterministic average anom score failed:" << avgAnom10.getCurrentAvg() << " should be: " << goldAnAvg;
+	if(useSPglobal) { NTA_CHECK(0.0f == spGlobal.anomaly) << "Deterministic spatial anomaly mismatch!" << spGlobal.anomaly; }
       }
 
       // check runtime speed

src/htm/algorithms/SpatialPooler.cpp

@@ -459,7 +459,11 @@ void SpatialPooler::initialize(
 }
 
 
-void SpatialPooler::compute(const SDR &input, const bool learn, SDR &active) {
+void SpatialPooler::compute(const SDR &input, 
+		            const bool learn, 
+			    SDR &active,
+			    const Real spatialAnomalyInputValue) {
+
   input.reshape(  inputDimensions_ );
   active.reshape( columnDimensions_ );
   updateBookeepingVars_(learn);
@@ -485,6 +489,11 @@ void SpatialPooler::compute(const SDR &input, const bool learn, SDR &active) {
       updateMinDutyCycles_();
     }
   }
+
+  //update spatial anomaly
+  if(this->spAnomaly.enabled) {
+    spAnomaly.compute(spatialAnomalyInputValue);
+  }
 }
 
 

src/htm/algorithms/SpatialPooler.hpp

@@ -233,8 +233,15 @@ class SpatialPooler : public Serializable
   @param active An SDR representing the winning columns after
         inhibition. The size of the SDR is equal to the number of
         columns (also returned by the method getNumColumns).
+
+  @param  spatialAnomalyInputValue (optional) `Real` used for computing "spatial anomaly", see @ref `this.spatial_anomaly.compute()`,
+          obtained by @ref `SP.anomaly`
+
    */
-  virtual void compute(const SDR &input, const bool learn, SDR &active);
+  virtual void compute(const SDR &input, 
+		       const bool learn, 
+		       SDR &active, 
+		       const Real spatialAnomalyInputValue = std::numeric_limits<Real>::min());
 
 
   /**
@@ -1209,6 +1216,74 @@ class SpatialPooler : public Serializable
   Random rng_;
 
 public:
+  /**
+   *  holds together functionality for computing
+   *  `spatial anomaly`. This is used in NAB.
+   */
+  struct spatial_anomaly {
+
+    /** Fraction outside of the range of values seen so far that will be considered
+     * a spatial anomaly regardless of the anomaly likelihood calculation. 
+     * This accounts for the human labelling bias for spatial values larger than what
+     * has been seen so far.
+     */
+    Real SPATIAL_TOLERANCE = 0.05;
+
+    /**
+     * toggle if we should compute spatial anomaly
+     */
+    bool enabled = true;
+
+    /**
+     *  compute if the current input `value` is considered spatial-anomaly.
+     *  update internal variables.
+     *
+     *  @param value Real, input #TODO currently handles only 1 variable input, and requires value passed to compute! 
+     *  # later remove, and implement using SP's internal state. But for now we are compatible with NAB's implementation. 
+     *
+     * @return nothing, but updates internal variable `anomalyScore_`, which is either `NO_ANOMALY` (0.0f) , 
+     *   or `SPATIAL_ANOMALY` (exactly 0.9995947141f). Accessed by public @ref `SP.anomaly`.
+     *
+     */
+    void compute(const Real value) {
+      anomalyScore_ = NO_ANOMALY;
+      if(not enabled) return;
+      NTA_ASSERT(SPATIAL_TOLERANCE > 0.0f and SPATIAL_TOLERANCE <= 1.0f); 
+
+      if(minVal_ != std::numeric_limits<Real>::max() and
+	 maxVal_ != std::numeric_limits<Real>::min() ) {
+        const Real tolerance = (maxVal_ - minVal_) * SPATIAL_TOLERANCE;
+        const Real maxExpected = maxVal_ + tolerance;
+        const Real minExpected = minVal_ - tolerance;
+
+        if(value > maxExpected or value < minExpected) { //spatial anomaly
+          anomalyScore_ = SPATIAL_ANOMALY; 
+	}
+      }
+      if(value > maxVal_) maxVal_ = value;
+      if(value < minVal_) minVal_ = value;
+      NTA_ASSERT(anomalyScore_ == NO_ANOMALY or anomalyScore_ == SPATIAL_ANOMALY) << "Out of bounds of acceptable values";
+    }
+
+    private:
+      Real minVal_ = std::numeric_limits<Real>::max(); //TODO fix serialization
+      Real maxVal_ = std::numeric_limits<Real>::min();
+
+    public:
+      const Real NO_ANOMALY = 0.0f;
+      const Real SPATIAL_ANOMALY = 0.9995947141f; //almost 1.0 = max anomaly. Encodes value specific to spatial anomaly (so this can be recognized on results),
+        // "5947141" would translate in l33t speech to "spatial" :)
+      Real anomalyScore_ = NO_ANOMALY; //default score = no anomaly
+  } spAnomaly;
+
+  /**
+   *  spatial anomaly
+   *
+   *  updated on each @ref `compute()`. 
+   *
+   *  @return either 0.0f (no anomaly), or exactly 0.9995947141f (spatial anomaly). This specific value can be recognized in results. 
+   */
+  const Real& anomaly = spAnomaly.anomalyScore_;
   const Connections &connections = connections_;
 };
 

src/test/unit/algorithms/SpatialPoolerTest.cpp

@@ -32,6 +32,7 @@
 #include <htm/types/Types.hpp>
 #include <htm/utils/Log.hpp>
 #include <htm/os/Timer.hpp>
+#include <htm/encoders/ScalarEncoder.hpp>
 
 namespace testing {
 
@@ -2103,5 +2104,68 @@ TEST(SpatialPoolerTest, ExactOutput) {
   ASSERT_EQ( columns, gold_sdr );
 }
 
+TEST(SpatialPoolerTest, spatialAnomaly) {
+  SDR inputs({ 1000 });
+  SDR columns({ 200 });
+  SpatialPooler sp({inputs.dimensions}, {columns.dimensions},
+    /*potentialRadius*/ 99999,
+    /*potentialPct*/ 0.5f,
+    /*globalInhibition*/ true,
+    /*localAreaDensity*/ 0.05f,
+    /*stimulusThreshold*/ 3u,
+    /*synPermInactiveDec*/ 0.008f,
+    /*synPermActiveInc*/ 0.05f,
+    /*synPermConnected*/ 0.1f,
+    /*minPctOverlapDutyCycles*/ 0.001f,
+    /*dutyCyclePeriod*/ 200,
+    /*boostStrength*/ 10.0f,
+    /*seed*/ 42,
+    /*spVerbosity*/ 0,
+    /*wrapAround*/ true);
+
+
+  // test too large threshold
+#ifdef NTA_ASSERTIONS_ON  //only for Debug
+  sp.spAnomaly.SPATIAL_TOLERANCE = 1.2345f; //out of bounds, will crash!
+  EXPECT_ANY_THROW(sp.compute(inputs, false, columns, 0.1f /*whatever, fails on TOLERANCE */)) << "Spatial anomaly should fail if SPATIAL_TOLERANCE is out of bounds!";
+  sp.spAnomaly.SPATIAL_TOLERANCE = 0.01f; //within bounds, OK
+  EXPECT_NO_THROW(sp.compute(inputs, false, columns, 0.1f /*whatever */));
+#endif
+
+  //test spatial anomaly computation
+  sp.spAnomaly.SPATIAL_TOLERANCE = 0.2f; //threshold 20%
+  ScalarEncoderParameters params;
+  params.minimum = 0.0f;
+  params.maximum = 100.0f;
+  params.size = 1000;
+  params.sparsity = 0.3f;
+  ScalarEncoder enc(params);
+
+  Real val;
+
+  val = 0.0f;
+  enc.encode(val, inputs); //TODO can SDR hold .origValue = Real which encoders would set? Classifier,Predictor, and spatia_anomaly would use that
+  sp.compute(inputs, true, columns, val);
+  EXPECT_EQ(0.0f, sp.anomaly);
+  EXPECT_EQ(sp.spAnomaly.NO_ANOMALY, sp.anomaly) << "should be the same as above";
+
+  val = 10.0f;
+  enc.encode(val, inputs); 
+  sp.compute(inputs, true, columns, val);
+  EXPECT_EQ(sp.spAnomaly.NO_ANOMALY, sp.anomaly);
+
+  val = 11.99f; //(10-0) * 0.2 == 2 -> <-2, +12> is not anomalous
+  enc.encode(val, inputs); 
+  sp.compute(inputs, true, columns, val);
+  EXPECT_EQ(0.0f, sp.anomaly) << "This shouldn't be an anomaly!";
+
+  val = 100.0f; //(12-0) * 0.2 == ~2.2 -> <-2.2, +14.2> is not anomalous, but 100 is!
+  enc.encode(val, inputs);
+  sp.compute(inputs, true, columns, val);
+  EXPECT_EQ(0.9995947141f, sp.anomaly) << "This should be an anomaly!";
+  EXPECT_EQ(sp.spAnomaly.SPATIAL_ANOMALY, sp.anomaly) << "Should be same as above!";
+
+
+}
 
 } // end anonymous namespace