1 Описание проблемы

Один из казахстанских банков предлагает задачу для соискателей на должности аналитиков-программистов. Во второй задаче имеются данные о клиентах компании “Рога и Копыта”. Набор по клиентам содержит 30 000 обучающие наблюдений по 132 характеристикам: порядковый номер, 130 интервальных числовых признаков и 1 номинальный признак, а также искомый бинарный признак TARGET - является ли клиент мошенником (“1” - мошенник / “0”" – обычный клиент). Также в отдельном листе книги MS Excel содержится 10 000 наблюдений с незаполненным полем TARGET. Задача - классифицировать тестовый набор и проставить в ней аналогичный флаг мошенника.

Какого-либо описания признаков нет, но F120 содержит указание на образование клиента компании.

Загрузим необходимые библиотеки R и данные из книги MS Excel.

## Load libraries
library('tidyverse')

## -- Attaching packages --------------------------------------------------------------------------------- tidyverse 1.2.1 --

## v ggplot2 2.2.1     v purrr   0.2.4
## v tibble  1.4.2     v dplyr   0.7.4
## v tidyr   0.7.2     v stringr 1.2.0
## v readr   1.1.1     v forcats 0.2.0

## -- Conflicts ------------------------------------------------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library('caret')

## Loading required package: lattice

## 
## Attaching package: 'caret'

## The following object is masked from 'package:purrr':
## 
##     lift

library('caretEnsemble')

## 
## Attaching package: 'caretEnsemble'

## The following object is masked from 'package:ggplot2':
## 
##     autoplot

##  A data frame containing 130 features: 129 numeric and 1 ordered predictors & dichatomic nominal outcome.
DF1 <- readxl::read_excel("C:/Soft/R/Examples/caret/task2.xlsx", col_names=TRUE, sheet = "Мошенники_train")
DF2 <- readxl::read_excel("C:/Soft/R/Examples/caret/task2.xlsx", col_names=TRUE, sheet = "Мошенники_test")

DF <- bind_rows(DF1, DF2)
DF$F120 <- as.factor(DF$F120) # Transformation character variable into factor
str(DF)

## Classes 'tbl_df', 'tbl' and 'data.frame':    40000 obs. of  132 variables:
##  $ NUM   : num  0 1 2 3 4 5 6 7 8 9 ...
##  $ F0    : num  0 0 0 2 0 0 0 0 0 1 ...
##  $ F1    : num  1 0 0 10 5 ...
##  $ F2    : num  1 0 0 1 4 ...
##  $ F3    : num  199900 0 0 1187750 318771 ...
##  $ F4    : num  0 0 0 130990 0 ...
##  $ F5    : num  1 0 0 6 4 ...
##  $ F6    : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ F7    : num  1 0 0 2 1 ...
##  $ F8    : num  0 0 2 0 0 ...
##  $ F9    : num  0 0 2 0 0 ...
##  $ F10   : num  3 0 0 0 0 0 0 0 0 0 ...
##  $ F11   : num  2 0 0 0 0 0 0 0 0 0 ...
##  $ F12   : num  49725 2121 2121 2121 57831 ...
##  $ F13   : num  5 0 0 0 5 ...
##  $ F14   : num  1 0 0 0 1 ...
##  $ F15   : num  120 NA NA NA 179 ...
##  $ F16   : num  0 0 0 0 0 ...
##  $ F17   : num  1 0 0 0 1 ...
##  $ F18   : num  0 1 1 1 0 0 1 0 0 1 ...
##  $ F19   : num  0 1 1 1 0 0 1 1 0 1 ...
##  $ F20   : num  0 17 4 5 6 ...
##  $ F21   : num  0 1 4 5 6 ...
##  $ F22   : num  276677 142918 276677 276677 156218 ...
##  $ F23   : num  0 0 0 0 0 4 0 0 0 0 ...
##  $ F24   : num  0 3 0 0 1 ...
##  $ F25   : num  1 3 1 6 6 ...
##  $ F26   : num  34 111 5 6 10 ...
##  $ F27   : num  5 5 4.43 4.43 4.43 ...
##  $ F28   : num  1 0 0 0 0 ...
##  $ F29   : num  660000 125000 120005 144500 125000 ...
##  $ F30   : num  34 32 23 21 24 ...
##  $ F31   : num  6 876 499 88 32 ...
##  $ F32   : num  0 5 0 1 0 ...
##  $ F33   : num  0 0 0 0 0 ...
##  $ F34   : num  0 3 0 1 0 ...
##  $ F35   : num  0 3 0 0 0 ...
##  $ F36   : num  1 19 0 45 1 ...
##  $ F37   : num  1 0 0 29 7 ...
##  $ F38   : num  220275 92879 -2120 62879 27169 ...
##  $ F39   : num  610275 122879 117884 142379 67169 ...
##  $ F40   : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ F41   : num  9 0 0 12 1 ...
##  $ F42   : num  1 0 0 6 1 ...
##  $ F43   : num  0 4 1 2 0 ...
##  $ F44   : num  0 100 33.3 100 0 ...
##  $ F45   : num  100 0 0 85.7 100 ...
##  $ F46   : num  2.46 42616 42491 2.46 2.46 ...
##  $ F47   : num  199900 0 0 198990 154000 ...
##  $ F48   : num  0 0 2 0 0 ...
##  $ F49   : num  0 0 2 0 0 ...
##  $ F50   : num  0e+00 0e+00 5e+04 0e+00 0e+00 0e+00 8e+04 0e+00 0e+00 4e+05 ...
##  $ F51   : num  1 0 0 4.01 42406 ...
##  $ F52   : num  199900 0 0 65135 67756 ...
##  $ F53   : num  12 6 9 42500 22 ...
##  $ F54   : num  1 0 0 2 5 ...
##  $ F55   : num  0 2 0 1 0 ...
##  $ F56   : num  0 0 0 20 5 ...
##  $ F57   : num  1 0 0 1 0 0 0 0 0 24 ...
##  $ F58   : num  1.56 5.84 9.67e+05 6.58 7.18 ...
##  $ F59   : num  5 1977 160 182 238 ...
##  $ F60   : num  0e+00 0e+00 3e+04 0e+00 0e+00 0e+00 8e+04 0e+00 0e+00 4e+05 ...
##  $ F61   : num  199900 0 0 23990 31200 ...
##  $ F62   : num  1 1 1 1 2 ...
##  $ F63   : num  0 83.3 0 16.7 0 ...
##  $ F64   : num  1 0 0 6 4 ...
##  $ F65   : num  1 0 0 9 11 ...
##  $ F66   : num  12.7 16 6 2 2 ...
##  $ F67   : num  0 1 0 0 0 ...
##  $ F68   : num  92020 92020 112700 73100 48000 ...
##  $ F69   : num  0 1 0 0 0 0 0 0 0 0 ...
##  $ F70   : num  0 0 1 0 0 ...
##  $ F71   : num  0 1 1 1 0 ...
##  $ F72   : num  1 1 1 3 0 ...
##  $ F73   : num  1 0 0 0 0 1 0 1 0 0 ...
##  $ F74   : num  0 0 0 4 0 ...
##  $ F75   : num  100 NA 0 100 100 ...
##  $ F76   : num  1 1 0 1 0 1 0 1 1 0 ...
##  $ F77   : num  0 1 0 3 0 1 0 1 1 0 ...
##  $ F78   : num  0 1 0 0 0 ...
##  $ F79   : num  0 1 0 0 0 ...
##  $ F80   : num  0 1 0 0 0 0 0 0 0 0 ...
##  $ F81   : num  0 0 1 1 0 ...
##  $ F82   : num  0 1 1 11 1 ...
##  $ F83   : num  1 1 1 11 1 ...
##  $ F84   : num  0 0 0 12 0 ...
##  $ F85   : num  100 64147 0 130 100 ...
##  $ F86   : num  1 1 0 1 2 ...
##  $ F87   : num  0 1 0 7 0 ...
##  $ F88   : num  0 1 0 0 0 1 0 0 0 0 ...
##  $ F89   : num  0 1 0 0 0 1 0 0 0 0 ...
##  $ F90   : num  0 1 0 0 0 0 0 0 0 0 ...
##  $ F91   : num  0 0 0 0 0 ...
##  $ F92   : num  0 1 1 1 0 ...
##  $ F93   : num  0 1 1 3 0 ...
##  $ F94   : num  0 0 0 4 0 ...
##  $ F95   : num  1 1 1 1 0 1 1 1 1 0 ...
##  $ F96   : num  0 1 0 3 0 1 0 1 1 0 ...
##  $ F97   : num  0 1 1 0 0 ...
##   [list output truncated]

## To list all columns with missing values you might write: 
is.na(DF) %>% colSums()

##    NUM     F0     F1     F2     F3     F4     F5     F6     F7     F8 
##      0      0      0      0      0      0      0      0      0      0 
##     F9    F10    F11    F12    F13    F14    F15    F16    F17    F18 
##      0      0      0   6873      0      0  20822      0      0      0 
##    F19    F20    F21    F22    F23    F24    F25    F26    F27    F28 
##      0      0      0      0      0      0      0      0      0      0 
##    F29    F30    F31    F32    F33    F34    F35    F36    F37    F38 
##      0      0      0      0      0      0      0      0      0      0 
##    F39    F40    F41    F42    F43    F44    F45    F46    F47    F48 
##      0      0      0      0      0      0      0      0      0      0 
##    F49    F50    F51    F52    F53    F54    F55    F56    F57    F58 
##      0      0      0      0      0      0      0      0      0      0 
##    F59    F60    F61    F62    F63    F64    F65    F66    F67    F68 
##      0      0      0      0      0      0      0      0      0      0 
##    F69    F70    F71    F72    F73    F74    F75    F76    F77    F78 
##      0      0      0      0      0      0  17648      0      0      0 
##    F79    F80    F81    F82    F83    F84    F85    F86    F87    F88 
##      0      0      0      0      0      0      0      0      0      0 
##    F89    F90    F91    F92    F93    F94    F95    F96    F97    F98 
##      0      0      0      0      0      0      0      0      0      0 
##    F99   F100   F101   F102   F103   F104   F105   F106   F107   F108 
##      0      0      0      0      0      0      0      0      0      0 
##   F109   F110   F111   F112   F113   F114   F115   F116   F117   F118 
##      0      0      0      0      0      0      0      0      0      0 
##   F119   F120   F121   F122   F123   F124   F125   F126   F127   F128 
##      0      0      0      0      0      0      0      0      0      0 
##   F129 TARGET 
##      0  10000

X <- select(DF, -one_of(c("NUM", "F12", "F15", "F75", "TARGET")))
Y <- factor( ifelse( DF[, "TARGET"] == 1, 'Yes', 'No') )

## Split the data into two sets: training (inTrain) and testing (inTest) 
seed <- 1991
set.seed(seed)

inTrain1 <- runif(nrow(DF1)) <= 9/10
inTrain = c( inTrain1, rep( FALSE, times = dim(DF2)[1] ) )
inTest = c( !inTrain1, rep( FALSE, times = dim(DF2)[1] ) )
inProblem = c( rep( FALSE, times = dim(DF1)[1] ), rep( TRUE, times = dim(DF2)[1] ) )

table(Y[inTrain]) ## significant Class Imbalance Problem

## 
##    No   Yes 
## 17420  9580

# Weights of cases to resolve Class Imbalances Problem
weight.cases <- as.numeric( Y[inTrain]  )
for(val in unique(weight.cases)) {weight.cases[weight.cases==val]=
    1/sum(weight.cases==val)*length(weight.cases)/2} # normalized to sum to length(samples)

Посмотрим на распределение классов в целевом факторе. Как видим, классы не сбалансированы. Поэтому необходимо при расчетам прибегать к паре способов, которые исправляют проблему разбалансирования классов (англ. Imbalanced Classification Problem). Поскольку количество мошенников значительного меньше числа обычные клиентов. В этой задаче мы будем в методах классификации допускающих применение весов наблюдений прибегать к этому способу, а для прочих алгоритмов будем задавать сэмплинг на увеличение количества наблюдений малочисленного класса до уровня наибольшего класса.

При построении любой модели важно избежать переобучения (англ. "overfitted"), при котором переусложненная модель дает хорошие результаты подгонки на обучающем наборе данных, но часто не способна корректно предсказывать, выйдя за рамки обучения, что не преминет снизить аккуратность в дальнейшем. Построение модели сбалансированной сложности, в которой точность уравновешена с устойчивостью предсказанных значений, является важнейшей задачей развития алгоритмов Машинного обучения.

Поэтому рекомендуется использовать контрольный набор данных (англ. "Testing Set"). Это образец данных, которые мы изымаем при построении модели. Однако мы используем его непосредственно в конце нашего анализа, чтобы подтвердить точность нашей финальной модели. Это тест, который позволит нам понять, испорчены ли наша модель и дать нам уверенность в наших оценках точности (англ. accuracy) по неприменявшемся для её построения данным.

Обучающий набор (inTrain) составит 90% полных исходных данных. Кроме этого мы оставляем собственно проблемный набор данных (inProblem), для которые необходимо предсказать значения TARGET.

2 Преобразование данных и Разведочный Анализ

2.1 Преобразование данных

Заполняем пропущенные значения медианами, приводим исходные интервальные признаки к одинаковым единицам измерения посредством стандартизации.
Номинальные признаки были развернуты в дихотомические, а затем преобразованы в интервальные числовые признаки для согласованности со всеми другими числовыми характеристиками. При этом сами номинальные признаки удаляются.
После трансформации или преобразовании признаков могут появится предикторы с единственным уникальным значением или малым числом таких значений, происходящими с очень низкой частотой. Для многих алгоритмов это может привести к нестабильной работе. Эти «почти-нулевой-дисперсии» (англ. "near-zero-variance") характеристики необходимо выявить и удалить до моделирования.
Во избежании мультиколлинеарности убираем все признаки, которые имеют между собой уровень корреляции выше 90%.

##
##       Filtering predictors - Remove Redundant Variables
##

## 1. Centering, Scaling and Transformation, for example "YeoJohnson"
X <- cbind(X %>% select_if(., is.numeric) %>%
             predict(preProcess( ., method = c( "medianImpute", "center", "scale", "YeoJohnson" ) ), newdata = . ), #
             #  %>% spatialSign %>% data.frame

## 2. Convert factor (nominal and ordered) variables into a full set of dummy integer variables without linear dependencies induced between these predictors
             X %>% select_if(., is.factor) %>%
             data.frame(predict(dummyVars(" ~ .", data = ., fullRank = TRUE), newdata = .)) %>%
               select_if(., is.numeric)
      )

## 3. Dropping zero variance predictors: the cutoff for the ratio of the first most common value to the second value. See https://www.mql5.com/ru/articles/2029
nzv <- caret::nearZeroVar( X[inTrain, ], freqCut = 99/1 , saveMetrics= TRUE )
nzv[nzv$nzv,]

##              freqRatio percentUnique zeroVar  nzv
## F4            180.8034   8.922222222   FALSE TRUE
## F120.ПРОФИЛЬ 2249.0000   0.007407407   FALSE TRUE
## F120.СРЕДН    674.0000   0.007407407   FALSE TRUE
## F120.СРЕПЕНЬ    0.0000   0.003703704    TRUE TRUE
## F120.ШКОЛ     120.0762   0.007407407   FALSE TRUE

zv_cols = caret::nearZeroVar( X[inTrain, ], freqCut = 99/1, saveMetrics = FALSE )
print( sprintf("Dropping %d zero variance predictors from %d (fraction=%10.6f)",
        length(zv_cols), dim(X[inTrain, ])[2], length(zv_cols)/dim(X[inTrain, ])[2]) )

## [1] "Dropping 5 zero variance predictors from 134 (fraction=  0.037313)"

zv_cols

## [1]   5 129 131 132 134

if ( length(zv_cols) != 0 )  {
      X = X[, -zv_cols]
      }

# # 4. Remove NUMERIC variables with high correlation (> .90) to others (multicollinearity)
cor.matrix <- cor( sapply( X[inTrain, ], function(x)
                  { as.numeric(x) } ) )
cor.high   <- findCorrelation(cor.matrix, 0.90)

high.cor.remove <- row.names(cor.matrix)[cor.high]
print( sprintf("Dropping %d predictors due to high correlation to others (multicollinearity) %d (fraction=%10.6f)",
               length(high.cor.remove), dim(X[inTrain, ])[2],  length(high.cor.remove)/dim(X[inTrain, ])[2]) );

## [1] "Dropping 28 predictors due to high correlation to others (multicollinearity) 129 (fraction=  0.217054)"

high.cor.remove

##  [1] "F2"   "F3"   "F42"  "F48"  "F49"  "F54"  "F64"  "F65"  "F94"  "F1"  
## [11] "F5"   "F8"   "F10"  "F18"  "F32"  "F7"   "F37"  "F47"  "F50"  "F52" 
## [21] "F34"  "F78"  "F88"  "F41"  "F97"  "F84"  "F61"  "F106"

if (length(high.cor.remove) != 0) {
      X <- X[, -cor.high]
   }

## Remove some other features that do not add useful information for Machine Learning

# names( X )

В ходе преобразования после развертывания в дихотомические признаки количество характеристик удалось существенно сократить.

2.2 Описательная статистика

Первый шаг в этом процессе лучше понять проблему, ведь после таких преобразований получено 101 характеристик. Кроме размерности набора данных, т.е количества столбцов и строк в нем, полезно оценить основные свойства признаков, которые отобраны предикторами (минимум, максимум, среднее значение и квартили).

# dimensions of dataset
dim(X[inTrain, ])

## [1] 27000   101

# take a peek at the first 10 rows of the data
head(X[inTrain, ], n = 10)

##            F0         F6          F9        F11         F13        F14
## 1  -0.3750648 -0.2284258 -0.31831607  2.7204685  0.96521891  0.4504196
## 2  -0.3750648 -0.2284258 -0.31831607 -0.2327712 -1.23094841 -1.2444212
## 3  -0.3750648 -0.2284258  1.89357213 -0.2327712 -1.23094841 -1.2444212
## 4   1.0364339 -0.2284258 -0.31831607 -0.2327712 -1.23094841 -1.2444212
## 5  -0.3750648 -0.2284258 -0.31831607 -0.2327712  0.96521891  0.4504196
## 6  -0.3750648 -0.2284258 -0.31831607 -0.2327712  0.34758968  0.2965521
## 7  -0.3750648 -0.2284258 -0.03739404 -0.2327712 -1.23094841 -1.2444212
## 8  -0.3750648 -0.2284258 -0.31831607 -0.2327712  0.02675067  0.4504196
## 9  -0.3750648 -0.2284258 -0.31831607 -0.2327712  0.59322308  0.6163794
## 10  0.3306846 -0.2284258  0.78762803 -0.2327712 -1.23094841 -1.2444212
##            F16        F17        F19        F20         F21        F22
## 1  -0.44816644  0.6937819 -0.8249551 -1.5781085 -1.43366633  0.3672956
## 2  -0.44816644 -1.0094700  1.2389567  0.3547613 -1.04804555 -0.4950095
## 3  -0.44816644 -1.0094700  1.2389567 -0.6495012 -0.26126211  0.3672956
## 4  -0.44816644 -1.0094700  1.2389567 -0.5227584 -0.05333639  0.3672956
## 5  -0.44816644  0.6937819 -0.8249551 -0.4116836  0.13896654 -0.3835732
## 6   0.01458063  0.5339933 -0.8249551  1.3323777  0.64957139  2.1262730
## 7  -0.44816644 -1.0094700  1.2389567 -1.5781085 -1.43366633  0.3672956
## 8   1.20225336  1.2958858  1.2389567 -0.6495012 -0.26126211  0.3672956
## 9   0.30265448  1.1088741 -0.8249551  0.5987956  0.61595194  3.1907504
## 10 -0.44816644 -1.0094700  1.2389567 -1.5781085 -1.43366633 -1.5675889
##           F23        F24         F25         F26       F27        F28
## 1  -0.2104201 -1.2509932 -1.06690021 -0.03570962 0.4600757  1.4137814
## 2  -0.2104201  0.7638708 -0.01034756  1.45000560 0.4600757 -0.7145697
## 3  -0.2104201 -1.2509932 -1.06690021 -1.35921382 0.2390391 -0.7145697
## 4  -0.2104201 -1.2509932  0.82508981 -1.27483259 0.2390391 -0.7145697
## 5  -0.2104201 -0.1615936  0.82508981 -0.99775246 0.2390391 -0.7145697
## 6   6.7366735  0.7638708  0.32464147  0.35469219 0.2390391  0.9899409
## 7  -0.2104201 -1.2509932 -1.06690021  0.35469219 0.2390391 -0.7145697
## 8  -0.2104201 -1.2509932 -0.44590476 -1.45299371 0.2390391 -0.7145697
## 9  -0.2104201  0.6764689  0.38692693  0.60828238 1.2263520  0.9075532
## 10 -0.2104201 -0.1615936 -0.44590476 -1.83933866 0.2390391 -0.7145697
##           F29        F30        F31        F33        F35        F36
## 1   1.6465634  0.5395362 -1.9293055 -0.4881577 -0.3019280 -0.4162380
## 2  -0.2389806  0.3272802  1.0974745 -0.4881577  1.9658229  0.9809022
## 3  -0.2713973 -1.1063806  0.2905676 -0.4881577 -0.3019280 -1.0572083
## 4  -0.1194256 -1.5980238 -1.1589754 -0.4881577 -0.3019280  1.2950885
## 5  -0.2389806 -0.8920836 -1.5762025 -0.4881577 -0.3019280 -0.4162380
## 6  -2.2482148  1.3577210  0.2905676 -0.4881577 -0.3019280 -1.0572083
## 7  -0.2384664 -1.2227624  0.2905676 -0.4881577 -0.3019280 -1.0572083
## 8   0.4822639 -1.5980238 -1.1179408 -0.4881577 -0.3019280  0.5341745
## 9  -0.2418754  0.4904756  0.1180242  1.9521585 -0.1185779  0.7577309
## 10  0.3572593 -0.1922343 -1.4502733 -0.4881577 -0.3019280  1.2104387
##            F38         F39        F40         F43        F44        F45
## 1   0.88222266  1.89375840 -0.1692653 -1.14710485 -1.2449802  1.2941496
## 2   0.10887116 -0.30403376 -0.1692653  1.01113866  0.9601963 -0.7940411
## 3  -0.50946207 -0.32751085 -0.1692653 -0.04542373  0.2879074 -0.7940411
## 4  -0.07752372 -0.21273723 -0.1692653  0.47029489  0.9601963  1.2630986
## 5  -0.30379407 -0.56856111 -0.1692653 -1.14710485 -1.2449802  1.2941496
## 6  -0.01375086 -0.05789479 -0.1692653 -0.04542373  0.9601963 -0.7940411
## 7   0.20025456 -0.30365703 -0.1692653 -1.14710485 -1.2449802 -0.7940411
## 8  -0.08588711  0.10739444 -0.1692653 -0.04542373  0.2879074  1.2941496
## 9  -0.02364675 -0.12441400 -0.1692653  1.51298943  0.9121865  1.2371206
## 10 -0.14028246  0.22808935 -0.1692653  0.47029489  0.9601963  1.2857218
##           F46         F51         F53        F55        F56        F57
## 1  -0.3901315  0.04298515 -0.26000841 -0.6625497 -0.6917855  1.0476847
## 2   1.7294415 -0.93499952 -0.65532895  1.3400524 -0.6917855 -0.6093206
## 3   1.7293630 -0.93499952 -0.42029015 -0.6625497 -0.6917855 -0.6093206
## 4  -0.3901315  0.89691629  1.82671172  0.9290168  1.8289660  1.0476847
## 5  -0.3901315  2.26648486  0.05504271 -0.6625497  1.6046620 -0.6093206
## 6   1.7293630 -0.93499952  0.36377677 -0.6625497 -0.6917855 -0.6093206
## 7  -0.4814283 -0.93499952 -0.47076775 -0.6625497 -0.6917855 -0.6093206
## 8  -0.2790204  0.35555688 -0.42029015  0.9290168  0.7886259 -0.6093206
## 9  -0.2621606 -0.11356681  0.13687632  1.8458717  1.0589170 -0.6093206
## 10 -0.5058924  1.02856371  1.82670593 -0.6625497 -0.6917855  1.9520750
##            F58        F59       F60        F62        F63        F66
## 1  -1.69032993 -1.6554893 -0.459881 -0.1999468 -0.6078729  0.5320560
## 2  -0.54893403  1.2729538 -0.459881 -0.1999468  1.7063139  0.7573221
## 3   1.99926138 -1.0350505  2.161576 -0.1999468 -0.6078729 -0.1455321
## 4  -0.45248757 -0.9797745 -0.459881 -0.1999468  1.5444048 -0.9280756
## 5  -0.38232091 -0.8508140 -0.459881  0.7936514 -0.6078729 -0.9280756
## 6   0.61292844  0.9408396 -0.459881  0.7936514 -0.6078729 -1.8385049
## 7  -0.30773776 -1.6155081  2.170144 -1.8029839 -0.6078729  0.5320560
## 8  -0.33511905 -0.6496110 -0.459881 -0.1999468  1.5444048  0.2066794
## 9  -0.08804534  1.5360651 -0.459881  0.9890872  1.1404621  1.3155136
## 10  0.45747615  0.4414223  2.178030 -0.1999468 -0.6078729 -0.4651792
##           F67         F68        F69        F70        F71        F72
## 1  -1.0124843  0.38292282 -0.4982955 -0.2484266 -1.1919466  0.9017929
## 2   0.4869147  0.38292282  2.0359734 -0.2484266  0.6210658  0.9017929
## 3  -1.0124843  0.63538838 -0.4982955  2.9322940  0.6210658  0.9017929
## 4  -1.0124843  0.10703120 -0.4982955 -0.2484266  0.6210658  1.8787681
## 5  -1.0124843 -0.36905265 -0.4982955 -0.2484266 -1.1919466 -0.9661258
## 6   1.0569691  0.38292282 -0.4982955 -0.2484266  0.6210658  0.9017929
## 7   0.4869147  0.38292282 -0.4982955 -0.2484266  0.6210658  0.2829302
## 8   0.4869147 -0.20871239 -0.4982955 -0.2484266  0.6210658  0.9017929
## 9   1.1473628 -0.15728237 -0.4982955  2.1607999  1.3434130  0.9017929
## 10 -1.0124843 -0.08913272 -0.4982955 -0.2484266 -1.1919466 -0.9661258
##           F73        F74        F76        F77         F79        F80
## 1   1.1070522 -0.2911677  0.8126094 -0.7693181 -0.87109453 -0.4265556
## 2  -0.5717034 -0.2911677  0.8126094  1.3198370  1.15498336  2.2453090
## 3  -0.5717034 -0.2911677 -1.0775210 -0.7693181 -0.87109453 -0.4265556
## 4  -0.5717034 11.1669701  0.8126094  1.7122802 -0.87109453 -0.4265556
## 5  -0.5717034 -0.2911677 -1.0775210 -0.7693181 -0.87109453 -0.4265556
## 6   1.1070522 -0.2911677  0.8126094  1.3198370  1.15498336 -0.4265556
## 7  -0.5717034 -0.2911677 -1.0775210 -0.7693181 -0.87109453 -0.4265556
## 8   1.1070522 -0.2911677  0.8126094  1.3198370 -0.87109453 -0.4265556
## 9  -0.5717034  1.8785647  0.8126094  1.3198370 -0.04266895 -0.4265556
## 10 -0.5717034 -0.2911677 -1.0775210 -0.7693181 -0.87109453 -0.4265556
##           F81        F82       F83        F85        F86        F87
## 1  -0.3248362 -1.1278916 0.7969679 -0.6813253  0.9147271 -0.6580809
## 2  -0.3248362  0.5681946 0.7969679  1.3706240  0.9147271  1.3019975
## 3   2.2854222  0.5681946 0.7969679 -4.3097633 -0.9577701 -0.6580809
## 4   2.2854222  2.4872841 2.1574481 -0.5488329  0.9147271 13.0624685
## 5  -0.3248362  0.5681946 0.7969679 -0.6813253  1.5366166 -0.6580809
## 6  -0.3248362  0.5681946 0.7969679  1.3706240  0.9147271  1.3019975
## 7  -0.3248362  0.5681946 0.2349388 -0.5375694 -0.9577701 -0.6580809
## 8  -0.3248362  1.6047051 1.6318923 -0.3429502  0.9147271  1.3019975
## 9   1.6522956  1.1119688 0.5605401 -0.3993642  0.6568204  0.8265742
## 10 -0.3248362 -1.1278916 0.7969679 -0.6813253 -0.9577701 -0.6580809
##           F89        F90        F91        F92        F93        F95
## 1  -0.6805641 -0.5624468 -0.1976263 -1.2419994 -0.9719803  0.6803102
## 2   1.4704452  1.8282091 -0.1976263  0.6234199  0.9427618  0.6803102
## 3  -0.6805641 -0.5624468 -0.1976263  0.6234199  0.9427618  0.6803102
## 4  -0.6805641 -0.5624468 -0.1976263  0.6234199  2.0192261  0.6803102
## 5  -0.6805641 -0.5624468 -0.1976263 -1.2419994 -0.9719803 -1.2101620
## 6   1.4704452 -0.5624468 -0.1976263  0.6234199  0.9427618  0.6803102
## 7  -0.6805641 -0.5624468 -0.1976263  0.6234199 -0.2222013  0.6803102
## 8  -0.6805641 -0.5624468 -0.1976263  0.6234199  0.9427618  0.6803102
## 9  -0.6805641 -0.5624468  2.6617653  1.4543727  0.9427618  0.6803102
## 10 -0.6805641 -0.5624468 -0.1976263 -1.2419994 -0.9719803 -1.2101620
##           F96        F98        F99       F100       F101       F102
## 1  -0.8173491 -1.0476443 -0.5209571 -0.2495719 -0.9451842 -0.3346865
## 2   1.2421769  0.9593085  1.8866760 -0.2495719  0.9169835 -0.3346865
## 3  -0.8173491  0.9593085 -0.5209571 -0.2495719  0.9169835  0.4299355
## 4   1.8028418 -1.0476443 -0.5209571  2.8946903  2.3162958  0.4299355
## 5  -0.8173491 -1.0476443 -0.5209571 -0.2495719 -0.9451842 -0.3346865
## 6   1.2421769  0.9593085 -0.5209571 -0.2495719  0.9169835  0.4299355
## 7  -0.8173491  0.9593085 -0.5209571 -0.2495719 -0.1756040 -0.3346865
## 8   1.2421769 -1.0476443 -0.5209571 -0.2495719  1.5023936  0.4299355
## 9   1.2421769 -0.4745416 -0.5209571  2.1320394  0.9169835 -0.3346865
## 10 -0.8173491 -1.0476443 -0.5209571 -0.2495719  0.9169835  0.4299355
##          F103       F104       F105       F107       F108       F109
## 1  -0.2590429 -0.9398693 -0.7671477 -0.9160352  1.1458719  1.2395271
## 2  -0.2590429  1.0530408  1.3189088  1.0961748 -0.9018194 -0.7497482
## 3  -0.2590429  1.0530408 -0.7671477  1.0961748 -0.9018194 -0.7497482
## 4  24.7654321 -0.8506739  1.7823868 -0.9160352 -0.9018194 -0.7497482
## 5  -0.2590429 -0.9398693 -0.7671477 -0.9160352  1.1458719  1.3767632
## 6  -0.2590429  1.0530408  1.3189088  1.0961748 -0.9018194 -0.7497482
## 7  -0.2590429  1.0530408 -0.7671477  1.0961748 -0.9018194 -0.7497482
## 8  -0.2590429 -0.7061216  1.3189088 -0.9160352  1.1458719  1.3767632
## 9   1.3205161  0.6825469  1.3189088 -0.1620377  0.8680676  1.2187020
## 10 -0.2590429 -0.9398693 -0.7671477 -0.9160352 -0.9018194 -0.7497482
##          F110       F111         F112       F113       F114       F115
## 1   0.4918618  0.7432707 -0.530292862 -1.2123811  0.5506399  0.6768136
## 2   0.4918618 -1.4122374  0.349523992  0.4895526 -1.8346367 -1.4979969
## 3   0.4918618  0.7432707 -1.201414501 -1.2123811  0.5506399  0.6768136
## 4   0.4918618  0.7432707  0.152434007  1.7258636  0.5506399  0.6768136
## 5   0.4918618 -0.4079186  0.006662063  1.5259050  0.5506399  0.6768136
## 6  -2.0571459  0.7432707  1.533796858 -1.2123811  0.5506399  0.6768136
## 7   0.4918618  0.7432707 -1.143036531 -1.2123811  0.5506399  0.6768136
## 8   0.4918618  0.7432707 -1.003756106  0.1457934  0.5506399  0.6768136
## 9   0.4918618 -0.6661775 -1.404050697 -0.4190670  0.5506399 -0.9704896
## 10  0.4918618  0.7432707 -1.928667415 -1.2123811  0.5506399 -1.4979969
##          F116       F117        F118       F119       F121       F122
## 1   0.8338207 -0.5189758  1.32967890 -0.7514207  1.0056885  2.3799892
## 2  -1.2155774 -0.5189758 -0.01080949 -0.7514207  1.0056885 -0.4369685
## 3  -1.2155774 -0.5189758 -1.82327976  0.9400330 -0.3399709 -0.4369685
## 4   0.8338207 -0.5189758 -0.98821752 -0.7514207 -0.3399709 -0.4369685
## 5   0.8338207 -0.5189758  0.16603839  0.9400330 -0.3399709 -0.4369685
## 6   0.8338207  1.9363993  0.44278853  0.9400330  1.0056885 -0.4369685
## 7  -1.2155774 -0.5189758  0.37368829  0.9400330 -0.3399709 -0.4369685
## 8  -1.2155774 -0.5189758  0.31415777  0.9400330 -0.3399709  2.3799892
## 9  -0.8044835 -0.5189758 -0.14961067  0.9400330  0.6575783 -0.4369685
## 10 -1.2155774 -0.5189758  0.24298042 -1.1998540  1.0056885 -0.4369685
##          F123       F124      F125       F126       F127      F128
## 1  -0.2694714 -0.1535941  2.340352  1.3958285  1.0352899 -0.226562
## 2  -0.2694714 -0.1535941 -0.448080  1.3958285 -0.9973472 -0.226562
## 3  -0.2694714  6.9301945 -0.448080 -0.7369079  1.0352899 -0.226562
## 4  -0.2694714 -0.1535941 -0.448080  1.3958285  1.0352899  4.619047
## 5  -0.2694714 -0.1535941 -0.448080 -0.7369079  1.0352899 -0.226562
## 6  -0.2694714 -0.1535941 -0.448080 -0.7369079 -0.9973472 -0.226562
## 7  -0.2694714 -0.1535941 -0.448080  1.3958285  1.0352899 -0.226562
## 8   3.8100812 -0.1535941 -0.448080 -0.7369079  1.0352899 -0.226562
## 9  -0.2694714 -0.1535941 -0.448080  1.2459157  0.6725517 -0.226562
## 10 -0.2694714 -0.1535941 -0.448080  1.3958285 -0.9973472 -0.226562
##          F129 F120.НЕОК.ВЫСШ F120.НЕТ F120.СПЕЦ F120.СТЕПЕНЬ
## 1  -0.2061257              0        0         0            1
## 2  -0.2061257              0        0         0            1
## 3  -0.2061257              0        0         1            0
## 4  -0.2061257              0        0         0            1
## 5  -0.2061257              0        0         0            1
## 6  -0.2061257              0        0         1            0
## 7  -0.2061257              0        0         1            0
## 8  -0.2061257              0        0         0            1
## 9  -0.2061257              0        0         0            0
## 10 -0.2061257              0        0         0            1

# summarize attribute distributions of the first 10 rows of the data
summary(X[inTrain, 1:10])

##        F0                  F6                  F9          
##  Min.   :-0.375065   Min.   :-0.228426   Min.   :-0.31832  
##  1st Qu.:-0.375065   1st Qu.:-0.228426   1st Qu.:-0.31832  
##  Median :-0.375065   Median :-0.228426   Median :-0.31832  
##  Mean   : 0.008101   Mean   :-0.000701   Mean   :-0.00325  
##  3rd Qu.:-0.375065   3rd Qu.:-0.228426   3rd Qu.:-0.31832  
##  Max.   :10.211176   Max.   :17.117248   Max.   :14.05896  
##       F11                 F13                 F14           
##  Min.   :-0.232771   Min.   :-1.230948   Min.   :-1.244421  
##  1st Qu.:-0.232771   1st Qu.:-1.230948   1st Qu.:-1.244421  
##  Median :-0.232771   Median : 0.347590   Median : 0.296552  
##  Mean   :-0.001695   Mean   : 0.002232   Mean   :-0.001015  
##  3rd Qu.:-0.232771   3rd Qu.: 0.965219   3rd Qu.: 0.450420  
##  Max.   :18.963287   Max.   : 4.729664   Max.   : 3.702115  
##       F16                 F17                 F19           
##  Min.   :-0.448166   Min.   :-1.009470   Min.   :-0.824955  
##  1st Qu.:-0.448166   1st Qu.:-1.009470   1st Qu.:-0.824955  
##  Median :-0.448166   Median : 0.533993   Median :-0.824955  
##  Mean   :-0.000319   Mean   : 0.002362   Mean   :-0.001062  
##  3rd Qu.: 0.014581   3rd Qu.: 0.693782   3rd Qu.: 1.238957  
##  Max.   :16.056032   Max.   : 2.272201   Max.   : 1.238957  
##       F20          
##  Min.   :-1.57811  
##  1st Qu.:-0.79806  
##  Median : 0.07306  
##  Mean   : 0.00498  
##  3rd Qu.: 0.83837  
##  Max.   : 2.08950

2.3 Correlation Matrix

Перед построением модели полезно изучить корреляционную матрицу отобранных предикторов, в которой они сгруппированы во взаимосвязанные блоки.

## Correlation Matrix Heatmap for Feartures Visualization
## http://www.sthda.com/english/wiki/ggplot2-quick-correlation-matrix-heatmap-r-software-and-data-visualization
library('reshape2')

## 
## Attaching package: 'reshape2'

## The following object is masked from 'package:tidyr':
## 
##     smiths

library('scales')

## 
## Attaching package: 'scales'

## The following object is masked from 'package:purrr':
## 
##     discard

## The following object is masked from 'package:readr':
## 
##     col_factor

## Get lower triangle of the correlation matrix
  get_lower_tri<-function(cormat){
    cormat[upper.tri(cormat)] <- NA
    return(cormat)
  }

## Get upper triangle of the correlation matrix
  get_upper_tri <- function(cormat){
    cormat[lower.tri(cormat)] <- NA
    return(cormat)
  }

## Use correlation between variables as distance
reorder_cormat <- function(cormat){
dd <- as.dist((1 - cormat) / 2)
hc <- hclust(dd)
cormat <-cormat[hc$order, hc$order]
}

cormat <- round(cor( cbind(Y = as.numeric(Y[inTrain]), X[inTrain, ]) ), 2)
## Reorder the correlation matrix
cormat <- reorder_cormat(cormat)
upper_tri <- get_upper_tri(cormat)
# Melt the correlation matrix
melted_cormat <- melt(upper_tri, na.rm = TRUE)

plot <- ggplot(melted_cormat, aes(x = Var1, y = Var2, fill = value)) +
          geom_tile() +
          scale_x_discrete() +
          scale_y_discrete() +
          coord_fixed() +
          theme(legend.title = element_text(size = 8),
                legend.position = "top", legend.direction = "horizontal",
                legend.key.width = unit(1.25, "cm"), legend.key.height = unit(0.25, "cm"),
                legend.spacing = unit(0, "cm"), # panel.margin = element_blank(), 
                axis.text.x = element_text(angle = -90, vjust = 0.5, hjust = 0), 
                axis.title.y = element_blank(), axis.title.x = element_blank(),
                plot.title = element_text(hjust = 1, size = 12)) +
          scale_fill_gradient2(high = "#2ecc71", low = "#d1338e", mid = "white",
                midpoint = 0, limit = c(-1.0, 1.0),
                name = "Pearson\nCorrelation", breaks = pretty_breaks(8))  +
          labs(title = sprintf("Correlation Matrix of Target Factor (Y) and %s Features", dim(X)[2]))

## Output Plot as PNG file
png(paste("C:/Soft/R/Examples/caret/Correlation.Matrix2", "png", sep="."), units = "px", width = 910, height = 970)
plot(plot)
garbage <- dev.off()

## Positive and negative predictors associated with `AVG_COST`
melted_cormat %>%  filter(Var1 == 'Y') %>% arrange(desc(value)) %>% slice(2:8)

## # A tibble: 7 x 3
##   Var1  Var2     value
##   <fct> <fct>    <dbl>
## 1 Y     F45      0.290
## 2 Y     F19      0.280
## 3 Y     F120.НЕТ 0.270
## 4 Y     F126     0.260
## 5 Y     F56      0.250
## 6 Y     F51      0.220
## 7 Y     F111     0.210

melted_cormat %>%  filter(Var1 == 'Y') %>% arrange(value) %>% slice(2:8)

## # A tibble: 7 x 3
##   Var1  Var2   value
##   <fct> <fct>  <dbl>
## 1 Y     F55   -0.300
## 2 Y     F46   -0.280
## 3 Y     F31   -0.270
## 4 Y     F63   -0.270
## 5 Y     F20   -0.260
## 6 Y     F96   -0.250
## 7 Y     F24   -0.250

На графике заметно, что предикторы и целевой признак довольно сильно взаимосвязаны между собой - основной цвет корреляционной матрицы не белый (т.е. нейтральный), но широко представлены оттенки малинового (отрицательные связи) и зеленого цвета (положительные связи). Это заметно отличается от корреляционной матрицы из первой задачи, где слабая мультиколлинеарность признаков. Обращает внимание, что значения Y - флаг мошенничества положительно коррелирует с отсутствием сведений об образовании F120.НЕТ и характеристиками F81, F125 F83, но демонстрирует отрицательную сязь с признаками F63, F31, F59, F34, а также F96. В связи с этим открываются большие перспективы для применения методов снижения размерности, например, метода главных компонент.

2.4 Расчет расстояния между Классами

Пакет caret для организации машинного обучения содержит функции генерации новых предикторов на основе расстояний до центроидов класса (подобно тому, как работает линейный дискриминантный анализ). Для каждого набора предикторов вычисляется матрица центроидов классов целевого фактора Y и ковариаций. Для новых образцов расстояние Махаланобиса до каждого из центроидов класса вычисляется и может использоваться как дополнительный предиктор. Это может быть полезно для нелинейных моделей, когда истинная граница решения фактически линейна.

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

На графиках мы видим, что конфигурация наблюдений в обучающей и контрольной выборке по обоим классам довольно схожа. Поэтому мы можем полагать, что эти выборки части одной генеральной совокупности.

В тех случаях, когда в классе больше предикторов, чем в выборках, функция classDist имеет аргументы получения главных компонент и поддерживает аргументы, которые позволяют анализировать основные компоненты в каждом классе, чтобы избежать проблем с сингулярными ковариационными матрицами.

Затем для генерации расстояний класса используется pred.classDist. По умолчанию расстояния регистрируются, но это можно изменить с помощью аргумента trans для pred.classDist.

3 Построение модели

3.1 Формирование объекта управления машинным обучением

Для запуска машинного обучения требуется объект trainControl, в котором прописываются все необходимые параметры, в частности, возможности генерации повторных выборок или численный ресемпрлинг (англ. resampling), т.е. построение из полной обучающей выборки собственно обучающей (англ. train set) и проверочных выборок (англ. validation set), способ проверки, например, однократная проверка, бутстреп либо множественная перекрестная проверка, порядок вычислений при проблеме несбалансированных классов и т.д.

Однако провести лишь однократную проверку моделей для полноценного отбора модели недостаточно. Зачастую отлично зарекомендовавшая модель на одном массиве данных и удачно проверенная на последующем, на других данных, уже при проявлении новых явлений, станет предсказывать заметно хуже, чем модели ранее не столь успешные на однократной валидации. Так проявляется переподгонка параметров метода.

Но если в первой задаче, где регрессионая задача решена путем Глубокого обучения нейронный сетей для устранения переобученности модели производился отсев некоторой части нейронных связей, то применяемые нами алгоритмы требует других подходов. Поскольку наблюдений несколько десятков тысяч, то для этой классификационной задачи применим множественную перекрестную проверку.

The Cross-Validation Procedure with 5 folds or experiments

Source: Kaggle

Мы проведем 5 экспериментов (англ. Fold) в два подхода (англ. repeats) перепирая все наблюдения полной обучающей выборки.

## Do 2 repeats of 5-Fold CV for the data Set
folds = 5 
occasions = 2
set.seed(seed)

## For the last model:
set.seed(seed)

control <- trainControl(method = "repeatedcv",      ## cross-validation resampling method with 'repeats' sets
                        number = folds,             ## the number of k-fold cross-validation
                        repeats = occasions,        ## the number of complete sets of folds to compute
                        verboseIter = FALSE,        ## for printing a training log
                        returnResamp = "final",     ## for Ensembling with package 'caretEnsemble'
                        returnData = TRUE,          ## for Ensembling with package 'caretEnsemble'
                        savePredictions ="final",   ## for Ensembling with package 'caretEnsemble'
                        classProbs = TRUE,          ## for computing for classification models in each resample
                        index = createMultiFolds(Y[inTrain], k = folds, times = occasions), # for multi-folds
                        sampling = NULL)            ##  sampling using weights

controlImbalance <- control
controlImbalance$sampling <- 'up' ## Weighting from a smaller class to a larger one to resolve Class Imbalances
metrum  <- "Kappa"

В качестве меры точности модели мы изберем коэффициент Cohen’s Kappa, так как в задаче наблюдается проблема несбалансированности классов. Коэффициент аккуратности из-за этой проблемы в данной задаче не применим.Коэффициент Kappa - мера того, насколько близко наблюдения, классифицированные входе машинного обучения, т.е. предсказанные метки классов, соответствуют истинному разбиению на классы, контролируя точность обученного классификатора. Более того, что показатель Kappa проливает свет на то, как обучился сам классификатор, коэффициент для одной модели напрямую сопоставим со коэффициентами для любой другой модели, применяемой для одной и той же задачи классификации, включая множественнeую классификацию.

3.2 Преимущества и недостатки Ассамблирования Моделей

Из теории известно, что одним из основных принципов моделирования сложных систем является принцип множественности моделей, сформулированный В. В. Налимовым (1910-1997) и заключающийся в возможности представления одной и той же системы множеством различных моделей в зависимости от целей исследования. В этой связи встает задача формирования ансамблей статистических моделей, т.е. ассамблирование.

3.2.1 Преимущества

• Ассамблирование является проверенным методом повышения точности модели и работает в большинстве случаев.

• Это подход дает выигрыш почти для всех методов машинного обучения.

• Ансамбль делает модель более надежной и стабильной, обеспечивая при этом приемлимую производительность в тестовых случаях.

• Вы можете использовать ансамбли для решения относительно как просты, так и сложных нелинейных проблем. Его можно применять как для объединения моделей, полученными одним методом (например, созданных отдельно для каждого класса), так сформированными разными методами.

3.2.2 Недостатки

• Ассамблирование уменьшает интерпретационную способность модели и усложняет получения с ее помощью критических сведений для бизнеса.

• Требует много времени и, следовательно, может быть не лучшей идеей для приложений реального времени.

• Выбор моделей для создания ансамбля - это искусство, которое действительно требуется осваивать.

3.3 Обзор методов ассамблирования

Существуют различные методы для построения ансамблевых моделей. Некоторые наиболее популярных алгоритмов перечислены ниже:

Простая средная
Средневзвешенный по весам отдельным методов
Голосование большинства
Усиление через пересечение результатов “слабых” алгоритмов машинного обучения (Boosting algorithms)
Усиление через объединение результатов “слабых” алгоритмов машинного обучения (Bagging algorithms)
Укладка результатов “слабых” алгоритмов через другой способ машинного обучения (Ensemble Stacking aka Blending)

Мы воспользумся несколькими классификационными методами бустинга, баггинга, а также создадим несколько вариантов стекинга, включая объединение всех созданные моделей. Для этих алгоритмов применяли гиперпараметры обучения по умолчания и не оптимизировали специально под данную задачу.

Для комбинирования предсказаний моделей, созданных при помощи пакета caret, мы применим пакет ассамблирования caretEnsemble[last month’s downloads]http://cranlogs.r-pkg.org/badges/caretEnsemble. Он способен работать с регрессионными задачами, а также с бинарной (двухклассной) задачей. Для решения множественной классификации можно рекомендовать пакет машинного обучения mlr, способного провдить кластерный анализ и анализ выживаемости.

3.4 Алгоритмы Бустинга

Построим классификаторы двумя самыми удачными методами бустинга:

• Boosted Tree Model

• Stochastic Gradient Boosting

Ниже приведен пример сравнения моделей Boosted Tree и Gradient Boosting Model. Каждый из них допускает использование весов наблюдений для решения проблемы несбалансированных классов.

## Example of Boosting Algorithms
start_time <- Sys.time()
set.seed(seed)
boost.models <- caretList( x = X[inTrain, ], y = Y[inTrain], trControl = control, 
    tuneList=list(
## Boosted Tree Model with 2 Parameters
      {
        writeLines("\n\n`Boosted Tree Model` ...\n")
        set.seed(seed)
        caretModelSpec(method = "blackboost", metric = metrum, weights = weight.cases)
      },
## Stochastic Gradient Boosting with 4 Parameters
      {
        writeLines("\n\nGradient Boosting Model ...\n")    
        set.seed(seed)
        caretModelSpec(method = "gbm",  metric = metrum, verbose = FALSE, weights = weight.cases)
      }
    ) )

## 
## 
## `Boosted Tree Model` ...
## 
## 
## 
## Gradient Boosting Model ...

## Measure running time of R code for caretEnsemble
writeLines("\n\n Measure running time of R code = ")

## 
## 
##  Measure running time of R code =

print(Sys.time() - start_time)

## Time difference of 30.57896 mins

## Summarize results
boosting_results <- resamples(boost.models)
summary(boosting_results)

## 
## Call:
## summary.resamples(object = boosting_results)
## 
## Models: blackboost, gbm 
## Number of resamples: 10 
## 
## Accuracy 
##                 Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
## blackboost 0.8446296 0.8520833 0.8539815 0.8540556 0.8562037 0.8625926
## gbm        0.8503704 0.8633796 0.8645370 0.8635556 0.8665741 0.8690741
##            NA's
## blackboost    0
## gbm           0
## 
## Kappa 
##                 Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
## blackboost 0.6693078 0.6851577 0.6886987 0.6881149 0.6917181 0.7058233
## gbm        0.6808324 0.7065274 0.7095979 0.7078269 0.7140694 0.7191499
##            NA's
## blackboost    0
## gbm           0

dotplot(boosting_results)

Мы видим, что Градиентный бустинг (gbm) построил модель с большей точностью - коэффициент Kappa равен 0.7078.

Обычно величина Kappa большая 0.80 соответствует весьма точному результату классификации, диапазон от 0.40 до 0.80 можно интрепретировать от удовлетворительного до хорошего качества классификации, а менее 0.40 - как слабый результат классификатора.

Больше информации о применении в пакете caret методов бустинга смотрите на странице Boosting Models.

3.5 Алгоритмы баггинга

Построим классификаторы двумя наиболее эффективными для нашей задачи алгоритмами баггинга:

• Bagged CART (treebag)

• Random Forest (ranger)

Ниже приведен пример сравнения моделей Bagged CART and Random Forest (применена “Быстрая реализация Случайного Леса”). Каждый из них также допускает использование весов наблюдений для решения проблемы несбалансированных классов.

# Example of Bagging algorithms
start_time <- Sys.time()
set.seed(seed)
bag.models <- caretList( x = X[inTrain, ], y = Y[inTrain], trControl=control, 
    tuneList=list(
## Bagged CART
      {
        writeLines("\n\nBagged CART Models ...\n")    
        set.seed(seed)
        caretModelSpec(method = "treebag", metric = metrum, weights = weight.cases)
      }, 
## Random Forest with 2 Parameters
      {
        writeLines("\n\nRandom Forest Model ...\n")    
        set.seed(seed)
        caretModelSpec(method = "ranger",  metric = metrum, weights = weight.cases)
      }
    ) )

## 
## 
## Bagged CART Models ...
## 
## 
## 
## Random Forest Model ...
## 
## Growing trees.. Progress: 44%. Estimated remaining time: 40 seconds.
## Growing trees.. Progress: 91%. Estimated remaining time: 5 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 40 seconds.
## Growing trees.. Progress: 48%. Estimated remaining time: 1 minute, 7 seconds.
## Growing trees.. Progress: 72%. Estimated remaining time: 35 seconds.
## Growing trees.. Progress: 97%. Estimated remaining time: 4 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 10 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 50 seconds.
## Growing trees.. Progress: 78%. Estimated remaining time: 17 seconds.
## Growing trees.. Progress: 45%. Estimated remaining time: 37 seconds.
## Growing trees.. Progress: 93%. Estimated remaining time: 4 seconds.
## Growing trees.. Progress: 23%. Estimated remaining time: 1 minute, 43 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 1 minute, 9 seconds.
## Growing trees.. Progress: 70%. Estimated remaining time: 39 seconds.
## Growing trees.. Progress: 94%. Estimated remaining time: 7 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 11 seconds.
## Growing trees.. Progress: 37%. Estimated remaining time: 53 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 21 seconds.
## Growing trees.. Progress: 45%. Estimated remaining time: 37 seconds.
## Growing trees.. Progress: 92%. Estimated remaining time: 5 seconds.
## Growing trees.. Progress: 22%. Estimated remaining time: 1 minute, 52 seconds.
## Growing trees.. Progress: 45%. Estimated remaining time: 1 minute, 15 seconds.
## Growing trees.. Progress: 70%. Estimated remaining time: 40 seconds.
## Growing trees.. Progress: 94%. Estimated remaining time: 7 seconds.
## Growing trees.. Progress: 71%. Estimated remaining time: 12 seconds.
## Growing trees.. Progress: 37%. Estimated remaining time: 53 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 21 seconds.
## Growing trees.. Progress: 45%. Estimated remaining time: 38 seconds.
## Growing trees.. Progress: 92%. Estimated remaining time: 5 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 39 seconds.
## Growing trees.. Progress: 48%. Estimated remaining time: 1 minute, 6 seconds.
## Growing trees.. Progress: 73%. Estimated remaining time: 34 seconds.
## Growing trees.. Progress: 98%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 75%. Estimated remaining time: 10 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 51 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 19 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 35 seconds.
## Growing trees.. Progress: 94%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 40 seconds.
## Growing trees.. Progress: 48%. Estimated remaining time: 1 minute, 8 seconds.
## Growing trees.. Progress: 73%. Estimated remaining time: 35 seconds.
## Growing trees.. Progress: 97%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 75%. Estimated remaining time: 10 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 49 seconds.
## Growing trees.. Progress: 78%. Estimated remaining time: 17 seconds.
## Growing trees.. Progress: 48%. Estimated remaining time: 33 seconds.
## Growing trees.. Progress: 96%. Estimated remaining time: 2 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 37 seconds.
## Growing trees.. Progress: 49%. Estimated remaining time: 1 minute, 3 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 32 seconds.
## Growing trees.. Progress: 99%. Estimated remaining time: 1 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 9 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 49 seconds.
## Growing trees.. Progress: 78%. Estimated remaining time: 17 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 34 seconds.
## Growing trees.. Progress: 95%. Estimated remaining time: 2 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 37 seconds.
## Growing trees.. Progress: 49%. Estimated remaining time: 1 minute, 4 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 31 seconds.
## Growing trees.. Progress: 100%. Estimated remaining time: 0 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 9 seconds.
## Growing trees.. Progress: 37%. Estimated remaining time: 52 seconds.
## Growing trees.. Progress: 74%. Estimated remaining time: 21 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 34 seconds.
## Growing trees.. Progress: 95%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 24%. Estimated remaining time: 1 minute, 37 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 1 minute, 8 seconds.
## Growing trees.. Progress: 72%. Estimated remaining time: 36 seconds.
## Growing trees.. Progress: 96%. Estimated remaining time: 4 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 9 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 49 seconds.
## Growing trees.. Progress: 78%. Estimated remaining time: 17 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 34 seconds.
## Growing trees.. Progress: 94%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 22%. Estimated remaining time: 1 minute, 52 seconds.
## Growing trees.. Progress: 45%. Estimated remaining time: 1 minute, 14 seconds.
## Growing trees.. Progress: 69%. Estimated remaining time: 42 seconds.
## Growing trees.. Progress: 93%. Estimated remaining time: 8 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 10 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 49 seconds.
## Growing trees.. Progress: 78%. Estimated remaining time: 17 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 34 seconds.
## Growing trees.. Progress: 95%. Estimated remaining time: 3 seconds.
## Growing trees.. Progress: 23%. Estimated remaining time: 1 minute, 42 seconds.
## Growing trees.. Progress: 47%. Estimated remaining time: 1 minute, 8 seconds.
## Growing trees.. Progress: 71%. Estimated remaining time: 38 seconds.
## Growing trees.. Progress: 95%. Estimated remaining time: 5 seconds.
## Growing trees.. Progress: 76%. Estimated remaining time: 9 seconds.
## Growing trees.. Progress: 38%. Estimated remaining time: 50 seconds.
## Growing trees.. Progress: 77%. Estimated remaining time: 18 seconds.
## Growing trees.. Progress: 34%. Estimated remaining time: 1 minute, 0 seconds.
## Growing trees.. Progress: 68%. Estimated remaining time: 29 seconds.

## Measure running time of R code for caretEnsemble
writeLines("\n\n Measure running time of R code = ")

## 
## 
##  Measure running time of R code =

print(Sys.time() - start_time)

## Time difference of 1.180772 hours

## Summarize results
bagging_results <- resamples(bag.models)
summary(bagging_results)

## 
## Call:
## summary.resamples(object = bagging_results)
## 
## Models: treebag, ranger 
## Number of resamples: 10 
## 
## Accuracy 
##              Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## treebag 0.9250000 0.9297685 0.9328704 0.9315185 0.9335185 0.9348148    0
## ranger  0.9287037 0.9328241 0.9365741 0.9357222 0.9387963 0.9407407    0
## 
## Kappa 
##              Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## treebag 0.8342217 0.8448403 0.8519083 0.8486719 0.8531029 0.8560359    0
## ranger  0.8434873 0.8526755 0.8610680 0.8590075 0.8657373 0.8696807    0

dotplot(bagging_results)

Мы видим, что “Случайный лес” обучил модель с большей точностью, чем бустинговые алгоритмы, поскольку коэффициент Kappa составил 0.859.

Больше информации о применении в пакете caret методов баггинга смотрите на странице Bagging Models.

3.6 Построение моделей машинного обучения для стекинга (Stacking Algorithms)

Ниже приведен список моделей, доступных в пакете caret, которые функцией caretStack () будут совместно уложены для построения ансамбля - модели более высокого порядка.

Первоначально рассмотрим настройку пяти моделей для стекинга. Все эти алгоритмы для разрешения проблемы несбалансированных классов могут применять сэмплинг на увеличение наблюдений, относящихся к малочисленным классам до уровня наибольшего класса.

• Flexible Discriminant Analysis (fda)

• Lasso and Ridge (Elastic-net) regularized Generalized Linear Model (glmnet)

• Recursive Partitioning and Regression Trees (rpart)

• Neural Networks with Feature Extraction (pcaNNet)

• Nearest Shrunken Centroids (pam)

# Example of Stacking algorithms
# create submodels
start_time <- Sys.time()

## Merge two list of models into one list
library("RCurl")

## Loading required package: bitops

## 
## Attaching package: 'RCurl'

## The following object is masked from 'package:tidyr':
## 
##     complete

writeLines("\n\nFlexible Discriminant Analysis (fda) ...\n")

## 
## 
## Flexible Discriminant Analysis (fda) ...

set.seed(seed)
stack.models <- caretList( x = X[inTrain, ], y = Y[inTrain], trControl = controlImbalance,  ## sampling='up'
    tuneList=list(
## Flexible Discriminant Analysis with 1 Parameter
        caretModelSpec(method = 'fda', metric = metrum#,
                       # tuneGrid = expand.grid(.degree = 1:2,
                       # ## Product Degree - argument which allows polynomial regression of the required total degree
                       #                        .nprune = 2:14)
                       # ## Terms
                       ) ) )

## Loading required package: earth

## Loading required package: plotmo

## Loading required package: plotrix

## 
## Attaching package: 'plotrix'

## The following object is masked from 'package:scales':
## 
##     rescale

## Loading required package: TeachingDemos

writeLines("\n\nLasso and Ridge (Elastic-net) regularized GLM (glmnet) ...\n")

## 
## 
## Lasso and Ridge (Elastic-net) regularized GLM (glmnet) ...

set.seed(seed)
stack.models <- merge.list(stack.models,
  caretList( x = X[inTrain, ], y = Y[inTrain], trControl = controlImbalance,  ## sampling='up'
    tuneList=list(
## Generalized Linear Model with 2 Parameters
    caretModelSpec(method = 'glmnet', metric = metrum#,
                   # tuneGrid = expand.grid(.alpha = seq(0, 1, length=11),
                   # ## The elasticnet mixing parameter: alpha=1 is the lasso penalty, and alpha=0 the ridge penalty
                   #                        .lambda = seq(0.001, 2, length=11))
                   # ##  LASSO (least absolute shrinkage and selection operator)
                   ) ) ) )

writeLines("\n\nRecursive Partitioning and Regression Trees ...\n")

## 
## 
## Recursive Partitioning and Regression Trees ...

set.seed(seed)
stack.models <- merge.list(stack.models,
  caretList( x = X[inTrain, ], y = Y[inTrain], trControl = controlImbalance,  ## sampling='up'
    tuneList=list(
## Recursive Partitioning and Regression Trees with 1 Parameter
      caretModelSpec(method = "rpart", metric = metrum#,
                     # tuneGrid = expand.grid(.cp =  seq(0, 0.25, length = 11))
                     # ## Complexity Parameter (split even if it does not improve the tree)
                     ) ) ) )

writeLines("\n\nNeural Networks with Feature Extraction ...\n")

## 
## 
## Neural Networks with Feature Extraction ...

set.seed(seed)
stack.models <- merge.list(stack.models,
  caretList( x = X[inTrain, ], y = Y[inTrain], trControl = controlImbalance,  ## sampling='up'
    tuneList=list(
## Neural Networks with Feature Extraction with 2 Parameters
      caretModelSpec(method = "pcaNNet", metric = metrum, trace = FALSE
                     ) ) ) )

writeLines("\n\nNearest Shrunken Centroids (pam) ...\n")

## 
## 
## Nearest Shrunken Centroids (pam) ...

set.seed(seed)
stack.models <- merge.list(stack.models,
  caretList( x = X[inTrain, ], y = Y[inTrain], trControl = controlImbalance,  ## sampling='up'
    tuneList=list(
## Nearest Shrunken Centroids with 1 Parameter
      caretModelSpec(method = "pam", metric = metrum
                     ) ) ) )

## 12345678910111213141516171819202122232425262728293011111111111

## Measure running time of R code for caretEnsemble
writeLines("\n\n Measure running time of R code = ")

## 
## 
##  Measure running time of R code =

print(Sys.time() - start_time)

## Time difference of 30.59671 mins

## Summarize results
stack_results <- resamples(stack.models)
summary(stack_results)

## 
## Call:
## summary.resamples(object = stack_results)
## 
## Models: fda, glmnet, rpart, pcaNNet, pam 
## Number of resamples: 10 
## 
## Accuracy 
##              Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## fda     0.7950000 0.7995370 0.8018519 0.8027593 0.8069907 0.8122222    0
## glmnet  0.8064815 0.8131019 0.8154630 0.8153889 0.8187037 0.8212963    0
## rpart   0.6879630 0.7150463 0.7206481 0.7185370 0.7245370 0.7355556    0
## pcaNNet 0.8153704 0.8277778 0.8313889 0.8316111 0.8336111 0.8498148    0
## pam     0.7150000 0.7202315 0.7265741 0.7247222 0.7293056 0.7314815    0
## 
## Kappa 
##              Min.   1st Qu.    Median      Mean   3rd Qu.      Max. NA's
## fda     0.5708488 0.5765538 0.5857036 0.5851763 0.5939457 0.6038228    0
## glmnet  0.5917572 0.6043654 0.6099749 0.6106337 0.6189457 0.6219811    0
## rpart   0.3779996 0.3829559 0.4026566 0.4085715 0.4343819 0.4444649    0
## pcaNNet 0.6092680 0.6363531 0.6443212 0.6438821 0.6476877 0.6788000    0
## pam     0.4144456 0.4239527 0.4370120 0.4348243 0.4454725 0.4497200    0

dotplot(stack_results)

Мы видим, что модель pcaNNet (Neural Networks with Feature Extraction) содержат классификатор с наибольшим коэффициентом точности Kappa 0.6439 среди всех методов для стекинга.

## Merging Models for Ensembling
all.models <- merge.list( boost.models, bag.models )
all.models <- merge.list( all.models, stack.models )
all.models[['rpart']]$modelInfo$label <- 'Classification And Regression Trees (CART)'
all.models[['glmnet']]$modelInfo$label <- 'Lasso and Ridge (Elastic-net) regularized GLM'

## Summarize results
all_results <- resamples(all.models)
summary(all_results)

## 
## Call:
## summary.resamples(object = all_results)
## 
## Models: blackboost, gbm, treebag, ranger, fda, glmnet, rpart, pcaNNet, pam 
## Number of resamples: 10 
## 
## Accuracy 
##                 Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
## blackboost 0.8446296 0.8520833 0.8539815 0.8540556 0.8562037 0.8625926
## gbm        0.8503704 0.8633796 0.8645370 0.8635556 0.8665741 0.8690741
## treebag    0.9250000 0.9297685 0.9328704 0.9315185 0.9335185 0.9348148
## ranger     0.9287037 0.9328241 0.9365741 0.9357222 0.9387963 0.9407407
## fda        0.7950000 0.7995370 0.8018519 0.8027593 0.8069907 0.8122222
## glmnet     0.8064815 0.8131019 0.8154630 0.8153889 0.8187037 0.8212963
## rpart      0.6879630 0.7150463 0.7206481 0.7185370 0.7245370 0.7355556
## pcaNNet    0.8153704 0.8277778 0.8313889 0.8316111 0.8336111 0.8498148
## pam        0.7150000 0.7202315 0.7265741 0.7247222 0.7293056 0.7314815
##            NA's
## blackboost    0
## gbm           0
## treebag       0
## ranger        0
## fda           0
## glmnet        0
## rpart         0
## pcaNNet       0
## pam           0
## 
## Kappa 
##                 Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
## blackboost 0.6693078 0.6851577 0.6886987 0.6881149 0.6917181 0.7058233
## gbm        0.6808324 0.7065274 0.7095979 0.7078269 0.7140694 0.7191499
## treebag    0.8342217 0.8448403 0.8519083 0.8486719 0.8531029 0.8560359
## ranger     0.8434873 0.8526755 0.8610680 0.8590075 0.8657373 0.8696807
## fda        0.5708488 0.5765538 0.5857036 0.5851763 0.5939457 0.6038228
## glmnet     0.5917572 0.6043654 0.6099749 0.6106337 0.6189457 0.6219811
## rpart      0.3779996 0.3829559 0.4026566 0.4085715 0.4343819 0.4444649
## pcaNNet    0.6092680 0.6363531 0.6443212 0.6438821 0.6476877 0.6788000
## pam        0.4144456 0.4239527 0.4370120 0.4348243 0.4454725 0.4497200
##            NA's
## blackboost    0
## gbm           0
## treebag       0
## ranger        0
## fda           0
## glmnet        0
## rpart         0
## pcaNNet       0
## pam           0

dotplot(all_results)

Теперь объединим все полученные модели в один комлекс. Очевидно, что точность любой из моделей стекинга (на графике ниже они завершают список алгоритмов) уступает ансамблям, сформированным бустингом и баггингом. Но даже эти “слабые” классификаторы способны привнести свою доля точности в обобщающий ансамбль.

3.7 Корреляция предсказаний между моделями собранными для ассамблирования

Когда мы комбинируем предсказания разных моделей с использованием стекинга, желательно, чтобы выданными ими результаты, имели низкую корреляцию, т.е не были мультиколлинеарными. Это предполагает, что модели умеют решать задачу, но по-разному, позволяя создаваемому ассамблированием новому классификатору выяснить, как получить от каждой модели лучшее.

Если предсказания для разных модели были сильно взаимосвязаны (> 0,75), то они будут делать одни и те же или очень похожие предсказания в большинстве случаев, не привнося какие-либо преимущества от объединения методов.

# correlation between results
# modelCor(all_results, metric = metrum)
# splom(all_results)

## www.cmap.polytechnique.fr/~lepennec/enseignement
library("GGally")

## 
## Attaching package: 'GGally'

## The following object is masked from 'package:dplyr':
## 
##     nasa

all_results$values %>% 
  dplyr::select(., ends_with(metrum)) %>% 
    GGally::ggpairs(., xlab = metrum, ylab = metrum,
                title = "Correlations Between Results from Sub-Models in Stacking Ensemble",
                lower = list(continuous = wrap("points", alpha = 0.3), combo = wrap("dot_no_facet", alpha = 0.4)),
                columnLabels = all_results$models)

Корреляционая матрица демонстрирует, что многие пары предсказаний тесно взаимосвязаны. Например, два метода с высокой корреляцией между предсказаниями между алгоритами - баггинга Bagged Classification and Regression Trees (treebag) и Nearest Shrunken Centroids (pam).

Обнаруженная мультиколлинеарность при обращении к ассамблированию путем линейной регресии часто может приводить к неустойчивости, т. е. сильно зависящей от обучающих данных, что обычно является проявлением тенденции к переобучению. Избежать такой переподгонки помогает регуляризация – общий метод, заключающийся в наложении дополнительных ограничений (пенализации) на искомые параметры, которые могут предотвратить излишнюю сложность модели. Больше информации о применении регуляризации в пакете caret смотрите на странице Boosting Models.

Ниже сформированы три ансамбля моделей посредством стекинга с различными способами регуляризации:

• Averaged Neural Network (pcaNNet) с применением метода “главных компонент”

• Flexible Discriminant Analysis (knn) с пенализация

• Elastic-net (glmnet) с регуляризацией по Lasso и Ridge

От использования самого распространенного классификатора на основе множественной линейной регресии, именуемого в caretEnsemble “жадным” (англ. Greedy) мы из-за мультиколлиниарности отобранных методов вынуждены отказаться.

3.8 Ансамбль моделей cтекингом Нейронной Сети (pcaNNet)

Применена усредненная Нейронная сеть на главных компонентах, сами по себе являющиеся сильным средством борьбы с мультиколлинеарностью. Вообще нейронная сеть - довольно изощренный алгоритм комбинирования предсказаний, который в ряде задач дает отличные результаты модерирования.

## Stack using Model Averaged Neural Network with 2 Parameters
controlEnsemble <- controlImbalance
controlEnsemble$index <- NULL # parameter 'index' for caretEnsemble must be 'NULL'

set.seed(seed)
stack.net <- caretStack(all.models, method = "pcaNNet", metric = metrum, trControl = controlEnsemble,
                        trace = FALSE)

# plot(stack.net)
plot(varImp(stack.net$ens_model), length(stack.net$models))

getTrainPerf( stack.net$ens_model )

##   TrainAccuracy TrainKappa  method
## 1     0.9346019  0.8571419 pcaNNet

Мы видим, что в результате ассамблирования Нейронной сетью точность достигла лишь 0.8571, что меньше коэффициента Kappa у модели “Случайного леса”, т.е. одного из способов баггинга (0.859). Наибольший вклад в этот ансамбль внесли данные классификации методов: ranger и treebag.

3.9 Ансамбль моделей cтекингом Гибкий дискриминантный анализ (fda)

Для построения очередного ансамбля обратимся к Гибкому дискриминантному анализу, который не столь подвержен мультиколлиниарности как его линейный предшественник, поскольку в нем производится пенализация.

## Stack using Flexible Discriminant Analysis with 1 Parameter

set.seed(seed)
stack.fda <- caretStack(all.models, method = "fda", metric = metrum, trControl = controlEnsemble)

# plot(stack.fda)
plot(varImp(stack.fda$ens_model), length(stack.fda$models))

getTrainPerf( stack.fda$ens_model )

##   TrainAccuracy TrainKappa method
## 1       0.93475  0.8575158    fda

Точность ансамбля при гибком дискриминантном анализе (кэффициент Kappa) ставил 0.8575, что является небольшим улучшением по сравнению с применением pcaNNet для ассамблирования. Преобладающий вклад в ансамбль, обученного классификатором Гибкого дискриминантного анализа, внесли методы ranger, pam, а также gbm.

3.10 Ансамбль моделей cтекингом Эластичной сети (glmnet)

Еще одним способом избежать влияния сильной взаимосвязи предсказаний разных моделей может стать классификатор эластичная сеть — модель обобщенной регрессии с двумя регуляризаторами L1, L2.

## Stack using Regularized GLM with 2 Parameters
set.seed(seed)
stack.glm <- caretStack(all.models, method = "glmnet", metric = metrum, trControl = controlEnsemble)

# plot(stack.glm)
plot(varImp(stack.glm$ens_model), length(stack.glm$models))

getTrainPerf( stack.glm$ens_model )

##   TrainAccuracy TrainKappa method
## 1     0.9351852  0.8587493 glmnet

Ансамбль, полученный методом обобщенной регрессионной модели с регуляцизацией, получил еще большую точность - 0.8587, чем при ассамблировании алгоритмом fda. Но он по-прежнему меньше точности ансамбля, скомбинированного баггингом, т.е. “Случайным лесом”. Максимальный вклад в этот ансамбль внесли методы уже отличившиеся ранее.

3.11 Апробирование модели

После построения набора ансамблей проведем их апробирование на контрольном наборе данных (inTest), который не принимал участие в построение ансамблей. Точность на нём будем сравнивать с проверочными наблюдениями для данной задачи.

## Predict for Testing Set:
preds <- data.frame(sapply(all.models, function(x){predict(x, X[inTest, ], type = 'prob')[, 2]}))

preds$`Stack NNet` <- predict(stack.net, newdata = X[inTest, ], type = 'prob')
preds$`Stack FDA` <- predict(stack.fda, newdata = X[inTest, ], type = 'prob')
preds$`Stack GLM` <- predict(stack.glm,  newdata = X[inTest, ], type = 'prob')

# Preallocate data types
i = 1
Name <- character()          # Name
Accuracy <- numeric()        # Accuracy
Kappa <- numeric()           # Cohen's Kappa
AccuracyLearn <- numeric()   # Accuracy of Learning
KappaLearn <- numeric()      # Cohen's Kappa of Learning
Note <- character()          # long model name

# output Confusion Matrix of Training and Testing Sets
for (i in 1:dim(preds)[2] ) {
  testCM <- as.factor(ifelse(preds[i] > .5, 'Yes', 'No')[, 1]) %>%
      caret::confusionMatrix( data = ., reference = Y[inTest], positive ="Yes", mode = "everything" )

  Name[i] <- colnames(preds[i])
  AccuracyLearn[i] <- ifelse(i <= length(all.models),
                       as.numeric( round( getTrainPerf( all.models[[i]] )$TrainAccuracy, 4 ) ),
                       case_when(
                          (i == dim(preds)[2] - 2) ~ round( getTrainPerf( stack.net$ens_model )$TrainAccuracy, 4 ),
                          (i == dim(preds)[2] - 1) ~ round( getTrainPerf( stack.fda$ens_model )$TrainAccuracy, 4 ),
                          (i == dim(preds)[2]) ~     round( getTrainPerf( stack.glm$ens_model )$TrainAccuracy, 4 )
                          ) )
  KappaLearn[i] <- ifelse(i <= length(all.models),
                       as.numeric( round( getTrainPerf( all.models[[i]] )$TrainKappa, 4 ) ),
                       case_when(
                          (i == dim(preds)[2] - 2) ~ round( getTrainPerf( stack.net$ens_model )$TrainKappa, 4 ),
                          (i == dim(preds)[2] - 1) ~ round( getTrainPerf( stack.fda$ens_model )$TrainKappa, 4 ),
                          (i == dim(preds)[2]) ~     round( getTrainPerf( stack.glm$ens_model )$TrainKappa, 4 )  
                          ) )
  
  Accuracy[i] <- round( as.numeric(testCM$overall["Accuracy"]),4 )
  Kappa[i] <- round( as.numeric(testCM$overall["Kappa"]), 4 )
  Note[i] <- ifelse(i <= length(all.models), all.models[[i]]$modelInfo$label, colnames(preds[i]))
}

## coerce to data frame
result.df <- data_frame(Name, AccuracyLearn, KappaLearn, Accuracy, Kappa, Note)

library('kableExtra')

## Warning: package 'kableExtra' was built under R version 3.4.4

library('formattable')

## Warning: package 'formattable' was built under R version 3.4.4

## 
## Attaching package: 'formattable'

## The following object is masked from 'package:TeachingDemos':
## 
##     digits

## The following objects are masked from 'package:scales':
## 
##     comma, percent, scientific

## Classification results by VALIDATION and TESTING Sets from `caret` models

result.df %>% 
  mutate(Name = cell_spec(Name, "html", color = ifelse(Kappa >= (arrange(result.df, desc(Kappa))[3, 'Kappa'] %>%  as.numeric()), "red", "black")),
         AccuracyLearn, KappaLearn,  Accuracy,
         Kappa = color_bar("lightgreen")(Kappa),
         Note = Note) %>% 
knitr::kable(format = "html", digits = 4, longtable = TRUE, booktabs = TRUE, escape = F, 
             caption = "Classification results by VALIDATION and TESTING Sets from `caret` models") %>% 
         kableExtra::kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive",
                                                         full_width = FALSE)) %>% 
           column_spec(5, width = "3cm") %>%
             add_header_above(c(" ", "Validation Sets" = 2, "Test Set" = 2, " ")) %>% 
               group_rows(index = c("Boosting Ensembles" = 2, "Bagging Ensembles" = 2,
                                    "weak methods (not Еnsembles)" = 5, "Stacking Ensembles" = 3))

