Linux, Performance, PostGIS, PostgreSQL

Speeding Up Spatial Queries: GiST Indexes and EXPLAIN ANALYZE with PostGIS on PostgreSQL 18

by jb

In the previous articles we built a production-ready PostGIS environment on PostgreSQL 18 and created an interactive Leaflet map that displays parks from OpenStreetMap data. While the application works, spatial queries can quickly become slow as the amount of data grows.

In this article we will investigate why certain spatial queries are slow, create our first GiST index, compare execution plans with EXPLAIN ANALYZE, and measure the real impact on map loading times.

The Problem: Slow Spatial Queries

Let’s start by looking at the query used in our Leaflet application:

SELECT
  json_build_object(
    'type', 'FeatureCollection',
    'features', json_agg(feature)
  )
FROM (
  SELECT
    json_build_object(
      'type', 'Feature',
      'geometry', ST_AsGeoJSON(ST_Transform(way, 4326))::json,
      'properties', json_build_object('name', name, 'type', leisure)
    ) AS feature
  FROM planet_osm_polygon
  WHERE leisure = 'park'
    AND ST_DWithin(
          ST_Transform(way, 4326)::geography,
          ST_SetSRID(ST_MakePoint(11.576124, 48.137154), 4326)::geography,
          10000
        )
) t;

This query finds all parks within a 10 km radius of Munich. Without proper indexing, PostgreSQL has to examine a large portion of the table.

Analyzing the Query Without an Index

Run the following statement to see the current execution plan:

EXPLAIN (ANALYZE, BUFFERS, FORMAT TEXT)
SELECT count(*)
FROM planet_osm_polygon
WHERE leisure = 'park'
  AND ST_DWithin(
        ST_Transform(way, 4326)::geography,
        ST_SetSRID(ST_MakePoint(11.576124, 48.137154), 4326)::geography,
        10000
      );
Sample Output (click to expand): 
gisdemo=# EXPLAIN (ANALYZE, BUFFERS, FORMAT TEXT)
SELECT count(*)
FROM planet_osm_polygon
WHERE leisure = 'park'
  AND ST_DWithin(
        ST_Transform(way, 4326)::geography,
        ST_SetSRID(ST_MakePoint(11.576124, 48.137154), 4326)::geography,
        10000
      );
                                                                                               QUERY PLAN

------------------------------------------------------------------------------------------------------------------------------------------------------------
---------------------------------------------
 Aggregate  (cost=82703580.95..82703580.96 rows=1 width=8) (actual time=817.697..823.312 rows=1.00 loops=1)
   Buffers: shared hit=1483 read=252556
   ->  Gather  (cost=1000.00..82703580.95 rows=1 width=0) (actual time=434.015..823.182 rows=987.00 loops=1)
         Workers Planned: 2
         Workers Launched: 2
         Buffers: shared hit=1483 read=252556
         ->  Parallel Seq Scan on planet_osm_polygon  (cost=0.00..82702580.85 rows=1 width=0) (actual time=428.027..804.879 rows=329.00 loops=3)
               Filter: ((leisure = 'park'::text) AND st_dwithin((st_transform(way, 4326))::geography, '0101000020E6100000A4E194B9F9262740FF4124438E114840'::
geography, '10000'::double precision, true))
               Rows Removed by Filter: 2636472
               Buffers: shared hit=1483 read=252556
 Planning:
   Buffers: shared hit=3 read=3
 Planning Time: 0.689 ms
 Execution Time: 823.348 ms
(14 rows)

gisdemo=#

You will typically see a sequential scan (Seq Scan) on planet_osm_polygon. On larger OpenStreetMap extracts this can take several seconds.

Creating a (functional) GiST Index

The most important index for spatial queries in PostGIS is a GiST index on the geometry column. In this case we use a functional GiST index:

CREATE INDEX idx_planet_osm_polygon_geog
ON planet_osm_polygon
USING GIST ((ST_Transform(way, 4326)::geography));

After creating the index, update the table statistics:

ANALYZE planet_osm_polygon;

Comparing Execution Plans

Now run the same EXPLAIN ANALYZE statement again. You should see a significant change:

  • Instead of a sequential scan, PostgreSQL uses an Index Scan using the new GiST index.
  • The number of rows examined drops dramatically.
  • Execution time is usually reduced by a factor of 10–100, depending on the size of your dataset.
Sample Output (click to expand): 
gisdemo=# EXPLAIN (ANALYZE, BUFFERS, FORMAT TEXT)
SELECT count(*)
FROM planet_osm_polygon
WHERE leisure = 'park'
  AND ST_DWithin(
        ST_Transform(way, 4326)::geography,
        ST_SetSRID(ST_MakePoint(11.576124, 48.137154), 4326)::geography,
        10000
      );
                                                                                            QUERY PLAN

