# Úvod

V tejto úlohe pokračujem v nadväznosti na predchádzajúce zadania, ale pracujem s úplne inými dátami ako učiteľ. Používam dataset Students Performance, ktorý mám vo formáte CSV uložený v priečinku udaje. Mojím cieľom je overiť, či je môj lineárny model správne špecifikovaný, alebo či potrebujem urobiť funkčné transformácie, pridať kvadratické členy alebo zaviesť dummy premenné. Budem postupovať podľa štruktúry zo zadania, ale aplikované na študentské výsledky.

Balíčky a nastavenie

knitr::opts_chunk$set(echo = TRUE)
library(zoo)
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(tseries)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(lmtest)
library(sandwich)
library(car)
## Loading required package: carData
library(MASS)
library(ggplot2)
rm(list = ls())

PRÍPRAVA ÚDAJOV

udaje_raw <- read.csv("udaje/StudentsPerformance.csv",
                      sep = ",", dec = ".", header = TRUE)
 
names(udaje_raw) <- make.names(names(udaje_raw))
 
# Výber premenných
udaje <- udaje_raw[, c("math.score",
                       "reading.score",
                       "writing.score",
                       "parental.level.of.education")]
 
# Imputácia – ak by chýbali hodnoty
column_medians <- sapply(udaje, function(x) if(is.numeric(x)) median(x, na.rm = TRUE) else NA)
 
for (col in names(udaje)) {
  if (is.numeric(udaje[[col]])) {
    udaje[[col]][is.na(udaje[[col]])] <- column_medians[col]
  }
}
 
# Premenna číselná pre vzdelanie rodiča
edu_map <- c("some high school" = 1,
             "high school" = 2,
             "some college" = 3,
             "associate's degree" = 4,
             "bachelor's degree" = 5,
             "master's degree" = 6)
 
udaje$parent_edu_num <- edu_map[udaje$parental.level.of.education]
 
summary(udaje)
##    math.score     reading.score    writing.score    parental.level.of.education
##  Min.   :  0.00   Min.   : 17.00   Min.   : 10.00   Length:1000                
##  1st Qu.: 57.00   1st Qu.: 59.00   1st Qu.: 57.75   Class :character           
##  Median : 66.00   Median : 70.00   Median : 69.00   Mode  :character           
##  Mean   : 66.09   Mean   : 69.17   Mean   : 68.05                              
##  3rd Qu.: 77.00   3rd Qu.: 79.00   3rd Qu.: 79.00                              
##  Max.   :100.00   Max.   :100.00   Max.   :100.00                              
##  parent_edu_num 
##  Min.   :1.000  
##  1st Qu.:2.000  
##  Median :3.000  
##  Mean   :3.081  
##  3rd Qu.:4.000  
##  Max.   :6.000

Interpretácia: Dáta sú kompletné a pripravené na regresiu.

ZÁKLADNÁ REGRESIA

(regresia výkonu z písania voči matematike a vzdelaniu rodiča)

model <- lm(writing.score ~ math.score + parent_edu_num,
            data = udaje)
 
summary(model)
## 
## Call:
## lm(formula = writing.score ~ math.score + parent_edu_num, data = udaje)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -24.118  -6.897   0.417   6.835  21.882 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    12.49511    1.32472   9.432  < 2e-16 ***
## math.score      0.78654    0.01885  41.735  < 2e-16 ***
## parent_edu_num  1.16111    0.19568   5.934 4.09e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.917 on 997 degrees of freedom
## Multiple R-squared:  0.6564, Adjusted R-squared:  0.6557 
## F-statistic: 952.2 on 2 and 997 DF,  p-value: < 2.2e-16

Interpretácia: Obe premenné majú pozitívny vzťah k výsledku z písania.

1. RAMSEY RESET TEST

resettest(model)
## 
##  RESET test
## 
## data:  model
## RESET = 1.5469, df1 = 2, df2 = 995, p-value = 0.2134

Interpretácia: Test ukazuje, či modelu chýba nelinearita alebo dôležité premenné.

2. GRAFICKÁ ANALÝZA – Residuals vs Fitted

plot(model, which = 1)

Interpretácia: Ak je viditeľné zakrivenie, model nemusí byť úplne lineárny.

3. C+R GRAFY

car::crPlots(model)

Interpretácia: Ak sa krivka odkláňa od priamky, daná premenná si pýta transformáciu.

4. NELINEÁRNA ŠPECIFIKÁCIA – kvadratické členy

