-
Notifications
You must be signed in to change notification settings - Fork 36
/
counts.Rmd
151 lines (115 loc) · 5.33 KB
/
counts.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
144
145
146
147
148
149
150
151
# Unbounded Count Models
Unbounded count models are models in which the response is a natural number with no upper bound,
$$
y_i = 0, 1, \dots, \infty .
$$
The two distributions most commonly used to model this are
- Poisson
- Negative Binomial
## Poisson
The [Poisson distribution](https://en.wikipedia.org/wiki/Poisson_distribution):
$$
y_i \sim \dpois(\lambda_i) ,
$$
In Stan, the Poisson distribution has two implementations:
- `r stanfunc("poisson_lpdf")`
- `r stanfunc("poisson_log_lpdf")`
The Poisson with a log link. This is for numeric stability.
In **rstanarm** use `r rdoc("rstanarm", "stan_glm")` for a Poisson GLM.
- **rstanarm** vignette on [count models](https://cran.r-project.org/web/packages/rstanarm/vignettes/count.html).
## Negative Binomial
The [Negative Binomial](https://en.wikipedia.org/wiki/Negative_binomial_distribution) model is also used for unbounded count data,
$$
Y = 0, 1, \dots, \infty
$$
The Poisson distribution has the restriction that the mean is equal to the variance, $\E(X) = \Var(X) = \lambda$.
The Negative Binomial distribution has an additional parameter that allows the variance to vary (though it is always larger than the mean).
The outcome is modeled as a negative binomial distribution,
$$
y_i \sim \dBinom(\alpha_i, \beta)
$$
with shape $\alpha \in \R^{+}$ and inverse scale $\beta \in \R^{+}$, and $\E(y) = \alpha_i / \beta$ and $\Var(Y) = \frac{\alpha_i}{\beta^2}(\beta + 1)$.
Then the mean can be modeled and transformed to the
$$
\begin{aligned}[t]
\mu_i &= \log( \vec{x}_i \vec{\gamma} ) \\
\alpha_i &= \mu_i / \beta
\end{aligned}
$$
**Important** The negative binomial distribution has many different parameterizations.
An alternative parameterization of the negative binomial uses the mean and a over-dispersion parameter.
$$
y_i \sim \dnbinalt(\mu_i, \phi)
$$
with location parameter $\mu \in \R^{+}$ and over-dispersion parameter $\phi \in \R^{+}$, and $\E(y) = \mu_i$ and $\Var(Y) = \mu_i + \frac{\mu_i^2}{\phi}$.
Then the mean can be modeled and transformed to the
$$
\begin{aligned}[t]
\mu_i &= \log( \vec{x}_i \vec{\gamma} ) \\
\end{aligned}
$$
### Stan
In Stan, there are multiple parameterizations of the negative binomial distribution:
- `neg_binomial_lpdf(y | alpha, beta)`with shape parameter `alpha` and inverse scale parameter `beta`.
- `neg_binomial_2_lpdf(y | mu, phi)` with mean `mu` and over-dispersion parameter `phi`.
- `neg_binomial_2_log_lpdf(y | eta, phi)` with log-mean `eta` and over-dispersion parameter `phi`
Also, `rstanarm` supports Poisson and [negative binomial models](https://cran.r-project.org/web/packages/rstanarm/vignettes/count.html).
### Link functions
In many applications, $\lambda$, is modeled as some function of covariates or other parameters.
Since $\lambda_i$ must be positive, the most common link function is the log,
$$
\log(\lambda_i) = \exp(\vec{x}_i' \vec{\beta})
$$
which has the inverse,
$$
\lambda_i = \log(\vec{x}_i \vec{\beta})
$$
### Example: Bilateral Sanctions
A regression model of bilateral sanctions for the period 1939 to 1983 [@Martin1992a].
The outcome variable is the number of countries imposing sanctions.
```{r results='hide'}
mod_poisson1 <- stan_model("stan/poisson1.stan")
```
```{r}
data("sanction", package = "Zelig")
```
## Negative Binomial
The [Negative Binomial](https://en.wikipedia.org/wiki/Negative_binomial_distribution) model is also used for unbounded count data,
$$
Y = 0, 1, \dots, \infty
$$
The Poisson distribution has the restriction that the mean is equal to the variance, $\E(X) = \Var(X) = \lambda$.
The Negative Binomial distribution has an additional parameter that allows the variance to vary (though it is always larger than the mean).
The outcome is modeled as a negative binomial distribution,
$$
y_i \sim \dnegbinom(\alpha_i, \beta)
$$
with shape $\alpha \in \R^{+}$ and inverse scale $\beta \in \R^{+}$, and $\E(y) = \alpha_i / \beta$ and $\Var(Y) = \frac{\alpha_i}{\beta^2}(\beta + 1)$.
Then the mean can be modeled and transformed to the
$$
\begin{aligned}[t]
\mu_i &= \log( \vec{x}_i \vec{\gamma} ) \\
\alpha_i &= \mu_i / \beta
\end{aligned}
$$
**Important** The negative binomial distribution has many different parameterizations.
An alternative parameterization of the negative binomial uses the mean and a over-dispersion parameter.
$$
y_i \sim \dnbinomalt(\mu_i, \phi)
$$
with location parameter $\mu \in \R^{+}$ and over-dispersion parameter $\phi \in \R^{+}$, and $\E(y) = \mu_i$ and $\Var(Y) = \mu_i + \frac{\mu_i^2}{\phi}$.
Then the mean can be modeled and transformed to the
$$
\begin{aligned}[t]
\mu_i &= \log( \vec{x}_i \vec{\gamma} ) \\
\end{aligned}
$$
In Stan, there are multiple parameterizations of the
- `neg_binomial_lpdf(y | alpha, beta)`with shape parameter `alpha` and inverse scale parameter `beta`.
- `neg_binomial_2_lpdf(y | mu, phi)` with mean `mu` and over-dispersion parameter `phi`.
- `neg_binomial_2_log_lpdf(y | eta, phi)` with log-mean `eta` and over-dispersion parameter `phi`
Also, `rstanarm` supports Poisson and [negative binomial models](https://cran.r-project.org/web/packages/rstanarm/vignettes/count.html).
### Example: Economic Sanctions II {ex-econ-sanctions-2}
Continuing the [economic sanctions example](#ex-econ-sanctions-2) of @Martin1992a.
### References
For general references on count models see @GelmanHill2007a [p. 109-116], @McElreath2016a [Ch 10], @Fox2016a [Ch. 14], and @BDA3 [Ch. 16].