비선형 모형과 추세 분석(아파트 매매 데이터 기반)

지난번 포스팅을 통해서 얻은 데이터와 그 이전 포스팅의 아파트 매매가에 미치는 층수, 크기, 년도 효과에 대한 분석의 후속 분석으로 같은 데이터를 기반을 하는 분석이지만 2015년 Q1의 데이터가 어느정도 모였으니 이의 가격동향과 더불어 다중 회귀모형의 비선형적인 효과를 좀더 다른 방식의 알고리즘으로 모델링 해보고 시각화 해보는 과정을 거쳐보도록 하겠다.

참고로, 분석을 위해 2015년 2월까지의 아파트 매매 데이터를 다시 크롤링했다.

노파심에서 이야기 하지만, 이 분석은 부동산에 전혀 지식이 없는 한 개인이 내가 살고 있는 아파트 단지와 근처 몇몇 단지를 분석한 결과이므로 확대 해석은 매우 위험한 생각이라 말씀드리고 싶다. 그저 생활 데이터 분석의 하나일 뿐이다.

suppressPackageStartupMessages(
{
library(data.table)
library(mgcv)
library(party)
library(randomForest)
library(ggplot2)
library(dplyr)
library(extrafont)
library(randtests)
})

#http://freesearch.pe.kr/archives/4298 에서 만든 데이터를 로딩한다. 
load('result_sales_dt.RData')

setnames(result_sales_dt, 1:10,c('si_gun_gu', 'm_bun', 's_bun', 'dangi', 'area','cont_date', 'price', 'floor', 'year_of_construct', 'road_nm'))

result_sales_dt[,price:=as.numeric(price)]
result_sales_dt[,floor:=as.numeric(floor)]
result_sales_dt[,mm:=substr(yyyymm, 5,6)]
result_sales_dt[between(as.numeric(mm), 1, 3), qrt:='Q1']
result_sales_dt[between(as.numeric(mm), 4, 6), qrt:='Q2']
result_sales_dt[between(as.numeric(mm), 7, 9), qrt:='Q3']
result_sales_dt[between(as.numeric(mm), 10, 12), qrt:='Q4']
result_sales_dt[,yyyyqrt:=paste0(substr(yyyymm, 1,4), qrt)]
result_sales_dt[,yyyy:=factor(substr(yyyymm, 1,4))]
#이전 포스팅과 다르게 쿼터로 날짜는 자른다. 
result_sales_dt[,yyyyqrt:=factor(yyyyqrt)]


my_apt <- result_sales_dt[si_gun_gu %chin% " 서울특별시 관악구 신림동" & m_bun == 1730]

ggplot(my_apt, aes(floor, price)) + geom_point(aes(colour=yyyy)) 

plot of chunk unnamed-chunk-1

일단 층별 매매가를 보면 위와 같은 비선형적인 패턴이 보인다. 익히 알고 있는 로열층이 가격대가 높은것을 실제 데이터를 통해 확인이 가능하다.

이런 사전지식 하에서 이전 포스팅에서 했던 층변수에 2차 함수를 할당한 회귀모형 시각화를 확인해 보자.

#아웃라이어 레코드 제거 (잔차 분석 참고 ...)
my_apt <- my_apt[-c(54,67,202)]

lm_mdl <- lm(price ~ yyyyqrt + I(floor^2) + floor  + area, data=my_apt)

#지난번 하지 않았던 오차상관성 검정... 상관성은 없다고 봐도 무방하다. 
runs.test(resid(lm_mdl))
## 
##  Runs Test
## 
## data:  resid(lm_mdl)
## statistic = 1.0647, runs = 226, n1 = 214, n2 = 214, n = 428,
## p-value = 0.287
## alternative hypothesis: nonrandomness
levelss <- levels(factor(my_apt$yyyyqrt))

aggrs_fun <- function(floor, area, mdls){
  predict(mdls, data.table(yyyyqrt=factor(rep("2015Q1"),levels=levelss), floor=floor, area=area))
}


floor_seq <- seq(1,24, length=30)
area_seq <- seq(59,140, length=30)

persp(floor_seq,area_seq, z=outer(floor_seq,area_seq, aggrs_fun, lm_mdl), col="green", ticktype="detailed",shade=TRUE, theta=-25,phi = 5, xlab="floor", ylab="area", zlab = "price")

plot of chunk _lm

나쁘지 않아 보인다. 사실 다중 회귀 모형에서 2차 함수 형태를 직접 데이터 시각화를 통해서 발견해 채용한 그 과정이 흥미로웠다고 생각한다.

regression tree를 동일한 방식으로 살펴보자!

ctree_mdl <- ctree(price ~ yyyyqrt + floor + area, data=my_apt)

persp(floor_seq,area_seq, z=outer(floor_seq,area_seq, aggrs_fun, ctree_mdl), col="green", ticktype="detailed",shade=TRUE, theta=-25,phi = 5, xlab="floor", ylab="area", zlab = "price")

plot of chunk tree

물론 예상했지만 거의 비선형적인 패턴을 추출하지 못한것으로 보인다. 왜냐면 6층부터 24층까지 동일한 층수 효과를 적용했기 때문이다.

