Success is a long continuous journey.

Google Calendar에 추가한 Goal Alarm이 오질 않네요 여튼 그래서 수동으로 다시 시작합니다. 5일차때 다룬 Rossmann Store Sales 풀이 이어서 나가도록 하겠습니다.  저번에 제공 데이터를 str. summary 등의 함수를 가지고 한번 훑어보았데 train데이터만 보여줬는데 마저 test 데이터도 조회해보죠

summary(test)

##        Id            Store          DayOfWeek          Date           
##  Min.   :    1   Min.   :   1.0   Min.   :1.000   Min.   :2015-08-01  
##  1st Qu.:10273   1st Qu.: 279.8   1st Qu.:2.000   1st Qu.:2015-08-12  
##  Median :20544   Median : 553.5   Median :4.000   Median :2015-08-24  
##  Mean   :20544   Mean   : 555.9   Mean   :3.979   Mean   :2015-08-24  
##  3rd Qu.:30816   3rd Qu.: 832.2   3rd Qu.:6.000   3rd Qu.:2015-09-05  
##  Max.   :41088   Max.   :1115.0   Max.   :7.000   Max.   :2015-09-17  
##                                                                       
##       Open            Promo        StateHoliday       SchoolHoliday     
##  Min.   :0.0000   Min.   :0.0000   Length:41088       Length:41088      
##  1st Qu.:1.0000   1st Qu.:0.0000   Class :character   Class :character  
##  Median :1.0000   Median :0.0000   Mode  :character   Mode  :character  
##  Mean   :0.8543   Mean   :0.3958                                        
##  3rd Qu.:1.0000   3rd Qu.:1.0000                                        
##  Max.   :1.0000   Max.   :1.0000                                        
##  NA's   :11

test[is.na(test$Open), ] # Only store 622

##        Id Store DayOfWeek       Date Open Promo StateHoliday SchoolHoliday
##  1: 10752   622         6 2015-09-05   NA     0            0             0
##  2:  9040   622         1 2015-09-07   NA     0            0             0
##  3:  8184   622         2 2015-09-08   NA     0            0             0
##  4:  7328   622         3 2015-09-09   NA     0            0             0
##  5:  6472   622         4 2015-09-10   NA     0            0             0
##  6:  5616   622         5 2015-09-11   NA     0            0             0
##  7:  4760   622         6 2015-09-12   NA     0            0             0
##  8:  3048   622         1 2015-09-14   NA     1            0             0
##  9:  2192   622         2 2015-09-15   NA     1            0             0
## 10:  1336   622         3 2015-09-16   NA     1            0             0
## 11:   480   622         4 2015-09-17   NA     1            0             0

test$Open[test$Store == 622]

##  [1]  1  0  1  1  1  1  1  1  0  1  1  1  1  1  1  0  1  1  1  1  1  1  0
## [24]  1  1  1  1  1  1  0  1  1  1  1  1 NA  0 NA NA NA NA NA NA  0 NA NA
## [47] NA NA

여기에서 보면 위에 제가 볼드체로 표시한 부분 즉 실제 rows수를 보면 train데이터가 월등히 많네요. 1017209rows에 비해 test데이터는 턱없이 부족하고 그리고 dashboard에 issue로 올라온 것보면 622개의 데이터가 Open 컬럼이 NA로 되어있어 못쓴다라고 올라와있는데요. 하지만 이 row를 못쓰는건 아닙니다.  여기를 보면 아시다시피 Sale라는 컬럼이 닫혀있다면 0일수가 없으니 1로 간주하면 될거 같습니다. 이런걸 data munging이라고도 하는데 저희가 도메인을 파악하고 어떤 이유에서인지 알수 없거나 정리가 안된걸 직접 수정하여 trim처리를 하는 센스는 있어야 합니다.  게다가 테스트 데이터에는 Customers컬럼이 아예 빠져있는걸 알수 있습니다. 사후 데이터이기 때문이라네요 ?

여튼  그래서 아래와 같은 처리 합니다. NA인 녀석들을 싹다!!!!  1로 만들어버립시다

test[is.na(test)] <- 1

위의 문법은 test의 컬럼중 NA를 다 1로 채우는걸로 이해하면 됩니다

그다음에 다시

test[is.na(test$Open), ]

# Unique values per column
train[, lapply(.SD, function(x) length(unique(x)))]

##    Store DayOfWeek Date Sales Customers Open Promo StateHoliday
## 1:  1115         7  942 21734      4086    2     2            4
##    SchoolHoliday
## 1:             2

test[, lapply(.SD, function(x) length(unique(x)))]

##       Id Store DayOfWeek Date Open Promo StateHoliday SchoolHoliday
## 1: 41088   856         7   48    2     2            2             2

# All test stores are also in the train data
sum(unique(test$Store) %in% unique(train$Store)) 

## [1] 856

# 259 train stores are not in the test data
sum(!(unique(train$Store) %in% unique(test$Store))) 

## [1] 259

table(train$Open) / nrow(train) # Percent Open Train

## 
##         0         1 
## 0.1698933 0.8301067

table(test$Open) / nrow(test) # Percent Open Test 

## 
##         0         1 
## 0.1456386 0.8543614

