Step-by-step plotting choropleth map of Nepal

map.choropleth.nepal-updated.R

library(rgdal)
library(ggplot2)
library(dplyr)

# clone NepalMaps from https://github.com/anjesh/NepalMaps 

# read shapefile
nepal_shp <- readOGR(dsn="NepalMaps/baselayers/NPL_adm", layer="NPL_adm3", stringsAsFactors = FALSE)
# fortify shapefile data to data frame
shp_df <- fortify(nepal_shp, region = "NAME_3")

# update the HPI data in districts.csv as per http://data.opennepal.net/content/human-poverty-index-value-districts-2011
hpi.csvdata <- read.csv("https://github.com/opennepal/odp-poverty/raw/master/Human%20Poverty%20Index%20Value%20by%20Districts%20(2011)/data.csv", stringsAsFactors = FALSE)
hpi.data <- hpi.csvdata %>% 
  filter(Sub.Group == "HPI") %>% 
  select(id = District, HPI = Value)

# identify the mismatched districts by uncommenting 2 lines below
# unique(hpi.data$id[!hpi.data$id %in% unique(nepal.adm3.shp.df$id)])
# unique(unique(nepal.adm3.shp.df$id)[!unique(nepal.adm3.shp.df$id) %in% hpi.data$id])

# fix the mismatched districts
hpi.data$id[hpi.data$id == "Darchaula"] <- "Darchula"
hpi.data$id[hpi.data$id == "Kavre"] <- "Kavrepalanchok"

# merge hpi data with shapefile dataframe
shp_df <- shp_df %>% 
  left_join(hpi.data, by="id")

# finding the centers of all the district polygons. 
centroids <- rgeos::gCentroid(nepal_shp, byid = TRUE, id = unique(nepal_shp$NAME_3))
centroids_df <- centroids %>% 
  as.data.frame() %>% 
  mutate(label = row.names(.))

# showing district names for which HPI>40
centroids.selected <- centroids_df[centroids_df$label %in% (hpi.data[hpi.data$HPI>40,]$id),]

# setting up bare theme
theme_bare <- theme(
  axis.line = element_blank(), 
  axis.text.x = element_blank(), 
  axis.text.y = element_blank(),
  axis.ticks = element_blank(), 
  axis.title.x = element_blank(), 
  axis.title.y = element_blank(),
  legend.text=element_text(size=7),
  legend.title=element_text(size=8),
  panel.background = element_blank(),
  panel.border = element_rect(colour = "gray", fill=NA, size=0.5)
)

# plot map
ggplot(data = shp_df,aes(x = long, y = lat, group = group)) +
  geom_polygon(aes(fill = HPI), color = 'gray', size = 0.1) +
  ggtitle("Human Poverty Index Map") +
  guides(fill=guide_colorbar(title="HP Index")) + 
  scale_fill_gradient(high = "#e34a33", low = "#fee8c8", guide = "colorbar") +
  coord_fixed(1.3) +
  theme(legend.justification=c(0,-0.1), legend.position=c(0.05,0)) +
  with(centroids.selected, annotate(geom="text", x = x, y = y, label=label, size=2)) +
  theme_bare

map.choropleth.nepal.R

library(rgdal)
library(ggplot2)
library(dplyr)

# clone NepalMaps from https://github.com/anjesh/NepalMaps 

# read shapefile
nepal.adm3.shp <- readOGR(dsn="./NepalMaps/baselayers/NPL_adm", layer="NPL_adm3", stringsAsFactors = FALSE)

# fortify shapefile data to data frame
nepal.adm3.shp.df <- fortify(nepal.adm3.shp, region = "NAME_3")

# write districts to csv file
nepal.adm3.shp.df %>%
  distinct(id) %>%
  write.csv("districts.csv", row.names = FALSE)

# update the HPI data in districts.csv as per http://data.opennepal.net/content/human-poverty-index-value-districts-2011
hpi.data <- read.csv("districts.csv")
colnames(hpi.data) <- c("id","HPI")

# merge hpi data with shapefile dataframe
nepal.adm3.shp.df <- merge(nepal.adm3.shp.df, hpi.data, by ="id")

# finding the centers of all the district polygons. Thanks http://stackoverflow.com/questions/28962453/how-can-i-add-labels-to-a-choropleth-map-created-using-ggplot2
centroids <- setNames(do.call("rbind.data.frame", by(nepal.adm3.shp.df, nepal.adm3.shp.df$group, function(x) {Polygon(x[c('long', 'lat')])@labpt})), c('long', 'lat')) 
centroids$label <- nepal.adm3.shp.df$id[match(rownames(centroids), nepal.adm3.shp.df$group)]

# showing district names for which HPI>40
centroids.selected <- centroids[centroids$label %in% (hpi.data[hpi.data$HPI>40,]$id),]

# setting up bare theme
theme_bare <- theme(
  axis.line = element_blank(), 
  axis.text.x = element_blank(), 
  axis.text.y = element_blank(),
  axis.ticks = element_blank(), 
  axis.title.x = element_blank(), 
  axis.title.y = element_blank(),
  legend.text=element_text(size=7),
  legend.title=element_text(size=8),
  panel.background = element_blank(),
  panel.border = element_rect(colour = "gray", fill=NA, size=0.5)
)

# plot map
ggplot(data = nepal.adm3.shp.df,aes(x = long, y = lat, group = group)) +
  geom_polygon(aes(fill = HPI), color = 'gray', size = 0.1) +
  ggtitle("Human Poverty Index Map") +
  guides(fill=guide_colorbar(title="HP Index")) + 
  scale_fill_gradient(high = "#e34a33", low = "#fee8c8", guide = "colorbar") +
  coord_fixed(1.3) +
  theme(legend.justification=c(0,-0.1), legend.position=c(0.05,0)) +
  with(centroids.selected, annotate(geom="text", x = long, y = lat, label=label, size=2)) +
  theme_bare

while using line 26 and 27 produces errors message and deleting 75 districts records but did indirectly
dd=read.csv(("C:/Users/yagyarimal/Desktop/publication/NepalMaps-master/districts.csv"))
centroids=data.frame(centroids,dd)
colnames(centroids)=c("long","lat","label")
head(centroids)

I have simplified the script https://gist.github.com/anjesh/64cd61e3f94352ce5249ef368082df0c#file-map-choropleth-nepal-updated-r You should be able to make changes to the code to make the centroids work as per your need

Thankx a lot Anjesh ! I am a beginner and am self-learning and this helps me a lot.

Thanks so much Anjesh. I worked through the steps and I got it right to a point where I got stocked.. I am working on Africa and I am only selecting 20 countries out of Africa.. I used a percapita GDP of these countries for a particular year and I set the other countries which is not among my study to zero.. How ever when I tried to plot the cloroplethr with label name.. It couldn't imprint names on the countries .. Showing geom_text deleted 768 missing and it returned to the normal map without labels..

You might need to install extra libraries as well.

install.packages('maptools')
install.packages('rgeos')
install.packages('rgdal')