### Install histolytics from source with uv

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/installation.md

Install the latest development version of histolytics from its GitHub repository using uv.

```bash
uv pip install git+https://github.com/HautaniemiLab/histolytics.git
```

--------------------------------

### Install Required Libraries

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/finetuning.ipynb

Installs necessary libraries for the workflow, including torchdata and tables.

```python
# Some installations
# !pip install torchdata
# !pip install tables
```

--------------------------------

### Install Histolytics Package

Source: https://github.com/hautaniemilab/histolytics/blob/main/README.md

Use pip to install the histolytics package. This is the primary method for setting up the library.

```shell
pip install histolytics
```

--------------------------------

### Clone Repository and Install Dependencies

Source: https://github.com/hautaniemilab/histolytics/blob/main/CONTRIBUTING.md

Steps to clone the Histolytics repository, set up a virtual environment, and install development dependencies.

```bash
git clone https://github.com/HautaniemiLab/histolytics.git
cd histolytics

python -m venv venv
source venv/bin/activate  # On Windows: venv\Scripts\activate

pip install -e ".[dev]"
```

--------------------------------

### Install Seaborn

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/nuclear_morphology.ipynb

Install the seaborn library for data visualization. This is a prerequisite for visualizing the results of the morphological analysis.

```python
# Install seaborn for visualizing the results.
# !pip install seaborn
```

--------------------------------

### Install Dependencies

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/panoptic_segmentation.ipynb

Installs necessary libraries for the workflow. Ensure you have the correct CUDA version for cucim.

```python
# !pip install albumentations
# !pip install cucim-cu12
```

--------------------------------

### ApplyEach Transform Example

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/api/transforms/apply_each.md

Demonstrates the basic usage of the ApplyEach transform to apply a function to each item in a list.

```python
from histolytics.transforms import ApplyEach

def add_one(x):
    return x + 1

apply_each_transform = ApplyEach(add_one)

# Example usage:
input_list = [1, 2, 3, 4]
output_list = apply_each_transform(input_list)
print(output_list)  # Expected output: [2, 3, 4, 5]
```

--------------------------------

### Load and Plot Example Data

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/tls_lymphoid_aggregate.ipynb

Loads nuclei and tissue WSI data and plots them for visualization. Ensures necessary data modules are imported.

```python
from histolytics.data import hgsc_nuclei_wsi, hgsc_tissue_wsi

nuc = hgsc_nuclei_wsi()
tis = hgsc_tissue_wsi()

ax = tis.plot(figsize=(10, 10), column="class_name", alpha=0.2, aspect=1, legend=True)
nuc.plot(ax=ax, column="class_name", aspect=1)
ax.set_axis_off()
```

--------------------------------

### Initialize WSIPatchIterator for Feature Extraction

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/chromatin_patterns.ipynb

Sets up the WSIPatchIterator to process WSI patches. It requires a SlideReader, grid, nuclei segmentation map, and a pipeline function for feature extraction. This example uses a pipeline to compute chromatin area and proportion.

```python
from histolytics.wsi.wsi_iterator import WSIPatchIterator
from histolytics.nuc_feats.chromatin import chromatin_feats
from functools import partial

# define our pipeline
pipeline = partial(chromatin_feats, metrics=("chrom_area", "chrom_nuc_prop"))

crop_loader = WSIPatchIterator(
    slide_reader=reader,
    grid=gr,
    nuclei=nuc[nuc["class_name"] == "neoplastic"],  # use only neoplastic nuclei
    pipeline_func=pipeline,
    batch_size=8,
    num_workers=8,
    pin_memory=False,
    shuffle=False,
    drop_last=False,
)
```

--------------------------------

### Initialize Panoptic Cellpose with UNI Backbone

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Initialize the CellposePanoptic model using the UNI backbone from Hugging Face Hub. This example demonstrates passing custom keyword arguments for the timm encoder initialization and specifying output layer indices for U-net skip connections.

```python
from histolytics.models.cellpose_panoptic import CellposePanoptic

cpose_panoptic = CellposePanoptic(
    n_nuc_classes=6,
    n_tissue_classes=6,
    enc_name="hf_hub:MahmoodLab/uni",
    enc_pretrain=True,
    model_kwargs={
        "encoder_kws": {"init_values": 1e-5, "dynamic_img_size": True},
        "enc_out_indices": (2, 4, 6, 8),  # using layers 2, 4, 6, 8 from UNI
    },
)

# Print only a summary of the model instead of the full output
print(str(cpose_panoptic.model)[:500] + "\n...")
```

--------------------------------

### Check Library Versions

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/panoptic_segmentation.ipynb

Verifies the installed versions of PyTorch, Albumentations, cellseg_models_pytorch, and Python. This is useful for ensuring compatibility.

```python
from platform import python_version

import torch
import albumentations as A
import cellseg_models_pytorch

print("torch version:", torch.__version__)
print("albumentations version:", A.__version__)
print("cellseg_models_pytorch version:", cellseg_models_pytorch.__version__)
print("python version:", python_version())
```

--------------------------------

### Load Example Data and Plot Nuclei

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/nuclear_features.ipynb

Loads example nuclei segmentation and H&E image data. Displays the nuclei with their class names. This is a prerequisite for computing nuclear features.