table(train$Promo) / nrow(train) # Percent of the time promo in train

## 
##         0         1 
## 0.6184855 0.3815145

table(test$Promo) / nrow(test) # Percent of the time promo in test

## 
##         0         1 
## 0.6041667 0.3958333

table(train$StateHoliday) / nrow(train) # Percent of the time holiday in train

## 
##           0           a           b           c 
## 0.969475300 0.019917244 0.006576820 0.004030637

table(test$StateHoliday) / nrow(test) # no b and c = no easter holiday and no christmas

## 
##           0           a 
## 0.995619159 0.004380841

table(train$SchoolHoliday) / nrow(train) # Percent of the time school holiday in train

## 
##         0         1 
## 0.8213533 0.1786467

table(test$SchoolHoliday) / nrow(test) # Percent of the time school holiday in test

## 
##         0         1 
## 0.5565129 0.4434871

plot(train$Date, type = "l")

plot(test$Date, type = "l")

# As expected all 856 stores to be predicted daily
all(table(test$Date) == 856) 

## [1] TRUE

hist(train$Sales, 100)

hist(aggregate(train[Sales != 0]$Sales, 
               by = list(train[Sales != 0]$Store), mean)$x, 100, 
     main = "Mean sales per store when store was not closed")

hist(train$Customers, 100)

hist(aggregate(train[Sales != 0]$Customers, 
               by = list(train[Sales != 0]$Store), mean)$x, 100,
     main = "Mean customers per store when store was not closed")

ggplot(train[Sales != 0], aes(x = factor(SchoolHoliday), y = Sales)) +
    geom_jitter(alpha = 0.1) +
    geom_boxplot(color = "yellow", outlier.colour = NA, fill = NA)

벌써 5번째 Kaggler Day이군요 쨌든 오늘은 가장 kaggle script중에 가장 많은 관심을 받은 script 하나를 들고 왔습니다. 물론 이문제 또한 kaggle competition이 이루어진 문제구요. 

문제 설명을 하자면,, 가게, 프로모션, 그리고 경쟁자 정보를 가지고 판매량을 예측하는 겁니다. Rossman이라는 브랜드가 있나봅니다 여튼 여기서 아마 competition 문제를 의뢰한거 같구요.    이 회사가 머하는 곳이냐면 7개의 유럽 국가에 3000개의 약국을 운영하고 있는 회사입니다. 현재 Rossmann 가게 매니저들은 미리 6주되기까지 그들의 daily 판매량을 예상하는 일을 받았다고 합니다. 가게 판매는 많은 인자 즉 프로모션, 경쟁 업체, 학교, 그리고 주의 휴일, 계절, 지역특성에 등등 많은 인자에 의해 영향을 받고  수천명의 개별 매니저들과 그들의 특정 환경을 토대로 판매량을 예측하는함에 있어 결과의 정확성은 꽤 변동이 클것이라고 합니다.

그래서!! 캐글에서 Rossman은 독일에 있는 1,115개의 약국의 6개월간의 판매량을 예측하는걸 목표로 한다고 하네요? 신뢰할만한 판매 예측은 가게 매니저에게 생산성과 동기부여를 증가시킬수 있는 효율적인 직원 스케쥴을 만들수 있게 해준다고 하네요. Rossmann이 견고한 예측 모델을 만들도록 도와주면서 우리는 그들에게 고객과 그들의 팀에게 무엇이 가장 중요한지 알려주고자 합니다. 

자 문제 설명은 마쳤구요.. 여기 가보면 데이터셋을 확보할수 있습니다.

데이터셋을 보니 

손님 수, 가게 오픈여부, 주휴일, 학교휴일, 가게타입. assortment(무슨 factor중의 하나겠죠?) 그리고 경쟁업체와의 거리, 그리고 가장가까운 경쟁업체의 오픈 년월일. 프로모션, 프로모션2, 프로모션 주기가 주어져있습니다..

과연 이걸로 얻을수 있을까... 궁금하시죠? 그럼 award를 흠 그러니까 1등한분거는 공개가 안되는거 같으니.. 35k달러의 후보가 될뻔 한 소스를 한번 살펴볼까요.....!!!!

탐색전 분석(Exploratory Analysis)을 해보죠!!

사용한 라이브러리 들이 나오고요 data.table을 사용하여 속도를 높였다고 합니다. 특정 해석을 가능케 하기 위해 unmasked data가 중요하다합니다.

library(data.table)
library(zoo)

## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library(forecast)

## Loading required package: timeDate
## Loading required package: methods
## This is forecast 6.2

저번에는 ggplot을 보았는데 이번엔 ggplot2가 나왔습니다. 그리고 문제에서 주어지는 데이터가 이번엔 아예 test / train 으로 나뉘어져 있군요. 그리고 별도로 store라는 데이터셋이 있어서 가게 대장 정보가 별도로 잇네요.

library(ggplot2)
test <- fread("../input/test.csv")
train <- fread("../input/train.csv")
store <- fread("../input/store.csv")

아래는 입력받은 데이터셋을 간략하게 보기 좋은 함수죠~~

str(train)

