nlme 具有这样的能力,因为它可以拟合非线性混合模型。通过允许随机效应和相关误差,您可以将其视为nls 的扩展(仅具有固定效应的非线性回归)。
nlme 可以处理 ARMA 相关性,例如 correlation = corARMA(0.2, ~ 1, p = 0, q = 1, fixed = TRUE)。这意味着,残差是 MA(1) 过程,初始猜测系数为 0.2,但要在模型拟合期间更新。 ~ 1 表明 MA(1) 处于拦截状态,没有进一步的分组结构。
我不是nlme 方面的专家,但我知道nlme 是您所需要的。我制作了以下示例,但由于我不是专家,目前我无法获得nlme 的工作。我把它贴在这里给一个开始/味道。
set.seed(0)
x1 <- runif(100)
x2 <- runif(100)
## MA(1) correlated error, with innovation standard deviation 0.1
e <- arima.sim(model = list(ma = 0.5), n = 100, sd = 0.1)
## a true model, with `a = 0.2, g = 0.5`
y0 <- 0.2 + log(x1 ^ 0.5 + x2 ^ 0.5)
## observations
y <- y0 + e
## no need to install; it comes with R; just `library()` it
library(nlme)
fit <- nlme(y ~ a + log(x1 ^ g + x2 ^ g), fixed = a + g ~ 1,
start = list(a = 0.5, g = 1),
correlation = corARMA(0.2, form = ~ 1, p = 0, q = 1, fixed = FALSE))
类似于nls,我们有一个整体模型公式y ~ a + log(x1 ^ g + x2 ^ g),迭代过程需要起始值。我选择了start = list(a = 0.5, g = 1)。 correlation 位已经在开头解释过了。
nlme 中的fixed 和random 参数指定在整个公式中应该被视为固定效应和随机效应。由于我们没有随机效应,我们不指定它。我们想要a 和g 作为固定效果,所以我尝试了fixed = a + g ~ 1 之类的东西。不幸的是,由于某种我不知道的原因,它不太有效。我阅读了?nlme,并认为这个公式意味着我们想要一个共同的a 和g 用于所有观察,但后来nlme 报告了一个错误,说这不是一个有效的组公式。
我也在这方面进行投资;正如我所说,上面给了我们一个开始。我们已经非常接近最终答案了。
感谢user20650 指出我尴尬的错误。我应该使用gnls 函数而不是nlme。根据nlme 包的设计性质,函数lme 和nlme 必须采用random 参数才能工作。幸运的是,nlme 包中还有其他几个例程用于扩展线性模型和非线性模型。
-
gls 和 gnls 通过允许非对角方差函数扩展 lm 和 nls。
所以,我真的应该改用gnls:
## no `fixed` argument as `gnls` is a fixed-effect only
fit <- gnls(y ~ a + log(x1 ^ g + x2 ^ g), start = list(a = 0.5, g = 1),
correlation = corARMA(0.2, form = ~ 1, p = 0, q = 1, fixed = FALSE))
#Generalized nonlinear least squares fit
# Model: y ~ a + log(x1^g + x2^g)
# Data: NULL
# Log-likelihood: 92.44078
#
#Coefficients:
# a g
#0.1915396 0.5007640
#
#Correlation Structure: ARMA(0,1)
# Formula: ~1
# Parameter estimate(s):
# Theta1
#0.4184961
#Degrees of freedom: 100 total; 98 residual
#Residual standard error: 0.1050295