fit_mvn_smooth: Spatial multivariate normal mixture model clustering

In carter-allen/spruce: Spatial Random Effects Clustering of Single Cell Data

Spatial multivariate normal mixture model clustering

Description

Implement Gibbs sampling for MVN model with spatial smoothing

Usage

fit_mvn_smooth(
  Y,
  coords_df,
  K,
  r,
  nsim = 2000,
  burn = 1000,
  z_init = NULL,
  verbose = FALSE
)

Arguments

Y	An n x g matrix of gene expression values. n is the number of cell spots and g is the number of features.
coords_df	An n x 2 data frame or matrix of 2d spot coordinates.
K	The number of mixture components to fit.
r	Empirical spatial smoothing
nsim	Number of total MCMC iterations to run.
burn	Number of MCMC iterations to discard as burn in. The number of saved samples is nsim - burn.
z_init	Optional initialized allocation vector. Randomly initialized if NULL.
verbose	Logical for printing cluster allocations at each iteration.

Value

a list of posterior samples

Examples

## Not run: 
# parameters
data(coords_df_sim)
coords_df <- coords_df_sim[,1:2]
z <- remap_canonical2(coords_df_sim$z)
W_nn <- scran::buildKNNGraph(as.matrix(coords_df),
                             k = 4,
                             transposed = TRUE)
A <- igraph::as_adjacency_matrix(W_nn,
                                 type = "both",
                                 sparse = FALSE)
                                 
n <- nrow(coords_df) # number of observations
g <- 3 # number of features
K <- length(unique(coords_df_sim$z)) # number of clusters (mixture components)
pi <- table(z)/length(z) # cluster membership probability

# Cluster Specific Parameters
# cluster specific means
Mu <- list(
  Mu1 = rnorm(g,-2,1),
  Mu2 = rnorm(g,-1,1),
  Mu3 = rnorm(g,1,1),
  Mu4 = rnorm(g,2,1)
)
# cluster specific variance-covariance
S <- matrix(0.5,nrow = g,ncol = g) # y covariance matrix
diag(S) <- 1
Sig <- list(
  Sig1 = S,
  Sig2 = S, 
  Sig3 = S,
  Sig4 = S
)

Y <- matrix(0, nrow = n, ncol = g)
for(i in 1:n)
{
  Y[i,] <- mvtnorm::rmvnorm(1,mean = Mu[[z[i]]],sigma = Sig[[z[i]]])
}

# sometimes helps to initialize using heuristic like kmeans
fitk <- stats::kmeans(Y,4)
z_km <- remap_canonical2(fitk$cluster)

# fit model
# use more iterations in practice
fit1 <- fit_mvn_smooth(Y,coords_df,4,2,10,0,z_km)
## End(Not run)

carter-allen/spruce documentation built on April 6, 2024, 8:34 p.m.