## Classes 'data.table' and 'data.frame':   1017209 obs. of  9 variables:
##  $ Store        : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ DayOfWeek    : int  5 5 5 5 5 5 5 5 5 5 ...
##  $ Date         : chr  "2015-07-31" "2015-07-31" "2015-07-31" "2015-07-31" ...
##  $ Sales        : int  5263 6064 8314 13995 4822 5651 15344 8492 8565 7185 ...
##  $ Customers    : int  555 625 821 1498 559 589 1414 833 687 681 ...
##  $ Open         : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Promo        : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ StateHoliday : chr  "0" "0" "0" "0" ...
##  $ SchoolHoliday: chr  "1" "1" "1" "1" ...
##  - attr(*, ".internal.selfref")=<externalptr>

str(test)

## Classes 'data.table' and 'data.frame':   41088 obs. of  8 variables:
##  $ Id           : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Store        : int  1 3 7 8 9 10 11 12 13 14 ...
##  $ DayOfWeek    : int  4 4 4 4 4 4 4 4 4 4 ...
##  $ Date         : chr  "2015-09-17" "2015-09-17" "2015-09-17" "2015-09-17" ...
##  $ Open         : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Promo        : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ StateHoliday : chr  "0" "0" "0" "0" ...
##  $ SchoolHoliday: chr  "0" "0" "0" "0" ...
##  - attr(*, ".internal.selfref")=<externalptr>

str(store)

## Classes 'data.table' and 'data.frame':   1115 obs. of  10 variables:
##  $ Store                    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ StoreType                : chr  "c" "a" "a" "c" ...
##  $ Assortment               : chr  "a" "a" "a" "c" ...
##  $ CompetitionDistance      : int  1270 570 14130 620 29910 310 24000 7520 2030 3160 ...
##  $ CompetitionOpenSinceMonth: int  9 11 12 9 4 12 4 10 8 9 ...
##  $ CompetitionOpenSinceYear : int  2008 2007 2006 2009 2015 2013 2013 2014 2000 2009 ...
##  $ Promo2                   : int  0 1 1 0 0 0 0 0 0 0 ...
##  $ Promo2SinceWeek          : int  NA 13 14 NA NA NA NA NA NA NA ...
##  $ Promo2SinceYear          : int  NA 2010 2011 NA NA NA NA NA NA NA ...
##  $ PromoInterval            : chr  "" "Jan,Apr,Jul,Oct" "Jan,Apr,Jul,Oct" "" ... 

##  - attr(*, ".internal.selfref")=<externalptr>

# head(train); tail(train)
# head(test); tail(test)
train[, Date := as.Date(Date)]
test[, Date := as.Date(Date)]
store

##       Store StoreType Assortment CompetitionDistance
##    1:     1         c          a                1270
##    2:     2         a          a                 570
##    3:     3         a          a               14130
##    4:     4         c          c                 620
##    5:     5         a          a               29910
##   ---                                               
## 1111:  1111         a          a                1900
## 1112:  1112         c          c                1880
## 1113:  1113         a          c                9260
## 1114:  1114         a          c                 870
## 1115:  1115         d          c                5350
##       CompetitionOpenSinceMonth CompetitionOpenSinceYear Promo2
##    1:                         9                     2008      0
##    2:                        11                     2007      1
##    3:                        12                     2006      1
##    4:                         9                     2009      0
##    5:                         4                     2015      0
##   ---                                                          
## 1111:                         6                     2014      1
## 1112:                         4                     2006      0
## 1113:                        NA                       NA      0
## 1114:                        NA                       NA      0
## 1115:                        NA                       NA      1
##       Promo2SinceWeek Promo2SinceYear    PromoInterval
##    1:              NA              NA                 
##    2:              13            2010  Jan,Apr,Jul,Oct
##    3:              14            2011  Jan,Apr,Jul,Oct
##    4:              NA              NA                 
##    5:              NA              NA                 
##   ---                                                 
## 1111:              31            2013  Jan,Apr,Jul,Oct
## 1112:              NA              NA                 
## 1113:              NA              NA                 
## 1114:              NA              NA                 
## 1115:              22            2012 Mar,Jun,Sept,Dec

train <- train[order(Date)]
test <- test[order(Date)]
summary(train)

##      Store          DayOfWeek          Date                Sales      
##  Min.   :   1.0   Min.   :1.000   Min.   :2013-01-01   Min.   :    0  
##  1st Qu.: 280.0   1st Qu.:2.000   1st Qu.:2013-08-17   1st Qu.: 3727  
##  Median : 558.0   Median :4.000   Median :2014-04-02   Median : 5744  
##  Mean   : 558.4   Mean   :3.998   Mean   :2014-04-11   Mean   : 5774  
##  3rd Qu.: 838.0   3rd Qu.:6.000   3rd Qu.:2014-12-12   3rd Qu.: 7856  
##  Max.   :1115.0   Max.   :7.000   Max.   :2015-07-31   Max.   :41551  
##    Customers           Open            Promo        StateHoliday      
##  Min.   :   0.0   Min.   :0.0000   Min.   :0.0000   Length:1017209    
##  1st Qu.: 405.0   1st Qu.:1.0000   1st Qu.:0.0000   Class :character  
##  Median : 609.0   Median :1.0000   Median :0.0000   Mode  :character  
##  Mean   : 633.1   Mean   :0.8301   Mean   :0.3815                     
##  3rd Qu.: 837.0   3rd Qu.:1.0000   3rd Qu.:1.0000                     
##  Max.   :7388.0   Max.   :1.0000   Max.   :1.0000                     
##  SchoolHoliday     
##  Length:1017209    
##  Class :character  
##  Mode  :character  
##                    
##                    
## 

