**Wiekvoet**, and kindly contributed to R-bloggers)

I am listening to the audiobook Sync: How Order Emerges from Chaos in the Universe, Nature, and Daily Lifeby Steven Strogatz which I got from Audible. Obviously a mathematical book is not ideal to listen to, but lacking illustrations I can make them myself. On top of that, you learn more from programming than listening.

He mentioned a simple system which synchronizes itself. As far as I understand a large number of items are considered. Each item has more or less the same properties. Small steps are taken. Every step has a small additive effect, a multiplicative decrease. When an item gets over a limit, then it resets to zero and causes a small increase in the other items.

What is supposed to happen is that at low values you get steep slopes and at high values small slopes. As time progresses the small increases due to resets cause a synchronization.

To program this was not difficult, except for the initialization. If you start the items at random levels you have too many at low levels. In the end I decided to add an item at each time point for an init period. What did surprise me, was that calculating what happens actually takes less computer time than plotting it. But maybe that is because I wanted to use ggplot and a low alpha to show the bunching of items. I also discovered it is also not wise to display too complex a figure in ggplot running in Eclipse. It can freeze the whole program. I took to saving all plots rather than displaying within Eclipse.

In the init phase every item has different start and ending times. Note the black line at level 0 which consists of not yet running items

After 10000 steps it is still chaos, with a few somewhat fatter lines.

Not that they are perfect after 300000 steps, but it gets more and more closer.

#### R code

library(ggplot2)

plotpart <- function(xmat,fname=’een.png’,xstart=1) {

fin <- ncol(xmat)

nitems <- nrow(xmat)

df <- data.frame(score=as.numeric(xmat),

time=rep((1:fin)+xstart,each=nitems),

item =rep(1:nitems,fin))

g<- ggplot(df,aes(x=time,y=score,group=item,alpha=score)) +

geom_line() +

scale_alpha_continuous(range=c(.05,.0501) )+

theme(legend.position=’none’)

png(fname)

print(g)

dev.off()

}

nextgen <- function(x,add,shrink=.99,limit=1,spil=.0105) {

x <- x*shrink+add

maxed <- x>limit

x[maxed]<- .01#add[maxed]

x[!maxed] <- x[!maxed] + sum(maxed)*spil

x[x>limit] <- limit

x

}

ngenmat <- function(xstart,niter,add,spil) {

nitems <- length(xstart)

xmat <- matrix(0,nrow=nitems,ncol=niter)

xmat[,1 ] <- xstart

for ( i in 2:niter) {

xmat[,i] <- nextgen(xmat[,i-1],add=add,spil=spil)

}

xmat

}

ngenmove <- function(xstart,niter,add,spil) {

for ( i in 1:niter) {

xstart <- nextgen(xstart,add=add,spil=spil)

}

xstart

}

# init a run

nitems <- 1000

niter <- 1000

xmat <- matrix(0,nrow=nitems,ncol=niter)

for (i in 1:nitems) xmat[i,i] <- runif(1,min=.01,max=.02)

add <- rnorm(nitems,.0103,.00005)

spil=.000025

for ( i in 2:niter) {

todo <- xmat[,i-1]>0

xmat[todo,i] <- nextgen(xmat[todo,i-1],add=add[todo],spil=spil)

}

#and show it

show <- sample(1:nrow(xmat),50)

plotpart(xmat[show,],fname=’een.png’,xstart=1)

# move 9000 steps

x1 <- ngenmove(xmat[,ncol(xmat)],niter=9000,add=add,spil=spil)

xm1 <- ngenmat(x1,niter=1000,add=add,spil=spil)

plotpart(xm1[show,],fname=’twee.png’,xstart=10000)

x2 <- ngenmove(xm1[,ncol(xm1)],niter=39000,add=add,spil=spil)

xm2 <- ngenmat(x2,niter=1000,add=add,spil=spil)

plotpart(xm2[show,],fname=’drie.png’,xstart=50000)

x3 <- ngenmove(xm2[,ncol(xm2)],niter=49000,add=add,spil=spil)

xm3 <- ngenmat(x3,niter=1000,add=add,spil=spil)

plotpart(xm3[show,],99000,fname=’vier.png’)

x4 <- ngenmove(xm3[,ncol(xm3)],niter=49000,add=add,spil=spil)

xm4 <- ngenmat(x4,niter=1000,add=add,spil=spil)

plotpart(xm4[show,],150000,fname=’vijf.png’)

x5 <- ngenmove(xm4[,ncol(xm4)],niter=49000,add=add,spil=spil)

xm5 <- ngenmat(x5,niter=1000,add=add,spil=spil)

plotpart(xm5[show,],200000,fname=’zes.png’)

x6 <- ngenmove(xm5[,ncol(xm5)],niter=99000,add=add,spil=spil)

xm6 <- ngenmat(x6,niter=1000,add=add,spil=spil)

plotpart(xm6[show,],300000,fname=’zeven.png’)

**leave a comment**for the author, please follow the link and comment on his blog:

**Wiekvoet**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...