A web-based graphics device for animated plots in R

class: left, middle, inverse, title-slide

# A web-based graphics device for animated plots in R
### Chun Fung (Jackson) Kwok
### St. Vincent’s Institute of Medical Research
### 2021/06/17 @ BCG seminar, SVI

---

## Outline

### 1. The R package 'animate'

.medium[
- Introduction
- Basic usage and examples
]

### 2. Some applications in dynamical systems and agent-based modelling

.medium[
- Lorenz system
- Predator-prey process
- Schelling's model of segregation
]

---

## The R package 'animate' - Introduction

### What is it?
.medium[
A web-based graphics device for animated plots in R.
]

### Why a new device for animated plots?
.medium[
- Existing work: `animation`, `gganimate`
  - String many images into a .gif or .mp4 file
  - `base` plot takes ~5ms per frame, `ggplot` takes ~150ms per frame
  - Results are available after all frames are generated ⇒ can be slow to work with
]

.medium[
- Live plotting requires 60 frames per second (= 16.67ms per frame)
  - Works okay for simple `base` plots, but in general, there is the "blinking" problem.
]

.medium[
- Advance the R language by __providing support to animated plotting__ and __improving its integration with Web technologies__.
]

---

## The R package 'animate' - Introduction

### How does it work?

.medium[
#### Communication between R and the browser

- Use `httpuv` to create a WebSocket server in R
- Launch the __graphics device__ at the browser
- Exchange JSON data between the server and the browser. Conversion is done by `jsonlite`.
]

.medium[
#### The graphics device

- Developed in R using `sketch`, i.e. it is transpiled to JavaScript and runs on the browser.
- The core engine is based on the d3 JavaScript library.
- Unconventional top-down approach for graphics device development
- The API mimics R base primitives.
]

---

## The R package 'animate' - Basic usage

### Initialising and using the device

```r
library(animate)
device <- animate$new(width = 600, height = 400)  # Initialisation takes ~0.5s
device$plot(1:10, 1:10)  # Main function
device$clear()           # Clear a plot 
device$off()             # Switch off the device when you are done
```

.medium[
Alternatively, use `attach` to access the main functions directly.
]

```r
# Common set-up
library(animate)
device <- animate$new(width = 600, height = 400)  # Initialisation takes ~0.5s
*attach(device)
```

```r
plot(1:10, 1:10)  # Main function
clear()           # Clear a plot 
off()             # Switch off the device when you are done
*detach(device)
```

---

## The R package 'animate' - Basic usage

### Creating animated plots

.medium[
The __most important__ idea of this package is 
- __every rendered object on the screen must have an ID__;
- these IDs are used to decide which objects need to be modified (to create the animation effect).
]

#### The `plot` function

```r
x <- 1:10
y <- 1:10
id <- new_id(x)   # Give each point an ID
*plot(x, y, id = id)

new_y <- 10:1
*plot(x, new_y, id = id, transition = TRUE)
```
---

## The R package 'animate' - Basic usage

#### The `points` function

```r
y2 <- 10 * runif(10)
id <- new_id(y2, prefix = "points")   # Give each point an ID
*points(x, y2, id = id, bg = "red")

new_y2 <- 10 * runif(10)
*points(x, new_y2, id = id, bg = "lightgreen", transition = list(duration = 2000))
```

Note that when a drawing function is called, the scale of the data is computed based on the data provided. Use the `xlim` and `ylim` arguments, if one wants to keep the scale (and axes) constant.

#### The `lines` function

```r
y3 <- 10 * runif(10)
id <- "line-1"   # a line needs 1 ID only (despite containing multiple points)
*lines(x, y3, id = id)

new_y3 <- 10 * runif(10)
*lines(x, new_y3, id = id, transition = list(delay = 1000))
```

---

## Applications - Lorenz system

.medium[
A system of ODEs with chaotic behaviour
]

`$$\begin{aligned}
\dfrac{dx}{dt} = \sigma (y - x), \quad
\dfrac{dy}{dt} = x (\rho - z) - y, \quad
\dfrac{dz}{dt} = xy - \beta z
\end{aligned}$$`

