Sur ce site

dynSBM : dynamic stochastic block model in R


dynSBM is a R/C++ implementation of a model that combines a stochastic block model (SBM) for its static part with independent Markov chains for the evolution of the nodes groups through time.
It deals with binary or weighted dynamic networks (with discrete or continuous edges).

Package availability

The package was initially developed in R. It has been refactored with most of the code in C++ to enable the analysis of larger networks (larger number of nodes, time steps and number of groups).

  • last version (recommended) : available on CRAN
  • version 0.3-beta : improvement of plot functions. Tested on Windows (zip archive below).
    GZ - 38.8 ko
    Zip - 629.8 ko
  • version 0.2 : implementation of more robust estimation procedure. Full support for directed networks and/or with self loops. New functions for plotting the results.
    GZ - 37.7 ko
  • version 0.1 : code refactoring in C++ (now requires the package Rcpp), speedup x100
    GZ - 22.8 ko
  • version 0.0 : fully written in R, used for the publication in JRSSB
    GZ - 15.2 ko

Full demo

A complete demo is available online here with the required data and source files available here

Fully developed example

R CMD INSTALL dynsbm_x.x.tar.gz

Loading a simulated dataset (T=5 time points, N=40 nodes and Q=4 groups used to simulate the data) :


Reading the help of the estimation procedure (in particular, understanding that the model can deal with different edge types - binary or weighted) :


Estimating the model with Q=1 to Q=6 groups (can take about 1 minute) :

list.dynsbm <- select.dynsbm(simdataT5Q4N40binary, Qmin=1, Qmax=6, edge.type="binary", nstart=10)

Selecting Q=4 groups (because it maximizes the ICL criteria - if not, please increase nstart above) :

dynsbm <- list.dynsbm[[4]]

Examining the group membership over time :


Plotting intra/inter connectivity patterns :

connectivity.plot(dynsbm, simdataT5Q4N40binary)

Plotting switches between groups :



If you use dynSBM in a published work, please cite one of the following reference :

Catherine Matias and Vincent Miele, Statistical clustering of temporal networks through a dynamic stochastic block model, Journal of the Royal Statistical Society : Series B (2017)

PDF - 802.1 ko

Vincent Miele and Catherine Matias, Revealing the hidden structure of dynamic ecological networks, Royal Society Open Science (2017)


For any bugs, information or feedback, please contact :

Vincent Miele - vincent.miele followed by


Can I use dynSBM if I have only one aggregated network (i.e. only one time step) ?
Yes. It will consist in estimating a classical SBM on your network. The only trick is that dynSBM expects a R array TxNxN and in this case T=1. A short example :

# random data
Y1 <- as.matrix(get.adjacency(make_full_graph(10) %du% make_full_graph(10)))

# converting to R array
Y <- array(dim=c(1,20,20))
Y[1,,] <- Y1

# dynSBM
list.sbm <- select.dynsbm(Y, Qmin=1, Qmax=3, edge.type="binary")
sbm <- list.sbm[[2]]
connectivity.plot(sbm, Y)

Can I use dynSBM on a list of igraph objects ?
No but yes. dynSBM expects an array TxNxN (the adjacency matrices), therefore it is necessary to convert the list into a R array. A short example :

# random data
g <- make_full_graph(20) %du% make_full_graph(20)
V(g)$name <- paste("N",1:40,sep="")
g1 <- induced_subgraph(g, sample(1:40,30))
g2 <- induced_subgraph(g, sample(1:40,30))
g3 <- induced_subgraph(g, sample(1:40,30))
list.g <- list(g1,g2,g3)

# converting to R array <- unique(sort(unlist(lapply(list.g, function(g) V(g)$name))))
T <- length(list.g)
N <- length(
Y <- array(dim=c(T,N,N))
for (t in 1:T){
   g.tmp <- graph.empty(n=N, directed=F)
   V(g.tmp)$name <-
   Y[t,,] <- as.matrix(get.adjacency(union(g.tmp,list.g[[t]])))

# dynSBM
list.dynsbm <- select.dynsbm(Y, Qmin=1, Qmax=3, edge.type="binary")
dynsbm <- list.dynsbm[[2]]
connectivity.plot(dynsbm, Y)

Ho can I discretize the weighted adjacency matrices ?
Assuming K categories. From Y, it is possible to use the quantiles like this :

T <- dim(Y)[1]
N <- dim(Y)[2]
Yd <- array(0L, c(T,N,N))
for (t in 1:T){
       Yt <- Y[t,,]
       Yt.res <- matrix(0,N,N)
       bounds <- quantile(Yt[Yt>0], probs=seq(0,1,by=1/K))
       Yt.res[which(Yt>0 & Yt<=bounds[2])] <- 1
       for(k in 2:K){
           Yt.res[which(Yt>bounds[k] & Yt<=bounds[k+1])] <- k
       Yd[t,,] <- Yt.res