torch_em.data.datasets.histopathology.camelyon

  1"""This dataset contains annotations for tumor region segmentation in
  2whole-slide histopathology images of breast cancer sentinel lymph node sections.
  3
  4The data is from the CAMELYON16 and CAMELYON17 challenges, described in
  5https://doi.org/10.1093/gigascience/giy065 ("1399 H&E-stained sentinel lymph
  6node sections of breast cancer patients: the CAMELYON dataset"). It is hosted
  7on the public AWS Open Data Registry (https://registry.opendata.aws/camelyon/)
  8under a CC0 license. Please cite the publication if you use this dataset.
  9
 10Two segmentation-compatible sources are exposed:
 11- CAMELYON16: all 399 whole-slide images (159 normal, 111 tumor, 129 test),
 12  each paired with a tumor segmentation mask.
 13- CAMELYON17: the subset of 100 training whole-slide images that ship with a
 14  segmentation mask. The remaining CAMELYON17 slides only carry a patient-level
 15  pN-stage classification label and are out of scope for this loader.
 16
 17Mask label values: 0 (background / non-tumor tissue), 1 (tumor), and for
 18CAMELYON16 tumor slides additionally 2 (non-tumor tissue excluded from a
 19tumor annotation).
 20
 21NOTE: The whole-slide images and masks are multi-resolution pyramidal TIFFs of
 22several gigabytes each. On the first use each requested slide is converted into
 23a chunked HDF5 file at the requested pyramid level, which requires some time
 24and disk space.
 25"""
 26
 27import os
 28from pathlib import Path
 29from typing import List, Literal, Optional, Tuple, Union
 30
 31from tqdm import tqdm
 32
 33import torch
 34from torch.utils.data import Dataset, DataLoader
 35
 36import torch_em
 37
 38from .. import util
 39
 40
 41BASE_URL = "https://camelyon-dataset.s3.us-west-2.amazonaws.com"
 42
 43CHECKSUM_URLS = {
 44    "CAMELYON16": f"{BASE_URL}/CAMELYON16/checksums.md5",
 45    "CAMELYON17": f"{BASE_URL}/CAMELYON17/checksums.md5",
 46}
 47
 48
 49def _get_checksums(path, version, download):
 50    checksum_path = os.path.join(path, f"{version.lower()}_checksums.md5")
 51    util.download_source(path=checksum_path, url=CHECKSUM_URLS[version], download=download, checksum=None)
 52
 53    checksums = {}
 54    with open(checksum_path) as f:
 55        for line in f:
 56            checksum, name = line.split(maxsplit=1)
 57            checksums[name.strip().lstrip("*")] = checksum
 58    return checksums
 59
 60
 61def _verify_md5(path, expected):
 62    import hashlib
 63
 64    actual = hashlib.md5(Path(path).read_bytes()).hexdigest()
 65    if actual != expected:
 66        raise RuntimeError(f"The checksum of '{path}' does not match the expected checksum: {expected} != {actual}")
 67
 68
 69def _download_slide(path, version, stem, checksums, download):
 70    raw_dir = os.path.join(path, "raw", version)
 71    os.makedirs(raw_dir, exist_ok=True)
 72
 73    image_path = os.path.join(raw_dir, f"{stem}.tif")
 74    mask_path = os.path.join(raw_dir, f"{stem}_mask.tif")
 75
 76    for out_path, key in [(image_path, f"images/{stem}.tif"), (mask_path, f"masks/{stem}_mask.tif")]:
 77        if os.path.exists(out_path):
 78            continue
 79        util.download_source(path=out_path, url=f"{BASE_URL}/{version}/{key}", download=download, checksum=None)
 80        _verify_md5(out_path, checksums[key])
 81
 82    return image_path, mask_path
 83
 84
 85def _open_level(series, level_index):
 86    import zarr
 87
 88    # The pyramidal TIFFs are natively tiled, so a zarr view reads only the requested tiles lazily.
 89    array = zarr.open(series.aszarr(), mode="r")
 90    return array if hasattr(array, "shape") else array[str(level_index)]
 91
 92
 93def _convert_slide(image_path, mask_path, output_path, resolution_level, tile=4096):
 94    import h5py
 95    import tifffile
 96
 97    image_series = tifffile.TiffFile(image_path).series[0]
 98    mask_series = tifffile.TiffFile(mask_path).series[0]
 99
