for_liz.R

In smorabit/scicat:

library(devtools)
devtools::install_github(repo='smorabit/scicat', force=T)
library(scicat)


# check the stability of a clustering of a particular resolution
#
# specify the seurat object, the number of bootstrap iterations, the number of cells to sample,
# and the cluster resolution.
#
# stability of cluster resolution 0.9:
n_bootstrap = 10 # how many repetitions
n_cells = round(ncol(seurat_object)*0.9) # how many cells to sample
bl <- sc_cluster_stability(seurat_object, n_cells=n_cells, n_bootstrap=n_bootstrap, cluster_res=0.9, verbose=T)

# plot this as a barplot
pdf("figures/cluster_stability_0.9res.pdf", width=8, height=3)
colors <- ifelse(bl >= 0.60, "blue", "red")
colors <- ifelse(bl >= 0.85, "green", colors)
stability <- ifelse(bl >= 0.60, "stable", "unstable")
stability <- ifelse(bl >= 0.85, "highly stable", stability)
# having problems using the scicat function directly
p <- cluster_stability_barplot(bl, stability=stability, clusters=cell_type_identities$Cluster.ID) + ggtitle("cluster res=0.9, n_bootstrap=10, n_cells=90%")
print(p)
dev.off()


# heatmap
# similar to above, but this time we are testing a range of cluster resolutions. This takes a long time to run!!
cluster_resolutions = seq(1:10)/10
n_bootstrap = 25
n_cells = round(ncol(NucSeq)*85/100)
resm <- sc_cluster_res_stability(NucSeq, cluster_resolutions, n_cells=n_cells, n_bootstrap=n_bootstrap, reduction="inmf", verbose=T)
pdf("figures/cluster_stability_heatmap.pdf", width=20, height=8)
print(cluster_res_stability_heatmap(resm))
dev.off()