---
title: "Imaging data analysis test 2025-05-28"
author: "Burnet Bioinformatics group"
date: "`r Sys.Date()`"
output:
  html_document:
    toc: true
    toc_float: true
    fig_width: 7
    fig_height: 5
theme: cosmo
---

## Introduction

1.	File (path to *sis)
2.	ImageSet (image batch *czi)
3.	ROI (blank)
4.	Type (Segment)
5.	Tags (Mito Child, Nucleus Child, PV1 child, RBC, RBC area Filter, RBC segmented, Touching Edge Filter RBC)
6.	Image Set
7.	Id (C3, Field P1 - Experiment-161.czi)
8.	Name (Segment #001 (RBC segmented))
9.	Position, Well (well C3, C4)
10.	Field, Well
11.	Area, 2D Oriented Bounds (¬µm¬≤)
12.	no. Children
13.	Children Ids
14.	Children Names
15.	Parent Ids
16.	Parent Names
17.	Min, Intensities #1
18.	Max, Intensities #1
19.	Mean, Intensities #1
20.	Sum, Intensities #1
21.	SD, Intensities #1
22.	SNR (Mean/SD), Intensities #1
23.	Min, Intensities #2
24.	Max, Intensities #2
25.	Mean, Intensities #2
26.	Sum, Intensities #2
27.	SD, Intensities #2
28.	SNR (Mean/SD), Intensities #2
29.	Min, Intensities #3
30.	Max, Intensities #3
31.	Mean, Intensities #3
32.	Sum, Intensities #3
33.	SD, Intensities #3
34.	SNR (Mean/SD), Intensities #3
35.	Min, Intensities #4
36.	Max, Intensities #4
37.	Mean, Intensities #4
38.	Sum, Intensities #4
39.	SD, Intensities #4
40.	SNR (Mean/SD), Intensities #4

Tags entries (E3):

* iRBC filter (child only), Mito Child, Nucleus Child, PV1 child, RBC, RBC area Filter, RBC segmented, Touching Edge Filter RBC (512)
* Mito Child, Nucleus Child, PV1 child, RBC, RBC area Filter, RBC segmented, Touching Edge Filter RBC (1440)
* Mito, Mito Child, Mito Segmenter, Mito size filter, Mito within RBC (594)
* Nuclei, Nuclei Segmenter, Nucleus Child, Nucleus size filter, Nucleus within RBC (421)
* PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC inside (521)
* PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC intersecting (49)

The desired format is:

1.	Treatment group (Well C3,D3)
2.	RBC area um^2 (#1, mean)
3.	RBC count (number of entries with (#1, no. entries)
4.	iRBC area um^2 (#2, mean)
5.	iRBC count (#2, no. entries)
6.	PV1 area um^2 (#2, mean)
7.	PV1 intensity RLU (#2, mean)
8.	PV1 count (#2, no. entries)
9.	Mitochondria area um^2 (#3)
10.	Mitochondria intensity RLU (#3)
11.	Mitochondria count (#3)
12.	Nucleus area um^2 (#4)
13.	Nucleus intensity RLU (#4)
14.	Nucleus count (#4)

```{r,setup}

library("kableExtra")
library("DT")

```

## Example with well C3

```{r,read1}

x <- read.csv("20250701_COMBINED_Drug_pair_2_measurements.tsv",header=TRUE,sep="\t")

myt <- as.data.frame(table(x$Tags) )

myt |> kbl(caption="Frequency of features across all wells") |> kable_paper("hover", full_width = F)

WELL="C3"
y <- subset(x,Position..Well==WELL)

myt <- table(y$Tags)

myt |> kbl(caption="Frequency of features in well C3") |> kable_paper("hover", full_width = F)

# here is the important code for one well (C3)
rbc <- y[which(y$Parent.Names==""),]
RBC_count <- nrow(rbc)
RBC_area <- mean(rbc$Area..2D.Oriented.Bounds...µm...)

irbc <- subset(y,Tags=="iRBC filter (child only), Mito Child, Nucleus Child, PV1 child, RBC, RBC area Filter, RBC segmenter, Touching Edge Filter RBC")
iRBC_count <- nrow(irbc)
iRBC_area <- mean(irbc$Area..2D.Oriented.Bounds...µm...)

pv1 <- subset(y, Tags == "PV1 child, RBC, RBC area Filter, RBC segmenter, Touching Edge Filter RBC" | Tags =="PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC inside" | Tags == "PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC intersecting" )
PV1_area <- mean(pv1$Area..2D.Oriented.Bounds...µm...)
PV1_intensity <- mean(pv1$Mean..Intensities..2)
PV1_count <- nrow(pv1)

mito <- subset(y, Tags == "Mito Child, Mito Segmenter, Mito size filter, Mito within RBC, Mito within RBC inside, Mitochondria" | Tags == "Mito Child, Mito Segmenter, Mito size filter, Mito within RBC, Mito within RBC intersecting, Mitochondria")
Mito_area <- mean(mito$Area..2D.Oriented.Bounds...µm...)
Mito_intensity <- mean(mito$Mean..Intensities..1)
Mito_count <- nrow(mito)

nucleus <- subset(y, Tags == "Nuclei, Nuclei Segmenter, Nucleus Child, Nucleus size filter, Nucleus within RBC, Nucleus within RBC inside" | Tags == "Nuclei, Nuclei Segmenter, Nucleus Child, Nucleus size filter, Nucleus within RBC, Nucleus within RBC intersecting" )
Nucleus_area <- mean(nucleus$Area..2D.Oriented.Bounds...µm...)
Nucleus_intensity <- mean(nucleus$Mean..Intensities..3)
Nucleus_count <- nrow(nucleus)

res <- c("RBC count"=RBC_count,
  "RBC area"=RBC_area,
  "iRBC count"=iRBC_count,
  "iRBC area"=iRBC_area,
  "PV1_area"=PV1_area,
  "PV1_intensity"=PV1_intensity,
  "PV1_count"=PV1_count,
  "Mito_area"=Mito_area,
  "Mito_intensity"=Mito_intensity,
  "Mito_count"=Mito_count,
  "Nucleus_area"=Nucleus_area,
  "Nucleus_intensity"=Nucleus_intensity,
  "Nucleus_count"=Nucleus_count)

res

WELLS <- unique(x$Position..Well)

```

## Make a function and run for all wells

Now put it in a function so we can iterate over all wells.

```{r,func1}

# define the function
smry <- function(y) {
  rbc <- y[which(y$Parent.Names==""),]
  RBC_count <- nrow(rbc)
  RBC_area <- mean(rbc$Area..2D.Oriented.Bounds...µm...)

  irbc <- subset(y,Tags=="iRBC filter (child only), Mito Child, Nucleus Child, PV1 child, RBC, RBC area Filter, RBC segmenter, Touching Edge Filter RBC")
  iRBC_count <- nrow(irbc)
  iRBC_area <- mean(irbc$Area..2D.Oriented.Bounds...µm...)

  pv1 <- subset(y, Tags == "PV1 child, RBC, RBC area Filter, RBC segmenter, Touching Edge Filter RBC" | Tags =="PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC inside" | Tags == "PV1, PV1 child, PV1 Segmenter, PV1 size filter, PV1 within RBC, PV1 within RBC intersecting" )
  PV1_area <- mean(pv1$Area..2D.Oriented.Bounds...µm...)
  PV1_intensity <- mean(pv1$Mean..Intensities..2)
  PV1_count <- nrow(pv1)

  mito <- subset(y, Tags == "Mito Child, Mito Segmenter, Mito size filter, Mito within RBC, Mito within RBC inside, Mitochondria" | Tags == "Mito Child, Mito Segmenter, Mito size filter, Mito within RBC, Mito within RBC intersecting, Mitochondria")
  Mito_area <- mean(mito$Area..2D.Oriented.Bounds...µm...)
  Mito_intensity <- mean(mito$Mean..Intensities..1)
  Mito_count <- nrow(mito)

  nucleus <- subset(y, Tags == "Nuclei, Nuclei Segmenter, Nucleus Child, Nucleus size filter, Nucleus within RBC, Nucleus within RBC inside" | Tags == "Nuclei, Nuclei Segmenter, Nucleus Child, Nucleus size filter, Nucleus within RBC, Nucleus within RBC intersecting" )
  Nucleus_area <- mean(nucleus$Area..2D.Oriented.Bounds...µm...)
  Nucleus_intensity <- mean(nucleus$Mean..Intensities..3)
  Nucleus_count <- nrow(nucleus)

  res <- c("RBC count"=RBC_count,
    "RBC area"=RBC_area,
    "iRBC count"=iRBC_count,
    "iRBC area"=iRBC_area,
    "PV1_area"=PV1_area,
    "PV1_intensity"=PV1_intensity,
    "PV1_count"=PV1_count,
    "Mito_area"=Mito_area,
    "Mito_intensity"=Mito_intensity,
    "Mito_count"=Mito_count,
    "Nucleus_area"=Nucleus_area,
    "Nucleus_intensity"=Nucleus_intensity,
    "Nucleus_count"=Nucleus_count)

  res
}

# get a vector of all wells
WELLS <- unique(x$Position..Well)

# now run the function for all wells
z <- lapply(WELLS, function(w) {
  y <- subset(x,Position..Well==w)
  smry(y)
} )

# z is a list of vectors, which we can convert to a matrix
z <- do.call(rbind,z)
rownames(z) <- WELLS
# peek at the ddata format
head(z,2)

```

Show the data as a kableExtra table.

```{r,table1}

z |> kbl(caption="Summarised data by well") |> kable_paper("hover", full_width = F)

```

Show as data.table. This function allows you to download as an Excel file.

```{r,table2}

DT::datatable(
  { head(z,1000) },

  extensions = 'Buttons',
  filter = list(position = 'top', clear = FALSE),

  options = list(
    paging = TRUE,
    searching = TRUE,
    fixedColumns = TRUE,
    autoWidth = TRUE,
    ordering = TRUE,
    dom = 'tB',
    buttons = c('copy', 'csv', 'excel'),
    pageLength = 50,
    search = list(regex = TRUE, caseInsensitive = TRUE)),

  rownames= TRUE)

```

## Principal component analysis

```{r,pca1,fig.width=9,fig.height=9}

pca <- prcomp(z, scale=TRUE)

scree <- pca$sdev
names(scree) <- 1:length(scree)
barplot(scree, main="Screeplot")

biplot(pca, main="PCA")

```

More variations on biplot clusstering here: https://www.geeksforgeeks.org/how-to-create-a-biplot-in-r/

## Session information

For reproucibility

```{r,session}

sessionInfo()

```
