UPDATE: 2024-01-05 14:14:28.09333
このノートは「StanとRでベイズ統計モデリング」の内容を写経することで、ベイズ統計への理解を深めていくために作成している。
基本的には気になった部分を写経しながら、ところどころ自分用の補足をメモすることで、「StanとRでベイズ統計モデリング」を読み進めるための自分用の補足資料になることを目指す。私の解釈がおかしく、メモが誤っている場合があるので注意。
今回は第11章「離散値をとるパラメータを使う」の後半から写経していく。
ある飲食店への来店回数を取得したデータをもとに来店回数が多い人の特徴を調べたい。来店回数はポアソン分布を仮定できそうではあるが、0が多いことから、ゼロ過剰ポアソン分布を利用するほうが望ましい。
Sex
: 0が男性、1が女性Sake
: 0が飲まない、1が飲むlibrary(dplyr)
library(rstan)
library(ggplot2)
library(GGally)
options(max.print = 999999)
rstan_options(auto_write=TRUE)
options(mc.cores=parallel::detectCores())
d <- read.csv('https://raw.githubusercontent.com/MatsuuraKentaro/RStanBook/master/chap11/input/data-ZIP.txt')
head(d)
## Sex Sake Age Y
## 1 0 1 18 5
## 2 1 0 18 2
## 3 1 1 18 1
## 4 0 0 19 3
## 5 0 0 19 5
## 6 1 0 19 7
解析の目的は「リピーターになりそうな人を知りたい」とする。また「説明変数でリピーターになるかをどれほど予測できるかを知りたい、影響度合いも知りたい」とする。
d$Sex <- as.factor(d$Sex)
d$Sake <- as.factor(d$Sake)
N_col <- ncol(d)
ggp <- ggpairs(d, upper = 'blank', diag = 'blank', lower = 'blank')
# 対角成分のヒストグラムを作成
for (i in 1:N_col) {
x <- d[,i]
p <- ggplot(data.frame(x), aes(x = x)) +
theme_bw(base_size = 14) +
theme(axis.text.x = element_text(angle = 40, vjust = 1, hjust = 1))
if (class(x) == 'factor') {
p <- p + geom_bar(color = 'grey5', fill = 'grey80')
} else {
bw <- ifelse(colnames(d)[i] == 'Y', 1, (max(x)-min(x))/10)
p <- p + geom_histogram(binwidth = bw, color = 'grey5', fill = 'grey80') +
geom_line(aes(y = after_stat(count)*bw), stat = 'density')
}
p <- p + geom_label(data = data.frame(x = -Inf, y = Inf, label = colnames(d)[i]), aes(x = x, y = y, label = label), hjust = 0, vjust = 1)
ggp <- putPlot(ggp, p, i, i)
}
# 上三角成分の相関係数を作成
zcolat <- seq(-1, 1, length = 81)
zcolre <- c(zcolat[1:40]+1, rev(zcolat[41:81]))
for (i in 1:(N_col-1)) {
for (j in (i+1):N_col) {
x <- as.numeric(d[,i])
y <- as.numeric(d[,j])
r <- cor(x, y, method = 'spearman', use = 'pairwise.complete.obs')
zcol <- lattice::level.colors(r, at = zcolat, col.regions = grey(zcolre))
textcol <- ifelse(abs(r) < 0.4, 'grey20', 'white')
ell <- ellipse::ellipse(r, level = 0.95, type = 'l', npoints = 50, scale = c(.2, .2), centre = c(.5, .5))
p <- ggplot(data.frame(ell), aes(x = x, y = y)) + theme_bw() + theme(
plot.background = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.border = element_blank(), axis.ticks = element_blank()) +
geom_polygon(fill = zcol, color = zcol) +
geom_text(data = NULL, x = .5, y = .5, label = 100*round(r, 2), size = 6, col = textcol)
ggp <- putPlot(ggp, p, i, j)
}
}
# 下三角成分の箱ヒゲ図を作成
for (j in 1:(N_col-1)) {
for (i in (j+1):N_col) {
x <- d[,j]
y <- d[,i]
if (class(x) == 'factor' && class(y) == 'factor') {
p <- ggplot(reshape2::melt(table(x,y)), aes(x = x, y = y)) +
theme_bw(base_size = 14) +
geom_point(aes(size = value), color = 'grey80') +
geom_text(aes(label = value)) +
scale_size_area(max_size = 8) +
scale_x_continuous(breaks = 0:1, limits = c(-0.5,1.5)) + scale_y_continuous(breaks = 0:1, limits = c(-0.5,1.5))
} else {
p <- ggplot(data.frame(x, y, Y = as.numeric(d$Y)), aes(x = x, y = y, color = Y)) +
theme_bw(base_size = 14) +
theme(axis.text.x = element_text(angle = 40, vjust = 1, hjust = 1))
if (class(x) == 'factor') {
p <- p + geom_boxplot(alpha = 3/6, outlier.size = 0, fill = 'white') +
geom_point(position = position_jitter(w = 0.4, h = 0), size = 1)
} else {
p <- p + geom_point(size = 1)
}
p <- p + scale_color_gradient(low = 'grey65', high = 'grey5')
}
ggp <- putPlot(ggp, p, i, j)
}
}
print(ggp, left = 0.3, bottom = 0.3)
来店回数\(Y\)の0が突出している点が特徴的なので、このメカニズムを考える。1人が来店するメカニズムとして、「とにかく1回来店する」場合と「気に入って複数回来店する」場合のパターンがある。そして、とにかく1回来店するのは確率\(q\)のベルヌイ分布に従うと考え、リピーターの来店回数は平均\(\lambda\)のポアソン分布に従うと考える。
つまり、確率\(q\)のコインを投げて、裏が出れば来店回数は0、表が出れば平均\(\lambda\)のポアソン分布に従って生成される。そして、\(q\)はロジスティック回帰、\(\lambda\)はポアソン回帰でパラメタの推定を分ける。コイン投げは離散パラメタとなるので、周辺化消去が必要になる。
\(y=0\)のときは「ベルヌイ分布で0が出る確率」と「ベルヌイ分布で1が出る確率」と「ポアソン分布で0が出る確率」の和となる。\(y=0\)のときはベルヌイ分布の可能性が高いが、\(y=1\)のときはどちらの分布か判断し難くく、\(y\ge2\)のときはポアソン分布に従うと考える。そのため、\(Bernoulli(0|q) + Bernoulli(1|q) × Poisson(y=0|\lambda)\)とする。
\[ \begin{eqnarray} ZIP(y|q,\lambda) = \begin{cases} Bernoulli(0|q) + Bernoulli(1|q) × Poisson(y=0|\lambda) \ if \ y = 0 \\ Bernoulli(1|q) × Poisson(y|\lambda) \ if \ y \ge 1 \\ \end{cases} \end{eqnarray} \]
\(\boldsymbol{ X }\)はN×4の行列で、\(\overrightarrow{ b}\)は長さ4のベクトル。
\[ \begin{eqnarray} q[n] &=& inv\_logit((\boldsymbol{ X } \overrightarrow{ b_{1}})[n]) \\ \lambda[n] &=& (\boldsymbol{ X } \overrightarrow{ b_{2}})[n] \\ Y[n] &\sim& ZIP(q[n], \lambda[n]) \end{eqnarray} \]
Stanのモデルは下記の通り。
functions {
real ZIP_lpmf(int Y, real q, real lambda) {
if (Y == 0) {
return log_sum_exp(
bernoulli_lpmf(0 | q),
bernoulli_lpmf(1 | q) + poisson_log_lpmf(0 | lambda)
);
} else {
return bernoulli_lpmf(1 | q) + poisson_log_lpmf(Y | lambda);
}
}
}
data {
int N;
int D;
int<lower=0> Y[N];
matrix[N,D] X;
}
parameters {
vector[D] b[2];
}
transformed parameters {
vector[N] q_x;
vector[N] q;
vector[N] lambda;
q_x = X*b[1];
lambda = X*b[2];
for (n in 1:N)
q[n] = inv_logit(q_x[n]);
}
model {
for (n in 1:N)
Y[n] ~ ZIP(q[n], lambda[n]);
}
データを用意する。
d$Age <- d$Age/10
X <- cbind(1, d[,-ncol(d)])
data <- list(N = nrow(d), D = ncol(X), Y = d$Y, X = X)
data
## $N
## [1] 200
##
## $D
## [1] 4
##
## $Y
## [1] 5 2 1 3 5 7 5 0 8 8 0 7 6 0 0 1 4 0 0 3 0 0 4 8 7
## [26] 3 7 0 2 4 0 9 0 7 0 0 9 4 6 1 13 5 7 1 6 2 10 0 2 1
## [51] 6 9 0 9 8 3 11 4 0 6 3 1 4 3 0 0 5 0 8 7 2 4 0 7 6
## [76] 5 0 0 0 4 2 11 10 5 0 0 0 0 0 2 7 11 0 3 6 4 1 0 0 5
## [101] 0 0 6 0 5 9 0 0 10 0 2 0 6 10 6 8 5 8 9 3 6 5 9 0 0
## [126] 11 3 0 5 6 6 7 9 2 0 12 0 0 5 4 0 0 7 0 9 6 0 4 8 0
## [151] 0 4 0 8 6 3 8 11 0 7 0 15 0 0 0 9 0 12 0 9 7 9 0 15 9
## [176] 6 3 0 0 0 12 2 9 5 17 6 0 0 10 15 0 0 5 4 13 6 0 0 12 6
##
## $X
## 1 Sex Sake Age
## 1 1 0 1 1.8
## 2 1 1 0 1.8
## 3 1 1 1 1.8
## 4 1 0 0 1.9
## 5 1 0 0 1.9
## 6 1 1 0 1.9
## 7 1 0 1 2.0
## 8 1 0 0 2.0
## 9 1 0 0 2.1
## 10 1 0 0 2.1
## 11 1 0 0 2.1
## 12 1 0 1 2.1
## 13 1 0 0 2.1
## 14 1 1 0 2.2
## 15 1 0 0 2.2
## 16 1 1 0 2.2
## 17 1 0 1 2.2
## 18 1 0 0 2.2
## 19 1 1 0 2.3
## 20 1 1 0 2.3
## 21 1 0 0 2.3
## 22 1 1 0 2.3
## 23 1 1 0 2.3
## 24 1 1 0 2.3
## 25 1 0 1 2.3
## 26 1 0 0 2.3
## 27 1 0 0 2.3
## 28 1 0 0 2.4
## 29 1 1 0 2.4
## 30 1 1 1 2.4
## 31 1 1 1 2.5
## 32 1 0 0 2.5
## 33 1 0 0 2.6
## 34 1 0 1 2.6
## 35 1 0 0 2.6
## 36 1 0 0 2.6
## 37 1 0 0 2.6
## 38 1 1 0 2.6
## 39 1 0 0 2.6
## 40 1 1 0 2.6
## 41 1 0 0 2.6
## 42 1 1 0 2.7
## 43 1 0 1 2.7
## 44 1 1 0 2.7
## 45 1 0 1 2.7
## 46 1 1 1 2.7
## 47 1 0 1 2.8
## 48 1 1 0 2.8
## 49 1 1 0 2.8
## 50 1 1 0 2.8
## 51 1 0 0 2.8
## 52 1 0 0 2.8
## 53 1 0 0 2.8
## 54 1 0 0 2.8
## 55 1 0 0 2.9
## 56 1 1 1 2.9
## 57 1 0 1 2.9
## 58 1 0 1 2.9
## 59 1 0 0 3.0
## 60 1 1 0 3.0
## 61 1 1 1 3.0
## 62 1 1 0 3.0
## 63 1 1 0 3.0
## 64 1 1 0 3.0
## 65 1 0 0 3.1
## 66 1 0 0 3.1
## 67 1 0 0 3.1
## 68 1 0 0 3.1
## 69 1 0 0 3.2
## 70 1 0 1 3.2
## 71 1 1 0 3.2
## 72 1 0 1 3.2
## 73 1 0 0 3.2
## 74 1 0 1 3.2
## 75 1 1 0 3.2
## 76 1 1 0 3.2
## 77 1 1 0 3.2
## 78 1 0 0 3.3
## 79 1 0 0 3.3
## 80 1 1 1 3.3
## 81 1 1 0 3.3
## 82 1 0 1 3.3
## 83 1 0 0 3.3
## 84 1 0 0 3.3
## 85 1 0 1 3.3
## 86 1 0 0 3.3
## 87 1 0 0 3.4
## 88 1 0 0 3.4
## 89 1 0 0 3.4
## 90 1 1 1 3.4
## 91 1 0 1 3.4
## 92 1 0 0 3.5
## 93 1 0 0 3.5
## 94 1 1 1 3.5
## 95 1 1 0 3.5
## 96 1 1 1 3.5
## 97 1 1 0 3.5
## 98 1 0 0 3.6
## 99 1 0 0 3.6
## 100 1 1 1 3.6
## 101 1 1 0 3.6
## 102 1 0 0 3.6
## 103 1 1 0 3.6
## 104 1 1 0 3.7
## 105 1 1 0 3.7
## 106 1 0 0 3.7
## 107 1 0 0 3.8
## 108 1 0 0 3.8
## 109 1 0 1 3.8
## 110 1 0 0 3.8
## 111 1 0 0 3.8
## 112 1 0 0 3.8
## 113 1 1 1 3.9
## 114 1 0 0 3.9
## 115 1 0 1 3.9
## 116 1 1 0 3.9
## 117 1 1 0 3.9
## 118 1 0 1 4.0
## 119 1 0 0 4.0
## 120 1 1 1 4.0
## 121 1 1 0 4.0
## 122 1 0 0 4.0
## 123 1 0 0 4.1
## 124 1 1 0 4.1
## 125 1 0 0 4.1
## 126 1 0 0 4.1
## 127 1 1 1 4.1
## 128 1 0 0 4.1
## 129 1 0 1 4.1
## 130 1 1 0 4.2
## 131 1 1 1 4.2
## 132 1 0 0 4.2
## 133 1 1 0 4.2
## 134 1 1 1 4.2
## 135 1 0 0 4.2
## 136 1 0 1 4.3
## 137 1 1 0 4.3
## 138 1 1 1 4.3
## 139 1 1 0 4.3
## 140 1 1 0 4.3
## 141 1 0 1 4.3
## 142 1 0 0 4.3
## 143 1 0 1 4.4
## 144 1 0 0 4.4
## 145 1 0 1 4.4
## 146 1 0 1 4.4
## 147 1 1 0 4.4
## 148 1 0 1 4.4
## 149 1 0 1 4.5
## 150 1 1 0 4.5
## 151 1 0 0 4.6
## 152 1 1 0 4.6
## 153 1 0 0 4.6
## 154 1 0 1 4.6
## 155 1 1 0 4.6
## 156 1 1 0 4.6
## 157 1 1 0 4.7
## 158 1 0 0 4.7
## 159 1 0 0 4.7
## 160 1 1 0 4.7
## 161 1 0 0 4.8
## 162 1 0 1 4.8
## 163 1 1 1 4.8
## 164 1 0 0 4.8
## 165 1 0 0 4.8
## 166 1 1 0 4.8
## 167 1 0 0 4.9
## 168 1 0 1 4.9
## 169 1 0 0 5.0
## 170 1 0 1 5.0
## 171 1 0 0 5.0
## 172 1 0 1 5.0
## 173 1 1 0 5.1
## 174 1 0 0 5.1
## 175 1 0 1 5.1
## 176 1 1 1 5.1
## 177 1 1 0 5.1
## 178 1 0 0 5.1
## 179 1 0 0 5.1
## 180 1 1 0 5.1
## 181 1 0 1 5.2
## 182 1 1 1 5.2
## 183 1 0 1 5.2
## 184 1 1 0 5.2
## 185 1 0 1 5.2
## 186 1 1 0 5.2
## 187 1 0 0 5.3
## 188 1 0 0 5.3
## 189 1 0 1 5.3
## 190 1 0 0 5.4
## 191 1 1 0 5.4
## 192 1 0 0 5.4
## 193 1 1 1 5.4
## 194 1 1 0 5.4
## 195 1 0 0 5.4
## 196 1 0 1 5.5
## 197 1 0 0 5.5
## 198 1 0 0 5.5
## 199 1 0 1 5.5
## 200 1 1 0 5.5
ここでは、stan_model()
関数で最初にコンパイルしておいてから、
model117 <- stan_model('note_ahirubayes15-117.stan')
sampling()
関数でサンプリングする。
fit <- sampling(object = model117, data = data, seed = 1989)
推定結果は下記の通り。b[1,]
は来店確率に関わる回帰係数で、b[2,]
はリピーターが何回来店するかに関わる回帰係数。b[,1−4]
は切片、Sex
、Sake
、Age
に対応する。Age
は10で割っているので、10倍する必要がある。1回来る来店確率が高いのは、「女性(b[1,2]=1.6
)でお酒を飲んでいる人(b[1,3]=3.4
)」である。リピーターになる確率が高いのは、「男性(b[2,2]=-0.7
)でお酒を飲んでいない人(b[2,3]=-0.2
)で、年齢が高い人(b[2,4]=0.2
)」である。
print(fit, prob = c(0.025, 0.5, 0.975), pars = c('b', 'q', 'lambda'), digits_summary = 1)
## Inference for Stan model: anon_model.
## 4 chains, each with iter=2000; warmup=1000; thin=1;
## post-warmup draws per chain=1000, total post-warmup draws=4000.
##
## mean se_mean sd 2.5% 50% 97.5% n_eff Rhat
## b[1,1] -4.0 0 1.3 -6.6 -3.9 -1.7 2089 1
## b[1,2] 1.6 0 0.4 0.8 1.6 2.5 2918 1
## b[1,3] 3.4 0 0.9 2.0 3.2 5.4 2121 1
## b[1,4] -0.4 0 0.2 -0.7 -0.4 0.0 3029 1
## b[2,1] 2.4 0 0.2 2.0 2.4 2.7 2065 1
## b[2,2] -0.7 0 0.1 -0.9 -0.7 -0.6 2728 1
## b[2,3] -0.2 0 0.1 -0.3 -0.2 0.0 2408 1
## b[2,4] 0.2 0 0.0 0.1 0.2 0.3 2618 1
## q[1] 1.0 0 0.0 0.9 1.0 1.0 3471 1
## q[2] 0.9 0 0.1 0.7 0.9 1.0 2949 1
## q[3] 1.0 0 0.0 1.0 1.0 1.0 2623 1
## q[4] 0.6 0 0.1 0.4 0.6 0.7 3198 1
## q[5] 0.6 0 0.1 0.4 0.6 0.7 3198 1
## q[6] 0.9 0 0.1 0.7 0.9 0.9 2964 1
## q[7] 1.0 0 0.0 0.9 1.0 1.0 3495 1
## q[8] 0.6 0 0.1 0.4 0.6 0.7 3213 1
## q[9] 0.5 0 0.1 0.4 0.5 0.7 3229 1
## q[10] 0.5 0 0.1 0.4 0.5 0.7 3229 1
## q[11] 0.5 0 0.1 0.4 0.5 0.7 3229 1
## q[12] 1.0 0 0.0 0.9 1.0 1.0 3505 1
## q[13] 0.5 0 0.1 0.4 0.5 0.7 3229 1
## q[14] 0.8 0 0.1 0.7 0.8 0.9 3011 1
## q[15] 0.5 0 0.1 0.4 0.5 0.7 3246 1
## q[16] 0.8 0 0.1 0.7 0.8 0.9 3011 1
## q[17] 1.0 0 0.0 0.9 1.0 1.0 3516 1
## q[18] 0.5 0 0.1 0.4 0.5 0.7 3246 1
## q[19] 0.8 0 0.1 0.7 0.8 0.9 3027 1
## q[20] 0.8 0 0.1 0.7 0.8 0.9 3027 1
## q[21] 0.5 0 0.1 0.4 0.5 0.7 3250 1
## q[22] 0.8 0 0.1 0.7 0.8 0.9 3027 1
## q[23] 0.8 0 0.1 0.7 0.8 0.9 3027 1
## q[24] 0.8 0 0.1 0.7 0.8 0.9 3027 1
## q[25] 1.0 0 0.0 0.9 1.0 1.0 3525 1
## q[26] 0.5 0 0.1 0.4 0.5 0.7 3250 1
## q[27] 0.5 0 0.1 0.4 0.5 0.7 3250 1
## q[28] 0.5 0 0.1 0.4 0.5 0.7 3254 1
## q[29] 0.8 0 0.1 0.7 0.8 0.9 3044 1
## q[30] 1.0 0 0.0 1.0 1.0 1.0 2597 1
## q[31] 1.0 0 0.0 1.0 1.0 1.0 2611 1
## q[32] 0.5 0 0.1 0.4 0.5 0.6 3259 1
## q[33] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[34] 1.0 0 0.0 0.9 1.0 1.0 3547 1
## q[35] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[36] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[37] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[38] 0.8 0 0.1 0.7 0.8 0.9 3079 1
## q[39] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[40] 0.8 0 0.1 0.7 0.8 0.9 3079 1
## q[41] 0.5 0 0.1 0.4 0.5 0.6 3264 1
## q[42] 0.8 0 0.1 0.7 0.8 0.9 3097 1
## q[43] 1.0 0 0.0 0.9 1.0 1.0 3552 1
## q[44] 0.8 0 0.1 0.7 0.8 0.9 3097 1
## q[45] 1.0 0 0.0 0.9 1.0 1.0 3552 1
## q[46] 1.0 0 0.0 1.0 1.0 1.0 2657 1
## q[47] 1.0 0 0.0 0.9 1.0 1.0 3556 1
## q[48] 0.8 0 0.1 0.7 0.8 0.9 3115 1
## q[49] 0.8 0 0.1 0.7 0.8 0.9 3115 1
## q[50] 0.8 0 0.1 0.7 0.8 0.9 3115 1
## q[51] 0.5 0 0.1 0.4 0.5 0.6 3276 1
## q[52] 0.5 0 0.1 0.4 0.5 0.6 3276 1
## q[53] 0.5 0 0.1 0.4 0.5 0.6 3276 1
## q[54] 0.5 0 0.1 0.4 0.5 0.6 3276 1
## q[55] 0.5 0 0.1 0.4 0.5 0.6 3283 1
## q[56] 1.0 0 0.0 1.0 1.0 1.0 2654 1
## q[57] 1.0 0 0.0 0.9 1.0 1.0 3559 1
## q[58] 1.0 0 0.0 0.9 1.0 1.0 3559 1
## q[59] 0.5 0 0.1 0.3 0.5 0.6 3289 1
## q[60] 0.8 0 0.1 0.7 0.8 0.9 3154 1
## q[61] 1.0 0 0.0 1.0 1.0 1.0 2647 1
## q[62] 0.8 0 0.1 0.7 0.8 0.9 3154 1
## q[63] 0.8 0 0.1 0.7 0.8 0.9 3154 1
## q[64] 0.8 0 0.1 0.7 0.8 0.9 3154 1
## q[65] 0.5 0 0.1 0.3 0.5 0.6 3297 1
## q[66] 0.5 0 0.1 0.3 0.5 0.6 3297 1
## q[67] 0.5 0 0.1 0.3 0.5 0.6 3297 1
## q[68] 0.5 0 0.1 0.3 0.5 0.6 3297 1
## q[69] 0.4 0 0.1 0.3 0.4 0.6 3304 1
## q[70] 0.9 0 0.0 0.9 1.0 1.0 3561 1
## q[71] 0.8 0 0.1 0.7 0.8 0.9 3193 1
## q[72] 0.9 0 0.0 0.9 1.0 1.0 3561 1
## q[73] 0.4 0 0.1 0.3 0.4 0.6 3304 1
## q[74] 0.9 0 0.0 0.9 1.0 1.0 3561 1
## q[75] 0.8 0 0.1 0.7 0.8 0.9 3193 1
## q[76] 0.8 0 0.1 0.7 0.8 0.9 3193 1
## q[77] 0.8 0 0.1 0.7 0.8 0.9 3193 1
## q[78] 0.4 0 0.1 0.3 0.4 0.5 3311 1
## q[79] 0.4 0 0.1 0.3 0.4 0.5 3311 1
## q[80] 1.0 0 0.0 1.0 1.0 1.0 2617 1
## q[81] 0.8 0 0.1 0.7 0.8 0.9 3213 1
## q[82] 0.9 0 0.0 0.9 1.0 1.0 3559 1
## q[83] 0.4 0 0.1 0.3 0.4 0.5 3311 1
## q[84] 0.4 0 0.1 0.3 0.4 0.5 3311 1
## q[85] 0.9 0 0.0 0.9 1.0 1.0 3559 1
## q[86] 0.4 0 0.1 0.3 0.4 0.5 3311 1
## q[87] 0.4 0 0.1 0.3 0.4 0.5 3318 1
## q[88] 0.4 0 0.1 0.3 0.4 0.5 3318 1
## q[89] 0.4 0 0.1 0.3 0.4 0.5 3318 1
## q[90] 1.0 0 0.0 1.0 1.0 1.0 2605 1
## q[91] 0.9 0 0.0 0.9 0.9 1.0 3556 1
## q[92] 0.4 0 0.1 0.3 0.4 0.5 3325 1
## q[93] 0.4 0 0.1 0.3 0.4 0.5 3325 1
## q[94] 1.0 0 0.0 1.0 1.0 1.0 2591 1
## q[95] 0.8 0 0.1 0.6 0.8 0.9 3251 1
## q[96] 1.0 0 0.0 1.0 1.0 1.0 2591 1
## q[97] 0.8 0 0.1 0.6 0.8 0.9 3251 1
## q[98] 0.4 0 0.1 0.3 0.4 0.5 3331 1
## q[99] 0.4 0 0.1 0.3 0.4 0.5 3331 1
## q[100] 1.0 0 0.0 1.0 1.0 1.0 2582 1
## q[101] 0.8 0 0.1 0.6 0.8 0.9 3269 1
## q[102] 0.4 0 0.1 0.3 0.4 0.5 3331 1
## q[103] 0.8 0 0.1 0.6 0.8 0.9 3269 1
## q[104] 0.8 0 0.1 0.6 0.8 0.9 3286 1
## q[105] 0.8 0 0.1 0.6 0.8 0.9 3286 1
## q[106] 0.4 0 0.1 0.3 0.4 0.5 3336 1
## q[107] 0.4 0 0.1 0.3 0.4 0.5 3340 1
## q[108] 0.4 0 0.1 0.3 0.4 0.5 3340 1
## q[109] 0.9 0 0.0 0.8 0.9 1.0 3478 1
## q[110] 0.4 0 0.1 0.3 0.4 0.5 3340 1
## q[111] 0.4 0 0.1 0.3 0.4 0.5 3340 1
## q[112] 0.4 0 0.1 0.3 0.4 0.5 3340 1
## q[113] 1.0 0 0.0 0.9 1.0 1.0 2572 1
## q[114] 0.4 0 0.1 0.3 0.4 0.5 3343 1
## q[115] 0.9 0 0.0 0.8 0.9 1.0 3451 1
## q[116] 0.7 0 0.1 0.6 0.7 0.9 3308 1
## q[117] 0.7 0 0.1 0.6 0.7 0.9 3308 1
## q[118] 0.9 0 0.0 0.8 0.9 1.0 3422 1
## q[119] 0.4 0 0.1 0.3 0.4 0.5 3346 1
## q[120] 1.0 0 0.0 0.9 1.0 1.0 2570 1
## q[121] 0.7 0 0.1 0.6 0.7 0.9 3320 1
## q[122] 0.4 0 0.1 0.3 0.4 0.5 3346 1
## q[123] 0.4 0 0.1 0.3 0.4 0.5 3347 1
## q[124] 0.7 0 0.1 0.6 0.7 0.9 3340 1
## q[125] 0.4 0 0.1 0.3 0.4 0.5 3347 1
## q[126] 0.4 0 0.1 0.3 0.4 0.5 3347 1
## q[127] 1.0 0 0.0 0.9 1.0 1.0 2569 1
## q[128] 0.4 0 0.1 0.3 0.4 0.5 3347 1
## q[129] 0.9 0 0.0 0.8 0.9 1.0 3390 1
## q[130] 0.7 0 0.1 0.6 0.7 0.8 3359 1
## q[131] 1.0 0 0.0 0.9 1.0 1.0 2567 1
## q[132] 0.4 0 0.1 0.3 0.4 0.5 3348 1
## q[133] 0.7 0 0.1 0.6 0.7 0.8 3359 1
## q[134] 1.0 0 0.0 0.9 1.0 1.0 2567 1
## q[135] 0.4 0 0.1 0.3 0.4 0.5 3348 1
## q[136] 0.9 0 0.0 0.8 0.9 1.0 3324 1
## q[137] 0.7 0 0.1 0.6 0.7 0.8 3377 1
## q[138] 1.0 0 0.0 0.9 1.0 1.0 2565 1
## q[139] 0.7 0 0.1 0.6 0.7 0.8 3377 1
## q[140] 0.7 0 0.1 0.6 0.7 0.8 3377 1
## q[141] 0.9 0 0.0 0.8 0.9 1.0 3324 1
## q[142] 0.4 0 0.1 0.2 0.3 0.5 3348 1
## q[143] 0.9 0 0.0 0.8 0.9 1.0 3288 1
## q[144] 0.3 0 0.1 0.2 0.3 0.5 3348 1
## q[145] 0.9 0 0.0 0.8 0.9 1.0 3288 1
## q[146] 0.9 0 0.0 0.8 0.9 1.0 3288 1
## q[147] 0.7 0 0.1 0.6 0.7 0.8 3394 1
## q[148] 0.9 0 0.0 0.8 0.9 1.0 3288 1
## q[149] 0.9 0 0.0 0.8 0.9 1.0 3252 1
## q[150] 0.7 0 0.1 0.6 0.7 0.8 3409 1
## q[151] 0.3 0 0.1 0.2 0.3 0.5 3309 1
## q[152] 0.7 0 0.1 0.5 0.7 0.8 3424 1
## q[153] 0.3 0 0.1 0.2 0.3 0.5 3309 1
## q[154] 0.9 0 0.1 0.8 0.9 1.0 3216 1
## q[155] 0.7 0 0.1 0.5 0.7 0.8 3424 1
## q[156] 0.7 0 0.1 0.5 0.7 0.8 3424 1
## q[157] 0.7 0 0.1 0.5 0.7 0.8 3437 1
## q[158] 0.3 0 0.1 0.2 0.3 0.5 3290 1
## q[159] 0.3 0 0.1 0.2 0.3 0.5 3290 1
## q[160] 0.7 0 0.1 0.5 0.7 0.8 3437 1
## q[161] 0.3 0 0.1 0.2 0.3 0.5 3272 1
## q[162] 0.9 0 0.1 0.8 0.9 1.0 3173 1
## q[163] 1.0 0 0.0 0.9 1.0 1.0 2546 1
## q[164] 0.3 0 0.1 0.2 0.3 0.5 3272 1
## q[165] 0.3 0 0.1 0.2 0.3 0.5 3272 1
## q[166] 0.7 0 0.1 0.5 0.7 0.8 3443 1
## q[167] 0.3 0 0.1 0.2 0.3 0.5 3259 1
## q[168] 0.9 0 0.1 0.8 0.9 1.0 3151 1
## q[169] 0.3 0 0.1 0.2 0.3 0.5 3251 1
## q[170] 0.9 0 0.1 0.8 0.9 1.0 3128 1
## q[171] 0.3 0 0.1 0.2 0.3 0.5 3251 1
## q[172] 0.9 0 0.1 0.8 0.9 1.0 3128 1
## q[173] 0.7 0 0.1 0.5 0.7 0.8 3449 1
## q[174] 0.3 0 0.1 0.2 0.3 0.4 3241 1
## q[175] 0.9 0 0.1 0.8 0.9 1.0 3105 1
## q[176] 1.0 0 0.0 0.9 1.0 1.0 2516 1
## q[177] 0.7 0 0.1 0.5 0.7 0.8 3449 1
## q[178] 0.3 0 0.1 0.2 0.3 0.4 3241 1
## q[179] 0.3 0 0.1 0.2 0.3 0.4 3241 1
## q[180] 0.7 0 0.1 0.5 0.7 0.8 3449 1
## q[181] 0.9 0 0.1 0.8 0.9 1.0 3083 1
## q[182] 1.0 0 0.0 0.9 1.0 1.0 2507 1
## q[183] 0.9 0 0.1 0.8 0.9 1.0 3083 1
## q[184] 0.6 0 0.1 0.5 0.6 0.8 3449 1
## q[185] 0.9 0 0.1 0.8 0.9 1.0 3083 1
## q[186] 0.6 0 0.1 0.5 0.6 0.8 3449 1
## q[187] 0.3 0 0.1 0.1 0.3 0.4 3221 1
## q[188] 0.3 0 0.1 0.1 0.3 0.4 3221 1
## q[189] 0.9 0 0.1 0.7 0.9 1.0 3060 1
## q[190] 0.3 0 0.1 0.1 0.3 0.4 3210 1
## q[191] 0.6 0 0.1 0.4 0.6 0.8 3443 1
## q[192] 0.3 0 0.1 0.1 0.3 0.4 3210 1
## q[193] 1.0 0 0.0 0.9 1.0 1.0 2489 1
## q[194] 0.6 0 0.1 0.4 0.6 0.8 3443 1
## q[195] 0.3 0 0.1 0.1 0.3 0.4 3210 1
## q[196] 0.9 0 0.1 0.7 0.9 1.0 3031 1
## q[197] 0.3 0 0.1 0.1 0.3 0.4 3199 1
## q[198] 0.3 0 0.1 0.1 0.3 0.4 3199 1
## q[199] 0.9 0 0.1 0.7 0.9 1.0 3031 1
## q[200] 0.6 0 0.1 0.4 0.6 0.8 3439 1
## lambda[1] 1.6 0 0.1 1.5 1.6 1.8 2895 1
## lambda[2] 1.1 0 0.1 0.9 1.1 1.2 3487 1
## lambda[3] 0.9 0 0.1 0.7 0.9 1.1 2818 1
## lambda[4] 1.8 0 0.1 1.7 1.8 2.0 2362 1
## lambda[5] 1.8 0 0.1 1.7 1.8 2.0 2362 1
## lambda[6] 1.1 0 0.1 0.9 1.1 1.3 3525 1
## lambda[7] 1.7 0 0.1 1.5 1.7 1.8 2942 1
## lambda[8] 1.8 0 0.1 1.7 1.8 2.0 2362 1
## lambda[9] 1.9 0 0.1 1.7 1.9 2.0 2363 1
## lambda[10] 1.9 0 0.1 1.7 1.9 2.0 2363 1
## lambda[11] 1.9 0 0.1 1.7 1.9 2.0 2363 1
## lambda[12] 1.7 0 0.1 1.5 1.7 1.9 2969 1
## lambda[13] 1.9 0 0.1 1.7 1.9 2.0 2363 1
## lambda[14] 1.1 0 0.1 1.0 1.1 1.3 3646 1
## lambda[15] 1.9 0 0.1 1.7 1.9 2.0 2364 1
## lambda[16] 1.1 0 0.1 1.0 1.1 1.3 3646 1
## lambda[17] 1.7 0 0.1 1.6 1.7 1.9 2999 1
## lambda[18] 1.9 0 0.1 1.7 1.9 2.0 2364 1
## lambda[19] 1.2 0 0.1 1.0 1.2 1.3 3687 1
## lambda[20] 1.2 0 0.1 1.0 1.2 1.3 3687 1
## lambda[21] 1.9 0 0.1 1.8 1.9 2.0 2367 1
## lambda[22] 1.2 0 0.1 1.0 1.2 1.3 3687 1
## lambda[23] 1.2 0 0.1 1.0 1.2 1.3 3687 1
## lambda[24] 1.2 0 0.1 1.0 1.2 1.3 3687 1
## lambda[25] 1.7 0 0.1 1.6 1.7 1.9 3031 1
## lambda[26] 1.9 0 0.1 1.8 1.9 2.0 2367 1
## lambda[27] 1.9 0 0.1 1.8 1.9 2.0 2367 1
## lambda[28] 1.9 0 0.1 1.8 1.9 2.1 2372 1
## lambda[29] 1.2 0 0.1 1.0 1.2 1.3 3729 1
## lambda[30] 1.0 0 0.1 0.8 1.0 1.2 2699 1
## lambda[31] 1.0 0 0.1 0.8 1.0 1.2 2675 1
## lambda[32] 1.9 0 0.1 1.8 1.9 2.1 2378 1
## lambda[33] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[34] 1.8 0 0.1 1.7 1.8 1.9 3151 1
## lambda[35] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[36] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[37] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[38] 1.2 0 0.1 1.1 1.2 1.4 3808 1
## lambda[39] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[40] 1.2 0 0.1 1.1 1.2 1.4 3808 1
## lambda[41] 2.0 0 0.1 1.8 2.0 2.1 2386 1
## lambda[42] 1.2 0 0.1 1.1 1.2 1.4 3844 1
## lambda[43] 1.8 0 0.1 1.7 1.8 2.0 3199 1
## lambda[44] 1.2 0 0.1 1.1 1.2 1.4 3844 1
## lambda[45] 1.8 0 0.1 1.7 1.8 2.0 3199 1
## lambda[46] 1.1 0 0.1 0.9 1.1 1.3 2627 1
## lambda[47] 1.8 0 0.1 1.7 1.8 2.0 3251 1
## lambda[48] 1.3 0 0.1 1.1 1.3 1.4 3876 1
## lambda[49] 1.3 0 0.1 1.1 1.3 1.4 3876 1
## lambda[50] 1.3 0 0.1 1.1 1.3 1.4 3876 1
## lambda[51] 2.0 0 0.1 1.9 2.0 2.1 2409 1
## lambda[52] 2.0 0 0.1 1.9 2.0 2.1 2409 1
## lambda[53] 2.0 0 0.1 1.9 2.0 2.1 2409 1
## lambda[54] 2.0 0 0.1 1.9 2.0 2.1 2409 1
## lambda[55] 2.0 0 0.1 1.9 2.0 2.1 2424 1
## lambda[56] 1.1 0 0.1 0.9 1.1 1.3 2576 1
## lambda[57] 1.9 0 0.1 1.7 1.9 2.0 3308 1
## lambda[58] 1.9 0 0.1 1.7 1.9 2.0 3308 1
## lambda[59] 2.0 0 0.1 1.9 2.0 2.2 2432 1
## lambda[60] 1.3 0 0.1 1.2 1.3 1.4 3923 1
## lambda[61] 1.1 0 0.1 1.0 1.1 1.3 2550 1
## lambda[62] 1.3 0 0.1 1.2 1.3 1.4 3923 1
## lambda[63] 1.3 0 0.1 1.2 1.3 1.4 3923 1
## lambda[64] 1.3 0 0.1 1.2 1.3 1.4 3923 1
## lambda[65] 2.1 0 0.1 2.0 2.1 2.2 2439 1
## lambda[66] 2.1 0 0.1 2.0 2.1 2.2 2439 1
## lambda[67] 2.1 0 0.1 2.0 2.1 2.2 2439 1
## lambda[68] 2.1 0 0.1 2.0 2.1 2.2 2439 1
## lambda[69] 2.1 0 0.1 2.0 2.1 2.2 2447 1
## lambda[70] 1.9 0 0.1 1.8 1.9 2.0 3441 1
## lambda[71] 1.3 0 0.1 1.2 1.3 1.5 3940 1
## lambda[72] 1.9 0 0.1 1.8 1.9 2.0 3441 1
## lambda[73] 2.1 0 0.1 2.0 2.1 2.2 2447 1
## lambda[74] 1.9 0 0.1 1.8 1.9 2.0 3441 1
## lambda[75] 1.3 0 0.1 1.2 1.3 1.5 3940 1
## lambda[76] 1.3 0 0.1 1.2 1.3 1.5 3940 1
## lambda[77] 1.3 0 0.1 1.2 1.3 1.5 3940 1
## lambda[78] 2.1 0 0.1 2.0 2.1 2.2 2455 1
## lambda[79] 2.1 0 0.1 2.0 2.1 2.2 2455 1
## lambda[80] 1.2 0 0.1 1.0 1.2 1.4 2474 1
## lambda[81] 1.4 0 0.1 1.2 1.4 1.5 3935 1
## lambda[82] 1.9 0 0.1 1.8 1.9 2.1 3483 1
## lambda[83] 2.1 0 0.1 2.0 2.1 2.2 2455 1
## lambda[84] 2.1 0 0.1 2.0 2.1 2.2 2455 1
## lambda[85] 1.9 0 0.1 1.8 1.9 2.1 3483 1
## lambda[86] 2.1 0 0.1 2.0 2.1 2.2 2455 1
## lambda[87] 2.1 0 0.1 2.0 2.1 2.2 2467 1
## lambda[88] 2.1 0 0.1 2.0 2.1 2.2 2467 1
## lambda[89] 2.1 0 0.1 2.0 2.1 2.2 2467 1
## lambda[90] 1.2 0 0.1 1.0 1.2 1.4 2450 1
## lambda[91] 2.0 0 0.1 1.8 2.0 2.1 3551 1
## lambda[92] 2.1 0 0.1 2.0 2.1 2.3 2483 1
## lambda[93] 2.1 0 0.1 2.0 2.1 2.3 2483 1
## lambda[94] 1.2 0 0.1 1.1 1.2 1.4 2426 1
## lambda[95] 1.4 0 0.1 1.3 1.4 1.5 3893 1
## lambda[96] 1.2 0 0.1 1.1 1.2 1.4 2426 1
## lambda[97] 1.4 0 0.1 1.3 1.4 1.5 3893 1
## lambda[98] 2.2 0 0.1 2.1 2.2 2.3 2501 1
## lambda[99] 2.2 0 0.1 2.1 2.2 2.3 2501 1
## lambda[100] 1.3 0 0.1 1.1 1.3 1.4 2403 1
## lambda[101] 1.4 0 0.1 1.3 1.4 1.5 3857 1
## lambda[102] 2.2 0 0.1 2.1 2.2 2.3 2501 1
## lambda[103] 1.4 0 0.1 1.3 1.4 1.5 3857 1
## lambda[104] 1.4 0 0.1 1.3 1.4 1.6 3813 1
## lambda[105] 1.4 0 0.1 1.3 1.4 1.6 3813 1
## lambda[106] 2.2 0 0.1 2.1 2.2 2.3 2523 1
## lambda[107] 2.2 0 0.1 2.1 2.2 2.3 2554 1
## lambda[108] 2.2 0 0.1 2.1 2.2 2.3 2554 1
## lambda[109] 2.0 0 0.1 1.9 2.0 2.1 3871 1
## lambda[110] 2.2 0 0.1 2.1 2.2 2.3 2554 1
## lambda[111] 2.2 0 0.1 2.1 2.2 2.3 2554 1
## lambda[112] 2.2 0 0.1 2.1 2.2 2.3 2554 1
## lambda[113] 1.3 0 0.1 1.1 1.3 1.5 2343 1
## lambda[114] 2.2 0 0.1 2.1 2.2 2.3 2589 1
## lambda[115] 2.1 0 0.1 1.9 2.1 2.2 3920 1
## lambda[116] 1.5 0 0.1 1.3 1.5 1.6 3701 1
## lambda[117] 1.5 0 0.1 1.3 1.5 1.6 3701 1
## lambda[118] 2.1 0 0.1 2.0 2.1 2.2 3960 1
## lambda[119] 2.2 0 0.1 2.1 2.2 2.4 2628 1
## lambda[120] 1.3 0 0.1 1.2 1.3 1.5 2325 1
## lambda[121] 1.5 0 0.1 1.4 1.5 1.6 3637 1
## lambda[122] 2.2 0 0.1 2.1 2.2 2.4 2628 1
## lambda[123] 2.3 0 0.1 2.1 2.3 2.4 2669 1
## lambda[124] 1.5 0 0.1 1.4 1.5 1.6 3570 1
## lambda[125] 2.3 0 0.1 2.1 2.3 2.4 2669 1
## lambda[126] 2.3 0 0.1 2.1 2.3 2.4 2669 1
## lambda[127] 1.4 0 0.1 1.2 1.4 1.5 2310 1
## lambda[128] 2.3 0 0.1 2.1 2.3 2.4 2669 1
## lambda[129] 2.1 0 0.1 2.0 2.1 2.2 3987 1
## lambda[130] 1.5 0 0.1 1.4 1.5 1.7 3502 1
## lambda[131] 1.4 0 0.1 1.2 1.4 1.5 2297 1
## lambda[132] 2.3 0 0.1 2.2 2.3 2.4 2709 1
## lambda[133] 1.5 0 0.1 1.4 1.5 1.7 3502 1
## lambda[134] 1.4 0 0.1 1.2 1.4 1.5 2297 1
## lambda[135] 2.3 0 0.1 2.2 2.3 2.4 2709 1
## lambda[136] 2.1 0 0.1 2.0 2.1 2.2 3969 1
## lambda[137] 1.6 0 0.1 1.4 1.6 1.7 3434 1
## lambda[138] 1.4 0 0.1 1.2 1.4 1.6 2285 1
## lambda[139] 1.6 0 0.1 1.4 1.6 1.7 3434 1
## lambda[140] 1.6 0 0.1 1.4 1.6 1.7 3434 1
## lambda[141] 2.1 0 0.1 2.0 2.1 2.2 3969 1
## lambda[142] 2.3 0 0.1 2.2 2.3 2.4 2751 1
## lambda[143] 2.2 0 0.1 2.0 2.2 2.3 3937 1
## lambda[144] 2.3 0 0.1 2.2 2.3 2.4 2792 1
## lambda[145] 2.2 0 0.1 2.0 2.2 2.3 3937 1
## lambda[146] 2.2 0 0.1 2.0 2.2 2.3 3937 1
## lambda[147] 1.6 0 0.1 1.4 1.6 1.7 3367 1
## lambda[148] 2.2 0 0.1 2.0 2.2 2.3 3937 1
## lambda[149] 2.2 0 0.1 2.1 2.2 2.3 3893 1
## lambda[150] 1.6 0 0.1 1.5 1.6 1.7 3303 1
## lambda[151] 2.4 0 0.1 2.2 2.4 2.5 2881 1
## lambda[152] 1.6 0 0.1 1.5 1.6 1.8 3241 1
## lambda[153] 2.4 0 0.1 2.2 2.4 2.5 2881 1
## lambda[154] 2.2 0 0.1 2.1 2.2 2.3 3840 1
## lambda[155] 1.6 0 0.1 1.5 1.6 1.8 3241 1
## lambda[156] 1.6 0 0.1 1.5 1.6 1.8 3241 1
## lambda[157] 1.6 0 0.1 1.5 1.6 1.8 3183 1
## lambda[158] 2.4 0 0.1 2.3 2.4 2.5 2894 1
## lambda[159] 2.4 0 0.1 2.3 2.4 2.5 2894 1
## lambda[160] 1.6 0 0.1 1.5 1.6 1.8 3183 1
## lambda[161] 2.4 0 0.1 2.3 2.4 2.5 2904 1
## lambda[162] 2.2 0 0.1 2.1 2.2 2.4 3718 1
## lambda[163] 1.5 0 0.1 1.3 1.5 1.7 2260 1
## lambda[164] 2.4 0 0.1 2.3 2.4 2.5 2904 1
## lambda[165] 2.4 0 0.1 2.3 2.4 2.5 2904 1
## lambda[166] 1.7 0 0.1 1.5 1.7 1.8 3128 1
## lambda[167] 2.4 0 0.1 2.3 2.4 2.6 2911 1
## lambda[168] 2.3 0 0.1 2.1 2.3 2.4 3654 1
## lambda[169] 2.4 0 0.1 2.3 2.4 2.6 2916 1
## lambda[170] 2.3 0 0.1 2.2 2.3 2.4 3591 1
## lambda[171] 2.4 0 0.1 2.3 2.4 2.6 2916 1
## lambda[172] 2.3 0 0.1 2.2 2.3 2.4 3591 1
## lambda[173] 1.7 0 0.1 1.6 1.7 1.9 2988 1
## lambda[174] 2.5 0 0.1 2.3 2.5 2.6 2919 1
## lambda[175] 2.3 0 0.1 2.2 2.3 2.4 3529 1
## lambda[176] 1.6 0 0.1 1.4 1.6 1.7 2264 1
## lambda[177] 1.7 0 0.1 1.6 1.7 1.9 2988 1
## lambda[178] 2.5 0 0.1 2.3 2.5 2.6 2919 1
## lambda[179] 2.5 0 0.1 2.3 2.5 2.6 2919 1
## lambda[180] 1.7 0 0.1 1.6 1.7 1.9 2988 1
## lambda[181] 2.3 0 0.1 2.2 2.3 2.4 3469 1
## lambda[182] 1.6 0 0.1 1.4 1.6 1.7 2256 1
## lambda[183] 2.3 0 0.1 2.2 2.3 2.4 3469 1
## lambda[184] 1.7 0 0.1 1.6 1.7 1.9 2948 1
## lambda[185] 2.3 0 0.1 2.2 2.3 2.4 3469 1
## lambda[186] 1.7 0 0.1 1.6 1.7 1.9 2948 1
## lambda[187] 2.5 0 0.1 2.3 2.5 2.7 2921 1
## lambda[188] 2.5 0 0.1 2.3 2.5 2.7 2921 1
## lambda[189] 2.3 0 0.1 2.2 2.3 2.5 3413 1
## lambda[190] 2.5 0 0.1 2.4 2.5 2.7 2920 1
## lambda[191] 1.8 0 0.1 1.6 1.8 1.9 2880 1
## lambda[192] 2.5 0 0.1 2.4 2.5 2.7 2920 1
## lambda[193] 1.6 0 0.1 1.4 1.6 1.8 2243 1
## lambda[194] 1.8 0 0.1 1.6 1.8 1.9 2880 1
## lambda[195] 2.5 0 0.1 2.4 2.5 2.7 2920 1
## lambda[196] 2.4 0 0.1 2.2 2.4 2.5 3310 1
## lambda[197] 2.5 0 0.1 2.4 2.5 2.7 2918 1
## lambda[198] 2.5 0 0.1 2.4 2.5 2.7 2918 1
## lambda[199] 2.4 0 0.1 2.2 2.4 2.5 3310 1
## lambda[200] 1.8 0 0.1 1.6 1.8 2.0 2850 1
##
## Samples were drawn using NUTS(diag_e) at Fri Jan 5 14:14:40 2024.
## 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).
\(q, \lambda\)の順位相関係数の中央値は負であるため、来店確率が高いことと、リピーターで多数来店することは、負の関係がある。
ms <- rstan::extract(fit)
N_mcmc <- length(ms$lp__)
r <- sapply(1:N_mcmc,
function(i) cor(ms$lambda[i,], ms$q[i,], method = 'spearman'))
quantile(r, prob = c(0.025, 0.25, 0.5, 0.75, 0.975))
## 2.5% 25% 50% 75% 97.5%
## -0.8081299 -0.6974237 -0.6496459 -0.6050023 -0.4827052