Added kNN_FEDSA_Classifier + tests.

Signed-off-by: Daniel Scheffler <danschef@gfz-potsdam.de>

gms_preprocessing/algorithms/classification.py

@@ -428,23 +428,17 @@ class FEDSA_Classifier(_ImageClassifier):
     def fedsa(self):
         return self._distance_metrics
 
-    def _predict(self, tilepos):
-        assert global_shared_endmembers is not None and global_shared_im2classify is not None
-
-        (rS, rE), (cS, cE) = tilepos
-        tileimdata = global_shared_im2classify[rS: rE + 1, cS: cE + 1, :]
-        endmembers = global_shared_endmembers
-
-        if not tileimdata.shape[2] == self.train_spectra.shape[1]:
+    def _calc_fedsa(self, image, endmembers):
+        if not image.shape[2] == self.train_spectra.shape[1]:
             raise RuntimeError('Matrix dimensions are not aligned. Input image has %d bands but input spectra '
-                               'have %d.' % (tileimdata.shape[2], self.train_spectra.shape[1]))
+                               'have %d.' % (image.shape[2], self.train_spectra.shape[1]))
 
         # normalize input data because SAM asserts only data between -1 and 1
-        train_spectra_norm, tileimdata_norm = normalize_endmembers_image(endmembers, tileimdata)
+        train_spectra_norm, tileimdata_norm = normalize_endmembers_image(endmembers, image)
 
-        angles = np.zeros((tileimdata.shape[0], tileimdata.shape[1], self.n_samples), np.float)
-        ed = np.zeros((tileimdata.shape[0], tileimdata.shape[1], self.n_samples), np.float)
-        tileimspectra = im2spectra(tileimdata)
+        angles = np.zeros((image.shape[0], image.shape[1], self.n_samples), np.float)
+        ed = np.zeros((image.shape[0], image.shape[1], self.n_samples), np.float)
+        tileimspectra = im2spectra(image)
         # if np.std(tileimdata) == 0:  # skip tiles that only contain the same value
 
         # loop over all training spectra and compute spectral angle for each pixel
@@ -453,12 +447,22 @@ class FEDSA_Classifier(_ImageClassifier):
             angles[:, :, n_sample] = calc_sam(tileimdata_norm, train_spectrum, axis=2)
             ed[:, :, n_sample] = np.sqrt(np.sum((tileimspectra.astype(np.float) -
                                                  train_spectrum.flatten().astype(np.float)) ** 2, axis=1))\
-                .reshape(tileimdata.shape[:2])
+                .reshape(image.shape[:2])
 
         angles_norm = angles / angles.max()
         ed_norm = ed / ed.max()
 
         fedsa = (angles_norm + ed_norm) / 2
+
+        return fedsa
+
+    def _predict(self, tilepos):
+        assert global_shared_endmembers is not None and global_shared_im2classify is not None
+
+        (rS, rE), (cS, cE) = tilepos
+        tileimdata = global_shared_im2classify[rS: rE + 1, cS: cE + 1, :]
+
+        fedsa = self._calc_fedsa(tileimdata, global_shared_endmembers)
         fedsa_min = np.min(fedsa, axis=2).astype(np.float32)
         cmap = np.argmin(fedsa, axis=2).astype(np.int16)
 
@@ -481,6 +485,44 @@ class FEDSA_Classifier(_ImageClassifier):
         self._show_distance_metrics(**kwargs)
 
 
