Skip to content

ratiopath.augmentations.estimate_stain_vectors

DAB = np.array([0.27, 0.57, 0.78]) module-attribute

EOSIN = np.array([0.2159, 0.8012, 0.5581]) module-attribute

HDAB = np.array([HEMATOXYLIN, DAB, make_residual_stain(HEMATOXYLIN, DAB)], dtype=np.float32) module-attribute

HE = np.array([HEMATOXYLIN, EOSIN, make_residual_stain(HEMATOXYLIN, EOSIN)], dtype=np.float32) module-attribute

HEMATOXYLIN = np.array([0.65, 0.7, 0.29]) module-attribute

discard_pixels(od, min_stain, max_stain, gray_threshold=np.cos(0.15))

Discard pixels based on optical density thresholds.

Parameters:

Name Type Description Default
od ndarray

A numpy array of optical densities for red, green, and blue channels.

 required
min_stain float

Minimum optical density threshold.

 required
max_stain float

Maximum optical density threshold.

 required
gray_threshold float

Threshold for excluding very gray pixels (default is cos(0.15)).

 cos(0.15)

Returns:

Type Description
ndarray

A numpy array containing the filtered optical densities for red, green, and blue channels.
Source code in ratiopath/augmentations/estimate_stain_vectors.py 
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
def discard_pixels(
    od: np.ndarray,
    min_stain: float,
    max_stain: float,
    gray_threshold: float = np.cos(0.15),
) -> np.ndarray:
    """Discard pixels based on optical density thresholds.

    Parameters:
        od: A numpy array of optical densities for red, green, and blue channels.
        min_stain: Minimum optical density threshold.
        max_stain: Maximum optical density threshold.
        gray_threshold: Threshold for excluding very gray pixels (default is cos(0.15)).

    Returns:
        A numpy array containing the filtered optical densities for red, green, and blue channels.
    """
    mag_squared = np.sum(od**2, axis=1)

    keep_mask = (
        (mag_squared <= max_stain**2)
        & np.all(od >= min_stain, axis=1)
        & (mag_squared > 0)
    )

    # Filtering first reduces the number of expensive operations (sqrt and division)
    od = od[keep_mask]
    mag_squared = mag_squared[keep_mask]

    sqrt3 = 1 / np.sqrt(3)
    gray_mask = (np.sum(od * sqrt3, axis=1) / np.sqrt(mag_squared)) < gray_threshold

    return od[gray_mask]

estimate_stain_vectors(image, default_stain_vectors, i0=256, min_stain=0.05, max_stain=1.0, alpha=0.01)

Estimate stain vectors from an image using optical density transformation.

Parameters:

Name Type Description Default
image ndarray

A numpy array representing the input image.

 required
default_stain_vectors ndarray

A numpy array of default unit stain vectors.

 required
i0 int

The intensity value for normalization.

 256
min_stain float

Minimum optical density threshold for discarding pixels.

 0.05
max_stain float

Maximum optical density threshold for discarding pixels.

 1.0
alpha float

The percentage of pixels to use for estimating the stain vectors (default is 0.01, which corresponds to 1%).

 0.01

Returns:

Type Description
ndarray

A numpy array of estimated stain vectors.
References

Paper: A method for normalizing histology slides for quantitative analysis, M Macenko, M Niethammer, JS Marron, D Borland, JT Woosley, G Xiaojun, C Schmitt, NE Thomas, IEEE ISBI, 2009. dx.doi.org/10.1109/ISBI.2009.5193250

Source code in ratiopath/augmentations/estimate_stain_vectors.py 
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
def estimate_stain_vectors(
    image: np.ndarray,
    default_stain_vectors: np.ndarray,
    i0: int = 256,
    min_stain: float = 0.05,
    max_stain: float = 1.0,
    alpha: float = 0.01,
) -> np.ndarray:
    """Estimate stain vectors from an image using optical density transformation.

    Parameters:
        image: A numpy array representing the input image.
        default_stain_vectors: A numpy array of default unit stain vectors.
        i0: The intensity value for normalization.
        min_stain: Minimum optical density threshold for discarding pixels.
        max_stain: Maximum optical density threshold for discarding pixels.
        alpha: The percentage of pixels to use for estimating the stain vectors
            (default is 0.01, which corresponds to 1%).

    Returns:
        A numpy array of estimated stain vectors.

    References:
        Paper: A method for normalizing histology slides for quantitative analysis,
            M Macenko, M Niethammer, JS Marron, D Borland, JT Woosley, G Xiaojun,
            C Schmitt, NE Thomas, IEEE ISBI, 2009. dx.doi.org/10.1109/ISBI.2009.5193250
    """
    od = -np.log10(np.maximum(image.reshape(-1, 3), 1) / i0)

    od_hat = discard_pixels(od, min_stain, max_stain)

    cov = np.cov(od_hat.T)

    _, eigvecs = np.linalg.eigh(cov)

    t_hat = od_hat.dot(eigvecs[:, 1:3])

    phi = np.arctan2(t_hat[:, 1], t_hat[:, 0])

    min_phi = np.percentile(phi, alpha)
    max_phi = np.percentile(phi, 100 - alpha)

    direction_min = np.array([np.cos(min_phi), np.sin(min_phi)], dtype=np.float32)
    stain1 = eigvecs[:, 1:3].dot(direction_min).reshape(-1)
    stain1 /= np.linalg.norm(stain1)

    direction_max = np.array([np.cos(max_phi), np.sin(max_phi)], dtype=np.float32)
    stain2 = eigvecs[:, 1:3].dot(direction_max).reshape(-1)
    stain2 /= np.linalg.norm(stain2)

    cos_angle11 = np.dot(stain1, default_stain_vectors[0])
    cos_angle12 = np.dot(stain1, default_stain_vectors[1])
    cos_angle21 = np.dot(stain2, default_stain_vectors[0])
    cos_angle22 = np.dot(stain2, default_stain_vectors[1])

    if max(cos_angle12, cos_angle21) > max(cos_angle11, cos_angle22):
        stain1, stain2 = stain2, stain1

    return np.array([stain1, stain2, make_residual_stain(stain1, stain2)])

make_residual_stain(stain1, stain2)

Create a residual stain vector from two stain vectors.

Parameters:

Name Type Description Default
stain1 ndarray

A numpy array representing the first stain vector.

 required
stain2 ndarray

A numpy array representing the second stain vector.

 required

Returns:

Type Description
ndarray

A numpy array representing the residual stain vector.
Source code in ratiopath/augmentations/estimate_stain_vectors.py 
10
11
12
13
14
15
16
17
18
19
20
21
def make_residual_stain(stain1: np.ndarray, stain2: np.ndarray) -> np.ndarray:
    """Create a residual stain vector from two stain vectors.

    Parameters:
        stain1: A numpy array representing the first stain vector.
        stain2: A numpy array representing the second stain vector.

    Returns:
        A numpy array representing the residual stain vector.
    """
    res = np.linalg.cross(stain1, stain2)
    return res / np.linalg.norm(res)