DuckDB Geospatial Analytics: Turning Coordinates into Business Insights

-- Load the spatial extension (first time requires download)
LOAD spatial;
INSTALL spatial;

-- Verify installation
SELECT version();

LOAD spatial;

-- Create your store table
CREATE TABLE my_stores AS
SELECT * FROM (VALUES
    ('Beijing Sanlitun',   116.4510, 39.9320),
    ('Shanghai Nanjing Rd', 121.4730, 31.2300),
    ('Guangzhou Tianhe',   113.3250, 23.1300),
    ('Shenzhen Huaqiang',   114.0650, 22.5400),
    ('Chengdu Chunxi Rd',  104.0800, 30.6600),
    ('Hangzhou West Lake',  120.1500, 30.2800),
    ('Wuhan Jianghan Rd',   114.2800, 30.5800),
    ('Chongqing Jiefangbei', 106.5800, 29.5600)
) AS t(name, longitude, latitude);

-- Create competitor store table
CREATE TABLE competitor_stores AS
SELECT * FROM (VALUES
    ('Starbucks Sanlitun',  116.4530, 39.9340),
    ('Luckin Sanlitun',     116.4490, 39.9310),
    ('Starbucks Nanjing Rd', 121.4750, 31.2320),
    ('Luckin Nanjing Rd',    121.4710, 31.2280),
    ('Manner Tianhe',        113.3270, 23.1320),
    ('Luckin Huaqiang',      114.0670, 22.5420)
) AS t(name, longitude, latitude);

-- Convert coordinates to WKT geometry objects
ALTER TABLE my_stores ADD COLUMN geom GEOMETRY;
UPDATE my_stores SET geom = ST_SetSRID(ST_MakePoint(longitude, latitude), 4326);

ALTER TABLE competitor_stores ADD COLUMN geom GEOMETRY;
UPDATE competitor_stores SET geom = ST_SetSRID(ST_MakePoint(longitude, latitude), 4326);

-- Calculate distance between every pair of stores (in meters)
SELECT 
    m.name AS my_store,
    c.name AS competitor,
    ROUND(
        ST_Distance(
            ST_Transform(m.geom, 3857),  -- Web Mercator (metric CRS)
            ST_Transform(c.geom, 3857)
        )::INTEGER
    ) AS distance_meters,
    CASE 
        WHEN ST_Distance(
            ST_Transform(m.geom, 3857),
            ST_Transform(c.geom, 3857)
        ) < 500 THEN 'Direct Competition'
        WHEN ST_Distance(
            ST_Transform(m.geom, 3857),
            ST_Transform(c.geom, 3857)
        ) < 1000 THEN 'Potential Competition'
        ELSE 'Safe Distance'
    END AS competition_level
FROM my_stores m
CROSS JOIN competitor_stores c
WHERE ST_DWithin(
    ST_Transform(m.geom, 3857),
    ST_Transform(c.geom, 3857),
    1000  -- Within 1000 meters
)
ORDER BY distance_meters;

-- Create 500-meter service buffers for each store
CREATE TABLE store_buffers AS
SELECT 
    name,
    ST_Buffer(ST_Transform(geom, 3857), 500) AS buffer_geom
FROM my_stores;

-- Count how many competitors fall within each buffer
SELECT 
    s.name AS store,
    COUNT(c.name) AS nearby_competitors
FROM store_buffers s
LEFT JOIN competitor_stores c 
    ON ST_Contains(s.buffer_geom, ST_Transform(c.geom, 3857))
GROUP BY s.name
ORDER BY nearby_competitors DESC;

-- Given a new order location, find the 3 nearest stores
WITH order_location AS (
    SELECT ST_SetSRID(ST_MakePoint(116.4550, 39.9350), 4326)::GEOMETRY AS geom
)
SELECT 
    m.name,
    ROUND(ST_Distance(
        ST_Transform(m.geom, 3857),
        ST_Transform(o.geom, 3857)
    )) AS distance_m
FROM my_stores m, order_location o
ORDER BY distance_m
LIMIT 3;

-- Load city district GeoJSON
CREATE TABLE districts AS
SELECT 
    feature->>'name' AS district_name,
    feature->>'population' AS population,
    st_geomfromjson(feature->'geometry') AS geom
FROM read_geojson('https://example.com/china-districts.geojson');

-- Count stores per district and calculate coverage ratio
SELECT 
    d.district_name,
    d.population,
    COUNT(m.name) AS store_count,
    ROUND(d.population::DOUBLE / NULLIF(COUNT(m.name), 0), 0) AS per_store_population
FROM districts d
LEFT JOIN my_stores m ON ST_Contains(d.geom, m.geom)
GROUP BY d.district_name, d.population
ORDER BY per_store_population;