Density2d Plot using another variable for the fill (similar to geom_tile)?
I am trying to plot a map for my final project, and I am trying to do a heat map of crime by BLock in the US.
For each block, I have Lat, Lon, and a prediction of the crime rate. It follows this structure:
Lat / Lon / Prediction
-76.0 / 40.0 / 125
-76.120 / 40.5 / 145
-75.98 / 41.001 / 95
And so on.
Is there a way to plot a heat map showing the Prediction as the fill?
I think this is what geom_tiles do, but that geom is not working (maybe because the points are not evenly spaced)
Any help would be more than welcome. Please!
EDIT
This is what I have tried so far:
-geom_density2d:
ggplot(ny2,aes(x=GEO_CENTROID_LON,y=GEO_CENTROID_LON,fill=prediction))+geom_density2d()
Gives me the error: "Error in unit(tic_pos.c, "mm") : 'x' and 'units' must have length > 0"
-geom_tiles:
ggplot(ny2,aes(x=GEO_CENTROID_LON,y=GEO_CENTROID_LON,fill=prediction))+geom_tile()
Produces a plot with the proper scale, but not data shown on the map.
Regarding chloropeth, it would work if I happend to have block level information for the whole US, but I can't find such data.
SUBSAMPLE of data can be found here
First, let's load the data:
data<-read.csv(file = "NY subsample.csv")
Data points
Then, let's try just plotting the basic locations and values of the data:
require('ggplot2')
# start with points
pred.points <- ggplot(data = data,
aes(x = GEO_CENTROID_LON,
y = GEO_CENTROID_LAT,
colour = prediction)) +
geom_point()
print(pred.points)
ggsave(filename = "NYSubsamplePredPoints.png",
plot = p2,
scale = 1,
width = 5, height = 3,
dpi = 300)
which gives us this:
Binned data
Then, you can try to plot the mean in a 2-D region using stat_summary2d()
:
pred.stat <- ggplot(data = data,
aes(x = GEO_CENTROID_LON,
y = GEO_CENTROID_LAT,
z = prediction)) +
stat_summary2d(fun = mean)
print(pred.stat)
ggsave(filename = "NYSubsamplePredStat.png",
plot = pred.stat,
scale = 1,
width = 5, height = 3,
dpi = 300)
which gives us this plot of the mean value of prediction in each box.
Binned and with custom colormap and correct projection
Next, we can set the bin size, color scales, and fix the projection:
# refine breaks and palette ----
require('RColorBrewer')
YlOrBr <- c("#FFFFD4", "#FED98E", "#FE9929", "#D95F0E", "#993404")
pred.stat.bin.width <- ggplot(data = data,
aes(x = GEO_CENTROID_LON,
y = GEO_CENTROID_LAT,
z = prediction)) +
stat_summary2d(fun = median, binwidth = c(.05, .05)) +
scale_fill_gradientn(name = "Median",
colours = YlOrBr,
space = "Lab") +
coord_map()
print(pred.stat.bin.width)
ggsave(filename = "NYSubsamplePredStatBinWidth.png",
plot = pred.stat.bin.width,
scale = 1,
width = 5, height = 3,
dpi = 300)
which gives us this:
Plotted over a base map
And last of all, here's the data overlain on a map.
require('ggmap')
map.in <- get_map(location = c(min(data$GEO_CENTROID_LON),
min(data$GEO_CENTROID_LAT),
max(data$GEO_CENTROID_LON),
max(data$GEO_CENTROID_LAT)),
source = "osm")
theme_set(theme_bw(base_size = 8))
pred.stat.map <- ggmap(map.in) %+% data +
aes(x = GEO_CENTROID_LON,
y = GEO_CENTROID_LAT,
z = prediction) +
stat_summary2d(fun = median,
binwidth = c(.05, .05),
alpha = 0.5) +
scale_fill_gradientn(name = "Median",
colours = YlOrBr,
space = "Lab") +
labs(x = "Longitude",
y = "Latitude") +
coord_map()
print(pred.stat.map)
ggsave(filename = "NYSubsamplePredStatMap.png",
plot = pred.stat.map,
scale = 1,
width = 5, height = 3,
dpi = 300)
Setting the colormap
And finally, to set the colormap to something like http://www.cadmaps.com/images/HeatMapImage.jpg, we can take a guess at the colormap:
colormap <- c("Violet","Blue","Green","Yellow","Red","White")
and do the plotting again:
pred.stat.map.final <- ggmap(map.in) %+% data +
aes(x = GEO_CENTROID_LON,
y = GEO_CENTROID_LAT,
z = prediction) +
stat_summary2d(fun = median,
binwidth = c(.05, .05),
alpha = 1.0) +
scale_fill_gradientn(name = "Median",
colours = colormap,
space = "Lab") +
labs(x = "Longitude",
y = "Latitude") +
coord_map()
print(pred.stat.map.final)
I'm still not really sure what you're aiming at, but if you want a map with some points and contours plotted on it, that is possible. The output looks like the screenshot below; obviously there are many ways that could be tweaked. (The state shown is Connecticut.)
Code follows:
library(ggplot2)
library(maps)
mytext <- "GEOGRAPHY_ID GEO_CENTROID_LAT GEO_CENTROID_LON prediction
90012003024 41.364744 -73.373055 90.1
90012052001 41.488312 -73.381211 97.5
90012452002 41.300794 -73.467077 70.8
90012452003 41.319154 -73.479132 95.7
90012453001 41.293613 -73.456416 94.9
90012453003 41.284093 -73.493029 84.6
90012454001 41.263554 -73.462508 133.7
90010716002 41.181796 -73.196349 264.1
90010714003 41.183947 -73.198804 157.9
90010714004 41.185397 -73.200341 133.1
90010716001 41.184585 -73.195373 245.4
90010719001 41.195616 -73.203213 241.3
90010719002 41.191719 -73.20067 201.4
90010810001 41.222401 -73.155751 99.8
90012053002 41.437178 -73.396647 88.9
90012053003 41.451281 -73.413168 87.2
90012101001 41.400968 -73.458019 95.8
90012101002 41.394998 -73.449385 94.9
90012455002 41.270781 -73.52177 130
90012456001 41.360282 -73.530252 76
90012456002 41.288121 -73.507822 72.7
90012456003 41.313295 -73.528298 99.5
90010720002 41.187461 -73.206022 183.8
90010721001 41.189995 -73.216666 122.1
90010721003 41.180815 -73.216566 146.6
90010812001 41.233187 -73.117965 104.3
90012102003 41.393535 -73.440432 157.3
90012456004 41.296988 -73.522458 79.8
90012571001 41.613474 -73.503408 100.3
90012572001 41.203247 -73.202377 137.5
90012572002 41.194794 -73.192015 140
90010722002 41.195945 -73.210519 151.4
90010722003 41.193973 -73.214259 136.6
90010723001 41.203358 -73.20714 143.9
90010723002 41.198932 -73.206447 142.4
90010723003 41.203037 -73.212283 186.3
90010724001 41.208966 -73.201209 165.8
90010813001 41.246832 -73.131348 72.6
90010813003 41.246515 -73.111748 155.4
90012103003 41.402606 -73.440486 167.3
90012104002 41.390956 -73.434375 94.8
90012104004 41.381899 -73.431034 102.2
90012105001 41.369942 -73.437917 123.2
90012572003 41.195887 -73.197341 173.9
90012572004 41.198997 -73.199825 162.1
90010724002 41.209585 -73.20547 148.2
90010725001 41.212546 -73.215288 91.6
90010725002 41.208281 -73.212059 98.6
90010901002 41.287263 -73.250588 102.3
90010902003 41.266822 -73.23491 92
90012105002 41.353974 -73.458436 99
90012106001 41.383904 -73.465779 73.8
90012106002 41.382667 -73.457351 83.2
90012106003 41.387105 -73.458921 189.7
90010603003 41.211193 -73.280762 90.6
90010604001 41.1995 -73.306179 92.5
90010604002 41.181699 -73.295321 88.5
90010604003 41.171775 -73.283657 74.7
90010726004 41.222752 -73.235936 95.1
90010726005 41.211959 -73.228826 161.5
90010727001 41.225968 -73.207215 131.2
90010727002 41.214759 -73.208042 132.3
90012106004 41.390604 -73.45726 162.3
90012107012 41.399524 -73.464468 133.7
90012107013 41.396556 -73.471967 99.7
90012107021 41.386497 -73.468571 79
90010604004 41.167694 -73.309524 84.1
90010604005 41.195134 -73.325438 83.8
90010605001 41.155009 -73.2852 84.2
90010606001 41.139112 -73.283955 110.6
90010607002 41.147077 -73.261239 94.2
90010728001 41.224439 -73.198242 99.1
90010728002 41.217875 -73.197045 202.6
90010728003 41.210812 -73.195182 190.9
90010729001 41.216106 -73.185495 102.2
90010729002 41.221654 -73.19115 120.5
90010904004 41.228466 -73.191144 129.8
90010905001 41.245708 -73.153402 94.4
90012107022 41.390499 -73.47572 122
90012108003 41.395869 -73.494013 169
90010607004 41.162216 -73.266692 108.4
90010729003 41.208007 -73.192118 119.6
90010731002 41.209897 -73.173502 200.1
90010732001 41.199521 -73.157307 104
90010907002 41.292969 -73.213286 97.2
90012109002 41.422006 -73.509898 130.3
90012109003 41.408571 -73.469713 96.7
90012110002 41.427622 -73.478071 71
90010732002 41.199779 -73.15968 129.9
90010733001 41.194838 -73.161245 220.4
90010733002 41.196553 -73.167201 223.9
90010734001 41.200568 -73.177405 167.6
90010734003 41.202783 -73.185394 132.6
90011002001 41.320343 -73.199107 84.5
90011002002 41.3086 -73.220255 88
90011002004 41.331716 -73.221468 78.4
90012113001 41.445838 -73.445831 79.8
90010612001 41.179606 -73.224896 81.8
90010612002 41.181751 -73.232609 87.5
90010613001 41.175294 -73.231124 101
90010613002 41.173877 -73.239615 81.5
90010735001 41.196816 -73.182355 201
90010735002 41.193029 -73.183216 175.8
90010735003 41.19529 -73.184959 231.1
90010736001 41.194455 -73.177717 198.3
90011003002 41.341963 -73.192365 77.7
90011051001 41.233168 -73.263024 123.6
90011052002 41.285673 -73.287781 111.5
90012113003 41.441636 -73.436266 116.9
90012114001 41.430746 -73.428211 78.7
90012114002 41.447087 -73.427341 108.3
90012114003 41.422319 -73.420856 65.4
90010614002 41.167362 -73.237612 95.1
90010615001 41.154024 -73.24035 146
90010616001 41.137264 -73.264706 72.9
90010737003 41.190196 -73.163451 174.4
90010737004 41.191454 -73.168643 136.9
90010737005 41.188779 -73.168872 125.4
90010738001 41.190175 -73.175898 149.3
90010738002 41.188112 -73.176438 104.1
90011102011 41.30204 -73.07854 93.2
90011102021 41.292759 -73.089266 221.6
90011102022 41.277148 -73.094239 209
90012201003 41.450671 -73.524399 95.3
90012201004 41.464062 -73.509055 91.4
90012202001 41.494752 -73.495944 71.2
90010616002 41.131092 -73.255091 114.8
90010616003 41.126782 -73.269031 116.1
90010701001 41.155218 -73.22345 222
90010701003 41.1547 -73.228366 121.2
90010701004 41.158021 -73.229301 119.4
90010738003 41.189804 -73.178993 120.1
90010739001 41.190496 -73.1826 120
90010739002 41.187738 -73.18267 173.6
90010739003 41.186232 -73.18633 218.8
90010739004 41.190055 -73.185757 115.7
90010740001 41.182438 -73.180944 187.2
90010740002 41.181587 -73.17741 145.9
90012203001 41.521653 -73.470862 103.1
90012203002 41.487018 -73.45847 95.8
90012301002 41.442759 -73.311749 103.8
90010702001 41.161921 -73.22484 130.6
90010702002 41.163873 -73.22144 137.8
90010702003 41.15867 -73.219191 181.1
90010703001 41.160754 -73.21476 215.2
90010705001 41.167282 -73.191851 175
90010743001 41.18634 -73.159732 166.9
90010743002 41.182201 -73.158335 109
90010743003 41.180505 -73.162622 152.2
90010743004 41.178543 -73.168002 157.3
90010743005 41.181 -73.168064 143.4
90010743006 41.182831 -73.16868 122.4
90011103022 41.269651 -73.130548 88.6
90011104001 41.294908 -73.161871 111.3
90012302001 41.413873 -73.304276 93.1
90012303003 41.387196 -73.321197 118.5
90010705002 41.167522 -73.196366 213.2
90010706001 41.177971 -73.193067 236.6
90010706002 41.167732 -73.187421 227.7
90010709001 41.170459 -73.203877 159.7
90010709002 41.17138 -73.198856 186.7
90010710001 41.171339 -73.218595 126.6
90010744001 41.179957 -73.157483 135.6
90010744002 41.177069 -73.162353 133.4
90010744003 41.172159 -73.168866 185.4
90010744004 41.176011 -73.16636 169.5
90010801001 41.200673 -73.150037 105.6
90010801003 41.193386 -73.151688 140.3
90011105004 41.309581 -73.141151 77.3
90012304001 41.384815 -73.294462 118.4
90012304002 41.36671 -73.293527 94.8
90012304003 41.368872 -73.348572 194.7
90012305011 41.406865 -73.223661 89.3
90010710002 41.172463 -73.214353 184.7
90010712001 41.180827 -73.207841 133.3
90010712002 41.179002 -73.205237 206.9
90010801004 41.199783 -73.15404 140.2
90010802001 41.190012 -73.142391 106.3
90010802002 41.189527 -73.147431 86.6
90010802003 41.188494 -73.152008 121.6
90010804001 41.182524 -73.140233 153.5
90010804002 41.177371 -73.137682 151.9
90012002001 41.381013 -73.410319 87.5
90012305021 41.375831 -73.251305 105.6
90012305023 41.379858 -73.219021 96.8
90010712003 41.177009 -73.202378 168.4
90010712004 41.173824 -73.203364 197.4
90010713001 41.18072 -73.201927 184
90010713002 41.178332 -73.198693 183.5
90010714001 41.188573 -73.196834 115.9
90010714002 41.185972 -73.201774 118.1
90010805002 41.152953 -73.123918 71.8
90010806001 41.185641 -73.129554 112.4
90012003021 41.36255 -73.396786 86.7
90012003022 41.351419 -73.408962 154
90012402002 41.317929 -73.400807 91.5
90012402003 41.300946 -73.348506 126.8
90010202001 41.134991 -73.600739 118.3
90010202002 41.114341 -73.58158 84.7
90010203001 41.154604 -73.545957 90.1
90010203003 41.15133 -73.581789 109.2
90010353001 41.126492 -73.477104 191.2
90010353002 41.129396 -73.513016 221
90010354001 41.192382 -73.482245 120.4
90010354002 41.160459 -73.474248 114.7
90010354003 41.144611 -73.464233 216.9
90010354004 41.186241 -73.495502 88.6
90010503004 41.149218 -73.336642 192.9
90010503005 41.169512 -73.351917 140.6
90010504002 41.116467 -73.375013 85.9
90010505002 41.116352 -73.352404 85.4
90010203004 41.136534 -73.573863 86.2
90010205002 41.066638 -73.562444 77.2
90010425001 41.160707 -73.404989 125.9
90010425003 41.144366 -73.403332 121.1
90010426003 41.133452 -73.385111 81.7
90010505003 41.111412 -73.358066 80.8
90010505004 41.129533 -73.36223 82.9
90010551001 41.233676 -73.362879 108.9
90010206001 41.103581 -73.554734 85.2
90010206003 41.085439 -73.552623 82.5
90010426004 41.12787 -73.392559 136.8
90010427001 41.149039 -73.424079 107
90010428002 41.131288 -73.416989 91
90010551003 41.21145 -73.381718 77.7
90010551004 41.247205 -73.406188 77.1
90010552001 41.213714 -73.338138 108.3
90010552003 41.201885 -73.365519 115.7
90010208003 41.087838 -73.544622 77.2
90010209001 41.098298 -73.515074 83.6
90010428004 41.127102 -73.399138 85.3
90010429001 41.149545 -73.437497 77
90010430002 41.12977 -73.433261 106.9
90010602003 41.193534 -73.252047 60.8
90010602004 41.187197 -73.257167 141.4
90010210001 41.086533 -73.521131 112.1
90010211001 41.072121 -73.517654 101.2
90010211003 41.075939 -73.52826 78.2
90010431001 41.102395 -73.446269 83.8
90010431002 41.1173 -73.459513 73.3
90010432002 41.110378 -73.442314 92.3
90010212002 41.070902 -73.543556 89.6
90010213001 41.070755 -73.553697 110.1
90010213003 41.057543 -73.550212 98.6
90010433002 41.120549 -73.433645 95.1
90010434001 41.124754 -73.414431 91.3
90010434002 41.122855 -73.419197 87.9
90010434003 41.129448 -73.422258 98.3
90010435001 41.125176 -73.386176 144.2
90010214001 41.053676 -73.554293 202.3
90010214002 41.050584 -73.555773 176.2
90010214003 41.04669 -73.557427 103
90010214004 41.053666 -73.562879 172.1
90010215001 41.053793 -73.54926 145.9
90010215002 41.050154 -73.548156 139.8
90010435002 41.120456 -73.394265 87.9
90010435003 41.117986 -73.388516 78.1
90010436002 41.111907 -73.404062 183.2
90010437001 41.117371 -73.414574 202.5
90010215003 41.046584 -73.550538 148.7
90010215004 41.050229 -73.551746 126
90010216001 41.063583 -73.538563 185.4
90010216002 41.060054 -73.543822 132.2
90010216004 41.062798 -73.535453 124.6
90010437002 41.11312 -73.413212 213.7
90010438001 41.114109 -73.420907 167.6
90010438002 41.114952 -73.424774 122.8
90010438003 41.109739 -73.423166 149.8
90010438004 41.10928 -73.42707 137.6
90010101023 41.117859 -73.646129 87.8
90010102013 41.061299 -73.62279 157
90010217001 41.056839 -73.527159 171.4
90010217003 41.055427 -73.533343 129.1
90010217004 41.057926 -73.534384 138.6
90010217005 41.057235 -73.529757 104.5
90010218011 41.067416 -73.5255 144.2
90010439001 41.097703 -73.435218 153.8
90010439004 41.08607 -73.443137 68.4
90010440001 41.10246 -73.422491 199.8
90010102022 41.05194 -73.591855 97.6
90010103002 41.03803 -73.628652 87.4
90010218012 41.063118 -73.529794 118.4
90010218013 41.062331 -73.523699 120.1
90010218023 41.059329 -73.520307 107.8
90010219001 41.057436 -73.508325 96.8
90010440002 41.102279 -73.426726 201.9
90010440004 41.096429 -73.424339 166.7
90010440005 41.091685 -73.428074 265.9
90010441001 41.09657 -73.419621 174.2
90010441002 41.095089 -73.417375 165.8
90010103005 41.029176 -73.648761 130.3
90010220002 41.053379 -73.518227 148.5
90010221001 41.051607 -73.522784 157.5
90010442002 41.104538 -73.408893 154.1
90010442003 41.100721 -73.407695 121.1
90010105004 41.018874 -73.641155 144.3
90010221002 41.046457 -73.522742 133.3
90010221003 41.039903 -73.52882 126.5
90010223001 41.038742 -73.54947 135.6
90010223002 41.0331 -73.548214 78.3"
ny2 <- read.table(textConnection(mytext), sep = " ",
check.names = FALSE,
strip.white = TRUE,
header = TRUE)
names(ny2) <- c('id','lat','long','prediction')
us.states <- map_data("state")
plot.states <- c('connecticut') # subset to Connecticut for ease of use in this example
p <- ggplot() +
geom_polygon(data = us.states[us.states$region %in% plot.states, ],
aes(x = long, y = lat, group = group),
colour = "white", fill = "lightblue") +
geom_point(data = ny2, aes(x = long, y = lat,
colour = prediction), group = NULL) +
stat_density2d(data = ny2, aes(x = long, y = lat), size = 0.02) +
coord_map() +
theme()
print(p)