Warning: package 'ggplot2' was built under R version 4.5.28 Сызықтық регрессия
Сызықтық регрессия
Мысалға қарайық: бізде екі айнымалы бар:
\(X\) -
sleep- кездейсоқ студенттің тәуліктік орташа ұйқы уақыты\(Y\) -
GPA- оның GPA
| sleep | GPA |
|---|---|
| 5.373546 | 2.541726 |
| 6.183643 | 2.925353 |
| 5.164371 | 2.418771 |
| 7.595281 | 3.386143 |
| 6.329508 | 2.831231 |
| 5.179532 | 2.959112 |
Визуализация
- Шашыраңқы диаграмманы қарастырайық
Визуалды талдау
- График трендті анық көрсетеді, оны елестету үшін қолданып көрейік
Warning in geom_abline(intercept = 0, qiyalik = 0.27, color = "tomato"):
Ignoring unknown parameters: `qiyalik`Warning in geom_abline(intercept = 1, горизонт = 0.27, color = "steelblue"):
Ignoring unknown parameters: `горизонт`Warning in geom_abline(intercept = 1, горизонт = 0.33, color = "wheat3"):
Ignoring unknown parameters: `горизонт`
- Сіздің ойыңызша, қай сызық трендті ең жақсы сипаттайды? Және неге?
- Қызыл
- Көк
- Бидай
Кішкене алгебра
- Түзу сызық - бұл негізінен деректер моделі: \(GPA = \beta_0 + \beta_1*ұйқы\)
- \(\beta_0\) - y-қиылысу нүктесі
- \(\beta_1\) - көлбеу (яғни көлбеу)
- Түзу сызық сызған кезде, біз \(\beta_1\) және \(\beta_0\) параметрлерін көрсетеміз
- \(\beta_1\), көлбеуді келесідей түсіндіре аламыз: әрбір қосымша \(1\) сағат ұйқы үшін студент орташа есеппен GPA-сын \(\beta_1\)-қа жақсартады
- Бұл жақсы, бірақ сіз қалағаныңызша көп түзу сызықтар сыза аласыз. Қайсысы деректерді ең жақсы сипаттайтынын қалай шешесіз?
- blue деректерге red қарағанда жақсырақ сәйкес келетіні интуитивті.
- Дегенмен, blue және wheat сызықтары арасындағы айырмашылық онша айқын емес.
- Сызықтың деректерге қаншалықты сәйкес келетінін сандық өлшеммен өлшеу жақсы болар еді.
- Ойланайық. Біз GPA мәндерін, \(y_i\), сондай-ақ болжамды мәндерді, \(\hat{y}_i\) байқадық.
- Олардың арасындағы айырмашылық болжам қателігі, \(\epsilon_i = y_i - \hat{y}_i\).
- Бұл қателіктер residuals[^chapter8-1] (residuals) деп аталады.
- Барлық қалдықтарды алып, оларды квадраттап, бірге қоссақ не болады? Бұл қосындыны \(J\) деп белгілейік.
- \(J(\beta_0, \beta_1) = \sum_{i=1}^{n}(\epsilon^2) = \sum_{i=1}^{n}(y_i - \hat{y}_i)^2 = \sum_{i=1}^{n}(y_i - \beta_0 - \beta_1 x_i)^2\)
- Бұл сызықтың деректерге жақындығын сипаттаудың жақсы тәсілі сияқты.
- \(J\)-тың кіші мәндері қателіктеріміздің қосындысының аз екенін білдіреді, яғни сызық деректерге жақсы сәйкес келеді.
- \(J\)-тың үлкен мәндері қателіктеріміздің қосындысының үлкен екенін білдіреді, яғни сызық деректерге нашар сәйкес келеді.
- \(J\) жоғалту функциясы деп аталады.
- Енді осы метриканы пайдаланып, үш сызықты салыстырайық.
Тексерейік
- Бізде үш модель бар
- 1-модель, қызыл сызық: \(\beta_0 = 0, \beta_1 = 0.27\) немесе \(\hat{Y}_{red} = 0 + 0.27X\)
- 2-модель, көк сызық: \(\beta_0 = 1, \beta_1 = 0.27\) немесе \(\hat{Y}_{blue} = 1 + 0.27X\)
- 3-модель, бидай сызығы: \(\beta_0 = 1, \beta_1 = 0.33\) немесе \(\hat{Y}_{бидай} = 1 + 0.33X\)
students_sleep_reg <-
students_sleep %>%
mutate(res_red = GPA - (0 + 0.27*sleep),
res_blue = GPA - (1 + 0.27*sleep),
res_wheat = GPA - (1 + 0.33*sleep))
students_sleep_reg %>%
summarise(j_red = sum(res_red^2),
j_blue = sum(res_blue^2),
j_wheat = sum(res_wheat^2))| j_red | j_blue | j_wheat |
|---|---|---|
| 156.7258 | 9.366998 | 5.349284 |
Ең жақсысын қалай табуға болады?
Жарайды, алғашқы үш модельден біреуін таңдай аламыз. Бірақ шексіз көп модель бар ма? Абсолютті ең жақсы сызық бар ма?
Бұл нені білдіреді? \(J(\beta_0, \beta_1)\) минимумы бар ма?
Бұл квадраттардың қосындысы, сондықтан иә, солай, және онда тек біреуі ғана бар.
Сондықтан математикада барлығын квадраттағанды ұнатады.
Біз \(J(\beta_0^{'}, \beta_1^{'})\) минималды мәнді қабылдайтындай \(\beta_0^{'}\) және \(\beta_1^{'}\) таба аламыз.
Немесе, “математикалық” терминдермен айтқанда, \(\arg\min_{\beta_0, \beta_1}(J(\beta_0, \beta_1)) = \arg\min_{\beta_0, \beta_1}(\sum_{i=1}^{n}(y_i - \beta_0 - \beta_1 x_i)^2)\).
Шығын функциясын минимумдайтын мәндерді табу оңай; тек жартылай туындыларды алып, оларды нөлге теңеп, сызықтық теңдеулер жүйесін шешіңіз. = 0 \ = 0
$$
R құтқаруға келеді
Call:
lm(formula = GPA ~ sleep, data = students_sleep)
Residuals:
Min 1Q Median 3Q Max
-0.3754 -0.1227 -0.0279 0.1079 0.4692
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.99373 0.13302 7.471 3.36e-11 ***
sleep 0.30979 0.02155 14.378 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1926 on 98 degrees of freedom
Multiple R-squared: 0.6784, Adjusted R-squared: 0.6751
F-statistic: 206.7 on 1 and 98 DF, p-value: < 2.2e-16- Енді есептелген коэффициенттерді қабылдай аламыз
Warning in geom_abline(intercept = 0, slop = 0.27, color = "tomato", alpha =
0.3): Ignoring unknown parameters: `slop`Warning in geom_abline(intercept = 1, slop = 0.27, color = "steelblue", :
Ignoring unknown parameters: `slop`Warning in geom_abline(intercept = 1, көлбеу = 0.33, түс = "бидай3", : Ignoring
unknown parameters: `көлбеу`, `түс`, and `альфа`Warning in geom_abline(intercept = the_beta0, көлбеу = the_beta1, : Ignoring
unknown parameters: `көлбеу` and `түсі`
- Сандар да растайды
| j_red | j_blue | j_wheat | j_best |
|---|---|---|---|
| 156.7258 | 9.366998 | 5.349284 | 3.633424 |
Қалдықтар туралы ойланайық
- Қалдықтар мен түсіндірме айнымалы мәндері үшін шашыраңқы графигін салайық.
- red моделі үшін орташа қалдықтар \(1.25\) құрайды
students_sleep_reg %>%
ggplot(aes(x = sleep, y = res_red)) +
geom_point() +
geom_hline(yintercept = mean(students_sleep_reg$res_red))
- “көк” моделі үшін орташа қалдықтар \(0.25\) құрайды
students_sleep_reg %>%
ggplot(aes(x = sleep, y = res_blue)) +
geom_point() +
geom_hline(yintercept = mean(students_sleep_reg$res_blue))
- “бидай” моделі үшін орташа қалдықтар \(-0.1\) құрайды
students_sleep_reg %>%
ggplot(aes(x = sleep, y = res_wheat)) +
geom_point() +
geom_hline(yintercept = mean(students_sleep_reg$res_wheat))
- Бірақ біздің ең жақсы моделіміз үшін бұл \(0\) болып шығады
students_sleep_reg %>%
ggplot(aes(x = sleep, y = res_best)) +
geom_point() +
geom_hline(yintercept = mean(students_sleep_reg$res_best))
tibble(R_calculated = my_linear_model$residuals,
Darkhan_calculated = students_sleep_reg$res_best) %>%
head()| R_calculated | Darkhan_calculated |
|---|---|
| -0.1166677 | -0.1166677 |
| 0.0160006 | 0.0160006 |
| -0.1748230 | -0.1748230 |
| 0.0394826 | 0.0394826 |
| -0.1233085 | -0.1233085 |
| 0.3608220 | 0.3608220 |
Ең жақсы модельге тағы бір көз жүгіртейік
- Коэффициенттерді түсінуге тырысайық. Модельге көз жүгіртейік
my_linear_model %>% summary()
Call:
lm(formula = GPA ~ sleep, data = students_sleep)
Residuals:
Min 1Q Median 3Q Max
-0.3754 -0.1227 -0.0279 0.1079 0.4692
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.99373 0.13302 7.471 3.36e-11 ***
sleep 0.30979 0.02155 14.378 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1926 on 98 degrees of freedom
Multiple R-squared: 0.6784, Adjusted R-squared: 0.6751
F-statistic: 206.7 on 1 and 98 DF, p-value: < 2.2e-16- Немесе оқылатын түрде. Мақалалардағы регрессия кестелері әдетте келесідей көрінеді:
| Model | |
|---|---|
| (Intercept) | 0.994 (0.133)*** |
| [0.730, 1.258] | |
| sleep | 0.310 (0.022)*** |
| [0.267, 0.353] | |
| Num.Obs. | 100 |
| R2 | 0.678 |
| R2 Adj. | 0.675 |
| F | 206.738 |
| RMSE | 0.19 |
- Коэффициенттерді алып, оларды сызықтық модельге ауыстырайық.
\[\begin{align*} GPA = 0.99373 + 0.30979*sleep \\ \textrm{if sleep} = 7 \\ GPA = 0.99373 + 0.30979*7 = 3.16226 \\ \textrm{if sleep} = 8 \\ GPA = 0.99373 + 0.30979*8 = 3.47205 \end{align*}\]
- ‘*’ жұлдызшалары коэффициенттің \(0\) екендігі туралы нөлдік гипотезаға қарсы маңыздылық деңгейлерін көрсетеді
- p-мәні неғұрлым кіші болса, біздікіне ұқсас деректерді байқау ықтималдығы соғұрлым аз болады, яғни коэффициент нөлге тең емес.
- Соңында, R-squared, бұл не? Алдымен регрессия теңдеуін қарастырайық.
\[\begin{align*} Y = \beta_0 + \beta_1X + \epsilon (\sim N(\mu_{\epsilon}, \sigma^2_{\epsilon})) \\ Var(Y) = Var(\hat{Y}) + \sigma^2_{\epsilon} \\ Var(Y) = \beta_{1}^{2}Var(\hat{X}) + \sigma^2_{\epsilon} \\ \end{align*}\]
- Мұны былай түсіндіруге болады:
- жауап айнымалысы дисперсиясын түсіндіру айнымалысының дисперсиясы және қалдықтардың дисперсиясы (яғни, шу)- R-квадрат1 болып табылады жауап айнымалысының дисперсиясы мен қалдықтардың дисперсиясы арасындағы қатынас, \((var(Y) - \sigma^2_{\epsilon})/var(Y)\). Басқаша айтқанда, бұл түсіндірілген дисперсия пайызы.
- Сіз қателіктер жібермейтін модельді орната алдыңыз деп елестетіп көріңіз, онда \(\sigma^2_{\epsilon} = 0\) және \(\textrm{R-squared} = 1\), бұл сіз жауап айнымалысындағы барлық дисперсияны түсіндірдіңіз дегенді білдіреді.
| var_y | var_resid | R-squared |
|---|---|---|
| 0.114125 | 0.0367013 | 0.6784119 |
Call:
lm(formula = GPA ~ sleep, data = students_sleep)
Residuals:
Min 1Q Median 3Q Max
-0.3754 -0.1227 -0.0279 0.1079 0.4692
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.99373 0.13302 7.471 3.36e-11 ***
sleep 0.30979 0.02155 14.378 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.1926 on 98 degrees of freedom
Multiple R-squared: 0.6784, Adjusted R-squared: 0.6751
F-statistic: 206.7 on 1 and 98 DF, p-value: < 2.2e-16R-squared туралы аздап көбірек
`geom_smooth()` using formula = 'y ~ x'
- Қалдықтар
Warning: `fortify(<lm>)` was deprecated in ggplot2 4.0.0.
ℹ Please use `broom::augment(<lm>)` instead.
ℹ The deprecated feature was likely used in the ggplot2 package.
Please report the issue at <https://github.com/tidyverse/ggplot2/issues>.
tibble(`Total Variance` = students_sleep_2$GPA_1 %>% var() %>% round(2),
`Explained Variance` = fitr1$fitted.values %>% var() %>% round(2),
`Residual Variance` = fitr1$residuals %>% var() %>% round(2))| Total Variance | Explained Variance | Residual Variance |
|---|---|---|
| 0.1 | 0.06 | 0.04 |
R-squared 2
`geom_smooth()` using formula = 'y ~ x'
- Қалдықтары
| Total Variance | Explained Variance | Residual Variance |
|---|---|---|
| 0.21 | 0.04 | 0.17 |
R-squared 3
`geom_smooth()` using formula = 'y ~ x'
- Қалдықтары
| Total Variance | Explained Variance | Residual Variance |
|---|---|---|
| 1.146 | 0.014 | 1.132 |
Сызықтық регрессия: Болжамдар
Жоғарыдағы мысалдарда біз сызықтық модельге сәйкестендіруге тырысқан деректерге ие болдық. Ол үшін біз бірнеше жасырын болжамдар жасадық. Оларды толығырақ қарастырайық.
Сызықтық: Біз жауап айнымалысы түсіндірме айнымалының сызықтық функциясына тәуелді деп анық болжадық.
Қалыпты қалдықтар: Біз сондай-ақ қалдықтар/қателер қалыпты түрде таралған деп болжадық. Яғни, біз модель қателіктер жібереді деп болжаймыз, бірақ бұл қателіктер симметриялы және әдетте бір-бірін жояды.
Тұрақты дисперсия[^chapter8-3]: Біз сондай-ақ жауап айнымалысының дисперсиясы түсіндірме айнымалының әртүрлі мәндері үшін өзгермейді деп болжадық (яғни, қателіктер әрқашан шамамен бірдей).
Тәуелсіз бақылаулар: Соңында, біз айнымалылардың бақылаулары бір-бірінен тәуелсіз деп болжадық. Бұл болжамның бұзылған мысалы - уақыт қатары. Мысалы, температура мен тәулік уақыты арасындағы байланыс. Егер температура сағат 15:00-де 3˚C болса, онда сағат 16:00-де ол 3˚C-қа жақын болады.
Шарттар орындалмаған кезде не болады
- Айталық, біз стресс, «стресс» және академиялық үлгерім, «өнімділік» арасындағы байланысты қарастырып жатырмыз.
set.seed(123)
temp_df <-
tibble(stress = rnorm(1000, mean = 6, 1),
performance = 13 - ((stress - 5)^2) + rnorm(1000, mean = 0, sd = 2)) %>%
mutate(performance =
10*(performance - min(performance))/(max(performance) - min(performance) + 0.5))
p_ex2 <-
temp_df %>%
ggplot(aes(x = stress, y = performance)) +
geom_point() +
labs(x = "Stress Level", y = "Academic Performance")
p_ex2
Сызықтық еместік
- Деректерді түзу сызықпен сәйкестендіруге тырысайық.
p_ex2 +
geom_smooth(method = "lm")
- Бірдеңе анық дұрыс емес!
Қалдықтар
- Күмәндануларға қарамастан, біз модельді деректерге сәйкестендіреміз делік.
- Нәтижелерді қарастырамыз, және бәрі жақсы болып көрінеді. R-squared $0.37. Уау! Біз дисперсияның үштен бір бөлігінен астамын түсіндірдік!
fit1 <-
lm(performance ~ stress, data = temp_df)
summary(fit1)
Call:
lm(formula = performance ~ stress, data = temp_df)
Residuals:
Min 1Q Median 3Q Max
-4.4317 -0.6240 0.0841 0.6980 2.9457
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 11.15337 0.19052 58.54 <2e-16 ***
stress -0.76947 0.03125 -24.62 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.9794 on 998 degrees of freedom
Multiple R-squared: 0.378, Adjusted R-squared: 0.3773
F-statistic: 606.4 on 1 and 998 DF, p-value: < 2.2e-16- Онша тез емес! Қалдықтарды қарастырайық!
ggplot(fit1, aes(x = .fitted, y = .resid)) +
geom_point() +
geom_hline(yintercept = 0, linetype = "dashed", color = "red")
- Сонымен, қалдықтар бізге не айтады?
- Біздің модель стресс деңгейі төмен немесе жоғары адамдардың өнімділігін жүйелі түрде асыра бағалайды.
- Сондай-ақ, ол орташа стресс деңгейі бар адамдардың өнімділігін де төмен бағалайды.
Қалдықтардың қалыптылығы
- Егер қалдықтардың таралуына қарасақ, оның солға қарай қисайғанын көреміз.
ggplot(data = fit1, aes(x = .resid)) +
geom_blank() +
geom_histogram(aes(y = ..density..), color = "grey50",
bins = 30) +
geom_density(color = "steelblue") +
xlab("Residuals") +
stat_function(fun = dnorm, args = c(mean = 0, sd = sd(fit1$residuals)),
color = "tomato")Warning: The dot-dot notation (`..density..`) was deprecated in ggplot2 3.4.0.
ℹ Please use `after_stat(density)` instead.
- Соңында, Q-Q графигі деп аталатын жаңа нәрсені қолданып көрейік.
Q-Q графигі - екі үлестірімді олардың квантильдерін бір-біріне қарсы қою арқылы салыстырудың графикалық әдісі. Графиктегі нүкте екінші үлестірімнің квантильдерінің біріне сәйкес келеді, бірінші үлестірімнің сол квантильіне қарсы қойылған.
ggplot(data = fit1, aes(sample = .resid)) +
stat_qq() +
stat_qq_line(color = "red") +
labs(x = "Қалыпты квантилдер", y = "Қалдық квантилдер")
Мысалы, бұл графикте қалдық үлестірімді центрі мен дисперсиясы бірдей теориялық қалыпты үлестіріммен салыстырамыз. Қалдықтардың ауырырақ “сол жақ құйрығы” бар екенін көреміз — жеке бақылаулар орташа мәннен 4 стандартты ауытқуларға дейін ауытқиды. Егер біздің қалдықтарымыз қалыпты үлестірілсе, нүктелер түзу сызықта жатады.
Тағы бір мысал
tibble(x = rweibull(100, shape = 1) - 1) %>%
ggplot(aes(x = x)) +
geom_histogram(bins = 30, color = "gray")
set.seed(1212)
df <- tibble(
x = runif(100, -1, 1),
error = rweibull(100, shape = 1) - 1,
y = 1*x + error)
fit <-
df %>%
lm(data = ., y ~ x)
summary(fit)
Call:
lm(formula = y ~ x, data = .)
Residuals:
Min 1Q Median 3Q Max
-1.1496 -0.7185 -0.1893 0.3811 4.2376
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.16465 0.09833 1.674 0.0972 .
x 1.23180 0.17224 7.152 1.56e-10 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.9755 on 98 degrees of freedom
Multiple R-squared: 0.3429, Adjusted R-squared: 0.3362
F-statistic: 51.15 on 1 and 98 DF, p-value: 1.559e-10df %>%
ggplot(aes(x = x, y = y)) +
geom_point() +
geom_smooth(method = "lm")`geom_smooth()` using formula = 'y ~ x'
Жаттығулар
- Жаттығайық. Мен сендерге регрессияларды беремін, ал сендер оларды түсіндіресіңдер.
pisa_kz <-
read_rds("data/pisa_2018.rds") %>%
filter(CNT == "Kazakhstan")
pisa_kz %>% head()| CNT | LANGTEST_COG | ST011Q01TA | ST011Q02TA | ST011Q03TA | ST011Q04TA | ST011Q05TA | ST011Q06TA | ST011Q07TA | ST019AQ01T | ST012Q01TA | ST012Q02TA | ST012Q03TA | ST012Q05NA | ST012Q06NA | ST012Q07NA | ST012Q08NA | ST012Q09NA | ST013Q01TA | ST097Q01TA | ST097Q02TA | ST097Q03TA | ST097Q04TA | ST097Q05TA | ST100Q01TA | ST100Q02TA | ST100Q03TA | ST100Q04TA | ST102Q01TA | ST102Q02TA | ST102Q03TA | ST102Q04TA | ST158Q01HA | ST158Q02HA | ST158Q03HA | ST158Q04HA | ST158Q05HA | ST158Q06HA | ST158Q07HA | ST160Q01IA | ST160Q02IA | ST160Q03IA | ST160Q04IA | ST160Q05IA | ST167Q01IA | ST167Q02IA | ST167Q03IA | ST167Q04IA | ST167Q05IA | ST175Q01IA | FISCED | MISCED | TMINS | ESCS | WEALTH | PV1MATH | PV1READ | PV1SCIE | Gender |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Kazakhstan | Kazakh | Yes | No | Yes | Yes | Yes | Yes | No | Country of test | One | One | None | None | One | None | None | None | 11-25 books | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Every lesson | No | Yes | No | No | No | Yes | No | Disagree | Agree | Disagree | Disagree | Disagree | A few times a year | A few times a year | Several times a month | About once a month | A few times a year | More than 30 minutes to less than 60 minutes a day | ISCED 3A, ISCED 4 | ISCED 3A, ISCED 4 | 1620 | -1.3482 | -1.9957 | 414.520 | 396.672 | 348.232 | Male |
| Kazakhstan | Kazakh | Yes | Yes | Yes | No | No | Yes | Yes | Country of test | One | None | None | None | None | Three or more | Three or more | One | 0-10 books | Never or hardly ever | Never or hardly ever | Never or hardly ever | Never or hardly ever | Never or hardly ever | Every lesson | Most lessons | Most lessons | Most lessons | Most lessons | Most lessons | Most lessons | Most lessons | No | No | Yes | No | Yes | No | No | Disagree | Agree | Disagree | Strongly disagree | Disagree | About once a month | About once a month | About once a month | Several times a month | About once a month | More than 30 minutes to less than 60 minutes a day | ISCED 3A, ISCED 4 | ISCED 3A, ISCED 4 | 2240 | -0.6018 | -0.8302 | 353.766 | 308.628 | 371.455 | Male |
| Kazakhstan | Kazakh | Yes | Yes | Yes | NA | No | NA | NA | Other country | One | NA | NA | NA | NA | NA | NA | NA | 0-10 books | Some lessons | Never or hardly ever | Never or hardly ever | Never or hardly ever | Never or hardly ever | Every lesson | Some lessons | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Yes | Yes | Yes | Yes | Yes | Yes | Yes | Disagree | Agree | Strongly agree | Disagree | Disagree | A few times a year | Several times a month | Several times a week | Several times a month | Never or almost never | More than 30 minutes to less than 60 minutes a day | ISCED 5B | ISCED 3A, ISCED 4 | NA | -1.1250 | NA | 350.357 | 414.048 | 335.833 | Female |
| Kazakhstan | Kazakh | Yes | No | Yes | No | No | Yes | Yes | Country of test | One | One | None | None | None | None | None | One | 11-25 books | Some lessons | Never or hardly ever | Never or hardly ever | Some lessons | Never or hardly ever | Never or hardly ever | Every lesson | Every lesson | Every lesson | Never or hardly ever | Never or hardly ever | Never or hardly ever | Never or hardly ever | Yes | Yes | Yes | No | No | No | No | Strongly disagree | Agree | Agree | Disagree | Disagree | Several times a week | About once a month | Several times a week | Several times a month | Several times a month | 1 to 2 hours a day | ISCED 3A, ISCED 4 | ISCED 3A, ISCED 4 | NA | -0.6677 | -1.8445 | 275.045 | 268.966 | 207.007 | Male |
| Kazakhstan | Kazakh | Yes | Yes | Yes | No | No | Yes | No | Country of test | One | One | None | One | None | Two | None | None | 11-25 books | Some lessons | Never or hardly ever | Never or hardly ever | Never or hardly ever | Never or hardly ever | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Every lesson | Yes | Yes | Yes | Yes | Yes | Yes | Yes | Disagree | Agree | Strongly disagree | Disagree | Disagree | Several times a month | Never or almost never | About once a month | About once a month | Several times a week | 1 to 2 hours a day | ISCED 5B | ISCED 3A, ISCED 4 | 1530 | -1.1238 | -1.9130 | 250.600 | 411.113 | 366.879 | Female |
| Kazakhstan | Kazakh | Yes | No | Yes | Yes | Yes | Yes | No | Country of test | Two | Two | None | None | One | None | None | None | 0-10 books | Never or hardly ever | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | Some lessons | NA | Some lessons | Every lesson | Every lesson | Every lesson | Every lesson | No | No | No | No | Yes | No | No | Agree | Agree | Disagree | Disagree | Agree | A few times a year | Never or almost never | A few times a year | A few times a year | Several times a week | I do not read for enjoyment | ISCED 3A, ISCED 4 | ISCED 3A, ISCED 4 | NA | -1.4306 | -1.2780 | 408.725 | 361.051 | 334.163 | Male |
- Нұсқаулыққа назар аударайық
Rows: 60 Columns: 2
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (2): name, value
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.| name | value |
|---|---|
| CNT | Country code 3-character |
| CNTSCHID | Intl. School ID |
| LANGTEST_COG | Language of Assessment |
| ST011Q01TA | In your home: A desk to study at |
| ST011Q02TA | In your home: A room of your own |
| ST011Q03TA | In your home: A quiet place to study |
| ST011Q04TA | In your home: A computer you can use for school work |
| ST011Q05TA | In your home: Educational software |
| ST011Q06TA | In your home: A link to the Internet |
| ST011Q07TA | In your home: Classic literature (e.g. |
| ST019AQ01T | In what country were you and your parents born? You |
| ST012Q01TA | How many in your home: Televisions |
| ST012Q02TA | How many in your home: Cars |
| ST012Q03TA | How many in your home: Rooms with a bath or shower |
| ST012Q05NA | How many in your home: |
| ST012Q06NA | How many in your home: Computers (desktop computer, portable laptop, or notebook) |
| ST012Q07NA | How many in your home: |
| ST012Q08NA | How many in your home: E-book readers (e.g. |
| ST012Q09NA | How many in your home: Musical instruments (e.g. guitar, piano) |
| ST013Q01TA | How many books are there in your home? |
| ST097Q01TA | How often during |
| ST097Q02TA | How often during |
| ST097Q03TA | How often during |
| ST097Q04TA | How often during |
| ST097Q05TA | How often during |
| ST100Q01TA | How often during |
| ST100Q02TA | How often during |
| ST100Q03TA | How often during |
| ST100Q04TA | How often during |
| ST102Q01TA | How often during |
| ST102Q02TA | How often during |
| ST102Q03TA | How often during |
| ST102Q04TA | How often during |
| ST158Q01HA | Taught at school: How to use keywords when using a search engine such as <Google©>, <Yahoo©>, etc. |
| ST158Q02HA | Taught at school: How to decide whether to trust information from the Internet |
| ST158Q03HA | Taught at school: How to compare different web pages and decide what information is more relevant for your school work |
| ST158Q04HA | Taught at school: To understand the consequences of making information publicly available online on |
| ST158Q05HA | Taught at school: How to use the short description below the links in the list of results of a search |
| ST158Q06HA | Taught at school: How to detect whether the information is subjective or biased |
| ST158Q07HA | Taught at school: How to detect phishing or spam emails |
| ST160Q01IA | How much do you agree or disagree? I read only if I have to. |
| ST160Q02IA | How much do you agree or disagree? Reading is one of my favourite hobbies. |
| ST160Q03IA | How much do you agree or disagree? I like talking about books with other people. |
| ST160Q04IA | How much do you agree or disagree? For me, reading is a waste of time. |
| ST160Q05IA | How much do you agree or disagree? I read only to get information that I need. |
| ST167Q01IA | How often do you read these materials because you want to? Magazines |
| ST167Q02IA | How often do you read these materials because you want to? Comic books |
| ST167Q03IA | How often do you read these materials because you want to? Fiction (novels, narratives, stories) |
| ST167Q04IA | How often do you read these materials because you want to? Non-fiction books (informational, documentary) |
| ST167Q05IA | How often do you read these materials because you want to? Newspapers |
| ST175Q01IA | About how much time do you usually spend reading for enjoyment? |
| FISCED | Father’s Education (ISCED) |
| MISCED | Mother’s Education (ISCED) |
| TMINS | Learning time (minutes per week) - in total |
| ESCS | Index of economic, social and cultural status |
| WEALTH | Family wealth (WLE) |
| PV1MATH | Plausible Value 1 in Mathematics |
| PV1READ | Plausible Value 1 in Reading |
| PV1SCIE | Plausible Value 1 in Science |
| ST004D01T | Student (Standardized) Gender |
Алдымен айнымалыларды қарастырайық
Warning: Removed 7152 rows containing non-finite outside the scale range
(`stat_bin()`).
Call:
lm(formula = PV1MATH ~ TMINS, data = pisa_kz)
Residuals:
Min 1Q Median 3Q Max
-359.30 -61.79 1.03 61.21 336.15
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.016e+02 2.840e+00 141.42 <2e-16 ***
TMINS 3.292e-02 1.627e-03 20.24 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 90.48 on 12353 degrees of freedom
(7152 observations deleted due to missingness)
Multiple R-squared: 0.03209, Adjusted R-squared: 0.03201
F-statistic: 409.5 on 1 and 12353 DF, p-value: < 2.2e-16Түсіндірейік:
Модель: \(\textrm{Math Scores} = (Қиылысу) + \beta_1*TMINS + \textrm{Noise}\)
Коэффициенттерді ауыстырыңыз: \(\textrm{Math Scores} = 401.6 + 0.033*TMINS + \textrm{Noise}\)
Енді аптасына \(2400\) минут (\(40\) сағат) оқитын оқушыдан қандай математикалық балл күтуге болатынын есептеп көрейік.
- \(\textrm{E[Math Score]} = 401.6 + 0.033*2400 + Noise = 480.8\)
Нақты деректерге қаншалықты жақын екенімізді көрейік.
pisa_kz %>%
ggplot(aes(x = TMINS, y = PV1MATH)) +
geom_point(alpha = 0.4) +
geom_smooth(method = "lm") +
geom_vline(xintercept = 2400, color = "tomato", size = 3, alpha = 0.4) +
geom_hline(yintercept = 480.8, color = "green", size = 3, alpha = 0.4)
Ал қалдықтарды қарастырайық.
Жалпы алғанда, өте жақсы.
Мұны Q-Q графигін пайдаланып тексерейік.
Дәл солай. Қалдықтар қалыпты таралған сияқты.
[1] 90.47834Көптік регрессия
Әзірге біз тек бір түсіндірме айнымалымен жұмыс істедік. Алайда, іс жүзінде біз көптеген басқа айнымалыларды қолдануға көбірек қызығушылық танытамыз. Сызықтық регрессия идеясын бірнеше түсіндірме айнымалылар жағдайына кеңейтейік.
\[ Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \dots + \beta_p X_p + Noise \]
PISA деректерімізді алып, ескі модельге тек бір айнымалыны, Gender қосып көрейік.
Формальды түрде, біздің моделіміз келесідей болады:
\[ \textrm{Math Score} = (Intercept) + \beta_1*TMINS + \beta_2*I(Gender = Male) + Noise \]
R тілінде біз осы синтаксисті пайдаланып, деректерге модельді сәйкестендіре аламыз.
Регрессия нәтижелерін қарастырайық.
Call:
lm(formula = PV1MATH ~ TMINS + Gender, data = .)
Residuals:
Min 1Q Median 3Q Max
-361.18 -61.82 0.99 61.32 334.16
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.994e+02 2.978e+00 134.124 <2e-16 ***
TMINS 3.305e-02 1.627e-03 20.309 <2e-16 ***
GenderMale 3.953e+00 1.629e+00 2.427 0.0152 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 90.46 on 12352 degrees of freedom
(7152 observations deleted due to missingness)
Multiple R-squared: 0.03255, Adjusted R-squared: 0.03239
F-statistic: 207.8 on 2 and 12352 DF, p-value: < 2.2e-16Сонымен, алдымен екі коэффициент те маңызды болды. Бұл деректер оқуға жұмсалған уақыт пен студенттің жынысы математикалық ұпай үшін статистикалық тұрғыдан маңызды екенін көрсетеді дегенді білдіреді. Не болып жатқанын жақсы түсіну үшін коэффициент бағаларын алып, оларды формулаға қоса аламыз.
\[ \textrm{Math Scores} = 399.4 + 0.033*TMINS + 3.953*I(Жынысы = Еркек) + Шу \]
Аптасына 40 сағат оқитын қыз бен ұл үшін күтілетін ұпайларды есептеп көрейік.
Қыз бала үшін:
\[ E[\textrm{Math Score}] = 399.4 + 0.033*2400 + 3.953*0 + Шу \апрокс 478.6 \]
Ұл бала үшін:
\[ E[\textrm{Math Score}] = 399.4 + 0.033*2400 + 3.953*1 + Шу \апрокс 482.6 \]
Сонымен, басқа барлық жағдайлар бірдей болғанда, ұл балалар PISA математика емтиханында қыздарға қарағанда шамамен 4 балл жақсы нәтиже көрсететінін анықтадық. Айырмашылық аз, бірақ статистикалық тұрғыдан маңызды.
Көбірек айнымалылар
Енді көбірек айнымалылар қосайық. Содан кейін сіз нәтижелерді түсіндіресіз.
Call:
lm(formula = PV1MATH ~ TMINS + Gender + ESCS, data = .)
Residuals:
Min 1Q Median 3Q Max
-349.20 -59.21 1.25 59.27 332.17
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.163e+02 2.944e+00 141.407 <2e-16 ***
TMINS 2.675e-02 1.591e-03 16.815 <2e-16 ***
GenderMale 4.518e+00 1.578e+00 2.863 0.0042 **
ESCS 2.674e+01 9.401e-01 28.439 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 87.61 on 12342 degrees of freedom
(7161 observations deleted due to missingness)
Multiple R-squared: 0.09187, Adjusted R-squared: 0.09165
F-statistic: 416.2 on 3 and 12342 DF, p-value: < 2.2e-16| Gender | mean(ESCS, na.rm = TRUE) |
|---|---|
| Female | -0.3198707 |
| Male | -0.3441408 |
Сіз модельдеріңізге өзара әрекеттесу әсерлерін қоса аласыз. Өзара әрекеттесу әсері бір айнымалының әсері басқа айнымалының мәніне байланысты болған кезде пайда болады.
Мысалы, отбасылық байлық ұлдар мен қыздарға әртүрлі әсер етуі мүмкін.
Call:
lm(formula = PV1MATH ~ TMINS + Gender * ESCS, data = .)
Residuals:
Min 1Q Median 3Q Max
-348.63 -59.05 1.30 59.37 331.80
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.161e+02 2.954e+00 140.855 < 2e-16 ***
TMINS 2.677e-02 1.591e-03 16.822 < 2e-16 ***
GenderMale 4.961e+00 1.645e+00 3.016 0.00257 **
ESCS 2.584e+01 1.329e+00 19.450 < 2e-16 ***
GenderMale:ESCS 1.774e+00 1.863e+00 0.952 0.34101
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 87.61 on 12341 degrees of freedom
(7161 observations deleted due to missingness)
Multiple R-squared: 0.09194, Adjusted R-squared: 0.09165
F-statistic: 312.4 on 4 and 12341 DF, p-value: < 2.2e-16\[ \textrm{Математикалық ұпай} = 416.1 + 0.0268*TMINS + 4.96*Жынысы + 25.84*ESCS + 0*ЖынысыЕр:ESCS \]
Бай отбасынан шыққан қыз: TMINS = 2400, ESCS = 1
\[ \textrm{Математикалық ұпай} = 416.1 + 0.0268*2400 + 25.84*1 = 506.26 \]
Орташа отбасынан шыққан қыз: TMINS = 2400, ESCS = 0
\[ \textrm{Математикалық ұпай} = 416.1 + 0.0268*2400 + 25.84*0 = 481 \]
Бай отбасынан шыққан бала:
\[ \textrm{Математикалық ұпай} = 416.1 + 0.0268*2400 + 4.96*1 + 25.84*1 = 511.26 \]
Орта тап отбасынан шыққан бала:
\[ \textrm{Математикалық ұпай} = 416.1 + 0.0268*2400 + 4.96*1 + 25.84*0 = 485.38 \]
Өзара әрекеттесу әсерлері
- Мысал ретінде, Қазақстан мен Жапония үшін PISA нәтижелерін салыстырайық
| CNT | Mean |
|---|---|
| Japan | 528.2104 |
| Kazakhstan | 439.9113 |
- Бір түсіндірме айнымалысы бар қарапайым модельді қарастырайық (ел)
| (1) | |
|---|---|
| (Intercept) | 528.210 |
| (1.174) | |
| CNTKazakhstan | -88.299 |
| (1.345) | |
| Num.Obs. | 25616 |
| R2 | 0.144 |
| R2 Adj. | 0.144 |
| AIC | 304222.3 |
| BIC | 304246.7 |
| Log.Lik. | -152108.137 |
| F | 4308.341 |
| RMSE | 91.75 |
Мұнда біз Қазақстандағы оқушылардың математикалық тесттерден Жапонияға қарағанда орташа есеппен 88 балл нашар жинағанын көреміз.
Модельді күрделендіріп, BAYLIQ айнымалысын (отбасылық байлық) қосайық.
| (1) | |
|---|---|
| (Intercept) | 534.606 |
| (1.198) | |
| CNTKazakhstan | -79.115 |
| (1.397) | |
| WEALTH | 14.850 |
| (0.697) | |
| Num.Obs. | 25487 |
| R2 | 0.159 |
| R2 Adj. | 0.159 |
| AIC | 302112.7 |
| BIC | 302145.3 |
| Log.Lik. | -151052.338 |
| RMSE | 90.71 |
Қазақстан үшін коэффициент өзгерді! Неліктен?
Себебі CNT айнымалысы Wealth айнымалысымен корреляцияланған. Қарап көріңіз.
Warning: Removed 129 rows containing non-finite outside the scale range
(`stat_bin()`).
| (1) | |
|---|---|
| (Intercept) | 490.772 |
| (2.812) | |
| CNTKazakhstan | -71.391 |
| (1.489) | |
| WEALTH | 12.564 |
| (0.837) | |
| THOURS | 1.776 |
| (0.089) | |
| Num.Obs. | 17872 |
| R2 | 0.174 |
| R2 Adj. | 0.174 |
| AIC | 210561.1 |
| BIC | 210600.1 |
| Log.Lik. | -105275.571 |
| RMSE | 87.50 |
Қазақстан үшін коэффициент қайтадан төмендеді (абсолютті мәнде)!
Бұл
THOURSжәнеCNTайнымалыларының да корреляцияланғанын білдіреді.
Warning: Removed 7715 rows containing non-finite outside the scale range
(`stat_bin()`).
- Енді
CNTелі менWEALTHбайлығы арасындағы өзара әрекеттесу әсері бар модель құруға тырысайық.
| (1) | |
|---|---|
| (Intercept) | 486.059 |
| (2.834) | |
| CNTKazakhstan | -59.828 |
| (1.810) | |
| WEALTH | -2.063 |
| (1.554) | |
| THOURS | 1.733 |
| (0.089) | |
| CNTKazakhstan × WEALTH | 20.525 |
| (1.839) | |
| Num.Obs. | 17872 |
| R2 | 0.180 |
| R2 Adj. | 0.180 |
| AIC | 210439.0 |
| BIC | 210485.8 |
| Log.Lik. | -105213.510 |
| RMSE | 87.19 |
- Мұны қалай түсінеміз? Формуланы қарастырайық:
\[\begin{align*} \textrm{Math Score} &= \beta_0 + \beta_1 \times I(CNT = Қазақстан) + \beta_2 \times BAYLIQ + \beta_3 \times THOURS + \beta_4 \times (BAYLIQ \times Қазақстан) \\ \textrm{Math Score} &= 486 − 59.8 \times I(CNT = Қазақстан) − 0 \times BAYLIQ + 1.733 \times THOURS + 20.525 \times (BAYLIQ \times Қазақстан) \end{align*}\]
- Енді Қазақстан мен Жапониядағы бай және кедей отбасылардан шыққан оқушылардың математикадан қалай нәтиже көрсететінін түсінуге тырысайық.
Кедей отбасынан шыққан еңбекқор жапон оқушысы: САҒАТ = 40, БАЙЛЫҚ = 0
\[ Математикадан ұпайлар = 486 - 0 - 0 + 1.733 \есеп 40 + 20.5 \есеп 0 \шамамен 555 \]
Кедей отбасынан шыққан еңбекқор қазақ жігіті: САҒАТ = 40, БАЙЛЫҚ = 0
\[ Математикадан ұпайлар = 486 - 59.8 - 0 + 1.733 \есеп 40 + 20.5 \есеп 0 \шамамен 495.5 \]
Бай отбасынан шыққан еңбекқор жапон жігіті: САҒАТ = 40, БАЙЛЫҚ = 3
\[ Математикадан ұпайлар = 486 - 0 - 0 + 1.733 \есеп 40 + 20.5 \есеп 0 \шамамен 555 \]
Еңбекқор Бай отбасынан шыққан қазақ: САҒАТ = 40, БАЙЛЫҚ = 3
\[ Математикадан ұпайлар = 486 - 59.8 - 0 + 1.733 \times 40 + 20.5 \times 3 \approx 557 \]
- Графиктерді салуға болады. Алдымен синтетикалық деректерді құрайық.
| CNT | THOURS | WEALTH |
|---|---|---|
| Kazakhstan | 30 | -3 |
| Kazakhstan | 30 | -2 |
| Kazakhstan | 30 | -1 |
| Kazakhstan | 30 | 0 |
| Kazakhstan | 30 | 1 |
| Kazakhstan | 30 | 2 |
| Kazakhstan | 30 | 3 |
| Japan | 30 | -3 |
| Japan | 30 | -2 |
| Japan | 30 | -1 |
| Japan | 30 | 0 |
| Japan | 30 | 1 |
| Japan | 30 | 2 |
| Japan | 30 | 3 |
- Содан кейін оны саламыз.
Ең жақсы модельді қалай таңдауға болады?
Жоғарыдағы мысалда бізде бірнеше модель болды. Қайсысы жақсырақ екенін қалай анықтауға болады?
Бір жолы - модельдердің түсіндірілген дисперсиясын (R-squared) салыстыру.
| model1 | model2 | model3 | model4 | |
|---|---|---|---|---|
| (Intercept) | 528.210 (1.174)*** | 534.606 (1.198)*** | 490.772 (2.812)*** | 486.059 (2.834)*** |
| CNTKazakhstan | -88.299 (1.345)*** | -79.115 (1.397)*** | -71.391 (1.489)*** | -59.828 (1.810)*** |
| WEALTH | 14.850 (0.697)*** | 12.564 (0.837)*** | -2.063 (1.554) | |
| THOURS | 1.776 (0.089)*** | 1.733 (0.089)*** | ||
| CNTKazakhstan × WEALTH | 20.525 (1.839)*** | |||
| Num.Obs. | 25616 | 25487 | 17872 | 17872 |
| R2 | 0.144 | 0.159 | 0.174 | 0.180 |
| R2 Adj. | 0.144 | 0.159 | 0.174 | 0.180 |
| F | 4308.341 | |||
| RMSE | 91.75 | 90.71 | 87.50 | 87.19 |
Бұл айырмашылық қаншалықты үлкен? Мысалы, 4-ші модель 3-ші модельден айтарлықтай жақсы деп айта аламыз ба?
Түсіндірілген дисперсия (яғни, R-squared) да кездейсоқ айнымалы болып шығады. Және тек кездейсоқ емес, белгілі үлестірімі бар айнымалы. Бұл түсіндірілген дисперсия айтарлықтай өзгерді ме деген сұраққа жауап беру үшін ANOVA сынағын пайдалана алатынымызды білдіреді.
| term | df.residual | rss | df | sumsq | statistic | p.value |
|---|---|---|---|---|---|---|
| PV1MATH ~ CNT * WEALTH + THOURS | 17867 | 135872904 | NA | NA | NA | NA |
| PV1MATH ~ CNT + WEALTH + THOURS | 17868 | 136819845 | -1 | -946940.8 | 124.5207 | 0 |
| Res.Df | RSS | Df | Sum of Sq | F | Pr(>F) |
|---|---|---|---|---|---|
| 17867 | 135872904 | NA | NA | NA | NA |
| 17868 | 136819845 | -1 | -946940.8 | 124.5207 | 0 |