Classification results by VALIDATION and TESTING Sets from `caret` models
	Validation Sets		Test Set
Name	AccuracyLearn	KappaLearn	Accuracy	Kappa	Note
Boosting Ensembles
blackboost	0.8541	0.6881	0.8557	0.6957	Boosted Tree
gbm	0.8636	0.7078	0.8607	0.7073	Stochastic Gradient Boosting
Bagging Ensembles
treebag	0.9315	0.8487	0.9397	0.8689	Bagged CART
ranger	0.9357	0.8590	0.9490	0.8899	Random Forest
weak methods (not Еnsembles)
fda	0.8028	0.5852	0.7963	0.5780	Flexible Discriminant Analysis
glmnet	0.8154	0.6106	0.8200	0.6235	Lasso and Ridge (Elastic-net) regularized GLM
rpart	0.7185	0.4086	0.7163	0.3943	Classification And Regression Trees (CART)
pcaNNet	0.8316	0.6439	0.8153	0.6176	Neural Networks with Feature Extraction
pam	0.7247	0.4348	0.7317	0.4537	Nearest Shrunken Centroids
Stacking Ensembles
Stack NNet	0.9346	0.8571	0.9487	0.8892	Stack NNet
Stack FDA	0.9348	0.8575	0.9477	0.8874	Stack FDA
Stack GLM	0.9352	0.8587	0.9483	0.8889	Stack GLM

