-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathc4.Rmd
143 lines (105 loc) · 12.1 KB
/
c4.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
---
title: "Rethinking. Поглавје 4. Белешки"
author: "discindo"
date: "February 2, 2020"
output:
html_document:
df_print: paged
toc: yes
toc_float: yes
toc_depth: 2
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, eval=TRUE)
```
```{r}
suppressPackageStartupMessages(library(rethinking))
suppressPackageStartupMessages(library(ggplot2))
```
# За
Целта на ова поглавје е очигледна. Да научиме како да правиме линеарна регресија од нула, во Баесова смисла. Ова е супер интересно, бидејќи нѐ охрабрува да мислиме за претпоставките директно и да го користиме ова сфаќање за да го креираме моделот од почеток. Речиси секој од нас има направено регресија претходно (и Excel се важи), но ретко кој од нас го има тоа направено преку рачно местење на релевантните параметри или дистрибуции. Добро разбирање на оваа процедура е неопходно за покомплицираните модели (машинско учење, вештачка интелегенција, итн.).
# Што направивме до сега
Материјалот во првите три поглавја од книгата е резимиран во претходниот документ (`c3.Rmd/c3.html`). Претходните поглавја го втемелуваат Баесовиот начин на размислување (кој преовладува во статистичките школи во светот) и воведуваат основни алатки за креирање, сумаризирање, и симулирање на примерици од постериорната дистрибуција. Иако може да се почне со учење тука (со линеарна регресија) овој документ претпоставува дека претходниот материјал е совладан.
# Јазик за опишување модели
Почнуваме со практичниот дел од поглавје 4, имено како да опишеме, креираме, и извришиме еден линеарен модел. Претходниот текст е супер интересен, како од историска и концептуална смисла, така и од аспект на изворите на нормалната (Гаусова) дистрибуција (собирање или множење на мали флуктуации). Читањето на тие страници ако ништо друго може да ве научи како да симулирате нормални дистрибуции без готови функции (пр. `rnorm` во `R`). Но исто така посочува на разликата помегу предвидувања и каузални (причинско-последични) врски: геоцентричниот модел е еквивалентен на машинско учење бидејќи вистинските причини и последици не се од примарен интерес, она што не интересира е предвидување на вредност што е можно подобро (координати на небото или кликање на вебсајт). Додека прашањето „зошто“ е второстепено.
Да се вратиме на терминологија.
1) Прво имаме некои мерења (опсервации) што сакаме да ги разбереме или предвидиме. Тоа е _резултантната_ променлива или променливи (избегнувам да го користам терминот „зависна“)
2) За секоја од овие _резултантни_ променливи дефинираме дистрибуција за веројатност, која ги дава можностите за индивидуалните мерења. За линеарна регресија оваа дистрибуција е секогаш _нормална_.
3) Второ, наоѓаме мерења (опсервации) кои сакаме да ги користиме за да ја разбереме или предвидиме резултантните променливи. Овие ги нарекуваме _предвидувачки_ или _објаснувачки_ променливи (пак избегнувам „независни“).
4) Ја поврзуваме точната форма на дистрибуцијата на веројатноста со променливите за предвидување. Точната форма подразбира прецизната локација, варијанца, и други компоненти, ако ги има. Со овој акт, приморани сме да ги дефинираме и именуваме сите параметри на моделот.
5) На крај, одбираме приорни дистрибуции за сите параметри на моделот. Приорните дистрибуции ги дефинираат првичните информации на моделот, пред да бидат земени во предвид видените податоци.
Тогаш го сумираме моделот на начин сличен на овој:
$$
outcome_i ∼ Normal(μ_i, σ)\\
μ_i = β × predictor_i\\
β ∼ Normal(0, 10)\\
σ ∼ HalfCauchy(0, 1)\\
$$
Она што следи е изкажување на овој математички јазик во форма на код, без да правиме математички манипулации. Во основа, овие модели дефинираат начини на кои може да дојдеме до обзервациите на _резултантните_ променливи имајќи ги во предвид _предвидувачките_ променливи.
# Линеарна регресија рачно!
Тука ќе ја дефинираме и извршиме нашата прва линеарна регресија со решеткаво приближување. Не интересира висина на популација на некоја етничка група на луѓе од Африка. Она што се обидуваме да го процениме е средната вредност ($\mu$) и стандардната девијација ($\sigma$) на популацијата користејќи нормална (Гаусова) дистрибуција. Со нашиот нов статистички јазик овој модел можеме лесно да го опишеме:
$$h_i ∼ Normal(\mu, \sigma)$$
```{r}
data(Howell1)
d <- Howell1
# d$height
# за сега не интересира само висината на возрасни
d2 <- d[ d$age >= 18, ]
# tidy equivalent
# d2 <- d %>% filter(age >= 18)
```
Нормално во Баесова средина ни требаат приорни дистрибуции за сите параметри, па нашиот модел го доопишуваме со:
$$
h_i ∼ Normal(\mu, \sigma) \\
\mu ∼ Normal(178, 20) \\
\sigma ∼ Uniform(0, 50) \\
$$
Што значат овие приори? Дека средната вредност за висината на популацијата од интерес ја влечеме од нормална дистрибуција со средна вредност 178 cm и стандардна девијација од 20 cm. Дека стандардната девијација на висината на популацијата што ја моделираме ја влечеме од униформна дистрибуција со граници 0 и 50. Ајде да га визуелизираме овие приорни дистрибуции.
```{r}
par(mfrow=c(1,2))
curve(dnorm(x = x, mean = 178, sd = 20), from = 100, to = 250)
curve(dunif(x = x, min=0, max = 50), from = -10, to = 60)
```
Значи очекуваме средната висина да биде некаде околу 170-180 cm, но дозволуваме поприлично голем распон така да моделот ќе евалуира бајѓи можности. За стандардната девијација немаме посебни ограничувања (помалку интуиција за дисперзијата на висина?) освен тоа дека овој параметар мора да биде позитивен. ОК, сега може да влечеме вредности од приорните дистрибуции за да видиме какви вредности поточно ќе евалуира моделот
```{r}
sample_mu <- rnorm( 1e4 , 178 , 20 ) # можни вредности за просек
sample_sigma <- runif( 1e4 , 0 , 50 ) # можни врености за девијација
prior_h <- rnorm( 1e4 , sample_mu , sample_sigma ) # можни вредности за висина
par(mfrow=c(1,1))
dens( prior_h )
```
Поприлично рамен приор. Дури и вредности од 100 и 250 cm имаат некаква осетлива приорна веројатност. Може ќе треба да пробаме малку потесна дистрибуција, но за сега одиме со ова.
Сега правиме решеткаво приближување за нашата линеарна регресија. Треба да процениме два параметри, средната вредност и стандардната девијација. Така да кодот е малку повозбудлив.
```{r}
# решетка за просек
mu.list <- seq( from=140, to=160 , length.out=200 )
# решетка за дисперзија
sigma.list <- seq( from=4 , to=9 , length.out=200 )
# формираме табела со сите парови на можности од решетката
post <- expand.grid( mu=mu.list , sigma=sigma.list )
# head(post)
# nrow(post)
# калкулираме лог веројатност (LL = log likelihood)
# post$LL <- sapply(1:nrow(post) , function(i)
# sum(dnorm(d2$height , mean = post$mu[i] , sd = post$sigma[i] , log = TRUE)))
# горниот код со `sapply` е само по-компактен начин да го напишеме овој for циклус
for (i in 1:nrow(post)) {
# за секоја видена вредност за висина
# пресметај ја веројатноста на сите парови на параметри во решетката
# врати веројатност на логаритамска скала
post$LL[i] <- sum(dnorm(d2$height , mean = post$mu[i] , sd = post$sigma[i] , log = TRUE))
}
# множиме веројатност со приор. за логаритми множење е собирање
post$prod <- post$LL +
dnorm( post$mu , 178 , 20 , log = TRUE ) +
dunif( post$sigma , 0 , 50 , TRUE )
# тука прави трик за да избегне грешка со заокружување мали броеви (веројатност)
# но овие вредности се (пропорционални на) нашата постериорна дистрибуција
post$prob <- exp( post$prod - max(post$prod) )
```
Конечно ја визуелизираме постериорната дистрибуција.
```{r}
contour_xyz( post$mu , post$sigma , post$prob )
image_xyz( post$mu , post$sigma , post$prob )
ggplot(post, aes(x=mu, y=prob)) + geom_line()
```