LAB REPORT 1

LOW POH SWEAN

SD23018- 02G

1. Extract data from the following excel file: gmp.txt.

file.choose()

library(readxl)
data<- read_excel("C:\\Users\\joone\\OneDrive\\Desktop\\Y3 SMS\\gmp.xlsx")
data

library(Hmisc)

# Check for missing values
summary(data$x3)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
   81.0   171.2   243.0   217.9   258.8   366.0       2 

# Impute and convert to numeric
data$x3 <- as.numeric(impute(data$x3, mean))

# Confirm no missing values remain
summary(data$x3)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   81.0   173.8   231.5   217.9   256.2   366.0 

data

2. Develop a multiple linear regression model by considering all variables.

Part A

1. Construct a normal probability plot of the residuals. Does there seem to be any problems with the normality assumption?

# Fit the multiple linear regression model
model2 <- lm(y~., data =data)

# Model summary
summary_model2 <- summary(model2)
print(summary_model2)

Call:
lm(formula = y ~ ., data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9134 -1.8769 -0.2887  1.7538  4.7129 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept)  9.755693  28.797397   0.339   0.7383  
x1          -0.026803   0.034760  -0.771   0.4497  
x2          -0.003633   0.056401  -0.064   0.9493  
x3           0.026598   0.029589   0.899   0.3794  
x4           1.620045   2.356393   0.688   0.4997  
x5           5.003852   2.998637   1.669   0.1108  
x6          -0.136788   1.182675  -0.116   0.9091  
x7          -2.612686   2.881391  -0.907   0.3753  
x8           0.210631   0.109101   1.931   0.0678 .
x9          -0.318245   0.303778  -1.048   0.3073  
x10         -0.006969   0.004370  -1.595   0.1265  
x11          0.617416   2.940955   0.210   0.8358  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.169 on 20 degrees of freedom
Multiple R-squared:  0.8377,    Adjusted R-squared:  0.7484 
F-statistic: 9.382 on 11 and 20 DF,  p-value: 1.094e-05

residuals <- residuals(model2)
qqnorm(residuals)
qqline(residuals)

Interpret:

The points mostly follow the straight line in the middle section.

However, there is some deviation at both tails such as the extreme left and right.

There is no serious violation of the normality assumption.

The residuals are approximately normally distributed, though there is slight deviation in the tails, which suggests the presence of minor outliers or mild non-normality.

Conclusion, the Q-Q plot shows that the residuals generally lie along the straight line, indicating that the normality assumption is reasonably satisfied. Minor deviations at the tails suggest a few outliers or slight non-normality, but not enough to be a major concern.

plot(model2,which=2)  # Q-Q plot

residuals2 <- resid(model2)
hist(residuals2, main = "Histogram of Residuals", col = "skyblue")

Interpret:

The histogram is roughly bell-shaped, centered around zero.

There is some slight skewness, the right tail seems a bit longer than the left.

However, the overall pattern looks reasonably symmetric, which supports the normality assumption.

Conclusion, the histogram of residuals appears approximately symmetric and bell-shaped, with most residuals clustering around zero. This suggests that the residuals are roughly normally distributed. Although there is slight skewness in the right tail, the deviation from normality is not severe. Therefore, the normality assumption for the multiple linear regression model is reasonably satisfied.

2. Support your answer in (i) by using appropriate normality test.

#Two-Sample K–S Test
ks.test(residuals2, "pnorm", mean(residuals2), sd(residuals2))

    Exact one-sample Kolmogorov-Smirnov test

data:  residuals2
D = 0.095955, p-value = 0.9027
alternative hypothesis: two-sided

Interpret:

\(H_{0}\): Residuals are normally distributed.

\(H_{1}\): Residuals are not normally distributed.

\(p-value=0.9027\)

Since \((p-value=0.9027)\) >\((\alpha=0.05)\), do not reject \(H_{0}\).

At \(\alpha=0.05\), normality assumption holds.

