Narco killings map of Mexico

narco-killings.R

#Download the narco killings map compiled by Walter Mckay available
#at LEAP Mexico https://sites.google.com/site/policereform/leap-mexico/narco-killings
library(rgdal)
library(gtools)
library(maptools)
library(ggplot2)
library(maps)
library(Cairo)
library(stringr)
library(MASS)
library(gpclib)
gpclibPermit()

urlMap <- "http://maps.google.com/maps/ms?ie=UTF8&hl=en&source=embed&msa=0&output=kml&msid=101594476481747955479.000486b8afd1e2a371d55"
download.file(urlMap, "Narco Killings.kml")

narco <- readOGR("Narco Killings.kml","Narco Killings")
narco <- (data.frame(narco))

names(narco) <- c("deaths", "Description", "long", "lat", "coords")
cleanData <- defmacro(str1, str2,
                      expr = {
                          narco$deaths <- gsub(str1, str2,
                                               narco$deaths,
                                            ignore.case = TRUE)
                          })
#Clean the data
cleanData("un", "1")
cleanData("doble", "2")
cleanData("seis", "6")
cleanData("ex comandante y su esposa", "2")
cleanData("Director of corrections murdered", "2")
cleanData("[0-9]+ injured", "0")
cleanData("[0-9]+ lesionadas", "0")
cleanData("Placemark [0-9]+", "0")

narco$involved <- "civilians"
narco$involved[grep("poli", narco$deaths, ignore.case = TRUE)] <- "police"
narco$involved[grep("milit", narco$deaths, ignore.case = TRUE)] <- "military"
narco$involved[grep("milit", narco$Description, ignore.case = TRUE)] <- "military"

#Add the numbers present in the data
x <- sapply(narco$deaths, function(x) strsplit(x, "[a-zA-Z ]"))
narco$value <- sapply(x, function(x) sum(as.numeric(x), na.rm = TRUE))

#Fix the narco fosa in Monterrey and the 72 immigrantes
narco$involved[which(narco$value == 75)] <- "civilians"
narco$involved[which(narco$value == 50)] <- "civilians"

narco <- narco[order(narco$value),]

write.csv(narco[,c(3,4,7)], "narco.csv", row.names = FALSE)
write.csv(narco, "narco-all.csv", row.names = FALSE)

narco$date <- as.Date(str_sub(narco$Description, 16, 23),
                      format = "%d-%m-%y")

#summary statistics
summary(narco$value)
qplot(narco$value)
ggplot(narco, aes(value, group = involved, fill = involved)) +
  geom_histogram()
ggplot(subset(narco, date > as.Date("2010-05-11")),
              aes(date, value, group = involved,
                  color = involved)) +
  geom_point() +
  scale_x_date()
sum(narco$value)

#A simple map
world.map <- map_data("world")
mexico <- subset(world.map, region == "Mexico")

#Uncomment if you want to use the INEGI map
#mexico <- readOGR("ESTADOS-90.shp", "ESTADOS-90")
#Create a ESTADOS-90.prj with the following txt (proj.4):
#PROJCS["North_America_Lambert_Conformal_Conic",GEOGCS["ITRF92",DATUM["D_ITRF_1992",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",2500000.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-102.0],PARAMETER["Standard_Parallel_1",17.5],PARAMETER["Standard_Parallel_2",29.5],PARAMETER["Latitude_Of_Origin",12.0],UNIT["Meter",1.0]]
#mexico <- spTransform(mexico, CRS("+proj=longlat +ellps=WGS84"))
#mexico <- fortify(mexico, region = "CVE_ENT")

Cairo("png", file = "mexico-drug-killings.png",
      width = 800, height = 600 )
ggplot(mexico, aes(long, lat)) +
         geom_polygon(aes(group = group), fill = "white",
                      color = "gray40", size = .2) +
         geom_jitter(data = narco, aes(long, lat, size = value,
                    color = involved),
                    alpha = .6) +
         scale_y_continuous(breaks = NA) +
         scale_x_continuous(breaks = NA) + xlab("") + ylab("") +
         opts(panel.background = theme_rect(fill = "white",
                                            colour = "white")) +
         opts(title = "Location of Drug-Related Homicides in Mexico (May-Oct 2010)") +
         coord_map(projection = "mercator") +
         scale_size("number\nkilled")
dev.off()