Из таблицы видно, что при сравнении предсказаний на контрольной выборке с реальным разбиением на классы тройку лидеров по аккуратности составили:

• ансамбль , обученный классификатором (ranger) - точность 0.8899

• ансамбль стэкинга, оптимизированный классификатором (Stack NNet) - точность 0.8892

• ансамбль стэкинга, построенный классификатором (Stack GLM) - точность 0.8889.

Выберем точнейший ансамбль на этом контрольном массиве - классификатор ranger.

4 Завершение модели

Полученные результаты бинарной классификации для предсказаний мошенничества являются довольно точными. Завершим моделирование, создав целевой бинарный признак TARGET на обученной модели, с использованием набор данных inProblem. Расчитанные результаты выведем в MS Excel файл.

## Solution of a Classification Problem
problem.df <- bind_cols(NUM = DF2$NUM, # c(30000:39999)
                        predict(bag.models$ranger, newdata = X[inProblem, ], type = 'prob') %>%
                          as.data.frame %>% 
                            rename("Probability" = !!names(.[1])) %>%
                              mutate( TARGET = if_else(Probability > 0.5 , 1, 0) ) )

## Open workbook into temporary file
openxlsx::addWorksheet(wb0 <- openxlsx::createWorkbook(), sheetName = "Ensembles", gridLines = FALSE)
openxlsx::writeData(wb0, sheet = 1, x = result.df, withFilter = TRUE)

