fixed some bugs in jensen shannon and kullback leibler

main.cpp

@@ -30,7 +30,7 @@ int main(int argc, char** argv) {
 
 
 	//::testing::GTEST_FLAG(filter) = "*Offline.readWrite*";
-	::testing::GTEST_FLAG(filter) = "*Angle*";
+    ::testing::GTEST_FLAG(filter) = "*Jensen*";
 
 
     //::testing::GTEST_FLAG(filter) = "*Barometer*";

math/Random.h

@@ -18,10 +18,10 @@ class RandomGenerator : public std::minstd_rand {
 public:
 
 	/** ctor with default seed */
-	RandomGenerator() : std::minstd_rand(RANDOM_SEED) {;}
+    RandomGenerator() : std::minstd_rand(RANDOM_SEED) {;}
 
 	/** ctor with custom seed */
-	RandomGenerator(result_type seed) : std::minstd_rand(seed) {;}
+    RandomGenerator(result_type seed) : std::minstd_rand(seed) {;}
 
 };
 

math/divergence/JensenShannon.h

@@ -18,6 +18,13 @@ public:
      * @param Q is a vector containg the densities of the same samples set then P
     */
     static inline Scalar getGeneralFromSamples(Eigen::VectorXd P, Eigen::VectorXd Q, LOGMODE mode){
+        // normalize
+        P /= P.sum();
+        Q /= Q.sum();
+
+        Assert::isNear((double)P.sum(), 1.0, 0.01,"Normalization failed.. this shouldn't happen");
+        Assert::isNear((double)Q.sum(), 1.0, 0.01, "Normalization failed.. this shouldn't happen");
+
         Eigen::VectorXd M = 0.5 * (P + Q);
 
         return (0.5 * KullbackLeibler<Scalar>::getGeneralFromSamples(P, M, mode)) + (0.5 * KullbackLeibler<Scalar>::getGeneralFromSamples(Q, M, mode));

math/divergence/KullbackLeibler.h

@@ -56,8 +56,8 @@ namespace Divergence {
 
                         // calc PQratio
                         if(Q[i] == 0.0){
-                            Assert::doThrow("Division by zero is not allowed ;).");
-                            //PQratio = P[i] / 0.00001;
+                            //Assert::doThrow("Division by zero is not allowed ;).");
+                            PQratio = P[i] / 0.00001;
                         } else {
                             PQratio = P[i] / Q[i];
                         }

sensors/radio/VAPGrouper.h

@@ -89,7 +89,7 @@ public:
 
 		}
 
-		Log::add(name, "grouped " + std::to_string(original.entries.size()) + " measurements into " + std::to_string(result.entries.size()), true);
+        //Log::add(name, "grouped " + std::to_string(original.entries.size()) + " measurements into " + std::to_string(result.entries.size()), true);
 
 
 		// done

tests/math/divergence/TestJensenShannon.cpp

@@ -0,0 +1,94 @@
+#ifdef WITH_TESTS
+
+#include "../../Tests.h"
+#include "../../../math/divergence/JensenShannon.h"
+#include "../../../math/Distributions.h"
+
+#include <random>
+
+
+TEST(JensenShannon, generalFromSamples) {
+    //ge cov
+    Eigen::VectorXd mu1(2);
+    mu1[0] = 1.0;
+    mu1[1] = 1.0;
+
+    Eigen::VectorXd mu2(2);
+    mu2[0] = 1.0;
+    mu2[1] = 1.0;
+
+    Eigen::VectorXd mu3(2);
+    mu3[0] = 1.0;
+    mu3[1] = 1.0;
+
+    Eigen::VectorXd mu4(2);
+    mu4[0] = 9.0;
+    mu4[1] = 9.0;
+
+    Eigen::MatrixXd cov1(2,2);
+    cov1(0,0) = 1.0;
+    cov1(0,1) = 0.0;
+    cov1(1,0) = 0.0;
+    cov1(1,1) = 1.0;
+
+    Eigen::MatrixXd cov2(2,2);
+    cov2(0,0) = 1.0;
+    cov2(0,1) = 0.0;
+    cov2(1,0) = 0.0;
+    cov2(1,1) = 1.0;
+
+    Eigen::MatrixXd cov3(2,2);
+    cov3(0,0) = 1.0;
+    cov3(0,1) = 0.0;
+    cov3(1,0) = 0.0;
+    cov3(1,1) = 1.0;
+
+    Eigen::MatrixXd cov4(2,2);
+    cov4(0,0) = 3.0;
+    cov4(0,1) = 0.0;
+    cov4(1,0) = 0.0;
+    cov4(1,1) = 13.0;
+
+    Distribution::NormalDistributionN norm1(mu1, cov1);
+    Distribution::NormalDistributionN norm2(mu2, cov2);
+    Distribution::NormalDistributionN norm3(mu3, cov3);
+    Distribution::NormalDistributionN norm4(mu4, cov4);
+
+    int size = 10000;
+    Eigen::VectorXd samples1(size);
+    Eigen::VectorXd samples2(size);
+    Eigen::VectorXd samples3(size);
+    Eigen::VectorXd samples4(size);
+
+    //random numbers
+    std::mt19937_64 rng;
+   // initialize the random number generator with time-dependent seed
+   uint64_t timeSeed = std::chrono::high_resolution_clock::now().time_since_epoch().count();
+   std::seed_seq ss{uint32_t(timeSeed & 0xffffffff), uint32_t(timeSeed>>32)};
+   rng.seed(ss);
+   // initialize a uniform distribution between 0 and 1
+   std::uniform_real_distribution<double> unif(-9, 10);
+
+    //generate samples
+    for(int i = 0; i < size; ++i){
+
+        double r1 = unif(rng);
+        double r2 = unif(rng);
+        Eigen::VectorXd v(2);
+        v << r1, r2;
+
+        samples1[i] = norm1.getProbability(v);
+        samples2[i] = norm2.getProbability(v);
+        samples3[i] = norm3.getProbability(v);
+        samples4[i] = norm4.getProbability(v);
+    }
+
+    double kld12 = Divergence::JensenShannon<float>::getGeneralFromSamples(samples1, samples2, Divergence::LOGMODE::NATURALIS);
+    double kld34 = Divergence::JensenShannon<float>::getGeneralFromSamples(samples3, samples4, Divergence::LOGMODE::NATURALIS);
+    std::cout << kld34 << " > " << kld12 << std::endl;
+
+    ASSERT_GE(kld34, kld12);
+
+}
+
+#endif