The most circular lakes in the US

Created in November 21, 2024

2024

At CUGOS this month we decided to do the 30-day map challenge collaboratively. Everyone picked a map to do and we all came together to show off what we came up with. My task was “only circular shapes”. After some thinking, I had an idea: find the most circular lake in each US state.

Lake geometry is pretty cool. And wouldn’t you know, the LAGOS team has curated a fantastic dataset of lake polygons for the conterminous United States. I wanted to give myself an additional challenge for this map: don’t use more than 1 GB of RAM. Let us begin!

Calculating circle-ness

There are a few methods out there for calculating polygon roundness, but for this map I want a direct measure of how circular a lake is. To do so, we construct a circle with the same area as a basin, center it on the lake centroid, then calculate the intersection of the lake and the circle. The ratio of the intersection area to lake area is our measure of circle-ness. Fortunately the shapely library makes this easy to do in Python.

import numpy as np
import shapely

def circle_of_same_area(shape: shapely.geometry.MultiPolygon):
    radius = np.sqrt(shape.area / np.pi)
    circle = shapely.centroid(shape).buffer(radius)
    return circle

def overlap_with_circle(shape: shapely.geometry.MultiPolygon):
    circle = circle_of_same_area(shape)
    
    return circle.intersection(shape).area / circle.area

For example, here is the result for Trout Pond in Aroostook county, Maine (LAGOS lake ID 55038).

Reading in chunks

The full uncompressed LAGOS geopackage is a little over 18 GB on disk. So, we have to read the data in chunks to stay within our 1 GB limit. The function gpd.read_file has a parameter rows that lets us read a subset of a file. We select a small chunk size and iterate over the whole geopackage, only retaining the overlap index we need.

from osgeo import ogr
ogr.UseExceptions()
import pandas as pd
import geopandas as gpd

dataset = "gis_locus_v1.0.gpkg"

# Get number of features without reading
# entire layer into memory
ds = ogr.Open(dataset)
n_feats = ds.GetLayer(0).GetFeatureCount()
ds = None # close connection

# Arbitrarily pick 20 chunks
chunk_size = n_feats // 20

# Generator to create slice objects we pass to the
# rows parameter.
def enumerate_chunks(chunk_size, stop):
    upto = 0
    while upto < stop:
        # Add one because slice object are right-exclusive
        next = min(stop+1, upto+chunk_size)
        yield upto, next
        upto = next

def read_chunk_get_summary(chunk_tuple):
    this_chunk = gpd.read_file(dataset, rows=slice(*chunk_tuple), layer="lake")
    this_chunk["overlap"] = this_chunk.geometry.apply(overlap_with_circle)
    return this_chunk[["lagoslakeid", "Shape_Area", "overlap"]]

summary_df = pd.concat((read_chunk_get_summary(ch) for ch in enumerate_chunks(chunk_size, n_feats)))

Winners by state

Now we attach a state attribute to our summary data frame, pick the one with the highest circle score, and extract only those rows from the original dataset.

Since this isn’t meant to be a code post, you can find all the details on the github repo linked at the bottom of the page :P

Now, let’s see some circles!

Yeah those are some circles all right. If you check out where all these lakes are on the below folium map, you’ll see that they are a) pretty small and b) usually made by people. The size part makes sense - lake shorelines become more complex as they grow larger. Since a circle has a very undeveloped shoreline, it follows that the most circular lakes will all be quite small.

