Correct H&E stain ordering heuristic in ExtractHEStains (#8551)

### Description

I noticed when attempting to extract the hematoxylin stain from an H&E
stained image that the previous heuristic for ordering hematoxylin first
and eosin second incorrectly assumed hematoxylin has a higher red
channel value, when the opposite is typically true. I've corrected this
and implemented a more robust heuristic that compares red-blue ratios
instead, since:

- Hematoxylin (nuclear, blue): lower red/blue ratio
- Eosin (cytoplasm, pink): higher red/blue ratio

I've updated the tests to reflect the corrected stain ordering, and the
output now accurately reflects the documented behaviour (first column
hematoxylin, second eosin).

### Types of changes

- [x] Non-breaking change (fix or new feature that would not break
existing functionality).
- [x] Integration tests passed locally by running `./runtests.sh -f -u
--net --coverage`.
- [x] Quick tests passed locally by running `./runtests.sh --quick
--unittests --disttests`.

---------

Signed-off-by: Iyassou Shimels <s.iyassou@gmail.com>
Signed-off-by: Bruce Hashemian <3968947+bhashemian@users.noreply.github.com>
Co-authored-by: Bruce Hashemian <3968947+bhashemian@users.noreply.github.com>

Author: iyassou

monai/apps/pathology/transforms/stain/array.py

@@ -85,7 +85,12 @@ def _deconvolution_extract_stain(self, image: np.ndarray) -> np.ndarray:
         v_max = eigvecs[:, 1:3].dot(np.array([(np.cos(max_phi), np.sin(max_phi))], dtype=np.float32).T)
 
         # a heuristic to make the vector corresponding to hematoxylin first and the one corresponding to eosin second
-        if v_min[0] > v_max[0]:
+        # Hematoxylin: high blue, lower red (low R/B ratio)
+        # Eosin: high red, lower blue (high R/B ratio)
+        eps = np.finfo(np.float32).eps
+        v_min_rb_ratio = v_min[0] / (v_min[2] + eps)
+        v_max_rb_ratio = v_max[0] / (v_max[2] + eps)
+        if v_min_rb_ratio < v_max_rb_ratio:
             he = np.array((v_min[:, 0], v_max[:, 0]), dtype=np.float32).T
         else:
             he = np.array((v_max[:, 0], v_min[:, 0]), dtype=np.float32).T

tests/apps/pathology/transforms/test_pathology_he_stain.py

@@ -48,7 +48,7 @@
 # input pixels not uniformly filled, leading to two different stains extracted
 EXTRACT_STAINS_TEST_CASE_5 = [
     np.array([[[100, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]]]),
-    np.array([[0.70710677, 0.18696113], [0.0, 0.0], [0.70710677, 0.98236734]]),
+    np.array([[0.18696113, 0.70710677], [0.0, 0.0], [0.98236734, 0.70710677]]),
 ]
 
 # input pixels all transparent and below the beta absorbance threshold
@@ -68,7 +68,7 @@
 NORMALIZE_STAINS_TEST_CASE_4 = [
     {"target_he": np.full((3, 2), 1)},
     np.array([[[100, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]]]),
-    np.array([[[87, 87, 87], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]]]),
+    np.array([[[31, 31, 31], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]]]),
 ]
 
 
@@ -135,7 +135,7 @@ def test_result_value(self, image, expected_data):
           [[0.18696113],[0],[0.98236734]] and
           [[0.70710677],[0],[0.70710677]] respectively
         - the resulting extracted stain should be
