Карта административных регионов страны с помощью ggmap и ggplot2

Я могу сделать график уровня безработицы в США с помощью следующего кода.

library(XML)
library(ggplot2)
library(plyr)
library(maps)

unemp <-
  readHTMLTable('http://www.bls.gov/web/laus/laumstrk.htm',
    colClasses = c('character', 'character', 'numeric'))[[2]]

names(unemp) <- c('rank', 'region', 'rate')
unemp$region <- tolower(unemp$region)

us_state_map <- map_data('state')
map_data <- merge(unemp, us_state_map, by = 'region')

map_data <- arrange(map_data, order)

states <- data.frame(state.center, state.abb)

p1 <- ggplot(data = map_data, aes(x = long, y = lat, group = group))
p1 <- p1 + geom_polygon(aes(fill = cut_number(rate, 5)))
p1 <- p1 + geom_path(colour = 'gray', linestyle = 2)
p1 <- p1 + scale_fill_brewer('Unemployment Rate (Jan 2011)', palette  = 'PuRd')
p1 <- p1 + coord_map()
p1 <- p1 + geom_text(data = states, aes(x = x, y = y, label = state.abb, group = NULL), size = 2)
p1 <- p1 + theme_bw()
p1

Теперь я хочу подобный вид графика для Пакистана. Результаты моих нескольких попыток приведены ниже:

data(world.cities)
Pakistan <- data.frame(map("world", "Pakistan", plot=FALSE)[c("x","y")])

p <- ggplot(Pakistan, aes(x=x, y=y)) +
     geom_path(colour = 'green', linestyle = 2) +
     coord_map() + theme_bw()
p <- p + labs(x=" ", y=" ")
p <- p + theme(panel.grid.minor=element_blank(), panel.grid.major=element_blank())
p <- p + theme(axis.ticks = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank())
p <- p + theme(panel.border = element_blank())
print(p)

а также

library(mapproj)

Country <- "Pakistan"

Get_Map_Country <-
  get_map(
      location = Country
    , zoom = 5
    , scale = "auto"
    , maptype = "roadmap"
    , messaging = FALSE
    , urlonly = FALSE
    , filename = "ggmapTemp"
    , crop = TRUE
    , color = "color"
    , source = "google"
    , api_key
    )

Country1 <-
  ggmap(
      ggmap = Get_Map_Country
    , extent = "panel"
  #  , base_layer
    , maprange = FALSE
    , legend = "right"
    , padding = 0.02
    , darken = c(0, "black")
    )

Country1 <- Country1 + labs(x="Longitude", y="Latitude")
print(Country1)

Country2 <- Country1 + geom_polygon(data = Pakistan
                    , aes(x=x, y=y)
                    , color = 'white', alpha = .75, size = .2)

print(Country2)

Вопросы

Интересно, как получить карту административных районов Пакистана как США? Я знаю, что для этого нам нужны долгота и широта административных границ. Мне интересно, как получить долготу и широту административных границ для страны. Я пробовал глобальные административные зоны, но безуспешно. Любая помощь в этом отношении будет высоко оценена. Спасибо

2 ответа

Решение

Я не знаю необходимого вам пространственного уровня административных областей, но вот два способа считывания данных шейп-файлов и форматов.RData из глобальных административных областей (gadm.org) и преобразования их в фреймы данных для использования в ggplot2. Кроме того, для того, чтобы скопировать карту США, вам нужно будет построить названия административных областей, расположенных на центроидах многоугольника.

library(ggplot2)
library(rgdal)

Способ 1. SpatialPolygonDataFrames хранится в формате.RData

# Data from the Global Administrative Areas
# 1) Read in administrative area level 2 data

load("/Users/jmuirhead/Downloads/PAK_adm2.RData")
pakistan.adm2.spdf <- get("gadm")

Method2. Формат шейп-файла читается с помощью rgdal::readOGR

pakistan.adm2.spdf <- readOGR("/Users/jmuirhead/Downloads/PAK_adm", "PAK_adm2", 
 verbose = TRUE, stringsAsFactors = FALSE)

Например, создание data.frame из пространственного PolygonDataframes и объединение с data.frame, содержащим информацию о безработице.

pakistan.adm2.df <- fortify(pakistan.adm2.spdf, region = "NAME_2")

# Sample dataframe of unemployment info
unemployment.df <- data.frame(id= unique(pakistan.adm2.df[,'id']),
  unemployment = runif(n = length(unique(pakistan.adm2.df[,'id'])), min = 0, max = 25))

pakistan.adm2.df <- merge(pakistan.adm2.df, unemployment.df, by.y = 'id', all.x = TRUE)

Извлечение названий и центоидов административных областей для построения

# Get centroids of spatialPolygonDataFrame and convert to dataframe
# for use in plotting  area names. 

pakistan.adm2.centroids.df <- data.frame(long = coordinates(pakistan.adm2.spdf)[, 1], 
   lat = coordinates(pakistan.adm2.spdf)[, 2]) 

# Get names and id numbers corresponding to administrative areas
pakistan.adm2.centroids.df[, 'ID_2'] <- pakistan.adm2.spdf@data[,'ID_2']
pakistan.adm2.centroids.df[, 'NAME_2'] <- pakistan.adm2.spdf@data[,'NAME_2']

Создайте ggplot с метками для административных областей

p <- ggplot(pakistan.adm2.df, aes(x = long, y = lat, group = group)) + geom_polygon(aes(fill = cut(unemployment,5))) +
geom_text(data = pakistan.adm2.centroids.df, aes(label = NAME_2, x = long, y = lat, group = NAME_2), size = 3) + 
labs(x=" ", y=" ") + 
theme_bw() + scale_fill_brewer('Unemployment Rate (Jan 2011)', palette  = 'PuRd') + 
coord_map() + 
theme(panel.grid.minor=element_blank(), panel.grid.major=element_blank()) + 
theme(axis.ticks = element_blank(), axis.text.x = element_blank(), axis.text.y = element_blank()) + 
theme(panel.border = element_blank())

print(p)

Пример выходных данных по безработице в Пакистане

Вы пробовали http://www.diva-gis.org/gdata? Это должно дать вам шейп-файлы с границами штата / округа для Пакистана.

В пакете maptools есть функции, которые позволяют вам читать и преобразовывать шейп- файлы в фрейм данных, который можно использовать для разметки граничных полигонов.

Другие вопросы по тегам