100    # Some scanners (e.g. Philips) round the declared level shape slightly differently between the raw
101    # slide and the independently rasterized mask, so tolerate a small mismatch and crop to the overlap.
102    image_height, image_width = image_series.levels[resolution_level].shape[:2]
103    mask_height, mask_width = mask_series.levels[resolution_level].shape[:2]
104    height, width = min(image_height, mask_height), min(image_width, mask_width)
105    tolerance = 0.02
106    if abs(image_height - mask_height) > tolerance * image_height or \
107            abs(image_width - mask_width) > tolerance * image_width:
108        raise RuntimeError(
109            f"The mask '{mask_path}' does not match the raw shape ({image_height}, {image_width}) "
110            f"at level {resolution_level}: got ({mask_height}, {mask_width})."
111        )
112
113    image = _open_level(image_series, resolution_level)
114    mask = _open_level(mask_series, resolution_level)
115
116    tmp_path = output_path + ".tmp"
117    with h5py.File(tmp_path, "w") as f:
118        raw = f.create_dataset(
119            "images/raw", shape=(3, height, width), dtype="uint8", compression="gzip", chunks=(1, 512, 512)
120        )
121        labels = f.create_dataset(
122            "labels/mask", shape=(height, width), dtype="uint8", compression="gzip", chunks=(512, 512)
123        )
124        for y in tqdm(range(0, height, tile), desc=f"Converting {Path(image_path).stem}"):
125            for x in range(0, width, tile):
126                th, tw = min(tile, height - y), min(tile, width - x)
127                raw[:, y:y + th, x:x + tw] = image[y:y + th, x:x + tw].transpose(2, 0, 1)
128                labels[y:y + th, x:x + tw] = mask[y:y + th, x:x + tw]
129
130    os.replace(tmp_path, output_path)
131
132
133def _restrict_to_sample_ids(stems, sample_ids):
134    if sample_ids is None:
135        return stems
136    missing = sorted(set(sample_ids) - set(stems))
137    if missing:
138        raise ValueError(f"The following sample ids are not part of this dataset: {missing}")
139    return [stem for stem in stems if stem in sample_ids]
140
141
142def _resolve_camelyon16_stems(checksums, split, sample_ids):
143    stems = sorted(Path(name).stem for name in checksums if name.startswith("images/"))
144    if split is not None:
145        assert split in ("train", "test"), "Please choose from the available `train` / `test` splits"
146        prefixes = ("normal_", "tumor_") if split == "train" else ("test_",)
147        stems = [stem for stem in stems if stem.startswith(prefixes)]
148    return _restrict_to_sample_ids(stems, sample_ids)
149
150
151def _resolve_camelyon17_stems(checksums, sample_ids):
152    stems = sorted(Path(name).stem[:-len("_mask")] for name in checksums if name.startswith("masks/"))
153    return _restrict_to_sample_ids(stems, sample_ids)
154
155
156def get_camelyon16_data(
157    path: Union[os.PathLike, str],
158    split: Optional[Literal["train", "test"]] = None,
159    sample_ids: Optional[List[str]] = None,
160    resolution_level: int = 0,
161    download: bool = False,
162) -> str:
163    """Download and preprocess the CAMELYON16 tumor segmentation data.
164
165    Args:
166        path: Filepath to a folder where the data will be saved.
167        split: The split to use. Either 'train' (normal and tumor slides) or 'test'. By default all slides are used.
168        sample_ids: The slide names to restrict the data to, e.g. ['tumor_091', 'normal_108'].
169            By default all slides (matching `split`) are used.
170        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
171            reduce the size of the preprocessed data.
172        download: Whether to download the data if it is not present.
173
174    Returns:
175        Filepath to the folder where the preprocessed data is stored.
176    """
177    checksums = _get_checksums(path, "CAMELYON16", download)
178    stems = _resolve_camelyon16_stems(checksums, split, sample_ids)
179
180    preprocessed_dir = os.path.join(path, "preprocessed", "CAMELYON16")
181    os.makedirs(preprocessed_dir, exist_ok=True)
182
183    for stem in stems:
184        output_path = os.path.join(preprocessed_dir, f"{stem}_level{resolution_level}.h5")
185        if os.path.exists(output_path):
186            continue
187        image_path, mask_path = _download_slide(path, "CAMELYON16", stem, checksums, download)
188        _convert_slide(image_path, mask_path, output_path, resolution_level)
189
190    return preprocessed_dir
191
192
193def get_camelyon16_paths(
194    path: Union[os.PathLike, str],
195    split: Optional[Literal["train", "test"]] = None,
196    sample_ids: Optional[List[str]] = None,
197    resolution_level: int = 0,
198    download: bool = False,
199) -> List[str]:
200    """Get paths to the CAMELYON16 tumor segmentation data.
201
202    Args:
203        path: Filepath to a folder where the data will be saved.
204        split: The split to use. Either 'train' (normal and tumor slides) or 'test'. By default all slides are used.
205        sample_ids: The slide names to restrict the data to, e.g. ['tumor_091', 'normal_108'].
206            By default all slides (matching `split`) are used.
207        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
208            reduce the size of the preprocessed data.
209        download: Whether to download the data if it is not present.
210
211    Returns:
212        List of filepaths to the preprocessed HDF5 files.
213    """
214    preprocessed_dir = get_camelyon16_data(path, split, sample_ids, resolution_level, download)
215    stems = _resolve_camelyon16_stems(_get_checksums(path, "CAMELYON16", download=False), split, sample_ids)
216    return [os.path.join(preprocessed_dir, f"{stem}_level{resolution_level}.h5") for stem in stems]
217
218
219def get_camelyon17_data(
220    path: Union[os.PathLike, str],
221    sample_ids: Optional[List[str]] = None,
222    resolution_level: int = 0,
223    download: bool = False,
224) -> str:
225    """Download and preprocess the CAMELYON17 tumor segmentation subset.
226
227    This is the subset of CAMELYON17 training slides that ship with a segmentation mask.
228    The remaining slides only carry a patient-level classification label and are not included.
229
230    Args:
231        path: Filepath to a folder where the data will be saved.
232        sample_ids: The slide names to restrict the data to, e.g. ['patient_000_node_4'].
233            By default all slides with a segmentation mask are used.
234        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
235            reduce the size of the preprocessed data.
236        download: Whether to download the data if it is not present.
237
238    Returns:
239        Filepath to the folder where the preprocessed data is stored.
240    """
241    checksums = _get_checksums(path, "CAMELYON17", download)
242    stems = _resolve_camelyon17_stems(checksums, sample_ids)
243
244    preprocessed_dir = os.path.join(path, "preprocessed", "CAMELYON17")
245    os.makedirs(preprocessed_dir, exist_ok=True)
246
247    for stem in stems:
248        output_path = os.path.join(preprocessed_dir, f"{stem}_level{resolution_level}.h5")
249        if os.path.exists(output_path):
250            continue
251        image_path, mask_path = _download_slide(path, "CAMELYON17", stem, checksums, download)
252        _convert_slide(image_path, mask_path, output_path, resolution_level)
253
254    return preprocessed_dir
255
256
257def get_camelyon17_paths(
258    path: Union[os.PathLike, str],
259    sample_ids: Optional[List[str]] = None,
260    resolution_level: int = 0,
261    download: bool = False,
262) -> List[str]:
263    """Get paths to the CAMELYON17 tumor segmentation subset.
264
265    Args:
266        path: Filepath to a folder where the data will be saved.
267        sample_ids: The slide names to restrict the data to, e.g. ['patient_000_node_4'].
268            By default all slides with a segmentation mask are used.
269        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
270            reduce the size of the preprocessed data.
271        download: Whether to download the data if it is not present.
272
273    Returns:
274        List of filepaths to the preprocessed HDF5 files.
275    """
276    preprocessed_dir = get_camelyon17_data(path, sample_ids, resolution_level, download)
277    stems = _resolve_camelyon17_stems(_get_checksums(path, "CAMELYON17", download=False), sample_ids)
278    return [os.path.join(preprocessed_dir, f"{stem}_level{resolution_level}.h5") for stem in stems]
279
280
281def _get_dataset(
282    volume_paths, patch_shape, label_dtype, resize_inputs, kwargs,
283) -> Dataset:
284    if resize_inputs:
285        resize_kwargs = {"patch_shape": patch_shape, "is_rgb": True}
286        kwargs, patch_shape = util.update_kwargs_for_resize_trafo(
287            kwargs=kwargs, patch_shape=patch_shape, resize_inputs=resize_inputs, resize_kwargs=resize_kwargs
288        )
289
290    return torch_em.default_segmentation_dataset(
291        raw_paths=volume_paths,
292        raw_key="images/raw",
293        label_paths=volume_paths,
294        label_key="labels/mask",
295        patch_shape=patch_shape,
296        label_dtype=label_dtype,
297        is_seg_dataset=True,
298        with_channels=True,
299        ndim=2,
300        **kwargs
301    )
302
303
304def get_camelyon16_dataset(
305    path: Union[os.PathLike, str],
306    patch_shape: Tuple[int, int],
307    split: Optional[Literal["train", "test"]] = None,
308    sample_ids: Optional[List[str]] = None,
309    resolution_level: int = 0,
310    download: bool = False,
311    label_dtype: torch.dtype = torch.int64,
312    resize_inputs: bool = False,
313    **kwargs
314) -> Dataset:
315    """Get the CAMELYON16 dataset for tumor segmentation in whole-slide histopathology images.
316
317    Args:
318        path: Filepath to a folder where the data will be saved.
319        patch_shape: The patch shape to use for training.
320        split: The split to use. Either 'train' (normal and tumor slides) or 'test'. By default all slides are used.
321        sample_ids: The slide names to restrict the data to, e.g. ['tumor_091', 'normal_108'].
322            By default all slides (matching `split`) are used.
323        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
324            reduce the size of the preprocessed data.
325        download: Whether to download the data if it is not present.
326        label_dtype: The datatype of the labels.
327        resize_inputs: Whether to resize the input images.
328        kwargs: Additional keyword arguments for `torch_em.default_segmentation_dataset`.
329
330    Returns:
331        The segmentation dataset.
332    """
333    volume_paths = get_camelyon16_paths(path, split, sample_ids, resolution_level, download)
334    return _get_dataset(volume_paths, patch_shape, label_dtype, resize_inputs, kwargs)
335
336
337def get_camelyon16_loader(
338    path: Union[os.PathLike, str],
339    patch_shape: Tuple[int, int],
340    batch_size: int,
341    split: Optional[Literal["train", "test"]] = None,
342    sample_ids: Optional[List[str]] = None,
343    resolution_level: int = 0,
344    download: bool = False,
345    label_dtype: torch.dtype = torch.int64,
346    resize_inputs: bool = False,
347    **kwargs
348) -> DataLoader:
349    """Get the CAMELYON16 dataloader for tumor segmentation in whole-slide histopathology images.
350
351    Args:
352        path: Filepath to a folder where the data will be saved.
353        patch_shape: The patch shape to use for training.
354        batch_size: The batch size for training.
355        split: The split to use. Either 'train' (normal and tumor slides) or 'test'. By default all slides are used.
356        sample_ids: The slide names to restrict the data to, e.g. ['tumor_091', 'normal_108'].
357            By default all slides (matching `split`) are used.
358        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
359            reduce the size of the preprocessed data.
360        download: Whether to download the data if it is not present.
361        label_dtype: The datatype of the labels.
362        resize_inputs: Whether to resize the input images.
363        kwargs: Additional keyword arguments for `torch_em.default_segmentation_dataset` or for the PyTorch DataLoader.
364
365    Returns:
366        The DataLoader.
367    """
368    ds_kwargs, loader_kwargs = util.split_kwargs(torch_em.default_segmentation_dataset, **kwargs)
369    dataset = get_camelyon16_dataset(
370        path=path, patch_shape=patch_shape, split=split, sample_ids=sample_ids, resolution_level=resolution_level,
371        download=download, label_dtype=label_dtype, resize_inputs=resize_inputs, **ds_kwargs
372    )
373    return torch_em.get_data_loader(dataset, batch_size, **loader_kwargs)
374
375
376def get_camelyon17_dataset(
377    path: Union[os.PathLike, str],
378    patch_shape: Tuple[int, int],
379    sample_ids: Optional[List[str]] = None,
380    resolution_level: int = 0,
381    download: bool = False,
382    label_dtype: torch.dtype = torch.int64,
383    resize_inputs: bool = False,
384    **kwargs
385) -> Dataset:
386    """Get the CAMELYON17 tumor segmentation subset for whole-slide histopathology images.
387
388    Args:
389        path: Filepath to a folder where the data will be saved.
390        patch_shape: The patch shape to use for training.
391        sample_ids: The slide names to restrict the data to, e.g. ['patient_000_node_4'].
392            By default all slides with a segmentation mask are used.
393        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
394            reduce the size of the preprocessed data.
395        download: Whether to download the data if it is not present.
396        label_dtype: The datatype of the labels.
397        resize_inputs: Whether to resize the input images.
398        kwargs: Additional keyword arguments for `torch_em.default_segmentation_dataset`.
399
400    Returns:
401        The segmentation dataset.
402    """
403    volume_paths = get_camelyon17_paths(path, sample_ids, resolution_level, download)
404    return _get_dataset(volume_paths, patch_shape, label_dtype, resize_inputs, kwargs)
405
406
407def get_camelyon17_loader(
408    path: Union[os.PathLike, str],
409    patch_shape: Tuple[int, int],
410    batch_size: int,
411    sample_ids: Optional[List[str]] = None,
412    resolution_level: int = 0,
413    download: bool = False,
414    label_dtype: torch.dtype = torch.int64,
415    resize_inputs: bool = False,
416    **kwargs
417) -> DataLoader:
418    """Get the CAMELYON17 tumor segmentation subset dataloader for whole-slide histopathology images.
419
420    Args:
421        path: Filepath to a folder where the data will be saved.
422        patch_shape: The patch shape to use for training.
423        batch_size: The batch size for training.
424        sample_ids: The slide names to restrict the data to, e.g. ['patient_000_node_4'].
425            By default all slides with a segmentation mask are used.
426        resolution_level: The pyramid level to convert. 0 is the native resolution; use a higher level to
427            reduce the size of the preprocessed data.
428        download: Whether to download the data if it is not present.
429        label_dtype: The datatype of the labels.
430        resize_inputs: Whether to resize the input images.
431        kwargs: Additional keyword arguments for `torch_em.default_segmentation_dataset` or for the PyTorch DataLoader.
432
433    Returns:
434        The DataLoader.
435    """
436    ds_kwargs, loader_kwargs = util.split_kwargs(torch_em.default_segmentation_dataset, **kwargs)
437    dataset = get_camelyon17_dataset(
438        path=path, patch_shape=patch_shape, sample_ids=sample_ids, resolution_level=resolution_level,
439        download=download, label_dtype=label_dtype, resize_inputs=resize_inputs, **ds_kwargs
440    )
441    return torch_em.get_data_loader(dataset, batch_size, **loader_kwargs)