library(nortest)
#Anderson-Darling Test to test for normality
ad.test(residuals2)

    Anderson-Darling normality test

data:  residuals2
A = 0.37293, p-value = 0.3982

Interpret:

\(H_{0}\): Residuals are normally distributed.

\(H_{1}\): Residuals are not normally distributed.

\(p-value=0.3982\)

Since \((p-value=0.3982)\) >\((\alpha=0.05)\), do not reject \(H_{0}\).

At \(\alpha=0.05\), normality assumption holds.

#conduct shapirowilk Test to test for normality
shapiro.test(residuals2)

    Shapiro-Wilk normality test

data:  residuals2
W = 0.96269, p-value = 0.3249

Interpretation:

\(H_{0}\): Residuals are normally distributed.

\(H_{1}\): Residuals are not normally distributed.

\(p-value=0.3249\)

Since \((p-value=0.3249)\) >\((\alpha=0.05)\), do not reject \(H_{0}\).

At \(\alpha=0.05\), normality assumption holds.

Based on the Shapiro–Wilk, Anderson–Darling, and Kolmogorov–Smirnov tests, the p-values are greater than 0.05, indicating that the residuals are not significantly different from a normal distribution. Therefore, the normality assumption of the multiple linear regression model is supported.

3. Construct and interpret a plot of residuals versus the predicted response.

residuals <- resid(model2)
predicted_response <- fitted(model2)
plot(predicted_response, residuals,
         xlab = "Predicted Response",
         ylab = "Residuals",
         main = "Residuals vs. Predicted Response Plot")
abline(h = 0, col = "red", lty = 2) # Add a horizontal line at 0

Interpret:

The residuals are randomly scattered around the horizontal line at zero.

There is no clear pattern or systematic curvature, meaning the linearity assumption is satisfied.

The spread which is variance of residuals appears roughly constant across predicted values, so it is no obvious funnel or cone shape.

A few points deviate slightly from zero, but there are no extreme outliers.

Conclusion, the residuals vs. predicted response plot shows that the residuals are randomly distributed around zero with no clear trend or pattern. This indicates that the assumptions of linearity and constant variance (homoscedasticity) are reasonably met. Therefore, the model appears appropriate, and there is no evidence of heteroscedasticity or non-linearity.

Part B

Then, detect any outliers occurs using Cook’s Distance method.

1. Plot the influential observation by Cook’s Distance and comments.

#Cook distance
model2 <- lm(y ~., data =data)
cooksd<-cooks.distance(model2)
plot(cooksd, pch="*", cex=2, main="Influential Obsby Cooks distance") # plot cook's distance
abline(h = 4*mean(cooksd, na.rm=T), col="red") # add cutoffline

text(x=1:length(cooksd)+1, y=cooksd, labels=ifelse(cooksd>4*mean(cooksd, na.rm=T),names(cooksd),""), col="red") # add labels

Interpret:

The red horizontal line represents the cutoff value (4 * mean(Cook’s distance)).

Most points are below this line, meaning they are not influential.

Observations 14 and 17 are above the red line, indicating that they are potentially influential observations.

This means these two data points have a disproportionate effect on the regression coefficients and may influence the model’s fit.

Conclusion, the Cook’s Distance plot shows that most observations have small Cook’s distance values, indicating they do not have an undue influence on the regression model. However, observations 14 and 17 exceed the cutoff line (4 × mean Cook’s Distance), suggesting that they are potentially influential. These points should be further investigated to determine whether they are valid data points or outliers that may distort the regression results.

2. Examine and list the point observation(s) which consider as outlier(s). Justify your answer.

# influential row numbers
influential <-as.numeric(names(cooksd)[(cooksd> 4*mean(cooksd, na.rm=T))])
influential 

[1] 14 17

head(data[influential,]) # influential observations.

NA

Interpret:

In these 2 row have very high weight(x10), and the type of transmission both is automatic.