------------------------------------------------------------------------------------------------------------------------------------------------------------
---------------------------------------
 Aggregate  (cost=8238451.19..8238451.20 rows=1 width=8) (actual time=167.990..167.991 rows=1.00 loops=1)
   Buffers: shared hit=82358
   ->  Index Scan using idx_planet_osm_polygon_geog on planet_osm_polygon  (cost=0.54..8238451.19 rows=1 width=0) (actual time=0.285..167.586 rows=987.00 lo
ops=1)
         Index Cond: ((st_transform(way, 4326))::geography && _st_expand('0101000020E6100000A4E194B9F9262740FF4124438E114840'::geography, '10000'::double pr
ecision))
         Filter: ((leisure = 'park'::text) AND st_dwithin((st_transform(way, 4326))::geography, '0101000020E6100000A4E194B9F9262740FF4124438E114840'::geogra
phy, '10000'::double precision, true))
         Rows Removed by Filter: 298751
         Index Searches: 1
         Buffers: shared hit=82358
 Planning:
   Buffers: shared hit=5 dirtied=1
 Planning Time: 1.059 ms
 Execution Time: 168.036 ms
(12 rows)

gisdemo=#

Measuring the Impact on the Leaflet Map

After creating the index, reload the Leaflet application in your browser. The parks should appear noticeably faster.

You can also measure the response time of the PHP endpoint directly:

time curl -s http://lin6.fritz.box/api/parks.php > /dev/null

Compare the results before and after the index creation.

Sample Output (click to expand): 
# without the index:
[root@lin6 ~]# time curl -s http://lin6.fritz.box/api/parks.php > /dev/null

real    0m0.855s
user    0m0.001s
sys     0m0.009s
[root@lin6 ~]#

# with the index
[root@lin6 ~]# time curl -s http://lin6.fritz.box/api/parks.php > /dev/null

real    0m0.239s
user    0m0.004s
sys     0m0.004s
[root@lin6 ~]#

Creating an even better index

If we mostly query parks we can also create an index limited to parks:

drop index drop index idx_planet_osm_polygon_geog;
CREATE INDEX idx_planet_osm_polygon_geog_parks
ON planet_osm_polygon
USING GIST ((ST_Transform(way, 4326)::geography))
WHERE leisure = 'park';
ANALYZE planet_osm_polygon;

This gives us an even more efficient plan and runtimes:

gisdemo=# CREATE INDEX idx_planet_osm_polygon_geog_parks
ON planet_osm_polygon
USING GIST ((ST_Transform(way, 4326)::geography))
WHERE leisure = 'park';
CREATE INDEX
gisdemo=# ANALYZE planet_osm_polygon;
ANALYZE
gisdemo=# EXPLAIN (ANALYZE, BUFFERS, FORMAT TEXT)
SELECT count(*)
FROM planet_osm_polygon
WHERE leisure = 'park'
  AND ST_DWithin(
        ST_Transform(way, 4326)::geography,
        ST_SetSRID(ST_MakePoint(11.576124, 48.137154), 4326)::geography,
        10000
      );
                                                                              QUERY PLAN

------------------------------------------------------------------------------------------------------------------------------------------------------------
----------
 Aggregate  (cost=33.30..33.30 rows=1 width=8) (actual time=24.092..24.093 rows=1.00 loops=1)
   Buffers: shared hit=972
   ->  Index Scan using idx_planet_osm_polygon_geog_parks on planet_osm_polygon  (cost=0.27..33.29 rows=1 width=0) (actual time=0.089..23.994 rows=987.00 lo
ops=1)
         Index Cond: ((st_transform(way, 4326))::geography && _st_expand('0101000020E6100000A4E194B9F9262740FF4124438E114840'::geography, '10000'::double pr
ecision))
         Filter: st_dwithin((st_transform(way, 4326))::geography, '0101000020E6100000A4E194B9F9262740FF4124438E114840'::geography, '10000'::double precision
, true)
         Rows Removed by Filter: 99
         Index Searches: 1
         Buffers: shared hit=972
 Planning:
   Buffers: shared hit=26
 Planning Time: 0.238 ms
 Execution Time: 24.119 ms
(12 rows)

gisdemo=#

[root@lin6 ~]# time curl -s http://lin6.fritz.box/api/parks.php > /dev/null

real    0m0.084s
user    0m0.005s
sys     0m0.004s
[root@lin6 ~]#

Best Practices

  • Always run ANALYZE after creating spatial indexes.
  • Prefer GiST indexes for most PostGIS use cases (SP-GiST is useful only in specific scenarios).
  • Use partial indexes when you repeatedly filter on the same attribute (e.g. leisure = ‘park’).
  • Monitor index usage with pg_stat_user_indexes.
  • Rebuild bloated indexes periodically with REINDEX INDEX CONCURRENTLY.

Conclusion

Adding a GiST index is one of the highest-impact optimizations you can make in a PostGIS environment. In our example, a single index transformed a slow sequential scan into a fast index scan and significantly improved the responsiveness of the Leaflet map.

In the next article we will dig deeper into spatial query patterns and explore the differences between ST_DWithin, ST_Intersects, and bounding-box filters — and when to use each of them.

0