期中小組競賽, Ta-Feng

建立模型 集群分析 交叉驗證 結論

library(dplyr)
library(ggplot2)
library(caTools)
library(Matrix)
library(rpart)
library(rpart.plot)
library(caret)
library(doParallel)

老師原本模型 取 log 分成 TR2 & TS2

rm(list=ls(all=TRUE))
load("data/tf2.rdata")
n = nrow(A2)
set.seed(2018); spl2 = 1:n %in% sample(1:n, round(0.7*n))
A2 = subset(A, A$buy) %>% mutate_at(c("m","rev","amount"), log10)
TR2 = subset(A2, spl2)
TS2 = subset(A2, !spl2)
lm1 = lm(amount ~ ., TR2[,c(2:6,8:10)])
pred =  predict(lm1, TS2)
r2.tr = summary(lm1)$r.sq
SST = sum((TS2$amount - mean(TR2$amount))^ 2)
SSE = sum((predict(lm1, TS2) -  TS2$amount)^2)
r2.ts = 1 - (SSE/SST)
c(r2.tr, r2.ts)

[1] 0.2909908 0.2575966

嘗試增加變數,再分TR2 & TS2

A2$aa = A2$f*A2$m^2*A2$r^2
A2$uu = A2$f*A2$r
A2$dd = A2$m*A2$r^2
A2$cc = A2$s*A2$f
A2$bb = A2$s*A2$m^0.5
A2$vv = A2$m*A2$f^4
A2$ii = A2$f*A2$f^0.5
A2$pp = A2$r*A2$r^2
A2$jj = (A2$s-A2$r)^2*A2$m^2
A2$FF = (A2$s-A2$r)^2*A2$f*A2$s^3
TR2 = subset(A2, spl2)
TS2 = subset(A2, !spl2)

cx=c(2:10, 12:21)
colnames(TR2[,cx])

 [1] "r"      "s"      "f"      "m"      "rev"    "raw"    "age"   
 [8] "area"   "amount" "aa"     "uu"     "dd"     "cc"     "bb"    
[15] "vv"     "ii"     "pp"     "jj"     "FF"    

lm1 = lm(amount ~ ., TR2[,cx])
summary(lm1)

Call:
lm(formula = amount ~ ., data = TR2[, cx])

Residuals:
     Min       1Q   Median       3Q      Max 
-1.83854 -0.22761  0.04917  0.27798  1.51672 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.816e+00  1.064e-01  17.060  < 2e-16 ***
r            5.197e-03  1.141e-03   4.556 5.29e-06 ***
s           -2.449e-02  3.093e-03  -7.917 2.72e-15 ***
f            2.419e-02  3.619e-02   0.668  0.50387    
m            2.716e-01  1.224e-01   2.220  0.02647 *  
rev          3.381e-02  1.195e-01   0.283  0.77729    
raw          4.281e-05  9.955e-06   4.301 1.72e-05 ***
ageB         7.397e-02  2.490e-02   2.971  0.00298 ** 
ageC         1.186e-01  2.285e-02   5.191 2.14e-07 ***
ageD         1.235e-01  2.254e-02   5.478 4.43e-08 ***
ageE         1.323e-01  2.305e-02   5.738 9.86e-09 ***
ageF         1.065e-01  2.405e-02   4.429 9.56e-06 ***
ageG         7.922e-02  2.625e-02   3.018  0.00255 ** 
ageH         7.100e-02  3.096e-02   2.293  0.02185 *  
ageI         7.207e-02  3.181e-02   2.266  0.02348 *  
ageJ        -2.037e-02  2.797e-02  -0.728  0.46651    
ageK         1.126e-01  3.926e-02   2.867  0.00415 ** 
areaB        8.270e-02  4.312e-02   1.918  0.05513 .  
areaC        4.041e-02  3.502e-02   1.154  0.24859    
areaD       -9.158e-03  3.682e-02  -0.249  0.80359    
areaE        6.240e-03  3.228e-02   0.193  0.84670    
areaF        1.429e-02  3.250e-02   0.440  0.66027    
areaG        2.310e-02  3.463e-02   0.667  0.50479    
areaH        1.376e-02  3.854e-02   0.357  0.72100    
aa           6.068e-07  5.277e-07   1.150  0.25023    
uu          -1.264e-03  2.936e-04  -4.306 1.68e-05 ***
dd          -3.612e-05  7.385e-06  -4.891 1.02e-06 ***
cc           7.126e-04  3.405e-04   2.093  0.03639 *  
bb           1.417e-02  1.976e-03   7.173 7.92e-13 ***
vv           2.259e-08  1.662e-08   1.360  0.17401    
ii          -8.650e-03  3.033e-03  -2.852  0.00436 ** 
pp           7.185e-07  1.653e-07   4.345 1.41e-05 ***
jj          -1.670e-06  1.089e-06  -1.534  0.12502    
FF          -2.082e-12  1.674e-12  -1.244  0.21351    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.4174 on 9235 degrees of freedom
Multiple R-squared:  0.3058,    Adjusted R-squared:  0.3033 
F-statistic: 123.3 on 33 and 9235 DF,  p-value: < 2.2e-16

