r语言 - 如何向 rjags 模型添加样条曲线

Question

我很难找到使用 rjag 拟合样条的信息（我的动机是尝试在 jag 中重新创建 glm 以插补缺失的相关值）。无论如何，我只能找到有关此查找的信息，只有交叉验证的答案：https://stats.stackexchange.com/questions/79973/how-to-analyze-this-data-using-rjags-or-any-other-way/80650#80650。但是，我无法理解那里的样条代码（并且没有提出问题的声誉！首先，我不明白为什么该代码同时循环S&G。

因此，我制作了一个锯齿状的玩具线性模型：

library(datasets)
library(rjags)
library(ggplot2)

# Specify a JAGS linear model
mk_jags_lin_mod <- function(prior.a, prior.b){
    sink(paste("lin_reg_jags.mod.txt", sep=""))
    cat(paste0("
               model {
               for (i in 1:N){
               y[i] ~ dnorm(y.hat[i], tau)
               y.hat[i] <- a + b * x[i]
               }
               a ~ ",prior.a,
               "tb ~ ",prior.b,
               "ttau <- pow(sigma, -2)
               sigma ~ dunif(0, 100)
               }
               "))
    sink()
}
# Define a default vague prior
default <- "dnorm(0, .0001)n"
mk_jags_lin_mod(default, default)
# Initialise
jags.cars <- jags.model('lin_reg_jags.mod.txt',
                              data = list('x' = mtcars$hp,
                                          'y' = mtcars$mpg,
                                          'N' = nrow(mtcars)),
                              n.chains = 2,
                              n.adapt = 1000)
# Burn-in
update(jags.cars, 5000) 
# Sample
coda.cars <- coda.samples(jags.cars, variable.names = c('a', 'b', 'y.hat','tau'), n.iter = 1000)
# Extract posterior estimates
coda.sum <- summary(coda.cars)
q <- coda.sum$quantiles
mtcars$fit <- q[4:35 , 3]
ggplot(data=mtcars, aes(x=hp, y=mpg)) + geom_point() +
    geom_line(aes(y=fit, col="red"))

我的问题是：如何将受限三次样条添加到线性模型中的"b"估计中？

Answer 1

好吧，

我花了大半天的时间，但我终于想通了......我认为。。。。我修改了模型，如下所示：

 model {
   for (i in 1:N){
     y[i] ~ dnorm(y.hat[i], tau)
     y.hat[i] <- a + beta[1] + beta[2]*x[i] + beta[3]*pow(x[i], 2) + beta[4]*pow(x[i], 3)
   }
   a ~ dnorm(0, .0001)
   # Specify priors for spline terms
   for (k in 1:4) {
     beta.mu[k] ~ dnorm(0, 100)
     beta.tau[k] ~ dgamma(0.01, 10)
     beta[k] ~ dnorm(beta.mu[k], beta.tau[k])
   }
     tau <- pow(sigma, -2)
     sigma ~ dunif(0, 100)
   }

k

的值必须至少为 4，但此后增加 k 会增加平滑。否则我认为我的数学是正确的（尽管我愿意接受更正！