+class kNN_FEDSA_Classifier(FEDSA_Classifier):
+    def __init__(self, train_spectra, n_neighbors=3, CPUs=1):
+        # type: (np.ndarray, int, Union[int, None]) -> None
+        super(kNN_FEDSA_Classifier, self).__init__(train_spectra, CPUs=CPUs)
+
+        self.clf_name = 'k-nearest neighbour spectral angle mapper (kNN_SAM; k=%d)' % n_neighbors
+        self.n_neighbors = n_neighbors
+
+    def _predict(self, tilepos):
+        assert global_shared_endmembers is not None and global_shared_im2classify is not None
+
+        (rS, rE), (cS, cE) = tilepos
+        tileimdata = global_shared_im2classify[rS: rE + 1, cS: cE + 1, :]
+
+        fedsa = self._calc_fedsa(tileimdata, global_shared_endmembers)
+        k = self.n_neighbors if self.n_neighbors <= fedsa.shape[2] else fedsa.shape[2]
+
+        if self.n_neighbors < fedsa.shape[2]:
+            cmap = np.argpartition(fedsa, k, axis=2)[:, :, :k].astype(np.int16)
+            fedsa_min_k = np.partition(fedsa, k, axis=2)[:, :, :k].astype(np.float32)
+
+            # sort cmap by ascending angles
+            idx_2D = np.argsort(fedsa_min_k, axis=2).reshape(-1, cmap.shape[2])
+            cmap = \
+                cmap.reshape(-1, cmap.shape[2])[np.arange(cmap.shape[0] * cmap.shape[1])[:, np.newaxis], idx_2D]\
+                .reshape(*cmap.shape)
+            fedsa_min_k = np.sort(fedsa_min_k, axis=2)
+        else:
+            cmap = np.tile(np.arange(fedsa.shape[2]).reshape(1, 1, -1), (*fedsa.shape[:2], 1))
+            fedsa_min_k = fedsa
+
+        if global_shared_im2classify.nodata is not None and self._cmap_nodataVal is not None:
+            cmap = self.overwrite_cmap_at_nodata_positions(cmap, tileimdata,
+                                                           self._cmap_nodataVal, global_shared_im2classify.nodata)
+
+        return tilepos, cmap.astype(np.int16), fedsa_min_k
+
+
 class SID_Classifier(_ImageClassifier):
     def __init__(self, train_spectra, CPUs=1):
         # type: (np.ndarray, Union[int, None]) -> None
@@ -633,6 +675,11 @@ def classify_image(image, train_spectra, train_labels, classif_alg, in_nodataVal
             train_spectra,
             CPUs=CPUs)
 
