Voting Districts Day 9: Back to Python for Mapping

Because of the large size of the resulting JSON file, I don't think I'm going to be able to use the leafletjs javascript module to create a map.  I think that Python is my best bet.  

Mapping in Python might not be so bad anyway.  Matplotlib has a bunch of tools available, particularly modules Basemap and pyplot.  

This page has a lot of great examples that I can use to better understand what's going on.  The trick, though, is getting my shapefile read into it...

Looking at the documentation, I have to figure out some way of dealing with projections again.  My voter district shapefile isn't the standard lat/long values. 
Pyproj has a great method where you can set the inverse property and spit out lat long coordinates from the weird NC coordinates.  So, I'll create a big list of all the coordinates of all the shapes after running them through this method.

import shapefile
from pyproj import Proj

vote = shapefile.Reader('ncsbe\\Precincts.shp') #creates an instance that has the lists of data we want.

shapes = vote.shapes() #create lists of coordinates making up the shape for each voting district.
#for each shape in the shapefile...
for x in range(0,len(shapes)):
    #for each point in the shape...
    for y in range(0,len(shapes[x].points)):

        lon, lat = nc(shapes[x].points[y][0], shapes[x].points[y][1], inverse=True) #convert the shape file points into traditional lat/long coords.

and that's going to give me some data that I can plot according to the documentation in matplotlib (I think).

Voting Districts Day 8 Changing Projections

Well, I can't map my json file because the coordinates are not the traditional latitude/longitude values, but numbers in the 800Ks.  Here's an example of a coordinate:  [1884558.6496061385, 851226.325625971]

Why's this?  The shapefile I downloaded is using a different, less popular, projection to describe where the shapes go in relation to the earth.  So, I need to change the projection of the voting district shapefile to have traditional latitude longitude coordinates.

There's a python package called pyproj that'll hopefully help me.  I was able to install the binary version with pip.  

What projection do I use?  Well, to figure this out, I looked at the prj file that came with the voting district's shapefile and it said:  NAD_1983_StatePlane_North_Carolina_FIPS_3200_Feet.

Thanks to spatialreference.org, I was able to get the projection the pyproj needed here (clicked the Proj4 link).  I plugged the value in the Proj4 link into the right parameter using pyproj's Proj class:

from pyproj import Proj

myProj = Proj("+proj=lcc +lat_1=34.33333333333334 +lat_2=36.16666666666666 +lat_0=33.75 +lon_0=-79 +x_0=609601.2199999999 +y_0=0 +ellps=GRS80 +datum=NAD83 +to_meter=0.3048006096012192 +no_defs")

and then test it out using a coordinate from the voting district shapefile using the corresponding method for Proj and the myProj object I created.

 lat = my_json['features'][1000]['geometry']['coordinates'][0][0][1]

 lon = my_json['features'][1000]['geometry']['coordinates'][0][0][0]

 longi, latti = myProj(lon, lat, inverse=True)

From that, the long and lat values looked to be in the right ballpark (NOT!).  

It turns out that the projection requires another parameter according to this forum post due to the measurements being different (feet vs meters).  I need to set preserve_units to True.    

nc = Proj("+proj=lcc +lat_1=34.33333333333334 +lat_2=36.16666666666666 +lat_0=33.75 +lon_0=-79 +x_0=609601.2192024384 +y_0=0 +datum=NAD83 +to_meter=0.3048006096012192 +no_defs ", preserve_units=True)

Now when I use Proj with the new parameters, I get something much more realistic:  (-81.23183299105735, 36.13090600029268).

So, now I've got to plug the use of pyproj into my overall script, figure out how to convert all of my shapes' coordinates so I can finally create a map!

Voting Districts Day 7: Post a VERY basic map

So, I'm starting to figure out how to create a map with a leaflet package in javascript.

Leafletjs' tutorial has been very helpful so far in figuring out how to make a map appear and adding a layer to it. I can then add a layer to the map by pasting the coordinates via a json file. I was able to find a shape representing NC on Github from a gentleman named Johan.

