Visualizing parameters for clustering PBMCs using gganimate

gganimate_test.R

##################### Visualizing parameters for clustering PBMCs using gganimate
library(ggplot2)
library(gganimate)

## PBMC dataset
library(MUDAN)
data("pbmcA")

## Downsample to run faster
cd <- MUDAN::cleanCounts(pbmcA)[, 1:2000]
## CPM normalization
mat <- MUDAN::normalizeCounts(cd, 
                              verbose=FALSE) 
## variance normalize, identify overdispersed genes
matnorm.info <- MUDAN::normalizeVariance(mat, 
                                         details=TRUE, 
                                         verbose=FALSE) 
## log transform
matnorm <- log10(matnorm.info$mat+1) 
## 30 PCs on overdispersed genes
pcs <- MUDAN::getPcs(matnorm[matnorm.info$ods,], 
                     nGenes=length(matnorm.info$ods), 
                     nPcs=30, 
                     verbose=FALSE) 
## Community detection
com <- getComMembership(pcs, 
                        k=10, method=igraph::cluster_louvain, 
                        verbose=FALSE) 

##################### TSNE
library(Rtsne)

## Loop through different perplexities
perplexities <- c(3, 5, 10, 15, 30, 50, 100, 200)
embs <- do.call(rbind, lapply(perplexities, function(p) {
  print(p)
  
  ## TSNE embedding with PCs
  emb <- Rtsne::Rtsne(pcs, perplexity=p)$Y
  rownames(emb) <- rownames(pcs)
  
  ## keep track of p
  emb <- cbind(emb, p, com)
  
  return(emb)
}))

## Create animation
df <- data.frame(embs)
colnames(df) <- c('X', 'Y', 'p', 'com')
df$com <- as.factor(df$com)
head(df)

## Plot one
vi <- embs[, 3] == 30 ## perplexity = 30
ggplot(df[vi,], aes(X, Y, color = com)) + 
  geom_point(alpha=0.5) +
  theme_void() +
  theme(legend.position="none")

## Animation
ggplot(df, aes(X, Y, color = com)) + 
  geom_point(alpha=0.5) +
  theme_void() + 
  theme(legend.position="none") + 
  transition_states(p) + 
  labs(title = 'Perplexity: {closest_state}')


##################### UMAP
library(umap)

## Loop through different min_dists
mindists <- c(0.01, 0.1, 0.3, 0.5, 0.8, 1, 1.5, 2)
embs.umap <- do.call(rbind, lapply(mindists, function(md) {
  print(md)
  
  ## UMAP embedding with PCs
  config <- umap.defaults
  config$min_dist = md
  emb <- umap::umap(pcs, config=config)$layout
  
  ## keep track of p
  emb <- cbind(emb, md, com)
  
  return(emb)
}))

## Create animation
df <- data.frame(embs.umap)
colnames(df) <- c('X', 'Y', 'md', 'com')
df$com <- as.factor(df$com)
head(df)

## Plot one
vi <- embs.umap[, 3] == 1 ## min_dist = 1
ggplot(df[vi,], aes(X, Y, color = com)) + 
  geom_point(alpha=0.5) +
  theme_void() +
  theme(legend.position="none")

## Animation
ggplot(df, aes(X, Y, color = com)) + 
  geom_point(alpha=0.5) +
  theme_void() + 
  theme(legend.position="none") + 
  transition_states(md) + 
  labs(title = 'min_dist: {closest_state}')