```python
from histolytics.data import hgsc_cancer_he, hgsc_cancer_nuclei

# Load some example data
nuc = hgsc_cancer_nuclei()  # Nuclei vector segmentation
he = hgsc_cancer_he()  # Corresponding H&E image

ax = nuc.plot(figsize=(10, 10), column="class_name")
ax.set_axis_off()
nuc.head(5)
```

--------------------------------

### Visualize Grid on Slide Thumbnail

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/collagen_orientation.ipynb

This is a commented-out example showing how to visualize the processed grid on top of a slide thumbnail. Uncomment and adapt the 'reader' object to use this functionality.

```python
# Visualize the grid on top of the slide
# thumbnail = reader.read_level(-2)
# reader.get_annotated_thumbnail(thumbnail, gr, linewidth=1)
```

--------------------------------

### Initialize SlideReader and Process WSI

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/panoptic_segmentation.ipynb

Initializes the SlideReader with a specified backend (e.g., 'CUCIM'), generates a tissue mask, obtains tile coordinates, and filters them to get contiguous sub-grids. This is useful for focusing analysis on relevant tissue areas.

```python
from histolytics.wsi.slide_reader import SlideReader
from histolytics.wsi.utils import get_sub_grids
from pathlib import Path

sl_path = Path().home() / "sample.tiff"  # Replace with your WSI path

# Initialize SlideReader with CUCIM backend
reader = SlideReader(sl_path, backend="CUCIM")

# Get the tissue mask and coordinates
thresh, tissue_mask = reader.get_tissue_mask(level=-2)
coordinates = reader.get_tile_coordinates(width=1024, tissue_mask=tissue_mask)

# Get the contiguous sub-grids
sub_grids = get_sub_grids(coordinates)
filtered_sub_grids = [grid for grid in sub_grids if len(grid) > 100]

# pick the first sub-grid
sub_coords = filtered_sub_grids[0]

# visualize the thumbnail and the sub-coordinates
# thumbnail = reader.read_level(-2)
# reader.get_annotated_thumbnail(thumbnail, sub_coords, linewidth=1)
```

--------------------------------

### Python Setup for Legendgram Analysis

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/legendgram.ipynb

Imports necessary libraries and defines helper functions for calculating immune-connective cell ratios and performing quantile normalization. This code prepares the data for spatial aggregation and visualization.

```python
import pandas as pd
from histolytics.data import cervix_nuclei, cervix_tissue
from histolytics.spatial_ops.ops import get_interfaces
from histolytics.spatial_ops.h3 import h3_grid
from histolytics.spatial_agg import grid_aggregate
from histolytics.utils.plot import legendgram
from scipy.stats import rankdata


# This function will compute the ratio of immune to connective cells within each grid cell
# In general, any function that takes a GeoDataFrame and returns a scalar
# can be used here. Typically, this will be a function that calculates
# a count, sum, mean, or other statistic of interest out of the nuclei.
def immune_connective_ratio(nuclei):
    """Calculate the immune to connective cell ratio in a grid cell."""
    if "inflammatory" in nuclei.value_counts("class_name"):
        immune_cnt = nuclei.value_counts("class_name", normalize=True)["inflammatory"]
    else:
        immune_cnt = 0

    if "connective" in nuclei.value_counts("class_name"):
        connective_cnt = nuclei.value_counts("class_name", normalize=True)["connective"]
    else:
        connective_cnt = 0

    if connective_cnt > 0:
        return float(immune_cnt / connective_cnt)
    else:
        return 0


# quantile normalization function
def qnorm(values: pd.Series):
    """Quantile normalize a pandas Series."""
    ranks = rankdata(values, method="average")
    quantile_normalized = (ranks - 1) / (len(ranks) - 1)
    return quantile_normalized


# Load data
tis = cervix_tissue()
nuc = cervix_nuclei()

# get the stroma and lesion (cin) tissue segmentations
stroma = tis[tis["class_name"] == "stroma"]
cin_tissue = tis[tis["class_name"] == "cin"]

# get the lesion-stroma-interface
interface = get_interfaces(cin_tissue, stroma, buffer_dist=300)
interface = interface.assign(class_name="lesion-stroma-interface")

# fit a hexagonal grid to the interface
h3_res10_inter = h3_grid(interface, resolution=10)
# fit a hexagonal grid to the stroma
h3_res10_stroma = h3_grid(stroma, resolution=10)

# immune-connective-ratio at the interface
h3_res10_inter = grid_aggregate(
    objs=nuc,
    grid=h3_res10_inter,
    metric_func=immune_connective_ratio,
    new_col_names=["immune_connective_ratio"],
    predicate="intersects",
    num_processes=2,
)

# immune-connective-ratio at the stroma
h3_res10_stroma = grid_aggregate(
    objs=nuc,
    grid=h3_res10_stroma,
    metric_func=immune_connective_ratio,
    new_col_names=["immune_connective_ratio"],
    predicate="intersects",
    num_processes=2,
)

h3_res10_inter = h3_res10_inter.assign(
    immune_connective_ratio_qnorm=qnorm(h3_res10_inter["immune_connective_ratio"])
)
h3_res10_stroma = h3_res10_stroma.assign(
    immune_connective_ratio_qnorm=qnorm(h3_res10_stroma["immune_connective_ratio"])
)
```