.medium[
First-order solution can be numerically solved using Euler's method.
]

.height-20pc[

```r
sigma <- 10
beta <- 8 / 3
rho <- 28
x <- y <- z <- 1
xs <- x
ys <- y
zs <- z
dt <- 0.015

for (i in 1:2000) {
  # Euler's method
  dx <- sigma * (y - x) * dt
  dy <- (x * (rho - z) - y) * dt
  dz <- (x * y - beta * z) * dt
  x <- x + dx
  y <- y + dy
  z <- z + dz
  # Keep a record for line trace plot
  xs <- c(xs, x)
  ys <- c(ys, y)
  zs <- c(zs, z)
  # Plot the x-y solution plane
* plot(x, y, id = "ID-1", xlim = c(-30, 30), ylim = c(-30, 40))
* lines(xs, ys, id = "lines-1", xlim = c(-30, 30), ylim = c(-30, 40))
  Sys.sleep(0.025)
}
# Change to x-z solution plane
*plot(x, z, id = "ID-1", xlim = c(-30, 30), ylim = range(zs), transition = TRUE)
*lines(xs, zs, id = "lines-1", xlim = c(-30, 30), ylim = range(zs), transition = TRUE)
```
]

---

## Applications - Lorenz system

---

## Applications - Predator-Prey process

---

## Applications - Predator-Prey process

### The model

.medium[
- Modelled as a 2d space
- Foxes
  - reproduce every 300 time steps
  - have energy that lasts for ~500 time steps and die when energy goes to 0.
  - eat rabbits, each of which refills 200 energy
  - die of old age after 1000 time steps (w.p. 0.1 for each time step)
  - disappear when they die
- Rabbits
  - reproduces every 150 time steps
  - have infinite energy and outrun fox
  - die of old age after 500 time steps (w.p. 0.1 for each time step)
  - become corpse when they die
]

---

## Applications - Predator-Prey process

### Code skeleton

```r
n <- 100
animal <- sample(c("rabbit", "fox"), size = n, replace = T, prob = c(0.8, 0.2))
```

Use dots to represent foxes and rabbits:

```r
points(x = runif(n), y = runif(n), xlim = c(0, 1), ylim = c(0, 1),
       bg = ifelse(animal == "rabbit", "pink", "brown"))
```

Use mini-icons:

```r
rabbit_icon <- "https://bit.ly/2UwGxxj"
fox_icon <- "https://bit.ly/3D2sL6Q"
image(x = runif(n), y = runif(n), xlim = c(0, 1), ylim = c(0, 1),
      width = 20, height = 20, href = ifelse(animal == "rabbit", rabbit_icon, fox_icon))
```

---

## Applications - Predator-Prey process

---

## Applications - Schelling's model of segregation

### The model

.medium[
- Originally a model of society in social science (1971, 1978)  
(later also adapted to physics and biology)

- Assume a `$n \times n$` grid and `$m$` agents split into two groups; each agent occupies one space

- Stepping forward in time,
  - if an agent is surrounded by more than `$x$`% (e.g. 50%) of its own kind, it will stay;
  - otherwise, it will move to a random empty space.
  
- The question of interest is: __how segregated would this society be?__
]

---

## Applications - Schelling's model of segregation

### Code skeleton

```r
# Use canvas of size 500 x 500
n <- 30
id <- new_id(1:n^2, prefix = "agent")
race <- sample(c("red", "blue", "white"), prob = c(0.4, 0.4, 0.1), size = n^2, replace = TRUE)
grid <- data.frame(expand.grid(x = 1:n, y = 1:n), race = race, id = id)

*points(x = grid$x, y = grid$y, id = grid$id, pch = "square", cex = 140, bg = grid$race)
```

.medium[
Nicky Case has a playable implementation @ [https://ncase.me/polygons/](https://ncase.me/polygons/)
]

---

## Applications - Schelling's model of segregation

---

## Planned work / milestones

.medium[
1. Make plots detachable from R after the first run

2. Integration with R Markdown Document

3. Integration with Shiny

4. Provide a more complete support to `base` primitives
]