Berechnet man eine Posteriori-Verteilung mit stan_glm(), so kann man entweder die schwach informativen Prioriwerte der Standardeinstellung verwenden, oder selber Prioriwerte definieren.
Grundsätzlich hat ein Regressionsmodell die folgenden Parameter:
einen Parameter für den Intercept, \(\beta_0\)
pro UV ein weiterer Parameter, \(\beta_1, \beta_2, \ldots\)
für sigma (\(\sigma\)) noch ein zusätzlicher Parameter
Zu beachten ist aber, dass bei einer nominalen Variablen mit zwei Stufen nur ein Regressionsgewicht (\(\beta_1\)) berechnet wird. Allgemein gilt bei nominalen also, dass bei \(k\) Stufen nur \(k-1\) Regressionsgewichte berechnet werden.
Im vorliegenden Fall hat die Variable cut 5 Stufen, also werden 4 Regressiongewiche berechnet.
In Summe werden also 6 Parameter berechnet.
Berechnet man das Modell, so kann man sich auch Infos über die Prioris ausgeben lassen:
Priors for model 'm1'
------
Intercept (after predictors centered)
~ normal(location = 3000, scale = 500)
Coefficients
~ normal(location = [100,100,100,...], scale = [100,100,100,...])
Auxiliary (sigma)
~ exponential(rate = 1)
------
See help('prior_summary.stanreg') for more details
Wie man sieht, wird für die Streuung im Standard eine Exponentialverteilung verwendet von stan_glm(). Gibt man also nicht an - wie im Beispiel m1 oben, so wird stan_glm() für die Streuung, d.h. prior_aux eine Exponentialverteilung verwenden. Zu beachten ist, dass stan_glm() ein automatische Skalierung vornimmt.
Priors for model 'm1'
------
Intercept (after predictors centered)
~ normal(location = 3000, scale = 500)
Coefficients
~ normal(location = [100,100,100,...], scale = [100,100,100,...])
Auxiliary (sigma)
~ exponential(rate = 1)
------
See help('prior_summary.stanreg') for more details
Zu beachten ist beim selber definieren der Prioris, dass dann keine Auto-Skalierung von stan_glm() vorgenommen wird, es sei denn, man weist es explizit an:
m3 <-stan_glm(price ~ cut, data = diamonds, prior =normal(location =c(100, 100, 100, 100),scale =c(100, 100, 100, 100),autoscale =TRUE),prior_intercept =normal(3000, 500, autoscale =TRUE),prior_aux =exponential(1, autoscale =TRUE),chain =1, # nur 1 mal Stichproben ziehen, um Zeit zu sparenrefresh =0)prior_summary(m3)
Grundsätzlich ist es nützlich für die numerische Stabilität, dass die Zahlen (hier die Parameterwerte) etwa die gleiche Größenordnung haben, am besten um die 0-1 herum. Daher bietet sich oft eine z-Standardisierung an.
Unabhängig von der der Art der Parameter ist die Anzahl immer gleich.
Die Anzahl der geschätzten Parameter werden im Modell-Summary unter Estimates gezeigt:
Model Info:
function: stan_glm
family: gaussian [identity]
formula: price ~ cut
algorithm: sampling
sample: 4000 (posterior sample size)
priors: see help('prior_summary')
observations: 53940
predictors: 5
Estimates:
mean sd 10% 50% 90%
(Intercept) 4027.0 22.1 3998.3 4026.8 4055.8
cut.L -248.8 52.7 -316.9 -249.5 -180.9
cut.Q -283.6 46.9 -342.3 -283.4 -223.2
cut.C -580.5 41.8 -633.6 -580.1 -527.2
cut^4 -266.4 36.8 -314.6 -266.5 -218.6
sigma 3830.6 11.1 3816.4 3830.4 3845.4
Fit Diagnostics:
mean sd 10% 50% 90%
mean_PPD 3931.9 23.1 3902.3 3931.8 3962.0
The mean_ppd is the sample average posterior predictive distribution of the outcome variable (for details see help('summary.stanreg')).
MCMC diagnostics
mcse Rhat n_eff
(Intercept) 0.4 1.0 3025
cut.L 1.1 1.0 2393
cut.Q 1.0 1.0 2135
cut.C 0.8 1.0 2550
cut^4 0.7 1.0 3132
sigma 0.1 1.0 7241
mean_PPD 0.4 1.0 3673
log-posterior 0.0 1.0 1883
For each parameter, mcse is Monte Carlo standard error, n_eff is a crude measure of effective sample size, and Rhat is the potential scale reduction factor on split chains (at convergence Rhat=1).
Das sind:
1 Intercept (Achsenabschnitt) - prior_intercept
4 Gruppen (zusätzlich zur Referenzgruppe, die mit dem Achsenabschnitt dargestellt ist) - prior_normal
1 Sigma (Ungewissheit “innerhalb des Modells”) - prior_aux
Categories:
~
Source Code
---exname: stan-glm-parameterzahlextype: numexsolution: r solexshuffle: noextol: 0expoints: 1date: '2022-12-13'slug: stan_glm_parameterzahltitle: stan_glm_parameterzahlcategories:- bayes- regression- parameters---```{r global-knitr-options, include=FALSE}knitr::opts_chunk$set(fig.pos = 'H', fig.asp = 0.618, fig.width = 4, fig.cap = "", fig.path = "", echo = FALSE, message = FALSE, fig.show = "hold")```# ExerciseBerechnet man eine Posteriori-Verteilung mit `stan_glm()`, so kann man entweder die schwach informativen Prioriwerte der Standardeinstellung verwenden, oder selber Prioriwerte definieren.Betrachten Sie dazu dieses Modell:```stan_glm(price ~ cut, data = diamonds, prior = normal(location = c(100, 100, 100, 100), scale = c(100, 100, 100, 100)), prior_aux = exponential(1), prior_intercept = normal(3000, 500))```Wie viele Parameter gibt es in diesem Modell?Hinweise:- Geben Sie nur eine (ganze) Zahl ein.</br></br></br></br></br></br></br></br></br></br># SolutionGrundsätzlich hat ein Regressionsmodell die folgenden Parameter:1. einen Parameter für den Intercept, $\beta_0$2. pro UV ein weiterer Parameter, $\beta_1, \beta_2, \ldots$3. für sigma ($\sigma$) noch ein zusätzlicher ParameterZu beachten ist aber, dass bei einer *nominalen* Variablen mit zwei Stufen nur *ein* Regressionsgewicht ($\beta_1$) berechnet wird. Allgemein gilt bei nominalen also, dass bei $k$ Stufen nur $k-1$ Regressionsgewichte berechnet werden.Im vorliegenden Fall hat die Variable `cut` 5 Stufen, also werden 4 Regressiongewiche berechnet.In Summe werden also *6* Parameter berechnet.```{r libs, message = FALSE}library(tidyverse)library(rstanarm)```Berechnet man das Modell, so kann man sich auch Infos über die Prioris ausgeben lassen:```{r echo = TRUE}m1 <- stan_glm(price ~ cut, data = diamonds, prior = normal(location = c(100, 100, 100, 100), scale = c(100, 100, 100, 100)), prior_intercept = normal(3000, 500), prior_aux = exponential(1), refresh = 0)prior_summary(m1)```Wie man sieht, wird für die Streuung im Standard eine Exponentialverteilung verwendet von `stan_glm()`. Gibt man also nicht an - wie im Beispiel `m1` oben, so wird `stan_glm()` für die Streuung, d.h. `prior_aux` eine Exponentialverteilung verwenden.Zu beachten ist, dass `stan_glm()` ein automatische Skalierung vornimmt.S. [hier](http://mc-stan.org/rstanarm/articles/priors.html#auxiliary-parameters) für weitere Erläuterung.Möchte man den Prior für die Streuung direkt ansprechen, so kann man das so formulieren:```{r echo = TRUE}m2 <- stan_glm(price ~ cut, data = diamonds, prior = normal(location = c(100, 100, 100, 100), scale = c(100, 100, 100, 100)), prior_intercept = normal(3000, 500), prior_aux = exponential(1), refresh = 0)prior_summary(m1)```Zu beachten ist beim selber definieren der Prioris,dass dann keine Auto-Skalierung von `stan_glm()` vorgenommen wird,es sei denn, man weist es explizit an:```{r echo = TRUE}m3 <- stan_glm(price ~ cut, data = diamonds, prior = normal(location = c(100, 100, 100, 100), scale = c(100, 100, 100, 100), autoscale = TRUE), prior_intercept = normal(3000, 500, autoscale = TRUE), prior_aux = exponential(1, autoscale = TRUE), chain = 1, # nur 1 mal Stichproben ziehen, um Zeit zu sparen refresh = 0)prior_summary(m3)```Grundsätzlich ist es nützlich für die numerische Stabilität,dass die Zahlen (hier die Parameterwerte) etwa die gleiche Größenordnung haben,am besten um die 0-1 herum.Daher bietet sich oft eine z-Standardisierung an.Unabhängig von der der Art der Parameter ist die Anzahl immer gleich.Die Anzahl der *geschätzten Parameter* werden im Modell-Summary unter `Estimates` gezeigt:```{r}summary(m1)```Das sind:- 1 Intercept (Achsenabschnitt) - `prior_intercept`- 4 Gruppen (zusätzlich zur Referenzgruppe, die mit dem Achsenabschnitt dargestellt ist) - `prior_normal`- 1 Sigma (Ungewissheit "innerhalb des Modells") - `prior_aux`---Categories: ~