Site icon R-bloggers

Plotting Watts-Strogatz model

[This article was first published on R snippets, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
Recently I wanted to reproduce Figure 2 from Watts and Strogatz (1998). The task using igraph is simple but an interesting task was annotation of the resulting plot.

Watts-Strogatz model generates graphs that have so called small-world network property. Such networks should have low average path length and high clustering coefficient. The algorithm has three parameters: number of nodes in the graph, initial number of neighbors of each node distributed on a ring and rewiring probability.

Interestingly in Watts-Strogatz model having small but positive values of rewiring probability generates graphs having desired properties – and this is exactly depicted on Figure 2 in their article.

I decided to replicate it. To enhance it I wanted to plot median and 5 and 95 percentile of distribution of average path length and clustering coefficient as a function of rewiring probability.

Here you have the code that generates the graph (warning: it takes about 1 minute to run):

library(igraph)< o:p>
set.seed(1)< o:p>

avg.stat <- function(nei, p) {< o:p>
  result <- replicate(1000, {< o:p>
    wsg <- watts.strogatz.game(1, 100, nei, p)< o:p>
    c(average.path.length(wsg),< o:p>
      transitivity(wsg))< o:p>
  })< o:p>
  apply(result, 1, quantile, probs = c(0.5, 0.05, 0.95))< o:p>
}< o:p>

nei <- 6< o:p>
p <- 2 ^ seq(0, 10, len = 21)< o:p>
result <- sapply(p, avg.stat, nei = nei)< o:p>
result <- t(result / rep(avg.stat(nei, 0)[1,], each = 3))< o:p>
par(mar=c(3.220.20.2), mgp=c(210))< o:p>
matplot(p, result, type = “l”, log = “x”, xaxt = “n”, ylab = “”,< o:p>
        lty = rep(c(1,2,2),2), col=rep(c(1,2), each=3))< o:p>
axis(1, at = 2 ^ -(0:10),< o:p>
     labels =  c(1, parse(text = paste(2, 1:10, sep = “^-“,< o:p>
                                       collapse “;”))))< o:p>
legend(“bottomleft”, c(“average path length”, “clustering coefficient”),< o:p>
       lty = 1, col = c(1, 2))

The result of the procedure is the following picture:


It looks very similar to what is shown in the article (apart from adding lines depicting 5 and 95 percentile of distributions of both graph characteristics).

However, the interesting part was to properly annotate X-axis on the plot. Of course you can use expression function to get it but then the problem is that you have to do it ten times. Interestingly parsing a string containing those ten expressions separated by semicolons works just as needed.

To leave a comment for the author, please follow the link and comment on their blog: R snippets.

R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.