Анализ экспортированных аннотаций .geojson из Qupath в SpatStat в R

Я новичок в программировании, имею опыт работы в области наук о жизни/биологии.

Я стремлюсь проанализировать совместную локализацию лейкозных клеток в областях кровеносных сосудов и адипоцитов.

Я экспортировал местоположения ячеек в виде точек с отметками, обозначающими либо положительный результат CKIT, либо положительный результат только DAPI.

Я экспортировал точечные узоры клеток в Qupath и аннотации расположения структур, окрашенных кровеносными сосудами/эндомуцином, и структур, окрашенных адипоцитами/перилипином, в виде файлов .geojson из Qupath.

Файл .geojson с экспортированными аннотациями Qupath

Окно должно представлять собой квадратную коробку, из которой я взял образец ткани.

Я приложил изображение примера того, как я получил точки и аннотации от qupath. Аннотации, окрашенные эндомуцином и перилипином, имеют белый и зеленый цвета соответственно и идентифицируются с помощью классификатора пикселей qupath, который я обучил. Точки — это ячейки.

Точки ячеек Qupath и аннотации

Я хотел бы, чтобы окно представляло собой квадратную внешнюю рамку, точки - в виде многотипного шаблона точек, а аннотации для окрашивания перилипином и эндомуцином - выступали в качестве пространственных ковариат для анализа шаблона точек. Однако как мне лучше всего смоделировать эти последние аннотации? После импорта их как .geojson в R и преобразования в формат SF некоторые края соприкасаются при проверке с помощью st_is_valid.

Должен ли я преобразовать это в пиксельное изображение или бинарную маску? Могут ли они остаться полигонами и быть преобразованы в объекты OWIN? (например, данные о золоте Мерчисона в замечательной книге Бэддели по анализу пространственных точечных структур). Мне не удалось преобразовать эти полигоны в объекты OWIN.

Другая ключевая часть этого заключается в том, что аннотации должны сохранять свой класс как окрашенные перилипином или эндомуцином... Возможно, это было бы проще, если бы это были пиксельные изображения с уровнями факторов?

Я попытался импортировать файлы .geojson и преобразовать их в окна на основе spatstat с помощью следующего кода:

library(sf)
library(spatstat)
library(spatstat.geom)

setwd("~/RStudio/ImageAnalysis")

pt_geojson <- "Veh Tibia - Proximal Tibia.geojson"
annotations_sf <- st_read(pt_geojson)

# Check resulting sf object

head(annotations_sf)

annotations_transform <- st_transform(annotations_sf, crs = 26717)

annotations_owin <- as.owin(annotations_transform)

plot(annotations_owin)

К сожалению, это просто генерирует следующую ошибку:

Error in if (!all(st_dimension(W) == 2)) stop("as.owin.sfc needs polygonal geometries") : 
  missing value where TRUE/FALSE needed

Я добавил следующий код:

# Check for invalid geometries
invalid_geometries <- st_is_valid(annotations_sf, reason = TRUE)
print(invalid_geometries)

# Fix invalid geometries
annotations_sf <- st_make_valid(annotations_sf)

# Verify geometries are now valid
valid_geometries <- st_is_valid(annotations_sf, reason = TRUE)
print(valid_geometries)

Это привело к следующему:

> # Check for invalid geometries
> invalid_geometries <- st_is_valid(annotations_sf, reason = TRUE)
Warning messages:
1: In st_is_longlat(x) :
  bounding box has potentially an invalid value range for longlat data
2: In st_is_longlat(x) :
  bounding box has potentially an invalid value range for longlat data
> print(invalid_geometries)
 [1] "Valid Geometry"          "Edge 2 crosses edge 4"  
 [3] "Edge 16 crosses edge 20" "Valid Geometry"         
 [5] "Valid Geometry"          "Edge 14 crosses edge 20"
 [7] "Valid Geometry"          "Valid Geometry"         
 [9] "Edge 0 crosses edge 8"   "Edge 16 crosses edge 18"
