Task 1

a) Finding and explaining data set

I have found my data set in R with “data(Wong)” function, then I have used “write.csv(Wong,”Wong.csv”” function to save Wong data set into current working directory.

mydata <- read.table("./Wong.csv", header=TRUE, sep=",")

head(mydata)
##     id days duration sex      age piq viq
## 1 3358   30        4   0 20.67077  87  89
## 2 3535   16       17   0 55.28816  95  77
## 3 3547   40        1   0 55.91513  95 116
## 4 3592   13       10   0 61.66461  59  73
## 5 3728   19        6   0 30.12731  67  73
## 6 3790   13        3   0 57.06229  76  69

Explanation of variables:

  • id: Patient ID number

  • days: Number of days post coma at which IQs were measured

  • duration: Duration of the coma in days

  • sex: 0-male, 1-female

  • age: Age in years at the time of injury

  • piq: Performance in points - (i.e., mathematical) IQ

  • viq: Performance in points - verbal IQ

b) Data manipulation

# Changing sex variable into categorical factor
mydata$sex <- factor(mydata$sex, 
                         levels = c(0, 1), 
                         labels = c("Male", "Female"))
head(mydata)
##     id days duration  sex      age piq viq
## 1 3358   30        4 Male 20.67077  87  89
## 2 3535   16       17 Male 55.28816  95  77
## 3 3547   40        1 Male 55.91513  95 116
## 4 3592   13       10 Male 61.66461  59  73
## 5 3728   19        6 Male 30.12731  67  73
## 6 3790   13        3 Male 57.06229  76  69
# Checking if any n/a in the dataset
anyNA(mydata)
## [1] FALSE

Explanation: There are no n/a. Now we can perform further analysis. *any() checks if at least one of those is TRUE.

c) Describtive Statistics

library(psych)
summary(mydata[,-1]) 
##       days            duration         sex           age        
##  Min.   :   13.0   Min.   :  0.0   Male  :260   Min.   : 6.513  
##  1st Qu.:   59.0   1st Qu.:  1.0   Female: 71   1st Qu.:21.737  
##  Median :  150.0   Median :  7.0                Median :26.877  
##  Mean   :  481.9   Mean   : 14.3                Mean   :31.853  
##  3rd Qu.:  416.0   3rd Qu.: 16.0                3rd Qu.:40.923  
##  Max.   :11628.0   Max.   :255.0                Max.   :80.033  
##       piq              viq        
##  Min.   : 50.00   Min.   : 64.00  
##  1st Qu.: 77.00   1st Qu.: 85.00  
##  Median : 87.00   Median : 94.00  
##  Mean   : 87.56   Mean   : 94.96  
##  3rd Qu.: 97.00   3rd Qu.:105.00  
##  Max.   :133.00   Max.   :132.00

Explanation:

  • Mean for “duration”: On average patients were in coma for 14.3 days.

  • Median for “days”: Median is 150 days, it means that 50% of the patients were assessed on 150th day post-coma or before, while the other 50% were assessed after 150 days.

  • Maximum for “PIQ”: The maximum PIQ score in the sample was 133 points.

d) Graphical Distribution

i) Graphical Distribtion of PIQ using histogram

hist(mydata$piq,
     main = "Distribution of Performance IQ (PIQ)",
     xlab = "Performance IQ Score",
     ylab = "Number of Patients",
     col  = "lightblue",
     border = "white")

Explanation: As we can see from the histogram, the distribution of PIQ scores is slightly asymmetrical to the right or positively skewed.

ii) Graphical Distribution of PIQ using box plot

library(ggplot2)

ggplot(mydata, aes(x = sex, y = piq, fill = sex)) +
  geom_boxplot() +
  labs(title = "Boxplot of PIQ by Gender",
       x = "Gender",
       y = "Performance IQ") +
  theme_minimal() 

Explanation of boxplot: The boxplots suggest that the PIQ distribution is roughly symmetric for both genders, both having high value outliers. We can say that both boxplots are fairly simillar, only what I can see is that there is less variability in points for men. I will further use additional describtive statistics to analyze graph better.

Using describtive statistics for better interpretation of boxplot: Median:

median(mydata$piq[mydata$sex == "Female"])
## [1] 87
median(mydata$piq[mydata$sex == "Male"])
## [1] 87

