Data624-Homework 4
2026-02-24
Homework 4
3.1 and 3.2 in the Kuhn and Johnson book Applied Predictive Modeling.
Loading packages
library(ggplot2)
library(tidyr)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union3.1 The UC Irvine Machine Learning Repository6 contains a data set related to glass identification. The data consist of 214 glass samples labeled as one of seven class categories. There are nine predictors, including the refractive index and percentages of eight elements: Na, Mg, Al, Si, K, Ca, Ba, and Fe
Loading Data
library(mlbench)## Warning: package 'mlbench' was built under R version 4.4.3data(Glass)
str(Glass) ## 'data.frame': 214 obs. of 10 variables:
## $ RI : num 1.52 1.52 1.52 1.52 1.52 ...
## $ Na : num 13.6 13.9 13.5 13.2 13.3 ...
## $ Mg : num 4.49 3.6 3.55 3.69 3.62 3.61 3.6 3.61 3.58 3.6 ...
## $ Al : num 1.1 1.36 1.54 1.29 1.24 1.62 1.14 1.05 1.37 1.36 ...
## $ Si : num 71.8 72.7 73 72.6 73.1 ...
## $ K : num 0.06 0.48 0.39 0.57 0.55 0.64 0.58 0.57 0.56 0.57 ...
## $ Ca : num 8.75 7.83 7.78 8.22 8.07 8.07 8.17 8.24 8.3 8.4 ...
## $ Ba : num 0 0 0 0 0 0 0 0 0 0 ...
## $ Fe : num 0 0 0 0 0 0.26 0 0 0 0.11 ...
## $ Type: Factor w/ 6 levels "1","2","3","5",..: 1 1 1 1 1 1 1 1 1 1 ...(a) Using visualizations, explore tok ahe predictor variables to understand their distributions as well as the relationships between predictors.
names(Glass)## [1] "RI" "Na" "Mg" "Al" "Si" "K" "Ca" "Ba" "Fe" "Type"summary(Glass)## RI Na Mg Al
## Min. :1.511 Min. :10.73 Min. :0.000 Min. :0.290
## 1st Qu.:1.517 1st Qu.:12.91 1st Qu.:2.115 1st Qu.:1.190
## Median :1.518 Median :13.30 Median :3.480 Median :1.360
## Mean :1.518 Mean :13.41 Mean :2.685 Mean :1.445
## 3rd Qu.:1.519 3rd Qu.:13.82 3rd Qu.:3.600 3rd Qu.:1.630
## Max. :1.534 Max. :17.38 Max. :4.490 Max. :3.500
## Si K Ca Ba
## Min. :69.81 Min. :0.0000 Min. : 5.430 Min. :0.000
## 1st Qu.:72.28 1st Qu.:0.1225 1st Qu.: 8.240 1st Qu.:0.000
## Median :72.79 Median :0.5550 Median : 8.600 Median :0.000
## Mean :72.65 Mean :0.4971 Mean : 8.957 Mean :0.175
## 3rd Qu.:73.09 3rd Qu.:0.6100 3rd Qu.: 9.172 3rd Qu.:0.000
## Max. :75.41 Max. :6.2100 Max. :16.190 Max. :3.150
## Fe Type
## Min. :0.00000 1:70
## 1st Qu.:0.00000 2:76
## Median :0.00000 3:17
## Mean :0.05701 5:13
## 3rd Qu.:0.10000 6: 9
## Max. :0.51000 7:29Pivoting to long format
glass_piv <- pivot_longer(Glass,
cols = RI:Fe,
names_to = "Variable",
values_to = "Value") head(glass_piv)## # A tibble: 6 × 3
## Type Variable Value
## <fct> <chr> <dbl>
## 1 1 RI 1.52
## 2 1 Na 13.6
## 3 1 Mg 4.49
## 4 1 Al 1.1
## 5 1 Si 71.8
## 6 1 K 0.06Visualization - Histogram
ggplot(glass_piv, aes(x = Value)) +
geom_histogram(bins = 30, fill = "skyblue", color = "black") +
facet_wrap(~ Variable, scales = "free") +
theme_minimal() 
According to the above plots, there are many variations to the
distributions. The closest distributions to normal are for Al, Na, RI,
and Si, although they do show some skewness.
Visualization - Box Plot
ggplot(glass_piv, aes(x = Variable, y = Value)) +
geom_boxplot(fill = "skyblue") +
theme_light() +
facet_wrap(~Variable, scales = "free_y") +
labs(
title = "Boxplots of Glass Predictors",
x = "Predictor",
y = "Value"
)
Ba and Ca show the most amount of outliers. Al, Na, Si, and RI show the most normal distributions confirming earler plots. To best explore the real relationships I will create a heatmap.
Visualization - Correlation Heatmap
Libararies
library(reshape2)##
## Attaching package: 'reshape2'## The following object is masked from 'package:tidyr':
##
## smithsCorrelation matrix
predictor_cols <- c("RI", "Na", "Mg", "Al", "Si", "K", "Ca", "Ba", "Fe")
corr_matrix <- cor(Glass[, predictor_cols])
corr_gg <- melt(corr_matrix)head(corr_gg)## Var1 Var2 value
## 1 RI RI 1.0000000
## 2 Na RI -0.1918854
## 3 Mg RI -0.1222740
## 4 Al RI -0.4073260
## 5 Si RI -0.5420522
## 6 K RI -0.2898327Heatmap
ggplot(corr_gg, aes(Var1, Var2, fill = value)) +
geom_tile() +
scale_fill_gradient2(low = "purple", mid = "white", high = "red") +
labs(
x = "Predictor",
y = "Predictor",
fill = " Correlation"
)
RI and Si have a strong negative correlation, which means that as the
silicon content is increasing, the refractive index is decreasing.
Calcium and refractive index, on the other hand, have a strong positive
correlation indicated in the heat mapp.
(b) Do there appear to be any outliers in the data? Are any predictors skewed?
Al and Ca are moderately right skewed. Ba, Fe, and K are strongly right skewed. Mg is strongly left skewed. The rest are relatively following the shape of a normal bell curve. Ba, Fe, K, and Ca have the most outliers according to the box plot.
(c) Are there any relevant transformations of one or more predictors that might improve the classification model? I would need to work on the ones that are strongly skewed and have lots of zeros, so Ba, Fe, and K.
Transformations
Glass$Ba_log <- log(Glass$Ba + 1)
Glass$Fe_log <- log(Glass$Fe + 1)
Glass$K_log <- log(Glass$K + 1)Plotting to see if it helped
library(gridExtra)##
## Attaching package: 'gridExtra'## The following object is masked from 'package:dplyr':
##
## combineba <- ggplot(Glass, aes(x = Ba)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="Ba"
) +
theme_light()
ba_log <- ggplot(Glass, aes(x = Ba_log)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="Ba After Log Transformation"
) +
theme_light()
fe <- ggplot(Glass, aes(x = Fe)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="Fe"
) +
theme_light()
fe_log <- ggplot(Glass, aes(x = Fe_log)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="Fe After Log Transformation"
) +
theme_light()
k <- ggplot(Glass, aes(x = K)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="K"
) +
theme_light()
k_log <- ggplot(Glass, aes(x = K_log)) +
geom_histogram(fill = "skyblue", color = "black") +
labs(
title ="K After Log Transformation"
) +
theme_light()
grid.arrange(ba, ba_log, fe, fe_log, k, k_log) ## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
It looks like the log transformation here is not effective at all and the distributions of these predictors are inherently not following any normal distribution no matter the adjustments. Most values are just zero and that could just be the way the data is and that they are inherently skewed. These values just represent the minority of the glass samples that contain a meaningful amount of either Ba, Fe, or K.
The soybean data can also be found at the UC Irvine Machine Learning Repository. Data were collected to predict disease in 683 soybeans. The 35 predictors are mostly categorical and include information on the environmental conditions (e.g., temperature, precipitation) and plant conditions (e.g., left spots, mold growth). The outcome labels consist of 19 distinct classes.
Loading Libraries
library(mlbench)
data(Soybean)(a) Investigate the frequency distributions for the categorical predictors. Are any of the distributions degenerate in the ways discussed earlier in this chapter?
str(Soybean)## 'data.frame': 683 obs. of 36 variables:
## $ Class : Factor w/ 19 levels "2-4-d-injury",..: 11 11 11 11 11 11 11 11 11 11 ...
## $ date : Factor w/ 7 levels "0","1","2","3",..: 7 5 4 4 7 6 6 5 7 5 ...
## $ plant.stand : Ord.factor w/ 2 levels "0"<"1": 1 1 1 1 1 1 1 1 1 1 ...
## $ precip : Ord.factor w/ 3 levels "0"<"1"<"2": 3 3 3 3 3 3 3 3 3 3 ...
## $ temp : Ord.factor w/ 3 levels "0"<"1"<"2": 2 2 2 2 2 2 2 2 2 2 ...
## $ hail : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 2 1 1 ...
## $ crop.hist : Factor w/ 4 levels "0","1","2","3": 2 3 2 2 3 4 3 2 4 3 ...
## $ area.dam : Factor w/ 4 levels "0","1","2","3": 2 1 1 1 1 1 1 1 1 1 ...
## $ sever : Factor w/ 3 levels "0","1","2": 2 3 3 3 2 2 2 2 2 3 ...
## $ seed.tmt : Factor w/ 3 levels "0","1","2": 1 2 2 1 1 1 2 1 2 1 ...
## $ germ : Ord.factor w/ 3 levels "0"<"1"<"2": 1 2 3 2 3 2 1 3 2 3 ...
## $ plant.growth : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ leaves : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ leaf.halo : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
## $ leaf.marg : Factor w/ 3 levels "0","1","2": 3 3 3 3 3 3 3 3 3 3 ...
## $ leaf.size : Ord.factor w/ 3 levels "0"<"1"<"2": 3 3 3 3 3 3 3 3 3 3 ...
## $ leaf.shread : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ leaf.malf : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ leaf.mild : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
## $ stem : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ lodging : Factor w/ 2 levels "0","1": 2 1 1 1 1 1 2 1 1 1 ...
## $ stem.cankers : Factor w/ 4 levels "0","1","2","3": 4 4 4 4 4 4 4 4 4 4 ...
## $ canker.lesion : Factor w/ 4 levels "0","1","2","3": 2 2 1 1 2 1 2 2 2 2 ...
## $ fruiting.bodies: Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ ext.decay : Factor w/ 3 levels "0","1","2": 2 2 2 2 2 2 2 2 2 2 ...
## $ mycelium : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ int.discolor : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...
## $ sclerotia : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ fruit.pods : Factor w/ 4 levels "0","1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ fruit.spots : Factor w/ 4 levels "0","1","2","4": 4 4 4 4 4 4 4 4 4 4 ...
## $ seed : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ mold.growth : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ seed.discolor : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ seed.size : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ shriveling : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
## $ roots : Factor w/ 3 levels "0","1","2": 1 1 1 1 1 1 1 1 1 1 ...summary(Soybean)## Class date plant.stand precip temp
## brown-spot : 92 5 :149 0 :354 0 : 74 0 : 80
## alternarialeaf-spot: 91 4 :131 1 :293 1 :112 1 :374
## frog-eye-leaf-spot : 91 3 :118 NA's: 36 2 :459 2 :199
## phytophthora-rot : 88 2 : 93 NA's: 38 NA's: 30
## anthracnose : 44 6 : 90
## brown-stem-rot : 44 (Other):101
## (Other) :233 NA's : 1
## hail crop.hist area.dam sever seed.tmt germ plant.growth
## 0 :435 0 : 65 0 :123 0 :195 0 :305 0 :165 0 :441
## 1 :127 1 :165 1 :227 1 :322 1 :222 1 :213 1 :226
## NA's:121 2 :219 2 :145 2 : 45 2 : 35 2 :193 NA's: 16
## 3 :218 3 :187 NA's:121 NA's:121 NA's:112
## NA's: 16 NA's: 1
##
##
## leaves leaf.halo leaf.marg leaf.size leaf.shread leaf.malf leaf.mild
## 0: 77 0 :221 0 :357 0 : 51 0 :487 0 :554 0 :535
## 1:606 1 : 36 1 : 21 1 :327 1 : 96 1 : 45 1 : 20
## 2 :342 2 :221 2 :221 NA's:100 NA's: 84 2 : 20
## NA's: 84 NA's: 84 NA's: 84 NA's:108
##
##
##
## stem lodging stem.cankers canker.lesion fruiting.bodies ext.decay
## 0 :296 0 :520 0 :379 0 :320 0 :473 0 :497
## 1 :371 1 : 42 1 : 39 1 : 83 1 :104 1 :135
## NA's: 16 NA's:121 2 : 36 2 :177 NA's:106 2 : 13
## 3 :191 3 : 65 NA's: 38
## NA's: 38 NA's: 38
##
##
## mycelium int.discolor sclerotia fruit.pods fruit.spots seed
## 0 :639 0 :581 0 :625 0 :407 0 :345 0 :476
## 1 : 6 1 : 44 1 : 20 1 :130 1 : 75 1 :115
## NA's: 38 2 : 20 NA's: 38 2 : 14 2 : 57 NA's: 92
## NA's: 38 3 : 48 4 :100
## NA's: 84 NA's:106
##
##
## mold.growth seed.discolor seed.size shriveling roots
## 0 :524 0 :513 0 :532 0 :539 0 :551
## 1 : 67 1 : 64 1 : 59 1 : 38 1 : 86
## NA's: 92 NA's:106 NA's: 92 NA's:106 2 : 15
## NA's: 31
##
##
## library(caret)## Loading required package: latticenzv <- nearZeroVar(Soybean, saveMetrics = TRUE)
nzv## freqRatio percentUnique zeroVar nzv
## Class 1.010989 2.7818448 FALSE FALSE
## date 1.137405 1.0248902 FALSE FALSE
## plant.stand 1.208191 0.2928258 FALSE FALSE
## precip 4.098214 0.4392387 FALSE FALSE
## temp 1.879397 0.4392387 FALSE FALSE
## hail 3.425197 0.2928258 FALSE FALSE
## crop.hist 1.004587 0.5856515 FALSE FALSE
## area.dam 1.213904 0.5856515 FALSE FALSE
## sever 1.651282 0.4392387 FALSE FALSE
## seed.tmt 1.373874 0.4392387 FALSE FALSE
## germ 1.103627 0.4392387 FALSE FALSE
## plant.growth 1.951327 0.2928258 FALSE FALSE
## leaves 7.870130 0.2928258 FALSE FALSE
## leaf.halo 1.547511 0.4392387 FALSE FALSE
## leaf.marg 1.615385 0.4392387 FALSE FALSE
## leaf.size 1.479638 0.4392387 FALSE FALSE
## leaf.shread 5.072917 0.2928258 FALSE FALSE
## leaf.malf 12.311111 0.2928258 FALSE FALSE
## leaf.mild 26.750000 0.4392387 FALSE TRUE
## stem 1.253378 0.2928258 FALSE FALSE
## lodging 12.380952 0.2928258 FALSE FALSE
## stem.cankers 1.984293 0.5856515 FALSE FALSE
## canker.lesion 1.807910 0.5856515 FALSE FALSE
## fruiting.bodies 4.548077 0.2928258 FALSE FALSE
## ext.decay 3.681481 0.4392387 FALSE FALSE
## mycelium 106.500000 0.2928258 FALSE TRUE
## int.discolor 13.204545 0.4392387 FALSE FALSE
## sclerotia 31.250000 0.2928258 FALSE TRUE
## fruit.pods 3.130769 0.5856515 FALSE FALSE
## fruit.spots 3.450000 0.5856515 FALSE FALSE
## seed 4.139130 0.2928258 FALSE FALSE
## mold.growth 7.820896 0.2928258 FALSE FALSE
## seed.discolor 8.015625 0.2928258 FALSE FALSE
## seed.size 9.016949 0.2928258 FALSE FALSE
## shriveling 14.184211 0.2928258 FALSE FALSE
## roots 6.406977 0.4392387 FALSE FALSEFinding the degenerate classes
nzv[nzv$nzv == TRUE, ]## freqRatio percentUnique zeroVar nzv
## leaf.mild 26.75 0.4392387 FALSE TRUE
## mycelium 106.50 0.2928258 FALSE TRUE
## sclerotia 31.25 0.2928258 FALSE TRUEThis proves that leaf.mild, mycelium, and sclerotia are the degenrate classes so they have mostly zeros and not many ones. These three predictors a pretty much dominated by just one value (0) and can make it hard to model.
(b) Roughly 18% of the data are missing. Are there particular predictors that are more likely to be missing? Is the pattern of missing data related to the classes?
Finding Missing Values
missing_values <- colSums(is.na(Soybean))
missing_values## Class date plant.stand precip temp
## 0 1 36 38 30
## hail crop.hist area.dam sever seed.tmt
## 121 16 1 121 121
## germ plant.growth leaves leaf.halo leaf.marg
## 112 16 0 84 84
## leaf.size leaf.shread leaf.malf leaf.mild stem
## 84 100 84 108 16
## lodging stem.cankers canker.lesion fruiting.bodies ext.decay
## 121 38 38 106 38
## mycelium int.discolor sclerotia fruit.pods fruit.spots
## 38 38 38 84 106
## seed mold.growth seed.discolor seed.size shriveling
## 92 92 106 92 106
## roots
## 31miss_val <- missing_values[missing_values > 0]
miss_val## date plant.stand precip temp hail
## 1 36 38 30 121
## crop.hist area.dam sever seed.tmt germ
## 16 1 121 121 112
## plant.growth leaf.halo leaf.marg leaf.size leaf.shread
## 16 84 84 84 100
## leaf.malf leaf.mild stem lodging stem.cankers
## 84 108 16 121 38
## canker.lesion fruiting.bodies ext.decay mycelium int.discolor
## 38 106 38 38 38
## sclerotia fruit.pods fruit.spots seed mold.growth
## 38 84 106 92 92
## seed.discolor seed.size shriveling roots
## 106 92 106 31Visialization
missing_df <- data.frame(
predictor = names(miss_val),
count = miss_val
)
ggplot(missing_df, aes(x = reorder(predictor, count), y = count)) +
geom_bar(stat = "identity", fill = "skyblue", color = "black") +
coord_flip() +
labs(title = "Missing Values By Predictor",
x = "Predictor",
y = "Count Missing") +
theme_light() 
Yes, there are predictors that are likely to be missing values the
classes sever, seed.tmt, lodging, hail. germ, leaf.mild as the top ones.
I also see that for a lot of them they have the same missing count,
which makes me think that a lot of the missing observations are tied to
specific disease classes. Maybe if some diseases affect certain parts of
the plant but not others, so then the observer might just record the
absence of effects as missing and not zero.
Checking if there is a class connection
Soybean |>
group_by(Class) |>
summarise(missing = sum(is.na(across(everything())))) |>
arrange(desc(missing))## # A tibble: 19 × 2
## Class missing
## <fct> <int>
## 1 phytophthora-rot 1214
## 2 2-4-d-injury 450
## 3 cyst-nematode 336
## 4 diaporthe-pod-&-stem-blight 177
## 5 herbicide-injury 160
## 6 alternarialeaf-spot 0
## 7 anthracnose 0
## 8 bacterial-blight 0
## 9 bacterial-pustule 0
## 10 brown-spot 0
## 11 brown-stem-rot 0
## 12 charcoal-rot 0
## 13 diaporthe-stem-canker 0
## 14 downy-mildew 0
## 15 frog-eye-leaf-spot 0
## 16 phyllosticta-leaf-spot 0
## 17 powdery-mildew 0
## 18 purple-seed-stain 0
## 19 rhizoctonia-root-rot 0Phytophthora has the most missing cases, which can suggest that it could be harder to asses and determine the presence of the disease. The missing data is clearly strongly related to specific disease classes. There are many disease classes that have no missing data at all, which makes me believe that the missing values are not random at all and are definitely tied to the class.
(c) Develop a strategy for handling missing data, either by eliminating predictors or imputation.
I will encode the missing data because I feel like removing all missing values will remove too much of the data.
I will encode the missing data into it’s own category.
soybean_encoded <- Soybean
for(col in names(soybean_encoded)) {
if(is.factor(soybean_encoded[[col]])) {
soybean_encoded[[col]] <- addNA(soybean_encoded[[col]])
}
}levels(soybean_encoded$hail)## [1] "0" "1" NAsummary(soybean_encoded$hail)## 0 1 <NA>
## 435 127 121Since the missing values are systematic, I added a new category level for that value because I think it is important to note that there are missing values in the first place. It could indicate a pattern and that is valuable information, especially that it depends on specific disease classes.