Overlaying Shapefile datapoints on Density Map

Question

I am new to shapefiles and mapping in python so I was hoping to get some help with overlaying data points from a shapefile on a density map.

To be honest, I am a beginner with mapping and reading in shapefiles so what I have so far not much.

I have started off using pyshp but if there are better packages out there to do this then I would love any feedback.

The following code is to create the base map of the LA area:

def get_base_map(rides_clean):
    return folium.Map(locations=[rides_clean.start_lat.mean(),                                     
                                 rides_clean.start_lon.mean()],
                      zoom_start = 20, tiles = 'cartodbpositron')

The following code is to create the density/heat map:

from folium import plugins
stationArr = rides_clean[['start_lat', 'start_lon']][:40000].as_matrix()
get_base_map(rides_clean).add_child(plugins.HeatMap(stationArr, 
                                    radius=40, max_val=300))

The following code is the same heat map but with route lines added:

(draw_route_lines(get_base_map(rides_clean), 
        routedf_vol)).add_child(plugins.HeatMap(stationArr, radius=40, 
        max_val=300))

I want to see data points from the shapefile shown as markers on top of the density plot.

Answer 1

It is possible to do this with pyshp. I've only ever used Matplotlib to plot shapefile points on a map, but this method will create two arrays which will be the x and y coordinates of each point you'd like to plot. The first snippet is used if you have multiple shapes in your shapefile, while the second can be used if you only have one shape.

import shapefile
import numpy as np
sf = shapefile.Reader('/path/to/shapefile')
point_list = []

for shape in sf:
    temp = shape.points()
    point_list.append(temp)

point_list = np.array(point_list)
x = point_list[:,0]
y = point_list[:,1]

And for a shapefile with only a single shape:

import shapefile
import numpy as np
sf = shapefile.Reader('/path/to/shapefile')

point_list = np.array(sf.shape(0).points)
x = point_list[:,0]
y = point_list[:,1]

You can tell how many shapes are in your shapefile using sf.shapes() and it will print a list detailing all the different shapes. From your question it appeared you were wanting to plot it as points on the marker rather than lines, sorry if this is not the case.