r2.tr = summary(lm1)$r.sq
SST = sum((TS2$amount - mean(TR2$amount))^ 2)
SSE = sum((predict(lm1, TS2) -  TS2$amount)^2)
r2.ts = 1 - (SSE/SST)
c(r2.tr, r2.ts)

[1] 0.3057826 0.2775772

#0.2772139

加入變數

library(lubridate)
#顧客在1月的消費
Jan = filter(X, month(date)==1 ) %>% 
  group_by(cust) %>% 
  summarise(
    amount_m1 = sum(total),
    items_m1=sum(items),
    pieces_m1=sum(pieces),
    gross_m1=sum(gross),
    price_m1=sum(gross)
  ) 
Jan

加入變數到A2

A2 = merge(A2, Jan, by="cust", all.x=T)
A2

---

填補NA

##### 用平均值填補NA
for(i in 22:26){
  mean_col <- mean(A2[, i], na.rm = T)  # mean of col ith
  na.rows <- is.na(A2[, i])   #col ith na data
  A2[na.rows, i] <- mean_col
}
A2

A2新增變數加入Training Testing

TR2 = subset(A2, spl2)
TS2 = subset(A2, !spl2)

cx=c(2:10, 12:23)
colnames(TR2[,cx])

 [1] "r"         "s"         "f"         "m"         "rev"      
 [6] "raw"       "age"       "area"      "amount"    "aa"       
[11] "uu"        "dd"        "cc"        "bb"        "vv"       
[16] "ii"        "pp"        "jj"        "FF"        "amount_m1"
[21] "items_m1" 

lm1 = lm(amount ~ ., TR2[,cx])
summary(lm1)

Call:
lm(formula = amount ~ ., data = TR2[, cx])

Residuals:
     Min       1Q   Median       3Q      Max 
-1.85417 -0.22684  0.04689  0.27720  1.51675 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.799e+00  1.088e-01  16.530  < 2e-16 ***
r            5.168e-03  1.195e-03   4.325 1.54e-05 ***
s           -2.379e-02  3.155e-03  -7.542 5.08e-14 ***
f            1.625e-02  3.622e-02   0.449 0.653685    
m            2.622e-01  1.223e-01   2.144 0.032080 *  
rev          5.165e-02  1.199e-01   0.431 0.666572    
raw          4.835e-05  1.179e-05   4.099 4.18e-05 ***
ageB         7.350e-02  2.488e-02   2.955 0.003138 ** 
ageC         1.181e-01  2.283e-02   5.171 2.38e-07 ***
ageD         1.223e-01  2.253e-02   5.428 5.84e-08 ***
ageE         1.297e-01  2.304e-02   5.631 1.85e-08 ***
ageF         1.038e-01  2.404e-02   4.317 1.60e-05 ***
ageG         7.767e-02  2.623e-02   2.962 0.003068 ** 
ageH         6.980e-02  3.093e-02   2.256 0.024073 *  
ageI         7.282e-02  3.179e-02   2.291 0.021989 *  
ageJ        -2.116e-02  2.795e-02  -0.757 0.449003    
ageK         1.127e-01  3.923e-02   2.873 0.004072 ** 
areaB        8.299e-02  4.308e-02   1.927 0.054070 .  
areaC        3.967e-02  3.499e-02   1.134 0.257031    
areaD       -8.264e-03  3.679e-02  -0.225 0.822295    
areaE        4.430e-03  3.225e-02   0.137 0.890744    
areaF        1.264e-02  3.247e-02   0.389 0.697022    
areaG        2.361e-02  3.460e-02   0.682 0.494948    
areaH        1.404e-02  3.851e-02   0.365 0.715413    
aa           5.373e-07  5.281e-07   1.018 0.308913    
uu          -1.225e-03  2.957e-04  -4.144 3.45e-05 ***
dd          -3.543e-05  7.381e-06  -4.800 1.61e-06 ***
cc           7.364e-04  3.409e-04   2.161 0.030759 *  
bb           1.366e-02  2.033e-03   6.721 1.92e-11 ***
vv           1.286e-08  1.676e-08   0.767 0.442955    
ii          -7.657e-03  3.044e-03  -2.515 0.011919 *  
pp           7.111e-07  1.655e-07   4.296 1.76e-05 ***
jj          -1.459e-06  1.115e-06  -1.309 0.190518    
FF          -2.299e-12  1.673e-12  -1.374 0.169327    
amount_m1   -1.679e-05  4.897e-06  -3.428 0.000611 ***
items_m1     2.571e-03  6.205e-04   4.143 3.46e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.417 on 9233 degrees of freedom
Multiple R-squared:  0.3072,    Adjusted R-squared:  0.3045 
F-statistic:   117 on 35 and 9233 DF,  p-value: < 2.2e-16

