Wine Quality Predictor
Wine Quality Predictor
INTRODUCTION
Team Member : Shradhit Subudhi
Professor’s Name : Dr. Douglas Jones
Course Period : Fall - 2018 (September - December)
Date : “December 3, 2018”
Inspiration for Project
“In wine, there’s truth” - Pliny the Elder.
Question
Can we predict the quality of wine using data analysis techniques?
Wine Data Set - Info
The two datasets are related to red and white variants of the Portuguese “Vinho Verde” wine. Due to privacy and logistic issues, only physicochemical (inputs) and sensory (the output) variables are available (e.g. there is no data about grape types, wine brand, wine selling price, etc.).
These datasets can be viewed as classification or regression tasks. The classes are ordered and not balanced (e.g. there are much more normal wines than excellent or poor ones). Outlier detection algorithms could be used to detect the few excellent or poor wines. Also, we are not sure if all input variables are relevant. So it could be interesting to test feature selection methods.
We worked on predicting the white wine quality for this project. Worst wine is represented with 0 and the best wine is represented as 10.
Packages Installed
library("rmdformats")
library("corrgram")
library("MASS")
library("ggplot2")
library("naniar")
library("e1071")
library("lattice")##
## Attaching package: 'lattice'## The following object is masked from 'package:corrgram':
##
## panel.filllibrary("caret")
library("car")## Loading required package: carDatalibrary("caTools")
library("corrplot")## corrplot 0.84 loadedlibrary("knitr")Inputting Data
Importing the data with the help of a url and assigning it to white.
white_url <- "https://archive.ics.uci.edu/ml/machine-learning-databases/wine-quality/winequality-white.csv"
whitewine <- read.csv(white_url, header = TRUE, sep = ";")
white <- whitewineData - First View
Viewing the loaded data.
str(white)## 'data.frame': 4898 obs. of 12 variables:
## $ fixed.acidity : num 7 6.3 8.1 7.2 7.2 8.1 6.2 7 6.3 8.1 ...
## $ volatile.acidity : num 0.27 0.3 0.28 0.23 0.23 0.28 0.32 0.27 0.3 0.22 ...
## $ citric.acid : num 0.36 0.34 0.4 0.32 0.32 0.4 0.16 0.36 0.34 0.43 ...
## $ residual.sugar : num 20.7 1.6 6.9 8.5 8.5 6.9 7 20.7 1.6 1.5 ...
## $ chlorides : num 0.045 0.049 0.05 0.058 0.058 0.05 0.045 0.045 0.049 0.044 ...
## $ free.sulfur.dioxide : num 45 14 30 47 47 30 30 45 14 28 ...
## $ total.sulfur.dioxide: num 170 132 97 186 186 97 136 170 132 129 ...
## $ density : num 1.001 0.994 0.995 0.996 0.996 ...
## $ pH : num 3 3.3 3.26 3.19 3.19 3.26 3.18 3 3.3 3.22 ...
## $ sulphates : num 0.45 0.49 0.44 0.4 0.4 0.44 0.47 0.45 0.49 0.45 ...
## $ alcohol : num 8.8 9.5 10.1 9.9 9.9 10.1 9.6 8.8 9.5 11 ...
## $ quality : int 6 6 6 6 6 6 6 6 6 6 ...Looking at the output we get to know that we have 4898 samples and 12 variables.
The response variable is quality.
The eleven predictor variables are of the numeric class and the response variable, quality, is of the integer class.
We don’t see any categorical variable in the dataset.
Columns
We then check the names of the variables present in the data.
colnames(white)## [1] "fixed.acidity" "volatile.acidity" "citric.acid"
## [4] "residual.sugar" "chlorides" "free.sulfur.dioxide"
## [7] "total.sulfur.dioxide" "density" "pH"
## [10] "sulphates" "alcohol" "quality"Look at the summary to better understand the dataset.
summary(white)## fixed.acidity volatile.acidity citric.acid residual.sugar
## Min. : 3.800 Min. :0.0800 Min. :0.0000 Min. : 0.600
## 1st Qu.: 6.300 1st Qu.:0.2100 1st Qu.:0.2700 1st Qu.: 1.700
## Median : 6.800 Median :0.2600 Median :0.3200 Median : 5.200
## Mean : 6.855 Mean :0.2782 Mean :0.3342 Mean : 6.391
## 3rd Qu.: 7.300 3rd Qu.:0.3200 3rd Qu.:0.3900 3rd Qu.: 9.900
## Max. :14.200 Max. :1.1000 Max. :1.6600 Max. :65.800
## chlorides free.sulfur.dioxide total.sulfur.dioxide
## Min. :0.00900 Min. : 2.00 Min. : 9.0
## 1st Qu.:0.03600 1st Qu.: 23.00 1st Qu.:108.0
## Median :0.04300 Median : 34.00 Median :134.0
## Mean :0.04577 Mean : 35.31 Mean :138.4
## 3rd Qu.:0.05000 3rd Qu.: 46.00 3rd Qu.:167.0
## Max. :0.34600 Max. :289.00 Max. :440.0
## density pH sulphates alcohol
## Min. :0.9871 Min. :2.720 Min. :0.2200 Min. : 8.00
## 1st Qu.:0.9917 1st Qu.:3.090 1st Qu.:0.4100 1st Qu.: 9.50
## Median :0.9937 Median :3.180 Median :0.4700 Median :10.40
## Mean :0.9940 Mean :3.188 Mean :0.4898 Mean :10.51
## 3rd Qu.:0.9961 3rd Qu.:3.280 3rd Qu.:0.5500 3rd Qu.:11.40
## Max. :1.0390 Max. :3.820 Max. :1.0800 Max. :14.20
## quality
## Min. :3.000
## 1st Qu.:5.000
## Median :6.000
## Mean :5.878
## 3rd Qu.:6.000
## Max. :9.000Removing the Duplicate Rows
#white <- unique(white)
white <- white[!duplicated(white), ]
dim(white)## [1] 3961 12The rows got reduced to 3961 after removing the duplicates.
Checking for NA
Graphically
Checking if there are any missing (NA) values.
vis_miss(white)
100% values of all columns are present.
Technically
sum(is.na(white))## [1] 0There is no NA values in the data set.
Response Count
table(white$quality)##
## 3 4 5 6 7 8 9
## 20 153 1175 1788 689 131 5We can see that the class is skewed i.e. it has imbalance.
There are 3961 samples but only 20 points are from 3rd class and only 5 points are from the 9th class.
Wines with lowest and highest quality are rare. Generally, wines have quality 5 and 6 which is in the medium range(not low, not high).
#rounds the values in its first argument to the specified number of decimal places (default 0).
round(cor(white, method = "pearson"), 2)## fixed.acidity volatile.acidity citric.acid
## fixed.acidity 1.00 -0.02 0.30
## volatile.acidity -0.02 1.00 -0.16
## citric.acid 0.30 -0.16 1.00
## residual.sugar 0.08 0.10 0.11
## chlorides 0.02 0.09 0.13
## free.sulfur.dioxide -0.06 -0.10 0.09
## total.sulfur.dioxide 0.08 0.10 0.12
## density 0.27 0.06 0.16
## pH -0.43 -0.05 -0.18
## sulphates -0.02 -0.02 0.05
## alcohol -0.11 0.05 -0.08
## quality -0.12 -0.19 0.01
## residual.sugar chlorides free.sulfur.dioxide
## fixed.acidity 0.08 0.02 -0.06
## volatile.acidity 0.10 0.09 -0.10
## citric.acid 0.11 0.13 0.09
## residual.sugar 1.00 0.08 0.31
## chlorides 0.08 1.00 0.10
## free.sulfur.dioxide 0.31 0.10 1.00
## total.sulfur.dioxide 0.41 0.19 0.62
## density 0.82 0.25 0.29
## pH -0.17 -0.09 -0.01
## sulphates -0.02 0.02 0.04
## alcohol -0.40 -0.36 -0.25
## quality -0.12 -0.22 0.01
## total.sulfur.dioxide density pH sulphates alcohol
## fixed.acidity 0.08 0.27 -0.43 -0.02 -0.11
## volatile.acidity 0.10 0.06 -0.05 -0.02 0.05
## citric.acid 0.12 0.16 -0.18 0.05 -0.08
## residual.sugar 0.41 0.82 -0.17 -0.02 -0.40
## chlorides 0.19 0.25 -0.09 0.02 -0.36
## free.sulfur.dioxide 0.62 0.29 -0.01 0.04 -0.25
## total.sulfur.dioxide 1.00 0.54 0.01 0.14 -0.45
## density 0.54 1.00 -0.06 0.08 -0.76
## pH 0.01 -0.06 1.00 0.14 0.09
## sulphates 0.14 0.08 0.14 1.00 -0.02
## alcohol -0.45 -0.76 0.09 -0.02 1.00
## quality -0.18 -0.34 0.12 0.05 0.46
## quality
## fixed.acidity -0.12
## volatile.acidity -0.19
## citric.acid 0.01
## residual.sugar -0.12
## chlorides -0.22
## free.sulfur.dioxide 0.01
## total.sulfur.dioxide -0.18
## density -0.34
## pH 0.12
## sulphates 0.05
## alcohol 0.46
## quality 1.00Correlation Matrix
We now check the correlation matrix, confidence interval with the help of corrplot package.It also helps us to do matrix reordering.
corrplot(cor(white))
corrgram(white, type="data", lower.panel=panel.conf,
upper.panel=panel.shade, main= "Corrgram for wine quality dataset", order=T, cex.labels=1.2)
Sugar,pH and citric acid factors are not playing a really big role in the quality of wine.
alcohol is kind of strongly correlated to the quality of wine.
Histogram
Plotting the histograms using hist() for the data.
attach(white)
par(mfrow=c(2,2), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
barplot((table(quality)), col=c("slateblue4", "slategray", "slategray1", "slategray2", "slategray3", "skyblue4"))
mtext("Quality", side=1, outer=F, line=2, cex=0.8)
truehist(fixed.acidity, h = 0.5, col="slategray3")
mtext("Fixed Acidity", side=1, outer=F, line=2, cex=0.8)
truehist(volatile.acidity, h = 0.05, col="slategray3")
mtext("Volatile Acidity", side=1, outer=F, line=2, cex=0.8)
truehist(citric.acid, h = 0.1, col="slategray3")
mtext("Citric Acid", side=1, outer=F, line=2, cex=0.8)
par(mfrow=c(2,2), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
truehist(residual.sugar, h = 5, col="slategray3")
mtext("Residual Sugar", side=1, outer=F, line=2, cex=0.8)
truehist(chlorides, h = 0.01, col="slategray3")
mtext("Chloride", side=1, outer=F, line=2, cex=0.8)
truehist(alcohol, h = 0.5, col="slategray3")
mtext("Alcohol", side=1, outer=F, line=2, cex=0.8)
truehist(density, h = 0.005, col="slategray3")
mtext("Density", side=1, outer=F, line=2, cex=0.8)
par(mfrow=c(2,2), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
truehist(free.sulfur.dioxide, h = 10, col="slategray3")
mtext("Free Sulfur Dioxide", side=1, outer=F, line=2, cex=0.8)
truehist(pH, h = 0.1, col="slategray3")
mtext("pH values", side=1, outer=F, line=2, cex=0.8)
truehist(sulphates, h = 0.1, col="slategray3")
mtext("sulphates", side=1, outer=F, line=2, cex=0.8)
truehist(total.sulfur.dioxide, h = 20, col="slategray3")
mtext("total.sulfur.dioxide", side=1, outer=F, line=2, cex=0.8)
Boxplot
Y ^(lambda) where Y would be the values in the dataset and lambda would be the exponent (ranges from -5 to 5).
Examples
If lambda = 2, all the values would be squared.
If lambda = 0.5, the square root of the values would be taken.
If lambda = 0, the log of the values would be taken.
Plotting the boxplots using boxplot() for the data.
par(mfrow=c(1,5), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
boxplot(fixed.acidity, col="slategray2", pch=19)
mtext("Fixed Acidity", cex=0.8, side=1, line=2)
boxplot(volatile.acidity, col="slategray2", pch=19)
mtext("Volatile Acidity", cex=0.8, side=1, line=2)
boxplot(citric.acid, col="slategray2", pch=19)
mtext("Citric Acid", cex=0.8, side=1, line=2)
boxplot(residual.sugar, col="slategray2", pch=19)
mtext("Residual Sugar", cex=0.8, side=1, line=2)
boxplot(chlorides, col="slategray2", pch=19)
mtext("Chlorides", cex=0.8, side=1, line=2)
par(mfrow=c(1,6), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
boxplot(alcohol, col="slategray2", pch=19)
mtext("Alcohol", cex=0.8, side=1, line=2)
boxplot(density, col="slategray2", pch=19)
mtext("density", cex=0.8, side=1, line=2)
boxplot(free.sulfur.dioxide, col="slategray2", pch=19)
mtext("free.sulfur.dioxide", cex=0.8, side=1, line=2)
boxplot( pH, col="slategray2", pch=19)
mtext("pH", cex=0.8, side=1, line=2)
boxplot(sulphates, col="slategray2", pch=19)
mtext("sulphates", cex=0.8, side=1, line=2)
boxplot(total.sulfur.dioxide, col="slategray2", pch=19)
mtext("total.sulfur.dioxide", cex=0.8, side=1, line=2)
detach(white)Observations regarding variables: All variables have outliers.
Quality has most values concentrated in the categories 5, 6 and 7. Only a small proportion is in the categories [3, 4] and [8, 9] and none in the categories [1, 2, 10].
Fixed acidity, volatile acidity and citric acid have outliers. If those outliers are eliminated then distribution of the variables may be considered to be symmetric.
Residual sugar has a positively skewed distribution; even after eliminating the outliers, distribution will remain skewed.
Some of the variables, e.g . free sulphur dioxide and density, have a few outliers but these are very different from the rest.Mostly outliers are on the larger side.
Alcohol has an irregular shaped distribution but it does not have pronounced outliers.
Skewness
Checking the skewness of the individual variables of the data to check whether the data is normally distributed or not.
attach(white)
skewness(quality)## [1] 0.1119192skewness(chlorides)## [1] 4.965313skewness(free.sulfur.dioxide)## [1] 1.565494skewness(residual.sugar)## [1] 1.332629skewness(alcohol)## [1] 0.4503553skewness(citric.acid)## [1] 1.309609skewness(density)## [1] 1.272354skewness(fixed.acidity)## [1] 0.6955731skewness(volatile.acidity)## [1] 1.639838skewness(total.sulfur.dioxide)## [1] 0.4564538skewness(sulphates)## [1] 0.9371431skewness(pH)## [1] 0.4551119colnames(white)## [1] "fixed.acidity" "volatile.acidity" "citric.acid"
## [4] "residual.sugar" "chlorides" "free.sulfur.dioxide"
## [7] "total.sulfur.dioxide" "density" "pH"
## [10] "sulphates" "alcohol" "quality"detach(white)Skewness Rule of Thumb
If skewness is 0, the data are perfectly symmetrical.
As a general rule of thumb: If skewness is less than -1 or greater than 1, the distribution is highly skewed.
If skewness is between -0.5 and 0.5, the distribution is approximately symmetric.
Data Transformation / Preparation
Box Cox Transformation
We use the Boxcox transformation and transform the data and then again check the skewness.
preprocess_white <- preProcess(white[,1:11], c("BoxCox", "center", "scale"))
new_white <- data.frame(trans = predict(preprocess_white, white))
colnames(new_white)## [1] "trans.fixed.acidity" "trans.volatile.acidity"
## [3] "trans.citric.acid" "trans.residual.sugar"
## [5] "trans.chlorides" "trans.free.sulfur.dioxide"
## [7] "trans.total.sulfur.dioxide" "trans.density"
## [9] "trans.pH" "trans.sulphates"
## [11] "trans.alcohol" "trans.quality"Skewness - Transformed Data
skewness(new_white$trans.quality)## [1] 0.1119192skewness(new_white$trans.chlorides)## [1] -0.1496158skewness(new_white$trans.free.sulfur.dioxide)## [1] 0.0950346skewness(new_white$trans.residual.sugar)## [1] -0.05826498skewness(new_white$trans.alcohol)## [1] 0.03929004skewness(new_white$trans.citric.acid)## [1] 1.309609skewness(new_white$trans.density)## [1] 1.099659skewness(new_white$trans.fixed.acidity)## [1] 0.1131388skewness(new_white$trans.volatile.acidity)## [1] 0.171864skewness(new_white$trans.total.sulfur.dioxide)## [1] 0.008931003skewness(new_white$trans.sulphates)## [1] -0.01873289skewness(new_white$trans.pH)## [1] 0.0002881265Removal of Outliers
Most parametric statistics, like means, standard deviations, and correlations, and every statistic based on these, are highly sensitive to outliers.The assumptions of common statistical procedures, like linear regression and ANOVA, are also based on these statistics, outliers can disrupt the analysis. Thus, we remove the outliers.
Possibly, the most important step in data preparation is to identify outliers. Since this is a multivariate data, we consider only those points which do not have any predictor variable value to be outside of limits constructed by boxplots. The following rule is applied:
A predictor value is considered to be an outlier only if it is greater than SD - 3. The rationale behind this rule is that the extreme outliers are all on the higher end of the values and the distributions are all positively skewed.
new_white <- new_white[!abs(new_white$trans.fixed.acidity) > 3,]
new_white <- new_white[!abs(new_white$trans.volatile.acidity) > 3,]
new_white <- new_white[!abs(new_white$trans.citric.acid) > 3,]
new_white <- new_white[!abs(new_white$trans.residual.sugar) > 3,]
new_white <- new_white[!abs(new_white$trans.chlorides) > 3,]
new_white <- new_white[!abs(new_white$trans.density) > 3,]
new_white <- new_white[!abs(new_white$trans.pH) > 3,]
new_white <- new_white[!abs(new_white$trans.sulphates) > 3,]
new_white <- new_white[!abs(new_white$trans.alcohol) > 3,]We now check the correlation matrix, confidence interval with the help of corrplot package.It also helps us to do matrix reordering.
corrplot(cor(new_white), type = "lower")
corrgram(new_white, type="data", lower.panel=panel.conf,
upper.panel=panel.shade, main= "Corrgram for wine quality dataset", order=T, cex.labels=1.1)
Linear Model
Assumption for linear model
Linear regression is an analysis that assesses whether one or more predictor variables explain the dependent (criterion) variable. The regression has five key assumptions:
Linear relationship
Multivariate normality
No or little multicollinearity
No auto-correlation
Homoscedasticity
Train - Test Set
set.seed(100)
trainingRowIndex <- sample(1:nrow(new_white), 0.8*nrow(new_white))
whitetrain <- new_white[trainingRowIndex, ]
whitetest <- new_white[-trainingRowIndex, ] Model Selection
Base Model - Trying to choose the best linear model..
linear_0 <- lm(trans.quality ~ . , whitetrain) # taking all the points.
summary(linear_0)##
## Call:
## lm(formula = trans.quality ~ ., data = whitetrain)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7154 -0.4737 -0.0370 0.4768 3.0699
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.86221 0.01377 425.612 < 2e-16 ***
## trans.fixed.acidity 0.06476 0.02128 3.044 0.002356 **
## trans.volatile.acidity -0.16696 0.01530 -10.909 < 2e-16 ***
## trans.citric.acid 0.06819 0.01749 3.898 9.91e-05 ***
## trans.residual.sugar 0.35957 0.03351 10.729 < 2e-16 ***
## trans.chlorides -0.06810 0.01921 -3.545 0.000398 ***
## trans.free.sulfur.dioxide 0.12545 0.01886 6.652 3.43e-11 ***
## trans.total.sulfur.dioxide -0.05308 0.02141 -2.480 0.013199 *
## trans.density -0.47660 0.05724 -8.327 < 2e-16 ***
## trans.pH 0.13383 0.01873 7.144 1.13e-12 ***
## trans.sulphates 0.07597 0.01463 5.192 2.22e-07 ***
## trans.alcohol 0.15887 0.03593 4.422 1.01e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.751 on 2988 degrees of freedom
## Multiple R-squared: 0.308, Adjusted R-squared: 0.3055
## F-statistic: 120.9 on 11 and 2988 DF, p-value: < 2.2e-16Checking Multicollinarilty
In multicollinearity,collinearity exists between three or more variables even if no pair of variables has a particularly high correlation. This means that there is redundancy between predictor variables.
In the presence of multicollinearity, the solution of the regression model becomes unstable.
Multicollinearity can assessed by computing a score called the variance inflation factor (or VIF), which measures how much the variance of a regression coefficient is inflated due to multicollinearity in the model.
The VIF of a predictor is a measure for how easily it is predicted from a linear regression using the other predictors.
# Variance inflation factor
vif(linear_0)## trans.fixed.acidity trans.volatile.acidity
## 2.230946 1.162012
## trans.citric.acid trans.residual.sugar
## 1.172273 5.903351
## trans.chlorides trans.free.sulfur.dioxide
## 1.452565 1.898081
## trans.total.sulfur.dioxide trans.density
## 2.415348 16.176165
## trans.pH trans.sulphates
## 1.710817 1.124369
## trans.alcohol
## 6.687425We can see multicollinarity over here. We remove trans.density cause it exibits multicollinarity.
Now fitting the new model.
linear_1 <- lm(trans.quality ~ . -trans.density , whitetrain)
summary(linear_1)##
## Call:
## lm(formula = trans.quality ~ . - trans.density, data = whitetrain)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.6568 -0.4768 -0.0324 0.4722 3.0679
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.86572 0.01392 421.278 < 2e-16 ***
## trans.fixed.acidity -0.04727 0.01667 -2.835 0.004607 **
## trans.volatile.acidity -0.15950 0.01545 -10.322 < 2e-16 ***
## trans.citric.acid 0.05733 0.01764 3.250 0.001167 **
## trans.residual.sugar 0.11454 0.01622 7.062 2.03e-12 ***
## trans.chlorides -0.08891 0.01926 -4.616 4.08e-06 ***
## trans.free.sulfur.dioxide 0.13794 0.01901 7.255 5.11e-13 ***
## trans.total.sulfur.dioxide -0.07259 0.02152 -3.373 0.000752 ***
## trans.pH 0.05698 0.01649 3.456 0.000556 ***
## trans.sulphates 0.05205 0.01451 3.587 0.000340 ***
## trans.alcohol 0.41616 0.01853 22.455 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7596 on 2989 degrees of freedom
## Multiple R-squared: 0.292, Adjusted R-squared: 0.2896
## F-statistic: 123.3 on 10 and 2989 DF, p-value: < 2.2e-16vif(linear_1)## trans.fixed.acidity trans.volatile.acidity
## 1.338853 1.158030
## trans.citric.acid trans.residual.sugar
## 1.165763 1.351630
## trans.chlorides trans.free.sulfur.dioxide
## 1.427982 1.886078
## trans.total.sulfur.dioxide trans.pH
## 2.386416 1.295562
## trans.sulphates trans.alcohol
## 1.081025 1.739920We fit another model.
linear_2 <- lm(trans.quality ~ . -trans.density - trans.fixed.acidity, whitetrain)
summary(linear_2)##
## Call:
## lm(formula = trans.quality ~ . - trans.density - trans.fixed.acidity,
## data = whitetrain)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.6697 -0.4717 -0.0311 0.4764 2.9648
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.86541 0.01394 420.774 < 2e-16 ***
## trans.volatile.acidity -0.15724 0.01545 -10.178 < 2e-16 ***
## trans.citric.acid 0.04561 0.01717 2.657 0.007937 **
## trans.residual.sugar 0.11487 0.01624 7.075 1.86e-12 ***
## trans.chlorides -0.09141 0.01926 -4.745 2.18e-06 ***
## trans.free.sulfur.dioxide 0.14654 0.01879 7.797 8.66e-15 ***
## trans.total.sulfur.dioxide -0.08141 0.02132 -3.819 0.000137 ***
## trans.pH 0.07510 0.01522 4.935 8.45e-07 ***
## trans.sulphates 0.05175 0.01453 3.562 0.000373 ***
## trans.alcohol 0.41542 0.01855 22.391 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7605 on 2990 degrees of freedom
## Multiple R-squared: 0.2901, Adjusted R-squared: 0.2879
## F-statistic: 135.7 on 9 and 2990 DF, p-value: < 2.2e-16We fit another model removing trans.citric.acid.
linear_4 <- lm(trans.quality ~ . -trans.density - trans.fixed.acidity - trans.citric.acid, whitetrain)
summary(linear_4)##
## Call:
## lm(formula = trans.quality ~ . - trans.density - trans.fixed.acidity -
## trans.citric.acid, data = whitetrain)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7122 -0.4640 -0.0368 0.4873 3.0119
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.86240 0.01391 421.531 < 2e-16 ***
## trans.volatile.acidity -0.16539 0.01516 -10.911 < 2e-16 ***
## trans.residual.sugar 0.11679 0.01624 7.193 7.99e-13 ***
## trans.chlorides -0.09122 0.01928 -4.731 2.34e-06 ***
## trans.free.sulfur.dioxide 0.14497 0.01880 7.710 1.70e-14 ***
## trans.total.sulfur.dioxide -0.07492 0.02120 -3.534 0.000415 ***
## trans.pH 0.06711 0.01493 4.494 7.25e-06 ***
## trans.sulphates 0.05462 0.01450 3.767 0.000169 ***
## trans.alcohol 0.41760 0.01855 22.508 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7612 on 2991 degrees of freedom
## Multiple R-squared: 0.2884, Adjusted R-squared: 0.2865
## F-statistic: 151.5 on 8 and 2991 DF, p-value: < 2.2e-16Trying to figure out the important predictor values
corrplot(cor(new_white), type = "lower")
vif(linear_4)## trans.volatile.acidity trans.residual.sugar
## 1.109472 1.348887
## trans.chlorides trans.free.sulfur.dioxide
## 1.424964 1.836333
## trans.total.sulfur.dioxide trans.pH
## 2.305832 1.058009
## trans.sulphates trans.alcohol
## 1.074973 1.736160ANOVA
Analysis of variance (ANOVA) is a collection of statistical models and their associated estimation procedures used to analyze the differences among group means in a sample. ANOVA provides a statistical test of whether the population means of several groups are equal, and therefore generalizes the t-test to more than two groups. It is useful for comparing (testing) three or more group means for statistical significance.
par(mfrow=c(2,2), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
plot(linear_4)
F - Test
Null Hypothesis - Ho : B1 = B2 = B3 = B4 = 0 Alternate Hypothesis - Ho : B1 or B2 ….B11 != 0 (Bj where j ranges from 1 to 11)
F- Statistic - p value <= 2.2e -16
We can see that p value is less than 0.05. Thus we can reject the null hypothesis and accept the alternate hypothesis. Thus,we can say that the model is significant.
Stepwise
steplinear <- step(linear_0)## Start: AIC=-1705.85
## trans.quality ~ trans.fixed.acidity + trans.volatile.acidity +
## trans.citric.acid + trans.residual.sugar + trans.chlorides +
## trans.free.sulfur.dioxide + trans.total.sulfur.dioxide +
## trans.density + trans.pH + trans.sulphates + trans.alcohol
##
## Df Sum of Sq RSS AIC
## <none> 1685.4 -1705.8
## - trans.total.sulfur.dioxide 1 3.469 1688.9 -1701.7
## - trans.fixed.acidity 1 5.226 1690.6 -1698.6
## - trans.chlorides 1 7.090 1692.5 -1695.3
## - trans.citric.acid 1 8.571 1694.0 -1692.6
## - trans.alcohol 1 11.030 1696.4 -1688.3
## - trans.sulphates 1 15.206 1700.6 -1680.9
## - trans.free.sulfur.dioxide 1 24.957 1710.3 -1663.8
## - trans.pH 1 28.788 1714.2 -1657.0
## - trans.density 1 39.107 1724.5 -1639.0
## - trans.residual.sugar 1 64.932 1750.3 -1594.4
## - trans.volatile.acidity 1 67.124 1752.5 -1590.7summary(steplinear)##
## Call:
## lm(formula = trans.quality ~ trans.fixed.acidity + trans.volatile.acidity +
## trans.citric.acid + trans.residual.sugar + trans.chlorides +
## trans.free.sulfur.dioxide + trans.total.sulfur.dioxide +
## trans.density + trans.pH + trans.sulphates + trans.alcohol,
## data = whitetrain)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7154 -0.4737 -0.0370 0.4768 3.0699
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.86221 0.01377 425.612 < 2e-16 ***
## trans.fixed.acidity 0.06476 0.02128 3.044 0.002356 **
## trans.volatile.acidity -0.16696 0.01530 -10.909 < 2e-16 ***
## trans.citric.acid 0.06819 0.01749 3.898 9.91e-05 ***
## trans.residual.sugar 0.35957 0.03351 10.729 < 2e-16 ***
## trans.chlorides -0.06810 0.01921 -3.545 0.000398 ***
## trans.free.sulfur.dioxide 0.12545 0.01886 6.652 3.43e-11 ***
## trans.total.sulfur.dioxide -0.05308 0.02141 -2.480 0.013199 *
## trans.density -0.47660 0.05724 -8.327 < 2e-16 ***
## trans.pH 0.13383 0.01873 7.144 1.13e-12 ***
## trans.sulphates 0.07597 0.01463 5.192 2.22e-07 ***
## trans.alcohol 0.15887 0.03593 4.422 1.01e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.751 on 2988 degrees of freedom
## Multiple R-squared: 0.308, Adjusted R-squared: 0.3055
## F-statistic: 120.9 on 11 and 2988 DF, p-value: < 2.2e-16ANOVA
Analysis of variance (ANOVA) is a collection of statistical models and their associated estimation procedures used to analyze the differences among group means in a sample. ANOVA provides a statistical test of whether the population means of several groups are equal, and therefore generalizes the t-test to more than two groups. It is useful for comparing (testing) three or more group means for statistical significance.
anova(linear_4,steplinear)## Analysis of Variance Table
##
## Model 1: trans.quality ~ (trans.fixed.acidity + trans.volatile.acidity +
## trans.citric.acid + trans.residual.sugar + trans.chlorides +
## trans.free.sulfur.dioxide + trans.total.sulfur.dioxide +
## trans.density + trans.pH + trans.sulphates + trans.alcohol) -
## trans.density - trans.fixed.acidity - trans.citric.acid
## Model 2: trans.quality ~ trans.fixed.acidity + trans.volatile.acidity +
## trans.citric.acid + trans.residual.sugar + trans.chlorides +
## trans.free.sulfur.dioxide + trans.total.sulfur.dioxide +
## trans.density + trans.pH + trans.sulphates + trans.alcohol
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 2991 1733.2
## 2 2988 1685.4 3 47.827 28.264 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Used the stepwise model because it has lower RSS and ANOVA tells that it’s significant.
par(mfrow=c(2,2), oma = c(1,1,0,0) + 0.1, mar = c(3,3,1,1) + 0.1)
plot(steplinear)
Residuals vs Fitted plot shows if residuals have non-linear patterns.Residuals around a horizontal line without distinct patterns, that is a good indication we don’t have non-linear relationships.
Normal QQ plot shows residuals fitting the line. Hence, can call it norlly distibuted residuals.
Scale-Location plot shows if residuals are spread equally along the ranges of predictors. This is how you can check the assumption of equal variance (homoscedasticity). It’s good if you see a horizontal line with equally (randomly) spread points.
Residuals vs Leverage plot has a typical look when there is some influential case. You can barely see Cook’s distance lines (a red dashed line) because all cases are well inside of the Cook’s distance.
Predicting - Trained Set
Predict is a generic function for predictions from the results of various model fitting functions. The function invokes particular methods which depend on the class of the first argument.
distPred <- predict(steplinear, whitetrain)
head(distPred)## 1426 1180 2658 259 2246 2323
## 5.880671 5.626071 6.517168 6.331008 5.915127 5.247180Converting the Number to a Whole Number
We use the ceiling operation for rounding the numeric values towards near integer values.
distPred1 <- ceiling(distPred)
head(distPred1)## 1426 1180 2658 259 2246 2323
## 6 6 7 7 6 6Training Data Confusion Matrix
Table function performs categorical tabulation of data with the variable and its frequency. table() function is also helpful in creating Frequency tables with condition and cross tabulations.
trn_tab <- table(predicted = distPred1, actual = whitetrain$trans.quality)
trn_tab## actual
## predicted 3 4 5 6 7 8 9
## 5 2 20 71 20 1 0 0
## 6 7 85 718 700 114 16 1
## 7 4 8 106 601 421 92 3
## 8 0 0 0 2 5 3 0Accuracy for the Linear Model
We check the accuracy of the linear model.
sum(diag(trn_tab))/length(whitetest$trans.quality)## [1] 0.26Accuracy Prediction over train set Linear Model is 26%.
Testing or Validating the Model
distPred <- predict(steplinear, whitetest)
head(distPred)## 1 4 13 25 28 30
## 4.900123 5.784360 6.037997 5.133616 5.952417 6.431214Round the numeric values to the nearest integer values.
distPred1 <- ceiling(distPred)
head(distPred1)## 1 4 13 25 28 30
## 5 6 7 6 6 7tst_tab <- table(predicted = distPred1, actual = whitetest$trans.quality)
tst_tab## actual
## predicted 3 4 5 6 7 8
## 5 1 6 17 6 1 0
## 6 1 20 166 204 20 1
## 7 1 2 25 162 103 12
## 8 0 0 0 0 2 0Accuracy - Test Data
Checking the accuracy of the test data.
sum(diag(tst_tab))/length(whitetest$trans.quality)## [1] 0.06133333Accuracy Prediction over test set Linear Model is 6.13%.
Ordinal Logistic Regression Model
Assumptions for Logistic Regression
First, binary logistic regression requires the dependent variable to be binary and ordinal logistic regression requires the dependent variable to be ordinal.
Second, logistic regression requires the observations to be independent of each other. In other words, the observations should not come from repeated measurements or matched data.
Third, logistic regression requires there to be little or no multicollinearity among the independent variables. This means that the independent variables should not be too highly correlated with each other.
Fourth, logistic regression assumes linearity of independent variables and log odds. although this analysis does not require the dependent and independent variables to be related linearly, it requires that the independent variables are linearly related to the log odds. (Did Not check for this! To improve model. We can try this in future.)
Finally, logistic regression typically requires a large sample size.
Ordinal logistic regression can be used to model a ordered factor response.
We use as.factor() to covert the variables into factors, that is, categorical variables.
#using the orignal and making changes specifically required for Logistic Regession.
white$quality2 <- as.factor(white$quality)Train - Test Set
set.seed(3000)
spl = sample.split(white$quality2, SplitRatio = 0.7)
whitetrain = subset(white, spl==TRUE)
whitetest = subset(white, spl==FALSE)
head(whitetrain)## fixed.acidity volatile.acidity citric.acid residual.sugar chlorides
## 2 6.3 0.30 0.34 1.60 0.049
## 7 6.2 0.32 0.16 7.00 0.045
## 10 8.1 0.22 0.43 1.50 0.044
## 11 8.1 0.27 0.41 1.45 0.033
## 12 8.6 0.23 0.40 4.20 0.035
## 13 7.9 0.18 0.37 1.20 0.040
## free.sulfur.dioxide total.sulfur.dioxide density pH sulphates alcohol
## 2 14 132 0.9940 3.30 0.49 9.5
## 7 30 136 0.9949 3.18 0.47 9.6
## 10 28 129 0.9938 3.22 0.45 11.0
## 11 11 63 0.9908 2.99 0.56 12.0
## 12 17 109 0.9947 3.14 0.53 9.7
## 13 16 75 0.9920 3.18 0.63 10.8
## quality quality2
## 2 6 6
## 7 6 6
## 10 6 6
## 11 5 5
## 12 5 5
## 13 5 5Fitting Model
polr() is the fuction used to fit the ordinal logistic regression model.
require(MASS)
require(reshape2)## Loading required package: reshape2# Hess=TRUE to let the model output show the observed information matrix from optimization which is used to get standard errors.
o_lrm <- polr(quality2 ~ . - quality, data = whitetrain, Hess=TRUE)# should not use vif to check for multicollinarilty in case of categorical veriable.
vif(o_lrm)## Warning in vif.default(o_lrm): No intercept: vifs may not be sensible.## fixed.acidity volatile.acidity citric.acid
## 9.437704e+08 3.828385e+08 6.424758e+05
## residual.sugar chlorides free.sulfur.dioxide
## 1.845434e+12 1.457669e+00 6.425308e+00
## total.sulfur.dioxide density pH
## 2.228460e+16 4.943074e+10 1.096931e+10
## sulphates alcohol
## 1.373230e+08 2.243821e+16summary(o_lrm)## Call:
## polr(formula = quality2 ~ . - quality, data = whitetrain, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## fixed.acidity 0.19109 0.049513 3.859
## volatile.acidity -4.16006 0.396986 -10.479
## citric.acid 0.47473 0.315331 1.506
## residual.sugar 0.18950 0.009254 20.477
## chlorides -1.02123 1.618406 -0.631
## free.sulfur.dioxide 0.01804 0.002996 6.022
## total.sulfur.dioxide -0.00179 0.001256 -1.425
## density -405.55303 0.611027 -663.724
## pH 2.37597 0.279420 8.503
## sulphates 1.83492 0.325575 5.636
## alcohol 0.51080 0.042104 12.132
##
## Intercepts:
## Value Std. Error t value
## 3|4 -393.5020 0.6212 -633.4946
## 4|5 -391.1957 0.6196 -631.3735
## 5|6 -388.2891 0.6262 -620.0709
## 6|7 -385.6073 0.6410 -601.5335
## 7|8 -383.3112 0.6550 -585.2349
## 8|9 -380.0386 0.8192 -463.9258
##
## Residual Deviance: 6145.137
## AIC: 6179.137Smaller the AIC better is the model. So let’s try step wise logistic regression.
Stepwise
o_lr = step(o_lrm)## Start: AIC=6179.14
## quality2 ~ (fixed.acidity + volatile.acidity + citric.acid +
## residual.sugar + chlorides + free.sulfur.dioxide + total.sulfur.dioxide +
## density + pH + sulphates + alcohol + quality) - quality
##
## Df AIC
## - chlorides 1 6177.5
## - total.sulfur.dioxide 1 6179.0
## <none> 6179.1
## - citric.acid 1 6179.4
## - fixed.acidity 1 6183.7
## - alcohol 1 6203.4
## - sulphates 1 6206.0
## - density 1 6207.4
## - free.sulfur.dioxide 1 6213.1
## - pH 1 6217.6
## - residual.sugar 1 6222.9
## - volatile.acidity 1 6287.0
##
## Step: AIC=6177.52
## quality2 ~ fixed.acidity + volatile.acidity + citric.acid + residual.sugar +
## free.sulfur.dioxide + total.sulfur.dioxide + density + pH +
## sulphates + alcohol
##
## Df AIC
## - total.sulfur.dioxide 1 6177.4
## <none> 6177.5
## - citric.acid 1 6177.6
## - fixed.acidity 1 6183.0
## - alcohol 1 6201.6
## - sulphates 1 6204.9
## - density 1 6208.0
## - free.sulfur.dioxide 1 6211.1
## - pH 1 6218.7
## - residual.sugar 1 6225.4
## - volatile.acidity 1 6288.9
##
## Step: AIC=6177.36
## quality2 ~ fixed.acidity + volatile.acidity + citric.acid + residual.sugar +
## free.sulfur.dioxide + density + pH + sulphates + alcohol
##
## Df AIC
## - citric.acid 1 6177.3
## <none> 6177.4
## - fixed.acidity 1 6184.0
## - alcohol 1 6200.2
## - sulphates 1 6204.0
## - density 1 6213.2
## - free.sulfur.dioxide 1 6216.0
## - pH 1 6219.7
## - residual.sugar 1 6229.7
## - volatile.acidity 1 6307.1
##
## Step: AIC=6177.31
## quality2 ~ fixed.acidity + volatile.acidity + residual.sugar +
## free.sulfur.dioxide + density + pH + sulphates + alcohol
##
## Df AIC
## <none> 6177.3
## - fixed.acidity 1 6184.9
## - alcohol 1 6201.8
## - density 1 6211.8
## - sulphates 1 6212.1
## - pH 1 6218.1
## - free.sulfur.dioxide 1 6225.9
## - residual.sugar 1 6228.4
## - volatile.acidity 1 6316.1We see some of the variables got eliminated.
fixed.acidity, alcohol, density, sulphates, pH, free.sulfur.dioxide, residual.sugar, volatile.acidity are the variables being considered.
head(fitted(o_lr))## 3 4 5 6 7 8
## 2 0.012733201 0.102476569 0.5891421 0.2675954 0.02514971 0.002792462
## 7 0.007370319 0.062364804 0.5085991 0.3739248 0.04270917 0.004839804
## 10 0.002796939 0.024740172 0.3137434 0.5415504 0.10395809 0.012702967
## 11 0.001083919 0.009752403 0.1561342 0.5775992 0.22198247 0.032135733
## 12 0.005402251 0.046583924 0.4488420 0.4350332 0.05727276 0.006603401
## 13 0.001270388 0.011408834 0.1775801 0.5833461 0.19772064 0.027554412
## 9
## 2 0.0001105250
## 7 0.0001919682
## 10 0.0005080325
## 11 0.0013121200
## 12 0.0002624039
## 13 0.0011195326Training Set Accuracy
predicting -
p <- predict(o_lr, type = "class")
head(p)## [1] 5 5 6 6 5 6
## Levels: 3 4 5 6 7 8 9Confusion Matrix Test
cm1 = as.matrix(table(Actual = whitetrain$quality2, Predicted = p))
cm1## Predicted
## Actual 3 4 5 6 7 8 9
## 3 0 1 5 7 0 1 0
## 4 0 4 64 38 1 0 0
## 5 0 0 425 393 4 0 0
## 6 1 0 239 936 76 0 0
## 7 0 0 19 355 107 1 0
## 8 0 0 2 59 31 0 0
## 9 0 0 0 1 3 0 0sum(diag(cm1))/length(whitetrain$quality2)## [1] 0.530833Training Set Accuracy is 53.08%
Test Set Accuracy
Accuray for the test set.
tst_pred <- predict(o_lr, newdata = whitetest, type = "class")Confusion Matrix Test
cm2 <- table(predicted = tst_pred, actual = whitetest$quality2)
cm2## actual
## predicted 3 4 5 6 7 8 9
## 3 0 0 0 0 0 0 0
## 4 0 0 0 0 0 0 0
## 5 2 32 177 93 6 1 0
## 6 3 14 171 410 156 26 0
## 7 1 0 5 33 45 12 1
## 8 0 0 0 0 0 0 0
## 9 0 0 0 0 0 0 0Test Set Accuracy
sum(diag(cm2))/length(whitetrain$quality2)## [1] 0.227912Test Set Accuracy - 22.8 %
Binomial Logistic Regression Model
The variable to be predicted is binary and hence we use binomial logistic regression.
white$category[white$quality <= 5] <- 0
white$category[white$quality > 5] <- 1
white$category <- as.factor(white$category)head(white)## fixed.acidity volatile.acidity citric.acid residual.sugar chlorides
## 1 7.0 0.27 0.36 20.7 0.045
## 2 6.3 0.30 0.34 1.6 0.049
## 3 8.1 0.28 0.40 6.9 0.050
## 4 7.2 0.23 0.32 8.5 0.058
## 7 6.2 0.32 0.16 7.0 0.045
## 10 8.1 0.22 0.43 1.5 0.044
## free.sulfur.dioxide total.sulfur.dioxide density pH sulphates alcohol
## 1 45 170 1.0010 3.00 0.45 8.8
## 2 14 132 0.9940 3.30 0.49 9.5
## 3 30 97 0.9951 3.26 0.44 10.1
## 4 47 186 0.9956 3.19 0.40 9.9
## 7 30 136 0.9949 3.18 0.47 9.6
## 10 28 129 0.9938 3.22 0.45 11.0
## quality quality2 category
## 1 6 6 1
## 2 6 6 1
## 3 6 6 1
## 4 6 6 1
## 7 6 6 1
## 10 6 6 1Train Test Split
set.seed(3000)
spl = sample.split(white$category, SplitRatio = 0.7)
whitetrain = subset(white, spl==TRUE)
whitetest = subset(white, spl==FALSE)
head(whitetrain)## fixed.acidity volatile.acidity citric.acid residual.sugar chlorides
## 2 6.3 0.30 0.34 1.60 0.049
## 7 6.2 0.32 0.16 7.00 0.045
## 10 8.1 0.22 0.43 1.50 0.044
## 11 8.1 0.27 0.41 1.45 0.033
## 13 7.9 0.18 0.37 1.20 0.040
## 15 8.3 0.42 0.62 19.25 0.040
## free.sulfur.dioxide total.sulfur.dioxide density pH sulphates alcohol
## 2 14 132 0.9940 3.30 0.49 9.5
## 7 30 136 0.9949 3.18 0.47 9.6
## 10 28 129 0.9938 3.22 0.45 11.0
## 11 11 63 0.9908 2.99 0.56 12.0
## 13 16 75 0.9920 3.18 0.63 10.8
## 15 41 172 1.0002 2.98 0.67 9.7
## quality quality2 category
## 2 6 6 1
## 7 6 6 1
## 10 6 6 1
## 11 5 5 0
## 13 5 5 0
## 15 5 5 0glm()= Generalized linear model
We will use the glm() - Generalized linear model command to run a logistic regression, regressing success on the numeracy and anxiety scores.
model_glm <- glm(category ~ . - quality - quality2, data = whitetrain, family=binomial(link = "logit"))
#backwards = step(fullmod,trace=0) would suppress step by step output.Stepwise
model_gl <- step(model_glm)## Start: AIC=2821.85
## category ~ (fixed.acidity + volatile.acidity + citric.acid +
## residual.sugar + chlorides + free.sulfur.dioxide + total.sulfur.dioxide +
## density + pH + sulphates + alcohol + quality + quality2) -
## quality - quality2
##
## Df Deviance AIC
## - total.sulfur.dioxide 1 2797.9 2819.9
## - chlorides 1 2797.9 2819.9
## - fixed.acidity 1 2798.0 2820.0
## <none> 2797.8 2821.8
## - citric.acid 1 2800.2 2822.2
## - free.sulfur.dioxide 1 2806.9 2828.9
## - pH 1 2809.7 2831.7
## - density 1 2811.5 2833.5
## - sulphates 1 2817.4 2839.4
## - residual.sugar 1 2823.2 2845.2
## - alcohol 1 2826.3 2848.3
## - volatile.acidity 1 2916.6 2938.6
##
## Step: AIC=2819.87
## category ~ fixed.acidity + volatile.acidity + citric.acid + residual.sugar +
## chlorides + free.sulfur.dioxide + density + pH + sulphates +
## alcohol
##
## Df Deviance AIC
## - chlorides 1 2798.0 2818.0
## - fixed.acidity 1 2798.0 2818.0
## <none> 2797.9 2819.9
## - citric.acid 1 2800.2 2820.2
## - pH 1 2809.7 2829.7
## - density 1 2811.8 2831.8
## - free.sulfur.dioxide 1 2813.1 2833.1
## - sulphates 1 2817.8 2837.8
## - residual.sugar 1 2823.4 2843.4
## - alcohol 1 2826.6 2846.6
## - volatile.acidity 1 2920.6 2940.6
##
## Step: AIC=2817.97
## category ~ fixed.acidity + volatile.acidity + citric.acid + residual.sugar +
## free.sulfur.dioxide + density + pH + sulphates + alcohol
##
## Df Deviance AIC
## - fixed.acidity 1 2798.2 2816.2
## <none> 2798.0 2818.0
## - citric.acid 1 2800.3 2818.3
## - pH 1 2810.3 2828.3
## - density 1 2812.6 2830.6
## - free.sulfur.dioxide 1 2813.1 2831.1
## - sulphates 1 2818.1 2836.1
## - residual.sugar 1 2825.1 2843.1
## - alcohol 1 2826.7 2844.7
## - volatile.acidity 1 2922.7 2940.7
##
## Step: AIC=2816.17
## category ~ volatile.acidity + citric.acid + residual.sugar +
## free.sulfur.dioxide + density + pH + sulphates + alcohol
##
## Df Deviance AIC
## <none> 2798.2 2816.2
## - citric.acid 1 2800.6 2816.6
## - free.sulfur.dioxide 1 2813.1 2829.1
## - pH 1 2817.6 2833.6
## - sulphates 1 2818.2 2834.2
## - density 1 2824.1 2840.1
## - residual.sugar 1 2841.8 2857.8
## - alcohol 1 2863.4 2879.4
## - volatile.acidity 1 2927.8 2943.8#model_gl = step(model_glm,trace=0) would suppress step by step output.Prediction - Train
The default is on the scale of the linear predictors; the alternative “response” is on the scale of the response variable. Thus for a default binomial model the default predictions are of log-odds (probabilities on logit scale) and type = “response” gives the predicted probabilities
head(fitted(model_gl))## 2 7 10 11 13 15
## 0.3085439 0.3942105 0.6825519 0.8508986 0.8124451 0.4915685head(predict(model_gl))## 2 7 10 11 13 15
## -0.80693550 -0.42964756 0.76552441 1.74166660 1.46597646 -0.03372931head(predict(model_gl, type = "response"))## 2 7 10 11 13 15
## 0.3085439 0.3942105 0.6825519 0.8508986 0.8124451 0.4915685Catagorizing Set
trn_pred <- ifelse(predict(model_gl, type = "response") > 0.5,"Good Wine", "Bad Wine")
head(trn_pred)## 2 7 10 11 13 15
## "Bad Wine" "Bad Wine" "Good Wine" "Good Wine" "Good Wine" "Bad Wine"Confusion Matrix - Training Set
Obtaining confusion matrix of the training data.
trn_tab <- table(predicted = trn_pred, actual = whitetrain$category)
trn_tab## actual
## predicted 0 1
## Bad Wine 504 254
## Good Wine 440 1575Training Set Accuracy
Checking accuracy of the training set.
sum(diag(trn_tab))/length(whitetrain$category)## [1] 0.7497295We can see that Binomial Logistic Regressio gives an training set accuracy of 74.97295%.
Confusion Matrix - Test Set
Confusion matrix for the test data.
# Making predictions on the test set.
tst_pred <- ifelse(predict(model_gl, newdata = whitetest, type = "response") > 0.5, "Good Wine", "Bad Wine")
tst_tab <- table(predicted = tst_pred, actual = whitetest$category)
tst_tab## actual
## predicted 0 1
## Bad Wine 212 110
## Good Wine 192 674Test Set Accuracy
Checking accuracy for the test data.
sum(diag(tst_tab))/length(whitetest$category)## [1] 0.7457912We can see that Binomial Logistic Regression gives an test set accuracy of 74.58%.
Conclusion
Please, use binomial logistic regression to predict the quality of wine.