Tuesday 23 February 2016

shapefile - How to plot points on maps using ggplot2 and R?


Thanks for help me in GIS 101 Problem #1, now I have geocoded a few hospital in Connecticut using google map, I have able to visualize them in QGIS, now I am trying to do it in R.


ct <- readShapeSpatial("housect_37800_0000_2010_s100_census_1_shp/wgs84/housect_37800_0000_2010_s100_census_1_shp_wgs84.shp")

ct_mod <- fortify(ct,region="SLDLST10")
chart <- ggplot(data=ct_mod,aes(long,lat,group=group))
chart <- chart + scale_x_continuous(limits=c(-73.8,-71.7),breaks=seq(-74,-71,0.1))
chart <- chart + scale_y_continuous(limits=c(40.9,42.1),breaks=seq(40,43,0.1))
chart <- chart + geom_polygon(fill="grey80")

chart <- chart + geom_path(color="white")
chart1 <- chart + geom_point(data=hosp.list,aes(x=coord_x,y=coord_y))
chart2 <- chart + geom_point(aes(x=-72,y=41))

My hosp.list is a dataframe with hospital name, address, lat, long.


Now I can show chart and chart2, but not chart1, because of the error Error in eval(expr, envir, enclos) : object 'group' not found, can anyone help?


update 01


Thanks for celenius's comment, I tried the following instead:


chart1 <- chart + geom_point(aes(x=hosp.list[,"coord_x"],y=hosp.list[,"coord_y"]))  


And this gives me an error below:


Error in data.frame(x = c(-72.4509975, -72.4509975, -92.059305, -73.16584,  : 
arguments imply differing number of rows: 58, 133341

What I am trying to do is as below:



  1. I plot the Connecticut 2010 Census State Legislative District as a layer

  2. I created a dataframe containing several hospitals in Connecticut, geocoded their addresses into long and lat, and put the long/lat into the data frame also (i.e. the dataframe has 3 column)

  3. I want to create a new layer on top of the "district layer", showing the hospitals by dots

  4. the next step I will do somethings like: display name of each dot, highlight the region with hospitals, etc... which should by me future question poster here.



Thanks again.



Answer



I fixed that:


ct <- readShapeSpatial("housect_37800_0000_2010_s100_census_1_shp/wgs84/housect_37800_0000_2010_s100_census_1_shp_wgs84.shp")

ct_mod <- fortify(ct,region="SLDLST10")
# chart <- ggplot(data=ct_mod,aes(long,lat,group=group)) # the group is the issue, should not be used here as the hosp.list will also be looked for group, which does not exist
chart <- ggplot(data=ct_mod,aes(long,lat))
chart <- chart + scale_x_continuous(limits=c(-73.8,-71.7),breaks=seq(-74,-71,0.1))

chart <- chart + scale_y_continuous(limits=c(40.9,42.1),breaks=seq(40,43,0.1))
# chart <- chart + geom_polygon(fill="grey80")
chart <- chart + geom_polygon(fill="grey80",aes(group=group)) # add group here
# chart <- chart + geom_path(color="white")
chart <- chart + geom_path(color="white",aes(group=group)) # add group here
chart1 <- chart + geom_point(data=hosp.list,aes(x=coord_x,y=coord_y))

The above script works for me.


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...