rf_mdl <- randomForest(price ~ yyyyqrt + floor  + area, data=my_apt)


persp(floor_seq,area_seq, z=outer(floor_seq,area_seq, aggrs_fun, rf_mdl), col="green", ticktype="detailed",shade=TRUE, theta=-25,phi = 5, xlab="floor", ylab="area", zlab = "price")

plot of chunk rf

랜덤포레스트는 일반 regression tree보다는 비선형적인 특징이 잘 설명되는것을 볼 수 있으나 사실 눈으로 확인한 패턴에 가깝지는 않다.

gam_mdl <- gam(price ~ yyyyqrt + s(floor) + area, data=my_apt)


persp(floor_seq,area_seq, z=outer(floor_seq,area_seq, aggrs_fun, gam_mdl), col="green", ticktype="detailed",shade=TRUE, theta=-25,phi = 5, xlab="floor", ylab="area", zlab = "price")

plot of chunk gam

개인적으로는 GAM(generalized additive model)에서 Spline형태로 표현된 비선형 패턴이 가장 실제 데이터에서 보였던 패턴과 가장 가깝다고 생각한다.

그럼 몇층이 가장 가치가 높다고 볼 수 있을까?

plot(gam_mdl)

plot of chunk unnamed-chunk-2

8층에서 20층 사이는 오차 범위가 겹치니 이중에 본인이 좋아하는 층을 선택하는게 좋을것 같다.


2015년 1Q가 다 지나가고 데이터도 좀 쌓였으니 시간에 따른 매매 가격 변동 추이를 자세하게 볼 필요가 있을 것 같다.

theme_set(theme_bw(base_family = "UnBatang"))

ptbl <- data.frame(summary(gam_mdl)$p.table)
ptbl$params <- rownames(ptbl)
ptbl <- data.table(ptbl)  

ptbl_yyyy <- ptbl[params %like% "^yyyy"]

ggplot(ptbl_yyyy, aes(x=params,y=Estimate)) + geom_errorbar(aes(ymin=Estimate - 2*`Std..Error`, ymax=Estimate + 2*`Std..Error`), linetype=2) + geom_point() + geom_line(aes(group=1)) + theme(axis.text.x=element_text(angle=90)) + geom_smooth(aes(group=1), method=loess) + ggtitle("서울특별시 관악구 신림동 1730 신림 푸르지오 1차 ")

plot of chunk unnamed-chunk-3

위 그래프는 해당 단지에서 면적과 층수에 대한 효과를 제외한 시간에 따른 가격 변동을 2006년을 기준(분양 이후)으로 표현한 도식화된 그래프이다. 과연 금년 3월에 실시된 역대 최저 주택담보대출 금리의 효과가 앞으로 어떻게 반영이 될지 매우 궁금하다.

비교를 위해 근처 두 아파트 단지의 추이를 같은 방법으로 비교했다. 2007년, 2008년이 매매가가 고점이였다는 걸 감안하고 보면 된다.

my_apt2 <- result_sales_dt[si_gun_gu %chin% " 서울특별시 관악구 신림동" & m_bun == 1736]


gam_mdl2 <- gam(price ~ yyyyqrt + s(floor) + area, data=my_apt2)


ptbl <- data.frame(summary(gam_mdl2)$p.table)
ptbl$params <- rownames(ptbl)
ptbl <- data.table(ptbl)  

ptbl_yyyy <- ptbl[params %like% "^yyyy"]

ggplot(ptbl_yyyy, aes(x=params,y=Estimate)) + geom_errorbar(aes(ymin=Estimate - 2*`Std..Error`, ymax=Estimate + 2*`Std..Error`), linetype=2) + geom_point() + geom_line(aes(group=1)) + theme(axis.text.x=element_text(angle=90)) + geom_smooth(aes(group=1), method=loess) + ggtitle("서울특별시 관악구 신림동 1736 신림 푸르지오 2차")

plot of chunk unnamed-chunk-4

my_apt3 <- result_sales_dt[si_gun_gu %chin% " 서울특별시 관악구 신림동" & m_bun == 1735]


gam_mdl3 <- gam(price ~ yyyyqrt + s(floor) + area, data=my_apt3)


ptbl <- data.frame(summary(gam_mdl3)$p.table)
ptbl$params <- rownames(ptbl)
ptbl <- data.table(ptbl)  

ptbl_yyyy <- ptbl[params %like% "^yyyy"]

ggplot(ptbl_yyyy, aes(x=params,y=Estimate)) + geom_errorbar(aes(ymin=Estimate - 2*`Std..Error`, ymax=Estimate + 2*`Std..Error`), linetype=2) + geom_point() + geom_line(aes(group=1)) + theme(axis.text.x=element_text(angle=90)) + geom_smooth(aes(group=1), method=loess) + ggtitle("서울특별시 관악구 신림동 1735 관악산 휴먼시아2차")

plot of chunk unnamed-chunk-4

세단지 모두 최근의 추세가 반등을 하느냐, 못하느냐 가장 중요한 시점임에는 분명한것 같다. 정부의 정책과 더불어 앞으로의 매매 추이가 매우 궁금해진다.

CC BY-NC 4.0 비선형 모형과 추세 분석(아파트 매매 데이터 기반) by from __future__ import dream is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License.