wonlab-cs / vesalius Goto Github PK

View Code? Open in Web Editor NEW

20.0 20.0 8.0 151.2 MB

Vesalius: Dissecting Tissue Anatomy from Spatial Transcriptomic Data

Home Page: https://patrickcnmartin.github.io/Vesalius/

License: GNU General Public License v3.0

R 98.89% C++ 1.11%

vesalius's People

Contributors

Stargazers

Watchers

Forkers

neocaleb fmz349 heidiottesen henriettaholze faragmamdouh danohdoh whistle-ch0i dataxplora

vesalius's Issues

Inserting my own log normed PCA data

Hello,

How can I replace this step

... 
bead <- new(Class = 'SlideSeq', assay = 'Spatial', coordinates = coords[,c("x","y")])
rownames(bead@coordinates) <- rownames(coords)
exp <- CreateSeuratObject(counts, assay ="Spatial")
bead <- bead[Cells(x = exp)]
DefaultAssay(object = bead) <- "Spatial"
exp[["slice1"]] <- bead

#reducedDims(exp, withDimnames=TRUE) <- list(logPCA=t(data.matrix(lpc, rownames.force=TRUE)))

exp <- NormalizeData(exp)
exp <- FindVariableFeatures(exp, nfeatures = 100)
exp <- ScaleData(exp)
#message(paste0("Seurat Pre-Processed Data: maynard_noNA_",slideTag[i]))
 
rgbu <- rgbUMAP(exp, pcs=30, conserveSparse = FALSE)
...

when I have a preprocessed matrix of log normalized PCs for my spots?

Data simulation in Vesalius paper

Dear Authors, I wonder if you could share the simulated data sets that were used in the paper or the code for generating them? Thank you!

How to extract predicted labels for barcodes/coordinates?

Hello,
I have ground truth labels for the barcodes/coordinates in a .csv file and want to compare to predicted labels of those coordinates from Vesalius after I've input the coordinate and gene expression/count matrix and formed rgb and isolated territory images. I believe the answer lies somewhere in the following example code from the tutorial but I'm not sure where in Vesalius how to actually enter in the labels and then how to read out the predicted labels?

What is the relationship between homo_cells and my labels per bead/barcode?


