Animation of the paleogeography of earth for the past 542 million years

2018-08_paleodem.Rmd

---
title: "Paleogeography of earth for the past 542 million years"
author: "Dewey Dunnington"
date: '2018-08-20'
output: github_document
---

```{r setup, include=FALSE}
library(rvest)
library(tidyverse)
library(stars)
library(gganimate)

# cleanup
unlink("2018-08_paleodem_PaleoDEMS_nc.zip")
unlink("2018-08_paleodem_PaleoDEMS_nc", recursive = TRUE)

# download data
# data page: http://www.earthbyte.org/paleodem-resource-scotese-and-wright-2018/
curl::curl_download("https://www.earthbyte.org/webdav/ftp/Data_Collections/Scotese_Wright_2018_PaleoDEM/Scotese_Wright_2018_Maps_1-88_1degX1deg_PaleoDEMS_nc.zip", "2018-08_paleodem_PaleoDEMS_nc.zip")
unzip("2018-08_paleodem_PaleoDEMS_nc.zip", exdir = "2018-08_paleodem_PaleoDEMS_nc")

knitr::opts_chunk$set(echo = TRUE, dpi = 300)
```

```{r, eval=FALSE}
library(tidyverse)
library(stars)
library(gganimate)

curl::curl_download(
  "https://www.earthbyte.org/webdav/ftp/Data_Collections/Scotese_Wright_2018_PaleoDEM/Scotese_Wright_2018_Maps_1-88_1degX1deg_PaleoDEMS_nc.zip", 
  "2018-08_paleodem_PaleoDEMS_nc.zip"
)
unzip(
  "2018-08_paleodem_PaleoDEMS_nc.zip", 
  exdir = "2018-08_paleodem_PaleoDEMS_nc"
)
```

```{r paleomap-present}
present <- read_stars("2018-08_paleodem_PaleoDEMS_nc/PaleoDEMS netcdf/Map01_PALEOMAP_1deg_Holocene_0.0Ma.nc")
plot(present)
```

```{r}
library(rvest)

time_scale_page <- read_html("https://en.wikipedia.org/wiki/Geologic_time_scale")
time_scale <- time_scale_page %>% 
  rvest::html_nodes("table.wikitable") %>%
  rvest::html_table(fill = TRUE) %>%
  first() %>%
  select(period = starts_with("Period"), start_ma = starts_with("Start,")) %>%
  extract(start_ma, "age_ma", ".*?([0-9.]+)", convert = TRUE) %>%
  mutate(period = if_else(str_detect(period, "^n/a|Phanero|Cenozo"), NA_character_, period) %>%
           str_replace("Carbon", "Carboniferous")) %>%
  fill(period) %>%
  extract(period, "period", "^([a-zA-Z]+)") %>%
  filter(age_ma <= 2500) %>%
  group_by(period) %>%
  summarise(max_age = max(age_ma)) %>%
  arrange(max_age) %>%
  mutate(min_age = lag(max_age, default = 0))

time_scale
```

```{r}
geo_period <- function(age_ma) {
  map_chr(
    age_ma,
    ~time_scale %>%
      filter(max_age > .x, min_age <= .x) %>%
      pull(period) %>%
      first()
  )
}

geo_period(c(-1:10))
```

```{r}
all <- list.files(
  "2018-08_paleodem_PaleoDEMS_nc/PaleoDEMS netcdf/", 
  ".nc$",
  full.names = TRUE
) %>%
  read_stars(quiet = TRUE) %>%
  st_set_crs(4326) %>% 
  st_transform(54030)

all_df <- all %>%
  as_tibble() %>%
  gather(map_name, z, starts_with("Map")) %>%
  extract(
    map_name,
    into = c("map_number", "age_ma"),
    regex = "Map([0-9.]+).*?([0-9.]+)Ma.nc",
    convert = TRUE
  )

head(all_df)
```

```{r paleo-dem-anim, cache=TRUE}
x_res <- all_df$x %>% unique() %>% sort() %>% diff() %>% unique() %>% first()
y_res <- all_df$y %>% unique() %>% sort() %>% diff() %>% unique() %>% first()

anim <- all_df %>%
  ggplot(aes(x, y, fill = z)) +
  geom_tile(width = x_res, height = y_res) +
  geom_hline(yintercept = 0, alpha = 0.5, size = 0.3) +
  scale_fill_gradient2(midpoint = 0, low = scales::muted("blue"), high = scales::muted("red")) +
  coord_fixed(expand = FALSE) +
  labs(
    title = "Paleo DEM at {sprintf('%0.0f', -frame_time)} Ma ({geo_period(-frame_time)})",
    fill = "Elevation",
    caption = "Data: earthbyte.org/paleodem-resource-scotese-and-wright-2018",
    x = NULL,
    y = NULL
  ) +
  theme(axis.text = element_blank(), axis.ticks = element_blank()) +
  transition_time(-age_ma)

animate(anim, nframes = 541, width = 600, height = 300, res = 96)
```