--------------------------------

### Create H5 Files for Training and Validation

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/finetuning.ipynb

This script defines a function to create H5 files from image and mask data. It applies augmentations and saves processed patches for training and validation sets. Ensure necessary libraries like `tables`, `albumentations`, and `histolytics` are installed.

```python
import tables as tb
import albumentations as A

from pathlib import Path
from histolytics.utils import H5Handler
from cellseg_models_pytorch.inference.predictor import Predictor
from cellseg_models_pytorch.wsi.tiles import _pad_tile
from tqdm import tqdm

warnings.filterwarnings("ignore")

# we will apply some simple augmentations to the patches
img_transforms = A.Compose(
    [
        A.HorizontalFlip(p=0.33),
        A.VerticalFlip(p=0.33),
    ]
)

# we'll save the h5 file in a directory in the home folder
h5handler = H5Handler()
save_dir = Path.home() / "panoptils_refined"
save_dir.mkdir(parents=True, exist_ok=True)


# Helper function to patch the images and to create H5 files from the patches
def create_h5(df: pd.DataFrame, fold: str, stride: int, patch_size: tuple):
    fname = save_dir / f"panoptils_{fold}_p{patch_size[0]}_{stride}.h5"
    h5 = tb.open_file(fname, "w")
    h5handler.init_img(h5, patch_size)
    h5handler.init_mask(h5, "inst_mask", patch_size)
    h5handler.init_mask(h5, "type_mask", patch_size)
    h5handler.init_mask(h5, "sem_mask", patch_size)

    try:
        for i, row in tqdm(
            df.iterrows(), total=len(df), desc=f"Processing {fold} patches"
        ):
            image = png_bytes_to_array(row["image"])
            inst_mask = png_bytes_to_array(row["inst"])
            type_mask = png_bytes_to_array(row["type"])
            sem_mask = png_bytes_to_array(row["sem"])

            data = np.concatenate(
                [
                    image,
                    inst_mask[..., None],
                    type_mask[..., None],
                    sem_mask[..., None],
                ],
                axis=2,
            )
            slices, pady, padx = Predictor._get_slices(
                stride, patch_size, image.shape[:2]
            )

            padx, modx = divmod(padx, 2)
            pady, mody = divmod(pady, 2)
            padx += modx
            pady += mody

            for yslice, xslice in slices:
                x_i = _pad_tile(data[yslice, xslice, ...], shape=patch_size, fill=0)
                x_i = img_transforms(image=x_i)("image")  # apply flips
                im_i = x_i[..., :3]
                inst_i = x_i[..., 3]
                type_i = x_i[..., 4]
                sem_i = x_i[..., 5]
                h5handler.append_array(h5, im_i[None, ...], "image")
                h5handler.append_array(h5, inst_i[None, ...], "inst_mask")
                h5handler.append_array(h5, type_i[None, ...], "type_mask")
                h5handler.append_array(h5, sem_i[None, ...], "sem_mask")

        h5.close()
    except Exception as e:
        h5.close()
        raise e


# Create H5 files for train and validation sets
create_h5(df_train, "train", stride=224, patch_size=(224, 224))
create_h5(df_val, "val", stride=224, patch_size=(224, 224))
```

--------------------------------

### Fit Rectangular Grid to Tissue Segmentation

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/collagen_orientation.ipynb

Fit a rectangular grid to the tissue segmentation map. This example uses 256x256 sized patches and focuses analysis on stromal regions. The grid is visualized on top of the tissue segmentation.

```python
from histolytics.spatial_ops.rect_grid import rect_grid

# fit the grid. We'll use 256, 256 sized patches, focus only on stroma
patch_size = (256, 256)
gr = rect_grid(tis, patch_size, 0, "intersects")
ax = tis.plot(column="class_name", figsize=(10, 10), alpha=0.5, legend=True)
ax = gr.plot(ax=ax, edgecolor="black", facecolor="none", lw=0.5)
ax.set_axis_off()
```

--------------------------------

### Initialize WsiPanopticSegmenter

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/panoptic_segmentation.ipynb

Initialize the WsiPanopticSegmenter with a SlideReader, model, and image transformations. The MinMaxNormalization transform is recommended as it was likely used during model training.

```python
import albumentations as A
from histolytics.transforms import MinMaxNormalization
from histolytics.wsi.wsi_segmenter import WsiPanopticSegmenter

# We will use the MinMaxNormalization transform to normalize the input images
# This was also used in the training of the model
img_transforms = A.Compose([MinMaxNormalization()])

wsi_segmenter = WsiPanopticSegmenter(
    reader=reader,
    model=model,
    level=0,  # 20x magnification
    coordinates=sub_coords,
    batch_size=8,
    transforms=img_transforms,
)
```

--------------------------------

### FocalLoss Class

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/api/losses/focal.md

The FocalLoss class is designed to address class imbalance in deep learning models by down-weighting the loss assigned to well-classified examples. This helps the model focus more on hard-to-classify examples.

