[Python+SQL] 도로 geometry 로 네트워크 분석하기

Task:

도로 geometry 를 바탕으로 Network Centrality 분석에 활용하기 위한

Adjacency Matrix 생성 및 결과 값을 기반으로 한 shp 파일을 만들자.

 

Open street map 과 NetworkX 를 바탕을 Boeing 아저씨의 지분이 100% 인 네트워크 분석도 물론 가능하다.

 

네트워크 분석은 도시 공간 분석에서 빼놓지 않고 등장하는 떡밥이다.

공간을 수치화하여 정량적으로 분석할 수 있다니, 군침이 줄줄 흐른다.

 

위와 같은 task를 위해 물론 Boeing 이 개발한 Open Street Map 과 NetworkX 를 활용할 수 있고,

충분한 결과물을 내주지만 해당 분석은

"도로의 결절점" 을 중심으로 포인트를 내어준다는 점에서 도로 기반의 분석을 하고픈 우리에게 썩 기분좋은 일은 아니다.

 

도로의 인접 여부를 바탕으로 Adjacency Matrix 를 구성하고,

Network centrality (betweenness 와 closensss) 를 계산하는 방안을 알아보자.

 

1. Adjacency Matrix (PostGRESQL)

우선, 도로 geometry 를 Postgres 로 이동하여 데이터베이스화 시켜준다.

데이터베이스로 이동하는 것은 (Ctrl + Shift + D)를 통해 금방 할 수 있거니와,

[PostGRES] 다른 소스의 공간 데이터 병합을 위한 text similarity 기반 전처리 :: 이언배 연구노트 (tistory.com)

여기 초장에 설명했으니 자세한 설명은 생략한다.

[PostGRES] 다른 소스의 공간 데이터 병합을 위한 text similarity 기반 전처리

왼쪽의 geometry shp 파일이 오른쪽의 공간 데이터베이스로 이동된다. SQL을 활용하면 더욱 빠르고 자유로운 분석을 할 수 있다.

 

Adjacency matrix 구성은 매우 복잡한 병렬 계산을 필요로 한다.

그리고 SQL은 그런 계산을 매우 빠르게 진행해준다.

Python 이었으면 iteration 을 몇 번을 돌렸겠지만, SQL은 이러한 연산에 특화되어있다.

 

CREATE TABLE buk_adjacency_matrix AS
SELECT
    a.id AS id1,
    b.id AS id2,
    CASE WHEN ST_Intersects(a.geom, b.geom) THEN 1 ELSE 0 END AS adjacent
FROM
    buk_street a
CROSS JOIN
    buk_street b;

 

10줄이 안되는 SQL 코드면 아래와 같은 adjacency matrix가 만들어진다.

도로 id 마다 다른 도로에 인접해있으면 1, 아니면 0을 준다.

 

이제 남은 일은, 네트워크 분석이다. 

SQL이 공간 심화 분석까지 해주지는 않으므로, python 의 힘을 빌리자.

 

2.Network Analysis

SQL을 Python 으로 불러오기 위해서는, sql 연동 모듈 - psycpg2 가 필요하다.

numpy 랑 pandas 는 대두분 필요하니 고민 말고 import 하도록 하자.

import pandas as pd
import psycopg2 as psql
import numpy as np

앞으로 psql 이라고 부를 psycpg2 은,

데이터베이스가 설치된 서버 컴퓨터 (내 경우에는 local 컴퓨터이므로, localhost 라는 이름을 쓴다) 

에 접속 (conn) 을 만들고, 

원하는 SQL 명령어를 텍스트 형태로 데이터베이스에 전달해주는 'SQL 빵셔틀' 이라고 생각하면 된다.

"내가 SQL 문을 텍스트로 줄테니까 가서 데이터베이스에 접속한 다음 빵사와" 라고 하면 갖다주는 느낌.

db_name = 'jongno'
conn = psql.connect(host = 'localhost', port = 5432, user = 'postgres', password = '****', database = db_name)

 

접속 개체 (connect) 는 한땀한땀 공들여 데이터베이스를 뒤적대다가 (cursor)

내가 준 SQL 문에 적합한 데이터가 있으면 그제서야 한 뭉터기를 가져다가 (Query) python 으로 돌아온다.

나는 매번 이해하기 귀찮아서 그냥 text와 connect 를 집어넣으면 Dataframe 을 갖다주는 함수를 만들었다.

db_name = 'jongno'
conn = psql.connect(host = 'localhost', port = 5432, user = 'postgres', password = 'user', database = db_name)

def fetch_table(query, connect):
    '''
    Fetch table from database
    Input: str, str / Output: Dataframe
    '''
    cur = connect.cursor()
    cur.execute(query)
    rows = cur.fetchall()
    df = pd.DataFrame(rows)
    cur.close()
    return df

 

 

이제 Query 를 시작한다.

 

아까 만든 adjacency matrix 를 불러와 데이터프레임으로 변환해야 한다.

이 때, 멍청하게도 psql은  column name 을 잊어먹으므로 이름을 다시 지정해주는 일을 잊지 말자.

indexing 은 삽질만이 답이다.

