[GIS] Computing radius of gyration in PostGIS

centroidspostgispostgresql

The radius of gyration is a measure of the standard deviation of distances between points (and their center of mass), i.e., how frequently and how far points are moving.
A formal definition is:

radius of gyration equation

I have a large table of localization data (x-y coordinates, not a geometry field – might be transformed if faster?). What is the best way to compute the radius of gyration with PostGreSQL / PostGIS?

For the center of mass, ST_Centroid(ST_Union(ST_MakePoint(x,y)))?
Then ST_Distance? How to aggregate all this? ST_Distance(MakePoint(x,y),ST_Centroid(ST_Union(ST_MakePoint(x,y)))) does not work (ERROR: column "x" must appear in the GROUP BY clause or be used in an aggregate function).

Sorry, maybe a stupid question, I'm not so familiar with PostGreSQL.

Best Answer

I think it is better to store your points as geometries because you can use a spatial index for speeding up queries.

The query:

WITH points AS
    (SELECT
        ST_SetSRID(ST_MakePoint(x,y), 4326) AS geom
    FROM your_xy_table
    ),
centroid AS
    (SELECT
        ST_Centroid(ST_Union(points.geom)) AS geom
    FROM points),
multiobject AS
    (SELECT
        geom,
        generate_series(1,600) AS n 
    FROM points),
objects AS
    (SELECT
        n, 
        ST_GeometryN(geom, n) AS geometries 
    FROM multiobject 
        WHERE n <= ST_NumGeometries(geom)),
distance AS
    (SELECT DISTINCT ON
        (ST_Distance(geometries, centroid.geom))
        ST_Distance(geometries, centroid.geom) as distance
    FROM objects, centroid)

    SELECT
        stddev(distance) AS gyration
    FROM distance

First subquery: creates the points

Second subquery:calculates the centroid of your points

Third subquery: generates a series for your point objects, you have to change the value, if your data is bigger than 1000 objects

Fourth subquery: selects the single geometries for the distance calculation

Fifth subquery: calculates the distances

The last step is the calculation of the standard deviation of the distance values.

UPDATE

A second approach is this shorter query. The first subquery creates the points with the xy data in your table (but it is better to store the points as geometries, like mentioned above). You have to put in the EPSG number of your projection (e.g. 4326) in the first query. Change it when you work with a different projection.

WITH points AS
    (SELECT
        ST_SetSRID(ST_MakePoint(x,y), 4326) AS geom
    FROM your_xy_table),
centroid AS
    (SELECT
        ST_Centroid(ST_Union(points.geom)) AS geom
    FROM points),
objects AS
    (SELECT
        (ST_Dump(points.geom)).geom AS geometries
    FROM points),
distance AS
    (SELECT DISTINCT ON
        (ST_Distance(geometries, centroid.geom))
        ST_Distance(geometries, centroid.geom) AS distance
    FROM objects, centroid)

SELECT stddev(distance) FROM distance

With the first query maybe you can make other detailed calculations because you handle with the single objects with their numbers (n).

Related Question