plotting urban population time series.R

library(gganimate)
library(magick)
library(scales)
library(sf) 
library(tidyverse)


# Importing the India states shapefile into a simple feature, and reducing
# number of vertices for quick screen drawing
india = st_read("India/States.shp") %>%
    st_simplify(preserveTopology = T, dTolerance = 0.05)


# custom function to clean variable names upon import
clean_var_names <- function(var_name) {
    return(str_replace_all(tolower(var_name), "\\s", "_"))
}


# Creates a cleaned up tidy data frame (tibble) from the raw data.
# Steps described inline below...
data <- 
    # read CSV file
    read_csv("WUP-data.csv", col_types = cols()) %>%
    # rename imported variable names by converting to lower case and
    # replacing spaces with underscores
    rename_all(clean_var_names) %>% 
    # select only rows pertaining to India
    filter(country_or_area == "India") %>% 
    # selects columns of interest and renames 3 of them...
    select(
        city = urban_agglomeration,
        lon = longitude,
        lat = latitude,
        starts_with("19"),
        starts_with("20")) %>% 
    # converts data from wide to long format, which is ideal for programmatic
    # plotting of a time series
    gather(
        key = year, 
        value = population, 
        -c(city, lon, lat)) %>% 
        mutate(
        year = as.integer(year),
        population = as.numeric(str_replace_all(population, " ", ""))
    )


# defining a set of object vectors for use with plotting below
# this is a vector of years for which plots are needed
years <- data %>% select(year) %>% distinct() %>% pull()  
# this is a vector of 'break' values for creating class intervals
breaks <- c(100, 250, 500, 750, 1000, 2500, 5000, 7500, 10000, 45000)
# this is a vector of labels that will appear on the map legend
labels <- c("< 1.0L", "1.0L - 2.5L", "2.5L - 5.0L", "5.0L - 7.5L", "7.5L - 1.0M", "1.0M - 2.5M", "2.5M - 5.0M", "5.0M - 7.5M", "7.5M - 10M", "> 10M")


# create a custom function to make a map for any given year
make_plot <- function(yr) {
p <- data %>%
    filter(year==yr) %>%
    ggplot() +
    geom_sf(data = india, fill="grey25", col="white", size=0.1) +
    geom_point(aes(x=lon, y=lat, size=population, col=population)) +
    scale_colour_distiller("Population", breaks = breaks, limits = c(0, 45000), labels=labels, palette = "Spectral", direction = -1, values = rescale(breaks)) +
    scale_size_area("Population", breaks = breaks, limits = c(0, 45000), labels=labels, max_size = 15) +
    theme_minimal(base_size = 18)+
    guides(colour= guide_legend(), size=guide_legend())+
    xlab("Longitude") + ylab("Latitude") +
    annotate("text", label = yr, x=90, y=35, size=10, col="dark red")
print(p)
}


# use function above to make a map for the year 2015;
make_plot(2015)


# This part is entirely optional, and relevant if you want to produce
# an animated GIF of the maps for each year, playing in sequence.
img = image_graph(width = 1280, height = 1024)
out = lapply(years, make_plot)
dev.off()
image_animate(img, 4) #this opens a viewer window from which you can either open in a browser or export GIF

# Alternate way of animating, using gganimate
p <- data %>%
    ggplot() +
    geom_sf(data = india, fill="grey25", col="white", size=0.1) +
    geom_point(aes(x=lon, y=lat, size=population, col=population)) +
    scale_colour_distiller("Population", breaks = breaks, limits = c(0, 45000), labels=labels, palette = "Spectral", direction = -1, values = rescale(breaks)) +
    scale_size_area("Population", breaks = breaks, limits = c(0, 45000), labels=labels, max_size = 15) +
    theme_minimal(base_size = 18)+
    guides(colour= guide_legend(), size=guide_legend())+
    labs(title = 'Year: {frame_time}', x = 'Latitude', y = 'Longitude') +
    transition_time(year) +
    ease_aes('linear')

animate(p, fps = 10)