table_name = 'buk_adjacency_matrix'
query1 = f'SELECT distinct on (id1) id1 FROM {table_name}'
df_ids = fetch_table(query1, conn)
df_ids = df_ids.reset_index()
df_ids = df_ids.rename(columns = {'index': 'df_id', 0: 'rd_id'})

 

우리가 얻는 df_ids(복선) 는 아래와 같다.

Dataframe 상의 index 와 SHP 파일 내의 Road ID 가 다르고,

나중에 numpy 분석 결과를 붙여야 하므로 이 키는 꼭 기억해두자.

 

이제 numpy 변환을 시작한다.

 

위 형태의 adjacency matrix 는 매우 불편-하다.

SQL 세계에서는 편하지만, python 세계에서는 어색한 형태다.

 

(물론 더 좋은 방법이 세상 어디엔가 있겠지만, 내 머리의 한계로는)

모든 도로 id 를 순서대로 iteration 돌면서 한 줄씩 array 에 붙일 예정이다.

 

우선, adjacent matrix 의 첫번째 row가 될 array 를 만든다.

첫번째 road id 를 SQL에서 query 해오고, array 로 변환하면 될 것 같다.

query = f'SELECT * FROM {table_name} WHERE id1 = {df_ids.iloc[0, 1]}'
first_df = fetch_table(query,conn)
first_row = np.array([first_df[2].values.T])
adj_mat = first_row

 

이 경우에는 첫번째 road id 인 5 의 adjacency matrix 결과만을 query 해와서

transpose 하고, 가로 1열로 만든 상태.

 

그 다음 모든 road id 를 돌아가며 동일한 방식으로 query 해오고, array 로 변환하고, transpose 하고,

 

갖다가 붙이면 되겠다.

for k, ids in df_ids[1:].iterrows():
    df_id = ids['df_id']
    rd_id = ids['rd_id']
    query = f'SELECT * FROM {table_name} WHERE id1 = {rd_id}'
    df_temp = fetch_table(query, conn)
    arr_vertical = df_temp[2].values
    arr_temp = np.array([arr_vertical.T])
    adj_mat = np.concatenate((adj_mat, arr_temp))

 

이제 파이썬 유저가 보아도 편-안한 adjacency matrix 가 만들어졌다.

 

사실 보기 편하려고 한 것도 있지만, Boeing 아저씨가 이렇게 만들라고 해서 그런거다.

 

이제 network 분석을 시작한다.

 

다시 한 번, 감사하게도 위와 같은 형태의 adjacency matrix 만 있으면 곧장 네트워크 분석이 가능하다.

 

import networkx as nx

G = nx.from_numpy_array(adj_mat)

# Calculate betweenness centrality
betweenness_centrality = nx.betweenness_centrality(G)
print("Betweenness Centrality:", betweenness_centrality)

# Calculate closeness centrality
closeness_centrality = nx.closeness_centrality(G)
print("Closeness Centrality:", closeness_centrality)

 

세상 감사하게도 네트워크 분석이 한줄 컷으로 끝난다.

 

결과물은 아래와 같은 dictionary 를 준다.

 

3. Join to SHP

우리는 결과물을 지도로 봐야만 알아먹을 수 있다.

아까 기억해뒀던 도로id, 데이터프레임 id 테이블을 가지고 결과물을 붙여야 한다.

 

원래의 도로 shp 에는 도로 id가 있고,

우리가 얻은 network dictionary 에는 데이터프레임 id 가 있다.

 

일단 데이터프레임id 와 도로id 와 centrality 값을 합쳐보자. 

아까 깔아둔 복선 df_ids 을 에 centrality 값을 join 해본다.

 

def attach_centrality(i, network_dict):
    '''
    Return centrality value based on dataframe id and dictionary
    Input: Int, dictionary / Output: Float
    '''
    return network_dict[i]

df_ids['betweenness'] = df_ids['df_id'].apply(lambda x: attach_centrality(x, betweenness_centrality))
df_ids['closeness'] = df_ids['df_id'].apply(lambda x: attach_centrality(x, closeness_centrality))

 

나는 apply 함수를 매우 좋아한다.

 

그리고 geometry 를 포함하고 있는 파일을 불러와 id에 맞춰서 JOIN 해주면 완료.

여기서 한번 더, SQL을 사용해서 geometry 를 포함하고 있는 공간 데이터베이스를 가져와

도로 id를 사용해 join 해주겠다.

table_name = 'buk_street'
gdf_query = f'SELECT id, ST_AsText(geom) geometry FROM {table_name}'
street_df = fetch_table(gdf_query, conn)
street_df = street_df.rename(columns = {0: 'rd_id', 1: 'geometry'})
street_df_merge = street_df.merge(df_ids, on = 'rd_id')

 

아름다운 마무리.

 

(JOIN 은 SQL에서 해도 좋지만, 우리 빵셔틀은 데이터프레임을 SQL에 갖다주는 건 매우 어려워하므로, python 에서 join 했다)

 

결과물은 csv도 좋고, geodataframe 으로 변환해서 shp 로 저장해도 좋다.

street_df_merge.to_csv('../Paper2/buk_street_networkx.csv')
street_gdf_merge = gpd.GeoDataFrame(street_df_merge, geometry = 'geometry')
street_gdf_merge.to_file('../Paper2/buk_street_networkx.shp')

아, 이쁘다.