+    elif classif_alg == 'kNN_FEDSA':
+        clf = kNN_FEDSA_Classifier(
+            train_spectra,
+            CPUs=CPUs)
+
     elif classif_alg == 'SID':
         clf = SID_Classifier(
             train_spectra,
@@ -646,15 +693,15 @@ def classify_image(image, train_spectra, train_labels, classif_alg, in_nodataVal
 
     else:
         raise NotImplementedError("Currently only the methods 'kNN', 'MinDist', 'kNN_MinDist', 'SAM', 'kNN_SAM', "
-                                  "'FEDSA', 'SID' and 'RF' are implemented.")
+                                  "'FEDSA', 'kNN_FEDSA', 'SID' and 'RF' are implemented.")
 
     cmap = clf.classify(image, in_nodataVal=in_nodataVal, cmap_nodataVal=cmap_nodataVal, tiledims=tiledims)
 
     # label unclassified pixels
     if unclassified_threshold is not None:
-        if classif_alg not in ['MinDist', 'kNN_MinDist', 'SAM', 'kNN_SAM', 'FEDSA', 'SID']:
-            raise RuntimeError("Only the methods 'MinDist', 'kNN_MinDist', 'SAM', 'kNN_SAM', 'FEDSA' and 'SID' "
-                               "can label unclassifed pixels.")
+        if classif_alg not in ['MinDist', 'kNN_MinDist', 'SAM', 'kNN_SAM', 'FEDSA', 'kNN_FEDSA', 'SID']:
+            raise RuntimeError("Only the methods 'MinDist', 'kNN_MinDist', 'SAM', 'kNN_SAM', 'FEDSA', 'kNN_FEDSA' and "
+                               "'SID' can label unclassifed pixels.")
 
         clf.label_unclassified_pixels(label_unclassified=unclassified_pixVal, threshold=unclassified_threshold)
 
@@ -666,7 +713,7 @@ def classify_image(image, train_spectra, train_labels, classif_alg, in_nodataVal
             dist = clf.euclidian_distance
         elif classif_alg in ['SAM', 'kNN_SAM']:
             dist = clf.angles_deg
-        elif classif_alg == 'FEDSA':
+        elif classif_alg in ['FEDSA', 'kNN_FEDSA']:
             dist = clf.fedsa
         elif classif_alg == 'SID':
             dist = clf.sid

tests/test_image_classifier.py

@@ -20,7 +20,7 @@ from time import time
 from gms_preprocessing import set_config
 from gms_preprocessing.algorithms.classification import \
     MinimumDistance_Classifier, kNN_MinimumDistance_Classifier, kNN_Classifier, SAM_Classifier, kNN_SAM_Classifier, \
-    FEDSA_Classifier, SID_Classifier, RF_Classifier
+    FEDSA_Classifier, kNN_FEDSA_Classifier, SID_Classifier, RF_Classifier
 
 from . import db_host
 
@@ -101,7 +101,7 @@ class Test_kNN_MinimumDistance_Classifier(unittest.TestCase):
     def test_label_unclassified_pixels_absolute_th(self):
         MDC = kNN_MinimumDistance_Classifier(cluster_centers, cluster_labels, n_neighbors=self.n_neighbors, CPUs=None)
         MDC.classify(test_gA, in_nodataVal=-9999, cmap_nodataVal=-9999, tiledims=(400, 200))
-        MDC.label_unclassified_pixels(label_unclassified=-1, threshold=10)
+        MDC.label_unclassified_pixels(label_unclassified=-1, threshold=6000)
 
     def test_label_unclassified_pixels_relative_th(self):
         MDC = kNN_MinimumDistance_Classifier(cluster_centers, cluster_labels, self.n_neighbors, CPUs=None)
@@ -210,6 +210,39 @@ class Test_FEDSA_Classifier(unittest.TestCase):
         self.assertTrue(np.array_equal(cmap_mp.flatten(), cluster_labels))
 
 
+class Test_kNN_FEDSA_Classifier(unittest.TestCase):
+    def setUp(self) -> None:
+        self.n_neighbors = 5
+
+    def test_classify(self):
+        FC = kNN_FEDSA_Classifier(cluster_centers, n_neighbors=self.n_neighbors, CPUs=1)
+        cmap_sp = FC.classify(test_gA, in_nodataVal=-9999, cmap_nodataVal=-9999, tiledims=(400, 200))
+        self.assertIsInstance(cmap_sp, GeoArray)
+        self.assertEqual(cmap_sp.shape, (1010, 1010, self.n_neighbors))
+
+        FC = kNN_FEDSA_Classifier(cluster_centers, n_neighbors=self.n_neighbors, CPUs=None)
+        cmap_mp = FC.classify(test_gA, in_nodataVal=-9999, cmap_nodataVal=-9999, tiledims=(400, 200))
+        self.assertIsInstance(cmap_mp, GeoArray)
+        self.assertEqual(cmap_mp.shape, (1010, 1010, self.n_neighbors))
+
+        self.assertTrue(np.array_equal(cmap_sp, cmap_mp))
+
+        FC = kNN_FEDSA_Classifier(cluster_centers, n_neighbors=self.n_neighbors, CPUs=None)
+        cmap_mp = FC.classify(test_gA_pure_endmembers, in_nodataVal=-9999, cmap_nodataVal=-9999)
+
+        self.assertTrue(np.array_equal(cmap_mp[:, :, 0].flatten(), cluster_labels))
+
+    def test_label_unclassified_pixels_absolute_th(self):
+        FC = kNN_FEDSA_Classifier(cluster_centers, n_neighbors=self.n_neighbors, CPUs=None)
+        FC.classify(test_gA, in_nodataVal=-9999, cmap_nodataVal=-9999, tiledims=(400, 200))
+        FC.label_unclassified_pixels(label_unclassified=-1, threshold=10)
+
+    def test_label_unclassified_pixels_relative_th(self):
+        FC = kNN_FEDSA_Classifier(cluster_centers, self.n_neighbors, CPUs=None)
+        FC.classify(test_gA, in_nodataVal=-9999, cmap_nodataVal=-9999, tiledims=(400, 200))
+        FC.label_unclassified_pixels(label_unclassified=-1, threshold='10%')
+
+
 class Test_SID_Classifier(unittest.TestCase):
     def test_classify(self):
         SC = SID_Classifier(cluster_centers, CPUs=1)