Лабораторная №8
Косован Михаил
13 05 2020
Необходимо построить две модели для прогноза на основе дерева решений:
- Для непрерывной зависимой переменной;
- Для категориальной зависимой переменной.
Данные и переменные указаны в таблице с вариантами.
Ядро генератора случайных чисел – номер варианта.
Задания
Для каждой модели:
- Указать настроечные параметры метода из своего варианта (например: количество узлов, количество предикторов, скорость обучения).
- Подогнать модель на обучающей выборке (50% наблюдений). Рассчитать MSE на тестовой выборке.
- Перестроить модель с помощью метода, указанного в варианте.
- Сделать прогноз по модели с подобранными в п.3 параметрами на тестовой выборке, оценить его точность и построить график «прогноз-реализация».
Вариант - 10
Модели: бэггинг (количество предикторов).
Данные: `Wage {ISLR}’.
Деревья решений
# Загрузка пакетов
library('tree') # деревья tree()
library('GGally') # матричный график разброса ggpairs()
library('ISLR') # набор данных Auto
library('randomForest') # бэггинг
# Загрузка данных Auto
data('Auto')
# Название столбцов переменных
names(Auto)## [1] "mpg" "cylinders" "displacement" "horsepower"
## [5] "weight" "acceleration" "year" "origin"
## [9] "name"# Размерность данных
dim(Auto)## [1] 392 9# Ядро генератора случайных чисел
my.seed <- 10Модель 1 (для непрерывной зависимой переменной mpg)
# Избавляемся от Name
Auto <- Auto[, -9]
# ?Auto
head(Auto)# Матричные графики разброса переменных
p <- ggpairs(Auto[, c(1, 2:3)])
suppressMessages(print(p))
p <- ggpairs(Auto[, c(1, 4:5)])
suppressMessages(print(p))
p <- ggpairs(Auto[, c(1, 6:8)])
suppressMessages(print(p))
# Обучающая выборка
set.seed(my.seed)
# Обучающая выборка - 50%
train <- sample(1:nrow(Auto), nrow(Auto)/2)Построим дерево регрессии для зависимой переменной mpg: миль на галлон.
# Обучаем модель
tree.auto <- tree(mpg ~ ., Auto, subset = train)
summary(tree.auto)##
## Regression tree:
## tree(formula = mpg ~ ., data = Auto, subset = train)
## Variables actually used in tree construction:
## [1] "displacement" "year" "weight"
## Number of terminal nodes: 9
## Residual mean deviance: 8.357 = 1563 / 187
## Distribution of residuals:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -7.9050 -1.5020 -0.3897 0.0000 1.3200 14.2000# Визуализация
plot(tree.auto)
text(tree.auto, pretty = 0)
tree.auto # Посмотреть всё дерево в консоли## node), split, n, deviance, yval
## * denotes terminal node
##
## 1) root 196 11130.00 22.83
## 2) displacement < 198.5 111 3263.00 27.90
## 4) year < 78.5 67 1041.00 25.51
## 8) weight < 2123.5 17 113.30 30.05 *
## 9) weight > 2123.5 50 459.30 23.97
## 18) weight < 2315 20 145.80 25.90 *
## 19) weight > 2315 30 188.70 22.68 *
## 5) year > 78.5 44 1257.00 31.54
## 10) weight < 2462 21 367.90 35.53 *
## 11) weight > 2462 23 247.80 27.89 *
## 3) displacement > 198.5 85 1301.00 16.22
## 6) displacement < 284.5 31 493.50 19.60
## 12) weight < 3018 6 256.40 23.80 *
## 13) weight > 3018 25 105.70 18.59 *
## 7) displacement > 284.5 54 250.20 14.28
## 14) year < 76.5 41 115.70 13.46 *
## 15) year > 76.5 13 21.53 16.85 *# Прогноз по модели
yhat <- predict(tree.auto, newdata = Auto[-train, ])
auto.test <- Auto[-train, "mpg"]
# MSE на тестовой выборке
mse.test <- mean((yhat - auto.test)^2)
names(mse.test)[length(mse.test)] <- 'Auto.regr.tree.all'
mse.test## Auto.regr.tree.all
## 11.40431# Точность прогноза на тестовой выборке
acc.test <- sum(abs(yhat-auto.test))/sum(auto.test)
names(acc.test)[length(acc.test)] <- 'Auto.regr.tree.all'
acc.test## Auto.regr.tree.all
## 0.1051395Бэггинг (модель 1)
Используем бэггинг, причем возбмем все 7 предикторов на каждом шаге.
# бэггинг с 7 предикторами
set.seed(my.seed)
bag.auto <- randomForest(mpg ~ ., data = Auto, subset = train,
mtry = 7, importance = TRUE)
bag.auto##
## Call:
## randomForest(formula = mpg ~ ., data = Auto, mtry = 7, importance = TRUE, subset = train)
## Type of random forest: regression
## Number of trees: 500
## No. of variables tried at each split: 7
##
## Mean of squared residuals: 9.389023
## % Var explained: 83.47# прогноз
yhat.bag = predict(bag.auto, newdata = Auto[-train, ])
# MSE на тестовой
mse.test <- c(mse.test, mean((yhat.bag - auto.test)^2))
names(mse.test)[length(mse.test)] <- 'Auto.bag.model.1.7'
mse.test## Auto.regr.tree.all Auto.bag.model.1.7
## 11.404309 7.973838# Точность прогноза на тестовой выборке
acc.test <- sum(abs(yhat.bag-auto.test))/sum(auto.test)
names(acc.test)[length(acc.test)] <- 'Auto.regr.tree.model.1.7'
acc.test## Auto.regr.tree.model.1.7
## 0.08442945Ошибка на тестовой выборке равна 7.974. Можно изменить число деревьев с помощью аргумента.
# Бэггинг с 7 предикторами и 25 деревьями
bag.auto <- randomForest(mpg ~ ., data = Auto, subset = train,
mtry = 7, ntree = 25)
# прогноз
yhat.bag <- predict(bag.auto, newdata = Auto[-train, ])
# MSE на тестовой
mse.test <- c(mse.test, mean((yhat.bag - auto.test)^2))
names(mse.test)[length(mse.test)] <- 'Auto.bag.model.1.7.25'
mse.test## Auto.regr.tree.all Auto.bag.model.1.7 Auto.bag.model.1.7.25
## 11.404309 7.973838 8.404608# Точность прогноза на тестовой выборке
acc.test <- c(acc.test, sum(abs(yhat.bag-auto.test))/sum(auto.test))
names(acc.test)[length(acc.test)] <- 'Auto.regr.tree.model.1.7.25'
acc.test## Auto.regr.tree.model.1.7 Auto.regr.tree.model.1.7.25
## 0.08442945 0.08618288# График прогноз - реализация
plot(yhat.bag, auto.test)
# линия идеального прогноза
abline(0, 1)
Судя по полученным результатам наименьшая MSE наблюдается у модели с использованием бэггинга с 7 предикторами. Минимальная MSE на тестовой выборке равна 7.97, точность прогноза составила 0.08.
Модель 2 (для категориальной зависимой переменной high.medv)
Загрузим таблицу с данными по расходу бензина, лошадиной силе и другая информации для автомобилей и добавим к ней переменную high.mpg - миль на галлон:
1, если миля на галлон >= 29;
0- в противном случае.
# Новая переменная
high.mpg <- ifelse(Auto$mpg < 29, '0', '1')
# Присоединяем к таблице данных
Auto <- cbind(Auto, high.mpg)
# Название столбцов переменных
names(Auto)## [1] "mpg" "cylinders" "displacement" "horsepower"
## [5] "weight" "acceleration" "year" "origin"
## [9] "high.mpg"# Размерность данных
dim(Auto)## [1] 392 9# Матричные графики разброса переменных
p <- ggpairs(Auto[, c(9, 1:2)], aes(color = high.mpg))
suppressMessages(print(p))
p <- ggpairs(Auto[, c(9, 3:5)], aes(color = high.mpg))
suppressMessages(print(p))
p <- ggpairs(Auto[, c(9, 6:8)], aes(color = high.mpg))
suppressMessages(print(p))
Судя по графикам, класс 0 превосходит по размеру класс 1 по переменной high.mpg приблизительно в 3 раза. Классы на графиках разброса объясняющих переменных сильно смешаны, поэтому модели с непрерывной разрешающей границей вряд ли работают хорошо. Построим дерево для категориального отклика high.mpg, отбросив непрерывный отклик mpg (мы оставили его на первом графике, чтобы проверить, как сработало разделение по значению mpg = 29).
# Модель бинарного дерева
tree.auto <- tree(high.mpg ~ . -mpg, Auto)
summary(tree.auto)##
## Classification tree:
## tree(formula = high.mpg ~ . - mpg, data = Auto)
## Variables actually used in tree construction:
## [1] "horsepower" "year" "weight" "displacement"
## Number of terminal nodes: 10
## Residual mean deviance: 0.2449 = 93.56 / 382
## Misclassification error rate: 0.05102 = 20 / 392# График результата
plot(tree.auto) # Ветви
text(tree.auto, pretty = 0) # Подписи
tree.auto # Посмотреть всё дерево в консоли## node), split, n, deviance, yval, (yprob)
## * denotes terminal node
##
## 1) root 392 451.500 0 ( 0.737245 0.262755 )
## 2) horsepower < 84.5 132 161.900 1 ( 0.303030 0.696970 )
## 4) year < 76.5 51 67.350 0 ( 0.627451 0.372549 )
## 8) weight < 2091 22 23.580 1 ( 0.227273 0.772727 )
## 16) displacement < 90.5 15 7.348 1 ( 0.066667 0.933333 ) *
## 17) displacement > 90.5 7 9.561 0 ( 0.571429 0.428571 ) *
## 9) weight > 2091 29 14.560 0 ( 0.931034 0.068966 ) *
## 5) year > 76.5 81 52.220 1 ( 0.098765 0.901235 )
## 10) weight < 2255 58 0.000 1 ( 0.000000 1.000000 ) *
## 11) weight > 2255 23 29.720 1 ( 0.347826 0.652174 ) *
## 3) horsepower > 84.5 260 91.110 0 ( 0.957692 0.042308 )
## 6) year < 79.5 227 12.850 0 ( 0.995595 0.004405 )
## 12) year < 78.5 208 0.000 0 ( 1.000000 0.000000 ) *
## 13) year > 78.5 19 7.835 0 ( 0.947368 0.052632 ) *
## 7) year > 79.5 33 40.490 0 ( 0.696970 0.303030 )
## 14) displacement < 137 10 12.220 1 ( 0.300000 0.700000 )
## 28) displacement < 115.5 5 6.730 0 ( 0.600000 0.400000 ) *
## 29) displacement > 115.5 5 0.000 1 ( 0.000000 1.000000 ) *
## 15) displacement > 137 23 17.810 0 ( 0.869565 0.130435 ) *Теперь построим дерево на обучающей выборке и оценим ошибку на тестовой.
# Тестовая выборка
Auto.test <- Auto[-train,]
high.mpg.test <- high.mpg[-train]
# Строим дерево на обучающей выборке
tree.auto <- tree(high.mpg ~ . -mpg, Auto, subset = train)
# Делаем прогноз
tree.pred <- predict(tree.auto, Auto.test, type = "class")
# Матрица неточностей
tbl <- table(tree.pred, high.mpg.test)
tbl## high.mpg.test
## tree.pred 0 1
## 0 131 8
## 1 5 52# ACC на тестовой
acc.test.2 <- sum(diag(tbl))/sum(tbl)
names(acc.test.2)[length(acc.test.2)] <- 'Auto.class.tree.all.model.2'
acc.test.2## Auto.class.tree.all.model.2
## 0.9336735Обобщённая характеристика точности: доля верных прогнозов: 0.93.
Бэггинг (модель 2)
set.seed(my.seed)
bag.auto <- randomForest(high.mpg ~ . -mpg, data = Auto, subset = train,
mtry = 7, importance = TRUE)
# График и таблица относительной важности переменных
summary(bag.auto)## Length Class Mode
## call 6 -none- call
## type 1 -none- character
## predicted 196 factor numeric
## err.rate 1500 -none- numeric
## confusion 6 -none- numeric
## votes 392 matrix numeric
## oob.times 196 -none- numeric
## classes 2 -none- character
## importance 28 -none- numeric
## importanceSD 21 -none- numeric
## localImportance 0 -none- NULL
## proximity 0 -none- NULL
## ntree 1 -none- numeric
## mtry 1 -none- numeric
## forest 14 -none- list
## y 196 factor numeric
## test 0 -none- NULL
## inbag 0 -none- NULL
## terms 3 terms call# прогноз
yhat.bag <- predict(bag.auto, newdata = Auto[-train, ])
# Матрица неточностей
tbl <- table(yhat.bag, high.mpg.test)
tbl## high.mpg.test
## yhat.bag 0 1
## 0 133 9
## 1 3 51# Точность прогноза на тестовой выборке
acc.test.2 <- c(acc.test.2, sum(diag(tbl))/sum(tbl))
names(acc.test.2)[length(acc.test.2)] <- 'Auto.class.tree.model.2.7'
acc.test.2## Auto.class.tree.all.model.2 Auto.class.tree.model.2.7
## 0.9336735 0.9387755# бэггинг с 7 предикторами и 25 деревьями
bag.auto <- randomForest(high.mpg ~ .-mpg, data = Auto, subset = train,
mtry = 7, ntree = 25)
# прогноз
yhat.bag <- predict(bag.auto, newdata = Auto[-train, ])
# Матрица неточностей
tbl <- table(yhat.bag, high.mpg.test)
tbl## high.mpg.test
## yhat.bag 0 1
## 0 133 8
## 1 3 52# Точность прогноза на тестовой выборке
acc.test.2 <- c(acc.test.2, sum(diag(tbl))/sum(tbl))
names(acc.test.2)[length(acc.test.2)] <- 'Auto.class.tree.model.2.7.25'
acc.test.2## Auto.class.tree.all.model.2 Auto.class.tree.model.2.7
## 0.9336735 0.9387755
## Auto.class.tree.model.2.7.25
## 0.9438776# График "прогноз - реализация"
plot(yhat.bag, Auto$high.mpg[-train])
Точность модели на тестовой выборке с применением бэггинга с 7 предикторами и 25 деревьями является самой высокой и равна 0.94.