summary(test)

##        Id            Store          DayOfWeek          Date           
##  Min.   :    1   Min.   :   1.0   Min.   :1.000   Min.   :2015-08-01  
##  1st Qu.:10273   1st Qu.: 279.8   1st Qu.:2.000   1st Qu.:2015-08-12  
##  Median :20544   Median : 553.5   Median :4.000   Median :2015-08-24  
##  Mean   :20544   Mean   : 555.9   Mean   :3.979   Mean   :2015-08-24  
##  3rd Qu.:30816   3rd Qu.: 832.2   3rd Qu.:6.000   3rd Qu.:2015-09-05  
##  Max.   :41088   Max.   :1115.0   Max.   :7.000   Max.   :2015-09-17  
##                                                                       
##       Open            Promo        StateHoliday       SchoolHoliday     
##  Min.   :0.0000   Min.   :0.0000   Length:41088       Length:41088      
##  1st Qu.:1.0000   1st Qu.:0.0000   Class :character   Class :character  
##  Median :1.0000   Median :0.0000   Mode  :character   Mode  :character  
##  Mean   :0.8543   Mean   :0.3958                                        
##  3rd Qu.:1.0000   3rd Qu.:1.0000                                        
##  Max.   :1.0000   Max.   :1.0000                                        
##  NA's   :11

데이터가 흠 가게 Id는 정수형이고 다 주 5일제고 뭐...손님수는 다양하네요 500~1500명까지.... 휴일여부는 0,1로 나뉘어져있고 네..

여튼 데이터 구경은 이정도면 충분하고. 넘어갑니다.   그리고 summary함수를 통해서 또 최대/최소/최빈값등등을 보면서 대략적이나마 data distribution 형태를 좀 예측할 수 가 있지요... 

앗!!! 다음에 이어하겠습니다.

좀 더 해야하는데.. 다음 6일차때 추가하겠습니다.

Kaggler 4번째 시도입니다. 이번에는 RF(Random Forest)를 가지고 Facebook competition에서 Awards를 받은 팀의 스크립트에 대해서 정리하고자 합니다. 

간단하게 문제는 페이스북데이터(평방 10km 그리드안에 약 30만개의 check-ins 정보)가 주어지고 시공간 정보에 의존하여 유저가 무슨 비지니스에 체크인하는지를 파악하는 문제이라고 합니다. 웃긴 부분은 place id를 기준으로 분류하면 100k 개의 종류가 있어서 대부분의 supervised learning 기술은 전체 dataset에서 안먹힐 거라고 합니다. 하지만 다행히 대부분의 종류들이 그리드의 특정 지역에서 cluster가 되어서 여기 팀은 그리드안의 작은 직사각형을 선택해서 그 범위안에서 검증을 해보았다고 합니다. 그래서 더 작은 스퀘어에서 약간의 실험적인 데이터 분석을 하고나면 예측을 위해서 random forest를 사용하여 결과 분석할 것이다라고 합니다.

아래는 사용한 라이브러리

library(data.table) #reading in the data
library(dplyr) #dataframe manipulation
library(ggplot2) #viz
library(ranger) #the random forest implementation
library(plotly) #3D plotting
library(tidyr) #dataframe manipulation

fb에 training dataset을 집어넣습니다.

fb <- fread("../input/train.csv", integer64 = "character", showProgress = FALSE)

그리고 나서 아까 말한 부분집하을 데이터에 선택합니다. 그 기준은  250 x 250meters 크기의 특정 직사각형의 데이터를 뽑고

fb %>% filter(x >1, x <1.25, y >2.5, y < 2.75) -> fb
head(fb, 3)

##    row_id      x      y accuracy   time   place_id
## 1:    600 1.2214 2.7023       17  65380 6683426742
## 2:    957 1.1832 2.6891       58 785470 6683426742
## 3:   4345 1.1935 2.6550       11 400082 6889790653

place_id와 시간, x,y는 좌표일테고... 그리고 accuracy는 뭘까요? 

small_train = fb[fb$time < 7.3e5,]
small_val = fb[fb$time >= 7.3e5,] 

여튼 kaggle에서 test 데이터로 돌려서 scoring을 할텐데 좀 더 유사하기 위해서 최근 데이터를 validation으로 두었다고 합니다.

자 분석을 시작해봅시다!!!

ggplot(small_train, aes(x, y )) +
    geom_point(aes(color = place_id)) + 
    theme_minimal() +
    theme(legend.position = "none") +
    ggtitle("Check-ins colored by place_id")

