Estimated probability of each response by question and group
Copy Right @ Prof Dr Suhartono M.Kom @ Magister Informatika UIN Maulana Malik Ibrahim@ (http://suhartono.lecturer.uin-malang.ac.id/riwayat-hidup/)
21 Januari 2021
Sukses: program Stan pertama
Misalkan kita memiliki N responden survei yang masing-masing menjawab pertanyaan Q, dimana setiap pertanyaan memiliki A jawaban yang memungkinkan.
Ada dua kelompok responden yang ada: (1) membandingkan bagaimana probabilitas tanggapan mereka berbeda, dan (2) diberikan serangkaian tanggapan baru, memprediksi kelompok mana mereka termasuk.
Sebelum Anda menulis kode apa pun, saya sangat menyarankan Anda untuk menuliskan model Anda di atas kertas. Mencoret-coret sesuatu seperti berikut ini selalu menjelaskan pemikiran saya.
Misalkan Ri, j adalah jawaban responden ke-i untuk pertanyaan ke-j, dan misalkan yi 0, {0,1} adalah keanggotaan kelompoknya. Misalkan
Ri,j|θj∼Categorical(θj) θj|yi∼Dirichlet(5)
Jadi kami membiarkan setiap kelompok memiliki distribusi respons yang berbeda untuk setiap pertanyaan, dan mengecilkan distribusi ini satu sama lain dengan Dirichlet prior dengan 5 pseudo-counts di setiap kategori respons. Regularisasi ini masuk akal jika secara a-priori Anda mengharapkan kedua grup merespons dengan cara yang serupa.
Mari kita mulai dengan membuat beberapa data palsu dari proses pembuatan data yang kami usulkan. Pertama kita memasukkan beberapa parameter ukuran masalah dan Dirichlet sebelumnya:
Kita akan coba 24 pertanyaan,dengan 5 kategori jawaban.
library(tidyverse)## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──## ✓ ggplot2 3.3.3 ✓ purrr 0.3.4
## ✓ tibble 3.0.5 ✓ dplyr 1.0.3
## ✓ tidyr 1.1.2 ✓ stringr 1.4.0
## ✓ readr 1.4.0 ✓ forcats 0.5.0## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()set.seed(27)
Q <- 24
A <- 5
alpha <- rep(5, A)
theta_0 <- gtools::rdirichlet(Q, alpha)
theta_1 <- gtools::rdirichlet(Q, alpha)Kita mengambil sampel sekali:
sample_one <- function(theta) {
R <- numeric(Q)
for (q in 1:Q)
R[q] <- sample(1:A, 1, prob = theta[q, ])
names(R) <- paste0("q", 1:Q)
as.list(R)
}Kita lakukan pengambilan sampel sebanyak n kali
sample_n <- function(theta, n, id) {
samples <- map_dfr(1:n, ~sample_one(theta))
samples <- add_column(samples, y = id, .before = TRUE)
mutate_all(samples, as.integer)
}Kemudian kita mengambil sampel 55 kali dari kelompok 0 dan 35 kali dari kelompok 1, dan melihat data yang dihasilkan:
df <- sample_n(theta_0, 55, id = 0) %>%
bind_rows(sample_n(theta_1, 55, id = 1))
df## # A tibble: 110 x 25
## y q1 q2 q3 q4 q5 q6 q7 q8 q9 q10 q11 q12
## <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
## 1 0 4 2 1 3 3 3 5 3 4 3 5 1
## 2 0 4 3 2 5 4 5 2 1 5 4 2 5
## 3 0 1 2 5 3 1 3 2 1 4 1 2 1
## 4 0 4 2 1 5 2 1 3 1 4 2 4 3
## 5 0 5 2 5 5 2 2 4 2 4 5 1 4
## 6 0 1 2 2 4 5 2 1 2 2 3 4 2
## 7 0 1 2 4 3 1 1 1 2 3 5 2 1
## 8 0 1 1 5 4 3 5 3 1 4 1 1 3
## 9 0 1 2 4 3 5 3 2 1 3 5 4 3
## 10 0 1 4 4 5 5 2 2 1 3 5 4 4
## # … with 100 more rows, and 12 more variables: q13 <int>, q14 <int>, q15 <int>,
## # q16 <int>, q17 <int>, q18 <int>, q19 <int>, q20 <int>, q21 <int>,
## # q22 <int>, q23 <int>, q24 <int>library(rstan)## Loading required package: StanHeaders## rstan (Version 2.21.2, GitRev: 2e1f913d3ca3)## For execution on a local, multicore CPU with excess RAM we recommend calling
## options(mc.cores = parallel::detectCores()).
## To avoid recompilation of unchanged Stan programs, we recommend calling
## rstan_options(auto_write = TRUE)##
## Attaching package: 'rstan'## The following object is masked from 'package:tidyr':
##
## extractrstan_options(auto_write = TRUE)model <- stan_model("survey_1.stan")data <- list(
R = as.matrix(select(df, -y)),
N = nrow(df),
y = df$y,
Q = Q,
A = A
)
data## $R
## q1 q2 q3 q4 q5 q6 q7 q8 q9 q10 q11 q12 q13 q14 q15 q16 q17 q18 q19 q20
## [1,] 4 2 1 3 3 3 5 3 4 3 5 1 5 4 3 2 3 2 1 3
## [2,] 4 3 2 5 4 5 2 1 5 4 2 5 5 2 2 1 4 2 4 4
## [3,] 1 2 5 3 1 3 2 1 4 1 2 1 2 2 3 2 1 2 5 4
## [4,] 4 2 1 5 2 1 3 1 4 2 4 3 5 2 3 3 5 1 1 4
## [5,] 5 2 5 5 2 2 4 2 4 5 1 4 1 4 3 4 1 2 2 4
## [6,] 1 2 2 4 5 2 1 2 2 3 4 2 4 4 2 5 3 1 3 1
## [7,] 1 2 4 3 1 1 1 2 3 5 2 1 1 3 5 4 2 4 5 3
## [8,] 1 1 5 4 3 5 3 1 4 1 1 3 2 5 3 5 5 1 1 4
## [9,] 1 2 4 3 5 3 2 1 3 5 4 3 1 2 5 4 4 4 2 5
## [10,] 1 4 4 5 5 2 2 1 3 5 4 4 1 2 1 3 4 3 2 4
## [11,] 1 4 5 1 1 2 2 1 4 5 2 2 3 5 3 1 2 1 4 4
## [12,] 1 5 5 2 2 1 2 3 4 1 2 3 3 2 5 1 3 1 1 5
## [13,] 5 2 1 3 3 5 2 5 3 5 2 1 1 2 1 5 2 4 3 5
## [14,] 5 1 3 1 2 2 2 5 4 5 1 2 5 5 5 1 1 2 2 4
## [15,] 4 2 4 5 1 3 5 2 4 5 2 1 3 5 5 4 3 4 2 2
## [16,] 4 4 5 4 2 1 2 2 5 3 2 2 4 4 2 5 2 3 4 3
## [17,] 1 1 2 3 2 5 4 4 4 5 1 2 3 3 4 4 2 3 2 3
## [18,] 2 3 2 1 2 1 2 2 5 1 2 2 2 4 3 4 5 4 5 2
## [19,] 3 1 3 2 2 4 2 3 5 2 4 5 4 1 2 3 3 2 3 5
## [20,] 2 5 4 2 2 4 2 2 3 3 1 2 3 1 1 1 3 2 1 4
## [21,] 1 1 5 3 2 2 3 5 2 3 1 3 4 1 2 4 2 2 1 2
## [22,] 2 1 3 5 2 1 2 1 3 1 1 5 3 4 3 4 5 2 3 5
## [23,] 3 2 5 3 3 2 2 2 4 1 2 3 3 4 5 4 1 3 2 5
## [24,] 5 1 2 5 3 4 3 3 5 1 2 1 3 1 3 5 1 1 1 5
## [25,] 1 5 1 5 2 5 4 1 5 4 2 3 4 4 5 3 3 5 4 3
## [26,] 2 4 2 3 2 2 5 5 3 1 2 4 3 4 5 3 1 3 1 4
## [27,] 1 5 5 5 1 5 4 2 3 5 1 2 4 5 4 5 3 1 1 4
## [28,] 2 2 1 1 2 3 2 1 4 1 2 1 4 3 3 2 1 3 1 4
## [29,] 4 3 4 3 4 5 3 2 5 5 3 2 3 1 3 3 4 3 4 2
## [30,] 3 2 1 3 1 3 3 2 3 5 5 2 4 2 3 4 4 2 3 3
## [31,] 4 4 1 3 5 5 5 1 5 3 2 3 4 3 5 5 5 3 5 3
## [32,] 2 1 3 1 2 5 2 2 5 5 1 3 1 2 2 5 4 1 5 4
## [33,] 1 4 5 3 2 5 2 3 5 1 4 3 2 3 4 5 3 5 3 5
## [34,] 2 1 2 3 4 3 4 1 1 4 2 3 4 1 3 2 2 5 2 5
## [35,] 4 2 5 3 1 2 4 3 3 1 2 2 3 3 2 2 2 2 1 3
## [36,] 1 5 1 4 1 3 3 1 4 5 4 1 4 3 5 5 2 5 1 5
## [37,] 2 1 3 4 2 5 3 2 2 1 5 5 2 3 3 3 4 5 4 1
## [38,] 1 4 5 4 2 3 2 5 3 4 1 3 3 4 4 3 3 2 2 5
## [39,] 4 2 3 1 1 3 2 2 3 1 4 3 3 4 4 2 2 3 3 5
## [40,] 5 1 3 1 2 2 3 1 1 1 1 3 3 1 3 1 3 1 4 3
## [41,] 5 5 3 1 5 5 2 2 5 4 4 3 3 2 5 4 4 1 4 4
## [42,] 1 1 5 2 2 3 5 2 1 2 2 5 4 2 2 5 1 4 1 1
## [43,] 1 2 2 1 2 2 4 2 3 1 3 1 4 3 3 2 4 3 5 3
## [44,] 3 5 3 1 5 4 3 1 4 1 2 2 3 2 3 4 5 2 2 3
## [45,] 5 3 2 1 5 5 3 2 1 5 4 3 2 3 2 5 3 5 3 5
## [46,] 1 1 1 4 2 3 1 4 3 5 2 5 5 2 5 4 5 1 5 3
## [47,] 2 1 2 3 2 5 4 4 1 4 1 2 2 4 4 5 4 3 1 5
## [48,] 2 3 1 3 1 1 2 2 1 1 2 3 3 2 2 5 4 2 2 2
## [49,] 5 2 1 3 2 1 4 5 1 1 5 3 3 2 3 3 2 2 2 4
## [50,] 4 2 1 4 2 3 1 2 5 5 5 1 4 3 3 3 3 4 3 3
## [51,] 2 2 1 3 3 5 5 3 4 1 1 4 5 1 4 2 5 4 4 4
## [52,] 1 2 5 1 4 1 2 2 3 5 4 5 2 2 3 4 3 4 1 3
## [53,] 4 2 1 3 3 4 1 2 3 5 2 3 4 4 1 3 3 5 1 4
## [54,] 3 3 2 5 1 2 1 2 3 5 4 2 3 2 4 5 3 1 2 3
## [55,] 4 2 2 3 1 3 5 4 3 1 1 5 1 4 3 1 4 3 4 3
## [56,] 4 1 2 3 4 4 1 2 4 4 1 4 4 1 4 3 3 2 3 5
## [57,] 1 1 2 1 5 3 4 2 3 3 4 3 1 3 4 2 2 5 1 4
## [58,] 2 3 2 3 5 4 4 3 5 4 1 2 3 2 4 2 1 1 4 5
## [59,] 1 5 1 4 5 1 1 2 4 3 5 1 2 2 2 1 2 3 2 2
## [60,] 2 2 4 4 5 4 3 3 5 3 1 1 2 2 5 2 3 5 2 2
## [61,] 1 5 4 5 1 2 3 3 2 5 4 5 3 2 3 5 1 3 3 4
## [62,] 1 3 5 4 1 3 1 4 2 2 4 3 2 1 2 2 3 4 2 4
## [63,] 3 3 1 5 5 1 1 3 5 2 3 3 3 4 3 2 5 4 3 3
## [64,] 4 5 3 4 1 2 5 2 5 4 3 5 4 1 2 3 3 4 1 4
## [65,] 5 2 2 1 5 4 5 4 2 5 3 3 3 2 3 1 5 2 3 1
## [66,] 4 2 2 3 5 4 4 3 2 3 1 2 4 4 5 2 4 5 4 3
## [67,] 4 2 2 4 2 5 3 5 5 4 3 3 2 2 4 2 3 3 4 5
## [68,] 5 2 1 1 5 4 5 3 1 1 1 5 4 1 2 1 1 1 3 4
## [69,] 5 1 5 4 1 4 2 1 3 5 3 3 2 3 3 2 4 2 5 3
## [70,] 4 4 3 3 1 1 1 3 5 2 1 5 4 2 3 2 1 4 4 1
## [71,] 4 1 2 4 1 1 3 2 5 1 2 5 1 1 3 5 3 1 4 1
## [72,] 5 5 5 4 5 5 5 4 2 1 5 2 4 4 4 4 3 3 3 1
## [73,] 2 1 3 4 1 4 4 4 3 3 3 1 1 3 3 3 5 4 3 4
## [74,] 2 3 2 4 1 5 4 3 5 3 2 5 3 4 5 1 3 2 4 1
## [75,] 2 3 5 2 4 3 3 3 5 3 3 1 3 2 4 4 2 4 3 2
## [76,] 5 2 1 3 4 4 4 1 5 1 1 2 5 1 5 3 3 2 1 3
## [77,] 1 3 5 2 3 5 5 5 4 2 2 5 4 4 3 5 2 3 4 3
## [78,] 1 2 4 1 2 4 5 1 3 3 4 5 3 2 4 4 1 2 5 1
## [79,] 4 2 4 2 4 4 3 3 2 3 3 5 3 1 5 2 3 3 4 5
## [80,] 1 5 2 1 5 1 3 4 5 4 1 5 5 4 1 4 1 2 3 3
## [81,] 2 3 2 4 5 4 4 3 5 3 4 3 3 2 3 2 5 1 4 4
## [82,] 5 1 2 4 3 5 2 5 2 4 4 2 1 1 3 2 3 3 5 1
## [83,] 3 4 2 5 4 4 4 5 5 4 1 5 2 4 3 4 5 4 5 4
## [84,] 1 3 2 3 2 1 1 3 2 1 2 3 1 3 3 2 3 1 3 5
## [85,] 2 1 2 2 5 4 1 5 3 2 2 3 1 2 4 2 3 5 2 5
## [86,] 2 4 1 5 1 1 4 5 5 5 1 5 3 4 2 1 3 1 3 1
## [87,] 3 3 5 4 5 1 5 3 5 4 3 3 3 4 2 2 2 4 4 3
## [88,] 5 3 2 4 4 3 5 1 4 4 3 5 1 2 4 5 4 5 3 4
## [89,] 1 1 4 2 2 5 3 3 3 3 4 2 2 1 1 4 1 1 4 1
## [90,] 3 4 2 4 2 2 1 2 4 4 4 1 4 2 1 3 3 3 5 4
## [91,] 1 1 2 4 2 4 1 5 5 3 4 5 1 4 5 4 4 4 4 3
## [92,] 2 4 2 2 1 4 2 2 1 4 3 1 1 4 2 1 1 4 4 5
## [93,] 1 4 2 4 5 3 5 2 5 5 3 3 5 3 5 1 3 1 4 5
## [94,] 4 1 2 4 5 4 5 3 3 3 4 2 1 3 4 1 5 1 3 4
## [95,] 1 4 3 2 5 5 3 3 1 5 5 2 1 4 4 5 1 5 4 1
## [96,] 1 2 5 4 2 2 5 1 5 3 3 2 1 5 3 4 2 3 4 3
## [97,] 5 2 4 3 4 4 4 1 1 3 3 3 1 4 4 4 2 5 4 2
## [98,] 2 1 5 2 4 2 5 1 3 4 4 4 1 1 2 2 4 4 3 4
## [99,] 5 2 3 2 4 4 1 5 3 2 5 2 2 1 4 3 5 2 1 3
## [100,] 1 2 1 1 2 3 5 1 5 4 5 1 1 4 4 3 3 1 2 1
## [101,] 2 1 4 3 1 5 1 5 5 2 1 3 1 3 5 2 3 1 2 3
## [102,] 2 4 3 4 1 4 4 5 5 3 1 3 2 2 2 3 1 1 5 3
## [103,] 5 2 1 4 1 5 4 1 5 3 3 1 5 1 3 1 5 4 4 5
## [104,] 3 2 2 2 4 5 5 2 4 3 1 1 3 4 5 2 5 4 4 3
## [105,] 2 1 2 4 3 4 1 5 3 2 4 3 4 2 4 3 5 4 4 3
## [106,] 4 2 1 4 4 2 4 3 3 1 5 1 5 2 4 2 2 2 1 1
## [107,] 4 5 4 4 5 4 5 3 1 3 5 2 1 4 5 3 5 2 4 5
## [108,] 2 5 2 1 1 5 5 3 1 5 4 2 4 1 5 1 5 5 4 4
## [109,] 1 5 2 2 1 4 1 4 3 5 2 3 4 3 4 5 2 5 3 1
## [110,] 2 3 2 4 4 4 3 3 5 1 3 3 4 4 4 2 1 1 4 3
## q21 q22 q23 q24
## [1,] 5 3 1 1
## [2,] 4 4 1 2
## [3,] 5 3 2 2
## [4,] 2 4 4 3
## [5,] 5 2 2 1
## [6,] 3 4 3 2
## [7,] 4 1 3 1
## [8,] 5 1 2 3
## [9,] 5 1 4 4
## [10,] 1 3 3 3
## [11,] 1 4 1 1
## [12,] 5 5 2 2
## [13,] 4 3 1 2
## [14,] 3 1 4 5
## [15,] 5 2 4 5
## [16,] 1 2 3 2
## [17,] 5 4 2 4
## [18,] 3 4 1 3
## [19,] 4 4 4 3
## [20,] 1 2 1 4
## [21,] 4 4 5 3
## [22,] 1 3 3 3
## [23,] 1 1 4 5
## [24,] 5 4 5 2
## [25,] 4 4 5 1
## [26,] 4 3 4 2
## [27,] 5 2 3 5
## [28,] 1 4 4 1
## [29,] 2 2 3 2
## [30,] 3 5 1 3
## [31,] 2 1 1 5
## [32,] 1 2 4 1
## [33,] 5 1 3 1
## [34,] 4 5 4 1
## [35,] 4 1 4 5
## [36,] 4 5 1 2
## [37,] 4 2 3 2
## [38,] 2 4 1 5
## [39,] 1 5 3 2
## [40,] 5 4 3 4
## [41,] 4 1 3 5
## [42,] 2 1 1 5
## [43,] 4 1 4 5
## [44,] 2 3 4 1
## [45,] 1 2 5 3
## [46,] 2 5 4 2
## [47,] 5 1 5 3
## [48,] 5 2 3 3
## [49,] 3 3 5 5
## [50,] 4 4 4 4
## [51,] 2 1 1 1
## [52,] 2 4 4 5
## [53,] 1 3 4 5
## [54,] 1 1 4 4
## [55,] 5 2 5 1
## [56,] 4 5 1 3
## [57,] 3 5 5 5
## [58,] 4 5 4 5
## [59,] 5 5 2 3
## [60,] 3 5 1 3
## [61,] 4 5 2 1
## [62,] 1 2 3 4
## [63,] 3 3 4 3
## [64,] 2 5 5 4
## [65,] 5 3 2 2
## [66,] 5 2 4 2
## [67,] 5 5 5 3
## [68,] 3 2 4 2
## [69,] 3 3 4 2
## [70,] 5 5 3 3
## [71,] 5 1 1 1
## [72,] 1 1 4 4
## [73,] 1 5 3 2
## [74,] 4 5 3 2
## [75,] 1 1 3 5
## [76,] 2 1 3 2
## [77,] 3 3 1 2
## [78,] 1 4 3 5
## [79,] 3 3 3 2
## [80,] 5 1 2 5
## [81,] 5 1 2 1
## [82,] 2 5 3 2
## [83,] 1 5 3 3
## [84,] 2 5 3 1
## [85,] 5 4 3 5
## [86,] 1 3 3 1
## [87,] 1 1 3 5
## [88,] 1 5 2 2
## [89,] 2 4 2 2
## [90,] 4 5 2 3
## [91,] 1 1 5 4
## [92,] 5 5 3 3
## [93,] 4 1 1 1
## [94,] 5 1 3 1
## [95,] 1 4 5 5
## [96,] 2 4 1 3
## [97,] 1 4 3 3
## [98,] 5 5 5 5
## [99,] 3 1 2 3
## [100,] 5 5 4 3
## [101,] 3 1 4 5
## [102,] 2 1 3 5
## [103,] 5 3 4 1
## [104,] 1 4 3 2
## [105,] 1 1 3 4
## [106,] 5 4 1 4
## [107,] 2 5 3 5
## [108,] 5 5 4 3
## [109,] 2 3 1 3
## [110,] 1 1 3 4
##
## $N
## [1] 110
##
## $y
## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [38] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##
## $Q
## [1] 24
##
## $A
## [1] 5fit1 <- sampling(model, data = data, chains = 2, iter = 1000, refresh = 0)print(fit1, pars = "alpha", probs = c(0.025, 0.5, 0.975))## Inference for Stan model: survey_1.
## 2 chains, each with iter=1000; warmup=500; thin=1;
## post-warmup draws per chain=500, total post-warmup draws=1000.
##
## mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
## alpha[1] 12.96 0.19 3.21 7.94 12.52 20.07 271 1
## alpha[2] 13.24 0.20 3.36 7.89 12.85 20.78 288 1
## alpha[3] 13.95 0.22 3.56 8.51 13.47 21.90 260 1
## alpha[4] 12.96 0.21 3.26 8.04 12.57 20.53 253 1
## alpha[5] 12.00 0.18 3.04 7.38 11.53 18.88 279 1
##
## Samples were drawn using NUTS(diag_e) at Fri Jan 22 07:35:53 2021.
## For each parameter, 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).shinystan::launch_shinystan(fit1)##
## Launching ShinyStan interface... for large models this may take some time.## Loading required package: shiny##
## Listening on http://127.0.0.1:4901model1 <- stan_model("survey_2.stan")data <- list(
R = as.matrix(select(df, -y)),
N = nrow(df),
y = df$y,
Q = Q,
A = A,
alpha = alpha
)
data## $R
## q1 q2 q3 q4 q5 q6 q7 q8 q9 q10 q11 q12 q13 q14 q15 q16 q17 q18 q19 q20
## [1,] 4 2 1 3 3 3 5 3 4 3 5 1 5 4 3 2 3 2 1 3
## [2,] 4 3 2 5 4 5 2 1 5 4 2 5 5 2 2 1 4 2 4 4
## [3,] 1 2 5 3 1 3 2 1 4 1 2 1 2 2 3 2 1 2 5 4
## [4,] 4 2 1 5 2 1 3 1 4 2 4 3 5 2 3 3 5 1 1 4
## [5,] 5 2 5 5 2 2 4 2 4 5 1 4 1 4 3 4 1 2 2 4
## [6,] 1 2 2 4 5 2 1 2 2 3 4 2 4 4 2 5 3 1 3 1
## [7,] 1 2 4 3 1 1 1 2 3 5 2 1 1 3 5 4 2 4 5 3
## [8,] 1 1 5 4 3 5 3 1 4 1 1 3 2 5 3 5 5 1 1 4
## [9,] 1 2 4 3 5 3 2 1 3 5 4 3 1 2 5 4 4 4 2 5
## [10,] 1 4 4 5 5 2 2 1 3 5 4 4 1 2 1 3 4 3 2 4
## [11,] 1 4 5 1 1 2 2 1 4 5 2 2 3 5 3 1 2 1 4 4
## [12,] 1 5 5 2 2 1 2 3 4 1 2 3 3 2 5 1 3 1 1 5
## [13,] 5 2 1 3 3 5 2 5 3 5 2 1 1 2 1 5 2 4 3 5
## [14,] 5 1 3 1 2 2 2 5 4 5 1 2 5 5 5 1 1 2 2 4
## [15,] 4 2 4 5 1 3 5 2 4 5 2 1 3 5 5 4 3 4 2 2
## [16,] 4 4 5 4 2 1 2 2 5 3 2 2 4 4 2 5 2 3 4 3
## [17,] 1 1 2 3 2 5 4 4 4 5 1 2 3 3 4 4 2 3 2 3
## [18,] 2 3 2 1 2 1 2 2 5 1 2 2 2 4 3 4 5 4 5 2
## [19,] 3 1 3 2 2 4 2 3 5 2 4 5 4 1 2 3 3 2 3 5
## [20,] 2 5 4 2 2 4 2 2 3 3 1 2 3 1 1 1 3 2 1 4
## [21,] 1 1 5 3 2 2 3 5 2 3 1 3 4 1 2 4 2 2 1 2
## [22,] 2 1 3 5 2 1 2 1 3 1 1 5 3 4 3 4 5 2 3 5
## [23,] 3 2 5 3 3 2 2 2 4 1 2 3 3 4 5 4 1 3 2 5
## [24,] 5 1 2 5 3 4 3 3 5 1 2 1 3 1 3 5 1 1 1 5
## [25,] 1 5 1 5 2 5 4 1 5 4 2 3 4 4 5 3 3 5 4 3
## [26,] 2 4 2 3 2 2 5 5 3 1 2 4 3 4 5 3 1 3 1 4
## [27,] 1 5 5 5 1 5 4 2 3 5 1 2 4 5 4 5 3 1 1 4
## [28,] 2 2 1 1 2 3 2 1 4 1 2 1 4 3 3 2 1 3 1 4
## [29,] 4 3 4 3 4 5 3 2 5 5 3 2 3 1 3 3 4 3 4 2
## [30,] 3 2 1 3 1 3 3 2 3 5 5 2 4 2 3 4 4 2 3 3
## [31,] 4 4 1 3 5 5 5 1 5 3 2 3 4 3 5 5 5 3 5 3
## [32,] 2 1 3 1 2 5 2 2 5 5 1 3 1 2 2 5 4 1 5 4
## [33,] 1 4 5 3 2 5 2 3 5 1 4 3 2 3 4 5 3 5 3 5
## [34,] 2 1 2 3 4 3 4 1 1 4 2 3 4 1 3 2 2 5 2 5
## [35,] 4 2 5 3 1 2 4 3 3 1 2 2 3 3 2 2 2 2 1 3
## [36,] 1 5 1 4 1 3 3 1 4 5 4 1 4 3 5 5 2 5 1 5
## [37,] 2 1 3 4 2 5 3 2 2 1 5 5 2 3 3 3 4 5 4 1
## [38,] 1 4 5 4 2 3 2 5 3 4 1 3 3 4 4 3 3 2 2 5
## [39,] 4 2 3 1 1 3 2 2 3 1 4 3 3 4 4 2 2 3 3 5
## [40,] 5 1 3 1 2 2 3 1 1 1 1 3 3 1 3 1 3 1 4 3
## [41,] 5 5 3 1 5 5 2 2 5 4 4 3 3 2 5 4 4 1 4 4
## [42,] 1 1 5 2 2 3 5 2 1 2 2 5 4 2 2 5 1 4 1 1
## [43,] 1 2 2 1 2 2 4 2 3 1 3 1 4 3 3 2 4 3 5 3
## [44,] 3 5 3 1 5 4 3 1 4 1 2 2 3 2 3 4 5 2 2 3
## [45,] 5 3 2 1 5 5 3 2 1 5 4 3 2 3 2 5 3 5 3 5
## [46,] 1 1 1 4 2 3 1 4 3 5 2 5 5 2 5 4 5 1 5 3
## [47,] 2 1 2 3 2 5 4 4 1 4 1 2 2 4 4 5 4 3 1 5
## [48,] 2 3 1 3 1 1 2 2 1 1 2 3 3 2 2 5 4 2 2 2
## [49,] 5 2 1 3 2 1 4 5 1 1 5 3 3 2 3 3 2 2 2 4
## [50,] 4 2 1 4 2 3 1 2 5 5 5 1 4 3 3 3 3 4 3 3
## [51,] 2 2 1 3 3 5 5 3 4 1 1 4 5 1 4 2 5 4 4 4
## [52,] 1 2 5 1 4 1 2 2 3 5 4 5 2 2 3 4 3 4 1 3
## [53,] 4 2 1 3 3 4 1 2 3 5 2 3 4 4 1 3 3 5 1 4
## [54,] 3 3 2 5 1 2 1 2 3 5 4 2 3 2 4 5 3 1 2 3
## [55,] 4 2 2 3 1 3 5 4 3 1 1 5 1 4 3 1 4 3 4 3
## [56,] 4 1 2 3 4 4 1 2 4 4 1 4 4 1 4 3 3 2 3 5
## [57,] 1 1 2 1 5 3 4 2 3 3 4 3 1 3 4 2 2 5 1 4
## [58,] 2 3 2 3 5 4 4 3 5 4 1 2 3 2 4 2 1 1 4 5
## [59,] 1 5 1 4 5 1 1 2 4 3 5 1 2 2 2 1 2 3 2 2
## [60,] 2 2 4 4 5 4 3 3 5 3 1 1 2 2 5 2 3 5 2 2
## [61,] 1 5 4 5 1 2 3 3 2 5 4 5 3 2 3 5 1 3 3 4
## [62,] 1 3 5 4 1 3 1 4 2 2 4 3 2 1 2 2 3 4 2 4
## [63,] 3 3 1 5 5 1 1 3 5 2 3 3 3 4 3 2 5 4 3 3
## [64,] 4 5 3 4 1 2 5 2 5 4 3 5 4 1 2 3 3 4 1 4
## [65,] 5 2 2 1 5 4 5 4 2 5 3 3 3 2 3 1 5 2 3 1
## [66,] 4 2 2 3 5 4 4 3 2 3 1 2 4 4 5 2 4 5 4 3
## [67,] 4 2 2 4 2 5 3 5 5 4 3 3 2 2 4 2 3 3 4 5
## [68,] 5 2 1 1 5 4 5 3 1 1 1 5 4 1 2 1 1 1 3 4
## [69,] 5 1 5 4 1 4 2 1 3 5 3 3 2 3 3 2 4 2 5 3
## [70,] 4 4 3 3 1 1 1 3 5 2 1 5 4 2 3 2 1 4 4 1
## [71,] 4 1 2 4 1 1 3 2 5 1 2 5 1 1 3 5 3 1 4 1
## [72,] 5 5 5 4 5 5 5 4 2 1 5 2 4 4 4 4 3 3 3 1
## [73,] 2 1 3 4 1 4 4 4 3 3 3 1 1 3 3 3 5 4 3 4
## [74,] 2 3 2 4 1 5 4 3 5 3 2 5 3 4 5 1 3 2 4 1
## [75,] 2 3 5 2 4 3 3 3 5 3 3 1 3 2 4 4 2 4 3 2
## [76,] 5 2 1 3 4 4 4 1 5 1 1 2 5 1 5 3 3 2 1 3
## [77,] 1 3 5 2 3 5 5 5 4 2 2 5 4 4 3 5 2 3 4 3
## [78,] 1 2 4 1 2 4 5 1 3 3 4 5 3 2 4 4 1 2 5 1
## [79,] 4 2 4 2 4 4 3 3 2 3 3 5 3 1 5 2 3 3 4 5
## [80,] 1 5 2 1 5 1 3 4 5 4 1 5 5 4 1 4 1 2 3 3
## [81,] 2 3 2 4 5 4 4 3 5 3 4 3 3 2 3 2 5 1 4 4
## [82,] 5 1 2 4 3 5 2 5 2 4 4 2 1 1 3 2 3 3 5 1
## [83,] 3 4 2 5 4 4 4 5 5 4 1 5 2 4 3 4 5 4 5 4
## [84,] 1 3 2 3 2 1 1 3 2 1 2 3 1 3 3 2 3 1 3 5
## [85,] 2 1 2 2 5 4 1 5 3 2 2 3 1 2 4 2 3 5 2 5
## [86,] 2 4 1 5 1 1 4 5 5 5 1 5 3 4 2 1 3 1 3 1
## [87,] 3 3 5 4 5 1 5 3 5 4 3 3 3 4 2 2 2 4 4 3
## [88,] 5 3 2 4 4 3 5 1 4 4 3 5 1 2 4 5 4 5 3 4
## [89,] 1 1 4 2 2 5 3 3 3 3 4 2 2 1 1 4 1 1 4 1
## [90,] 3 4 2 4 2 2 1 2 4 4 4 1 4 2 1 3 3 3 5 4
## [91,] 1 1 2 4 2 4 1 5 5 3 4 5 1 4 5 4 4 4 4 3
## [92,] 2 4 2 2 1 4 2 2 1 4 3 1 1 4 2 1 1 4 4 5
## [93,] 1 4 2 4 5 3 5 2 5 5 3 3 5 3 5 1 3 1 4 5
## [94,] 4 1 2 4 5 4 5 3 3 3 4 2 1 3 4 1 5 1 3 4
## [95,] 1 4 3 2 5 5 3 3 1 5 5 2 1 4 4 5 1 5 4 1
## [96,] 1 2 5 4 2 2 5 1 5 3 3 2 1 5 3 4 2 3 4 3
## [97,] 5 2 4 3 4 4 4 1 1 3 3 3 1 4 4 4 2 5 4 2
## [98,] 2 1 5 2 4 2 5 1 3 4 4 4 1 1 2 2 4 4 3 4
## [99,] 5 2 3 2 4 4 1 5 3 2 5 2 2 1 4 3 5 2 1 3
## [100,] 1 2 1 1 2 3 5 1 5 4 5 1 1 4 4 3 3 1 2 1
## [101,] 2 1 4 3 1 5 1 5 5 2 1 3 1 3 5 2 3 1 2 3
## [102,] 2 4 3 4 1 4 4 5 5 3 1 3 2 2 2 3 1 1 5 3
## [103,] 5 2 1 4 1 5 4 1 5 3 3 1 5 1 3 1 5 4 4 5
## [104,] 3 2 2 2 4 5 5 2 4 3 1 1 3 4 5 2 5 4 4 3
## [105,] 2 1 2 4 3 4 1 5 3 2 4 3 4 2 4 3 5 4 4 3
## [106,] 4 2 1 4 4 2 4 3 3 1 5 1 5 2 4 2 2 2 1 1
## [107,] 4 5 4 4 5 4 5 3 1 3 5 2 1 4 5 3 5 2 4 5
## [108,] 2 5 2 1 1 5 5 3 1 5 4 2 4 1 5 1 5 5 4 4
## [109,] 1 5 2 2 1 4 1 4 3 5 2 3 4 3 4 5 2 5 3 1
## [110,] 2 3 2 4 4 4 3 3 5 1 3 3 4 4 4 2 1 1 4 3
## q21 q22 q23 q24
## [1,] 5 3 1 1
## [2,] 4 4 1 2
## [3,] 5 3 2 2
## [4,] 2 4 4 3
## [5,] 5 2 2 1
## [6,] 3 4 3 2
## [7,] 4 1 3 1
## [8,] 5 1 2 3
## [9,] 5 1 4 4
## [10,] 1 3 3 3
## [11,] 1 4 1 1
## [12,] 5 5 2 2
## [13,] 4 3 1 2
## [14,] 3 1 4 5
## [15,] 5 2 4 5
## [16,] 1 2 3 2
## [17,] 5 4 2 4
## [18,] 3 4 1 3
## [19,] 4 4 4 3
## [20,] 1 2 1 4
## [21,] 4 4 5 3
## [22,] 1 3 3 3
## [23,] 1 1 4 5
## [24,] 5 4 5 2
## [25,] 4 4 5 1
## [26,] 4 3 4 2
## [27,] 5 2 3 5
## [28,] 1 4 4 1
## [29,] 2 2 3 2
## [30,] 3 5 1 3
## [31,] 2 1 1 5
## [32,] 1 2 4 1
## [33,] 5 1 3 1
## [34,] 4 5 4 1
## [35,] 4 1 4 5
## [36,] 4 5 1 2
## [37,] 4 2 3 2
## [38,] 2 4 1 5
## [39,] 1 5 3 2
## [40,] 5 4 3 4
## [41,] 4 1 3 5
## [42,] 2 1 1 5
## [43,] 4 1 4 5
## [44,] 2 3 4 1
## [45,] 1 2 5 3
## [46,] 2 5 4 2
## [47,] 5 1 5 3
## [48,] 5 2 3 3
## [49,] 3 3 5 5
## [50,] 4 4 4 4
## [51,] 2 1 1 1
## [52,] 2 4 4 5
## [53,] 1 3 4 5
## [54,] 1 1 4 4
## [55,] 5 2 5 1
## [56,] 4 5 1 3
## [57,] 3 5 5 5
## [58,] 4 5 4 5
## [59,] 5 5 2 3
## [60,] 3 5 1 3
## [61,] 4 5 2 1
## [62,] 1 2 3 4
## [63,] 3 3 4 3
## [64,] 2 5 5 4
## [65,] 5 3 2 2
## [66,] 5 2 4 2
## [67,] 5 5 5 3
## [68,] 3 2 4 2
## [69,] 3 3 4 2
## [70,] 5 5 3 3
## [71,] 5 1 1 1
## [72,] 1 1 4 4
## [73,] 1 5 3 2
## [74,] 4 5 3 2
## [75,] 1 1 3 5
## [76,] 2 1 3 2
## [77,] 3 3 1 2
## [78,] 1 4 3 5
## [79,] 3 3 3 2
## [80,] 5 1 2 5
## [81,] 5 1 2 1
## [82,] 2 5 3 2
## [83,] 1 5 3 3
## [84,] 2 5 3 1
## [85,] 5 4 3 5
## [86,] 1 3 3 1
## [87,] 1 1 3 5
## [88,] 1 5 2 2
## [89,] 2 4 2 2
## [90,] 4 5 2 3
## [91,] 1 1 5 4
## [92,] 5 5 3 3
## [93,] 4 1 1 1
## [94,] 5 1 3 1
## [95,] 1 4 5 5
## [96,] 2 4 1 3
## [97,] 1 4 3 3
## [98,] 5 5 5 5
## [99,] 3 1 2 3
## [100,] 5 5 4 3
## [101,] 3 1 4 5
## [102,] 2 1 3 5
## [103,] 5 3 4 1
## [104,] 1 4 3 2
## [105,] 1 1 3 4
## [106,] 5 4 1 4
## [107,] 2 5 3 5
## [108,] 5 5 4 3
## [109,] 2 3 1 3
## [110,] 1 1 3 4
##
## $N
## [1] 110
##
## $y
## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [38] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##
## $Q
## [1] 24
##
## $A
## [1] 5
##
## $alpha
## [1] 5 5 5 5 5Kita akan mengkompilasi ini dan memberikannya alfa, lalu melihat theta_0 yang dihasilkan:
fit2 <- sampling(model1, data = data, chains = 2, iter = 1000, refresh = 0)
print(fit2, pars = "theta_0", probs = c(0.025, 0.5, 0.975))## Inference for Stan model: survey_2.
## 2 chains, each with iter=1000; warmup=500; thin=1;
## post-warmup draws per chain=500, total post-warmup draws=1000.
##
## mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
## theta_0[1,1] 0.30 0 0.05 0.21 0.30 0.40 3000 1
## theta_0[1,2] 0.20 0 0.04 0.12 0.20 0.29 2224 1
## theta_0[1,3] 0.12 0 0.04 0.06 0.12 0.21 2354 1
## theta_0[1,4] 0.21 0 0.04 0.13 0.21 0.30 2364 1
## theta_0[1,5] 0.16 0 0.04 0.09 0.16 0.25 2278 1
## theta_0[2,1] 0.24 0 0.05 0.16 0.23 0.33 1789 1
## theta_0[2,2] 0.32 0 0.05 0.23 0.32 0.43 2333 1
## theta_0[2,3] 0.14 0 0.04 0.07 0.14 0.22 1851 1
## theta_0[2,4] 0.15 0 0.04 0.09 0.15 0.23 1773 1
## theta_0[2,5] 0.15 0 0.04 0.09 0.15 0.24 1712 1
## theta_0[3,1] 0.24 0 0.05 0.15 0.23 0.34 2208 1
## theta_0[3,2] 0.21 0 0.05 0.13 0.21 0.31 2403 1
## theta_0[3,3] 0.18 0 0.04 0.10 0.17 0.26 2981 1
## theta_0[3,4] 0.14 0 0.04 0.08 0.13 0.22 2094 1
## theta_0[3,5] 0.24 0 0.05 0.15 0.23 0.34 1832 1
## theta_0[4,1] 0.21 0 0.05 0.13 0.21 0.31 2562 1
## theta_0[4,2] 0.11 0 0.03 0.06 0.11 0.19 1071 1
## theta_0[4,3] 0.32 0 0.05 0.24 0.32 0.43 1969 1
## theta_0[4,4] 0.16 0 0.04 0.09 0.16 0.26 1780 1
## theta_0[4,5] 0.19 0 0.04 0.11 0.18 0.26 1666 1
## theta_0[5,1] 0.21 0 0.05 0.13 0.21 0.32 2503 1
## theta_0[5,2] 0.38 0 0.05 0.27 0.37 0.49 2031 1
## theta_0[5,3] 0.15 0 0.04 0.08 0.15 0.23 2209 1
## theta_0[5,4] 0.11 0 0.04 0.05 0.11 0.19 1764 1
## theta_0[5,5] 0.15 0 0.04 0.09 0.15 0.24 1787 1
## theta_0[6,1] 0.17 0 0.04 0.11 0.17 0.25 2204 1
## theta_0[6,2] 0.21 0 0.05 0.13 0.21 0.31 2211 1
## theta_0[6,3] 0.24 0 0.05 0.15 0.24 0.34 2795 1
## theta_0[6,4] 0.12 0 0.04 0.06 0.12 0.20 2036 1
## theta_0[6,5] 0.25 0 0.05 0.16 0.25 0.36 2060 1
## theta_0[7,1] 0.14 0 0.04 0.07 0.13 0.22 2331 1
## theta_0[7,2] 0.34 0 0.05 0.24 0.34 0.45 1834 1
## theta_0[7,3] 0.20 0 0.04 0.12 0.20 0.29 1835 1
## theta_0[7,4] 0.18 0 0.04 0.10 0.17 0.27 1904 1
## theta_0[7,5] 0.15 0 0.04 0.08 0.15 0.23 1634 1
## theta_0[8,1] 0.25 0 0.05 0.17 0.25 0.35 2475 1
## theta_0[8,2] 0.35 0 0.05 0.25 0.35 0.46 2087 1
## theta_0[8,3] 0.15 0 0.04 0.08 0.15 0.24 3000 1
## theta_0[8,4] 0.11 0 0.04 0.05 0.11 0.19 1753 1
## theta_0[8,5] 0.14 0 0.04 0.07 0.14 0.22 2704 1
## theta_0[9,1] 0.15 0 0.04 0.08 0.15 0.24 2282 1
## theta_0[9,2] 0.10 0 0.03 0.04 0.10 0.18 2068 1
## theta_0[9,3] 0.29 0 0.05 0.19 0.28 0.40 2101 1
## theta_0[9,4] 0.25 0 0.05 0.17 0.25 0.35 1911 1
## theta_0[9,5] 0.21 0 0.05 0.13 0.21 0.31 2180 1
## theta_0[10,1] 0.32 0 0.05 0.23 0.32 0.41 2134 1
## theta_0[10,2] 0.10 0 0.03 0.05 0.10 0.17 1825 1
## theta_0[10,3] 0.14 0 0.04 0.08 0.13 0.22 2292 1
## theta_0[10,4] 0.14 0 0.04 0.07 0.13 0.21 2179 1
## theta_0[10,5] 0.31 0 0.05 0.22 0.31 0.41 2525 1
## theta_0[11,1] 0.24 0 0.04 0.16 0.24 0.33 2343 1
## theta_0[11,2] 0.34 0 0.05 0.24 0.34 0.44 2153 1
## theta_0[11,3] 0.09 0 0.03 0.04 0.08 0.17 1859 1
## theta_0[11,4] 0.21 0 0.05 0.13 0.21 0.31 2253 1
## theta_0[11,5] 0.13 0 0.04 0.07 0.12 0.20 2096 1
## theta_0[12,1] 0.19 0 0.04 0.11 0.18 0.28 2799 1
## theta_0[12,2] 0.24 0 0.05 0.15 0.24 0.34 2642 1
## theta_0[12,3] 0.30 0 0.06 0.20 0.30 0.41 2405 1
## theta_0[12,4] 0.11 0 0.03 0.05 0.11 0.20 2213 1
## theta_0[12,5] 0.16 0 0.04 0.09 0.16 0.24 2863 1
## theta_0[13,1] 0.15 0 0.04 0.09 0.15 0.23 2300 1
## theta_0[13,2] 0.16 0 0.04 0.09 0.16 0.25 1953 1
## theta_0[13,3] 0.30 0 0.05 0.21 0.30 0.40 2124 1
## theta_0[13,4] 0.25 0 0.05 0.16 0.25 0.36 2183 1
## theta_0[13,5] 0.14 0 0.04 0.07 0.13 0.23 2154 1
## theta_0[14,1] 0.16 0 0.04 0.09 0.16 0.25 2494 1
## theta_0[14,2] 0.27 0 0.05 0.19 0.27 0.37 2551 1
## theta_0[14,3] 0.20 0 0.05 0.12 0.20 0.30 2290 1
## theta_0[14,4] 0.24 0 0.05 0.15 0.23 0.33 2976 1
## theta_0[14,5] 0.12 0 0.04 0.06 0.12 0.20 1891 1
## theta_0[15,1] 0.11 0 0.03 0.06 0.11 0.18 3000 1
## theta_0[15,2] 0.19 0 0.04 0.11 0.18 0.28 1911 1
## theta_0[15,3] 0.33 0 0.06 0.22 0.32 0.44 2953 1
## theta_0[15,4] 0.16 0 0.04 0.09 0.16 0.25 2031 1
## theta_0[15,5] 0.21 0 0.05 0.13 0.21 0.30 1842 1
## theta_0[16,1] 0.15 0 0.04 0.08 0.15 0.24 2099 1
## theta_0[16,2] 0.16 0 0.04 0.09 0.16 0.25 1933 1
## theta_0[16,3] 0.20 0 0.05 0.12 0.20 0.30 1679 1
## theta_0[16,4] 0.24 0 0.05 0.16 0.24 0.33 2915 1
## theta_0[16,5] 0.25 0 0.05 0.16 0.25 0.35 2185 1
## theta_0[17,1] 0.16 0 0.04 0.09 0.16 0.24 2098 1
## theta_0[17,2] 0.20 0 0.04 0.13 0.20 0.29 2510 1
## theta_0[17,3] 0.26 0 0.05 0.17 0.26 0.37 2479 1
## theta_0[17,4] 0.21 0 0.04 0.13 0.21 0.31 2277 1
## theta_0[17,5] 0.16 0 0.04 0.09 0.16 0.25 2114 1
## theta_0[18,1] 0.21 0 0.05 0.12 0.21 0.32 1857 1
## theta_0[18,2] 0.25 0 0.05 0.17 0.25 0.34 1477 1
## theta_0[18,3] 0.21 0 0.04 0.13 0.21 0.31 2325 1
## theta_0[18,4] 0.17 0 0.04 0.10 0.17 0.26 1905 1
## theta_0[18,5] 0.15 0 0.04 0.08 0.15 0.24 1638 1
## theta_0[19,1] 0.26 0 0.05 0.17 0.26 0.37 2872 1
## theta_0[19,2] 0.22 0 0.05 0.14 0.22 0.33 2099 1
## theta_0[19,3] 0.17 0 0.04 0.11 0.17 0.26 2518 1
## theta_0[19,4] 0.19 0 0.04 0.11 0.18 0.28 2036 1
## theta_0[19,5] 0.15 0 0.04 0.08 0.15 0.24 2730 1
## theta_0[20,1] 0.10 0 0.03 0.05 0.10 0.18 2890 1
## theta_0[20,2] 0.12 0 0.04 0.07 0.12 0.20 1716 1
## theta_0[20,3] 0.26 0 0.05 0.18 0.26 0.36 1816 1
## theta_0[20,4] 0.27 0 0.05 0.18 0.27 0.38 2153 1
## theta_0[20,5] 0.24 0 0.05 0.15 0.23 0.34 1673 1
## theta_0[21,1] 0.21 0 0.04 0.14 0.21 0.31 1728 1
## theta_0[21,2] 0.17 0 0.04 0.10 0.17 0.26 2511 1
## theta_0[21,3] 0.12 0 0.04 0.06 0.12 0.20 2924 1
## theta_0[21,4] 0.24 0 0.05 0.15 0.24 0.33 1275 1
## theta_0[21,5] 0.25 0 0.05 0.16 0.25 0.35 1639 1
## theta_0[22,1] 0.24 0 0.05 0.16 0.23 0.34 2175 1
## theta_0[22,2] 0.20 0 0.04 0.13 0.20 0.28 1949 1
## theta_0[22,3] 0.17 0 0.04 0.10 0.17 0.27 2325 1
## theta_0[22,4] 0.25 0 0.05 0.17 0.25 0.34 2671 1
## theta_0[22,5] 0.14 0 0.04 0.07 0.13 0.22 3000 1
## theta_0[23,1] 0.21 0 0.05 0.13 0.21 0.31 3000 1
## theta_0[23,2] 0.13 0 0.04 0.06 0.12 0.21 3000 1
## theta_0[23,3] 0.23 0 0.05 0.14 0.22 0.33 2574 1
## theta_0[23,4] 0.29 0 0.05 0.19 0.28 0.39 2383 1
## theta_0[23,5] 0.15 0 0.04 0.08 0.15 0.23 3000 1
## theta_0[24,1] 0.21 0 0.04 0.13 0.21 0.30 2030 1
## theta_0[24,2] 0.23 0 0.05 0.13 0.22 0.33 3000 1
## theta_0[24,3] 0.20 0 0.05 0.12 0.20 0.29 2449 1
## theta_0[24,4] 0.14 0 0.04 0.07 0.14 0.22 2023 1
## theta_0[24,5] 0.23 0 0.05 0.15 0.23 0.32 2096 1
##
## Samples were drawn using NUTS(diag_e) at Fri Jan 22 07:36:49 2021.
## For each parameter, 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).Pada titik ini saya senang dengan logika inti model, jadi saya menggunakan ShinyStan untuk memeriksa diagnostik MCMC dan memastikan bahwa semua rantai tercampur, dll:
shinystan::launch_shinystan(fit2)##
## Launching ShinyStan interface... for large models this may take some time.##
## Listening on http://127.0.0.1:4901Itu semua diperiksa, yang tidak mengherankan, karena kami memperkirakan banyak proporsi secara efektif, jadi akan aneh jika semuanya mulai meledak.
Tapi sebenarnya, ide awal saya adalah untuk melihat seberapa baik model ini akan memprediksi kelompok responden, jadi kita perlu melakukan beberapa pekerjaan lagi.
Saya ingin P (y = 0) diberi tanggapan seseorang Ri1, …, RiA. Ini adalah salah satu contoh ketika saya diintimidasi oleh bagaimana mengetahuinya, tetapi kemudian semuanya menjadi jelas dengan satu atau dua menit dengan pena dan kertas kerja.
Aturan Bayes memberi kita bahwa P (yi = 0 | Ri1, …, RiA) sama
model2 <- stan_model("survey_3.stan")fit3 <- sampling(model2, data = data, chains = 2, iter = 1000)##
## SAMPLING FOR MODEL 'survey_3' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0.001252 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 12.52 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 1000 [ 0%] (Warmup)
## Chain 1: Iteration: 100 / 1000 [ 10%] (Warmup)
## Chain 1: Iteration: 200 / 1000 [ 20%] (Warmup)
## Chain 1: Iteration: 300 / 1000 [ 30%] (Warmup)
## Chain 1: Iteration: 400 / 1000 [ 40%] (Warmup)
## Chain 1: Iteration: 500 / 1000 [ 50%] (Warmup)
## Chain 1: Iteration: 501 / 1000 [ 50%] (Sampling)
## Chain 1: Iteration: 600 / 1000 [ 60%] (Sampling)
## Chain 1: Iteration: 700 / 1000 [ 70%] (Sampling)
## Chain 1: Iteration: 800 / 1000 [ 80%] (Sampling)
## Chain 1: Iteration: 900 / 1000 [ 90%] (Sampling)
## Chain 1: Iteration: 1000 / 1000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 6.7115 seconds (Warm-up)
## Chain 1: 7.14035 seconds (Sampling)
## Chain 1: 13.8518 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'survey_3' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 0.000929 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 9.29 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 1000 [ 0%] (Warmup)
## Chain 2: Iteration: 100 / 1000 [ 10%] (Warmup)
## Chain 2: Iteration: 200 / 1000 [ 20%] (Warmup)
## Chain 2: Iteration: 300 / 1000 [ 30%] (Warmup)
## Chain 2: Iteration: 400 / 1000 [ 40%] (Warmup)
## Chain 2: Iteration: 500 / 1000 [ 50%] (Warmup)
## Chain 2: Iteration: 501 / 1000 [ 50%] (Sampling)
## Chain 2: Iteration: 600 / 1000 [ 60%] (Sampling)
## Chain 2: Iteration: 700 / 1000 [ 70%] (Sampling)
## Chain 2: Iteration: 800 / 1000 [ 80%] (Sampling)
## Chain 2: Iteration: 900 / 1000 [ 90%] (Sampling)
## Chain 2: Iteration: 1000 / 1000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 6.90926 seconds (Warm-up)
## Chain 2: 7.30582 seconds (Sampling)
## Chain 2: 14.2151 seconds (Total)
## Chain 2:print(fit3, pars = "p", probs = c(0.025, 0.5, 0.975))## Inference for Stan model: survey_3.
## 2 chains, each with iter=1000; warmup=500; thin=1;
## post-warmup draws per chain=500, total post-warmup draws=1000.
##
## mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
## p[1] 0.81 0.01 0.19 0.28 0.88 0.99 1384 1.00
## p[2] 0.80 0.01 0.21 0.24 0.89 0.99 1304 1.00
## p[3] 0.99 0.00 0.02 0.92 1.00 1.00 712 1.00
## p[4] 0.97 0.00 0.06 0.79 0.99 1.00 691 1.00
## p[5] 0.99 0.00 0.02 0.93 1.00 1.00 533 1.01
## p[6] 0.09 0.00 0.12 0.00 0.04 0.45 772 1.00
## p[7] 0.92 0.00 0.12 0.58 0.97 1.00 747 1.00
## p[8] 0.94 0.00 0.10 0.64 0.98 1.00 895 1.00
## p[9] 0.97 0.00 0.06 0.82 0.99 1.00 644 1.00
## p[10] 0.91 0.00 0.12 0.51 0.97 1.00 965 1.00
## p[11] 0.99 0.00 0.02 0.96 1.00 1.00 727 1.00
## p[12] 0.93 0.00 0.10 0.67 0.97 1.00 998 1.00
## p[13] 0.99 0.00 0.03 0.92 1.00 1.00 946 1.00
## p[14] 0.99 0.00 0.04 0.89 1.00 1.00 777 1.00
## p[15] 0.99 0.00 0.02 0.95 1.00 1.00 766 1.00
## p[16] 0.88 0.00 0.15 0.40 0.94 1.00 921 1.00
## p[17] 0.93 0.00 0.10 0.61 0.97 1.00 672 1.00
## p[18] 0.92 0.00 0.13 0.48 0.97 1.00 683 1.00
## p[19] 0.47 0.01 0.29 0.03 0.45 0.95 2277 1.00
## p[20] 0.81 0.01 0.21 0.23 0.90 1.00 1200 1.00
## p[21] 0.93 0.00 0.10 0.64 0.97 1.00 1009 1.00
## p[22] 0.95 0.00 0.09 0.67 0.98 1.00 546 1.00
## p[23] 1.00 0.00 0.00 1.00 1.00 1.00 1059 1.00
## p[24] 0.78 0.01 0.21 0.22 0.85 0.99 1551 1.00
## p[25] 0.89 0.00 0.14 0.45 0.95 1.00 949 1.00
## p[26] 0.97 0.00 0.05 0.84 0.99 1.00 939 1.00
## p[27] 0.97 0.00 0.05 0.82 0.99 1.00 785 1.00
## p[28] 1.00 0.00 0.00 0.99 1.00 1.00 560 1.00
## p[29] 0.61 0.01 0.28 0.06 0.65 0.98 1711 1.00
## p[30] 0.93 0.00 0.11 0.58 0.97 1.00 765 1.00
## p[31] 0.76 0.01 0.22 0.18 0.84 0.99 1279 1.00
## p[32] 0.99 0.00 0.02 0.95 1.00 1.00 928 1.00
## p[33] 0.86 0.00 0.16 0.38 0.93 1.00 1181 1.00
## p[34] 0.69 0.01 0.26 0.10 0.77 0.99 1633 1.00
## p[35] 0.98 0.00 0.04 0.89 0.99 1.00 866 1.00
## p[36] 0.94 0.00 0.10 0.64 0.98 1.00 774 1.00
## p[37] 0.67 0.01 0.26 0.09 0.75 0.99 2215 1.00
## p[38] 0.92 0.00 0.12 0.56 0.96 1.00 669 1.00
## p[39] 0.81 0.01 0.21 0.22 0.90 0.99 1048 1.00
## p[40] 0.86 0.00 0.17 0.36 0.92 1.00 1253 1.00
## p[41] 0.80 0.01 0.21 0.24 0.88 0.99 1146 1.00
## p[42] 0.83 0.01 0.20 0.29 0.91 1.00 1219 1.00
## p[43] 0.95 0.00 0.09 0.67 0.98 1.00 865 1.00
## p[44] 0.95 0.00 0.08 0.73 0.98 1.00 827 1.00
## p[45] 0.70 0.01 0.26 0.13 0.78 0.99 1267 1.00
## p[46] 0.78 0.01 0.22 0.19 0.86 0.99 1089 1.00
## p[47] 0.74 0.01 0.25 0.15 0.84 0.99 1066 1.00
## p[48] 1.00 0.00 0.02 0.97 1.00 1.00 776 1.00
## p[49] 0.96 0.00 0.08 0.73 0.98 1.00 883 1.00
## p[50] 0.90 0.00 0.13 0.49 0.95 1.00 857 1.00
## p[51] 0.52 0.01 0.28 0.05 0.53 0.96 1818 1.00
## p[52] 0.98 0.00 0.04 0.88 0.99 1.00 953 1.00
## p[53] 0.97 0.00 0.06 0.80 0.99 1.00 624 1.00
## p[54] 0.81 0.01 0.21 0.26 0.89 0.99 1142 1.00
## p[55] 0.72 0.01 0.24 0.15 0.79 0.99 1435 1.00
## p[56] 0.23 0.01 0.23 0.01 0.15 0.81 1349 1.00
## p[57] 0.21 0.01 0.22 0.01 0.13 0.80 1100 1.00
## p[58] 0.01 0.00 0.03 0.00 0.00 0.07 882 1.00
## p[59] 0.21 0.01 0.21 0.01 0.14 0.76 1178 1.00
## p[60] 0.01 0.00 0.02 0.00 0.00 0.04 878 1.00
## p[61] 0.59 0.01 0.27 0.06 0.62 0.97 2193 1.00
## p[62] 0.13 0.01 0.16 0.00 0.07 0.68 951 1.00
## p[63] 0.04 0.00 0.07 0.00 0.01 0.28 864 1.00
## p[64] 0.10 0.01 0.14 0.00 0.05 0.59 718 1.00
## p[65] 0.03 0.00 0.06 0.00 0.01 0.20 722 1.00
## p[66] 0.10 0.00 0.13 0.00 0.05 0.55 808 1.00
## p[67] 0.02 0.00 0.03 0.00 0.01 0.10 798 1.00
## p[68] 0.25 0.01 0.24 0.01 0.16 0.85 991 1.00
## p[69] 0.61 0.01 0.26 0.08 0.66 0.98 1928 1.00
## p[70] 0.01 0.00 0.03 0.00 0.00 0.10 830 1.00
## p[71] 0.32 0.01 0.24 0.02 0.27 0.85 1802 1.00
## p[72] 0.11 0.00 0.15 0.00 0.05 0.55 1113 1.00
## p[73] 0.01 0.00 0.02 0.00 0.00 0.04 630 1.00
## p[74] 0.01 0.00 0.02 0.00 0.00 0.06 811 1.00
## p[75] 0.02 0.00 0.05 0.00 0.01 0.14 341 1.00
## p[76] 0.57 0.01 0.27 0.07 0.62 0.97 1935 1.00
## p[77] 0.61 0.01 0.27 0.09 0.65 0.98 2214 1.00
## p[78] 0.28 0.01 0.24 0.01 0.20 0.85 1584 1.00
## p[79] 0.00 0.00 0.01 0.00 0.00 0.02 572 1.01
## p[80] 0.15 0.01 0.18 0.00 0.08 0.68 956 1.00
## p[81] 0.00 0.00 0.01 0.00 0.00 0.02 885 1.00
## p[82] 0.04 0.00 0.08 0.00 0.01 0.27 845 1.00
## p[83] 0.01 0.00 0.03 0.00 0.00 0.10 833 1.00
## p[84] 0.30 0.01 0.25 0.01 0.23 0.84 1714 1.00
## p[85] 0.03 0.00 0.06 0.00 0.01 0.18 692 1.00
## p[86] 0.22 0.01 0.22 0.01 0.15 0.80 867 1.00
## p[87] 0.00 0.00 0.01 0.00 0.00 0.03 636 1.00
## p[88] 0.02 0.00 0.04 0.00 0.00 0.10 627 1.01
## p[89] 0.18 0.01 0.20 0.00 0.11 0.75 897 1.00
## p[90] 0.49 0.01 0.28 0.04 0.49 0.96 1873 1.00
## p[91] 0.02 0.00 0.03 0.00 0.01 0.10 725 1.00
## p[92] 0.01 0.00 0.03 0.00 0.00 0.08 575 1.00
## p[93] 0.26 0.01 0.23 0.01 0.19 0.83 1851 1.00
## p[94] 0.00 0.00 0.01 0.00 0.00 0.03 882 1.00
## p[95] 0.12 0.01 0.16 0.00 0.05 0.61 1016 1.00
## p[96] 0.43 0.01 0.29 0.02 0.39 0.95 2248 1.00
## p[97] 0.03 0.00 0.05 0.00 0.01 0.18 805 1.00
## p[98] 0.11 0.01 0.16 0.00 0.04 0.60 706 1.00
## p[99] 0.04 0.00 0.07 0.00 0.01 0.22 913 1.00
## p[100] 0.28 0.01 0.24 0.01 0.21 0.86 1034 1.00
## p[101] 0.22 0.01 0.22 0.01 0.14 0.77 1110 1.00
## p[102] 0.03 0.00 0.06 0.00 0.01 0.21 948 1.00
## p[103] 0.10 0.00 0.13 0.00 0.05 0.50 1176 1.00
## p[104] 0.10 0.00 0.14 0.00 0.05 0.56 877 1.00
## p[105] 0.01 0.00 0.02 0.00 0.00 0.07 845 1.00
## p[106] 0.39 0.01 0.27 0.02 0.34 0.91 1704 1.00
## p[107] 0.00 0.00 0.01 0.00 0.00 0.02 769 1.00
## p[108] 0.11 0.01 0.15 0.00 0.05 0.58 627 1.00
## p[109] 0.22 0.01 0.23 0.01 0.13 0.84 1470 1.00
## p[110] 0.00 0.00 0.00 0.00 0.00 0.01 522 1.00
##
## Samples were drawn using NUTS(diag_e) at Fri Jan 22 07:37:34 2021.
## For each parameter, 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).Ketika Q menjadi besar, kita dapat membuat p ini miring, berubah menjadi nol atau satu. Pada titik ini saya ingat bahwa mengalikan probabilitas bersama-sama dapat menyebabkan ledakan karena segala sesuatunya menjadi kecil dengan sangat cepat, dan komputer tidak menyukai angka kecil.
Kita bisa melakukan trik standar dan mengambil jumlah di ruang log, daripada mengalikan di ruang asli untuk memperbaikinya. Pada saat yang sama, mari tambahkan beberapa fungsi prediksi. Ada beberapa cara untuk melakukan ini saat ini.
Saya memutuskan untuk meneruskan data yang tidak berlabel ke Stan, dan menulis fungsi untuk mencegah kode duplikat. Fungsinya berjalan di bloknya sendiri. Butuh beberapa saat bagi saya untuk mencari tahu cara mendapatkan tanda tangan fungsi dengan benar, dan cara kerja pelingkupan variabel, tetapi saya akhirnya harus:
Sekarang kita perlu membuat set pengujian untuk diprediksi, dan melakukan prediksi.
new_df <- sample_n(theta_0, 55, id = 0) %>%
bind_rows(sample_n(theta_1, 55, id = 1))
new_data <- list(
R = as.matrix(select(df, -y)),
N = nrow(df),
y = df$y,
new_R = as.matrix(select(new_df, -y)),
new_N = nrow(new_df),
Q = Q,
A = A,
alpha = alpha
)model3 <- stan_model("survey_4.stan")fit4 <- sampling(model3, data = new_data, chains = 2, iter = 1000, refresh = 0)
print(fit4, pars = "new_pred", probs = c(0.025, 0.5, 0.975))## Inference for Stan model: survey_4.
## 2 chains, each with iter=1000; warmup=500; thin=1;
## post-warmup draws per chain=500, total post-warmup draws=1000.
##
## mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
## new_pred[1] 0.89 0.01 0.16 0.40 0.95 1.00 656 1.00
## new_pred[2] 0.74 0.01 0.24 0.17 0.82 0.99 1404 1.00
## new_pred[3] 0.93 0.00 0.11 0.59 0.97 1.00 717 1.00
## new_pred[4] 0.80 0.01 0.21 0.23 0.87 0.99 1241 1.00
## new_pred[5] 0.87 0.00 0.16 0.37 0.94 1.00 1233 1.00
## new_pred[6] 0.93 0.00 0.11 0.53 0.97 1.00 685 1.00
## new_pred[7] 0.89 0.00 0.14 0.49 0.95 1.00 1020 1.00
## new_pred[8] 0.57 0.01 0.29 0.04 0.60 0.98 2376 1.00
## new_pred[9] 0.62 0.01 0.28 0.06 0.69 0.98 2306 1.00
## new_pred[10] 0.72 0.01 0.24 0.13 0.79 0.99 1286 1.00
## new_pred[11] 0.44 0.01 0.29 0.02 0.41 0.95 2009 1.00
## new_pred[12] 0.79 0.01 0.22 0.19 0.86 0.99 1046 1.00
## new_pred[13] 0.12 0.00 0.16 0.00 0.05 0.62 1014 1.00
## new_pred[14] 0.85 0.00 0.17 0.34 0.92 1.00 1432 1.00
## new_pred[15] 0.58 0.01 0.27 0.07 0.61 0.97 1991 1.00
## new_pred[16] 0.93 0.00 0.11 0.58 0.98 1.00 592 1.00
## new_pred[17] 0.87 0.00 0.16 0.37 0.94 1.00 1029 1.00
## new_pred[18] 0.32 0.01 0.26 0.01 0.25 0.92 1502 1.00
## new_pred[19] 0.95 0.00 0.09 0.70 0.98 1.00 695 1.00
## new_pred[20] 0.74 0.01 0.24 0.17 0.82 0.99 1818 1.00
## new_pred[21] 0.55 0.01 0.28 0.06 0.56 0.97 1992 1.00
## new_pred[22] 0.99 0.00 0.02 0.96 1.00 1.00 689 1.00
## new_pred[23] 0.75 0.01 0.23 0.17 0.84 0.99 1383 1.00
## new_pred[24] 0.95 0.00 0.08 0.71 0.98 1.00 983 1.00
## new_pred[25] 0.98 0.00 0.05 0.87 0.99 1.00 1027 1.00
## new_pred[26] 0.88 0.00 0.16 0.41 0.94 1.00 1430 1.00
## new_pred[27] 0.91 0.00 0.13 0.49 0.96 1.00 1055 1.00
## new_pred[28] 0.94 0.00 0.10 0.65 0.97 1.00 793 1.00
## new_pred[29] 0.37 0.01 0.27 0.02 0.31 0.92 2132 1.00
## new_pred[30] 0.81 0.01 0.20 0.26 0.89 1.00 846 1.00
## new_pred[31] 0.84 0.01 0.19 0.29 0.92 1.00 1192 1.00
## new_pred[32] 0.64 0.01 0.28 0.07 0.70 0.99 2039 1.00
## new_pred[33] 0.96 0.00 0.06 0.77 0.98 1.00 896 1.00
## new_pred[34] 1.00 0.00 0.01 0.98 1.00 1.00 869 1.00
## new_pred[35] 0.94 0.00 0.10 0.64 0.98 1.00 714 1.00
## new_pred[36] 0.48 0.01 0.29 0.04 0.46 0.96 2642 1.00
## new_pred[37] 0.41 0.01 0.29 0.01 0.36 0.95 2539 1.00
## new_pred[38] 0.97 0.00 0.06 0.81 0.99 1.00 915 1.00
## new_pred[39] 0.73 0.01 0.24 0.17 0.82 0.99 1559 1.00
## new_pred[40] 0.94 0.00 0.09 0.69 0.98 1.00 751 1.01
## new_pred[41] 0.79 0.01 0.22 0.18 0.87 0.99 1553 1.00
## new_pred[42] 0.88 0.00 0.15 0.42 0.94 1.00 939 1.00
## new_pred[43] 0.95 0.00 0.09 0.68 0.98 1.00 613 1.00
## new_pred[44] 0.88 0.00 0.15 0.48 0.94 1.00 1092 1.00
## new_pred[45] 0.96 0.00 0.07 0.75 0.98 1.00 795 1.00
## new_pred[46] 0.59 0.01 0.27 0.08 0.64 0.97 2146 1.00
## new_pred[47] 0.90 0.00 0.14 0.48 0.96 1.00 951 1.00
## new_pred[48] 0.91 0.00 0.13 0.52 0.96 1.00 889 1.01
## new_pred[49] 0.27 0.01 0.23 0.01 0.18 0.81 1808 1.00
## new_pred[50] 0.94 0.00 0.09 0.63 0.98 1.00 726 1.00
## new_pred[51] 0.79 0.01 0.22 0.19 0.88 1.00 1288 1.00
## new_pred[52] 0.95 0.00 0.09 0.65 0.98 1.00 833 1.00
## new_pred[53] 0.78 0.01 0.22 0.22 0.87 0.99 1329 1.00
## new_pred[54] 0.67 0.01 0.27 0.11 0.76 0.99 1499 1.00
## new_pred[55] 0.98 0.00 0.04 0.90 0.99 1.00 922 1.00
## new_pred[56] 0.12 0.00 0.16 0.00 0.05 0.58 1161 1.00
## new_pred[57] 0.25 0.01 0.23 0.01 0.17 0.84 1269 1.00
## new_pred[58] 0.49 0.01 0.28 0.04 0.50 0.95 3000 1.00
## new_pred[59] 0.32 0.01 0.26 0.01 0.26 0.91 1899 1.00
## new_pred[60] 0.07 0.00 0.12 0.00 0.03 0.45 702 1.00
## new_pred[61] 0.03 0.00 0.06 0.00 0.01 0.18 707 1.00
## new_pred[62] 0.58 0.01 0.29 0.04 0.63 0.98 1960 1.00
## new_pred[63] 0.18 0.01 0.20 0.00 0.10 0.76 1275 1.00
## new_pred[64] 0.05 0.00 0.09 0.00 0.02 0.35 879 1.00
## new_pred[65] 0.46 0.01 0.29 0.04 0.41 0.94 2097 1.00
## new_pred[66] 0.48 0.01 0.29 0.03 0.47 0.96 2226 1.00
## new_pred[67] 0.18 0.01 0.19 0.01 0.11 0.71 859 1.00
## new_pred[68] 0.26 0.01 0.24 0.01 0.18 0.84 1668 1.00
## new_pred[69] 0.14 0.01 0.17 0.00 0.08 0.68 995 1.00
## new_pred[70] 0.01 0.00 0.01 0.00 0.00 0.05 688 1.00
## new_pred[71] 0.28 0.01 0.25 0.01 0.20 0.84 1739 1.00
## new_pred[72] 0.67 0.01 0.26 0.11 0.74 0.99 1532 1.00
## new_pred[73] 0.07 0.00 0.11 0.00 0.03 0.37 956 1.00
## new_pred[74] 0.03 0.00 0.05 0.00 0.01 0.15 781 1.00
## new_pred[75] 0.09 0.00 0.14 0.00 0.04 0.54 984 1.00
## new_pred[76] 0.07 0.00 0.11 0.00 0.03 0.38 938 1.00
## new_pred[77] 0.29 0.01 0.26 0.01 0.20 0.88 1772 1.00
## new_pred[78] 0.23 0.01 0.23 0.01 0.14 0.83 1231 1.00
## new_pred[79] 0.18 0.01 0.19 0.01 0.11 0.69 1230 1.00
## new_pred[80] 0.04 0.00 0.08 0.00 0.02 0.26 671 1.00
## new_pred[81] 0.18 0.01 0.20 0.00 0.10 0.72 1435 1.00
## new_pred[82] 0.29 0.01 0.24 0.01 0.21 0.86 1935 1.00
## new_pred[83] 0.04 0.00 0.07 0.00 0.01 0.25 765 1.00
## new_pred[84] 0.02 0.00 0.04 0.00 0.01 0.15 667 1.00
## new_pred[85] 0.41 0.01 0.27 0.02 0.37 0.93 2458 1.00
## new_pred[86] 0.01 0.00 0.02 0.00 0.00 0.05 896 1.00
## new_pred[87] 0.32 0.01 0.26 0.01 0.25 0.91 1917 1.00
## new_pred[88] 0.02 0.00 0.04 0.00 0.01 0.14 830 1.00
## new_pred[89] 0.27 0.01 0.25 0.01 0.18 0.87 1025 1.00
## new_pred[90] 0.36 0.01 0.26 0.02 0.31 0.89 1918 1.00
## new_pred[91] 0.66 0.01 0.26 0.10 0.73 0.98 2283 1.00
## new_pred[92] 0.08 0.00 0.12 0.00 0.03 0.47 861 1.00
## new_pred[93] 0.47 0.01 0.28 0.04 0.46 0.95 2383 1.00
## new_pred[94] 0.27 0.01 0.23 0.01 0.20 0.83 1708 1.00
## new_pred[95] 0.67 0.01 0.26 0.10 0.73 0.99 1886 1.00
## new_pred[96] 0.02 0.00 0.04 0.00 0.01 0.11 796 1.00
## new_pred[97] 0.78 0.01 0.22 0.21 0.86 0.99 1655 1.00
## new_pred[98] 0.13 0.01 0.17 0.00 0.07 0.68 927 1.00
## new_pred[99] 0.06 0.00 0.11 0.00 0.02 0.42 631 1.00
## new_pred[100] 0.35 0.01 0.29 0.01 0.28 0.93 2127 1.00
## new_pred[101] 0.52 0.01 0.29 0.03 0.52 0.97 2494 1.00
## new_pred[102] 0.01 0.00 0.03 0.00 0.00 0.10 630 1.00
## new_pred[103] 0.48 0.01 0.28 0.04 0.47 0.96 2059 1.00
## new_pred[104] 0.23 0.01 0.23 0.01 0.15 0.82 1096 1.00
## new_pred[105] 0.59 0.01 0.27 0.06 0.63 0.98 2040 1.00
## new_pred[106] 0.09 0.00 0.12 0.00 0.04 0.45 865 1.00
## new_pred[107] 0.25 0.01 0.24 0.01 0.17 0.84 1742 1.00
## new_pred[108] 0.18 0.01 0.20 0.00 0.10 0.75 1352 1.00
## new_pred[109] 0.34 0.01 0.25 0.02 0.28 0.87 2217 1.00
## new_pred[110] 0.26 0.01 0.25 0.01 0.17 0.85 1599 1.00
##
## Samples were drawn using NUTS(diag_e) at Fri Jan 22 07:38:12 2021.
## For each parameter, 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).Sekali lagi, Anda dapat mencoba shinystan :: launch_shinystan (fit4) untuk memeriksa diagnostik MCMC. Kami dapat dengan cepat menghitung akurasi pada set pengujian menggunakan ambang prediksi 0,5:
shinystan::launch_shinystan(fit4)##
## Launching ShinyStan interface... for large models this may take some time.##
## Listening on http://127.0.0.1:4901library(tidybayes)
# here `tidybayes::spread_draws()` and `tidybayes::median_qi()`
# save an immense amount of headache
pred <- fit4 %>%
spread_draws(new_pred[i]) %>%
median_qi(new_pred[i]) %>%
mutate(.pred = if_else(`new_pred[i]` > 0.5, 0, 1))
acc <- round(mean(pred$.pred == new_df$y), 2)
cat("The classification accuracy is:", acc)## The classification accuracy is: 0.51Awalnya saya berencana untuk membandingkan dengan LASSO, karena model ini terinspirasi oleh tweet yang menggunakan LASSO, tetapi saya sedikit kehabisan tenaga, jadi biarkan saja untuk saat ini (meskipun menurut saya ini akan keren untuk mencoba ExclusiveLasso untuk memilih paling banyak satu jawaban untuk setiap pertanyaan sebagai penting).
Akhirnya, kita juga dapat melihat bagaimana grup 0 dan grup 1 berbeda dengan memplot theta_0 dan theta_1. Kami munge sedikit dulu
theta_0_draws <- fit4 %>%
spread_samples(theta_0[i, j])## Warning: 'spread_samples' is deprecated.
## Use 'spread_draws' instead.
## See help("Deprecated") and help("tidybayes-deprecated").theta_1_draws <- fit4 %>%
spread_samples(theta_1[i, j])## Warning: 'spread_samples' is deprecated.
## Use 'spread_draws' instead.
## See help("Deprecated") and help("tidybayes-deprecated").theta_draws <- theta_0_draws %>%
left_join(theta_1_draws)## Joining, by = c("i", "j", ".chain", ".iteration", ".draw")dan kemudian dapat memvisualisasikan probabilitas respons dengan pertanyaan dan kelas
library(ggplot2)
theta_draws %>%
gather(group, theta, theta_0, theta_1) %>%
mutate(
group = if_else(group == "theta_0", "Group 0", "Group 1"),
question = i,
response = j
) %>%
ggplot(aes(theta, fill = group, color = group)) +
geom_density(alpha = 0.5) +
facet_grid(
rows = vars(question),
cols = vars(response)
) +
labs(
title = "Estimated probability of each response by question and group",
subtitle = "Columns correspond to response, rows to questions",
x = "Probability of response") +
theme(
legend.position = "none",
axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
Daftar pustaka :