[11] "Edge 6 crosses edge 8"   "Valid Geometry"         
[13] "Valid Geometry"          "Valid Geometry"         
[15] "Edge 0 crosses edge 26"  "Edge 4 crosses edge 10" 
[17] "Valid Geometry"          "Valid Geometry"         
[19] "Edge 12 crosses edge 18" "Valid Geometry"         
[21] "Edge 0 crosses edge 10"  "Valid Geometry"         
[23] "Valid Geometry"          "Valid Geometry"         
[25] "Edge 10 crosses edge 12" "Valid Geometry"         
[27] "Valid Geometry"          "Valid Geometry"         
[29] "Valid Geometry"          "Edge 4 crosses edge 24" 
> 
> # Fix invalid geometries
> annotations_sf <- st_make_valid(annotations_sf)
Warning message:
In st_is_longlat(x) :
  bounding box has potentially an invalid value range for longlat data
> 
> # Verify geometries are now valid
> valid_geometries <- st_is_valid(annotations_sf, reason = TRUE)
> print(valid_geometries)
 [1] "Valid Geometry"          "Edge 2 crosses edge 4"  
 [3] "Edge 16 crosses edge 20" "Valid Geometry"         
 [5] "Valid Geometry"          "Edge 14 crosses edge 20"
 [7] "Valid Geometry"          "Valid Geometry"         
 [9] "Edge 0 crosses edge 8"   "Edge 16 crosses edge 18"
[11] "Edge 6 crosses edge 8"   "Valid Geometry"         
[13] "Valid Geometry"          "Valid Geometry"         
[15] "Edge 0 crosses edge 26"  "Edge 4 crosses edge 10" 
[17] "Valid Geometry"          "Valid Geometry"         
[19] "Edge 12 crosses edge 18" "Valid Geometry"         
[21] "Edge 0 crosses edge 10"  "Valid Geometry"         
[23] "Valid Geometry"          "Valid Geometry"         
[25] "Edge 10 crosses edge 12" "Valid Geometry"         
[27] "Valid Geometry"          "Valid Geometry"         
[29] "Valid Geometry"          "Edge 4 crosses edge 24"

И последующий повторный запуск сюжета создал этот ужасно объединенный сюжет всех аннотаций.

Неудачный сюжет полигонального окна

Любая помощь будет очень признательна! Спасибо :))

Редактировать

Я думаю, что нашел решение, но не уверен, вызовет ли оно проблемы в дальнейшем. Я пришел к выводу, что экспорт аннотаций в виде отдельных объединенных полигонов избавляет от многих проблем и, похоже, работает. Я также изолировал ограничивающую рамку из исходного экспорта .geojson, но в последующих попытках я упрощу это, просто экспортировав ограничивающую рамку непосредственно из qupath.

#Load required packages

library(sf)
library(spatstat)
library(spatstat.geom)

#Set working directory (CHANGE AS REQUIRED)

setwd("~/RStudio/ImageAnalysis")

#Import .geojson files, generated by Qupath

pt_geojson <- "Veh Tibia - Proximal Tibia.geojson"
pt_peri <- "pt_peri.geojson"
pt_endo <- "pt_endo.geojson"

#Generate sf data

ann_sf <- st_read(pt_geojson)
ann_peri <- st_read(pt_peri)
ann_endo <- st_read(pt_endo)

#Generate outer bounding box

bbox <- st_bbox(ann_sf)
bbox <- st_as_sfc(bbox)
bbox <- st_bbox(bbox[[1]])

xmin <- bbox[[1]]
xmax <- bbox[[3]]
ymin <- bbox[[2]]
ymax <- bbox[[4]]

W <- owin(xrange = c(xmin, xmax), yrange = c(ymin, ymax))

# Generate adipocyte and vascular annotations

peri <- ann_peri$geometry[[1]]
peri <- as.owin(peri)

endo <- ann_endo$geometry[[1]]
endo <- as.owin(endo)

com <- solist(W, endo, peri, pt_ppp)

plot(W, main = "Adipocytes + Vasculature")
plot(peri, col = "green", add = TRUE)
plot(endo, col = "red", add = TRUE)
plot(pt_ppp, cols = c("yellow3","skyblue1"), chars = c(15,16), 
     use.marks = TRUE, legend = TRUE, cex = 0.5, 
     main = "Proximal Tibia Cell Point Pattern", add = TRUE)