PostGIS Clustering with MAX, MIN and SUM Limit of Attribute

Question

I have point data with number of peoples living in a house.

wkt111,members
SRID=4326;POINT(8.40771254932301 52.8935196103868),158
SRID=4326;POINT(8.35609033091261 52.8948388844689),122
SRID=4326;POINT(8.43633740614695 52.9290469721589),144
SRID=4326;POINT(8.43637376596481 52.9282851939852),144
SRID=4326;POINT(8.35612742505646 52.8953039259067),122
SRID=4326;POINT(8.40722533514059 52.8941621951891),158
SRID=4326;POINT(8.43974835467466 52.9295371710363),106
SRID=4326;POINT(8.4118531995535 52.893081792927),158
SRID=4326;POINT(8.4361907632209 52.9299321213647),144
SRID=4326;POINT(8.43860541290836 52.9287882893147),122
SRID=4326;POINT(8.43930763983385 52.929938063062),106
SRID=4326;POINT(8.43640308495477 52.92761686129),144
SRID=4326;POINT(8.43761323022107 52.9294136825416),122
SRID=4326;POINT(8.35565716705095 52.8947943141241),122
SRID=4326;POINT(8.40756478326444 52.8928834796089),116
SRID=4326;POINT(8.41197915975963 52.8922304811318),158
SRID=4326;POINT(8.35575141592804 52.8952636611962),122
SRID=4326;POINT(8.41199100509697 52.8928796332762),158
SRID=4326;POINT(8.35700296469407 52.8946903270974),122
SRID=4326;POINT(8.4081051389381 52.8927128902747),158
SRID=4326;POINT(8.40759193225009 52.8930786154292),116
SRID=4326;POINT(8.44019422628302 52.9295509228059),106
SRID=4326;POINT(8.41188170626939 52.8928475154659),158
SRID=4326;POINT(8.3562215746048 52.8948080186652),122
SRID=4326;POINT(8.43857527294205 52.929880156525),106
SRID=4326;POINT(8.4077242541519 52.8941146973073),158
SRID=4326;POINT(8.40775467492524 52.8922823085847),116
SRID=4326;POINT(8.43642552740359 52.927022173622),144
SRID=4326;POINT(8.35657574880137 52.8951443234909),122
SRID=4326;POINT(8.40848705845392 52.8920713315788),158
SRID=4326;POINT(8.40776048054732 52.893703405996),158
SRID=4326;POINT(8.43636411318751 52.9284504349415),144
SRID=4326;POINT(8.40741653275869 52.8934191858778),158
SRID=4326;POINT(8.35617545454214 52.8950776214254),122
SRID=4326;POINT(8.40823763848918 52.8923667872611),158
SRID=4326;POINT(8.43639461217891 52.9278189867993),144
SRID=4326;POINT(8.43853189619463 52.9295519446821),106
SRID=4326;POINT(8.40849838008488 52.8918852920387),158
SRID=4326;POINT(8.4389839529149 52.9295506923135),106
SRID=4326;POINT(8.40799108038239 52.8928973200872),158
SRID=4326;POINT(8.4381127659551 52.92906558007),122
SRID=4326;POINT(8.41209469177733 52.8922374887863),158
SRID=4326;POINT(8.40816140613239 52.8925378483298),158
SRID=4326;POINT(8.43985303594237 52.9299617157145),106
SRID=4326;POINT(8.3556053817816 52.8950242116499),122
SRID=4326;POINT(8.40793700015185 52.8930687863224),158
SRID=4326;POINT(8.38921628504428 52.8925246492037),112
SRID=4326;POINT(8.39473176195843 52.891683086891),112
SRID=4326;POINT(8.40739093896853 52.8944262688579),158
SRID=4326;POINT(8.40832773783288 52.8922089171594),158
SRID=4326;POINT(8.43937823057937 52.9295592639623),106
SRID=4326;POINT(8.44026594300409 52.9298270972058),106
SRID=4326;POINT(8.40735356098727 52.8935752290805),158
SRID=4326;POINT(8.43940727028956 52.9299369614778),106
SRID=4326;POINT(8.43634853799564 52.9286309912991),144
SRID=4326;POINT(8.35706977648743 52.8949976908418),122
SRID=4326;POINT(8.40775811103838 52.8938854805683),158
SRID=4326;POINT(8.39161917424377 52.8921613436314),112
SRID=4326;POINT(8.40780947454812 52.8944071474669),158
SRID=4326;POINT(8.43702951344901 52.9298314634418),122
SRID=4326;POINT(8.4375248181906 52.9295802736219),122
SRID=4326;POINT(8.43631696721841 52.9295696636839),144
SRID=4326;POINT(8.43628059188016 52.929764549172),144
SRID=4326;POINT(8.40727075744272 52.8939971673637),158
SRID=4326;POINT(8.4382991417383 52.9297268266898),106
SRID=4326;POINT(8.43633560151942 52.9288448256907),144
SRID=4326;POINT(8.43895687449178 52.9298073175027),106
SRID=4326;POINT(8.39273216869533 52.8920984174851),112
SRID=4326;POINT(8.43971473681796 52.9299635281786),106
SRID=4326;POINT(8.43799059571948 52.9298995849497),106
SRID=4326;POINT(8.43723205021406 52.9296960635931),122
SRID=4326;POINT(8.4075824477187 52.8928286205),116
SRID=4326;POINT(8.39351114085647 52.8916949848888),112
SRID=4326;POINT(8.41177354519619 52.8925962489103),158
SRID=4326;POINT(8.43834504310234 52.9289302986786),122
SRID=4326;POINT(8.44067887618606 52.9295576528102),106
SRID=4326;POINT(8.4363750768028 52.9280731928726),144
SRID=4326;POINT(8.40801362256092 52.8920202025446),116
SRID=4326;POINT(8.43632109984682 52.9292300058607),144
SRID=4326;POINT(8.40785381832656 52.8925270899493),116
SRID=4326;POINT(8.41171402780864 52.893658373826),158
SRID=4326;POINT(8.39014657911192 52.8922392583656),112
SRID=4326;POINT(8.35652068045738 52.8947780314362),122
SRID=4326;POINT(8.40701005464346 52.8943587912332),158
SRID=4326;POINT(8.40735082328216 52.8937746314017),158
SRID=4326;POINT(8.43640705432114 52.9274074791734),144
SRID=4326;POINT(8.41183287734992 52.8931357408673),158
SRID=4326;POINT(8.35682712053401 52.8944970640479),122
SRID=4326;POINT(8.4363201037938 52.9294231689194),144
SRID=4326;POINT(8.44057979522546 52.9297512515431),106
SRID=4326;POINT(8.39100985742862 52.8922600252125),112
SRID=4326;POINT(8.4378599136052 52.9293060812469),122
SRID=4326;POINT(8.41171179591219 52.8935838674255),158
SRID=4326;POINT(8.44105081209008 52.9298855598048),106
SRID=4326;POINT(8.43643577760023 52.9267984844396),144
SRID=4326;POINT(8.4120753326162 52.8926832890201),158
SRID=4326;POINT(8.40762272122407 52.8927048328018),116
SRID=4326;POINT(8.4119304894701 52.8924288327541),158
SRID=4326;POINT(8.43638879140181 52.9272268727532),144

