含類別預測變數的貝葉斯迴歸

使用 RJAGS 的貝氏建模

Alicia Johnson

Associate Professor, Macalester College

第 4 章目標

在貝葉斯模型中納入「類別型」預測變數
設計「多變量」貝葉斯迴歸模型
將常態迴歸的方法擴展到廣義線性模型：Poisson 迴歸

步道流量

目標：
探討麻州某自行車步道的日流量。

¹ Photo courtesy commons.wikimedia.org

建模流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數） $ \; $
$ \;$

模型
$Y_i \sim N(m_i, s^2)$

以平日狀態建模流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數）
$X_i$ = 平日為 1，週末為 0

模型
$Y_i \sim N(m_i, s^2)$

以平日狀態建模流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數）
$X_i$ = 平日為 1，週末為 0

模型
$Y_i \sim N(m_i, s^2)$

以平日狀態建模流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數）
$X_i$ = 平日為 1，週末為 0

模型
$Y_i \sim N(m_i, s^2)$
$m_i = a + bX_i$

以平日狀態建模流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數）
$X_i$ = 平日為 1，週末為 0

模型
$Y_i \sim N(m_i, s^2)$
$m_i = a + bX_i$

$a$ = 典型週末流量

$Y_i$ = 第 $i$ 天的步道流量（使用者數）
$X_i$ = 平日為 1，週末為 0

模型
$Y_i \sim N(m_i, s^2)$
$m_i = a + bX_i$

$a$ = 典型週末流量
$a + b$ = 典型平日流量

$b$ = 平日與週末典型流量的差異
$s$ = 殘差標準差

$a$ 與 $b$ 的先驗

典型「週末」流量最可能約為每天 400 位使用者，但也可能低至 100 或高至 700。

我們不確定平日與週末的流量高低關係。可能更多，也可能更少。

$s$ 的先驗

$\;$

每日流量的標準差（不論平日或週末）在 0 到 200 位使用者之間皆同等可能。

依平日狀態的貝葉斯流量模型

$Y_i \sim N(m_i, s^2)$
$m_i = a + b X_i$
$a \sim N(400, 100^2)$
$b \sim N(0, 200^2)$
$s \sim \text{Unif}(0, 200)$

在 RJAGS 中「定義」貝葉斯模型

$Y_i \sim N(m_i, s^2)$
$m_i = a + b X_i$
$a \sim N(400, 100^2)$
$b \sim N(0, 200^2)$
$s \sim \text{Unif}(0, 200)$

rail_model_1 <- "model{
    # Likelihood model for Y[i]





    # Prior models for a, b, s




}"

在 RJAGS 中「定義」貝葉斯模型

$Y_i \sim N(m_i, s^2)$
$m_i = a + b X_i$
$a \sim N(400, 100^2)$
$b \sim N(0, 200^2)$
$s \sim \text{Unif}(0, 200)$

rail_model_1 <- "model{
    # Likelihood model for Y[i]
    for(i in 1:length(Y)) {
        Y[i] ~ dnorm(m[i], s^(-2))

    }

    # Prior models for a, b, s
    a ~ dnorm(400, 100^(-2))
    s ~ dunif(0, 200)


}"

在 RJAGS 中「定義」貝葉斯模型

m[i] <- a + b[X[i]]

X[1] = 週末，X[2] = 平日
b 有 2 個層級：b[1]、b[2]
週末趨勢（$m_i = a$）
m[i] <- a + b[1]

rail_model_1 <- "model{
    # Likelihood model for Y[i]
    for(i in 1:length(Y)) {
        Y[i] ~ dnorm(m[i], s^(-2))
        m[i] <- a + b[X[i]]
    }

    # Prior models for a, b, s
    a ~ dnorm(400, 100^(-2))
    s ~ dunif(0, 200)


}"

在 RJAGS 中「定義」貝葉斯模型

m[i] <- a + b[X[i]]

X[1] = 週末，X[2] = 平日
b 有 2 個層級：b[1]、b[2]
週末趨勢（$m_i = a$）
m[i] <- a + b[1]
b[1] <- 0

rail_model_1 <- "model{
    # Likelihood model for Y[i]
    for(i in 1:length(Y)) {
        Y[i] ~ dnorm(m[i], s^(-2))
        m[i] <- a + b[X[i]]
    }

    # Prior models for a, b, s
    a ~ dnorm(400, 100^(-2))
    s ~ dunif(0, 200)
    b[1] <- 0

}"

在 RJAGS 中「定義」貝葉斯模型

m[i] <- a + b[X[i]]

X[1] = 週末，X[2] = 平日
b 有 2 個層級：b[1]、b[2]
週末趨勢（$m_i = a$）
m[i] <- a + b[1]
b[1] <- 0
平日（$m_i = a + b$）
m[i] <- a + b[2]

rail_model_1 <- "model{
    # Likelihood model for Y[i]
    for(i in 1:length(Y)) {
        Y[i] ~ dnorm(m[i], s^(-2))
        m[i] <- a + b[X[i]]
    }

    # Prior models for a, b, s
    a ~ dnorm(400, 100^(-2))
    s ~ dunif(0, 200)
    b[1] <- 0
    b[2] ~ dnorm(0, 200^(-2))
}"

b[2] ~ dnorm(0, 200^(-2))

一起來練習吧！

使用 RJAGS 的貝氏建模