그림이 이쁘게 plotting이 ggplot으로 그려지긴 하나 겹쳐지는 부분이 꽤 있어 보여서 place_id군들이 잘 구분이 안되어질거 같습니다. 그래서!!! 공간에 시간 개념을 넣어서 plot을 그려보자고 하는군요... 즉 x,y말고도 3차원을 만들기위해 z좌표로 hour을 쓰겠다는거군요.. 여튼 아래 소스로 직접 보면 우선 500개 이하인 place_id들은 날리고 place_id list를 ids에 담네요 그리고 small_trainz라는 subset을 만듭니다. 그리고 그걸로 그리는데...

small_train %>% count(place_id) %>% filter(n > 500) -> ids
small_trainz = small_train[small_train$place_id %in% ids$place_id,]

plot_ly(data = small_trainz, x = x , y = y, z = hour, color = place_id,  type = "scatter3d", mode = "markers") %>% layout(title = "Place_id's by position and Time of Day")

type을 scatter3d로 줘서 그런지 이쁘네요...  kaggle에서 plot을 못가져오겠네요. 동적이라서!!!  여기에서 직접 보세요..

그렇지만 캡쳐해서 올리면.. 이렇습니다. 캬..

 중간에 Place_id's by position and Time of Day라고.. 차트가 hovering하면 3차원의 plot을 동적으로 보여주네요.. ㄷㄷㄷ  시간단위는 먼지 모르겠지만 여튼... hour를 넣고 안겹치게 잘 그려지네요!! 와우!! 시간 dimension을 넣었더니 이렇게 clear가 되는군요..  특정 지점에서 체크인이 몇시간 동안 멈췄다가 다시 체크인이 이뤄지고 있는 모습도 보이고요. 어떤 비지니스들은 꽤 낮 내내 높은 사용률이 보이고 비지니스마다 꽤 서로 다른 peaks를 가지고 있는 경향이 있네요. 또한 z=0과 z=24인 시점은 같은 것이라는 것도 유의해서 3d plot을 보라고 하는군요.  

하지만 아직 우리는 매우 많은 클래스들이 있습니다. (아까 많이 제외를 했기때문에..)

다시 small_train set을 가지고 와서 unique로 종류를 보면 770개의 business들이 있고요. 

length(unique(small_train$place_id))

## [1] 770

여기에서 빈도가 낮은 녀석들은 날리고 ...  3개 이하로 일어나는것들은 날리면 종류가 그래도 많이 준다고 합니다.  우리는 validation set을 가지고 있기 때문에 우리는 항상 다시 돌아와서 filter level을 조절할수 있죠. 더 나은 정교한 모델을 만들기 위해서 말이지요~~

small_train %>% count(place_id) %>% filter(n > 3) -> ids
small_train = small_train[small_train$place_id %in% ids$place_id,]

이렇게 필터링하고 나면 15595 트레이닝 샘플과 224개의 클래스를 가지고 이제 ML을 할 준비를 했다고 합니다..... 자 숲으로!!

원래 RF보다 ranger라는모듈(별도의 RF를 구현)을 사용했다고 합니다. 약 5배의 빠른 성능과 메모리 효율성이 좋다고 합니다. 여튼 이건 중요한 얘긴 아니고..

seed는 넘어가고 place_id를 classification의 factor로 사용을 했고 ranger를 이용하여 model_rf라는 RF모듈을 만드네요.  x와 y, accuracy, hour를 넣고 돌립니다. 그리고 random sampling은 100개를 해서 decision tree를 만드네요.

set.seed(131L)
small_train$place_id <- as.factor(small_train$place_id) # ranger needs factors for classification
model_rf <- ranger(place_id ~ x + y + accuracy + hour,
                   small_train,
                   num.trees = 100,
                   write.forest = TRUE,
                   importance = "impurity")

## Growing trees.. Progress: 48%. Estimated remaining time: 33 seconds.

pred = predict(model_rf, small_val)
pred = pred$predictions
accuracy = mean(pred == small_val$place_id) 

오래 걸리는군요.  만든 모델을 가지고 실제 small_val에서 돌려보면 accuracy가 0.5240919라고 하네요. 나쁘지 않다고 합니다. 물론 처음치고 나쁜 건 아니죠.  이 competition의 평가 metric은 3개의 test set의 mean average precision이라고 하네요. 그래서 class에 의한 votes/probabilities을 예측하고 점수를 개선하기 위해 보증된 top 3개의 id를 카운팅한다?? 뭔말인지.. 여튼.. 넘어갑니다 ㅡ.ㅡ

small_val$Correct = (pred == small_val$place_id)

ggplot(small_val, aes(x, y)) +
    geom_point(aes(color = Correct)) + 
    theme_minimal() +
    scale_color_brewer(palette = "Set1") +
    ggtitle("RF Results")

위의 plot은 모델을 가지고 우리가 맞은것과 틀린것을 표시한건데 accuray가 50% 살짝 넘엇었죠?보기에 correct한 checkin이 더 클러스터링이 된거 같다고 합니다. uncorrect한건 보면 균일하게 분배되어있는거 같고. 

그래서.. 어떤 종류의 RF가 잘못되는지 보자고 합니다. 이를 위해서 우리는 id별 accuracy를 보고 또한 id기준의  plot, validation set에 그것들이 얼마나 자주 나타나는지를 기준으로 한 걸 나타낼려고 합니다. 이 모델이 더 인깅가 있을 수록 더 정확하고 아닐 경우는 accuracy가 줄어드는걸 확인할수 있네요.. 봅시다.