r2.tr = summary(lm1)$r.sq
SST = sum((TS2$amount - mean(TR2$amount))^ 2)
SSE = sum((predict(lm1, TS2) -  TS2$amount)^2)
r2.ts = 1 - (SSE/SST)
c(r2.tr, r2.ts)

[1] 0.3071653 0.2792819

#0.2793

2.嘗試先集群再LM

LTR = TR2[,c(2,4,5)]   ##r,f,m
LTS = TS2[,c(2,4,5)]
library(caret)
preproc = preProcess(LTR)
NTR = predict(preproc, LTR) ##標準化資料
NTS = predict(preproc, LTS)

km <- kmeans(NTR,5)
library(flexclust)
km.kcca = as.kcca(km,NTR)  ##拿TR2做分群
CTR = predict(km.kcca)     
CTS = predict(km.kcca, newdata=NTS)   ##預測TS2群

拿分群做回歸

apple = split(TR2, CTR)
M = lapply(1:5, function(x) 
  lm(amount ~ ., data = apple[[x]][,c(2:7,10,12:23)]))

預測log(amount)

Pred = lapply(1:5, function(i) 
  predict(M[[i]], TS2[CTS==i,]) )

prediction from a rank-deficient fit may be misleading

t = do.call(c, split(TS2$amount,CTS))
y = do.call(c, Pred)

計算R^2

SST = sum((TS2$amount - mean(TR2$amount))^ 2)
SSE = sum((y - t)^2)
r2.ts = 1 - (SSE/SST)
r2.ts

[1] 0.2616883

cor(TR2[,c(2:7,10,12:23)])

                    r           s           f            m         rev
r          1.00000000  0.09280689 -0.40129374 -0.081126640 -0.44027271
s          0.09280689  1.00000000  0.42269893 -0.038503261  0.37382516
f         -0.40129374  0.42269893  1.00000000 -0.052591356  0.53746879
m         -0.08112664 -0.03850326 -0.05259136  1.000000000  0.73023789
rev       -0.44027271  0.37382516  0.53746879  0.730237893  1.00000000
raw       -0.29568266  0.27730825  0.58067392  0.444250605  0.69106801
amount    -0.15318868  0.11573180  0.26075411  0.448832752  0.49819086
aa         0.73885222  0.23988889 -0.15141996  0.129895002 -0.01798996
uu         0.45421111  0.41233382  0.20091882 -0.093705354  0.10175458
dd         0.94090159  0.13572184 -0.32432534  0.055867743 -0.29336220
cc        -0.35813848  0.48208952  0.99335208 -0.053482978  0.52699089
bb         0.07156626  0.97702793  0.40478867  0.160789068  0.51390410
vv        -0.08120938  0.08234599  0.56792007 -0.008137503  0.16138921
ii        -0.29750436  0.31169698  0.96107444 -0.049063048  0.43032686
pp         0.89744937  0.16148585 -0.28441875 -0.076304783 -0.36618575
jj        -0.61988527  0.61834594  0.58949999  0.220160172  0.68107758
FF        -0.35556936  0.39572205  0.95126182 -0.042927934  0.45918225
amount_m1 -0.07099274  0.13250769  0.40762238  0.375026779  0.47183446
items_m1  -0.07773565  0.16806445  0.50023232  0.271614719  0.43343367
                  raw       amount          aa           uu          dd