I couldn't get it to work without just cutting and pasting the JSON text into the javascript portion of the file, but I'll work on that later. However, I was able to get a shape of NC to appear down below.

You need to reference the leaflet js in the head of the html.

Note: I stuck some hyphens at the beginning of the tag so it would show on the page below.

<--link href="docs/images/favicon.ico" rel="shortcut icon" type="image/x-icon">
<--link href="https://unpkg.com/leaflet@1.0.3/dist/leaflet.css" rel="stylesheet">
<--script src="https://unpkg.com/leaflet@1.0.3/dist/leaflet.js"><--/script>

You can then add the map and then the layer once you create it as a variable.

<--div id="mapid" style="height: 400px; width: 800px;">
<--/div>
<--script>
***Create Map***
var mymap = L.map('mapid').setView([35.505, -80.09], 7);
***Create layer variable***
var nc = ***JSON goes here***
***Add layer to map***
L.geoJSON(nc).addTo(mymap);
<--/script>

NOTE: To get the map to appear in Blogger, you have to use the HTML button, which is next to the compose button. I couldn't figure out to have my cake and eat it too. So, I typed up EVERYTHING manually, including the breaks and such, in html. Thank heavens for the preview button!

And that was the simplest way I could find to do this. This javascript is very new and foreign to me.

Next up, we are going to figure out how to create a json with our voting data (back to Python!) and shapefiles and how to read that in javascript as an external file.

Random Maps of Wake County Demographics

I've been playing around a lot with mapping data recently.  Here are some maps that represent voters in Wake County, North Carolina.  Each map represents a sample of 35000 voters collected from the North Carolina State Board of Elections in October 2013.  Each dot is a single voter and their residential location.

My data was collected from the following sources:

Voter registration information from the NC Board of Elections:  ftp://www.app.sboe.state.nc.us/

Mapping shapefiles from Wake county:  http://www.wakegov.com/gis/services/pages/data.aspx

and NC DENR:  http://portal.ncdenr.org/web/lr/geodetic/official-boundaries

Geocoding the addresses was done by Texas A&M's Geoservices:  http://geoservices.tamu.edu/

A Map of Registered Republicans and Democrats in Wake County

Here's my R code.  I utilized the R rgdal package for creating the maps.

This assumes that you've already got your voter data loaded into R.

roads <- readOGR("c:\\data\\poly\\wake_streets\\streets.shp","streets")

roadmap <- spTransform(roads, CRS("+proj=longlat +datum=WGS84"))

county <- readOGR("C:\\data\\poly\\nc_counties\\NC_State_County_Boundary_NAD83HARN.shp",'NC_State_County_Boundary_NAD83HARN')

countymap <- spTransform(county, CRS("+proj=longlat +datum=WGS84"))

plot(roadmap[roadmap$CLASSNAME == 'INT',],col='black',border='black', lwd=.5,axes=F,xlim=c(-79,-78.2),ylim=c(35.5,36.1))

plot(countymap[countymap$County == 'Wake',], add=T)

plot(roadmap[roadmap$CLASSNAME == 'USHWY',],col='black',border='black', lwd=.5, add=T)

points(vtx[vtx$party == 'REP','lng'],vtx[vtx$party == 'REP','lat'],col = rgb(255,0,0,50,maxColorValue=255),cex=.2,pch=20)

points(vtx[vtx$party == 'DEM','lng'],vtx[vtx$party == 'DEM','lat'],col = rgb(0,0,255,50,maxColorValue=255),cex=.2,pch=20)

title("Registered Wake County Republican or Democrats \n (sample of 35,000) - Oct 2013")

My data was collected from the following sources:

Voter registration information from the NC Board of Elections:  ftp://www.app.sboe.state.nc.us/

Mapping shapefiles from Wake county:  http://www.wakegov.com/gis/services/pages/data.aspx