#reordering the levels based on counts:
small_val$place_id <- factor(small_val$place_id,
                             levels = names(sort(table(small_val$place_id), decreasing = TRUE)))

small_val %>% 
    ggplot(aes(x = place_id)) + geom_bar(aes(fill = Correct)) + 
    theme_minimal() +
    theme(axis.text.x = element_blank()) +
    ggtitle("Prediction Accuracy by ID and Popularity") +
    scale_fill_brewer(palette = "Set1")

그럼 변수들의  중요성을.. 또한 보면...

data.frame(as.list(model_rf$variable.importance)) %>% gather() %>% 
    ggplot(aes(x = reorder(key, value), y = value)) +
    geom_bar(stat = "identity", width = 0.6, fill = "grey") +
    coord_flip() +
    theme_minimal() +
    ggtitle("Variable Importance (Gini Index)") +
    theme(axis.title.y = element_blank()) 

 흥미로운 점은 x보다 y가 중요하다는 것. 이건 다른 스크립트에서도 많이 관찰이 된다. place_id의 variaince는 y방향보다 x방향에서 더 크다. 이것은 y축이 더 place_id의 predictior로 좋고  RF가  그 스스로 이것을 알아낸다라는 것이다. hour또한 좋은 predictor이고 공간적인 지표보다는 못하다.  마지막으로 accuracy가 중요하다고 하빈다. 이것은 좀 미스테리하다고 합니다. 왜냐하면 우리는 실제 정의를 얻지 못하나 모델은 이게 중요하다고 한다. 마지막까지도... accuracy의 변수가 데이터에서 뭘 뜻하는지 모르고 종료합니다.

그리고.... 여기 아이디어에서 추가로 제안을 아래와 같이 하네요. 이것은 귀찮으니 걍 영어로... --a

한마디로 뭐.. grid사이즈를 다르게 주거나 다른 모델을 써보자는 얘기죠...

Here are a few sugesstions:

• play with how big the grid size is and see if you can get better validation results
• try different areas on the map
• try different models (maybe xgboost or somehing simple like logistic regression)
• split the grid into n*n little squares and run you algorithm in each of the squares - this might take a while though.

이상입니다.

Kaggler가 되어도 분석할 데이터, 도메인 지식이 없으면 아무 소용이 없다. 그래서 Kaggler가 되기 위한 시도로써 또 나는 분석할 대상, 도메인 지식도 쌓고 그리고 관련 데이터도 어디서 확보하면 좋을까를 목표로 또 하나의 습관을 들이기로 한다.

데이터 도메인은 부동산쪽이며 주제는 아파트가격 분석이 되겠다. 예측까지가면 좋을것 같지만 우선 과거 데이터를 가지고 분석정도를 목표로 한다. 

RAW 데이터 source는 다음이 되겠다.

아래는 데이터 분석 관련 통계청 데이터 목록 

1. 가계동향조사 : 가구의 생활수준실태와 그 변동사항을 파악하기 위해서 가계의 수입과 지출을 조사하여 - 국민소비수준 변화의 측정 및 분석 - 소비자물가지수 편제에 필요한 가중치 모집단 자료 등 각종 경제, 사회정책에 필요한 자료를 제공 - 공공사업 시행으로 인한 세입자의 주거대책비 산정 및 국민주택 공급대상의 기준설정 자료

2. 장래인구추계 : 2010년 인구주택총조사 결과를 기초로, 인구동태(출생․사망)와 국제인구이동통계를 활용, 코호트 요인법 ( Cohort components method)에 의해 2060년까지 향후 50년간의 인구규모 및 성․연령별 구조를 추계 ㅇ 2010년 인구주택총조사(11월 1일자)의 성․연령별 인구 보정 - 성․연령별 기준인구(2010년 7월 1일자) 확정 - 인구변동요인(출생․사망․국제이동) 분석 - 출산력․사망력․국제이동 모형 및 장래 변동 수준에 대한 가정설정 - 가정별 인구추계 시산 및 결과 분석

3. 국내인구이동통계 : 인구 이동 상황을 주기적으로 파악하여 인구정책 및 지역간 균형개발을 위한 기초자료로 활용

4. 한국도시통계 : - 전국 시·읍의 인구, 재정, 건설, 사회복지 및 문화 등 각 분야의 주요통계 자료를 수록하여 도시 변화양태 파악 및 도시행정을 위한 기초자료 제공 - 중앙부처, 지방자치단체, 유관기관 등에서 도시행정 수행에 활용

5. 노인(65세 이상)의 연간 총수입, 소득항목별 금액 및 구성비 

  - 일반특성별 (30개)

:전체, 지역별, 동부, 읍·면부, 성별, 남자, 여자, 연령별, 65~69세, 70~74세.....

- 소득항목별 (8개)

: 총수입, 근로소득, 사업소득, 재산소득, 사적이전소득, 공전이전소득, 기타수입, 기타소득

- 항목 (2개)

: 금액, 구성비