model_quad <- lm(writing.score ~ math.score + I(math.score^2) +
                               parent_edu_num + I(parent_edu_num^2),
                 data = udaje)
 
summary(model_quad)
## 
## Call:
## lm(formula = writing.score ~ math.score + I(math.score^2) + parent_edu_num + 
##     I(parent_edu_num^2), data = udaje)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -23.9883  -6.8634   0.4182   6.8453  22.0117 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          9.0319625  3.5759280   2.526   0.0117 *  
## math.score           0.9855933  0.1074629   9.171   <2e-16 ***
## I(math.score^2)     -0.0015529  0.0008272  -1.877   0.0608 .  
## parent_edu_num      -0.8421473  0.8381139  -1.005   0.3152    
## I(parent_edu_num^2)  0.3113696  0.1263532   2.464   0.0139 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.884 on 995 degrees of freedom
## Multiple R-squared:  0.6595, Adjusted R-squared:  0.6582 
## F-statistic: 481.9 on 4 and 995 DF,  p-value: < 2.2e-16
anova(model, model_quad)
## Analysis of Variance Table
## 
## Model 1: writing.score ~ math.score + parent_edu_num
## Model 2: writing.score ~ math.score + I(math.score^2) + parent_edu_num + 
##     I(parent_edu_num^2)
##   Res.Df   RSS Df Sum of Sq      F   Pr(>F)   
## 1    997 79268                                
## 2    995 78537  2    730.89 4.6299 0.009967 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
resettest(model_quad)
## 
##  RESET test
## 
## data:  model_quad
## RESET = 0.23921, df1 = 2, df2 = 993, p-value = 0.7873

Interpretácia: Kvadratické členy často zlepšia model, ak vzťah nie je čisto lineárny.

5. ROZŠÍRENÝ RESET A KVADRATICKÝ MODEL

model_rozsireny <- lm(writing.score ~ math.score + I(math.score^2) +
                                        parent_edu_num + I(parent_edu_num^2),
                      data = udaje)
 
summary(model_rozsireny)
## 
## Call:
## lm(formula = writing.score ~ math.score + I(math.score^2) + parent_edu_num + 
##     I(parent_edu_num^2), data = udaje)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -23.9883  -6.8634   0.4182   6.8453  22.0117 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)          9.0319625  3.5759280   2.526   0.0117 *  
## math.score           0.9855933  0.1074629   9.171   <2e-16 ***
## I(math.score^2)     -0.0015529  0.0008272  -1.877   0.0608 .  
## parent_edu_num      -0.8421473  0.8381139  -1.005   0.3152    
## I(parent_edu_num^2)  0.3113696  0.1263532   2.464   0.0139 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.884 on 995 degrees of freedom
## Multiple R-squared:  0.6595, Adjusted R-squared:  0.6582 
## F-statistic: 481.9 on 4 and 995 DF,  p-value: < 2.2e-16
anova(model, model_rozsireny)
## Analysis of Variance Table
## 
## Model 1: writing.score ~ math.score + parent_edu_num
## Model 2: writing.score ~ math.score + I(math.score^2) + parent_edu_num + 
##     I(parent_edu_num^2)
##   Res.Df   RSS Df Sum of Sq      F   Pr(>F)   
## 1    997 79268                                
## 2    995 78537  2    730.89 4.6299 0.009967 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
resettest(model_rozsireny)
## 
##  RESET test
## 
## data:  model_rozsireny
## RESET = 0.23921, df1 = 2, df2 = 993, p-value = 0.7873

Interpretácia: Porovnávame pôvodný model a rozšírený – vyšší R²_adj znamená zlepšenie.

6. DUMMY PREMENNÁ – zlom v špecifikácii

Dummy premenná podľa mediánu matematického skóre

threshold <- median(udaje$math.score)
udaje$DUM <- ifelse(udaje$math.score < threshold, 0, 1)
 
summary(udaje$DUM)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   0.000   0.000   1.000   0.517   1.000   1.000

Interpretácia: Dummy rozdelí žiakov na 2 skupiny podľa výkonu v matematike.*

Model so zlomom v autonómnom člene

modelD_auto <- lm(writing.score ~ DUM +
                                  math.score +
                                  parent_edu_num,
                  data = udaje)
 