Explanation: Both male and female patients have a median PIQ of 87, meaning that 50% of the patients in each group scored 87 points and below while the other 50% have scored above 87 points.

Mean:

round(mean(mydata$piq[mydata$sex == "Female"]),2)
## [1] 89.18
round(mean(mydata$piq[mydata$sex == "Male"]),2)
## [1] 87.11

Explanation: On average females have scored 2.07 points higher than males.

Standard Deviation:

round(sd(mydata$piq[mydata$sex == "Female"]),2)
## [1] 18
round(sd(mydata$piq[mydata$sex == "Male"]),2)
## [1] 14.26

Explanation: On the other hand there is difference in variability between Males and Females. Male PIQ scores are less spread out, indicating their scores are more clustered around the mean.

Conclusion: On average, females achieved higher PIQ scores compared to males; however, their scores exhibit greater variability.

Task 2

library(readxl)

# Importing dataset
mydata <- read_xlsx("~/Desktop/IMB/Bootcamp - R/IMB Bootcamp/R data/R - Task 2/Business School.xlsx")

mydata <- as.data.frame(mydata)

head(mydata)
##   Student ID Undergrad Degree Undergrad Grade MBA Grade
## 1          1         Business            68.4      90.2
## 2          2 Computer Science            70.2      68.7
## 3          3          Finance            76.4      83.3
## 4          4         Business            82.6      88.7
## 5          5          Finance            76.9      75.4
## 6          6 Computer Science            83.3      82.1
##   Work Experience Employability (Before) Employability (After) Status
## 1              No                    252                   276 Placed
## 2             Yes                    101                   119 Placed
## 3              No                    401                   462 Placed
## 4              No                    287                   342 Placed
## 5              No                    275                   347 Placed
## 6              No                    254                   313 Placed
##   Annual Salary
## 1        111000
## 2        107000
## 3        109000
## 4        148000
## 5        255500
## 6        103500

a) Graphical distribution of Undergrad degrees

library(ggplot2)

ggplot(mydata, aes(x = `Undergrad Degree`, fill = `Undergrad Degree`)) +
  geom_bar() +
  labs(title = "Distribution of Undergraduate Degrees",
       x = "Undergraduate Degree",
       y = "Number of Students") +
  theme_minimal() 

Answer: Graph is showing that Business Undergrad degree is the most common in the current generation of MBA students.

b) Descriptive statistics of the Annual Salary and its distribution with the histogram

library(ggplot2)

ggplot(mydata, aes(x = `Annual Salary`)) +
  geom_histogram(binwidth = 10000, fill = "skyblue", color = "black") +
  labs(title = "Distribution of Annual Salary",
       x = "Annual Salary",
       y = "Number of Students") +
  theme_minimal()

library(psych)
describe(mydata$`Annual Salary`)
##    vars   n   mean       sd median  trimmed     mad   min    max
## X1    1 100 109058 41501.49 103500 104600.2 25945.5 20000 340000
##     range skew kurtosis      se
## X1 320000 2.22     9.41 4150.15

Answer: The mean annual salary is about 109,058, while the median is slightly lower at 103,500. This gap between mean and median as well as graphical distribution tell us that we have a positively skewed distribution — a few students earn very high salaries (up to 340,000), pulling the mean upwards.

c) Testing the hypothesis with t-test

# t-test
t.test(mydata$`MBA Grade`, mu = 74, alternative = "two.sided")
## 
##  One Sample t-test
## 
## data:  mydata$`MBA Grade`
## t = 2.6587, df = 99, p-value = 0.00915
## alternative hypothesis: true mean is not equal to 74
## 95 percent confidence interval:
##  74.51764 77.56346
## sample estimates:
## mean of x 
##  76.04055

Answer: P-value is 0.00915, we can reject null hypothesis at 1% significance level and mean of this year class is 76.04. This means that this year average MBA grade is statistically significantly higher than previous year average grade of 74.

# Cohen's d - standardized measure of effect size
mean_val <- mean(mydata$`MBA Grade`, )
sd_val   <- sd(mydata$`MBA Grade`, )
mu0 <- 74


cohen_d <- (mean_val - mu0) / sd_val
round(cohen_d,2)
## [1] 0.27

Answer: Cohen’s* d is a standardized measure of effect size that expresses the mean difference in units of standard deviation. Cohen’s d = 0.27, indicates a small effect.