6.  각종 지수(통계청 발표)

    아파트 실거래가격지수(2006.1=100.0) 수록주기 월 2006.01~2016.03  아파트 실거래가격지수(2006.1=100.0) 파일 다운로드  아파트 실거래가격지수(2006.1=100.0) 구성정보

    규모별 아파트 실거래가격지수 전월대비 증감률(％) 수록주기 월 2006.02~2016.02  규모별 아파트 실거래가격지수 전월대비 증감률(％) 파일 다운로드  규모별 아파트 실거래가격지수 전월대비 증감률(％) 구성정보

    규모별 아파트 실거래가격지수(2006.1=100.0) 수록주기 월 2006.01~2016.02  규모별 아파트 실거래가격지수(2006.1=100.0) 파일 다운로드  규모별 아파트 실거래가격지수(2006.1=100.0) 구성정보

    아파트 실거래가격지수 전월대비 증감률(％) 수록주기 월 2006.02~2016.03  아파트 실거래가격지수 전월대비 증감률(％) 파일 다운로드  아파트 실거래가격지수 전월대비 증감률(％) 구성정보

 부동산시장 소비자심리  

    부동산시장 소비심리지수 수록주기 월 2011.07~2016.03  부동산시장 소비심리지수 파일 다운로드  부동산시장 소비심리지수 구성정보

    주택시장 소비심리지수 수록주기 월 2011.07~2016.03  주택시장 소비심리지수 파일 다운로드  주택시장 소비심리지수 구성정보

    주택매매시장 소비심리지수 수록주기 월 2011.07~2016.03  주택매매시장 소비심리지수 파일 다운로드  주택매매시장 소비심리지수 구성정보

    주택전세시장 소비심리지수 수록주기 월 2011.07~2016.03  주택전세시장 소비심리지수 파일 다운로드  주택전세시장 소비심리지수 구성정보

UPDATE ( 2016.05.27. "아는만큼 번다 돈이 되는 빅데이터" 책 내용 추가)

데이터를 다루는 기본 지식

비(ratio) : 분자 / 분모 독립적

비율(proportion)  : 분자가 분모에도 포함

율(rate) : 시간개념이 같이 있어 "속도"개념이 추가된 것.

기저효과(base effect) : 기준/비교 시점의 차이로 왜곡이 되는 현상.

%, %p : p가 붙으면 퍼센테이지의 차이를 얘기한다.

INDEX는 기준시점 + 가중치 + 지수식으로 결정

오늘은 Kaggler 되기 두번째 날입니다. 금일 다루는 내용은 어려운 내용이 될듯.. 털썩.  brief하게 더듬고 넘어가는 걸로 하겠습니다.  문제는 3달전 있었던 한 competition에서 우승한 팀의 how to solvec입니다. 1764개의 팀(1925명 참가)중에서 우승을 차지한 팀의 내용을 가져왔습니다.  해설로는 600개의 베이스 모델을 멋지게 혼합하여 우승을 할거라고 확신했다 하는데 봅시다.

문제 요약

Homesite라는 Home Insurance회사에서 어떤 고객이 특정 보험 상품을 구매할거냐를 분석하여 마케팅에 활용하고자 competition을 연거구요.  데이터는 여기에 트레이닝/테스트 데이터가 주어진답니다. 그리고 분석 후 제출할 샘플 양식까지 주어집니다.  제출하면 여기에서 점수가 주어지고 랭킹까지 되네요. 하지만 소스코드를 볼려고 했더니... 코드는 비공개인거 같습니다. 하긴 1등이 천만원의 상금이 걸려 있고 뭔가 이 회사에서 실제로 분석 모델로 사용할 거라 당연히 비공개일듯.. 

하지만 how to solve는 제가 보기에 나름(?) 정리가 잘 되어져 있습니다.  봅시다!!

Feature Engineering라는 것은 ML이 작동하기 위해서 도메인 지식을 가지고 features를 만드는 과정이라고 볼수 있습니다. 그럼 여기에서 feature는 머냐고 하면 관찰하고 있는 object/phenomenon에 대한 individual measurable property를 말하구요.  하지만 이 feature engineering을 사용하기에는 어려웠다고 합니다. feature를 잡아내기 위한 과정을 통하여 feature를 선별하긴 했으나 결국엔 극수의 engineered features를 얻었다고 합니다.  

TSNE과 마지막은은 이해가 잘 안가고. 나머지는 봅시다.  AUC가 처음부터 높게 나와서 신났다고 하네요. 

아래 그림은 PCA로  Red(부정)/Green(긍정)으로 Classification이 되서 말이죠.

그 다음 그림은 TSNE로 한거라고 합니다.

Dataset Modeling

Python Sklearn의 label encoder를 이용하여 categoriclal features를 인코딩하였다고 합니다. 간단하지만 막강한 이 녀석으로 트리 분류기에서도 작동이 꽤 잘 되었다고 합니다. 하지만 linear한 것들에게는 추천하지 안흔다고 합니다. 왜냐하면 우리 또한 one hot encoded features를 생성했기 때문이라고 합니다.. ( 이 말은.. 저도 이해가 안가네요) 결국 우리는 타겟에 영향을 주는 categorical features의 ratio를 구하기 위해서 target encoding도 시도했다고 합니다. 하지만 개선이 되지 않았다고 합니다.  <-- 결국 label encoding의 불편함을 얘기한거 같습니다 trial and error죠.

