1

The Guardian plotted homicide locations on top of city maps here:

https://www.theguardian.com/us-news/ng-interactive/2017/jan/09/special-report-fixing-gun-violence-in-america

Given lat/lon data for homicides, I'd like to create a census-block-group scatterplot of population density versus homicide prevalence.

I'm familiar with R, python, and linux commandline tools, and wouldn't mind learning some QGIS.

Having never done anything like this, there is probably a better way (?), but my thinking is to convert homicide prevalence to a heatmap, then calculate the average heatmap value over the area of each census block group.

Are there GIS utilities/libraries that would simplify the task?

  • I think you need to determine the block grp for each homicide and then join those totals to the pop for the block grp and form your scatter plot from that data set – Ian Turton Mar 21 '17 at 14:06
2

Here's a solution that I use often as an example. This has been adapted from Chris Brunsdon and Lex Comber. Definitely a useful resource for spatial analysis in R.

#Load the necessary libraries and packages
library(GISTools)
library(rgdal)
library(sp)

#Loading New Haven Data (as an example)
data(newhaven)

#Visualizing blocks
blocks

##Produce a map of each Census block in New Haven
choropleth(blocks, blocks$P_VACANT) #Plotting blocks and shading them by the percent vacant blocks

#You can change the shading as well as add a legend if you want
shades <- auto.shading(blocks$P_VACANT, cols= brewer.pal(5,"Greens")) #Suppose I want Green
choropleth(blocks,blocks$P_VACANT, shading = shades)
choro.legend(533000,161000,shades)

enter image description here

#Note you can convert the blocks to a dataframe and plot scatterplots etc
 blocks <- as.data.frame(blocks)

####For heatmaps, I am using an example from the meuse dataset
data("meuse.grid")

 #Convert it to a spatial Pixels data frame and use spplot to heatmap the dist data
 meuse.grid <- SpatialPixelsDataFrame(points=meuse.grid[c("x","y")], data = meuse.grid)
 spplot(meuse.grid, "dist", col.regions=terrain.colors(20))

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.