and NC DENR:  http://portal.ncdenr.org/web/lr/geodetic/official-boundaries

Geocoding the addresses was done by Texas A&M's Geoservices:  http://geoservices.tamu.edu/

Mapping Shapefiles from the State of North Carolina

I've been looking forever on how to use shapefiles originating from North Carolina for making maps.  Normally, I'll get crazy latitude and longitude coordinates if I plot the shapefiles using the default parameters through R's rgdal package.

Tonight, I stumbled upon this StackOverflow post, which shows that TWO functions are needed in order to fully utilize NC shapefiles: readOGR and spTransform.

I've always just used readOGR, which works fine for shapefile originating from other places, like the US Census bureau.  I've always had problems with files from my state, however, until I came across the aforementioned StackOverflow post.

Here's what I put into R to make a simple map of all roads in Wake County:

roads <- readOGR("c:\\data\\poly\\wake_streets\\streets.shp","streets", p4s = CRS("+proj=lcc +lat_1=34.33333333333334 +lat_2=36.16666666666666 +lat_0=33.75 +lon_0=-79 +x_0=609601.2199999997 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs"))
roadmap <- spTransform(roads, CRS("+proj=longlat +datum=WGS84"))
plot(roadmap,axes=T)
title("Roads in Wake County, North Carolina")

This'll produce a map like so:

If you want to try this at home, you can get a bunch of shapefiles from North Carolina's State Board of Elections FTP site that seem to work with this code.

Now that this mystery has finally been uncovered, my mapping options are much improved!

Getting Geocodes through R and Google's Web Service

Part of my new job as a Data Integration Analyst is learning how to study and manipulate data.  So far, I've really enjoyed this new challenge and I love having the opportunity to learn something new.  

I learned pretty quickly that R is a pretty popular programming language within the realm of data and analytics.  By itself, R can perform some complex data analysis. However, packages provided by other R enthusiasts can be loaded into the R interface to make it more powerful.  I've spent the last few months getting more familiar with the language and additional packages and learning to appreciate it.  Although I still have a lot to learn, I can already see that R can do a lot of really cool stuff.

One aspect of analytics that I've been particularly fascinated with involves analyzing data through geography.  R has a lot of packages that make this pretty straightforward.  The ones I've seen so far are great, but, in order to map a specific place, you need geocoordinates (latitude and longitude points).  Providing just an address to R and one of these mapping packages won't do.

I really want to map some data regarding voters in my home county, Wake county, North Carolina.  I think I figured out how to do it.

Google provides a free web service that allows you to collect geocoordinates for any address. All you have to do is provide Google with a residential address through a URL.

R has a function that allows you to collect data through the web.  It's as easy as this:

getweb <- url('http://maps.googleapis.com/maps/api/geocode/xml?address=1600+Pennsylvania+Avenue,+20500&sensor=true')
getaddress <- readLines(getweb)
close(getweb)

I've just requested the geocoordinates of the White House, placed the results in another object, then closed the connection with Google.

Google returns the data in an XML string, which is now in my 'getaddress' object.  Google can also return JSON, but R has a package that can interpret XML for you. Once you install the package, you can collect the coordinates from the XML like so:

lng <- xmlValue(getNodeSet(xmlParse(y),'//result//geometry//location//lng')[[1]])
lat <- xmlValue(getNodeSet(xmlParse(y),'//result//geometry//location//lat')[[1]])

You now have coordinates!  Using one of R's available mapping packages, you can plot it like so.  

This simple map was created using one of the easier of R's maps packages to create a map. Here's the process:

map('usa',bg='lightblue',col='tan',fill=T)
points(lng,lat,pch='*',cex=10,col='red')

This is really just a glimpse into the world of mapping through R.  There's a ton of resources out there that allow you to map all sorts of regions, locations, boundaries, and landmarks. 

The possibilities are endless.

NOTE:  Google is very generous to provide geocoordinates for free.  However, they do limit the number of daily queries for each person to 2500.