```APIDOC
## class histolytics.losses.FocalLoss

### Description

Focal Loss is a loss function that aims to solve class imbalance by down-weighting the loss assigned to easy examples, thereby focusing training on hard negatives.

### Parameters

*   **alpha** (float) - Default: 0.25 - Weighting factor for the loss.
*   **gamma** (float) - Default: 2.0 - Focusing parameter. Higher values give more weight to hard examples.
*   **reduction** (str) - Default: 'mean' - Specifies the reduction to apply to the output: 'none' | 'mean' | 'sum'.
*   **ignore_index** (int) - Default: -100 - Specifies a target value that is ignored and does not contribute to the input gradient.

### Methods

#### `forward`(input: torch.Tensor, target: torch.Tensor) -> torch.Tensor

Calculates the focal loss.

*   **input** (torch.Tensor) - The model predictions.
*   **target** (torch.Tensor) - The ground truth labels.
*   **Returns** (torch.Tensor) - The calculated focal loss.
```

--------------------------------

### Get Connected Components

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/api/spatial_graph/connected_components.md

Retrieves the connected components of a spatial graph. Ensure the graph object is properly initialized before calling this function.

```python
from histolytics.spatial_graph.utils import get_connected_components

# Assuming 'graph' is a pre-initialized spatial graph object
connected_components = get_connected_components(graph)
```

--------------------------------

### Initialize CellVit with SAM Huge Backbone

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Create a CellVitPanoptic model utilizing the SAM huge backbone. Customize the encoder's output layers via model_kwargs.

```python
from histolytics.models.cellvit_panoptic import CellVitPanoptic

cvit_panoptic = CellVitPanoptic(
    n_nuc_classes=6,
    n_tissue_classes=6,
    enc_name="samvit_huge_patch16",
    enc_pretrain=True,
    model_kwargs={
        "enc_out_indices": (2, 4, 6, 8),  # using layers 2, 4, 6, 8 from SAM Huge
    },
)

# Print only a summary of the model instead of the full output
print(str(cvit_panoptic.model)[:500] + "\n...")
```

--------------------------------

### Initialize Slide Reader

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/chromatin_patterns.ipynb

Initializes the SlideReader for accessing WSI data. Specify the correct slide path and choose a backend (e.g., 'CUCIM').

```python
from histolytics.wsi.slide_reader import SlideReader

sl_p = "/path/to/slide.tiff"  # <- modify path
reader = SlideReader(sl_p, backend="CUCIM")

# thumbnail = reader.read_level(-1)
# reader.get_annotated_thumbnail(thumbnail, gr)
```

--------------------------------

### Load and Visualize Cervix Data

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/querying.ipynb

Loads cervix tissue and nuclei data using Histolytics and visualizes them with Matplotlib. Ensure you have the necessary libraries installed.

```python
import matplotlib.pyplot as plt
from histolytics.data import cervix_nuclei, cervix_tissue

tis = cervix_tissue()
nuc = cervix_nuclei()

fig, ax = plt.subplots(figsize=(10, 10))
tis.plot(ax=ax, column="class_name", aspect=1, alpha=0.5, legend=True)
nuc.plot(ax=ax, column="class_name", aspect=1, legend=False)
ax.set_axis_off()
```

--------------------------------

### Get Interfaces Function Signature

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/api/spatial_ops/get_interfaces.md

Displays the function signature for get_interfaces, including type hints. This is useful for understanding the expected input parameters and return types.

```python
def get_interfaces(image: np.ndarray, mask: np.ndarray, radius: int = 1) -> list[tuple[tuple[int, int], tuple[int, int]]]
```

--------------------------------

### Visualize Tiles Selected by Annotation

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/selecting_tiles.ipynb

Visualizes the tiles selected by a bounding box annotation on a thumbnail of the WSI. This helps in verifying the annotation-based tile selection.

```python
reader.get_annotated_thumbnail(thumbnail, coordinates, linewidth=3)
```

--------------------------------

### Initialize StarDist Panoptic with Prov-GigaPath Backbone

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Initializes a StarDistPanoptic model using the Prov-GigaPath backbone from Hugging Face Hub. This snippet shows how to specify encoder output indices for the backbone.

```python
from histolytics.models.stardist_panoptic import StarDistPanoptic

sdist_panoptic = StarDistPanoptic(
    n_nuc_classes=6,
    n_tissue_classes=6,
    enc_name="hf_hub:prov-gigapath/prov-gigapath",
    enc_pretrain=True,
    model_kwargs={
        "enc_out_indices": (2, 4, 6, 8),  # using layers 2, 4, 6, 8 from Prov-GigaPath
    },
)

# Print only a summary of the model instead of the full output
print(str(sdist_panoptic.model)[:500] + "\n...")
```

--------------------------------

### Adjusting Tissue Mask Threshold

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/selecting_tiles.ipynb

Demonstrates adjusting the tissue mask threshold for a different slide (Philips-1.tiff). This is useful when the default Otsu threshold is not optimal, for example, with fatty tissue.

```python
slide = in_dir / "Philips-1.tiff"
reader = SlideReader(slide, backend="CUCIM")

thresh, tissue_mask = reader.get_tissue_mask(level=-2)
coordinates = reader.get_tile_coordinates(width=1024, tissue_mask=tissue_mask)

thumbnail = reader.read_level(-2)
anno_thumb = reader.get_annotated_thumbnail(thumbnail, coordinates, linewidth=2)
anno_thumb
```