Conclusion: Results suggests that although the grade improvement is statistically significant, the practical difference from the hypothesized mean of 74 is relatively modest.

*Cohen, J. (1988). Statistical power analysis for the behavioral sciences (2nd ed., pp. 19–25). Hillsdale, NJ: Lawrence Erlbaum Associates.

Task 3

a) Import the dataset Apartments.xlsx

library(readxl)

mydata <- read_xlsx("./Apartments.xlsx")

mydata <- as.data.frame(mydata)

mydata$ID <- 1:nrow(mydata) #An ID column was added to index the rows for easier referencing

head(mydata)
##   Age Distance Price Parking Balcony ID
## 1   7       28  1640       0       1  1
## 2  18        1  2800       1       0  2
## 3   7       28  1660       0       0  3
## 4  28       29  1850       0       1  4
## 5  18       18  1640       1       1  5
## 6  28       12  1770       0       1  6

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

  • ID: Apartment ID number

b) Change categorical variables into factors.

mydata$Parking <- factor(mydata$Parking, 
                         levels = c(0, 1), 
                         labels = c("No", "Yes"))

mydata$Balcony <- factor(mydata$Balcony, 
                         levels = c(0, 1), 
                         labels = c("No", "Yes"))
head(mydata)
##   Age Distance Price Parking Balcony ID
## 1   7       28  1640      No     Yes  1
## 2  18        1  2800     Yes      No  2
## 3   7       28  1660      No      No  3
## 4  28       29  1850      No     Yes  4
## 5  18       18  1640     Yes     Yes  5
## 6  28       12  1770      No     Yes  6

c) Test the hypothesis H0: Mu_Price = 1900 eur. What can you conclude?

t.test(mydata$Price, mu = 1900)
## 
##  One Sample t-test
## 
## data:  mydata$Price
## t = 2.9022, df = 84, p-value = 0.004731
## alternative hypothesis: true mean is not equal to 1900
## 95 percent confidence interval:
##  1937.443 2100.440
## sample estimates:
## mean of x 
##  2018.941

Answer: Since p < 0.01, we reject the null hypothesis at the 1% significance level. The test shows evidence that the mean apartment price per m2 is not equal to 1900€.

d) 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.

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.03401

Answers:

  • Regression coefficient: If the age of an apartment goes up by 1 year, the apartment price per m2 goes down by 8.98 €, on average.This negative coefficient is statistically significant at 5% level (p = 0.034).

  • Coefficient of determination: R2 = 0.053 means that only 5.3% of the variability in apartment prices is explained by the linear effect of age. The rest of variability is explained by other factors not included in this simple model (such as location, size, parking, balcony, etc.).

round(-sqrt(0.05302),2)
## [1] -0.23
  • Coefficient of correlation: This shows a weak and negative correlation between Age and Price.

e) Show the scateerplot matrix between Price, Age and Distance. Based on the matrix determine if there is potential problem with multicolinearity.

library(car, carData)

scatterplotMatrix(mydata[ , c(-4, -5, -6)], 
                  smooth = FALSE) 

Answer: Based on scatterplot matrix there is no problem with multicolinearity, because the scatterplot between distance and age shows points widely scattered with almost no slope.

f) Estimate the multiple regression function: Price = f(Age, Distance). Save it in object named fit2.

fit2 <- lm(Price ~ Age + Distance, data = mydata)

g) Check the multicolinearity with VIF statistics. Explain the findings.

vif(fit2) 
##      Age Distance 
## 1.001845 1.001845
mean(vif(fit2))
## [1] 1.001845

Answer: Since both predictors have VIF close to 1, there is no multicolinearity problem in this model. Multicolinearity would be a problem if VIF greater than 5.

h) Calculate standardized residuals and Cooks Distances for model fit2. Remove any potentially problematic units (outliers or units with high influence).

