{"word_count": 403, "character_count": 3044}
For some dissemination work, I want to create a nice graph to explain the exponential growth in pandemics, related to the value of [latex]R_0[/latex]. Recall that [latex]R_0[/latex] corresponds to the average number of people that a contagious person can infect. Hence, with [latex]R_0=1.5[/latex], 4 people will contaminate 6 people, and those 6 will contaminate 9, etc. After [latex]n[/latex] iteration, the number of contaminated people is simply [latex]R_0{}^n[/latex]. As explained by Daniel Kahnemanpeople, certainly including myself, don’t seem to be able to think straight about exponential growth. What we see today are infections that occurred 2 or 3 weeks ago and the deaths today are people who got infected 4 or 5 weeks ago. All of this is I think beyond intuitive human comprehensionFor different values of [latex]R_0[/latex] (on each row), I wanted to visualise the number of contaminated people after 3, 5 or 7 iterations, since graphs are usually the most simple way to give some intuition. The graph I had in mind was the following(to be honest, I am quite sure I had seen it somewhere, but I cannot find where). The main challenge here is pack optimally [latex]k[/latex] identical circles intro a unit circle: we need here the location of the points (center of the disks) and the radius. It seems to be a rather complicated mathematical problem. Nicely, on http://hydra.nat.uni-magdeburg.de/packing, it is possible to get the "best known packings of equal circles in a circle" (up to 5000, but many [latex]k[/latex]'s are missing). For instance, for [latex]k=37[/latex], we haveAnd interestingly, on the same website, we can get the coordinates of the centers, for example with 37 disks, so it is possible to recreate the R graph.k = 37base = read.table(paste("http://hydra.nat.uni-magdeburg.de/packing/cci/txt/cci",k,".txt",sep=""), header=FALSE)The problem, as discussed earlier, is that some cases are not solved, yes, for instance [latex]k=2^{12}=4096[/latex]: the next feasable case is 4105. To avoid that issue, one can useT = "Error"while(T == "Error"){    T = substr(try(base = read.table(paste("http://hydra.nat.uni-magdeburg.de/packing/cci/txt/cci",k,".txt",sep=""), header=FALSE),silent = TRUE),1,5)k=k+1} k=k-1Now we can almost plot it. The problem is that the radius of the circles is missing, here. But we can compute itD=as.matrix(dist(x = base[,2:3]))diag(D)=1e5i=which(D == min(D), arr.ind = TRUE)r = D[i[1,1],i[1,2]]Here the radius isr[1] 0.2959118To plot it, useplot(base$V2,base$V3,xlim=c(-1,1),ylim=c(-1,1))n=100theta=seq(0,pi,length=n+1)circ= function(x,y,r,h=1){  vu=x+r*cos(theta)  vv=r*sin(theta)  cbind(c(vu,rev(vu))*h,c(y+vv,y-rev(vv))*h)}for(i in 1:k) polygon(circ(base[i,2],base[i,3],r/2*.95),col=colr,border=NA)We can now use that code to create the graph above, with [latex]k=R_0{}^n[/latex] for various values of [latex]n[/latex]And we can also use it to visualize more subtile differences, like [latex]R_0=1.1[/latex], [latex]R_0=1.3[/latex], [latex]R_0=1.5[/latex] and [latex]R_0=1.7[/latex]