HW 3
YAREN YILMAZ
1. Import the dataset Apartments.xlsx
library(readxl)
mydata <- read_xlsx("./Apartments.xlsx")
mydata <- as.data.frame(mydata)Description:
- Age: Age of an apartment in years
- Distance: The distance from city center in km
- Price: Price per m2
- Parking: 0-No, 1-Yes
- Balcony: 0-No, 1-Yes
2. Change categorical variables into factors. (0.5 p)
mydata$ParkingF <- factor(mydata$Parking,
levels = c("0","1"),
labels = c("without_parking","with_parking"))
mydata$BalconyF <- factor(mydata$Balcony,
levels = c("0","1"),
labels = c("without_balcony","with_balcony"))3. Test the hypothesis H0: Mu_Price = 1900 eur. What can you conclude? (1 p)
- Evaluation of the assumptions of t-test:
- Price must be numeric: met.
- Price must be normally distributed: violated as histogram shows and Shapiro-Wilk normality test suggests below (p<0.001).
library(ggplot2)
ggplot(mydata, aes(x = Price)) +
geom_histogram(binwidth = 350, colour = "black") +
ylab("Frequency") +
xlab("Price")
shapiro.test(mydata$Price)##
## Shapiro-Wilk normality test
##
## data: mydata$Price
## W = 0.94017, p-value = 0.0006513- Therefore, using Wilcoxon Signed Rank Test,
median(mydata$Price)## [1] 1950wilcox.test(mydata$Price,
mu = 1900,
correct = FALSE)##
## Wilcoxon signed rank test
##
## data: mydata$Price
## V = 2328, p-value = 0.02828
## alternative hypothesis: true location is not equal to 1900library(effectsize)
effectsize(wilcox.test(mydata$Price,
mu = 1900,
correct = FALSE))## r (rank biserial) | 95% CI
## --------------------------------
## 0.27 | [0.04, 0.48]
##
## - Deviation from a difference of 1900.interpret_rank_biserial(0.27, rules = "funder2019")## [1] "medium"
## (Rules: funder2019)- Conclusion: Based on the sample data, we found that the median price per m2 for apartments in the Vienna region is not equal to 1900 (p=0.029, r=0.27 - medium effect), in fact, it was higher (1950).
4. Estimate the simple regression function: Price = f(Age). Save results in object fit1 and explain the estimate of regression coefficient, coefficient of correlation and coefficient of determination. (2 p)
fit1 <- lm(Price ~ Age,
data = mydata)
summary(fit1)##
## Call:
## lm(formula = Price ~ Age, data = mydata)
##
## Residuals:
## Min 1Q Median 3Q Max
## -623.9 -278.0 -69.8 243.5 776.1
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2185.455 87.043 25.108 <2e-16 ***
## Age -8.975 4.164 -2.156 0.034 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 369.9 on 83 degrees of freedom
## Multiple R-squared: 0.05302, Adjusted R-squared: 0.04161
## F-statistic: 4.647 on 1 and 83 DF, p-value: 0.03401cor(mydata$Price, mydata$Age,
method = "pearson")## [1] -0.230255Estimate of regression coefficient = -8.975: If age of the apartment increases by 1 year, the price per m2 of the apartment decreases on average by 8.975 units (euros, dollars, etc.) (p=0.034).
Coefficient of correlation = -0.230255: Linear relationship between price and age is negative and weak.
Coefficient of determination = 0.05302: 5.302% of variability of price can be explained by linear effect of age. All other factors not included in the model explain the remaining 94.698%.
5. Show the scatterplot matrix between Price, Age and Distance. Based on the matrix determine if there is potential problem with multicolinearity. (0.5 p)
library(car)## Zorunlu paket yükleniyor: carDatascatterplotMatrix(mydata[ , c(-4:-7)],
smooth = FALSE)
Based on the matrix above, age and distance have almost no correlation or correlation that is positive but very close to zero (very weak). Therefore, it is unlikely that there is a problem with multicolinearity.
6. Estimate the multiple regression function: Price = f(Age, Distance). Save it in object named fit2.
fit2 <- lm(Price ~ Age + Distance,
data = mydata)7. Check the multicolinearity with VIF statistics. Explain the findings. (1 p)
vif(fit2)## Age Distance
## 1.001845 1.001845mean(vif(fit2))## [1] 1.001845Since both VIF values are much less than 5 and also average VIF is very close to 1, there is no multicolinearity.
8. Calculate standardized residuals and Cooks Distances for model fit2. Remove any potentially problematic units (outliers or units with high influence). (2 p)
mydata$StdResid <- round(rstandard(fit2), 3) #Calculation of standardized residuals
mydata$CooksD <- round(cooks.distance(fit2), 3) #Calculation of Cooks distances
hist(mydata$StdResid, #to check outliers
xlab = "Standardized residuals",
ylab = "Frequency",
main = "Histogram of standardized residuals")
hist(mydata$CooksD, #to check units with high influence
xlab = "Cooks distance",
ylab = "Frequency",
main = "Histogram of Cooks distances")
Histogram of standardized residuals shows no unit falling outside of the interval [-3,+3]. Therefore, there is no outlier.
Histogram of Cooks distances shows a gap between units. The units on the very right of the histogram are units with high impact and should be removed.
head(mydata[order(-mydata$CooksD),], 6) #Six units with highest value of Cooks distance (to find the unit with high impact)## Age Distance Price Parking Balcony ParkingF BalconyF StdResid CooksD
## 38 5 45 2180 1 1 with_parking with_balcony 2.577 0.320
## 55 43 37 1740 0 0 without_parking without_balcony 1.445 0.104
## 33 2 11 2790 1 0 with_parking without_balcony 2.051 0.069
## 53 7 2 1760 0 1 without_parking with_balcony -2.152 0.066
## 22 37 3 2540 1 1 with_parking with_balcony 1.576 0.061
## 39 40 2 2400 0 1 without_parking with_balcony 1.091 0.038library(dplyr)##
## Attaching package: 'dplyr'## The following object is masked from 'package:car':
##
## recode## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionmydata_new <- mydata[c(-38,-55, -33, -53) , ] #Removing the units with the highest impact until there is no more gaps in the histogram of Cooks distances
newfit <- lm(Price ~ Age + Distance, #New regression without potentially problematic units
data = mydata_new)
mydata_new$StdResid <- round(rstandard(newfit), 3) #New standardized residuals
mydata_new$CooksD <- round(cooks.distance(newfit), 3) #New Cooks distances
hist(mydata_new$StdResid, #to check outliers
xlab = "Standardized residuals",
ylab = "Frequency",
main = "Histogram of standardized residuals")
hist(mydata_new$CooksD, #to check units with high influence
xlab = "Cooks distance",
ylab = "Frequency",
main = "Histogram of Cooks distances")
There is still no outliers and the histogram of Cooks distances shows no gaps between units anymore.
9. Check for potential heteroskedasticity with scatterplot between standarized residuals and standardized fitted values. Explain the findings. (1 p)
mydata_new$StdFitted <- scale(newfit$fitted.values)
library(car)
scatterplot(y = mydata_new$StdResid, x = mydata_new$StdFitted,
ylab = "Standardized residuals",
xlab = "Standardized fitted values",
boxplots = FALSE,
regLine = FALSE,
smooth = FALSE)
By looking at the scatterplot above, standardized residuals seem to have more or less constant variances across different standardized fitted values, suggesting homoskedasticity. However, to make sure that it is the case, I will use Breusch-Pagan heteroskedasticity test:
library(olsrr)##
## Attaching package: 'olsrr'## The following object is masked from 'package:datasets':
##
## riversols_test_breusch_pagan(newfit)##
## Breusch Pagan Test for Heteroskedasticity
## -----------------------------------------
## Ho: the variance is constant
## Ha: the variance is not constant
##
## Data
## ---------------------------------
## Response : Price
## Variables: fitted values of Price
##
## Test Summary
## ----------------------------
## DF = 1
## Chi2 = 1.64762
## Prob > Chi2 = 0.1992832The test suggests homoskedasticity (p=0.1992832).
10. Are standardized residuals distributed normally? Show the graph and formally test it. Explain the findings. (1 p)
hist(mydata_new$StdResid, #graph
xlab = "Standardized residuals",
ylab = "Frequency",
main = "Histogram of standardized residuals")
shapiro.test(mydata_new$StdResid) #formal test##
## Shapiro-Wilk normality test
##
## data: mydata_new$StdResid
## W = 0.94297, p-value = 0.001286Both the histogram and Shapiro-Wilk normality test (with p=0.002) suggest rejecting the normality of the distribution of standardized residuals.
11. Estimate the fit2 again without potentially excluded units and show the summary of the model. Explain all coefficients. (2 p)
fit2 <- lm(Price ~ Age + Distance,
data = mydata_new)
summary(fit2)##
## Call:
## lm(formula = Price ~ Age + Distance, data = mydata_new)
##
## Residuals:
## Min 1Q Median 3Q Max
## -408.72 -217.93 -38.86 220.60 506.18
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2492.715 75.544 32.997 < 2e-16 ***
## Age -7.496 3.169 -2.365 0.0205 *
## Distance -24.574 2.700 -9.100 6.9e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 259.4 on 78 degrees of freedom
## Multiple R-squared: 0.5324, Adjusted R-squared: 0.5204
## F-statistic: 44.4 on 2 and 78 DF, p-value: 1.335e-13sqrt(summary(fit2)$r.squared) #Coefficient of multiple correlation## [1] 0.7296492Estimate of regression coefficients:
- Regression coefficient estimate of Age = -7.496: If age of the apartment increases by 1 year, the price per m2 of the apartment decreases on average by 7.496 units (euros, dollars, etc.), holding distance constant (p=0.021).
- Regression coefficient estimate of Distance = -24.574: If the distance of the apartment from city center increases by 1 km, the price per m2 of the apartment decreases on average by 24.574 units (euros, dollars, etc.), holding the age of the apartment constant (p<0.001).
Coefficient of multiple correlation = 0.7296492: Linear relationship between price, age, and distance is strong.
Coefficient of determination = 0.5324: 53.24% of variability of price can be explained by linear effect of age and distance. All other factors not included in the model explain the remaining 46.76%.
12. Estimate the linear regression function Price = f(Age, Distance, Parking and Balcony). Be careful to correctly include categorical variables. Save the object named fit3.
fit3 <- lm(Price ~ Age + Distance + ParkingF + BalconyF,
data = mydata_new)13. With function anova check if model fit3 fits data better than model fit2. (1 p)
anova(fit2, fit3)## Analysis of Variance Table
##
## Model 1: Price ~ Age + Distance
## Model 2: Price ~ Age + Distance + ParkingF + BalconyF
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 78 5248925
## 2 76 4915944 2 332981 2.5739 0.08287 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1There is not enough statistically significant evidence at 5% to claim that model fit3 fits data better than model fit2 (p=0.083).
14. Show the results of fit3 and explain regression coefficient for both categorical variables. Can you write down the hypothesis which is being tested with F-statistics, shown at the bottom of the output? (2 p)
summary(fit3)##
## Call:
## lm(formula = Price ~ Age + Distance + ParkingF + BalconyF, data = mydata_new)
##
## Residuals:
## Min 1Q Median 3Q Max
## -403.16 -204.08 -41.24 239.52 528.62
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2372.693 93.344 25.419 < 2e-16 ***
## Age -6.906 3.118 -2.215 0.0298 *
## Distance -22.254 2.839 -7.837 2.25e-11 ***
## ParkingFwith_parking 137.479 60.854 2.259 0.0267 *
## BalconyFwith_balcony 17.234 57.099 0.302 0.7636
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 254.3 on 76 degrees of freedom
## Multiple R-squared: 0.5621, Adjusted R-squared: 0.539
## F-statistic: 24.38 on 4 and 76 DF, p-value: 5.29e-13Regression coefficient estimate of ParkingF = 137.479: Given the values of all other variables, price per m2 for apartments with parking space is on average 137.479 units (euros, dollars, etc.) more expensive compared to the price per m2 for apartments without parking space (p=0.027)
Regression coefficient estimate of BalconyF = 17.234: This estimate is not statistically significant (p=0.764), making it meaningless to interpret it. However, if it was statistically significant it would mean that given the values of all other variables, price per m2 for apartments with balcony is on average 17.234 units (euros, dollars, etc.) more expensive compared to the apartments without balcony.
The hypothesis which is being tested with F-statistics:
H0 (Null): ρ^2 = 0 (The explanatory variables which are age, distance, having parking space, and having balcony do not explain the differences in price.)
H1 (Alternative): ρ^2 > 0 (At least one of the explanatory variables which are age, distance, having parking space, and having balcony explains the differences in price.)
15. Save fitted values and calculate the residual for apartment ID2. (1 p)
fitted_values <- fit3$fitted.values
residual_ID2 <- mydata_new[2,3]-fitted_values[2]
residual_ID2## 2
## 436.3894