mydata$StdResid <- round(rstandard(fit2), 2) 
mydata$CooksD <- round(cooks.distance(fit2), 2) 
head(mydata[order(mydata$StdResid),], 10) 
##    Age Distance Price Parking Balcony ID StdResid CooksD
## 53   7        2  1760      No     Yes 53    -2.15   0.07
## 13  12       14  1650      No     Yes 13    -1.50   0.01
## 72  12       14  1650      No      No 72    -1.50   0.01
## 20  13        8  1800      No      No 20    -1.38   0.01
## 35  14       16  1660      No     Yes 35    -1.26   0.01
## 36  24        5  1830     Yes      No 36    -1.19   0.01
## 54  30       17  1560      No      No 54    -1.10   0.01
## 5   18       18  1640     Yes     Yes  5    -1.07   0.01
## 28  18       19  1620     Yes      No 28    -1.07   0.01
## 64  18       18  1640     Yes     Yes 64    -1.07   0.01
head(mydata[order(-mydata$StdResid),], 10) 
##    Age Distance Price Parking Balcony ID StdResid CooksD
## 38   5       45  2180     Yes     Yes 38     2.58   0.32
## 33   2       11  2790     Yes      No 33     2.05   0.07
## 2   18        1  2800     Yes      No  2     1.78   0.03
## 61  18        1  2800     Yes     Yes 61     1.78   0.03
## 58   8        2  2820     Yes      No 58     1.66   0.04
## 57  10        1  2810      No      No 57     1.60   0.03
## 22  37        3  2540     Yes     Yes 22     1.58   0.06
## 25   8       26  2300     Yes     Yes 25     1.57   0.03
## 11  12        5  2700     Yes      No 11     1.55   0.02
## 27  16        1  2750     Yes      No 27     1.55   0.02
head(mydata[order(-mydata$CooksD),], 10)
##    Age Distance Price Parking Balcony ID StdResid CooksD
## 38   5       45  2180     Yes     Yes 38     2.58   0.32
## 55  43       37  1740      No      No 55     1.44   0.10
## 33   2       11  2790     Yes      No 33     2.05   0.07
## 53   7        2  1760      No     Yes 53    -2.15   0.07
## 22  37        3  2540     Yes     Yes 22     1.58   0.06
## 39  40        2  2400      No     Yes 39     1.09   0.04
## 58   8        2  2820     Yes      No 58     1.66   0.04
## 2   18        1  2800     Yes      No  2     1.78   0.03
## 25   8       26  2300     Yes     Yes 25     1.57   0.03
## 31  45       21  1910      No     Yes 31     0.89   0.03
hist(mydata$StdResid, 
     xlab = "Standardized residuals", 
     ylab = "Frequency", 
     main = "Histogram of standardized residuals")

hist(mydata$CooksD, 
     xlab = "Cooks distance", 
     ylab = "Frequency", 
     main = "Histogram of Cooks distances")

Explanation: Since no standardized residuals exceed ±3, there are no clear outliers, but observation 38 has a Cook’s Distance of 0.32 indicating high influence on the model, so I have decided to remove it.

Removing ID 38:

library(dplyr)

mydata <- mydata %>%
filter(!ID == "38")

i) Check for potential heteroskedasticity with scatterplot between standarized residuals and standrdized fitted values. Explain the findings.

library(ggplot2)

fit2 <- lm(Price ~ Age + Distance, data = mydata)
mydata$StdFittedValues <- scale(fit2$fitted.values)

scatterplot(y = mydata$StdResid, x = mydata$StdFittedValues,
ylab = "Standardized residuals",
xlab = "Standardized fitted values",
boxplots = FALSE,
regLine = FALSE,
smooth = FALSE)

Explanation: The scatterplot of standardized residuals against standardized fitted values shows no clear pattern. The residuals appear randomly distributed around zero with roughly constant variance. This suggests that the assumption of homoskedasticity holds for our regression model.

j) Are standardized residuals ditributed normally? Show the graph and formally test it. Explain the findings.

hist(mydata$StdResid, 
     xlab = "Standardized residuals", 
     ylab = "Frequency", 
     main = "Histogram of standardized residuals")

shapiro.test(mydata$StdResid)
## 
##  Shapiro-Wilk normality test
## 
## data:  mydata$StdResid
## W = 0.9488, p-value = 0.00219

Answer: The histogram shows that standardized residuals are roughly normally distributed but not perfectly normal. The Shapiro–Wilk test (p = 0.00219 < 0.05) indicates that we reject the null hypothesis of normality, meaning the residuals are not normally distributed.

Conclusion: Although the assumption of normality is formally violated, the sample size is sufficiently large for the Central Limit Theorem to apply. Therefore, based on the theory we covered, this violation is not a major concern and the regression results remain reliable.

k) Estimate the fit2 again without potentially excluded units and show the summary of the model. Explain all coefficients.