cluster within 2km.
cluster should have min 100 or max 1000 of peoples/members. if exceeded the amount create a new cluster.
sub-group should be a cluster if meet number 2. Note: there may be more than 1 sub-cluster in the cluster.

I have used the clusterDBscan, but not sure how to make a cluster with max and min sum of peoples.

@MuhammadImranSiddique if members are below 100, do you want to make them part of neighboring clusters or NULL? — wondim, Jun 23 '20 at 14:46

wondim · Accepted Answer · 2020-06-23 15:08:29Z

To make your spatial table valid I suggest you to add an id column. In addition, I have added a cluster column which will store the clusters. Below is the table create query.

CREATE TABLE houses (
    id integer NOT NULL DEFAULT nextval('houses_id_seq'::regclass),
     members integer,
    geom GEOMETRY(POINT, 4326),
    CONSTRAINT houses_pkey PRIMARY KEY (id)
)

I added id column and added consecutive integers (1,2,3...), renamed wkt111 to geom. I made sure the order of columns in the DB table and csv table are the same. Then I imported your data using:

COPY houses FROM '/path/to/your/csvfile.csv';

Using this query below, the points are clustered by distance using ST_ClusterDBSCAN. And the geographic coordinate system (4326) is transformed to local projected coordinate system (4839) that uses meters, which you can change to any other.