--------------------------------

### Initialize HoverNet Panoptic with Virchow Backbone

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Initializes a HoverNetPanoptic model using the Virchow backbone from Hugging Face Hub. This snippet demonstrates setting encoder-specific keyword arguments for layers and activation functions.

```python
import torch
from timm.layers import SwiGLUPacked
from histolytics.models.hovernet_panoptic import HoverNetPanoptic

hnet_panoptic = HoverNetPanoptic(
    n_nuc_classes=6,
    n_tissue_classes=6,
    enc_name="hf-hub:paige-ai/Virchow",
    enc_pretrain=True,
    model_kwargs={
        "encoder_kws": {"mlp_layer": SwiGLUPacked, "act_layer": torch.nn.SiLU},
        "enc_out_indices": (2, 4, 6, 8),  # using layers 2, 4, 6, 8 from Virchow
    },
)

# Print only a summary of the model instead of the full output
print(str(hnet_panoptic.model)[:500] + "\n...")
```

--------------------------------

### Standalone Legendgram Plot

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/legendgram.ipynb

This snippet demonstrates how to generate a legendgram as a separate figure, independent of other Axes objects. It allows for customization of the colormap, using 'turbo' as an example. Ensure matplotlib is imported for plotting.

```python
import matplotlib.pyplot as plt

ax = legendgram(
    gdf=h3_res10_stroma,
    column="immune_connective_ratio_qnorm",
    n_bins=30,
    cmap="turbo",
)
ax
```

--------------------------------

### Initialize CellposePanoptic Model and Training Components

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/finetuning.ipynb

Initializes the CellposePanoptic model, defines class mappings, sets up training and validation datasets and dataloaders, configures multi-task loss functions with weights, and initializes the AdamW optimizer and CosineAnnealingLR scheduler.

```python
import torch.nn as nn
from torch.optim.lr_scheduler import CosineAnnealingLR

from histolytics.torch_datasets.h5dataset import DatasetH5
from histolytics.losses import CELoss, DiceLoss, JointLoss, MultiTaskLoss
from histolytics.models.cellpose_panoptic import CellposePanoptic
from histolytics.transforms import (
    AlbuStrongAugment,
    CellposeTransform,
    MinMaxNormalization,
    ApplyEach,
)

# panoptils (refined) classes
nuclei_classes = {
    "background": 0,
    "neoplastic": 1,
    "stromal": 2,
    "inflammatory": 3,
    "epithelial": 4,
    "other": 5,
    "unknown": 6,
}

tissue_classes = {
    "background": 0,
    "tumor": 1,
    "stroma": 2,
    "epithelium": 3,
    "junk/debris": 4,
    "blood": 5,
    "other": 6,
}


# Quick wrapper for MSE loss to make it fit the JointLoss API
class MSELoss(nn.Module):
    def __init__(self, **kwargs) -> None:
        super().__init__()

    def forward(
        self,
        yhat: torch.Tensor,
        target: torch.Tensor,
        **kwargs
    ) -> torch.Tensor:
        return F.mse_loss(yhat, target, reduction="mean").float()

batch_size = 8
in_keys = ("image", "sem_mask", "inst_mask", "type_mask")
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")

# init  CellposePanoptic model
cpose_panoptic = CellposePanoptic(
    n_nuc_classes=len(nuclei_classes),
    n_tissue_classes=len(tissue_classes),
    enc_name="efficientnet_b5",  # you can change this to any supported encoder
    enc_pretrain=True,
)
model = cpose_panoptic.model.to(device)

# init training transforms
img_transforms = A.Compose(
    [
        A.HorizontalFlip(p=0.33),
        A.VerticalFlip(p=0.33),
        AlbuStrongAugment(),
        MinMaxNormalization(),
    ]
)

inst_transforms = ApplyEach([CellposeTransform(deduplicate=False)], as_list=True)

# init val dataset and dataloader
ds_train = DatasetH5(
    Path().home() / "panoptils_refined" / "panoptils_train_p224_224.h5",
    img_key="image",
    inst_keys=["inst_mask"],
    mask_keys=["type_mask", "sem_mask"],
    transforms=img_transforms,
    inst_transforms=inst_transforms,
)

loader_train = NodesDataLoader(
    ds_train,
    batch_size,
    shuffle=False,
    num_workers=batch_size,
    pin_memory=True,
    drop_last=True,
)
# init train dataset and dataloader
ds_val = DatasetH5(
    Path().home() / "panoptils_refined" / "panoptils_val_p224_224.h5",
    img_key="image",
    inst_keys=["inst_mask"],
    mask_keys=["type_mask", "sem_mask"],
    transforms=img_transforms,
    inst_transforms=inst_transforms,
)

loader_val = NodesDataLoader(
    ds_val,
    batch_size,
    shuffle=True,
    num_workers=batch_size,
    pin_memory=True,
    drop_last=True,
)

# init training multi-task loss
# The losses are defined for each task, and the weights are set for each task.
losses = {
    "cellpose": JointLoss([MSELoss()]),
    "type_mask": JointLoss([CELoss(apply_sd=True), DiceLoss()]),
    "sem_mask": JointLoss([CELoss(apply_sd=True), DiceLoss()]),
}
loss_weights = {"cellpose": 1.0, "type_mask": 1.0, "sem_mask": 2.0}
multitask_loss = MultiTaskLoss(losses, loss_weights=loss_weights).to(device)

# optimizer and lr scheduler
optimizer = torch.optim.AdamW(params=model.parameters(), lr=0.0001)
lr_scheduler = CosineAnnealingLR(optimizer, T_max=1000, eta_min=0, last_epoch=-1)

```

