Skip to content

GitLab

Commit 15241bb6 authored by Stefano Beretta's avatar Stefano Beretta
Browse files

VisualZoneR analysis script

parent 43f07ff7
Show whitespace changes
Inline Side-by-side
scripts/VisualZoneR_analysis.R 0 → 100644
View file @ 15241bb6
### VisualZoneR ###
library(Seurat)
library(ggplot2)
library(patchwork)
library(dplyr)
library(cowplot)
library(OLIN)
library(fgsea)
library(gplots)
library(stringr)
# Directories
wdir <- "~/VisualZoneR"
out_dir <- paste(wdir, "results", sep = "/")
dir.create(path = out_dir, showWarnings = F)
# Read integrated Seurat object
crc.integrated <- readRDS(file = paste(out_dir, "crc.integrated.rds", sep = "/"))
for (sample in samples) {
crc.integrated <- SetIdent(object = crc.integrated, value = "orig.ident")
sub1 <- subset(crc.integrated, idents = sample)
sub1 <- SetIdent(object = sub1, value = "integrated_FinalClusters")
Coord <- merge(crc.integrated@images[[paste0(sample,"img")]]@coordinates, crc.integrated@meta.data[crc.integrated@meta.data$orig.ident == sample, "integrated_FinalClusters", drop = F], by = 0)
rownames(Coord) <- Coord$Row.names
Coord$Row.names <- NULL
Coord$ift <- ifelse(Coord$integrated_FinalClusters == 1 | Coord$integrated_FinalClusters == 6, 1, 0)
maxrow <- max(Coord$row)
maxcol <- max(Coord$col)
msub1 <- matrix(rep(NA, maxrow * maxcol), ncol = maxcol)
msub1 <- as.data.frame(msub1)
for (i in c(1:length(Coord$ift))) {
msub1[Coord$row[i], Coord$col[i]] <- Coord$ift[i]
}
msub1liv <- matrix(rep(NA, maxrow * maxcol), ncol = maxcol)
msub1liv <- as.data.frame(msub1liv)
for (i in c(1:length(Coord$ift))){
msub1liv[Coord$row[i], Coord$col[i]] <- (-Coord$ift[i]+1)
}
namesTokeep <- matrix(rep(NA, maxrow * maxcol), ncol = maxcol)
namesTokeep <- as.data.frame(namesTokeep)
for (i in c(1:length(Coord$ift))){
namesTokeep[Coord$row[i], Coord$col[i]] <- rownames(Coord)[i]
}
masub1 <- ma.matrix(as.matrix(msub1), av = "mean", delta = 2, edgeNA = FALSE )
masub1 <- ma.matrix(as.matrix(msub1), av = "mean", delta = 3, edgeNA = FALSE ) + (masub1*2)
indexes1 <- msub1
masub1 <- masub1 * indexes1
masub1liv <- ma.matrix(as.matrix(msub1liv), av = "mean", delta = 2, edgeNA = FALSE)
masub1liv <- ma.matrix(as.matrix(msub1liv), av = "mean", delta = 3, edgeNA = FALSE) * 0.1 + (masub1liv * 2)
indexes1 <- msub1
masub1liv <- masub1liv * (-1*(indexes1-1))
# Score of tumor
M1 = 2.96
M2 = 2.7
M3 = 2.3
masub1[masub1 > M1] <- (-4)
masub1[masub1 <= M1 & masub1 > M2] <- (-3)
masub1[masub1 <= M2 & masub1 > M3] <- (-2)
masub1[masub1 <= M3 & masub1 > 0] <- (-1)
# Score of liver
L1 = 2.095
L2 = 1.91
L3 = 1.7
masub1liv[masub1liv <= L3 & masub1liv > 0] <- 0
masub1liv[masub1liv <= L2 & masub1liv > L3] <- 1
masub1liv[masub1liv <= L1 & masub1liv > L2] <- 2
masub1liv[masub1liv > L1] <- 3
# Final Zones
finalScores <- masub1 + masub1liv
finalScores[finalScores == (-4)] <- "Zone_A"
finalScores[finalScores == (-3)] <- "Zone_B"
finalScores[finalScores == (-2)] <- "Zone_C"
finalScores[finalScores == (-1)] <- "Zone_D"
finalScores[finalScores == 0] <- "Zone_E"
finalScores[finalScores == 1] <- "Zone_F"
finalScores[finalScores == 2] <- "Zone_G"
finalScores[finalScores == 3] <- "Zone_H"
finalScores[finalScores == "NaN"] <- NA
FinalTable <- data.frame(SeuratID = as.vector(as.matrix(namesTokeep)), datasub1 = as.vector(as.matrix(finalScores)), drow = rep(1:nrow(masub1), ncol(masub1)), dcol = rep(1:ncol(masub1), nrow(masub1))[order(rep(1:ncol(masub1), nrow(masub1)))])
FinalTable <- na.omit(FinalTable)
SeuratPos <- names(crc.integrated@active.ident)
FinalTable <- FinalTable[order(FinalTable$SeuratID),]
if(!"Zones" %in% colnames(crc.integrated@meta.data)) {
crc.integrated@meta.data$Zones <- rep("Zone_A", length(SeuratPos))
}
crc.integrated@meta.data$Zones[crc.integrated@meta.data$orig.ident == sample] <- FinalTable$datasub1
}
full_md <- crc.integrated@meta.data
gz_out_md <- gzfile(paste(out_dir, "crc.integrated_metadata.csv.gz", sep = "/"), "w")
write.csv(x = full_md, gz_out_md)
close(gz_out_md)
saveRDS(crc.integrated, file = paste(out_dir, "crc.integrated.rds", sep = "/"))
# R version 4.0.3 (2020-10-10)
# Platform: x86_64-conda-linux-gnu (64-bit)
# Running under: Rocky Linux 8.9 (Green Obsidian)
#
# Matrix products: default
# BLAS/LAPACK: /opt/miniforge3/envs/singlecell4/lib/libopenblasp-r0.3.25.so; LAPACK version 3.11.0
#
# locale:
# [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 LC_MONETARY=en_US.UTF-8
# [6] LC_MESSAGES=en_US.UTF-8 LC_PAPER=en_US.UTF-8 LC_NAME=C LC_ADDRESS=C LC_TELEPHONE=C
# [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
#
# time zone: UTC
# tzcode source: system (glibc)
#
# attached base packages:
# [1] stats graphics grDevices utils datasets methods base
#
# other attached packages:
# [1] gplots_3.1.3 fgsea_1.28.0 reshape2_1.4.4 rayshader_0.37.1 OLIN_1.80.0 marray_1.80.0 limma_3.58.1
# [8] locfit_1.5-9.8 cowplot_1.1.2 dplyr_1.1.4 patchwork_1.2.0 ggplot2_3.4.4 SeuratObject_5.0.1 Seurat_4.1.1
#
# loaded via a namespace (and not attached):
# [1] RColorBrewer_1.1-3 rstudioapi_0.15.0 jsonlite_1.8.8 magrittr_2.0.3 spatstat.utils_3.0-4 farver_2.1.1
# [7] ragg_1.2.7 vctrs_0.6.5 ROCR_1.0-11 base64enc_0.1-3 htmltools_0.5.7 progress_1.2.3
# [13] sctransform_0.4.1 parallelly_1.36.0 KernSmooth_2.23-22 htmlwidgets_1.6.4 ica_1.0-3 plyr_1.8.9
# [19] plotly_4.10.4 zoo_1.8-12 igraph_1.6.0 mime_0.12 lifecycle_1.0.4 iterators_1.0.14
# [25] pkgconfig_2.0.3 Matrix_1.6-4 R6_2.5.1 fastmap_1.1.1 fitdistrplus_1.1-11 future_1.33.1
# [31] shiny_1.8.0 digest_0.6.34 colorspace_2.1-0 tensor_1.5 irlba_2.3.5.1 textshaping_0.3.7
# [37] labeling_0.4.3 progressr_0.14.0 fansi_1.0.6 spatstat.sparse_3.0-3 httr_1.4.7 polyclip_1.10-6
# [43] abind_1.4-5 mgcv_1.9-1 compiler_4.3.2 withr_3.0.0 doParallel_1.0.17 BiocParallel_1.36.0
# [49] MASS_7.3-60 caTools_1.18.2 gtools_3.9.5 tools_4.3.2 lmtest_0.9-40 httpuv_1.6.13
# [55] future.apply_1.11.1 goftest_1.2-3 glue_1.7.0 nlme_3.1-164 promises_1.2.1 grid_4.3.2
# [61] Rtsne_0.17 cluster_2.1.6 generics_0.1.3 gtable_0.3.4 spatstat.data_3.0-3 tidyr_1.3.0
# [67] spatstat.core_2.4-4 data.table_1.14.10 hms_1.1.3 sp_2.1-3 utf8_1.2.4 spatstat.geom_3.2-7
# [73] RcppAnnoy_0.0.21 ggrepel_0.9.5 RANN_2.6.1 foreach_1.5.2 pillar_1.9.0 stringr_1.5.1
# [79] spam_2.10-0 later_1.3.2 splines_4.3.2 lattice_0.22-5 survival_3.5-7 deldir_2.0-2
# [85] tidyselect_1.2.0 miniUI_0.1.1.1 pbapply_1.7-2 knitr_1.45 gridExtra_2.3 scattermore_1.2
# [91] xfun_0.42 statmod_1.5.0 matrixStats_1.2.0 stringi_1.8.3 lazyeval_0.2.2 codetools_0.2-19
# [97] tibble_3.2.1 cli_3.6.2 uwot_0.1.16 rpart_4.1.23 systemfonts_1.0.5 xtable_1.8-4
# [103] reticulate_1.34.0 munsell_0.5.0 Rcpp_1.0.12 globals_0.16.2 spatstat.random_3.2-2 png_0.1-8
# [109] parallel_4.3.2 ellipsis_0.3.2 rgl_1.2.13 prettyunits_1.2.0 dotCall64_1.1-1 bitops_1.0-7
# [115] listenv_0.9.0 viridisLite_0.4.2 scales_1.3.0 ggridges_0.5.5 leiden_0.4.3.1 purrr_1.0.2
# [121] crayon_1.5.2 rlang_1.1.3 fastmatch_1.1-4
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment