Downsampling

[karray]

I'm currently trying to produce patches. Is it possible to explicitly define a WSI level that will be used for the patches? Or the tile are always produced from the level 0?

I also would like to know what is the purpose of all *-downsample options. Can I define --output-downsample 1 to produce the patches with organal resolution?

[manuel-munoz-aguirre]

The tiles are produced from the level that is closest to the required downsample. For example, if level 0 contains the native resolution, level 1 contains 4x and level 2 contains 16x, if you specify --output-downsample 16 then level 2 will be used. Currently it's only possible to specify the downsampling factor and not the level. You can use openslide to figure out the downsampling factors:

>>> import openslide
>>> slide = openslide.OpenSlide("slide.svs")
>>> slide.level_count
3
>>> slide.level_dimensions
((47807, 38653), (11951, 9663), (2987, 2415))
>>> slide.level_downsamples
(1.0, 4.00017725627429, 16.005202391869254)

The *-downsample options control the downsampling level that will be used for the specific task: --output-downsample is the factor used to scale down the slide to generate the patches, so using a value of 1 will produce patches at the original resolution (level 0). The --mask-downsample flag is there as a tradeoff between segmentation quality and speed, if the value is large, the segmentation will be less detailed but it will also be faster.

[karray]

@manuel-munoz-aguirre thanks for your fast replay.

I'm currently trying to adapt your code for my project. It is not quite clear what you are doing in utility_functions.py:

# Get the best level to quickly downsample the image
# Add a pseudofactor of 0.1 to ensure getting the next best level
# (i.e. if 16x is chosen, avoid getting 4x instead of 16x)
best_downsampling_level = slide.get_best_level_for_downsample(downsampling_factor + 0.1)

# Get the image at the requested scale
svs_native_levelimg = slide.read_region((0, 0), best_downsampling_level, slide.level_dimensions[best_downsampling_level])
target_size = tuple([int(x//downsampling_factor) for x in slide.dimensions])
img = svs_native_levelimg.resize(target_size)

If I understand correctly, get_best_level_for_downsample already gives the best level for the requested downsample factor. Why do you resize this scaled version once again?

[manuel-munoz-aguirre]

This is done because a user might want to downscale the image at a factor that is not natively encoded in the slide. If you request a downsampling factor, for example, 4x, the final size will be the same returned by openslide:

>>> import openslide
>>> downsampling_factor = 4
>>> slide = openslide.OpenSlide("slide.svs")
>>> slide.level_dimensions
((47807, 38653), (11951, 9663), (2987, 2415))
>>> slide.level_downsamples
(1.0, 4.00017725627429, 16.005202391869254)
>>> best_downsampling_level = slide.get_best_level_for_downsample(downsampling_factor+0.1)
>>> best_downsampling_level
1
>>> svs_native_levelimg = slide.read_region((0,0), best_downsampling_level, slide.level_dimensions[best_downsampling_level])
>>> svs_native_levelimg.size
(11951, 9663)
>>> target_size = tuple([int(x//downsampling_factor) for x in slide.dimensions])
>>> target_size
(11951, 9663)
>>> img = svs_native_levelimg.resize(target_size)
>>> img.size
(11951, 9663)

Note that slide.dimensions is the size at level 0. If you wanted to do a downsample of 8x, then the target size would be (5975, 4831), but it is not natively encoded, so we do the last step to account for those cases. There's no need to do it from the original resolution (it's larger, so it would take more time), so we pick a level with closer dimensions to do the resizing. We can probably add a check to make this last step only when the dimension isn't natively encoded.

[karray]

Thank you for your support. I have one last question. Can you please give me a short hint on how I can reassemble the extracted patches back into WSI.

[manuel-munoz-aguirre]

To be able to reconstruct the original WSI, the call to PyHIST should include the --save-patches, --save-blank and --save-nonsquare flags (the second ensures that the non-selected tiles are also saved to disk, while the third flag accounts for the patches in the edges of the WSI that are not multiple of the tile size). For example:

python3 pyhist.py --output-downsample 16 --save-patches --save-blank --save-nonsquare --output output/ slide.svs

Then, you can reassemble the WSI by using the montage tool from Imagemagick. You can check in the generated metadata (tile_selection.tsv) for your slide the number of columns in the WSI. For example, if your WSI has 6 columns (indexed in the file from 0 to 5), it can be reconstructed to a file reconstructed.png with:

montage *.png -mode concatenate -tile 6x reconstructed.png

[karray]

Thank you for your support. The problem is that in order to train my model I want to remove blank patches. I ended up with the following script, which I am posting here in case anyone is looking for a similar solution.

import pandas as pd
from os.path import join as ospj
from PIL import Image as PImage
import pyvips
import numpy as np

def merge_patches(idx_file, input_path, output_path, ouput_name, patch_size=128, compression='jpeg', Q=80):
    df = pd.read_csv(idx_file, delimiter='\t')

    min_col = df.Column.min()
    max_col = df.Column.max()
    min_row = df.Row.min()
    max_row = df.Row.max()

    w = max_col*patch_size
    h = max_row *patch_size

    wsi= np.full((h, w, 3), 255, dtype=np.ubyte)
    print('wsi size: ', wsi.shape)

    df = df[df.Keep == 1].copy()
    i = 0
    for name, col, row in zip(df['Tile'], df['Column'], df['Row']):
        y_start = col*patch_size
        y_end = col*patch_size+patch_size

        x_start = row*patch_size
        x_end = row*patch_size+patch_size
        img = PImage.open(ospj(input_path,name+'.png'))
        wsi[x_start:x_end, y_start:y_end] = np.array(img, dtype=np.ubyte)

    slide = pyvips.Image.new_from_memory(wsi, wsi.shape[1], wsi.shape[0], wsi.shape[2], 'uchar')
    slide.tiffsave(ospj(output_path, ouput_name+'.tiff'), compression=compression, Q=Q, tile=True, tile_width=patch_size, tile_height=patch_size, pyramid=True)
    wsi = None

merge_patches('path/to/tile_selection.tsv', 'path/to/patches/',  'output/', 'slide_name')