How to use a train set and a test set to check if my model is over-parameterized in R?

Question

set.seed(1)
n<-11
x1<-1:n
x2<-rnorm(n)
x3<-1/20*rnorm(n,x1,x1/7)+x2/5
cor(x2,x3) ;cor(x1,x3)
y<-1+5*x1+6*x2+2*x3+rnorm(n)
x4<-runif(n) 
x5<-runif(n) ; x6<-runif(n) ; x7<-runif(n)
x8<-runif(n) ; x9<-runif(n) ; x10<-runif(n)
train<-data.frame(y=y,x1=x1,x2=x2,x3=x3,x4=x4,x5=x5,
                 x6=x6,x7=x7,x8=x8,x9=x9,x10=x10)

n<-11
x1<-1:n
x2<-rnorm(n)
x3<-1/20*rnorm(n,x1,x1/7)+x2/5
cor(x2,x3) ;cor(x1,x3)
y<-1+5*x1+6*x2+2*x3+rnorm(n)
x4<-runif(n) 
x5<-runif(n) ; x6<-runif(n) ; x7<-runif(n)
x8<-runif(n) ; x9<-runif(n) ; x10<-runif(n)
test<-data.frame(y=y,x1=x1,x2=x2,x3=x3,x4=x4,x5=x5,
                 x6=x6,x7=x7,x8=x8,x9=x9,x10=x10)

I just created a train set and a test set in order to understand how the comparison of them can help me to find which model to fit.

You can see from my construction that only $x_1$, $x_2$ and $x_3$ are important to $y$. As a matter of fact, I claim that at most two of them are important to $y$, since $x_3$ is highly correlated with $x_1$ and $x_2$.

So I claim for example that the model

m1<-lm(y~x1+x2,data=train)

is much better than the models

m2<-lm(y~x1+x2+x3,data=train)
m3<-lm(y~x1+x2+x3+x4+x5+x6+x7+x8+x9,data=train)
m4<-lm(y~.,data=train)

and that $m_3$ is over-parameterized.

How can I use the test set and the train set to confirm these two claims? I have read similar questions online and some answers which used the MSE (mean squared error) and the predictive power of a model, but I understood none of them and I don't know how to apply them to solve this problem. How can I use them in R to answer my question?

And lastly, what is actually the out-of-sample predictive power of a model? Can it be measured with some way? What is the predictive power of the models $m_1$ and $m_3$? Which of these two models has a greater predictive power?

Answer 1

Firstly, I would consider a fifth "baseline" model $M_0$ which basically is a model that always predicts using the average (In your example 33.165) and now we have the following 5:

m0<-lm(y~1,data=train)
m1<-lm(y~x1+x2,data=train)
m2<-lm(y~x1+x2+x3,data=train)
m3<-lm(y~x1+x2+x3+x4+x5+x6+x7+x8+x9,data=train)
m4<-lm(y~.,data=train)

Firstly, let's take a look at a the train data.

# Baseline model vs Model with X1 and X2 only
anova(m0, lm(y~x1,data=train) , test = "F")
anova(m0, lm(y~x2,data=train) , test = "F")

# Model with X1 vs m1 and model with X2 vs m1
anova(lm(y~x1,data=train), m1 , test = "F")
anova(lm(y~x2,data=train), m1 , test = "F")

# Model m1 vs m2
anova(m1, m2, test = "F")

From the training data, we can see that the simple regression models with X1 and X2 only yield statistically significant better fit vs the model with just the average (looking at the p-values). Similarly for the combination of X1 and X2 vs the models with X1 and X2 independently. However, we have no evidence than a model with more variables works better than m1. This is all from the training dataset.

Now, let's predict on the new test dataset and see how accurate each model is. We will then calculate the MSE to see how close are predicted values are to the real values (You can think of the MSE as the average difference between real values and the predicted ones squared. It basically measures how close the predicted to real are and then it squares it so that it penalises them the further they are)

pred_m0 <- predict(m0, newdata = test)
pred_m1 <- predict(m1, newdata = test)
pred_m2 <- predict(m2, newdata = test)
pred_m3 <- predict(m3, newdata = test)
pred_m4 <- predict(m4, newdata = test)

m0_MSE  <- mean((pred_m0 - test$y)^2)
m1_MSE  <- mean((pred_m1 - test$y)^2)
m2_MSE  <- mean((pred_m2 - test$y)^2)
m3_MSE  <- mean((pred_m3 - test$y)^2)
m4_MSE  <- mean((pred_m4 - test$y)^2)

It seems that $M_1$ has the lowest out-of-sample MSE and thus the best predictive power which in that case agrees with what we saw before in the training dataset.

To measure some kind of accuracy for models 1 and 3 I will compare them with the baseline model (predicting with just the average)

sum((pred_m1 - test$y)^2)/
sum((pred_m0 - test$y)^2)

sum((pred_m3 - test$y)^2)/
sum((pred_m0 - test$y)^2)

And to me this shows that $M1$ predicts with a less than 1% error compared to the baseline model whereas $M3$ predicts with 20% of the error of the baseline model. In other words, one predicts 260 times better than the average and the other one 5 times better than the average.