机器算法验证 - 通过 glmmTMB 包了解 AR1 - 吾爱随笔录

我一直在研究一个可重复的示例，以使用该glmmTMB包更好地理解 AR1 协方差矩阵。我有几个问题，即使只能回答一个，我将不胜感激。

这是一些合成数据，非常感谢 Ben Bolker 的代码：

library(glmmTMB)
library(MASS)    
simCor1 <- function(phi=0.8,sdgrp=2,sdres=1,
                    npergrp=20,ngrp=20,
                    seed=NULL,
                    ## set linkinv/simfun for GLMM sims
                    linkinv=identity,
                    simfun=identity) {
  if (!is.null(seed)) set.seed(seed)
  cmat <- sdres*phi^abs(outer(0:(npergrp-1),0:(npergrp-1),"-"))
  errs <- MASS::mvrnorm(ngrp,mu=rep(0,npergrp),Sigma=cmat)
  ranef <- rnorm(ngrp,mean=0,sd=sdgrp)
  d <- data.frame(f=rep(1:ngrp,each=npergrp))
  eta <- ranef[as.numeric(d$f)] + c(t(errs)) ## unpack errors by row
  mu <- linkinv(eta)
  d$y <- simfun(mu)
  d$tt <- factor(rep(1:npergrp,ngrp))
  return(d)
}
d <- simCor1(phi = 0.8, seed = 100)

这是一个基本的 AR1 模型glmmTMB：

mod <- glmmTMB(y~1 + (1|f) + ar1(tt-1|f), data=d,family=gaussian)

在分析模型拟合时，它很好地识别了估计值为 0.85 的 phi 参数：

summary(mod)
...    
Conditional model:
 Groups   Name        Variance Std.Dev. Corr      
 f        (Intercept) 4.15645  2.0387             
 f.1      tt1         1.06157  1.0303   0.85 (ar1)
 Residual             0.04183  0.2045             
Number of obs: 400, groups:  f, 20
...

以下是我的问题：

1）为什么包符号需要时间tt作为随机斜率ar1()才能工作？这样做的一个结果似乎是，当我使用该predict()函数时，它似乎具有追溯优势，即了解每个组在每个时间间隔内的表现，从而得出非常接近实际值的预测y。使用可重现的示例，其中phi=0：

d <- simCor1(phi = 0, npergrp = 100, seed = 50)
ar1mod <- glmmTMB(y~1 + (1|f) + ar1(tt-1|f), data=d, family=gaussian)
simpmod <- glmmTMB(y~1 + (1|f), data=d, family=gaussian)

尽管 phi 被估计为有效 0ar1mod并因此得出结论认为模型的 ar1 组件实际上是无用的，但两个模型之间的预测似乎大相径庭：

d$errar1 <- d$y - predict(ar1mod)
d$errsimp <- d$y - predict(simpmod)
sqrt(mean(d$errar1 ^ 2))
## 0.0005218347
sqrt(mean(d$errsimp ^ 2))
## 0.9983432

2）在试图理解为什么tt被用作随机斜率时，我试图重新创建glmmTMB必须进行的计算以确定它的估计phi。我认为，也许，通过取tt每个时间间隔的平均系数，然后确定滞后tt值之间的自相关，我们可以得出参数的glmmTMB估计值phi。使用modfrom 之前尝试重新创建估计的 .85 参数：

MeanTT <- colMeans(ranef(mod)$cond$f)[-1]
MeanTTlag1 <- MeanTT[-1]
MeanTTlag0 <- MeanTT[-length(MeanTT)]
cor(MeanTTlag1, MeanTTlag0)
## 0.430239

远低于 0.85。为什么这个计算没有达到 0.85 的值，glmmTMB而是用什么来计算phi？