-          [[0.70710677,0.18696113],[0,0],[0.70710677,0.98236734]]
+          [[0.18696113,0.70710677],[0,0],[0.98236734,0.70710677]]
         """
         if image is None:
             with self.assertRaises(TypeError):
@@ -206,17 +206,17 @@ def test_result_value(self, arguments, image, expected_data):
 
         For test case 4:
         - For this non-uniformly filled image, the stain extracted should be
-          [[0.70710677,0.18696113],[0,0],[0.70710677,0.98236734]], as validated for the
+          [[0.18696113,0.70710677],[0,0],[0.98236734,0.70710677]], as validated for the
           ExtractHEStains class. Solving the linear least squares problem (since
           absorbance matrix = stain matrix * concentration matrix), we obtain the concentration
-          matrix that should be [[-0.3101, 7.7508, 7.7508, 7.7508, 7.7508, 7.7508],
-          [5.8022, 0, 0, 0, 0, 0]]
+          matrix that should be [[5.8022, 0, 0, 0, 0, 0],
+          [-0.3101, 7.7508, 7.7508, 7.7508, 7.7508, 7.7508]]
         - Normalizing the concentration matrix, taking the matrix product of the
           target stain matrix and the concentration matrix, using the inverse
           Beer-Lambert transform to obtain the RGB image from the absorbance
           image, and finally converting to uint8, we get that the stain normalized
-          image should be [[[87, 87, 87], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]],
-          [[33, 33, 33], [33, 33, 33]]]
+          image should be [[[31, 31, 31], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]],
+          [[85, 85, 85], [85, 85, 85]]]
         """
         if image is None:
             with self.assertRaises(TypeError):

tests/apps/pathology/transforms/test_pathology_he_stain_dict.py

@@ -42,7 +42,7 @@
 # input pixels not uniformly filled, leading to two different stains extracted
 EXTRACT_STAINS_TEST_CASE_5 = [
     np.array([[[100, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]]]),
-    np.array([[0.70710677, 0.18696113], [0.0, 0.0], [0.70710677, 0.98236734]]),
+    np.array([[0.18696113, 0.70710677], [0.0, 0.0], [0.98236734, 0.70710677]]),
 ]
 
 # input pixels all transparent and below the beta absorbance threshold
@@ -62,7 +62,7 @@
 NORMALIZE_STAINS_TEST_CASE_4 = [
     {"target_he": np.full((3, 2), 1)},
     np.array([[[100, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]], [[0, 0, 0], [0, 0, 0]]]),
-    np.array([[[87, 87, 87], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]]]),
+    np.array([[[31, 31, 31], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]]]),
 ]
 
 
@@ -129,7 +129,7 @@ def test_result_value(self, image, expected_data):
           [[0.18696113],[0],[0.98236734]] and
           [[0.70710677],[0],[0.70710677]] respectively
         - the resulting extracted stain should be
-          [[0.70710677,0.18696113],[0,0],[0.70710677,0.98236734]]
+          [[0.18696113,0.70710677],[0,0],[0.98236734,0.70710677]]
         """
         key = "image"
         if image is None:
@@ -200,17 +200,17 @@ def test_result_value(self, arguments, image, expected_data):
 
         For test case 4:
         - For this non-uniformly filled image, the stain extracted should be
-          [[0.70710677,0.18696113],[0,0],[0.70710677,0.98236734]], as validated for the
+          [[0.18696113,0.70710677],[0,0],[0.98236734,0.70710677]], as validated for the
           ExtractHEStains class. Solving the linear least squares problem (since
           absorbance matrix = stain matrix * concentration matrix), we obtain the concentration
-          matrix that should be [[-0.3101, 7.7508, 7.7508, 7.7508, 7.7508, 7.7508],
-          [5.8022, 0, 0, 0, 0, 0]]
+          matrix that should be  [[5.8022, 0, 0, 0, 0, 0],
+          [-0.3101, 7.7508, 7.7508, 7.7508, 7.7508, 7.7508]]
         - Normalizing the concentration matrix, taking the matrix product of the
           target stain matrix and the concentration matrix, using the inverse
           Beer-Lambert transform to obtain the RGB image from the absorbance
           image, and finally converting to uint8, we get that the stain normalized
-          image should be [[[87, 87, 87], [33, 33, 33]], [[33, 33, 33], [33, 33, 33]],
-          [[33, 33, 33], [33, 33, 33]]]
+          image should be [[[31, 31, 31], [85, 85, 85]], [[85, 85, 85], [85, 85, 85]],
+          [[85, 85, 85], [85, 85, 85]]]
         """
         key = "image"
         if image is None: