[7.8][ML] Speed up the lat_long function

docs/CHANGELOG.asciidoc

@@ -32,6 +32,7 @@
 
 === Enhancements
 
+* Speedup anomaly detection for the lat_long function. (See {ml-pull}1102[#1102].)
 * Reduce CPU scheduling priority of native analysis processes to favor the ES JVM
   when CPU is constrained. (See {ml-pull}1109[#1109].)
 * Take `training_percent` into account when estimating memory usage for classification and regression. 

include/maths/CKMeans.h

@@ -7,6 +7,7 @@
 #ifndef INCLUDED_ml_maths_CKMeans_h
 #define INCLUDED_ml_maths_CKMeans_h
 
+#include <core/CContainerPrinter.h>
 #include <core/CLogger.h>
 
 #include <maths/CBasicStatistics.h>
@@ -23,6 +24,7 @@
 #include <algorithm>
 #include <cstddef>
 #include <cstdint>
+#include <limits>
 #include <sstream>
 #include <utility>
 #include <vector>
@@ -545,58 +547,41 @@ class CKMeansPlusPlusInitialization : private core::CNonCopyable {
     //! \param[out] result Filled in with the seed centres.
     template<typename ITR>
     void run(ITR beginPoints, ITR endPoints, std::size_t k, TPointVec& result) const {
+
         result.clear();
         if (beginPoints == endPoints || k == 0) {
             return;
         }
 
-        using TPointCRef = std::reference_wrapper<const POINT>;
-        using TPointCRefVec = std::vector<TPointCRef>;
-
-        std::size_t n = std::distance(beginPoints, endPoints);
+        std::size_t n(std::distance(beginPoints, endPoints));
         LOG_TRACE(<< "# points = " << n);
 
-        std::size_t select = CSampling::uniformSample(m_Rng, std::size_t(0), n);
-        LOG_TRACE(<< "select = " << select);
-
         result.reserve(k);
-        result.push_back(beginPoints[select]);
-        LOG_TRACE(<< "selected to date = " << core::CContainerPrinter::print(result));
-
-        TPointCRefVec selected{std::cref(result.back())};
-        selected.reserve(k);
-
-        CKdTree<TPointCRef> selectedLookup;
-        selectedLookup.reserve(k);
-
-        TDoubleVec distances(n, 0.0);
-        POINT distancesToHyperplanes{las::zero(result.back())};
+        result.push_back(beginPoints[CSampling::uniformSample(m_Rng, std::size_t(0), n)]);
 
+        m_Distances.assign(n, std::numeric_limits<double>::max());
         for (std::size_t i = 1; i < k; ++i) {
+            this->updateDistances(result.back(), beginPoints, endPoints);
+            m_Probabilities.assign(m_Distances.begin(), m_Distances.end());
+            result.push_back(beginPoints[CSampling::categoricalSample(m_Rng, m_Probabilities)]);
+        }
+        LOG_TRACE(<< "selected = " << core::CContainerPrinter::print(result));
+    }
 
-            selectedLookup.build(selected);
-
-            std::size_t j{0};
-            for (ITR point = beginPoints; point != endPoints; ++j, ++point) {
-                las::setZero(distancesToHyperplanes);
-                const auto* nn = selectedLookup.nearestNeighbour(*point, distancesToHyperplanes);
-                distances[j] = nn != nullptr
-                                   ? CTools::pow2(las::distance(*point, nn->get()))
-                                   : 0.0;
-            }
-
-            select = CSampling::categoricalSample(m_Rng, distances);
-            LOG_TRACE(<< "select = " << select);
-
-            result.push_back(beginPoints[select]);
-            selected.push_back(std::cref(result.back()));
-            LOG_TRACE(<< "selected to date = " << core::CContainerPrinter::print(result));
+private:
+    template<typename ITR>
+    void updateDistances(const POINT& selected, ITR beginPoints, ITR endPoints) const {
+        std::size_t j{0};
+        for (ITR point = beginPoints; point != endPoints; ++j, ++point) {
+            m_Distances[j] = std::min(
+                m_Distances[j], CTools::pow2(las::distance(*point, selected)));
         }
     }
 
 private:
-    //! The random number generator.
     RNG& m_Rng;
+    mutable TDoubleVec m_Distances;
+    mutable TDoubleVec m_Probabilities;
 };
 }
 }

include/maths/CKMeansOnline.h

@@ -126,11 +126,15 @@ class CKMeansOnline {
     }
 
     //! Construct with \p clusters.
-    CKMeansOnline(std::size_t k,
+    CKMeansOnline(TStoragePointMeanAccumulatorDoublePrVec& clusters,
+                  std::size_t k,
                   double decayRate,
-                  double minClusterSize,
-                  TStoragePointMeanAccumulatorDoublePrVec& clusters)
-        : CKMeansOnline{k, decayRate, minClusterSize} {
+                  double minClusterSize = MINIMUM_CATEGORY_COUNT,
+                  std::size_t bufferSize = BUFFER_SIZE,
+                  std::size_t numberSeeds = NUMBER_SEEDS,
+                  std::size_t maxIterations = MAX_ITERATIONS)
+        : CKMeansOnline{k,          decayRate,   minClusterSize,
+                        bufferSize, numberSeeds, maxIterations} {
         m_Clusters.swap(clusters);
         m_Clusters.reserve(m_K + m_BufferSize + 1);
     }
@@ -265,9 +269,9 @@ class CKMeansOnline {
         CBasicStatistics::SMin<double>::TAccumulator minCost;
         TSphericalClusterVec centres;
         TSphericalClusterVecVec candidates;
+        CKMeansPlusPlusInitialization<TSphericalCluster, RNG> seed{rng};
         for (std::size_t i = 0; i < numberSeeds; ++i) {
-            CKMeansPlusPlusInitialization<TSphericalCluster, RNG> seedCentres(rng);
-            seedCentres.run(kmeans.beginPoints(), kmeans.endPoints(), k, centres);
+            seed.run(kmeans.beginPoints(), kmeans.endPoints(), k, centres);
             kmeans.setCentres(centres);
             kmeans.run(maxIterations);
             kmeans.clusters(candidates);
@@ -309,7 +313,8 @@ class CKMeansOnline {
             for (std::size_t j = 0u; j < split[i].size(); ++j) {
                 clusters.push_back(m_Clusters[split[i][j]]);
             }
-            result.emplace_back(m_K, m_DecayRate, m_MinClusterSize, clusters);
+            result.emplace_back(clusters, m_K, m_DecayRate, m_MinClusterSize,
+                                m_BufferSize, m_NumberSeeds, m_MaxIterations);
         }
 
         return true;
@@ -561,11 +566,10 @@ class CKMeansOnline {
     //! \note We assume \p points is small so the bruteforce approach is fast.
     static void deduplicate(TStoragePointMeanAccumulatorDoublePrVec& clusters) {
         if (clusters.size() > 1) {
-            std::stable_sort(clusters.begin(), clusters.end(),
-                             [](const auto& lhs, const auto& rhs) {
-                                 return CBasicStatistics::mean(lhs.first) <
-                                        CBasicStatistics::mean(rhs.first);
-                             });
+            std::sort(clusters.begin(), clusters.end(), [](const auto& lhs, const auto& rhs) {
+                return CBasicStatistics::mean(lhs.first) <
+                       CBasicStatistics::mean(rhs.first);
+            });
             auto back = clusters.begin();
             for (auto i = back + 1; i != clusters.end(); ++i) {
                 if (CBasicStatistics::mean(back->first) == CBasicStatistics::mean(i->first)) {