{
"type": "FeatureCollection",
"name": "circles",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:OGC:1.3:CRS84" } },
"features": [
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -88.052978582739385, 30.744749829632678 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -111.663277706234354, 35.238665123551179 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -92.785771954407579, 35.019589690947797 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -120.655089268233297, 35.245524043081055 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -105.90595761034885, 37.485828100061063 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -72.013872216950972, 41.843845432139254 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -75.521415560417395, 39.756793179254636 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -77.110913065066455, 38.939000655761816 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -81.612496673935752, 28.37723727603387 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -84.709442449156626, 31.725740802429907 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -115.304509688456065, 47.108110571919021 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -89.481225611173443, 37.722487375452324 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -86.984682116194534, 38.37355829374944 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -90.435500845389555, 41.772685315106763 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -101.860077413811268, 39.257770882062687 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -86.196682979739364, 37.954913061408263 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -90.497379178653688, 30.705945741329355 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -67.886453411928969, 44.771150683726461 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -75.905774601169398, 38.351596358749184 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.482099657116166, 42.432173600041914 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -84.412192967715711, 45.183397569090417 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -94.394959110077949, 46.939573661473581 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -88.694574135425142, 30.735727505618176 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -89.292873324205715, 36.640427869758668 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -109.301019092621587, 45.926509366933644 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -102.5760036365091, 42.031061043304916 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -114.744357903867268, 37.979576817629955 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.457117895505718, 43.074919451899696 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -74.346084463562207, 40.370738175093663 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -104.346681768150745, 33.233270528122652 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -78.840111762656448, 42.90473122359348 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -80.197902923375793, 35.399128704137361 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -98.971676078885054, 46.638768322706262 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -82.776494214141778, 41.214599568072231 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -96.439613443518127, 34.085738126677292 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -120.162529037423738, 43.429186103422488 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -79.009994703121436, 41.830627900157886 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.611648184126125, 41.488113838888452 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -82.307613035968487, 34.88292192955906 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -98.099192091449041, 45.43429891672276 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -84.759007180771249, 35.29843411546554 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -95.326706014552173, 28.940246275547963 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -113.776011815437002, 37.156371211650587 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -72.896477818937228, 44.72112013260935 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -79.026934347372034, 37.99903138626523 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -122.283142347024381, 47.856041468376489 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -80.779537671678668, 37.963985484904946 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -88.133127229017248, 45.510319466500029 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -107.257940206575483, 41.74360131567795 ] } }
]
}

##

Weighting by area

Let’s see if we can make something a little more fun than a bunch of really tiny circles. To incentivize large, but still circular lakes, define a “circle score” as the overlap index we calculated above multiplied by the logarithm of lake area. Now who are the winners in each state?

And to complete the thought, where are these super-circles?

{
"type": "FeatureCollection",
"name": "circle_scores",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:OGC:1.3:CRS84" } },
"features": [
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -88.004739745528113, 30.759234754445082 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -111.453746896531698, 34.94484660551543 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -93.99318311534509, 33.748982432916428 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -119.014775717982062, 38.007607356350647 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -104.073689968861046, 40.391892913057163 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -73.445813270207722, 41.955695332371327 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -75.256907060303917, 38.859841868426024 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -77.091666452105684, 38.911177647560912 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -80.796598628241767, 26.948641852366585 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -82.856043200324478, 34.392250217305282 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -111.40362587826526, 44.641506570950099 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -88.636952633928871, 41.254665862319314 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -87.757842529339655, 38.350708329049233 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -95.095651420219639, 43.477604371228054 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -95.802898694329414, 38.251502145897248 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -87.8108766841316, 37.680422764420825 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -90.125539190672612, 30.184095828662578 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -70.548586413322667, 43.856179971162149 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -75.305475583668354, 38.108931695321886 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -70.897358457275089, 41.800503603131773 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -84.721626930797058, 44.336977383221964 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -93.646218081124601, 46.242505670070642 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -90.701952190502922, 33.79178111391478 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -93.170947429723583, 39.614039278967823 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -111.681990121260171, 45.430958184414884 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -101.134951619656974, 41.106980982026741 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -118.717762519764591, 38.695450279099461 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.156641717828961, 43.596659875204061 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -74.825583500532474, 40.615890602604026 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -103.987698526209783, 32.320293090866066 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -78.997526696720129, 43.148329472964562 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -76.148752397780726, 35.508800286614893 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -98.831995351275424, 48.2473957561538 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -84.233298596431439, 40.737509610617899 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -97.59543109822431, 35.566480174814288 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -122.106387909063372, 42.942009844684016 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -79.009994703121436, 41.830627900157886 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.574828524283689, 41.439122668005126 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -80.053563651657569, 33.310488747393158 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -97.608670285630453, 44.493294614509836 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -89.694894097691133, 35.793584980549575 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -94.684297544096111, 29.802553237039607 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -112.601554645334971, 40.874586028652892 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -71.715722343569809, 44.955706919654929 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -76.469727848999668, 36.602823876265425 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -122.720457014987645, 45.679084357261111 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -81.763204491582044, 38.824872361183857 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -89.253647164448083, 42.964593040811209 ] } },
{ "type": "Feature", "properties": { }, "geometry": { "type": "Point", "coordinates": [ -108.601954773085694, 43.184804199884212 ] } }
]
}

##

Code

Available here! This was a lot of fun. I wonder what happens if you look for the most square lakes… or the most spider-webby lakes.