我现在正在研究理查德·麦克雷思的优秀著作
统计反思
是的。我在第八章。我有两个麻烦的来源(解释如下)。注意你需要
stan
已安装此程序包才能正常工作。
首先创建一些数据来建模和查看它。这个简单的回归有截距10和斜率2。
df <- data.frame(x = 1:100,
noise = rnorm(100, 0, 20))
df$y <- 10 + df$x*2 + df$noise
ggplot(df, aes(x, y)) + geom_point() + geom_smooth(method = "lm")
现在是模特。这个
map2stan
函数是斯坦的包装器。下面的模型
library(rethinking)
mod <- map2stan(
alist(
y ~ dnorm(mu, sigma),
mu <- a + b*x,
a ~ dnorm(0, 100),
b ~ dnorm(0, 10),
sigma ~ dunif(0,10)
), data = df, chains = 2, iter = 4000, warmup = 1000)
检查模型是否返回了正确的参数
precis(mod)
# Mean StdDev lower 0.89 upper 0.89 n_eff Rhat
# a 7.68 2.00 4.41 10.79 2612 1
# b 1.94 0.03 1.88 1.99 2543 1
# sigma 9.96 0.04 9.90 10.00 3980 1
问题1
根据书上说如果我跑
plot(mod)
我应该为每个变量绘制一个可爱的跟踪图,绘制两条链
和
热身期。相反,我得到了错误
Error in as.double(y) :
cannot coerce type 'S4' to vector of type 'double'
问题1:
我做错什么了?为什么情节不成功?
问题2
如果我们打印
mod
对象
mod
它告诉我们
map2stan model fit
6000 samples from 2 chains
我们
可以
通过从
Map2stan地图
模型(这只提取
a
截距参数。
dfSamp1 <- data.frame(chain1 = extract.samples(mod, n=4e3)$a,
row = 1:4e3)
ggplot(dfSamp1, aes(x = row, y = chain1, group = 1)) + geom_line(colour = "skyblue")
但缺少的是多重链。从模型中采样两次
dfSamp2 <- data.frame(chain1 = extract.samples(mod, n=4e3)$a,
chain2 = extract.samples(mod, n=4e3)$a,
row = 1:4e3)
但两条链上的样本完全相同。
ggplot(dfSamp2) +
geom_line(aes(x = row, y = chain1, group = 1), colour = "skyblue") +
geom_line(aes(x = row, y = chain2, group = 1), colour = "black")
这个轨迹图上应该有两种颜色。但是只有一条,因为第二条黑线和第一条(蓝色)是一样的。
如果我们查看数据帧
dfSamp2[1:10,]
chain1 chain2 row
1 5.738354 5.738354 1
2 4.934791 4.934791 2
3 5.701355 5.701355 3
4 9.701134 9.701134 4
5 6.123426 6.123426 5
6 11.769974 11.769974 6
7 7.675946 7.675946 7
8 8.858765 8.858765 8
9 6.279705 6.279705 9
10 5.805226 5.805226 10
确切地
同样的:不是从两个独立的随机游动中得到的。
问题2:
那么我该如何从
Map2stan地图
使用
extract.samples()
功能?