I am reading a shapefile as geopandas DataFrame and them using pandas subset method to select a region.

geodata = gpd.read_file(bayshp)



1   06001   Alameda (POLYGON ((6065941.393835935 2104148.464510527...
2   06013   Contra Costa    (POLYGON ((6143913.640835938 2209458.230510532...
3   06041   Marin   (POLYGON ((5879149.417835938 2203020.920510533...
4   06055   Napa    POLYGON ((6075700.362835937 2441916.530510533,...
5   06075   San Francisco   (POLYGON ((5990480.312835939 2123810.13351053,..

Subset the df:

# Subset based on the index
geosub = geodata.iloc[0:2]

I've got a function that accepts geopandas DataFrame and number of points to sample as arguments.

def sample_random_geo(df, n):

    # Randomly sample geolocation data from defined polygon 
    points = np.random.sample(df, n)

    return points

However, the np.random.sample or for that matter any numpy random sampling doesn't support geopandas object type.

I am wondering if there is a way to randomly sample geocoordinates from the spatial region.

  • 1
    What you can do is sample from the integer range (0, 1, .., n) and then use the sampled integers to index the dataframe. But note that you also have a .sample() method on a (Geo)DataFrame that basically does that. – joris Aug 29 '18 at 22:05
  • The dataframe is a record of counties in the bay area with the geometry. What I am interested in sampling from different geometries. So for instance, in the above dataframe randomly sampling geocoordinates from the Alameda county. – kms Aug 29 '18 at 22:17
  • What do you mean by "sampling geocoordinates"? Random vertices of Alameda county polygon or coordinates of random points within that polygon or something else? – Kadir Şahbaz Aug 29 '18 at 22:24
  • @Kadir latter. Random points within the Alameda county polygon. – kms Aug 29 '18 at 22:55

GeoPandas uses Shapely geometries. As far as a know, there is no a function which gets random points within a polygon. So, you must write any like below. Add this script to yours.

from shapely.geometry import Point
import geopandas as gpd
import random

def random_points_in_polygon(number, polygon):
    points = []
    min_x, min_y, max_x, max_y = polygon.bounds
    i= 0
    while i < number:
        point = Point(random.uniform(min_x, max_x), random.uniform(min_y, max_y))
        if polygon.contains(point):
            i += 1
    return points  # returns list of shapely point

geodata = gpd.read_file("path/to/your/file.shp")

# generate 5 points within the first county polygon in geodata
points = random_points_in_polygon(5, geodata.iloc[0].geometry)

# Coordinates of the first point
# print(points[0].x, points[0].y)

# print coordinates of all points
for i, point in enumerate(points):
    print("Point {}: ({},{})".format(str(i+1), point.x, point.y))

Reference: How to generate random coordinates in a multipolygon in python


Here's another way to do it:

import geopandas as gpd
import numpy as np

# load an example polygons geodataframe
gdf_polys = gpd.read_file(gpd.datasets.get_path('nybb'))

It looks like the following:

enter image description here

# find the bounds of your geodataframe
x_min, y_min, x_max, y_max = gdf_polys.total_bounds

# set sample size
n = 100
# generate random data within the bounds
x = np.random.uniform(x_min, x_max, n)
y = np.random.uniform(y_min, y_max, n)

# convert them to a points GeoSeries
gdf_points = gpd.GeoSeries(gpd.points_from_xy(x, y))
# only keep those points within polygons
gdf_points = gdf_points[gdf_points.within(gdf_polys.unary_union)]

Now you have:

enter image description here


Here's a solution that takes advantage of MultiPoint and MultiPolygon to avoid loops.

import numpy as np
import geopandas as gpd
import shapely.geometry

def sample_geoseries(geoseries, size, overestimate=2):
    polygon = geoseries.unary_union
    min_x, min_y, max_x, max_y = polygon.bounds
    ratio = polygon.area / polygon.envelope.area
    samples = np.random.uniform((min_x, min_y), (max_x, max_y), (int(size / ratio * overestimate), 2))
    multipoint = shapely.geometry.MultiPoint(samples)
    multipoint = multipoint.intersection(polygon)
    samples = np.array(multipoint)
    while samples.shape[0] < size:
        # emergency catch in case by bad luck we didn't get enough within the polygon
        samples = np.concatenate([samples, random_points_in_polygon(polygon, size, overestimate=overestimate)])
    return samples[np.random.choice(len(samples), size)]

geodata = gpd.read_file(bayshp)
points = sample_geoseries(geodata['geometry'])

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.