Training & Ensembling

여기에서부터 흥미진진한데요.  데이터셋의 다른 버젼들에다가 아래와 같은 ML 모델을 돌렸다고 합니다. 정말 잘 알려져있는것들 잔뜩이네요. 저한텐 밑에 3개는 또 모르지만..

위의 3개는 알겠으나 그 밑에 3개는 모르겠습니다. 마지막껀 H20는 뭔가 RF의 새로운 버젼인거 같은데.. 좀 봐야겠고.. ^^ 여튼 약 600개의 모델로 1레벨이 이루어져 있다고 합니다. 100개는 수작업으로 위의 기술한 기술들의 different features와 hyperparameters를 가지고 만들어졌다고 합니다. 그리고 나서 다양성을 고려하여 500개의 나머지 모델을 만드는 로봇(?)을 만들었다고 합니다. 이 로봇은 자동으로 XGBoost, 로지스틱 회귀, NN가지고 임의로 선택된 features를 기반으로 하여 model을 훈련시킨다고 하네요...  -> 이부분은 신기하네요.. feature가 고정적인것도 아닌가봅니다.

여튼 우리의 모델은 5 fold stratified CV로 만들어졌으며(테스트한 train/test 데이터를 5개로 쪼개서 또 트레이닝 시킨다는 말입니다 CV= Cross Validation)  이로써 우리가 개선을 꾀하고 overfitting을 피할수 있게 도와준다 합니다. 게다가 CV로다가 앙상블 모델을 사용할 수 있었다 합니다(자세한 내용은 여기 뭔가 난해합니다. 여기까진 시간관계상 빠지는걸로..)

아래는 앙상블 예인거 같은데  RGF와 XGB의 관계도.. CORRELATION이 꽤 높죠?

그런데 전 이해가 잘안가네요... 그래서 어쩌라고?? 여튼  600개의 모델을 혼합하기 위해서 다양한 방법을 사용했고 몇번의 실패후에 3개의 멋진 조합을 발견했다고 합니다. 바로 고전적인 로지스틱 회귀와 XGB, NN의 조합이었다합니다.. 후아.. 끝....

그래서 요약하자면...

그리고 노하우를 알려주네요. 저는 아래 팁들이 상당히 공유가 꺼려질정도로 마음에 듭니다. 내 머리속에만 넣어두고 나만 활용하고 싶다는...

데이터 사이언스의 영역은 쉬운길이 아니고 어렵습니다. 하지만 아래와 같은 방식과 마인드셋으로 하면 언젠가.. 나도 FAMOUS DATA EXPERT가 될겁니다. 그럼 오늘 연재 끝..

난 일주일에 세번!!!! Kaggle에 올라온 Famous Script를 따라서 해보고 감상평과 내용요약을 하고자 한다.

1. 주제

Kaggle에서는 사용자들이 Public Dataset을 가지고 분석한 스크립트와 보고서를 보고 Voting을 한다. 그중 가장 많은 아마존의 Find Food의 리뷰 분석 실습이다.. 데이터는 스탠포드에서 제공하고 2012년 10월에 아마존 사용자들의 568,454 리뷰 데이터를 모은것이다.

데이터 구성은 csv로 아래와 같다.

• Id
• ProductId - unique identifier for the product
• UserId - unqiue identifier for the user
• ProfileName
• HelpfulnessNumerator - number of users who found the review helpful
• HelpfulnessDenominator - number of users who indicated whether they found the review helpful
• Score - rating between 1 and 5
• Time - timestamp for the review
• Summary - brief summary of the review
• Text - text of the review
• 
  

그리고 내가 선정한 스크립트는 바로 Ben Hamner이란 분이 올린 script가 되시겠다. 참고로 8개의 voting을 받았다. 

library(RSQLite)
library(tm)
library(wordcloud)

library(RSQLite)
db <- dbConnect(dbDriver("SQLite"), "../input/database.sqlite")

reviews <- dbGetQuery(db, "
SELECT *
FROM Reviews
LIMIT 10000")

make_word_cloud <- function(documents) {
  corpus = Corpus(VectorSource(tolower(documents)))
  corpus = tm_map(corpus, removePunctuation)
  corpus = tm_map(corpus, removeWords, stopwords("english"))
  
  frequencies = DocumentTermMatrix(corpus)
  word_frequencies = as.data.frame(as.matrix(frequencies))
  
  words <- colnames(word_frequencies)
  freq <- colSums(word_frequencies)
  wordcloud(words, freq,
            min.freq=sort(freq, decreasing=TRUE)[[400]],
            colors=brewer.pal(8, "Dark2"),
            random.color=TRUE)  
}

png("wordcloud.png")
make_word_cloud(reviews$Text) 
dev.off()

분명 데이터도 sqlite에서 잘 뽑아와서 corpus를 munging하고 wordcloud함수도 잘 생성되었으나... 쉬운예제임에도불구하고 그림이 안 나옴.. 

실패되시겠습니다. 그래도 여튼 원래장소에서 가져와서붙이고 내가 한걸로.. ㅡ.ㅡ