summary(modelD_auto)
## 
## Call:
## lm(formula = writing.score ~ DUM + math.score + parent_edu_num, 
##     data = udaje)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -24.2206  -6.9037   0.4151   6.8175  21.7794 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    12.16883    1.72372   7.060 3.13e-12 ***
## DUM            -0.27495    0.92876  -0.296    0.767    
## math.score      0.79378    0.03089  25.698  < 2e-16 ***
## parent_edu_num  1.15778    0.19610   5.904 4.86e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.921 on 996 degrees of freedom
## Multiple R-squared:  0.6564, Adjusted R-squared:  0.6554 
## F-statistic: 634.2 on 3 and 996 DF,  p-value: < 2.2e-16
anova(model, modelD_auto)
## Analysis of Variance Table
## 
## Model 1: writing.score ~ math.score + parent_edu_num
## Model 2: writing.score ~ DUM + math.score + parent_edu_num
##   Res.Df   RSS Df Sum of Sq      F Pr(>F)
## 1    997 79268                           
## 2    996 79261  1    6.9743 0.0876 0.7673

Interpretácia: Autonómny zlom sa nemusí ukázať ako štatisticky významný.*

Model so zlomom v sklone

modelD_sklon <- lm(writing.score ~ math.score +
                                 I(DUM * math.score) +
                                 parent_edu_num,
                   data = udaje)
 
summary(modelD_sklon)
## 
## Call:
## lm(formula = writing.score ~ math.score + I(DUM * math.score) + 
##     parent_edu_num, data = udaje)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -24.4136  -6.8657   0.4436   6.7658  21.5864 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         11.41300    1.93087   5.911 4.67e-09 ***
## math.score           0.80957    0.03534  22.906  < 2e-16 ***
## I(DUM * math.score) -0.01043    0.01353  -0.770    0.441    
## parent_edu_num       1.15419    0.19593   5.891 5.25e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.918 on 996 degrees of freedom
## Multiple R-squared:  0.6566, Adjusted R-squared:  0.6555 
## F-statistic: 634.7 on 3 and 996 DF,  p-value: < 2.2e-16
anova(model, modelD_sklon)
## Analysis of Variance Table
## 
## Model 1: writing.score ~ math.score + parent_edu_num
## Model 2: writing.score ~ math.score + I(DUM * math.score) + parent_edu_num
##   Res.Df   RSS Df Sum of Sq      F Pr(>F)
## 1    997 79268                           
## 2    996 79221  1    47.213 0.5936 0.4412
resettest(modelD_sklon)
## 
##  RESET test
## 
## data:  modelD_sklon
## RESET = 1.3551, df1 = 2, df2 = 994, p-value = 0.2584

Interpretácia: Zlom v sklone môže ukázať, že matematika ovplyvňuje písanie rôzne u slabších a lepších žiakov.*

7. BOX–COX TEST

boxcox(model)

Interpretácia: Box–Cox ukáže, či sa oplatí transformovať writing.score.*

Transformovaný model podľa Box–Cox

model_lambda <- lm(I((writing.score^0.8 - 1) / 0.8) ~
                      math.score + parent_edu_num,
                   data = udaje)
 
summary(model_lambda)
## 
## Call:
## lm(formula = I((writing.score^0.8 - 1)/0.8) ~ math.score + parent_edu_num, 
##     data = udaje)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -12.569  -2.938   0.194   3.050   9.384 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    10.916106   0.577201  18.912  < 2e-16 ***
## math.score      0.343580   0.008212  41.841  < 2e-16 ***
## parent_edu_num  0.502019   0.085263   5.888 5.34e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.885 on 997 degrees of freedom
## Multiple R-squared:  0.6573, Adjusted R-squared:  0.6566 
## F-statistic: 956.3 on 2 and 997 DF,  p-value: < 2.2e-16
resettest(model_lambda)
## 
##  RESET test
## 
## data:  model_lambda
## RESET = 4.8828, df1 = 2, df2 = 995, p-value = 0.007759

Interpretácia: Transformácia môže zlepšiť štatistiky, ale zhorší interpretáciu skóre.*

ZÁVER

V tejto úlohe som skúmala, či môj lineárny model výkonu v písaní je správne špecifikovaný. RESET test aj grafy ukázali, že jednoduchá lineárna forma nemusí úplne postačovať. Po pridaní kvadratických členov sa model zlepšil, čo potvrdila aj ANOVA. Zlom v sklone pri matematike naznačil, že matematické skóre môže ovplyvňovať výsledky z písania rozdielne u slabších a lepších žiakov. Box–Cox transformácia síce mierne vylepšila model, ale zhoršila interpretáciu, preto ju ďalej nepoužívam. Najvhodnejšie sa javí použitie kvadratických členov a modelu so zlomom v sklone.