SELECT *, ST_ClusterDBSCAN(ST_Transform(geom,4839), eps := 2000, minpoints := 1) OVER() AS cl_id FROM houses ORDER BY cl_id

The loop basically uses the above result to further cluster the points based on the members count.

The maximum - 1000 members is enforced. The minimum 100 members can be enforced when you suggest where they should go.

The following anonymous function updates the cluster column with cluster numbers:

DO $$
  DECLARE
    arow record;
    counter int := 0;
    curr_cluster int := 0;
    previous_cluster int := 0;
  BEGIN
    FOR arow IN 
      SELECT *, ST_ClusterDBSCAN(ST_Transform(geom,4839), eps := 2000, minpoints := 1)
        OVER() AS cl_id FROM houses ORDER BY cl_id 
    LOOP
     --- It means this is a new cluster 
    IF  previous_cluster <> arow.cl_id THEN 
        -- Increment by one as it is a new cluster
        curr_cluster = curr_cluster + 1;  
        counter := arow.members; -- reset members counter
     --- looping in the same original cluster group so increment members normally
     ELSIF previous_cluster = arow.cl_id THEN 
        counter := counter + arow.members;
    
    END IF;
       previous_cluster := arow.cl_id;

      IF counter < 1001 THEN
    
        UPDATE houses SET cluster = curr_cluster WHERE id=arow.id;
    
      ELSIF counter > 1000 THEN
          counter := arow.members; -- reset members counter
          -- increment current clusterr members should be below 1000
          curr_cluster = curr_cluster + 1; 
         UPDATE houses SET cluster = curr_cluster WHERE id=arow.id;
       END IF;
     
     END LOOP;
 
  END;
$$;

Finally, you can style the points based on the cluster column by setting different colors for each cluster.

Thanks for the effort but results should be following. the cluster should have min 100 or max 1000 of peoples/members not a number of rows. — Muhammad Imran Siddique, Jun 24 '20 at 18:20
Of course, that is how I wrote the code. It basically adds the members and when the count reaches 1000, it creates the next cluster. If the count is below 100, what do you want to do to the record? Add it to the previous record or make it NULL? — wondim, Jun 24 '20 at 19:05

Muhammad Imran Siddique · Accepted Answer · 2020-07-20 06:01:13Z

I have tried some different methods this one gives me near results but not fully accurate.1

DROP FUNCTION IF EXISTS public.I_Cluster(int4, int4, int4, int4);
CREATE OR REPLACE FUNCTION "public"."i_cluster"("distance" int4, "household_min" int4,
"min_threshold" int4, "max_threshold" int4)

