Изучить графический разведочный анализ данных, оценку информативности признаков, дискретизацию непрерывных переменных и автоматический выбор признаков. Используются пакеты caret, FSelector, arules и Boruta.
В работе четыре вычислительных задания и инструкция публикации отчёта.
Перед первым нажатием Knit выполнить в консоли RStudio следующий код. Этот блок показан в отчёте, но установка пакетов при сборке отчёта не запускается. Рекомендуется актуальная версия R; для новых версий arules может потребоваться обновление R. FSelector использует Java через свои зависимости: на Windows нужно установить JDK той же разрядности, что и R, затем перезапустить RStudio.
packages <- c("caret", "FSelector", "arules", "Boruta", "mlbench",
"ranger", "randomForest", "knitr", "rmarkdown")
new_packages <- packages[!vapply(packages, requireNamespace,
logical(1), quietly = TRUE)]
if (length(new_packages)) {
install.packages(new_packages, repos = "https://cloud.r-project.org")
}Фиксированные значения set.seed() позволяют повторить
генерацию данных и случайные процедуры в одинаковой программной среде.
Версии пакетов приведены в конце отчёта; при изменении версий отдельные
результаты могут отличаться. Каждый рисунок сохраняется в папку
graphs в формате JPG и включается в HTML.
## [1] "ada" "AdaBag" "AdaBoost.M1" "adaboost"
## [5] "amdai" "ANFIS" "avNNet" "awnb"
## [9] "awtan" "bag" "bagEarth" "bagEarthGCV"
## [13] "bagFDA" "bagFDAGCV" "bam" "bartMachine"
## [17] "bayesglm" "binda" "blackboost" "blasso"
## [21] "blassoAveraged" "bridge" "brnn" "BstLm"
## [25] "bstSm" "bstTree" "C5.0" "C5.0Cost"
## [29] "C5.0Rules" "C5.0Tree" "cforest" "chaid"
## [33] "CSimca" "ctree" "ctree2" "cubist"
## [37] "dda" "deepboost" "DENFIS" "dnn"
## [41] "dwdLinear" "dwdPoly" "dwdRadial" "earth"
## [45] "elm" "enet" "evtree" "extraTrees"
## [49] "fda" "FH.GBML" "FIR.DM" "foba"
## [53] "FRBCS.CHI" "FRBCS.W" "FS.HGD" "gam"
## [57] "gamboost" "gamLoess" "gamSpline" "gaussprLinear"
## [61] "gaussprPoly" "gaussprRadial" "gbm_h2o" "gbm"
## [65] "gcvEarth" "GFS.FR.MOGUL" "GFS.LT.RS" "GFS.THRIFT"
## [69] "glm.nb" "glm" "glmboost" "glmnet_h2o"
## [73] "glmnet" "glmStepAIC" "gpls" "hda"
## [77] "hdda" "hdrda" "HYFIS" "icr"
## [81] "J48" "JRip" "kernelpls" "kknn"
## [85] "knn" "krlsPoly" "krlsRadial" "lars"
## [89] "lars2" "lasso" "lda" "lda2"
## [93] "leapBackward" "leapForward" "leapSeq" "Linda"
## [97] "lm" "lmStepAIC" "LMT" "loclda"
## [101] "logicBag" "LogitBoost" "logreg" "lssvmLinear"
## [105] "lssvmPoly" "lssvmRadial" "lvq" "M5"
## [109] "M5Rules" "manb" "mda" "Mlda"
## [113] "mlp" "mlpKerasDecay" "mlpKerasDecayCost" "mlpKerasDropout"
## [117] "mlpKerasDropoutCost" "mlpML" "mlpSGD" "mlpWeightDecay"
## [121] "mlpWeightDecayML" "monmlp" "msaenet" "multinom"
## [125] "mxnet" "mxnetAdam" "naive_bayes" "nb"
## [129] "nbDiscrete" "nbSearch" "neuralnet" "nnet"
## [133] "nnls" "nodeHarvest" "null" "OneR"
## [137] "ordinalNet" "ordinalRF" "ORFlog" "ORFpls"
## [141] "ORFridge" "ORFsvm" "ownn" "pam"
## [145] "parRF" "PART" "partDSA" "pcaNNet"
## [149] "pcr" "pda" "pda2" "penalized"
## [153] "PenalizedLDA" "plr" "pls" "plsRglm"
## [157] "polr" "ppr" "pre" "PRIM"
## [161] "protoclass" "qda" "QdaCov" "qrf"
## [165] "qrnn" "randomGLM" "ranger" "rbf"
## [169] "rbfDDA" "Rborist" "rda" "regLogistic"
## [173] "relaxo" "rf" "rFerns" "RFlda"
## [177] "rfRules" "ridge" "rlda" "rlm"
## [181] "rmda" "rocc" "rotationForest" "rotationForestCp"
## [185] "rpart" "rpart1SE" "rpart2" "rpartCost"
## [189] "rpartScore" "rqlasso" "rqnc" "RRF"
## [193] "RRFglobal" "rrlda" "RSimca" "rvmLinear"
## [197] "rvmPoly" "rvmRadial" "SBC" "sda"
## [201] "sdwd" "simpls" "SLAVE" "slda"
## [205] "smda" "snn" "sparseLDA" "spikeslab"
## [209] "spls" "stepLDA" "stepQDA" "superpc"
## [213] "svmBoundrangeString" "svmExpoString" "svmLinear" "svmLinear2"
## [217] "svmLinear3" "svmLinearWeights" "svmLinearWeights2" "svmPoly"
## [221] "svmRadial" "svmRadialCost" "svmRadialSigma" "svmRadialWeights"
## [225] "svmSpectrumString" "tan" "tanSearch" "treebag"
## [229] "vbmpRadial" "vglmAdjCat" "vglmContRatio" "vglmCumulative"
## [233] "widekernelpls" "WM" "wsrf" "xgbDART"
## [237] "xgbLinear" "xgbTree" "xyf"
## [1] 239
Команда вывела 239 идентификаторов моделей,
доступных через train(). Это каталог моделей, а не
самостоятельный перечень методов выбора признаков. Наличие модели в
каталоге не означает, что установлены все её дополнительные
зависимости.
Для выбора признаков в CARET предусмотрены:
| Средство | Назначение |
|---|---|
rfe() |
Рекурсивное исключение признаков с оценкой качества подмножеств |
sbf() |
Отбор признаков фильтром внутри процедуры ресэмплинга |
gafs() |
Поиск подмножества признаков генетическим алгоритмом |
safs() |
Поиск подмножества методом имитации отжига |
varImp() |
Оценка важности; сама по себе не задаёт правило отбора |
Дополнительно выделим модели, отмеченные в каталоге тегом встроенного отбора:
model_info <- getModelInfo()
embedded_names <- names(model_info)[vapply(model_info, function(m) {
any(grepl("Feature Selection", m$tags, ignore.case = TRUE))
}, logical(1))]
print(embedded_names)## [1] "ada" "AdaBag" "AdaBoost.M1" "adaboost" "bagEarth"
## [6] "bagEarthGCV" "bagFDA" "bagFDAGCV" "bartMachine" "blasso"
## [11] "BstLm" "bstSm" "C5.0" "C5.0Cost" "C5.0Rules"
## [16] "C5.0Tree" "cforest" "chaid" "ctree" "ctree2"
## [21] "cubist" "deepboost" "earth" "enet" "evtree"
## [26] "extraTrees" "fda" "foba" "gamboost" "gbm_h2o"
## [31] "gbm" "gcvEarth" "glmnet_h2o" "glmnet" "glmStepAIC"
## [36] "J48" "JRip" "lars" "lars2" "lasso"
## [41] "leapBackward" "leapForward" "leapSeq" "lmStepAIC" "LMT"
## [46] "LogitBoost" "M5" "M5Rules" "msaenet" "nodeHarvest"
## [51] "OneR" "ordinalNet" "ordinalRF" "ORFlog" "ORFpls"
## [56] "ORFridge" "ORFsvm" "pam" "parRF" "PART"
## [61] "penalized" "PenalizedLDA" "qrf" "ranger" "Rborist"
## [66] "relaxo" "rf" "rFerns" "rfRules" "rotationForest"
## [71] "rotationForestCp" "rpart" "rpart1SE" "rpart2" "rpartCost"
## [76] "rpartScore" "rqlasso" "rqnc" "RRF" "RRFglobal"
## [81] "sdwd" "smda" "sparseLDA" "spikeslab" "stepLDA"
## [86] "stepQDA" "wsrf" "xgbDART" "xgbLinear" "xgbTree"
set.seed(233332)
x <- matrix(rnorm(50*5), ncol=5)
y <- factor(rep(c("A", "B"), 25))
colnames(x) <- paste0("X", 1:5)
dim(x)## [1] 50 5
## y
## A B
## 25 25
## X1 X2 X3 X4 X5
## Min. :-2.33681 Min. :-2.2586 Min. :-1.8520 Min. :-1.83286 Min. :-1.8567
## 1st Qu.:-0.63214 1st Qu.:-0.8816 1st Qu.:-0.9007 1st Qu.:-0.68846 1st Qu.:-0.3613
## Median : 0.09398 Median : 0.1130 Median :-0.2494 Median :-0.06545 Median : 0.1450
## Mean : 0.05918 Mean :-0.1163 Mean :-0.2176 Mean :-0.02486 Mean : 0.1461
## 3rd Qu.: 0.81534 3rd Qu.: 0.5570 3rd Qu.: 0.3765 3rd Qu.: 0.75455 3rd Qu.: 0.7775
## Max. : 2.42273 Max. : 1.6827 Max. : 1.4769 Max. : 2.39491 Max. : 2.3676
Матрица содержит 50 наблюдений и пять признаков. Каждый признак сгенерирован из стандартного нормального распределения. Классы A и B содержат по 25 наблюдений и назначены независимо от признаков.
p_box <- featurePlot(x = x, y = y, plot = "box",
auto.key = list(columns = 2),
scales = list(y = list(relation = "free")))
save_lattice(p_box, "caret_box.jpg")
knitr::include_graphics("graphs/caret_box.jpg")Диаграммы размаха позволяют сравнить медианы, межквартильные интервалы и крайние наблюдения. Перекрытие распределений показывает, что однозначно разделить классы по отдельному признаку трудно. Точки за пределами усов — потенциальные выбросы по правилу диаграммы, а не доказанные ошибки данных.
p_density <- featurePlot(x = x, y = y, plot = "density",
auto.key = list(columns = 2),
scales = list(x = list(relation = "free"),
y = list(relation = "free")))
save_lattice(p_density, "caret_density.jpg")
knitr::include_graphics("graphs/caret_density.jpg")Плотности классов могут немного различаться из-за небольшого объёма выборки. Такие различия не подтверждают устойчивую связь с классом: по способу генерации зависимость отсутствует.
p_pairs <- featurePlot(x = x, y = y, plot = "pairs",
auto.key = list(columns = 2))
save_lattice(p_pairs, "caret_pairs.jpg", width = 1800, height = 1800)
knitr::include_graphics("graphs/caret_pairs.jpg")Матрица рассеяния показывает попарные связи признаков и размещение классов. По конструкции данные не содержат заданных кластеров или линейной связи между признаками. Случайные визуальные закономерности возможны.
class_means <- aggregate(as.data.frame(x), by = list(Class = y), FUN = mean)
knitr::kable(class_means, digits = 3, caption = "Средние значения по классам")| Class | X1 | X2 | X3 | X4 | X5 |
|---|---|---|---|---|---|
| A | 0.324 | -0.216 | -0.169 | 0.031 | -0.059 |
| B | -0.206 | -0.016 | -0.266 | -0.080 | 0.351 |
Вывод. featurePlot() позволяет
рассмотреть данные несколькими способами. На этой искусственной выборке
наблюдаемые отличия классов обусловлены случайной генерацией. Сделать
вывод о полезности признака для новых данных по одним графикам нельзя.
Для оценки качества модели нужна независимая проверка. Получены файлы
caret_box.jpg, caret_density.jpg,
caret_pairs.jpg.
## 'data.frame': 150 obs. of 5 variables:
## $ Sepal.Length: num 5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
## $ Sepal.Width : num 3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
## $ Petal.Length: num 1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
## $ Petal.Width : num 0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
## $ Species : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...
##
## setosa versicolor virginica
## 50 50 50
Набор содержит 150 растений трёх видов, по 50 растений каждого вида.
Цель классификации — Species. Предикторы: длина и ширина
чашелистика (Sepal.Length, Sepal.Width) и
лепестка (Petal.Length, Petal.Width).
Используем три энтропийных фильтра FSelector: информационный выигрыш, отношение выигрыша и симметричную неопределённость. Все оценивают отдельные признаки относительно класса; их абсолютные шкалы отличаются.
ig <- FSelector::information.gain(Species ~ ., iris, unit = "log2")
gr <- FSelector::gain.ratio(Species ~ ., iris)
su <- FSelector::symmetrical.uncertainty(Species ~ ., iris)
features <- rownames(ig)
importance <- data.frame(
Feature = features,
InformationGain_bits = ig$attr_importance,
GainRatio = gr[features, "attr_importance"],
SymmetricalUncertainty = su[features, "attr_importance"]
)
importance <- importance[order(-importance$InformationGain_bits), ]
rownames(importance) <- NULL
knitr::kable(importance, digits = 4, caption = "Оценки важности признаков iris")| Feature | InformationGain_bits | GainRatio | SymmetricalUncertainty |
|---|---|---|---|
| Petal.Width | 1.3784 | 0.8714 | 0.8705 |
| Petal.Length | 1.3565 | 0.8585 | 0.8572 |
| Sepal.Length | 0.6523 | 0.4196 | 0.4156 |
| Sepal.Width | 0.3856 | 0.2473 | 0.2453 |
write.csv(importance, "results/iris_importance.csv", row.names = FALSE)
top2 <- FSelector::cutoff.k(ig, 2)
selected_formula <- FSelector::as.simple.formula(top2, "Species")
top2## [1] "Petal.Width" "Petal.Length"
## Species ~ Petal.Width + Petal.Length
## <environment: 0x5561124c16c8>
save_base(function() {
par(mfrow = c(1, 3), mar = c(8, 4, 3, 1))
for (column in names(importance)[-1]) {
barplot(importance[[column]], names.arg = importance$Feature,
las = 2, col = "steelblue", main = column, ylab = "Importance")
}
}, "iris_importance.jpg", width = 2160)
knitr::include_graphics("graphs/iris_importance.jpg")Максимальный информационный выигрыш получил
Petal.Width: 1.3784 бит. По
cutoff.k(ig, 2) выбраны Petal.Width,
Petal.Length. Полный порядок по информационному выигрышу:
Petal.Width → Petal.Length → Sepal.Length →
Sepal.Width.
Большая оценка означает более сильную индивидуальную связь признака с видом растения. Длина и ширина лепестка обычно хорошо различают виды iris, что следует сопоставить с полученной таблицей. Фильтр не оценивает все сочетания признаков; два высоко оценённых признака могут содержать общую информацию. Отбор здесь описательный, по полной выборке. При оценке классификатора отбор нужно выполнять только на обучающей части каждого разбиения, иначе возникает утечка информации. Повышение точности после отбора в данной работе не измерялось.
В качестве непрерывной переменной выбрана длина чашелистика
Sepal.Length. Во всех методах задаются три категории. Для
fixed используются заранее выбранные границы 5.5 и 6.5 см;
это учебные пороги, а не границы видов.
z <- iris$Sepal.Length
set.seed(233332)
disc <- list(
interval = arules::discretize(z, method = "interval", breaks = 3),
frequency = arules::discretize(z, method = "frequency", breaks = 3),
cluster = arules::discretize(z, method = "cluster", breaks = 3, nstart = 25),
fixed = arules::discretize(z, method = "fixed",
breaks = c(-Inf, 5.5, 6.5, Inf))
)
stopifnot(all(vapply(disc, function(v) length(v) == 150L &&
is.factor(v) && !anyNA(v), logical(1))))
disc_table <- do.call(rbind, lapply(names(disc), function(m) {
freq <- table(disc[[m]])
data.frame(Method = m, Interval = names(freq), Count = as.integer(freq))
}))
knitr::kable(disc_table, caption = "Интервалы и частоты после дискретизации")| Method | Interval | Count |
|---|---|---|
| interval | [4.3,5.5) | 52 |
| interval | [5.5,6.7) | 70 |
| interval | [6.7,7.9] | 28 |
| frequency | [4.3,5.4) | 46 |
| frequency | [5.4,6.3) | 53 |
| frequency | [6.3,7.9] | 51 |
| cluster | [4.3,5.45) | 52 |
| cluster | [5.45,6.46) | 63 |
| cluster | [6.46,7.9] | 35 |
| fixed | [-Inf,5.5) | 52 |
| fixed | [5.5,6.5) | 63 |
| fixed | [6.5, Inf] | 35 |
write.csv(disc_table, "results/discretization_counts.csv", row.names = FALSE)
disc_values <- data.frame(Sepal.Length = z, as.data.frame(disc),
Species = iris$Species)
knitr::kable(head(disc_values, 12), caption = "Первые 12 преобразованных наблюдений")| Sepal.Length | interval | frequency | cluster | fixed | Species |
|---|---|---|---|---|---|
| 5.1 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.9 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.7 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.6 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 5.0 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 5.4 | [4.3,5.5) | [5.4,6.3) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.6 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 5.0 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.4 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.9 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 5.4 | [4.3,5.5) | [5.4,6.3) | [4.3,5.45) | [-Inf,5.5) | setosa |
| 4.8 | [4.3,5.5) | [4.3,5.4) | [4.3,5.45) | [-Inf,5.5) | setosa |
breaks — актуальное имя аргумента: прежний
categories устарел. interval делит диапазон на
равные отрезки; frequency использует квантили;
cluster строит разбиение по k-means; fixed
применяет заданные пороги. У frequency одинаковые значения
не разделяются между категориями, поэтому численности могут отличаться.
Бесконечные крайние границы fixed охватывают все исходные
значения.
save_base(function() {
par(mfrow = c(2, 2), mar = c(4, 4, 3, 1))
for (m in names(disc)) {
hist(z, breaks = 15, col = "lightblue", border = "white",
main = m, xlab = "Sepal.Length (cm)", ylab = "Count")
bounds <- attr(disc[[m]], "discretized:breaks")
abline(v = bounds[is.finite(bounds)], col = "firebrick", lwd = 2, lty = 2)
}
}, "iris_discretization.jpg", height = 1440)
knitr::include_graphics("graphs/iris_discretization.jpg")Границы берутся из уже выполненной дискретизации, поэтому повторный запуск случайной кластеризации для рисунка не требуется.
for (m in names(disc)) {
cat("\nМетод:", m, "\n")
print(table(Category = disc[[m]], Species = iris$Species))
}##
## Метод: interval
## Species
## Category setosa versicolor virginica
## [4.3,5.5) 45 6 1
## [5.5,6.7) 5 38 27
## [6.7,7.9] 0 6 22
##
## Метод: frequency
## Species
## Category setosa versicolor virginica
## [4.3,5.4) 40 5 1
## [5.4,6.3) 10 31 12
## [6.3,7.9] 0 14 37
##
## Метод: cluster
## Species
## Category setosa versicolor virginica
## [4.3,5.45) 45 6 1
## [5.45,6.46) 5 35 23
## [6.46,7.9] 0 9 26
##
## Метод: fixed
## Species
## Category setosa versicolor virginica
## [-Inf,5.5) 45 6 1
## [5.5,6.5) 5 35 23
## [6.5, Inf] 0 9 26
Таблицы по видам нужны для интерпретации: вид растения не используется при расчёте границ этих четырёх методов.
| Метод | Что сохраняет | Ограничение |
|---|---|---|
| interval | Одинаковую ширину интервалов | Численность категорий может сильно различаться |
| frequency | Приблизительно равные численности | Повторы мешают получить строго равные группы |
| cluster | Группировку близких значений | Зависит от распределения и случайной инициализации |
| fixed | Заданные интерпретируемые пороги | Пороги требуют содержательного обоснования |
В данном запуске численности категорий составили:
Дискретизация преобразовала значения в факторы и уменьшила детализацию. Ни один из методов сам по себе не оптимизирует качество классификации. Для новых наблюдений нужно применять границы, рассчитанные по обучающей выборке, а не вычислять новые границы отдельно на тестовой выборке.
## [1] 366 13
## 'data.frame': 366 obs. of 13 variables:
## $ V1 : Factor w/ 12 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ V2 : Factor w/ 31 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ V3 : Factor w/ 7 levels "1","2","3","4",..: 4 5 6 7 1 2 3 4 5 6 ...
## $ V4 : num 3 3 3 5 5 6 4 4 6 7 ...
## $ V5 : num 5480 5660 5710 5700 5760 5720 5790 5790 5700 5700 ...
## $ V6 : num 8 6 4 3 3 4 6 3 3 3 ...
## $ V7 : num 20 NA 28 37 51 69 19 25 73 59 ...
## $ V8 : num NA 38 40 45 54 35 45 55 41 44 ...
## $ V9 : num NA NA NA NA 45.3 ...
## $ V10: num 5000 NA 2693 590 1450 ...
## $ V11: num -15 -14 -25 -24 25 15 -33 -28 23 -2 ...
## $ V12: num 30.6 NA 47.7 55 57 ...
## $ V13: num 200 300 250 100 60 60 100 250 120 120 ...
missing_counts <- colSums(is.na(Ozone))
knitr::kable(data.frame(Feature = names(missing_counts),
Missing = as.integer(missing_counts)),
caption = "Пропущенные значения Ozone")| Feature | Missing |
|---|---|
| V1 | 0 |
| V2 | 0 |
| V3 | 0 |
| V4 | 5 |
| V5 | 12 |
| V6 | 0 |
| V7 | 15 |
| V8 | 2 |
| V9 | 139 |
| V10 | 15 |
| V11 | 1 |
| V12 | 14 |
| V13 | 0 |
ozo <- na.omit(Ozone)
removed_n <- nrow(Ozone) - nrow(ozo)
c(original = nrow(Ozone), complete = nrow(ozo), removed = removed_n)## original complete removed
## 366 203 163
Ozone содержит 366 ежедневных наблюдений за загрязнением воздуха в
Лос-Анджелесе в 1976 году и 13 переменных. Целевая переменная
V4 — максимальное за день среднее часовое значение озона. В
отличие от предыдущей задачи это регрессия;
преобразовывать V4 в фактор без дополнительного условия не
нужно.
После удаления строк с пропусками осталось 203 наблюдений, удалено 163. Метод прост, но уменьшает выборку и может внести смещение, если пропуски неслучайны.
| Переменная | Смысл |
|---|---|
| V1–V3 | Месяц, день месяца, день недели |
| V4 | Целевое значение озона |
| V5 | Высота поверхности давления 500 мбар |
| V6–V7 | Скорость ветра и влажность |
| V8–V9 | Температура на двух станциях |
| V10 | Высота основания инверсии |
| V11 | Градиент давления |
| V12 | Температура основания инверсии |
| V13 | Видимость |
Сохраняется исходное числовое кодирование календарных признаков, как в справочном примере Boruta. Специальные календарные преобразования не выполняются.
set.seed(233332)
boruta_fit <- Boruta::Boruta(V4 ~ ., data = ozo,
maxRuns = 300, doTrace = 0,
holdHistory = TRUE,
getImp = Boruta::getImpRfZ, num.threads = 1)
print(boruta_fit)## Boruta performed 30 iterations in 3.657959 secs.
## 9 attributes confirmed important: V1, V10, V11, V12, V13 and 4 more;
## 3 attributes confirmed unimportant: V2, V3, V6;
initial_stats <- Boruta::attStats(boruta_fit)
knitr::kable(initial_stats, digits = 3,
caption = "Статистика важности и исходные решения Boruta")| meanImp | medianImp | minImp | maxImp | normHits | decision | |
|---|---|---|---|---|---|---|
| V1 | 9.591 | 9.497 | 8.266 | 11.353 | 1.000 | Confirmed |
| V2 | 1.018 | 0.827 | -0.743 | 3.845 | 0.200 | Rejected |
| V3 | -1.746 | -2.081 | -2.734 | -0.483 | 0.000 | Rejected |
| V5 | 8.923 | 8.914 | 7.222 | 10.419 | 1.000 | Confirmed |
| V6 | 1.426 | 1.373 | 0.069 | 2.937 | 0.067 | Rejected |
| V7 | 11.646 | 11.628 | 9.903 | 12.827 | 1.000 | Confirmed |
| V8 | 17.079 | 16.961 | 16.150 | 18.389 | 1.000 | Confirmed |
| V9 | 19.287 | 19.420 | 17.348 | 20.662 | 1.000 | Confirmed |
| V10 | 10.036 | 9.928 | 8.952 | 11.330 | 1.000 | Confirmed |
| V11 | 12.026 | 11.956 | 10.333 | 13.585 | 1.000 | Confirmed |
| V12 | 14.681 | 14.697 | 13.750 | 15.565 | 1.000 | Confirmed |
| V13 | 9.671 | 9.760 | 8.467 | 10.643 | 1.000 | Confirmed |
Boruta сравнивает важность признаков со случайными теневыми
признаками, которые получаются перестановкой значений. Источник важности
явно задан как getImpRfZ: он использует случайный лес из
ranger. Решения: Confirmed — подтверждён,
Rejected — отклонён, Tentative — пока не
определён.
save_base(function() {
par(mar = c(8, 4, 3, 1))
plot(boruta_fit, las = 2, cex.axis = 0.8,
main = "Boruta: Ozone (response V4)", xlab = "")
}, "ozone_boruta_boxplot.jpg", width = 2160, height = 1440)
knitr::include_graphics("graphs/ozone_boruta_boxplot.jpg")plot() для объекта Boruta строит именно boxplot
распределений важности по итерациям. Зелёный цвет обозначает
подтверждённые признаки, красный — отклонённые, жёлтый — неопределённые,
синий — теневые. Этот рисунок отражает исходное статистическое решение
до дополнительной обработки.
had_tentative <- any(boruta_fit$finalDecision == "Tentative")
boruta_final <- if (had_tentative) Boruta::TentativeRoughFix(boruta_fit) else boruta_fit
final_stats <- Boruta::attStats(boruta_final)
final_table <- data.frame(Feature = rownames(final_stats), final_stats,
InitialDecision = initial_stats[rownames(final_stats), "decision"],
row.names = NULL)
final_table <- final_table[order(-final_table$medianImp), ]
knitr::kable(final_table, digits = 3,
caption = "Итоговые оценки; исходные решения сохранены для сравнения")| Feature | meanImp | medianImp | minImp | maxImp | normHits | decision | InitialDecision | |
|---|---|---|---|---|---|---|---|---|
| 8 | V9 | 19.287 | 19.420 | 17.348 | 20.662 | 1.000 | Confirmed | Confirmed |
| 7 | V8 | 17.079 | 16.961 | 16.150 | 18.389 | 1.000 | Confirmed | Confirmed |
| 11 | V12 | 14.681 | 14.697 | 13.750 | 15.565 | 1.000 | Confirmed | Confirmed |
| 10 | V11 | 12.026 | 11.956 | 10.333 | 13.585 | 1.000 | Confirmed | Confirmed |
| 6 | V7 | 11.646 | 11.628 | 9.903 | 12.827 | 1.000 | Confirmed | Confirmed |
| 9 | V10 | 10.036 | 9.928 | 8.952 | 11.330 | 1.000 | Confirmed | Confirmed |
| 12 | V13 | 9.671 | 9.760 | 8.467 | 10.643 | 1.000 | Confirmed | Confirmed |
| 1 | V1 | 9.591 | 9.497 | 8.266 | 11.353 | 1.000 | Confirmed | Confirmed |
| 4 | V5 | 8.923 | 8.914 | 7.222 | 10.419 | 1.000 | Confirmed | Confirmed |
| 5 | V6 | 1.426 | 1.373 | 0.069 | 2.937 | 0.067 | Rejected | Rejected |
| 2 | V2 | 1.018 | 0.827 | -0.743 | 3.845 | 0.200 | Rejected | Rejected |
| 3 | V3 | -1.746 | -2.081 | -2.734 | -0.483 | 0.000 | Rejected | Rejected |
confirmed_initial <- Boruta::getSelectedAttributes(boruta_fit, withTentative = FALSE)
confirmed <- Boruta::getSelectedAttributes(boruta_final, withTentative = FALSE)
rejected <- names(boruta_final$finalDecision)[boruta_final$finalDecision == "Rejected"]
unresolved_initial <- names(boruta_fit$finalDecision)[boruta_fit$finalDecision == "Tentative"]
write.csv(final_table, "results/ozone_boruta_stats.csv", row.names = FALSE)
writeLines(confirmed, "results/ozone_selected_features.txt")
saveRDS(boruta_fit, "results/ozone_boruta_original.rds")
saveRDS(boruta_final, "results/ozone_boruta_final.rds")Вывод по выполненному запуску. Исходно подтверждены
V1, V5, V7, V8, V9, V10, V11, V12, V13. Неопределёнными
остались нет. TentativeRoughFix()
не потребовался. Это дополнительное эвристическое
решение по медианам важности, а не эквивалент продолжения основного
статистического теста.
Итоговое подмножество: V1, V5, V7, V8, V9, V10, V11, V12, V13 (9 признаков); отклонены V2, V3, V6. Лидирует по медиане важности V9. Нужно опираться на таблицу и boxplot, а не заранее назначать важные признаки. Алгоритм ищет все релевантные признаки, поэтому может сохранять несколько связанных между собой предикторов. Важность не доказывает причинную связь.
Выбор признаков выполнен на всём очищенном наборе для изучения алгоритма. Изменение ошибки прогноза не измерялось. При последующей проверке модели отбор следует повторять внутри обучающего разбиения; поскольку наблюдения ежедневные, для прогноза будущих дней разумно учитывать порядок времени.
В работе получены три графика разведочного анализа случайных данных CARET, таблица и график важности признаков iris, четыре разбиения длины чашелистика на категории и boxplot выбора признаков Ozone. Все рисунки сохранены в JPG.
Разведочный анализ показывает распределения и возможные связи, фильтры FSelector оценивают отдельные признаки, а Boruta сравнивает их важность со случайным ориентиром. Дискретизация меняет представление данных и не заменяет выбор признаков. Результаты отбора требуют проверки предсказательной полезности на новых данных.
| File | Size_KB |
|---|---|
| caret_box.jpg | 166.5 |
| caret_density.jpg | 233.7 |
| caret_pairs.jpg | 534.8 |
| iris_discretization.jpg | 255.1 |
| iris_importance.jpg | 198.1 |
| ozone_boruta_boxplot.jpg | 156.8 |
.Rmd в RStudio и установить пакеты из
раздела подготовки.HTML содержит встроенные изображения и подходит для публикации. Ссылка появляется после фактической публикации: локальный HTML её не заменяет.
versions <- data.frame(Package = required,
Version = vapply(required, function(p)
as.character(packageVersion(p)), character(1)))
knitr::kable(versions)| Package | Version | |
|---|---|---|
| caret | caret | 6.0.94 |
| FSelector | FSelector | 0.34 |
| arules | arules | 1.7.7 |
| Boruta | Boruta | 8.0.0 |
| mlbench | mlbench | 2.1.3.1 |
| ranger | ranger | 0.16.0 |
| randomForest | randomForest | 4.7.1.1 |
| knitr | knitr | 1.45 |
| rmarkdown | rmarkdown | 2.25 |
## R version 4.3.3 (2024-02-29)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 24.04.3 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.12.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0
##
## locale:
## [1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
## [4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
## [7] LC_PAPER=C.UTF-8 LC_NAME=C.UTF-8 LC_ADDRESS=C.UTF-8
## [10] LC_TELEPHONE=C.UTF-8 LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C.UTF-8
##
## time zone: Europe/Paris
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] caret_6.0-94 lattice_0.22-5 ggplot2_3.4.4
##
## loaded via a namespace (and not attached):
## [1] gtable_0.3.4 xfun_0.41 bslib_0.6.1 recipes_1.0.9
## [5] rJava_1.0-11 vctrs_0.6.5 tools_4.3.3 generics_0.1.3
## [9] stats4_4.3.3 parallel_4.3.3 tibble_3.2.1 fansi_1.0.5
## [13] highr_0.10 FSelector_0.34 pkgconfig_2.0.3 ModelMetrics_1.2.2.2
## [17] Matrix_1.6-5 data.table_1.14.10 lifecycle_1.0.4 compiler_4.3.3
## [21] stringr_1.5.1 mlbench_2.1-3.1 munsell_0.5.0 codetools_0.2-19
## [25] Boruta_8.0.0 arules_1.7-7 htmltools_0.5.7 class_7.3-22
## [29] sass_0.4.8 yaml_2.3.8 prodlim_2023.08.28 pillar_1.9.0
## [33] jquerylib_0.1.4 MASS_7.3-60.0.1 entropy_1.3.2 cachem_1.0.8
## [37] gower_1.0.1 iterators_1.0.14 rpart_4.1.23 foreach_1.5.2
## [41] nlme_3.1-164 parallelly_1.37.1 lava_1.7.3 tidyselect_1.2.0
## [45] digest_0.6.34 stringi_1.8.3 future_1.33.1 dplyr_1.1.4
## [49] reshape2_1.4.4 purrr_1.0.2 listenv_0.9.1 RWekajars_3.9.3-2
## [53] splines_4.3.3 fastmap_1.1.1 grid_4.3.3 colorspace_2.1-0
## [57] cli_3.6.2 magrittr_2.0.3 randomForest_4.7-1.1 survival_3.5-8
## [61] utf8_1.2.4 future.apply_1.11.1 withr_2.5.0 scales_1.3.0
## [65] lubridate_1.9.3 timechange_0.3.0 rmarkdown_2.25 globals_0.16.2
## [69] nnet_7.3-19 RWeka_0.4-50 timeDate_4032.109 ranger_0.16.0
## [73] evaluate_0.23 knitr_1.45 hardhat_1.3.1 rlang_1.1.3
## [77] Rcpp_1.0.12 glue_1.7.0 pROC_1.18.5 ipred_0.9-14
## [81] jsonlite_1.8.8 R6_2.5.1 plyr_1.8.9