• Overview
  • FASTQ+
  • PISA
  • Yano
  • My GitHub Page
  • Discussions

On this page

  • Introduction
  • INSTALL
  • Get started
  • Short cases
  • Changelog
  • Issues or questions

Yano

Spatial dissimilarity analysis in single cells

Introduction

Yano represents an R/C toolkit designed for conducting spatial dissimilarity analysis on single-cell RNA sequencing data. This method revolves around the core concept of examining the distinct expression patterns of a given feature (e.g. exon, snp allele) in relation to its associated binding feature (typically a gene or another allele at the same genomic locus) within the context of cell lineage (1D), spatial position (2D), or the multi-dimensional PCA space. The discernible differences in feature expression patterns and their binding features provide insights into a range of biological phenomena, including alternative splicing, cis-antisense RNA regulation, allele-specific gene expression, and more.

Yano is seamlessly integrated with Seurat, building upon the Seurat object’s framework. Users can perform conventional cell clustering analyses using the state-of-the-art Seurat pipeline and then incorporate exon, SNP counts as new “assays” within the Seurat objects. Subsequently, Yano facilitates the assessment of spatial dissimilarity between these two assays. For more details about the method, please refer to our manuscript.

INSTALL

if (!require("BiocManager")) install.packages('BiocManager') 
BiocManager::install("shiquan/Yano")

Notice: Multithread mode is disabled on macOS by default due to the lack of OpenMP support. However, data.table provides a useful tutorial on how to enable OpenMP on macOS. For more information, please refer to this guide: https://github.com/Rdatatable/data.table/wiki/Installation#enable-openmp-for-macos.

Get started

A typical workflow with Yano starts by using a built-in dataset.

require(Yano)
data("glbt_small")
DefaultAssay(glbt_small) <- "RNA"
glbt_small <- NormalizeData(glbt_small) %>% RunUMAP(dim = 1:20)

DimPlot(glbt_small, label = TRUE, label.size = 5)

DefaultAssay(glbt_small) <- "exon"
glbt_small <- NormalizeData(glbt_small)
Meta(glbt_small) %>% head
data frame with 0 columns and 6 rows
glbt_small <- ParseExonName(glbt_small)
Meta(glbt_small) %>% head
                                    chr     start       end gene_name strand
chr1:154169305-154169383/-/TPM3    chr1 154169305 154169383      TPM3      -
chr11:35197162-35197793/+/CD44    chr11  35197162  35197793      CD44      +
chr11:75421727-75422280/+/RPS3    chr11  75421727  75422280      RPS3      +
chr11:123060825-123061329/-/HSPA8 chr11 123060825 123061329     HSPA8      -
chr11:123061588-123061833/-/HSPA8 chr11 123061588 123061833     HSPA8      -
chr11:123061869-123062022/-/HSPA8 chr11 123061869 123062022     HSPA8      -
grep("_wm$",names(glbt_small), value=TRUE)
character(0)
glbt_small <- RunAutoCorr(glbt_small)
grep("_wm$",names(glbt_small), value=TRUE)
[1] "pca_wm"
# Perform spatial dissimilarity test
glbt_small <- RunSDT(glbt_small, bind.name = "gene_name", bind.assay = "RNA")
Meta(glbt_small) %>% head
                                    chr     start       end gene_name strand
chr1:154169305-154169383/-/TPM3    chr1 154169305 154169383      TPM3      -
chr11:35197162-35197793/+/CD44    chr11  35197162  35197793      CD44      +
chr11:75421727-75422280/+/RPS3    chr11  75421727  75422280      RPS3      +
chr11:123060825-123061329/-/HSPA8 chr11 123060825 123061329     HSPA8      -
chr11:123061588-123061833/-/HSPA8 chr11 123061588 123061833     HSPA8      -
chr11:123061869-123062022/-/HSPA8 chr11 123061869 123062022     HSPA8      -
                                   moransi.pval    moransi autocorr.variable
chr1:154169305-154169383/-/TPM3   2.144774e-108 0.09755130              TRUE
chr11:35197162-35197793/+/CD44     0.000000e+00 0.25567539              TRUE
chr11:75421727-75422280/+/RPS3     4.383662e-38 0.05692911              TRUE
chr11:123060825-123061329/-/HSPA8  1.099236e-63 0.07462144              TRUE
chr11:123061588-123061833/-/HSPA8 7.645339e-158 0.11891356              TRUE
chr11:123061869-123062022/-/HSPA8 1.466907e-150 0.11541286              TRUE
                                  gene_name.D gene_name.t gene_name.pval
chr1:154169305-154169383/-/TPM3     0.3436496   8.4268520   1.448736e-13
chr11:35197162-35197793/+/CD44      0.1966297   0.8886768   1.881655e-01
chr11:75421727-75422280/+/RPS3      0.3486643   8.1503513   5.710819e-13
chr11:123060825-123061329/-/HSPA8   0.2446247   3.7857530   1.314186e-04
chr11:123061588-123061833/-/HSPA8   0.3633522   8.8646031   1.631423e-14
chr11:123061869-123062022/-/HSPA8   0.3320379   8.4846775   1.086526e-13
                                  gene_name.padj
chr1:154169305-154169383/-/TPM3     5.695820e-12
chr11:35197162-35197793/+/CD44      1.000000e+00
chr11:75421727-75422280/+/RPS3      2.132991e-11
chr11:123060825-123061329/-/HSPA8   3.067802e-03
chr11:123061588-123061833/-/HSPA8   8.124484e-13
chr11:123061869-123062022/-/HSPA8   4.774324e-12
# Manhattan plot for spatial dissimilarity test result
FbtPlot(glbt_small, val = "gene_name.padj")

FeaturePlot(glbt_small, features = c("chr19:16095264-16095454/+/TPM4", "TPM4"), order=TRUE)

# Track plot for gene coverage at different cell types
db <- gtf2db("./gencode.v44.annotation.gtf.gz")
[2025-06-19 19:37:59] GTF loading..
[2025-06-19 19:38:20] Load 62700 genes.
TrackPlot(bamfile="./Parent_SC3v3_Human_Glioblastoma_possorted_genome_bam.bam", gtf =db, gene = "TPM4", junc = TRUE, cell.group = Idents(glbt_small), highlights = c(16095264,16095454))

See short cases for more details.

Short cases

  • Alternative splicing analysis for scRNA-seq
  • Allele-specific gene expression analysis for scRNA-seq
  • Annotating and prioritizing genetic variants for scRNA-seq
  • Perform alternative splicing analysis for multiple Visium samples
  • Select cells from reduction maps
  • Perform alternative splicing analysis for cell trajectory and user-defined embeddings

Changelog

1.2 2025/06/18
  • Rename RunBlockCorr to RunSDT.
  • Reimplement core algorithm with sparse matrix structure. Highly improved speed (>3x).
1.1 2025/06/02

Simplified parameters for RunCorrBlock() to streamlines the function and makes it more concise. In addition, new selector functions, *Selector(), introduced.

1.0 2025/02/19

First stable release.

0.0.0.9999 2023/03/22

Init version.

Issues or questions

  • https://github.com/shiquan/Yano/issues
  • https://github.com/shiquan/Yano/discussions
Back to top