# World Region Maps
# Michael Minn
# 18 January 2020

# Encapsulate map plotting into a reusable function

library(sp)
library(rgeos)
library(maptools)


worldmap = function(countries, datacolumn, categories, colors, legend, filename) {

	# NA as white boundaries only

	countries[is.na(countries@data[,datacolumn]), datacolumn] = ""

	categories = c(categories, "")

	colors = c(colors, "#ffffff")


	# Labels

	labels = gUnaryUnion(countries, id=countries@data[,datacolumn])

	labels = gCentroid(labels, byid=T) 


	# Robinson projection

	robinson = CRS("+proj=robin +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs")

	countries = spTransform(countries, robinson)

	labels = spTransform(labels, robinson)

	# Graticule

	wgs84 = CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")

	latlines= Lines(sapply(seq(-90, 90, 10), function(y) Line(t(sapply(seq(-180, 180, 10), function(x) c(x, y))))), 1)
	longlines = Lines(sapply(seq(-180, 180, 10), function(x) Line(t(sapply(seq(-90, 90, 10), function(y) c(x, y))))), 2)

	graticule = SpatialLines(c(latlines, longlines), proj4string = wgs84)

	graticule = spTransform(graticule, robinson)

	# Colors based on categories and palette

	col = colors[match(countries@data[,datacolumn], categories)]

	# Create output file and plot

	png(filename, width=800, height=450, pointsize=14)

	par(mar=c(0.3,0.3,0.3,0.3))

	plot(graticule, col="#c0c0c0")

	plot(countries, col=col, border="#a0a0a0", add=T)

	if (legend)
		legend(x="bottomleft", legend=categories[1:(length(categories) - 1)], 
			col=colors[1:(length(colors) - 1)], pch=15, pt.cex=2)
	else
		text(coordinates(labels), labels=rownames(labels@coords))

	dev.off()
} # worldmap()


# Color palettes

darkpalette = c("#000080", "#ffcc00", "#004802", "#800000", "#cccccc", "#660066", 
	"#fefbc6", "#666699", "#ff0000", "#969696", "#dcdcff", "#a5bf88")

names(darkpalette) = c("Dark Blue", "Yellow", "Dark Green", "Dark Red", "Light Grey", "Purple",
	"Light Yellow", "Slate", "Bright Red", "Dark Grey", "Light Blue", "Olive")

pastelpalette = c("#808080", "#ff8080", "#80ff80", "#ffff80", 
	"#8080ff", "#ff80ff", "#80ffff", "#ffffff",
	"#c08080", "#80c080", "#c0c080", "#8080c0",
	"#c080c0", "#80c0c0", "#c0c0c0")

names(pastelpalette) = c("Grey", "Red", "Green", "Yellow", 
	"Blue", "Magenta", "Cyan", "White",
	"Salmon", "Olive", "Tan", "Cornflower",
	"Purple", "Sky Blue", "Dark Gray")


# Load and merge data

data("wrld_simpl")

countries = read.csv("world-regions.csv", as.is=T, check.names=F)

countries = merge(wrld_simpl, countries, by="ISO3")


# UN Continent / Subregion Map 
# https://unstats.un.org/unsd/methodology/m49/

datacolumn = "UNSubRegionLabel"
categories = c(
"North America", "Latin America and the Caribbean",
"Sub-Saharan Africa", "North Africa",
"Northern Europe", "Eastern Europe", "Southern Europe", "Western Europe",
"W. Asia", "S. Asia", "Southeast Asia", "Eastern Asia", "Central Asia",
"Antarctica",
"Australia and New Zealand", "Polynesia", "Melanesia", "Micronesia")
colors = pastelpalette[c(
	"Blue", "Olive",
	"Yellow", "Tan",
	"Sky Blue", "Cornflower", "Sky Blue", "Cyan",
	"Red", "Purple", "Salmon", "Red", "Salmon",
	"Grey",
	"Magenta", "Purple", "Magenta", "Purple")]
filename = "world-regions-un-subregions.png"
worldmap(countries, datacolumn, categories, colors, F, filename)


# Continent maps

datacolumn = "FourContinent"
categories = c("Americas", "Afro-Eurasia", "Antarctica", "Australia")
colors = pastelpalette[c("Blue", "Red", "Grey", "Magenta")]
filename = "world-regions-four-continent.png"
worldmap(countries, datacolumn, categories, colors, F, filename)

datacolumn = "FiveContinent"
categories = c("Americas", "Africa", "Eurasia", "Antarctica", "Australia")
colors = pastelpalette[c("Blue", "Yellow", "Red", "Grey", "Magenta")]
filename = "world-regions-five-continent.png"
worldmap(countries, datacolumn, categories, colors, F, filename)

datacolumn = "SevenContinent"
categories = c("North America", "South America", "Africa", "Europe", "Asia", "Antarctica", "Oceana")
colors = pastelpalette[c("Blue", "Green", "Yellow", "Cyan", "Red", "Grey", "Magenta")]
filename = "world-regions-seven-continent.png"
worldmap(countries, datacolumn, categories, colors, F, filename)


# World Regions (Finlayson 2019)
# 

datacolumn = "Finlayson"
categories = c("North America", "Middle and South America", "North Africa and Southwest Asia",
	"Sub-Saharan Africa", "Europe", "Russia", "South Asia", "East and Southeast Asia",
	"Australia and the Pacific")
colors = pastelpalette[c("Blue", "Green", "Tan", "Yellow", "Cyan", "Sky Blue", "Red", "Salmon", "Magenta")]
filename = "world-regions-finlayson.png"
worldmap(countries, datacolumn, categories, colors, T, filename)

# World Bank Economic (Administrative) Regions

datacolumn = "WorldBankRegionLabels"
categories = c("North America", "Latin America & Caribbean", "MENA",
	"Sub-Saharan Africa", "Europe & Central Asia", "East Asia & Pacific", "South Asia")
colors = pastelpalette[c("Blue", "Green", "Tan", "Yellow", "Cyan", "Salmon", "Red")]
filename = "world-regions-world-bank.png"
worldmap(countries, datacolumn, categories, colors, F, filename)


# World Bank Economic Income Group Regions

datacolumn = "WorldBankIncomeGroup"
categories = c("High income", "Upper middle income", "Lower middle income", "Low income")
colors = c("#6060ff", "#4040c0", "#202080", "#000060")
filename = "world-regions-world-bank-income.png"
worldmap(countries, datacolumn, categories, colors, T, filename)


# Organisation for Economic Co-operation and Development (OECD) countries

countries$OECDAccession[!is.na(countries$OECDAccession)] = "OECD"
countries$OECDAccession[is.na(countries$OECDAccession)] = "Non-OECD"

datacolumn = "OECDAccession"
categories = c("OECD", "Non-OECD")
colors = pastelpalette[c("Blue", "Red")]
filename = "world-regions-oecd.png"
worldmap(countries, datacolumn, categories, colors, T, filename)

# US Military Unified Combatant Command
# https://www.defense.gov/Our-Story/Combatant-Commands/
# https://en.wikipedia.org/wiki/Unified_combatant_command

datacolumn = "CombatantCommand"
categories = c("NORTHCOM", "SOUTHCOM", "CENTCOM", "AFRICOM", "EUCOM", "INDOPACOM")
colors = pastelpalette[c("Blue", "Green", "Tan", "Yellow", "Cyan", "Red")]
filename = "world-regions-combatant-command.png"
worldmap(countries, datacolumn, categories, colors, F, filename)