# Multivariate regression in R

I found a very nice article on the internet that explained how to develop a multivariate regression model. In the article a dataset for several hundred 2005 used General Motors (GM) cars is supplied which is the basis for the step by step explanation. I think this is a good way to learn about regression models. I thought it was a good idea to try to reproduce the result, intercepts, slopes and plots in the article using R.

This (my) article will just show how to reproduce the results in the for-mentioned article using R. I have also added two other ways to get the best fit using the power of R.

The article uses 4 equations to develop the model into the best fit.

• Equation 1: Based on Mileage alone
• Equation 2: Based on Mileage, Cylinder, Doors, Cruise, Sound, Leather
• Equation 3: Based on Mileage, Buick, Cadillac, Chev, Pontiac, SAAB predicting the log of price
• Equation 4: Based on Mileage, Liter, Doors, Cruise, Sound, Leather, Buick, Cadillac, Chev, Pontiac, SAAB, Conv, Hatchback, Sedan predicting the log of price

Before starting modelling the regressions with R first we need to read the data and setup the required libraries.

```library("tidyr")
library("dplyr")

```

Equation 1

The first equation is simple and just to explore:

```model <- lm(Price ~ Mileage, data=kuiper)
```

Result:

``````## (Intercept) 24764.5590
## Mileage        -0.1725
## ---
## Multiple R-squared:  0.0205, Adjusted R-squared:  0.0192
``````

To create the plot of the model use this code:

```plot(model\$residuals, col="blue", type="l")
points(model\$residuals, col="red", type="p", pch=20)
abline(h=0)
``` Equation 2
Still a quite simple model, just with more predictors.

```model_acc <- lm(Price ~ Mileage + Cylinder +
Doors + Cruise + Sound + Leather, data=kuiper)
```
``````## (Intercept)  7323.1643
## Mileage        -0.1705
## Cylinder     3200.1246
## Doors       -1463.3991
## Cruise       6205.5113
## Sound       -2024.4007
## Leather      3327.1433
## ---
## Multiple R-squared:  0.446,  Adjusted R-squared:  0.442
``````

If we plot the model we get 4 plots which explain different aspects of the model.

```par(mfrow=c(2,2))
plot(model_acc)
par(mfrow=c(1,1))
``` Fot the histogram and distribution of the residuals we need different plots.

```hist(model_acc\$residuals, breaks=30)

plot(model_acc\$residuals, col="blue", type="l")
points(model_acc\$residuals, col="red", type="p", pch=20)
abline(h=0)
```  Equation 3

For this equation we first need to alter the dataset to add the factor variable Make as boolean variables and take the log of price.

```kuiper\$id <- 1:nrow(kuiper)
make <- as.data.frame(table(kuiper\$id, kuiper\$Make))
make <- make %>% spread(Var2, Freq) %>% rename(id = Var1)
kuiper <- merge(kuiper, make, by="id")
kuiper <- select(kuiper, -Make)

kuiper\$tPrice <- log10(kuiper\$Price)
model_make <- lm(tPrice ~ Mileage + Liter +
Buick + Cadillac + Chevrolet +
Pontiac + SAAB, data=kuiper)
```
``````## (Intercept)  3.979911512
## Mileage     -0.000003476
## Liter        0.099724553
## Buick        0.039968616
## Chevrolet   -0.009372280
## Pontiac      0.013613268
## SAAB         0.345305253
## ---
## Multiple R-squared:  0.917,  Adjusted R-squared:  0.916 ``````

Same plot as before.

```par(mfrow=c(2,2))
plot(model_make)
par(mfrow=c(1,1))
``` Histogram and residuals.

```hist(model_make\$residuals, breaks=30)
plot(model_make\$residuals, col="blue", type="l")
points(model_make\$residuals, col="red", type="p", pch=20)
abline(h=0)
```  Equation 4
For equation 4 the dataset needs to be altered for the Type variable.

```type <- as.data.frame(table(kuiper\$id, kuiper\$Type))
type <- type %>% spread(Var2, Freq) %>% rename(id = Var1)
kuiper <- merge(kuiper, type, by="id")
kuiper <- select(kuiper, -Type)
model_all <- lm(tPrice ~ Mileage + Liter + Doors +
Cruise + Sound + Leather + Buick +
Cadillac + Chevrolet + Pontiac + SAAB +
Convertible + Hatchback + Sedan, data=kuiper)
```
``````## (Intercept)  3.91810880
## Mileage     -0.00000358
## Liter        0.09576214
## Doors        0.03352658
## Cruise       0.00751699
## Sound        0.00522264
## Leather      0.00625998
## Buick        0.04165259
## Chevrolet   -0.01331520
## Pontiac     -0.00042065
## SAAB         0.28109769
## Convertible  0.13781907
## Hatchback   -0.08898857
## Sedan       -0.07114898
## ---
## Multiple R-squared:  0.952,  Adjusted R-squared:  0.951 ``````

Again the same plots.

```par(mfrow=c(2,2))
plot(model_all)
par(mfrow=c(1,1))

hist(model_all\$residuals, breaks=30)
plot(model_all\$residuals, col="blue", type="l")
points(model_all\$residuals, col="red", type="p", pch=20)
abline(h=0)
```   Let R find the best model

R has the Step function. You can supply two models and R steps through the model to get to the best fit. I’ve removed the variables I don’t want to use for predicting.

```kuiper <- kuiper %>% select(-Price, -id, -Model, -Trim)

model <- lm(tPrice ~ Mileage, data=kuiper)
model_all <- lm(tPrice ~ ., data=kuiper)
bestfit <- step(model, scope=list(lower=model, upper=model_all), direction="forward")
```

Result:

``````## (Intercept)  4.006511863
## Mileage     -0.000003579
## Liter        0.109135503
## SAAB         0.276731601
## Convertible  0.142481573
## Wagon        0.069803273
## Buick        0.041279545
## Hatchback   -0.013685317
## Chevrolet   -0.014330497
## Cylinder    -0.012359839
## Cruise       0.008355646
## Leather      0.005336807
## Sound        0.004812547
## ---
## Multiple R-squared:  0.952,  Adjusted R-squared:  0.952 ``````

This model is slightly different and produces almost the same R-squared.

Machine learning

Using the Caret package we can use machine learning to predict the best model on the dataset.

```library("caret")

ml_model <- train(tPrice ~ ., data=kuiper, method="lm")
```
``````## (Intercept)  3.93932861
## Mileage     -0.00000358
## Cylinder    -0.01224029
## Liter        0.10880574
## Doors        0.03377855
## Cruise       0.00808055
## Sound        0.00462693
## Leather      0.00517194
## Buick        0.04375913
## Chevrolet   -0.01264695
## Pontiac      0.00261046
## SAAB         0.27901919
## Saturn               NA
## Convertible  0.14078278
## Coupe                NA
## Hatchback   -0.08269889
## Sedan       -0.06948186
## Wagon                NA
## ---
## Multiple R-squared:  0.952,  Adjusted R-squared:  0.951 ``````

Remarks

I would like to add that normally you would not use the complete dataset to fit a model for prediction. You would want to split the dataset into a training and testset to be able to cross validate your model to prevent over fitting of the model.