--------------------------------

### Cross-Entropy Loss (CELoss) Initialization

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/api/losses/ce.md

Initializes the Cross-Entropy Loss function. Ensure that the input and target shapes are compatible for loss calculation.

```python
from histolytics.losses import CELoss

# Initialize CELoss with default parameters
celoss = CELoss()

# Initialize CELoss with specific parameters (example)
# celoss = CELoss(weight=torch.tensor([1.0, 2.0, 0.5]))
```

--------------------------------

### Initialize CPPNet with SAM Backbone

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Instantiate a CPPNetPanoptic model using a specified SAM encoder. Configure encoder output indices to control feature extraction layers.

```python
from histolytics.models.cppnet_panoptic import CPPNetPanoptic

cpp_panoptic = CPPNetPanoptic(
    n_nuc_classes=6,
    n_tissue_classes=6,
    enc_name="samvit_base_patch16",
    enc_pretrain=True,
    model_kwargs={
        "enc_out_indices": (2, 4, 6, 8),  # using layers 2, 4, 6, 8 from SAM
    },
)

# Print only a summary of the model instead of the full output
print(str(cpp_panoptic.model)[:500] + "\n...")
```

--------------------------------

### Query Inflammatory Nuclei within Stroma

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/querying.ipynb

Queries for inflammatory nuclei strictly contained within the stromal tissue region using the `get_objs` function. This showcases another example of spatial filtering based on containment.

```python
from histolytics.spatial_ops import get_objs

# get the stroma tissue and inflammatory nuclei
stroma_tissue = tis[tis["class_name"] == "stroma"]
infl_nuclei = nuc[nuc["class_name"] == "inflammatory"]

# select all the nuclei contained within stroma tissue
infl_within_stroma = get_objs(stroma_tissue, infl_nuclei, predicate="contains")

ax = tis.plot(figsize=(10, 10), column="class_name", aspect=1, alpha=0.5, legend=False)
infl_within_stroma.plot(ax=ax, column="class_name", aspect=1, legend=True)
ax.set_axis_off()
infl_within_stroma.head()
```

--------------------------------

### Initialize SlideReader with CUCIM Backend

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/selecting_tiles.ipynb

Initializes the SlideReader for a given slide file using the CUCIM backend. Ensure the slide file and directory exist.

```python
from histolytics.wsi.slide_reader import SlideReader
from pathlib import Path

in_dir = Path().home() / "openslide_test"
slide = in_dir / "CMU-1-JP2K-33005.svs"
reader = SlideReader(slide, backend="CUCIM")
```

--------------------------------

### Visualize Patches from H5 Train Dataset

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/finetuning.ipynb

Opens an H5 file containing training patches and visualizes a random sample of images along with their instance, type, and semantic masks. Requires matplotlib and PyTables.

```python
import matplotlib.pyplot as plt

# Open the H5 file for train patches
h5_train_path = save_dir / "panoptils_train_p224_224.h5"
with tb.open_file(h5_train_path, "r") as h5:
    n_patches = h5.root.image.shape[0]
    idxs = np.random.choice(n_patches, size=5, replace=False)
    fig, axes = plt.subplots(5, 4, figsize=(16, 16))
    for i, idx in enumerate(idxs):
        image = h5.root.image[idx]
        inst = h5.root.inst_mask[idx]
        type_ = h5.root.type_mask[idx]
        sem = h5.root.sem_mask[idx]
        axes[i, 0].imshow(image.astype(np.uint8))
        axes[i, 1].imshow(label2rgb(inst, bg_label=0))
        axes[i, 2].imshow(label2rgb(type_, bg_label=0))
        axes[i, 3].imshow(label2rgb(sem, bg_label=0))
        for j in range(4):
            axes[i, j].axis("off")
    plt.tight_layout()
    plt.show()
```

--------------------------------

### Login to Hugging Face Hub

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/backbones.ipynb

Use this snippet to log in to the Hugging Face Hub to load model weights. Ensure you have the necessary permissions for the weights before execution.

```python
from huggingface_hub import login

# login to huggingface to load the weights
# NOTE: You need to have granted permission for the weights before you can run this
login()
```

--------------------------------

### Benchmark Rectangular Grid Partitioning (512x512)

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/partitioning.ipynb

Measures the average time to partition stroma using a rectangular grid with 512x512 resolution and no overlap. This demonstrates performance with higher rectangular resolutions.

```python
%%timeit
rect_grid(stroma, resolution=(512, 512), overlap=0, predicate="intersects")
```

--------------------------------

### Extract Patch Images from WSI

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/chromatin_patterns.ipynb

Defines functions to convert polygon coordinates to xywh format and extract image regions from a Whole Slide Image (WSI) reader for given patches. This is used to get image data for selected low and high chromatin patches.

