## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
    collapse = TRUE,
    comment = "#>"
)

## ----installation, include = TRUE, eval=FALSE---------------------------------
# if (!requireNamespace("BiocManager")) {
#     install.packages("BiocManager")
# }
# BiocManager::install("sosta")

## ----setup, message=FALSE-----------------------------------------------------
library("sosta")
library("dplyr")
library("tidyr")
library("ggplot2")
library("ggspavis")
library("imcdatasets")

## ----loading, echo=FALSE, message=FALSE---------------------------------------
# load the data
spe <- imcdatasets::Damond_2019_Pancreas("spe", full_dataset = FALSE)

## -----------------------------------------------------------------------------
plotSpots(
    spe[, spe[["patient_id"]] == 6126 &
        spe[["image_name"]] %in% c("E02", "E03", "E04")],
    annotate = "cell_category",
    sample_id = "image_number",
    in_tissue = NULL,
    y_reverse = FALSE
) + facet_wrap(~image_name)

## -----------------------------------------------------------------------------
shapeIntensityImage(
    spe,
    marks = "cell_category",
    image_col = "image_name",
    image_id = "E04",
    mark_select = "islet"
)

## -----------------------------------------------------------------------------
n <- estimateReconstructionParametersSPE(
    spe,
    marks = "cell_category",
    image_col = "image_name",
    mark_select = "islet",
    nimages = 20,
    plot_hist = TRUE
)

## -----------------------------------------------------------------------------
(thres_spe <- mean(n$thres))
(bndw_spe <- mean(n$bndw))

## -----------------------------------------------------------------------------
islet <- reconstructShapeDensityImage(
    spe,
    marks = "cell_category",
    image_col = "image_name",
    image_id = "E04",
    mark_select = "islet",
    bndw = bndw_spe,
    dim = 500,
    thres = thres_spe
)

## -----------------------------------------------------------------------------
plotSpots(
    spe[, spe[["image_name"]] %in% c("E04")],
    annotate = "cell_category",
    sample_id = "image_number",
    in_tissue = NULL,
    y_reverse = FALSE,
) +
    geom_sf(
        data = islet,
        fill = NA,
        color = "darkblue",
        inherit.aes = FALSE, # this is important
        linewidth = 0.75
    )

## -----------------------------------------------------------------------------
islet_2 <- reconstructShapeDensityImage(
    spe,
    marks = "cell_category",
    image_col = "image_name",
    image_id = "E04",
    mark_select = "islet",
    dim = 500
)

## -----------------------------------------------------------------------------
plotSpots(
    spe[, spe[["image_name"]] %in% c("E04")],
    annotate = "cell_category",
    sample_id = "image_number",
    in_tissue = NULL,
    y_reverse = FALSE,
) +
    geom_sf(
        data = islet_2,
        fill = NA,
        color = "darkblue",
        inherit.aes = FALSE,
        linewidth = 0.75
    )

## ----eval=TRUE----------------------------------------------------------------
all_islets <- reconstructShapeDensitySPE(
    spe,
    marks = "cell_category",
    image_col = "image_name",
    mark_select = "islet",
    bndw = bndw_spe,
    thres = thres_spe,
    ncores = 2
)

## -----------------------------------------------------------------------------
islet_shape_metrics <- totalShapeMetrics(all_islets)

## -----------------------------------------------------------------------------
patient_data <- colData(spe) |>
    as_tibble() |>
    group_by(image_name) |>
    select(all_of(
        c(
            "patient_stage",
            "tissue_slide",
            "tissue_region",
            "patient_id",
            "patient_disease_duration",
            "patient_age",
            "patient_gender",
            "patient_ethnicity",
            "patient_BMI",
            "sample_id"
        )
    )) |>
    unique()

## -----------------------------------------------------------------------------
all_islets <- dplyr::left_join(all_islets, patient_data, by = "image_name")
all_islets <- cbind(all_islets, t(islet_shape_metrics))

## -----------------------------------------------------------------------------
library(ggfortify)

autoplot(
    prcomp(t(islet_shape_metrics), scale. = TRUE),
    x = 1,
    y = 2,
    data = all_islets,
    color = "patient_stage",
    size = 2,
    # shape = 'type',
    loadings = TRUE,
    loadings.colour = "steelblue3",
    loadings.label = TRUE,
    loadings.label.size = 3,
    loadings.label.repel = TRUE,
    loadings.label.colour = "black"
) +
    scale_color_brewer(palette = "Dark2") +
    theme_bw() +
    coord_fixed()

## ----fig.width=8, fig.height=6------------------------------------------------
all_islets |>
    sf::st_drop_geometry() |>
    select(patient_stage, rownames(islet_shape_metrics)) |>
    pivot_longer(-patient_stage) |>
    ggplot(aes(x = patient_stage, y = value, fill = patient_stage)) +
    geom_boxplot() +
    facet_wrap(~name, scales = "free") +
    scale_fill_brewer(palette = "Dark2") +
    scale_x_discrete(guide = guide_axis(n.dodge = 2)) +
    guides(fill = "none") +
    theme_bw()

## ----sessionInfo--------------------------------------------------------------
sessionInfo()