Based on the Cook’s Distance analysis, observations 14 and 17 exceed the threshold value (4 × mean Cook’s Distance). These observations are considered influential outliers, meaning they have a significant effect on the regression model’s coefficients. Their presence could distort the model fit and influence parameter estimates. It is advisable to examine these data points further to determine whether they are valid measurements or data entry errors. If they are genuine, they should be retained; if not, they may be removed or analyzed separately.

Part C

1. By considering only x1, x2, x3, x8, x9 and x10, construct the lack of fit test. Interpret and justify your answer.

#Develop model and consider the significant variables
model1 <-lm(y~x1+x2+x3+x8+x9+x10, data = data)
summary(model1)

Call:
lm(formula = y ~ x1 + x2 + x3 + x8 + x9 + x10, data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.8683 -1.8362 -0.2086  1.7344  6.2344 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept) 29.551590  16.249704   1.819    0.081 .
x1          -0.035411   0.025308  -1.399    0.174  
x2           0.019203   0.040693   0.472    0.641  
x3           0.001668   0.024560   0.068    0.946  
x8           0.140948   0.098906   1.425    0.167  
x9          -0.187418   0.251940  -0.744    0.464  
x10         -0.004439   0.003718  -1.194    0.244  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.089 on 25 degrees of freedom
Multiple R-squared:  0.8072,    Adjusted R-squared:  0.7609 
F-statistic: 17.44 on 6 and 25 DF,  p-value: 7.662e-08

model2 <-lm(y~., data = data)
summary(model2)

Call:
lm(formula = y ~ ., data = data)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9134 -1.8769 -0.2887  1.7538  4.7129 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept)  9.755693  28.797397   0.339   0.7383  
x1          -0.026803   0.034760  -0.771   0.4497  
x2          -0.003633   0.056401  -0.064   0.9493  
x3           0.026598   0.029589   0.899   0.3794  
x4           1.620045   2.356393   0.688   0.4997  
x5           5.003852   2.998637   1.669   0.1108  
x6          -0.136788   1.182675  -0.116   0.9091  
x7          -2.612686   2.881391  -0.907   0.3753  
x8           0.210631   0.109101   1.931   0.0678 .
x9          -0.318245   0.303778  -1.048   0.3073  
x10         -0.006969   0.004370  -1.595   0.1265  
x11          0.617416   2.940955   0.210   0.8358  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.169 on 20 degrees of freedom
Multiple R-squared:  0.8377,    Adjusted R-squared:  0.7484 
F-statistic: 9.382 on 11 and 20 DF,  p-value: 1.094e-05

#lack of fit test
anova(model1, model2)

Analysis of Variance Table

Model 1: y ~ x1 + x2 + x3 + x8 + x9 + x10
Model 2: y ~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 + x10 + x11
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1     25 238.62                           
2     20 200.90  5    37.723 0.7511  0.595

Interpretation:

From the table, we can see the table show the p-value is 0.5950, which is greater than 0.05.

\(H_{0}\):The reduced model (Model 1) fits the data adequately, there is no significant lack of fit.

\(H_{1}\):The reduced model does not fit adequately, the full model fits significantly better.

\(p-value=0.5950\)

Since \((p-value=0.5950)\) >\((\alpha=0.05)\), do not reject \(H_{0}\).

At \(\alpha=0.05\), the reduced model (Model 1) fits the data adequately; there is no significant lack of fit.

The simpler model (with x1, x2, x3, x8, x9, and x10) is adequate.

The additional predictors (x4, x5, x6, x7, x11) in the full model do not significantly improve the fit.

The lack of fit test was conducted by comparing the reduced model using variables x1, x2, x3, x8, x9, and x10 with the full model containing all predictors. The ANOVA result gives a p-value of 0.5950, which is greater than 0.05. Therefore, we fail to reject the null hypothesis, indicating that the reduced model provides an adequate fit to the data. The additional variables in the full model do not significantly improve the model performance. Hence, the simpler model is sufficient to explain the relationship between the response variable y and the selected predictors.