DATA 624: APM Exercises 3.1 and 3.2
Madina Kudanova
Exercise 3.1 - Glass data
data(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 ...x <- Glass[, -10] # predictors only (Type is column 10)
long <- pivot_longer(x, everything())(a) Distributions and relationships between predictors
ggplot(long, aes(value)) + geom_histogram(bins = 30) + facet_wrap(~name, scales = "free")
round(cor(x), 2)## RI Na Mg Al Si K Ca Ba Fe
## RI 1.00 -0.19 -0.12 -0.41 -0.54 -0.29 0.81 0.00 0.14
## Na -0.19 1.00 -0.27 0.16 -0.07 -0.27 -0.28 0.33 -0.24
## Mg -0.12 -0.27 1.00 -0.48 -0.17 0.01 -0.44 -0.49 0.08
## Al -0.41 0.16 -0.48 1.00 -0.01 0.33 -0.26 0.48 -0.07
## Si -0.54 -0.07 -0.17 -0.01 1.00 -0.19 -0.21 -0.10 -0.09
## K -0.29 -0.27 0.01 0.33 -0.19 1.00 -0.32 -0.04 -0.01
## Ca 0.81 -0.28 -0.44 -0.26 -0.21 -0.32 1.00 -0.11 0.12
## Ba 0.00 0.33 -0.49 0.48 -0.10 -0.04 -0.11 1.00 -0.06
## Fe 0.14 -0.24 0.08 -0.07 -0.09 -0.01 0.12 -0.06 1.00pairs(x)
Findings: (a): Ba (barium) and Fe (iron) are mostly zeros with a few large values, K (potassium) is concentrated below 1 with a few large values, and Mg (magnesium) is bimodal (a group at 0, the rest near 3.5). Ca (calcium), RI (refractive index), and Al (aluminum) have right tails, while Na (sodium) is the most symmetric. RI and Ca are strongly correlated (r = 0.81); the other pairs are weak to moderate.
(b) Outliers and skewness
ggplot(long, aes(y = value)) + geom_boxplot() + facet_wrap(~name, scales = "free")
sort(sapply(x, skew), decreasing = TRUE)## K Ba Ca Fe RI Al Na
## 6.4600889 3.3686800 2.0184463 1.7298107 1.6027151 0.8946104 0.4478343
## Si Mg
## -0.7202392 -1.1364523Findings: (b): Outliers are present in every predictor except Mg (magnesium). K (potassium), Ba (barium), and Fe (iron) have high outliers only, while Al (aluminum), Ca (calcium), Na (sodium), RI (refractive index), and Si (silicon) have outliers on both sides. K, Ba, Ca, Fe, and RI are strongly right-skewed (skewness 1.60 to 6.46), Al is moderately right-skewed (0.89), Na is roughly symmetric (0.45), and Si (-0.72) and Mg (-1.14) are left-skewed.
(c) Transformations
# skewness before and after log1p (log1p handles zeros)
round(cbind(raw = sapply(x, skew), log1p = sapply(log1p(x), skew)), 2)## raw log1p
## RI 1.60 1.60
## Na 0.45 0.06
## Mg -1.14 -1.31
## Al 0.89 0.00
## Si -0.72 -0.79
## K 6.46 1.95
## Ca 2.02 1.15
## Ba 3.37 2.69
## Fe 1.73 1.56# RI check: is log1p nearly a straight line over RI's range?
cor(x$RI, log1p(x$RI))## [1] 0.9999992Findings: (c): log1p helped the right-skewed predictors, though not all of them enough. Al (aluminum) went from 0.89 to 0.00, K (potassium) from 6.46 to 1.95, Ca (calcium) from 2.02 to 1.15, and Na (sodium) from 0.45 to 0.06. Ba (barium) and Fe (iron) improved only slightly (3.37 to 2.69 and 1.73 to 1.56) because of their many zeros. RI (refractive index) did not change (1.60) because its values are all close to 1.5, so log1p is nearly linear there (r = 0.9999992). Mg (magnesium) and Si (silicon) are left-skewed, so a log transform made them slightly worse (-1.14 to -1.31 and -0.72 to -0.79). Other options would be a Box-Cox or Yeo-Johnson transformation for the right-skewed predictors, a power greater than 1 for the left-skewed ones, and recoding the zero-heavy Ba, Fe, and Mg as zero/non-zero indicators.
Exercise 3.2 - Soybean data
data(Soybean)
x <- Soybean[, -1] # predictors only (Class is column 1)(a) Frequency distributions and degenerate predictors
ggplot(pivot_longer(data.frame(lapply(x, as.character)), everything()), aes(value)) +
geom_bar() + facet_wrap(~name, scales = "free")
freq_ratio <- sapply(x, function(v) { t <- sort(table(v), decreasing = TRUE); t[[1]] / t[[2]] })
pct_unique <- sapply(x, function(v) { v <- v[!is.na(v)]; length(unique(v)) / length(v) * 100 })
degenerate <- names(x)[freq_ratio > 20 & pct_unique < 10]
degenerate## [1] "leaf.mild" "mycelium" "sclerotia"Findings: (a):Most of the 35 predictors are categorical with 2 to 7 levels, and many are dominated by one level. Three predictors are degenerate: mycelium (frequency ratio 106.5), sclerotia (31.3), and leaf.mild (26.8). Each has a ratio above 20 and under 1% unique values, so nearly all observations fall in one level. The next highest, shriveling (14.2) and int.discolor (13.2), are skewed toward one level but fall below the cutoff. The three degenerate predictors are candidates for removal before modeling.
(b) Missing data by predictor and by class
sort(colSums(is.na(x)), decreasing = TRUE)## hail sever seed.tmt lodging germ
## 121 121 121 121 112
## leaf.mild fruiting.bodies fruit.spots seed.discolor shriveling
## 108 106 106 106 106
## leaf.shread seed mold.growth seed.size leaf.halo
## 100 92 92 92 84
## leaf.marg leaf.size leaf.malf fruit.pods precip
## 84 84 84 84 38
## stem.cankers canker.lesion ext.decay mycelium int.discolor
## 38 38 38 38 38
## sclerotia plant.stand roots temp crop.hist
## 38 36 31 30 16
## plant.growth stem date area.dam leaves
## 16 16 1 1 0sort(tapply(rowSums(is.na(x)), Soybean$Class, mean), decreasing = TRUE)## 2-4-d-injury cyst-nematode
## 28.12500 24.00000
## herbicide-injury phytophthora-rot
## 20.00000 13.79545
## diaporthe-pod-&-stem-blight alternarialeaf-spot
## 11.80000 0.00000
## anthracnose bacterial-blight
## 0.00000 0.00000
## bacterial-pustule brown-spot
## 0.00000 0.00000
## brown-stem-rot charcoal-rot
## 0.00000 0.00000
## diaporthe-stem-canker downy-mildew
## 0.00000 0.00000
## frog-eye-leaf-spot phyllosticta-leaf-spot
## 0.00000 0.00000
## powdery-mildew purple-seed-stain
## 0.00000 0.00000
## rhizoctonia-root-rot
## 0.00000Findings: (b): Missing values are concentrated in certain predictors. hail, sever, seed.tmt, and lodging are each missing 121 of 683 observations (about 17.7%), followed by germ (112), leaf.mild (108), and fruiting.bodies, fruit.spots, seed.discolor, and shriveling (106 each). leaves has none, and area.dam and date have only 1. Missingness is also strongly related to class. Only 5 of the 19 classes have any missing values: 2-4-d-injury (an average of about 28 of the 35 predictors missing per observation), cyst-nematode (24), herbicide-injury (20), phytophthora-rot (13.8), and diaporthe-pod-&-stem-blight (11.8). The other 14 classes have no missing values at all, so the missing data are not random.
(c) Strategy for missing data
x2 <- x[, !names(x) %in% degenerate]
x2[] <- lapply(x2, function(f) { f <- addNA(f); levels(f)[is.na(levels(f))] <- "missing"; f })
colSums(is.na(x2))## date plant.stand precip temp hail
## 0 0 0 0 0
## crop.hist area.dam sever seed.tmt germ
## 0 0 0 0 0
## plant.growth leaves leaf.halo leaf.marg leaf.size
## 0 0 0 0 0
## leaf.shread leaf.malf stem lodging stem.cankers
## 0 0 0 0 0
## canker.lesion fruiting.bodies ext.decay int.discolor fruit.pods
## 0 0 0 0 0
## fruit.spots seed mold.growth seed.discolor seed.size
## 0 0 0 0 0
## shriveling roots
## 0 0Strategy: I dropped the three predictors that are almost always the same value (leaf.mild, mycelium, sclerotia), because they tell us almost nothing. For the other 32 predictors, I made “missing” its own category instead of filling in guesses. Whether a value is missing is linked to the disease class (only 5 of the 19 classes have any missing values), so the missingness itself is useful information. Filling in the most common value would hide that pattern and could make those plants look like the more common diseases. After this change, no missing values are left.