openxlsx::addWorksheet(wb0, "Classification")
openxlsx::writeDataTable(wb0, sheet = "Classification", x = problem.df, colNames = TRUE,
                                 tableStyle = "TableStyleMedium4", withFilter = TRUE)
openxlsx::openXL(wb0)

5 Итоги

В этом задании мы провели решение бинарной классификационной задачи при помощи пакетов организации машинного обучения caret и формирования ансамблей caretEnsemble, используя язык программирования для статистической обработки данных R. В частности, описанные шаги были следующими:

Определение типа Проблемы (данные по клиентам компании).
Разведочный анализ данных (поиск неоднородных вариативности значений и высоких коррелляций ряда признаков).
Превращение номинальных признаков в набор числовых переменных.
Отбор предикторов для предсказания (удаление признаков с минимальной вариацией и избавление от пар тесно взаимосвязанных признаков).
Преобразование предикторов (заполнение пропусков, шкалирование и нормирование).
Построение ансамблей методом бустинга.
Построение ансамблей методом баггинга.
Построение ансамблей методом стэкинга.
Оценка полученной модели на контрольной выборке, отбор наилучшего ансамбля и вывод предсказанных значений.

Ensemble of Machine Learning Algorithms

Bank Task 2 - Classification Problem to identify Fraud in the Customer Base

Alexander Rodionov

02 Апрель 2018