r         -0.29568266 -0.153188684  0.73885222  0.454211113  0.94090159
s          0.27730825  0.115731805  0.23988889  0.412333822  0.13572184
f          0.58067392  0.260754113 -0.15141996  0.200918822 -0.32432534
m          0.44425061  0.448832752  0.12989500 -0.093705354  0.05586774
rev        0.69106801  0.498190855 -0.01798996  0.101754584 -0.29336220
raw        1.00000000  0.421625952 -0.06731293  0.063309711 -0.21378675
amount     0.42162595  1.000000000 -0.00180071 -0.009755721 -0.06830873
aa        -0.06731293 -0.001800710  1.00000000  0.753791502  0.75091862
uu         0.06330971 -0.009755721  0.75379150  1.000000000  0.35916408
dd        -0.21378675 -0.068308730  0.75091862  0.359164076  1.00000000
cc         0.57852509  0.259085976 -0.11901949  0.224913468 -0.28322191
bb         0.37057984  0.208395242  0.26543145  0.384065368  0.14421903
vv         0.36043920  0.132637142 -0.04356514  0.036113917 -0.05309565
ii         0.54955750  0.241391578 -0.12293414  0.155115137 -0.22846113
pp        -0.21288031 -0.108296885  0.64532873  0.303325579  0.95401175
jj         0.59298064  0.305639127 -0.37252335 -0.119131389 -0.51814039
FF         0.55855316  0.253255869 -0.19688190  0.066773206 -0.26928466
amount_m1  0.74927260  0.319771668 -0.01261175  0.007903888 -0.01865575
items_m1   0.67853640  0.311077080 -0.01669654  0.027278846 -0.02135074
                   cc         bb           vv          ii           pp
r         -0.35813848 0.07156626 -0.081209378 -0.29750436  0.897449375
s          0.48208952 0.97702793  0.082345995  0.31169698  0.161485854
f          0.99335208 0.40478867  0.567920071  0.96107444 -0.284418746
m         -0.05348298 0.16078907 -0.008137503 -0.04906305 -0.076304783
rev        0.52699089 0.51390410  0.161389208  0.43032686 -0.366185746
raw        0.57852509 0.37057984  0.360439203  0.54955750 -0.212880308
amount     0.25908598 0.20839524  0.132637142  0.24139158 -0.108296885
aa        -0.11901949 0.26543145 -0.043565137 -0.12293414  0.645328733
uu         0.22491347 0.38406537  0.036113917  0.15511514  0.303325579
dd        -0.28322191 0.14421903 -0.053095651 -0.22846113  0.954011750
cc         1.00000000 0.46279346  0.569576385  0.95804383 -0.244884522
bb         0.46279346 1.00000000  0.079370124  0.29576266  0.137911969
vv         0.56957638 0.07937012  1.000000000  0.73453682 -0.042572013
ii         0.95804383 0.29576266  0.734536822  1.00000000 -0.195318561
pp        -0.24488452 0.13791197 -0.042572013 -0.19531856  1.000000000
jj         0.61154559 0.66406621  0.142146860  0.45316346 -0.460465063
FF         0.96745798 0.37985761  0.617988624  0.95406941 -0.226752835
amount_m1  0.40196061 0.20261448  0.339104065  0.42714038 -0.008194315
items_m1   0.49531676 0.21831856  0.451695695  0.53121837 -0.008228728
                  jj          FF    amount_m1     items_m1
r         -0.6198853 -0.35556936 -0.070992738 -0.077735651
s          0.6183459  0.39572205  0.132507692  0.168064445
f          0.5895000  0.95126182  0.407622380  0.500232322
m          0.2201602 -0.04292793  0.375026779  0.271614719
rev        0.6810776  0.45918225  0.471834464  0.433433673
raw        0.5929806  0.55855316  0.749272604  0.678536397
amount     0.3056391  0.25325587  0.319771668  0.311077080
aa        -0.3725233 -0.19688190 -0.012611751 -0.016696538
uu        -0.1191314  0.06677321  0.007903888  0.027278846
dd        -0.5181404 -0.26928466 -0.018655749 -0.021350742
cc         0.6115456  0.96745798  0.401960609  0.495316755
bb         0.6640662  0.37985761  0.202614476  0.218318564
vv         0.1421469  0.61798862  0.339104065  0.451695695
ii         0.4531635  0.95406941  0.427140381  0.531218371
pp        -0.4604651 -0.22675284 -0.008194315 -0.008228728
jj         1.0000000  0.58762774  0.333491636  0.327896928
FF         0.5876277  1.00000000  0.411214247  0.507448063
amount_m1  0.3334916  0.41121425  1.000000000  0.832241010
items_m1   0.3278969  0.50744806  0.832241010  1.000000000

