← 메인 화면프로젝트 자료·코드
← 자료 목록

jochiwon_alley.py

프로그램·코드 · 16.5 KB원본 다운로드

원본 내용을 읽기 전용으로 표시합니다.

  1. # 조치원 상권 분석용 독립 모듈!
  3. import os
  5. import numpy as np
  6. import pandas as pd
  7. import matplotlib.pyplot as plt
  8. from sklearn.cluster import DBSCAN
  9. from sklearn.neighbors import BallTree
  10. from sklearn.model_selection import train_test_split
  11. from sklearn.ensemble import RandomForestClassifier
  12. from sklearn.inspection import permutation_importance
  13. from sklearn.metrics import roc_auc_score, classification_report
  15. from common import *
  18. # 4단계 : 조치원읍 골목상권 분석
  19. # 입지계수(LQ), 골목 덩어리 찾기, 업종별 소멸률, 골목 vs 대로, 핵심 골목 2곳 비교
  21. def 조치원분석(now, prev):
  22. # 조치원읍 골목상권을 자세히 분석 (보고서 4장, 활성화 방안의 근거 대부분이 여기서 나옴)
  23. # 입지계수 : 조치원에 많은/부족한 업종
  24. # DBSCAN : 골목 단위 덩어리 찾고 덩어리별 지표 계산
  25. # 업종별 소멸률, 골목 vs 대로, 핵심 골목 2곳 비교
  26. # 반환 파라미터: T(조치원 올해), Tp(조치원 작년), alley(골목 덩어리별 지표 표)
  27. log("\n" + "=" * 70)
  28. log(f"5. 대상 지역 정밀 분석: {TARGET}")
  29. log("=" * 70)
  31. T = now[now["행정동명"] == TARGET].copy() # 조치원 올해
  32. Tp = prev[prev["행정동명"] == TARGET].copy() # 조치원 작년
  33. new_city = now[now["권역"] == "동(신도심)"]
  34. log(f"{TARGET} 업소 {len(T):,}개 (세종 전체의 {len(T)/len(now):.1%}), "
  35. f"1년 소멸률 {Tp['소멸'].mean():.1%}, 체인 비율 {T['체인'].mean():.1%} "
  36. f"(신도심 체인 비율 {new_city['체인'].mean():.1%})")
  38. # '철물, 주점, 세탁, 식료품은 많고 교육, 의원은 부족'이 나옴 -> 방안에 작성 예정
  39. # 입지계수(LQ) = 조치원에서 그 업종 비율 / 세종 전체에서 그 업종 비율
  40. # 1보다 크면 조치원에 많은 업종, 작으면 적은 업종
  41. 조치원비율 = T["상권업종중분류명"].value_counts(normalize=True)
  42. 세종비율 = now["상권업종중분류명"].value_counts(normalize=True)
  43. lq = pd.DataFrame({"LQ": (조치원비율 / 세종비율).dropna()})
  44. lq["조치원업소수"] = T["상권업종중분류명"].value_counts()
  45. lq["세종업소수"] = now["상권업종중분류명"].value_counts()
  46. lq.index.name = "상권업종중분류명"
  47. lq = lq[lq["세종업소수"] >= 50] # 업소 수 너무 적은 업종은 빼기
  48. lq = lq.sort_values("LQ", ascending=False)
  49. log("\n[입지계수 상위(조치원 특화 업종)]")
  50. log(lq.head(12).round(2).to_string())
  51. log("\n[입지계수 하위(조치원에 상대적으로 부족한 업종)]")
  52. log(lq.tail(12).round(2).to_string())
  53. lq.round(3).to_csv(os.path.join(결과폴더, "table_조치원_입지계수.csv"), encoding="utf-8-sig")
  55. 보기 = pd.concat([lq.head(10), lq.tail(10)])
  56. 색 = []
  57. for v in 보기["LQ"]:
  58. if v >= 1:
  59. 색.append(COLORS[1])
  60. else:
  61. 색.append(COLORS[0])
  62. plt.figure(figsize=(7, 7))
  63. plt.barh(보기.index[::-1], 보기["LQ"][::-1], color=색[::-1])
  64. plt.axvline(1, color="gray")
  65. plt.title(f"{TARGET} 업종별 입지계수(LQ) 상·하위 10\n(1보다 크면 세종 평균보다 많음)")
  66. 그림저장("08_조치원_입지계수.png")
  68. # DBSCAN으로 골목상권 덩어리 찾기
  69. # 행정동 단위로는 골목별 차이가 안 보여서, 가게들이 붙어 있는 구역을 좌표로 찾음
  70. # K-Means와 달리 군집 수를 미리 안 정해도 되고, 띄엄띄엄 있는 가게는 -1(산재)로 빠져서 골목 찾기에 맞음
  71. # eps = 60m (걸어서 1분 정도), 최소 15개
  72. # (40, 80, 100m 등도 해봤는데 60m가 제일 적당했음 그냥 이걸로 하셈)
  73. coords = T[["x", "y"]].values
  74. db = DBSCAN(eps=60, min_samples=15)
  75. db.fit(coords)
  76. T["골목군집"] = db.labels_ # -1은 어느 덩어리에도 안 들어간 업소
  77. 골목수 = len(set(db.labels_)) - (1 if -1 in db.labels_ else 0)
  78. log(f"\nDBSCAN 결과: 골목상권 덩어리 {골목수}개, 군집 외(산재) 업소 {(db.labels_ == -1).sum()}개")
  80. # 작년 업소는 군집 번호가 없어서 제일 가까운 올해 업소(60m 안)의 번호를 따라감
  81. tree = BallTree(coords)
  82. dist, idx = tree.query(Tp[["x", "y"]].values, k=1)
  83. Tp["골목군집"] = np.where(dist[:, 0] <= 60, T["골목군집"].values[idx[:, 0]], -1)
  85. # 골목 덩어리마다 지표 계산
  86. # 덩어리 하나씩 돌면서 크기, 업종 비율, 소멸률, 신규비율, 청년소비형/생활밀착형 비율을 구해서 표로 만듦
  87. # 여기서 '행복길 = 음식 54%, 소멸 17%', '새내로 = 안정적이지만 신규 11%'가 나옴 -> 방안에 활용ㄱㄱ
  88. YOUTH = ["비알코올", "주점", "유원지·오락"] # 카페, 주점, 오락
  89. DAILY = ["이용·미용", "식료품 소매", "섬유·의복·신발 소매", "세탁", "종합 소매"] # 단순 소매점들
  90. 골목목록 = sorted(T.loc[T["골목군집"] >= 0, "골목군집"].unique())
  91. 행들 = []
  92. for k in 골목목록:
  93. a = T[T["골목군집"] == k] # 올해 이 골목 업소
  94. b = Tp[Tp["골목군집"] == k] # 작년 이 골목 업소
  95. 행 = {}
  96. 행["골목군집"] = k
  97. 행["업소수"] = len(a)
  98. 행["대표도로"] = ", ".join(a["도로"].value_counts().index[:3])
  99. 행["음식비율"] = (a["상권업종대분류명"] == "음식").mean()
  100. 행["소매비율"] = (a["상권업종대분류명"] == "소매").mean()
  101. 행["수리개인비율"] = (a["상권업종대분류명"] == "수리·개인").mean()
  102. 행["체인비율"] = a["체인"].mean()
  103. 행["신규비율"] = a["신규"].mean()
  104. 행["위도"] = a["위도"].mean()
  105. 행["경도"] = a["경도"].mean()
  106. 행["소멸률"] = b["소멸"].mean()
  107. 행["업종다양성"] = 다양성(a["상권업종중분류명"])
  108. 행["주요업종"] = ", ".join(a["상권업종중분류명"].value_counts().index[:4])
  109. 행["업소수_2025"] = len(b)
  110. 행["순증가율"] = len(a) / len(b) - 1
  111. 행["청년소비형비율"] = a["상권업종중분류명"].isin(YOUTH).mean()
  112. 행["생활밀착형비율"] = a["상권업종중분류명"].isin(DAILY).mean()
  113. 행들.append(행)
  114. alley = pd.DataFrame(행들).set_index("골목군집")
  115. alley = alley.sort_values("업소수", ascending=False)
  116. log("\n[조치원 골목상권 덩어리별 지표]")
  117. log(alley.drop(columns=["위도", "경도"]).round(3).to_string())
  118. alley.round(4).to_csv(os.path.join(결과폴더, "table_조치원_골목군집.csv"), encoding="utf-8-sig")
  120. # 조치원 업종(중분류)별 소멸률
  121. # 어떤 업종이 잘 망하고 어떤 업종이 오래 가는지 알아보는 코드
  122. mid = pd.DataFrame()
  123. mid["업소수_2025"] = Tp.groupby("상권업종중분류명").size()
  124. mid["소멸률"] = Tp.groupby("상권업종중분류명")["소멸"].mean()
  125. mid = mid.join(T.groupby("상권업종중분류명").size().rename("업소수_2026"), how="outer")
  126. mid = mid.join(T.groupby("상권업종중분류명")["신규"].sum().rename("신규수"), how="outer")
  127. mid = mid.fillna(0)
  128. mid["순증감"] = mid["업소수_2026"] - mid["업소수_2025"]
  129. mid = mid[mid["업소수_2025"] >= 20] # 작년에 20개 이상 있던 업종만
  130. mid = mid.sort_values("소멸률", ascending=False)
  131. log("\n[조치원 업종별 1년 소멸률(2025년 업소 20개 이상 업종)]")
  132. log(mid.round(3).to_string())
  133. mid.round(4).to_csv(os.path.join(결과폴더, "table_조치원_업종별_소멸률.csv"), encoding="utf-8-sig")
  135. # 골목(~길로 되어있는거) vs 대로(~로 단, 다솜로 같이 작게 되어 있는 것도 있어서 주의하세용) 비교
  136. T["도로유형"] = np.where(T["도로"].str.endswith("길"), "길(골목)", "로(대로)")
  137. Tp["도로유형"] = np.where(Tp["도로"].str.endswith("길"), "길(골목)", "로(대로)")
  138. road = pd.DataFrame()
  139. road["업소수"] = T.groupby("도로유형").size()
  140. road["체인비율"] = T.groupby("도로유형")["체인"].mean()
  141. road["신규비율"] = T.groupby("도로유형")["신규"].mean()
  142. road["소멸률"] = Tp.groupby("도로유형")["소멸"].mean()
  143. road["음식비율"] = (T["상권업종대분류명"] == "음식").groupby(T["도로유형"]).mean()
  144. log("\n[조치원 도로유형별 비교]")
  145. log(road.round(3).to_string())
  147. # 제일 큰 골목 2개 자세히 보기 (어떤 업종이 없어지고 생겼는지)
  148. # 행복길에서 없어진 업종과 생긴 업종이 거의 똑같다는 걸 여기서 확인함
  149. core = alley.index[:2].tolist()
  150. for k in core:
  151. a = T[T["골목군집"] == k]
  152. b = Tp[Tp["골목군집"] == k]
  153. 일층 = a.loc[a["1층여부"] >= 0, "1층여부"].mean()
  154. log(f"\n[핵심 골목 {k}: {alley.loc[k, '대표도로']}] 업소 {len(a)}개, "
  155. f"1년 소멸 {b['소멸'].sum()}개 / 신규 {a['신규'].sum()}개, 1층 비율 {일층:.0%}")
  156. 없어진업종 = b[b["소멸"] == 1]["상권업종중분류명"].value_counts().head(6)
  157. 생긴업종 = a[a["신규"] == 1]["상권업종중분류명"].value_counts().head(6)
  158. log(" - 소멸 상위 업종: " + ", ".join(f"{i}({v})" for i, v in 없어진업종.items()))
  159. log(" - 신규 상위 업종: " + ", ".join(f"{i}({v})" for i, v in 생긴업종.items()))
  161. # 핵심 골목 2개vs조치원 전체vs동지역 전체 비교표
  162. nc_prev = prev[~prev["행정동명"].str.endswith(("읍", "면"))] # 작년 동지역
  163. cmp = pd.DataFrame(columns=["1년 소멸률", "신규 비율", "체인 비율", "청년소비형 비율"])
  164. for k in core:
  165. 이름 = f"골목{k}\n({alley.loc[k, '대표도로'].split(',')[0]})"
  166. cmp.loc[이름] = [alley.loc[k, "소멸률"], alley.loc[k, "신규비율"], alley.loc[k, "체인비율"], alley.loc[k, "청년소비형비율"]]
  167. cmp.loc[f"{TARGET} 전체"] = [Tp["소멸"].mean(), T["신규"].mean(), T["체인"].mean(),
  168. T["상권업종중분류명"].isin(YOUTH).mean()]
  169. cmp.loc["신도심(동) 전체"] = [nc_prev["소멸"].mean(), new_city["신규"].mean(), new_city["체인"].mean(),
  170. new_city["상권업종중분류명"].isin(YOUTH).mean()]
  171. cmp = cmp.astype(float)
  172. log("\n[핵심 골목 비교표]")
  173. log(cmp.round(3).to_string())
  175. plt.figure(figsize=(12, 3.4))
  176. for n, col in enumerate(cmp.columns):
  177. plt.subplot(1, 4, n + 1)
  178. plt.barh(cmp.index[::-1], cmp[col][::-1], color=[COLORS[0], COLORS[0], COLORS[1], COLORS[1]])
  179. for i, v in enumerate(cmp[col][::-1]):
  180. plt.text(v, i, f" {v:.0%}", va="center", fontsize=8)
  181. plt.title(col)
  182. plt.xticks([])
  183. if n > 0:
  184. plt.yticks([])
  185. 그림저장("10_핵심골목_비교.png")
  187. # 조치원 업종별 소멸률 그래프
  188. mm = mid.sort_values("소멸률")
  189. plt.figure(figsize=(7, 6))
  190. plt.barh(mm.index, mm["소멸률"], color=COLORS[0])
  191. plt.axvline(Tp["소멸"].mean(), color="gray", linestyle=":")
  192. plt.title(f"{TARGET} 업종별 1년 소멸률 (점선 = 조치원 평균)")
  193. 그림저장("11_조치원_업종별_소멸률.png")
  195. return T, Tp, alley
  198. # 머신러닝 - 랜덤포레스트로 폐업(소멸) 위험 보기
  199. # 작년(2025.6) 업소가 1년 뒤에 없어졌는지(1) 남아있는지(0) 예측
  200. # 매출 같은 정보가 없어서 정확도는 낮을 것 같음
  202. def 소멸예측(prev, Tp, alley):
  203. # 랜덤포레스트로 1년 뒤 소멸 여부를 예측하고, 어떤 조건이 소멸과 관련 있는지 봄
  204. # 출력: 그림 09, table_조치원_골목군집.csv
  205. log("\n" + "=" * 70)
  206. log("6. 머신러닝 ②: 업소 소멸 위험 분류 (Random Forest)")
  207. log("=" * 70)
  209. P_ = prev.copy()
  211. # 모델 입력 변수 만들기
  212. # 작년 업소 기준이라 전처리에서 now 에 만든 것들을 prev 에도 똑같이 만들어 줌
  213. # 1층 여부 (작년 데이터도 같음!)
  214. 일층 = []
  215. for 층 in P_["층정보"]:
  216. if pd.isna(층):
  217. 일층.append(-1)
  218. elif str(층).strip() == "1":
  219. 일층.append(1)
  220. else:
  221. 일층.append(0)
  222. P_["1층여부"] = 일층
  224. 이름개수2 = P_["상호명"].value_counts()
  225. P_["체인"] = ((P_["지점명"].notna()) | (P_["상호명"].map(이름개수2) >= 3)).astype(int)
  226. P_["원도심"] = P_["행정동명"].str.endswith(("읍", "면")).astype(int)
  227. P_["골목길"] = P_["도로"].astype(str).str.endswith("길").astype(int)
  229. # 주변 100m 안에 업소가 몇 개 있는지 (자기 자신은 빼려고 -1)
  230. tree = BallTree(P_[["x", "y"]].values)
  231. P_["주변100m_업소수"] = tree.query_radius(P_[["x", "y"]].values, r=100, count_only=True) - 1
  233. # 주변 100m 안에 같은 업종(중분류)이 몇 개 있는지 => 경쟁 보려고
  234. same = np.zeros(len(P_), dtype=int)
  235. for 업종, 그룹 in P_.groupby("상권업종중분류명"):
  236. t = BallTree(그룹[["x", "y"]].values)
  237. 개수 = t.query_radius(그룹[["x", "y"]].values, r=100, count_only=True) - 1
  238. same[P_.index.get_indexer(그룹.index)] = 개수
  239. P_["주변100m_동종업소수"] = same
  240. P_["동종비중"] = P_["주변100m_동종업소수"] / (P_["주변100m_업소수"] + 1)
  242. # 업종은 글자라서 0, 1로 바꿈
  243. 업종더미 = pd.get_dummies(P_["상권업종대분류명"], prefix="업종").astype(int)
  244. num_cols = ["1층여부", "체인", "원도심", "골목길", "주변100m_업소수", "주변100m_동종업소수", "동종비중"]
  245. Xm = pd.concat([P_[num_cols], 업종더미], axis=1)
  246. y = P_["소멸"].values
  248. # 학습이랑 평가
  249. # 소멸이 13.7%뿐이라 그냥 정확도는 의미가 없어서 ROC-AUC로 평가
  250. # 로지스틱 회귀도 비교해 봤는데 AUC 비슷했음
  251. # 학습 75%, 검증 25%로 나누기
  252. X_tr, X_te, y_tr, y_te = train_test_split(Xm, y, test_size=0.25, stratify=y, random_state=SEED)
  254. # 소멸이 13%밖에 안 돼서 class_weight="balanced" 넣음 힘들다....
  255. rf = RandomForestClassifier(n_estimators=400, max_depth=10, min_samples_leaf=20,
  256. class_weight="balanced", random_state=SEED, n_jobs=-1)
  257. rf.fit(X_tr, y_tr)
  258. proba = rf.predict_proba(X_te)[:, 1] # 소멸할 확률
  259. auc = roc_auc_score(y_te, proba)
  260. log(f"학습 {len(X_tr):,}건 / 검증 {len(X_te):,}건, 소멸 비율 {y.mean():.1%}")
  261. log(f"검증 ROC-AUC = {auc:.3f} (0.5=무작위)")
  262. 예측 = (proba >= 0.5).astype(int)
  263. log(classification_report(y_te, 예측, target_names=["유지", "소멸"], digits=3))
  265. # 해석
  266. # 순열 중요도로 중요한 변수를 찾고, 그 변수별 실제 소멸률을 직접 계산해서 맞는지 다시 확인
  267. # 변수 중요도 (순열 중요도: 변수 하나를 섞었을 때 성능이 얼마나 떨어지는지)
  268. pi = permutation_importance(rf, X_te, y_te, scoring="roc_auc", n_repeats=10, random_state=SEED, n_jobs=-1)
  269. imp = pd.Series(pi.importances_mean, index=Xm.columns).sort_values(ascending=False)
  270. log("[순열 중요도 상위 10]")
  271. log(imp.head(10).round(4).to_string())
  273. top = imp.head(10)
  274. plt.figure(figsize=(7, 4.5))
  275. plt.barh(top.index[::-1], top.values[::-1], color=COLORS[0])
  276. plt.title(f"소멸 위험 모델 변수 중요도 (순열 중요도, AUC={auc:.2f})")
  277. 그림저장("09_소멸위험_변수중요도.png")
  279. # 중요하다고 나온 변수들이 실제로도 소멸률 차이가 있는지 직접 확인
  280. P_["주변밀도구간"] = pd.cut(P_["주변100m_업소수"], [-1, 5, 20, 50, 1000], labels=["0~5", "6~20", "21~50", "51+"])
  281. P_["동종비중구간"] = pd.cut(P_["동종비중"], [-0.01, 0.05, 0.15, 1], labels=["5%이하", "5~15%", "15%초과"])
  282. log("\n[조건별 실제 1년 소멸률]")
  283. for c in ["상권업종대분류명", "체인", "1층여부", "원도심", "골목길", "주변밀도구간", "동종비중구간"]:
  284. t = P_.groupby(c, observed=True)["소멸"].agg(["mean", "size"])
  285. 글 = []
  286. for k, v in t.iterrows():
  287. 글.append(f"{k}={v['mean']:.1%}(n={int(v['size'])})")
  288. log(f"- {c}: " + ", ".join(글))
  290. # 완료!
  291. P_["예측위험"] = rf.predict_proba(Xm)[:, 1]
  292. Tp["예측위험"] = P_.loc[Tp.index, "예측위험"]
  293. alley["평균예측위험"] = Tp[Tp["골목군집"] >= 0].groupby("골목군집")["예측위험"].mean()
  294. log("\n[조치원 골목군집별 평균 예측 위험도]")
  295. log(alley[["대표도로", "업소수", "소멸률", "평균예측위험"]].round(3).to_string())
  296. alley.round(4).to_csv(os.path.join(결과폴더, "table_조치원_골목군집.csv"), encoding="utf-8-sig")
파일 확인 정보

SHA-256: 280c92ee8a9884908229d1e8e51c31d6ed80c2e5ee88e57e9d1b4487e017c8f9