GSVA with a genelist on a "per-cluster" basis
1
0
Entering edit mode
9 weeks ago
bio_info ▴ 20

Hi everyone,

I have a scRNA-seq immune cell dataset (Seurat object) with 8 clusters and I have been asked to perform a GSVA or ssGSEA of a genelist on a "per-cluster" basis. I do not have much experience with enrichment analysis and I tried going through some tutorials but none of them explain how to do this for each cluster.

Could someone please explain how to run a GSVA or ssGSEA for this genelist on all clusters with some steps? Do I need to compute differentially expressed genes for each cluster? What packages can I use for the analysis?

GSVA scRNA-seq genelist • 407 views
ADD COMMENT
1
Entering edit mode
9 weeks ago

GSVA / ssGSEA were originally designed for bulk expression analysis and so some of the assumptions/caveats with single-cell data don't carry over well. For a very basic approach you could

log.mat <- GetAssayData(seurat.object, assay = "RNA", layer = "data")

And use that to run GSVA/ssGSEA using the GSVA package. IIRC there are some single-cell wrappers for this, but really doesn't extend functionallity very far.

After that you could calculate mean/median scores per cluster using something like

cluster.means <- t(apply(gsva.res, 1, \(x) tapply(x, seurat.object$your.clusters, mean)))

Originally, the advice was for people to use limma on top of GSVA/ssGSEA to do differential expression on the output, but this has caveats for use with single-cell data and I wouldn't advise it. Either way, for an exploratory analysis what you probably want to know more is what expression programmes are active in the clusters rather than have a specific p-value badge to pin on it.

ADD COMMENT
0
Entering edit mode

Hi! Thanks a lot for the explanation! How can I integrate my custom genelist in this analysis here? As an input I have my clustered scRNA-seq dataset and a genelist, but I am confused how to carry out a GSVA with these inputs.

ADD REPLY
1
Entering edit mode
library(GSVA)

log.mat <- GetAssayData(seurat.object, assay = "RNA", layer = "data")
genesets <- list(Genelist1 = some.list.of.genes1,
                 Genelist2 = some.list.of.genes2,
                 ... = ...)
gsva.par <- gsvaParam(log.mat,
                      genesets,
                      minSize = 10, 
                      maxSize = 500, ## tailor this to your genelist size

                      maxDiff = TRUE)
gsva.res <- gsva(gsva.par)

Something like this - I haven't tested it for obvious reasons. See the vignette for GSVA It's pretty self-explanatory

ADD REPLY
0
Entering edit mode

This is really helpful! Thanks a lot, I'll checkout the vignette and try to implement it on my data :)

ADD REPLY

Login before adding your answer.

Traffic: 2656 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6