#------------------------------------------------------------------------------#
# looping over all cells types in homotypic beads
# Keeping everything in a neat little data frame
#------------------------------------------------------------------------------#
by_cell <- [vector](https://rdrr.io/r/base/vector.html)("list",[length](https://rdrr.io/r/base/length.html)(homo_cells))
[names](https://rdrr.io/r/base/names.html)(by_cell) <- [names](https://rdrr.io/r/base/names.html)(homo_cells)
for(i in [seq_along](https://rdrr.io/r/base/seq.html)(homo_cells)){
    [print](https://rdrr.io/r/base/print.html)([length](https://rdrr.io/r/base/length.html)(homo_cells[[i]]))
    if([length](https://rdrr.io/r/base/length.html)(homo_cells[[i]]) <30){
      next()
    } else {
      deg <- [compareCells](https://patrickcnmartin.github.io/Vesalius/reference/compareCells.html)(imageBrain,
                          brainNorm,
                          homo_cells[[i]],
                          method ="wilcox",
                          seed=cortex,
                          logFC = 0.25,
                          pval = 0.01,
                          query=thalamus,
                          minBar = 30)
    }
    if([is.null](https://rdrr.io/r/base/NULL.html)(deg)){next()}
    deg$celltype <- [rep](https://rdrr.io/r/base/rep.html)([names](https://rdrr.io/r/base/names.html)(homo_cells)[i], [nrow](https://rdrr.io/r/base/nrow.html)(deg))
    deg$group_1 <- [rep](https://rdrr.io/r/base/rep.html)("Cortex", [nrow](https://rdrr.io/r/base/nrow.html)(deg))
    deg$group_2 <- [rep](https://rdrr.io/r/base/rep.html)("Thalamus", [nrow](https://rdrr.io/r/base/nrow.html)(deg))
    by_cell[[i]] <- deg
}
by_cell <- [do.call](https://rdrr.io/r/base/do.call.html)("rbind",by_cell)
[write.csv](https://rdrr.io/r/utils/write.table.html)(by_cell, file = "vesalius_terComp_CortexVSThalamus.csv")

Many thanks

Continuation of previous issue #10

Hello,

Thanks for your help!

I have ground truth labels for each coordinate, and I have 3,611 coordinates that I input into Vesalius. So I presume I want to do the first option.
But the territory column from isolateTerritories.array has 138,440 rows!

I'm using the processing steps in the first vignette leading up to following code where I extract territory as suggested.

imageBrain <- isolateTerritories.array(imageBrain,
                                       captureRadius = 0.008,
                                       minBar = 50)

write.csv(imageBrain$territory, "y_pred.csv")

"y_pred.csv" has the following with numbers that just increase incrementally to 138,440.

"","x"
"1","isolated"
"2","isolated"
"3","isolated"
"4","isolated"
"5","isolated"
"6","isolated"
"7","isolated"
"8","isolated"
"9","isolated"
"10","isolated"
"11","isolated"
"12","isolated"
"13","isolated"
"14","isolated"
"15","isolated"
"16","isolated"
"17","isolated"
"18","isolated"
"19","isolated"
...

My "ground_labels.csv" has the following 3,611 labels for each coordinate.

I need the 3,611 predictions for each coordinate.

Many thanks

Originally posted by @akshayMpatel in #10 (comment)

'arguments imply differing number of rows' Error during Converting Sparse Matrix to Matrix

Using my own data I get:

2022-09-12 08:28:38  Converting Sparse Matrix to Matrix 
Error in data.frame(rownames(coordinates), coordinates[, c("x", "y")],  : 
  arguments imply differing number of rows: 0, 33538

My UMImatrix has 33538 columns and 3566 rows. Given it doesn't have any rownames or header, just counts.

The error seems to eminate from

count <- rgbUMAP(count, pcs =30, conserveSparse = FALSE)
Trying to reproduce the issue with the vignette with Slide-seqV2 data:

Loading Bead file Puck_200115_08
Error in read.table(file = file, header = header, sep = sep, quote = quote,  : 
  more columns than column names

My coordinate .csv file has the same format for Slide-seqV2 coordinates file, including 'barcodes,xcoord,ycoord' header. However for my coordinates file the barcodes column consist of '1' to '3566' instead of the barcode letter sequences.

invalid territory.pdf file

My territory.pdf file seems to be corrupt after trying the vesalius vignette on a slices of the the Maynard et al 2020 dorsal lateral PFC data set. I've used the following code, mostly from the vesalius paper's vignette. The output is below, for which I get a couple of warnings.
Attached is also the UMAP_to_RGB file generated from rgbUMAP(). I appreciate you help!

library(vesalius)
library(imager)
library(imagerExtra)
library(Seurat)
library(ggplot2)
library(patchwork)
library(RColorBrewer)
library(dplyr)
library(tvR)
library(sp)
library(grid)

# Standard libraries - just in case you are a command line fiend
library(utils)
library(stats)
library(graphics)
library(grDevices)
### Set seed
set.seed(1)

#slideTag <- c("maynard_noNA_151673", "maynard_noNA_151675")
slideTag <- c("maynard_noNA_151675")

#slideBeads <-c("/gcm-lfs1/akshay/projects/Spatial_Benchmark/Vesalius/slide-seqV2_data_vesalius/Puck_190926_03_bead_locations.csv")

#slideBeads <-c("/data/hoan/scRNAseq_data/maynard/maynard_noNA_151673_coords.csv")
slideBeads <-c("/data/hoan/scRNAseq_data/maynard/maynard_noNA_151675_coords.csv")

#slideCounts <- c("/data/hoan/scRNAseq_data/maynard/maynard_noNA_151673_UMImatrix.csv")
slideCounts <- c("/data/hoan/scRNAseq_data/maynard/maynard_noNA_151675_UMImatrix.csv")

## number of variable features
nfeatures <- 100

## Plot output directory
plots <- "/gcm-lfs1/akshay/projects/Spatial_Benchmark/Vesalius/plots/"

## Output directory
IP <- "/gcm-lfs1/akshay/projects/Spatial_Benchmark/Vesalius/IP/"

#------------------------------------------------------------------------------#
# Looping over files and converting UMAP projection values to RGB colours
#------------------------------------------------------------------------------#
countList <- list()
for(i in seq_along(slideTag)){
  message(paste("Loading Bead file",slideTag[i]))
  #----------------------------------------------------------------------------#
  # Loading coordinates
  #----------------------------------------------------------------------------#
  #bead <- ReadSlideSeq(slideBeads[i])
  #bead <- read.csv(slideBeads[i], header = TRUE, row.names = 1)
  locs <- read.table(slideBeads[i], header = FALSE, sep = ',') # row names integers, coords
  colnames(locs) <- c('x', 'y')
  locs <- cbind(barcode = c(1:nrow(locs)), locs)

  bead <- new(Class = 'SlideSeq', assay = "Spatial", coordinates = locs[,c("x","y")])
  rownames(bead@coordinates) <- locs$barcode

  message(paste("Loading Count file",slideTag[i]))
  #----------------------------------------------------------------------------#
  # Unconventional loading - however required as some data sets
  # Fail to load  - This ensures all data sets can be loaded
  # Count matrix
  #----------------------------------------------------------------------------#
  # count <- read.table(slideCounts[i], header = TRUE )
  cnt <- read.table(slideCounts[i], header = FALSE, sep = ',') #maynard UMIMatrix has no header or rows, row.names=NULL forces integer row names
  cnt <- t(as.matrix(cnt)) # transpose to genes x barcodes
  colnames(cnt) <- c(1:ncol(cnt)) # give barcode column names
  #cnt <- cbind(Gene = c(1:nrow(cnt)), cnt) # add column for gene row names
  rownames(cnt) <- c(1:nrow(cnt))
  
  #----------------------------------------------------------------------------#
  # Creating seurat spatial object
  # NOTE this code is taken from the Seurat source code as it does not seem that
  # Slide seq loading function are all exported
  # If this has been updated - this section can be changed accordingly
  #----------------------------------------------------------------------------#
  count <- CreateSeuratObject(cnt, assay ="Spatial")
  bead <- bead[Cells(x = count)]
  DefaultAssay(object = bead) <- "Spatial"
  count[["slice1"]] <- bead
  
  #----------------------------------------------------------------------------#
  # Seuart pre-processing steps
  #----------------------------------------------------------------------------#
  count <- NormalizeData(count)
  count <- FindVariableFeatures(count, nfeatures = nfeatures)
  count <- ScaleData(count)
  
  message(paste("Seurat Pre-Processed Data for",slideTag[i]))
  #----------------------------------------------------------------------------#
  # Embed PCA loading into the RGB colour space.
  # NOTE that the output of this function only assigns colours to each
  # coordinates we do not have an image yet.
  #----------------------------------------------------------------------------#
  count <- rgbUMAP(count, pcs =30, conserveSparse = FALSE)
  countList[[i]] <- count
  #----------------------------------------------------------------------------#
  # exporting embeddings
  #----------------------------------------------------------------------------#
  
  filenames <- paste0(IP,"UMAP_to_RGB_log_",nfeatures,"_",slideTag[i],".csv")
  exportRGB.csv(count,file = filenames)
}

imageBrain <- buildImageArray(countList[[1]],
                              filterGrid = 0.01,
                              filterThreshold = 0.9975,
                              resolution = 40,
                              cores = 10)

imageBrain <- iterativeSegmentation.array(imageBrain, colDepth =11,
                                          smoothIter = 20,
                                          method = c("iso","box"),
                                          sigma=1,box = 15,
                                          useCenter = T)

imageBrain <- isolateTerritories.array(imageBrain,
                                       captureRadius = 0.008,
                                       minBar = 50)

imgSmoothed <- imagePlot(imageBrain, as.cimg = F,cex = 15) + theme_void()
imgSmoothed

imgTerritory <- territoryPlot(imageBrain, randomise = TRUE,cex =15 , cex.pt=0.25) +
  theme_void()+
  theme(legend.text = element_text(size = 12),
        legend.title = element_text(size=12),
        plot.title = element_text(size =15),
        legend.position = "right")

pdf("vesalius_maynard_Territory.pdf", width = 8, height =6)
imgTerritory
dev.off()

OUTPUT:

source('vesalius.R')
Loading required package: Seurat
Attaching SeuratObject
Attaching sp
Loading required package: imager
Loading required package: magrittr
Registered S3 method overwritten by 'imager':
  method      from
  plot.imlist     

Attaching package: ‘imager’

The following object is masked from ‘package:magrittr’:

    add

The following object is masked from ‘package:sp’:

    bbox

The following objects are masked from ‘package:stats’:

    convolve, spectrum

The following object is masked from ‘package:graphics’:

    frame

The following object is masked from ‘package:base’:

    save.image

Loading required package: dplyr

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

Loading required package: tidyr

Attaching package: ‘tidyr’

The following object is masked from ‘package:imager’:

    fill

The following object is masked from ‘package:magrittr’:

    extract

Loading required package: ggplot2
------------------------------------------------------------
** tvR
**  - Total Variation Regularization
** 
** Version    : 0.3.2      (2022)
** Maintainer : Kisung You ([email protected])
**
** Please share any bugs or suggestions to the maintainer.
------------------------------------------------------------

Attaching package: ‘grid’

The following object is masked from ‘package:imager’:

    depth

Loading Bead file maynard_noNA_151675
Loading Count file maynard_noNA_151675
Performing log-normalization
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating gene variances
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating feature variances of standardized and clipped values
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Centering and scaling data matrix
  |======================================================================| 100%
Seurat Pre-Processed Data for maynard_noNA_151675
#------------------------------------------------------------------------# 
2022-09-22 10:15:40  Converting Sparse Matrix to Matrix 
Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
This message will be shown once per session
Error in data.frame(rownames(coordinates), coordinates[, c("x", "y")],  : 
  arguments imply differing number of rows: 3566, 3567
In addition: Warning message:
replacing previous import ‘imager::fill’ by ‘tidyr::fill’ when loading ‘vesalius’ 
> source('vesalius.R')
Loading Bead file maynard_noNA_151675
Loading Count file maynard_noNA_151675
Performing log-normalization
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating gene variances
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Calculating feature variances of standardized and clipped values
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Centering and scaling data matrix
  |======================================================================| 100%
Seurat Pre-Processed Data for maynard_noNA_151675
#------------------------------------------------------------------------# 
2022-09-22 10:28:14  Converting Sparse Matrix to Matrix 
#------------------------------------------------------------------------# 
#------------------------------------------------------------------------# 
2022-09-22 10:28:26  Checking Vesalius Input - Using slice  1 
2022-09-22 10:28:26  Class coersion to numeric 
2022-09-22 10:28:26  Filtering outlier beads 
2022-09-22 10:28:26  Generating Voronoi Tesselation 
2022-09-22 10:28:26  Filtering Tiles
2022-09-22 10:28:26  Rasterising Tiles 
2022-09-22 10:28:38  Resizing Image 
#------------------------------------------------------------------------# 
#------------------------------------------------------------------------# 
2022-09-22 10:28:38  Smoothing Iteration  20 
2022-09-22 10:28:38  Segmentation Iteration  1 
#------------------------------------------------------------------------# 
#------------------------------------------------------------------------# 
2022-09-22 10:28:39  Pooling territory  11 
#------------------------------------------------------------------------# 
Warning messages:
1: In as.cimg.data.frame(select(image, c("x", "y", "cc", "value"))) :
  Guessing image dimensions from maximum coordinate values
2: In as.cimg.data.frame(select(image, c("x", "y", "cc", "value"))) :
  Guessing image dimensions from maximum coordinate values
> client_loop: send disconnect: Broken pipe
>

UMAP_to_RGB_log_100_maynard_noNA_151675.csv

Error in: vesalius <- rgbPCA(vesalius,slices = 1)

Hi,

I am trying to run the model on the Human Invasive Ductal Carcinoma breast tissue public data.
Loading the data with Load_10x_Spatial:

vesalius <- Load10X_Spatial(
data.dir,
filename = "filtered_feature_bc_matrix.h5",
assay = "Spatial",)

And after running the prepocessing functions:

vesalius <- NormalizeData(vesalius)
vesalius <- FindVariableFeatures(vesalius, nfeatures = 2000)
vesalius <- ScaleData(vesalius)

The function rgbPCA is showing this error:

> vesalius <- rgbPCA(vesalius,slices = 1)
#-------------------------------------------------------------------------------------------------------------------#
2021-08-17 12:17:07 Running Principle Component Analysis in 1 slices
2021-08-17 12:17:08 Converting Loading Values to RGB in slice 1
2021-08-17 12:21:30 Normalising RGB values
Error in [.data.frame(coordinates, , c("x", "y")) :
undefined columns selected

Any suggestion about how to solve this?
Thanks in advance.

Recommend Projects

React

A declarative, efficient, and flexible JavaScript library for building user interfaces.
Vue.js

🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
Typescript

TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
TensorFlow

An Open Source Machine Learning Framework for Everyone
Django

The Web framework for perfectionists with deadlines.
Laravel

A PHP framework for web artisans
D3

Bring data to life with SVG, Canvas and HTML. 📊📈🎉

Recommend Topics

javascript

JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
web

Some thing interesting about web. New door for the world.
server

A server is a program made to process requests and deliver data to clients.
Machine learning

Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Visualization

Some thing interesting about visualization, use data art
Game

Some thing interesting about game, make everyone happy.

Recommend Org

Facebook

We are working to build community through open source technology. NB: members must have two-factor auth.
Microsoft

Open source projects and samples from Microsoft.
Google

Google ❤️ Open Source for everyone.
Alibaba

Alibaba Open Source for everyone
D3

Data-Driven Documents codes.
Tencent

China tencent open source team.