```python
def polygon_to_xywh(polygon):
    """Convert polygon to xywh coordinates for slide reader."""
    minx, miny, maxx, maxy = polygon.bounds
    return (int(minx), int(miny), int(maxx - minx), int(maxy - miny))

def extract_patch_images(patches, reader):
    """Extract images from WSI for given patches."""
    images = []
    patch_info = []

    for idx, _ in patches.iterrows():
        original_geom = gr.loc[idx, "geometry"]
        x, y, w, h = polygon_to_xywh(original_geom)
        img = reader.read_region((int(x), int(y), int(w), int(h)), 0)
        images.append(img)

    return images, patch_info


# Extract images for low and high chromatin patches
low_images, low_info = extract_patch_images(low_patches, reader)
high_images, high_info = extract_patch_images(high_patches, reader)
```

--------------------------------

### Extract H&E Patch Images

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/collagen_orientation.ipynb

Defines helper functions to convert polygon coordinates to xywh format and extract images from a whole slide image (WSI) reader for given patches. This is used to get the actual H&E image data for visualization.

```python
def polygon_to_xywh(polygon):
    """Convert polygon to xywh coordinates for slide reader."""
    minx, miny, maxx, maxy = polygon.bounds
    return (int(minx), int(miny), int(maxx - minx), int(maxy - miny))


def extract_patch_images(patches, reader):
    """Extract images from WSI for given patches."""
    images = []
    patch_info = []

    for idx, _ in patches.iterrows():
        original_geom = gr.loc[idx, "geometry"]
        x, y, w, h = polygon_to_xywh(original_geom)
        img = reader.read_region((int(x), int(y), int(w), int(h)), 0)
        images.append(img)

    return images, patch_info


# Extract images for low and high chromatin patches
low_images, low_info = extract_patch_images(low_patches, reader)
high_images, high_info = extract_patch_images(high_patches, reader)
```

--------------------------------

### Benchmark Rectangular Grid Partitioning (256x256)

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/partitioning.ipynb

Measures the average time to partition stroma using a rectangular grid with 256x256 resolution and no overlap. This shows the speed of rectangular partitioning.

```python
%%timeit
rect_grid(stroma, resolution=(256, 256), overlap=0, predicate="intersects")
```

--------------------------------

### Compute RGB Intensity Features

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/nuclear_features.ipynb

Computes intensity metrics separately for each RGB channel (Red, Green, Blue) of an H&E image using a provided instance mask. This function also supports GPU acceleration. It requires the same setup as grayscale feature extraction.

```python
# compute the same features for every rgb channel separately
rgb_intensity_feats = rgb_intensity_feats(he, inst_mask, metrics=metrics, device="cuda")
rgb_intensity_feats.head(3)
```

--------------------------------

### Classify Immune Clusters by Area and Density

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/tls_lymphoid_aggregate.ipynb

This function classifies clusters into 'Potential TLS', 'Large Aggregate', 'Dense Aggregate', or 'Small Aggregate' based on user-defined area and density thresholds. It uses the median neighborhood distance as the default density threshold if none is provided. Ensure you have 'numpy' and 'matplotlib' installed.

```python
import matplotlib.pyplot as plt
import numpy as np


# Add a simple classification based on size and density
def classify_clusters(df, area_threshold=2500000, density_threshold=None):
    """Simple heuristic classification"""
    if density_threshold is None:
        density_threshold = df["nhood_dists_mean"].median()

    conditions = [
        (df["area"] > area_threshold) & (df["nhood_dists_mean"] < density_threshold),
        (df["area"] > area_threshold) & (df["nhood_dists_mean"] >= density_threshold),
        (df["area"] <= area_threshold) & (df["nhood_dists_mean"] < density_threshold),
        (df["area"] <= area_threshold) & (df["nhood_dists_mean"] >= density_threshold),
    ]

    choices = ["Potential TLS", "Large Aggregate", "Dense Aggregate", "Small Aggregate"]

    df["cluster_type"] = np.select(conditions, choices, default="Unclassified")
    return df


# Apply classification
clust_feats_classified = classify_clusters(clust_feats.copy())

# Plot with classification
plt.figure(figsize=(10, 6))
for cluster_type in clust_feats_classified["cluster_type"].unique():
    mask = clust_feats_classified["cluster_type"] == cluster_type
    plt.scatter(
        clust_feats_classified.loc[mask, "area"],
        clust_feats_classified.loc[mask, "nhood_dists_mean"],
        label=cluster_type,
        s=100,
        alpha=0.7,
    )

plt.xlabel("Cluster Area (μm²)")
plt.ylabel("Mean Neighborhood Distance")
plt.title("Immune Cluster Classification\n(Based on Area and Density)")
plt.legend()
plt.grid(True, alpha=0.3)

# Add threshold lines
plt.axvline(x=2500000, color="red", linestyle="--", alpha=0.5, label="Area Threshold")
plt.axhline(
    y=clust_feats["nhood_dists_mean"].median(),
    color="blue",
    linestyle="--",
    alpha=0.5,
    label="Density Threshold",
)

plt.tight_layout()
plt.show()

# Print classification results
print("\nCluster Classification Results:")
print(clust_feats_classified["cluster_type"].value_counts())
```