apple = split(TR2, CTR)
M = lapply(1:5, function(x) 
  lm(amount ~ ., data = apple[[x]][,c(2:6,10,13,14,17,18,20)]))

Pred = lapply(1:5, function(i) 
  predict(M[[i]], TS2[CTS==i,]) )

t = do.call(c, split(TS2$amount,CTS))
y = do.call(c, Pred)

SST = sum((TS2$amount - mean(TR2$amount))^ 2)
SSE = sum((y - t)^2)
r2.ts = 1 - (SSE/SST)
r2.ts

[1] 0.2638244

交叉驗證~~

ctrl = trainControl(
  method="repeatedcv", number=10,    # 10-fold, Repeated CV
  savePredictions = "final", classProbs=TRUE,
  summaryFunction=twoClassSummary)

ctrl2 = trainControl(
  method="repeatedcv", number=10,    # 10-fold, Repeated CV
  savePredictions = "final")

ctrl$repeats = 2
set.seed(2)
cv.lm2 = train(
  amount ~ ., data=TR2[,c(2:10, 12:23)], method="lm", 
  trControl=ctrl2, metric="Rsquared",
    tuneGrid = expand.grid( intercept = seq(1,2,0.01) ) 
  )
plot(cv.lm2)

cv.lm2$results

summary(cv.lm2)

Call:
lm(formula = .outcome ~ ., data = dat)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.85417 -0.22684  0.04689  0.27720  1.51675 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.799e+00  1.088e-01  16.530  < 2e-16 ***
r            5.168e-03  1.195e-03   4.325 1.54e-05 ***
s           -2.379e-02  3.155e-03  -7.542 5.08e-14 ***
f            1.625e-02  3.622e-02   0.449 0.653685    
m            2.622e-01  1.223e-01   2.144 0.032080 *  
rev          5.165e-02  1.199e-01   0.431 0.666572    
raw          4.835e-05  1.179e-05   4.099 4.18e-05 ***
ageB         7.350e-02  2.488e-02   2.955 0.003138 ** 
ageC         1.181e-01  2.283e-02   5.171 2.38e-07 ***
ageD         1.223e-01  2.253e-02   5.428 5.84e-08 ***
ageE         1.297e-01  2.304e-02   5.631 1.85e-08 ***
ageF         1.038e-01  2.404e-02   4.317 1.60e-05 ***
ageG         7.767e-02  2.623e-02   2.962 0.003068 ** 
ageH         6.980e-02  3.093e-02   2.256 0.024073 *  
ageI         7.282e-02  3.179e-02   2.291 0.021989 *  
ageJ        -2.116e-02  2.795e-02  -0.757 0.449003    
ageK         1.127e-01  3.923e-02   2.873 0.004072 ** 
areaB        8.299e-02  4.308e-02   1.927 0.054070 .  
areaC        3.967e-02  3.499e-02   1.134 0.257031    
areaD       -8.264e-03  3.679e-02  -0.225 0.822295    
areaE        4.430e-03  3.225e-02   0.137 0.890744    
areaF        1.264e-02  3.247e-02   0.389 0.697022    
areaG        2.361e-02  3.460e-02   0.682 0.494948    
areaH        1.404e-02  3.851e-02   0.365 0.715413    
aa           5.373e-07  5.281e-07   1.018 0.308913    
uu          -1.225e-03  2.957e-04  -4.144 3.45e-05 ***
dd          -3.543e-05  7.381e-06  -4.800 1.61e-06 ***
cc           7.364e-04  3.409e-04   2.161 0.030759 *  
bb           1.366e-02  2.033e-03   6.721 1.92e-11 ***
vv           1.286e-08  1.676e-08   0.767 0.442955    
ii          -7.657e-03  3.044e-03  -2.515 0.011919 *  
pp           7.111e-07  1.655e-07   4.296 1.76e-05 ***
jj          -1.459e-06  1.115e-06  -1.309 0.190518    
FF          -2.299e-12  1.673e-12  -1.374 0.169327    
amount_m1   -1.679e-05  4.897e-06  -3.428 0.000611 ***
items_m1     2.571e-03  6.205e-04   4.143 3.46e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.417 on 9233 degrees of freedom
Multiple R-squared:  0.3072,    Adjusted R-squared:  0.3045 
F-statistic:   117 on 35 and 9233 DF,  p-value: < 2.2e-16

結論

### 一般直接線性有時比群集再線性更有效
### 透過交叉驗證可以說明我們資料不只對TEST DATA 有效 ~~~

---