RETURNS TABLE("house" varchar, "cluster_id" int4) AS $BODY$ DECLARE
counter_cluster INTEGER := 1;
counter INTEGER := 0;
cluster_array CHARACTER VARYING [];
cluster_id INTEGER [];
cluster_next_loop_value CHARACTER VARYING;
total INTEGER := 0;
rec RECORD;
rec1 RECORD;
BEGIN
    FOR rec IN WITH res AS (
    SELECT
        A.haushalte,
        C.NAME,
        C.geom 
    FROM
        cluster_final A,
        house_polygons b,
        housecoordinates_final C 
    WHERE
        C."name" = b."coordname" 
        AND b."clustername" = A."name" 
        AND A.haushalte :: int4 > household_min
    ) SELECT A.*, st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) ),
    b.NAME name1, b.haushalte haushalte1 
FROM
    res A,
    res b 
WHERE
    st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) ) < distance 
ORDER BY
    A.NAME,
    st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) ),
    b.NAME
LOOP
IF counter = 0 THEN
        cluster_array := ARRAY_APPEND( cluster_array, rec.name1 );
        cluster_id := ARRAY_APPEND( cluster_id, counter_cluster );
        total := total + rec.haushalte1 :: INTEGER;
    
END IF;
IF rec.name1 = ANY ( cluster_array ) THEN
        ELSE
    IF total < max_threshold THEN
        
            cluster_id := ARRAY_APPEND( cluster_id, counter_cluster );
            cluster_array := ARRAY_APPEND( cluster_array, rec.name1 );
            total := total + rec.haushalte1 :: INTEGER;
            --RAISE NOTICE'Cluster Value %', rec.name1;
            --RAISE NOTICE'Total %', total;

    ELSIF total > min_threshold THEN
            total := 0;
            counter_cluster := counter_cluster + 1;
            cluster_next_loop_value := rec.name1;
            --RAISE NOTICE'Second Loop %', cluster_next_loop_value;
            --RAISE NOTICE'Total %', total;
            FOR rec1 IN WITH res AS (SELECT
                    A.haushalte, C.NAME, C.geom 
                FROM
                    cluster_final A, house_polygons b, housecoordinates_final C 
                WHERE
                    C."name" = b."coordname" 
                    AND b."clustername" = A."name" 
                    AND A.haushalte :: int4 > household_min 
                ) 
                SELECT A.*,
                st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) )::TEXT,
                b.NAME name1,
                b.haushalte haushalte1 
            FROM
                res A, res b 
            WHERE
                st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) ) < distance 
            ORDER BY A.NAME
                , st_distance ( st_transform ( A.geom, 3857 ), st_transform ( b.geom, 3857 ) )
                , b.name
        LOOP
            IF cluster_next_loop_value=rec1.name1 THEN
            IF rec1.name = ANY ( cluster_array ) THEN
                ELSE
                total := total + rec1.haushalte1::int8;

                IF total < max_threshold THEN
                      cluster_id := ARRAY_APPEND( cluster_id, counter_cluster );
                      cluster_array := ARRAY_APPEND( cluster_array, rec1.name1);
ELSIF total > max_threshold THEN
                    total := 0;
                    counter_cluster := counter_cluster + 1;
                    END IF;
                END IF;
            END IF;
        END LOOP;
ELSIF total <= min_threshold THEN
    total := 0;
    END IF;
    END IF;
    counter := counter + 1;
    
END LOOP;
RETURN query SELECT UNNEST( cluster_array ), UNNEST ( cluster_id );
END; $BODY$ LANGUAGE plpgsql VOLATILE COST 100 ROWS 1000;
select row_number() over() id, * from i_cluster(2000, 100, 100, 1000) as a, housecoordinates_final b where a.house=b.name;

Stack Exchange Network

current community

your communities

more stack exchange communities

PostGIS Clustering with MAX, MIN and SUM Limit of Attribute

2 Answers 2

Your Answer

Not the answer you're looking for? Browse other questions tagged postgis summarizing min or ask your own question.

Hot Network Questions

PostGIS Clustering with MAX, MIN and SUM Limit of Attribute

2 Answers 2

Your Answer

Sign up or log in

Post as a guest

Not the answer you're looking for? Browse other questions tagged postgis summarizing min or ask your own question.

Related

Hot Network Questions