--------------------------------

### Visualize Selected Sub-grid with Annotation

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/selecting_tiles.ipynb

Visualizes a selected sub-grid of tiles on a downsampled thumbnail of the WSI. Uses PIL for image manipulation and requires a SlideReader instance.

```python
from PIL import Image


#  helper downsampling function for visualization
def downsample_thumbnail(annotated_thumbnail, scale_factor=0.3):
    new_width = int(annotated_thumbnail.width * scale_factor)
    new_height = int(annotated_thumbnail.height * scale_factor)
    return annotated_thumbnail.resize((new_width, new_height), Image.LANCZOS)


# pick the first sub-grid
sub_coords = filtered_sub_grids[0]

thumbnail = reader.read_level(-1)
anno_thumb = reader.get_annotated_thumbnail(thumbnail, sub_coords, linewidth=3)
downsample_thumbnail(anno_thumb)
```

--------------------------------

### Define Desmoplastic Stroma Feature Extraction Pipeline

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/workflows/clustering_desmoplasia.ipynb

This Python function defines a pipeline to extract collagen and intensity features from WSI stromal patches. It calculates the standard deviation of collagen fiber orientation and the mean tortuosity, along with mean RGB intensities for hematoxylin and eosin stains. Ensure necessary libraries like numpy, pandas, and histolytics are installed.

```python
import numpy as np
import pandas as pd
from histolytics.wsi.wsi_iterator import WSIPatchIterator
from histolytics.stroma_feats.collagen import fiber_feats
from histolytics.stroma_feats.intensity import stromal_intensity_feats
from histolytics.wsi.slide_reader import SlideReader


def pipeline(img: np.ndarray, label: np.ndarray, mask: np.ndarray) -> pd.DataFrame:
    """A pipeline for extracting features from WSI stromal patches."""
    collagen_feats: gpd.GeoDataFrame = fiber_feats(
        img,
        label=label,
        mask=mask,
        metrics=("major_axis_angle", "tortuosity"),
        rm_bg=True,
        rm_fg=False,
        return_edges=False,
    )

    if len(collagen_feats) < 3:
        fiber_orientation_std = 0
        tortuosity_mean = 0
    else:
        fiber_orientation_std = collagen_feats["major_axis_angle"].std()
        tortuosity_mean = collagen_feats["tortuosity"].mean()

    collagen_feats = pd.Series(
        {
            "fiber_orientation_std": fiber_orientation_std,
            "tortuosity_mean": tortuosity_mean,
        }
    )

    stromal_feats: pd.Series = stromal_intensity_feats(
        img, label=label, mask=mask, metrics=["mean"]
    )
    combined_feats = pd.concat([collagen_feats, stromal_feats], axis=0)
    return combined_feats


sl_p = "/path/to/slide.mrxs"  # <- modify path
reader = SlideReader(sl_p, backend="OPENSLIDE")

crop_loader = WSIPatchIterator(
    slide_reader=reader,
    grid=gr,
    nuclei=nuc,
    pipeline_func=pipeline,
    batch_size=8,
    num_workers=8,
    pin_memory=False,
    shuffle=False,
    drop_last=False,
)
```

--------------------------------

### Benchmark Grid Aggregation (4 Cores)

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/partitioning.ipynb

Measures the time to aggregate grid data using four cores. This demonstrates the performance improvement gained by increasing the number of parallel processes.

```python
%%timeit
# four cores
grid_aggregate(
    objs=nuc,
    grid=h3_res10_inter,
    metric_func=immune_connective_ratio,
    new_col_names=["immune_connective_ratio"],
    predicate="intersects",
    num_processes=4,
)
```

--------------------------------

### Fit Spatial Graph and Visualize Neighborhoods

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/nhoods.ipynb

This snippet demonstrates how to fit a spatial graph using Delaunay triangulation and visualize the resulting neighborhood relationships. It requires the histolytics library and a GeoDataFrame with unique identifiers.

```python
from histolytics.spatial_graph.graph import fit_graph
from histolytics.utils.gdf import set_uid
from histolytics.data import hgsc_cancer_nuclei

nuc = hgsc_cancer_nuclei()
nuc = set_uid(nuc)  # Ensure the GeoDataFrame has a unique identifier

# fit spatial weights using Delaunay triangulation
w, w_gdf = fit_graph(nuc, "delaunay", id_col="uid", threshold=100)

ax = nuc.plot(figsize=(10, 10), column="class_name", aspect=1)
w_gdf.plot(ax=ax, linewidth=1, column="class_name", legend=True, aspect=1)
ax.set_axis_off()
w_gdf.head(5)
```

--------------------------------

### Print Python Version

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/seg/selecting_tiles.ipynb

Prints the current Python version. This is a basic utility snippet.

```python
from platform import python_version

print("python version:", python_version())
```

--------------------------------

### Benchmark H3 Grid Partitioning (Resolution 10)

Source: https://github.com/hautaniemilab/histolytics/blob/main/docs/user_guide/spatial/partitioning.ipynb

Measures the average time to partition stroma into a hexagonal grid at resolution 10. This demonstrates increased computation time at higher resolutions.

```python
%%timeit
h3_grid(stroma, resolution=10)
```