fit2 <- lm(Price ~ Age + Distance, data = mydata)

summary(fit2)
## 
## Call:
## lm(formula = Price ~ Age + Distance, data = mydata)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -604.92 -229.63  -56.49  192.97  599.35 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 2456.076     73.931  33.221  < 2e-16 ***
## Age           -6.464      3.159  -2.046    0.044 *  
## Distance     -22.955      2.786  -8.240 2.52e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 276.1 on 81 degrees of freedom
## Multiple R-squared:  0.4838, Adjusted R-squared:  0.4711 
## F-statistic: 37.96 on 2 and 81 DF,  p-value: 2.339e-12

Explanation:

  • Intercept: Intercept of 2456.1 is the predicted price of an apartment when both Age and Distance are equal to 0. It has no explanatory meaning in this case.

  • Age: Holding Distance constant, each additional year of apartment age decreases the price by 6.46 €/m2, on average. The effect is statistically significant at the 5% level (p=0.044)

  • Distance: Holding Age constant, each additional kilometer away from the city center decreases the price by about 22.96 €/m2, on average. This effect is statistically significant at 0.1% level, since p=2.52e-12.

  • Coefficient of determination: R2 = 0.484 means that 48.4% of the variability in apartment prices is explained by the linear effects of age and distance. Compared to the model with only age (R2 = 0.053), the inclusion of distance improves explanatory power, indicating that distance is an important predictor of apartment prices.

l) 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 + Parking + Balcony, data = mydata)

m) With function anova check if model fit3 fits data better than model fit2.

anova(fit2,fit3)
## Analysis of Variance Table
## 
## Model 1: Price ~ Age + Distance
## Model 2: Price ~ Age + Distance + Parking + Balcony
##   Res.Df     RSS Df Sum of Sq      F  Pr(>F)  
## 1     81 6176767                              
## 2     79 5654480  2    522287 3.6485 0.03051 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Answer: The p-value (0.03051) is below 0.05 significance level, which means that fit3 explains significantly more variance in apartment prices than fit2. Adding Parking and Balcony improves the model and model fit3 should be preferred over fit2.

n) 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?

summary(fit3)
## 
## Call:
## lm(formula = Price ~ Age + Distance + Parking + Balcony, data = mydata)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -473.21 -192.37  -28.89  204.17  558.77 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 2329.724     93.066  25.033  < 2e-16 ***
## Age           -5.821      3.074  -1.894  0.06190 .  
## Distance     -20.279      2.886  -7.026 6.66e-10 ***
## ParkingYes   167.531     62.864   2.665  0.00933 ** 
## BalconyYes   -15.207     59.201  -0.257  0.79795    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 267.5 on 79 degrees of freedom
## Multiple R-squared:  0.5275, Adjusted R-squared:  0.5035 
## F-statistic: 22.04 on 4 and 79 DF,  p-value: 3.018e-12

Explanation of regression coefficient for both categorical variables:

  • Parking: Apartments with parking are priced, on average, 167.53 €/m2 higher compared to apartments without parking, holding other variables constant. This effect is statistically significant at 1% level (p=0.00933).

  • Balcony: The coefficient for Balcony is not statistically significant (p = 0.80), which means we cannot conclude that having a balcony leads to any difference in apartment prices compared to apartments without a balcony.

F-Statistics part:

  • Null hypothesis (H0): All slope coefficients are equal to 0 (None of the predictors explain any variation in apartment prices)

  • Alternative hypothesis (H1): At least one slope coefficient is not equal to 0

Since the F-statistic is 22.04 with p < 0.001, we reject H0 at 0.1% significance level and conclude that at least one predictor explain a significant amount of variation in apartment prices.

o) Save fitted values and claculate the residual for apartment ID2.

mydata$fitted_fit3 <- fitted(fit3)
mydata$resid_fit3  <- resid(fit3)

mydata[mydata$ID == 2, c("Price", "fitted_fit3", "resid_fit3")]
##   Price fitted_fit3 resid_fit3
## 2  2800    2372.197   427.8029

Explanation: For apartment ID2, the actual price is 2800€, while the model fit3 predicts a price of 2372.20 €. The residual is therefore +427.80, meaning that the apartment is priced 427.80 € higher than what the model predicts. This positive residual indicates that the model underestimated the apartment’s price.