除了来自 的gls()函数外nlme，您还可以使用 MLE 包中的arima()函数。stats这是具有这两个功能的示例。

x <- 1:100
e <- 25*arima.sim(model=list(ar=0.3),n=100)
y <- 1 + 2*x + e

###Fit the model using gls()
require(nlme)
(fit1 <- gls(y~x, corr=corAR1(0.5,form=~1)))
Generalized least squares fit by REML
  Model: y ~ x 
  Data: NULL 
  Log-restricted-likelihood: -443.6371

Coefficients:
(Intercept)           x 
   4.379304    1.957357 

Correlation Structure: AR(1)
 Formula: ~1 
 Parameter estimate(s):
      Phi 
0.3637263 
Degrees of freedom: 100 total; 98 residual
Residual standard error: 22.32908 

###Fit the model using arima()
(fit2 <- arima(y, xreg=x, order=c(1,0,0)))

Call:
arima(x = y, order = c(1, 0, 0), xreg = x)

Coefficients:
         ar1  intercept       x
      0.3352     4.5052  1.9548
s.e.  0.0960     6.1743  0.1060

sigma^2 estimated as 423.7:  log likelihood = -444.4,  aic = 896.81 

arima() 函数的优点是您可以适应更多种类的 ARMA 错误过程。如果使用 forecast 包中的 auto.arima() 函数，可以自动识别 ARMA 错误：

require(forecast)    
fit3 <- auto.arima(y, xreg=x)

看看nlmegls包中的（广义最小二乘）

您可以为回归中的错误设置相关配置文件，例如 ARMA 等：

 gls(Y ~ X, correlation=corARMA(p=1,q=1))

对于 ARMA(1,1) 错误。

使用包nlme中的函数gls。这是示例。

##Generate data frame with regressor and AR(1) error. The error term is 
## \eps_t=0.3*\eps_{t-1}+v_t
df <- data.frame(x1=rnorm(100), err=filter(rnorm(100)/5,filter=0.3,method="recursive"))

##Create ther response
df$y <- 1 + 2*df$x + df$err

###Fit the model
gls(y~x, data=df, corr=corAR1(0.5,form=~1))

Generalized least squares fit by REML
  Model: y ~ x 
  Data: df 
  Log-restricted-likelihood: 9.986475

 Coefficients:
 (Intercept)           x 
   1.040129    2.001884 

 Correlation Structure: AR(1)
 Formula: ~1 
 Parameter estimate(s):
     Phi 
 0.2686271 
Degrees of freedom: 100 total; 98 residual
Residual standard error: 0.2172698 

由于模型是使用最大似然拟合的，因此您需要提供起始值。默认起始值​​为 0，但一如既往地尝试多个值以确保收敛是好的。

正如 G 博士指出的那样，您还可以使用其他相关结构，即 ARMA。

请注意，如果回归误差的协方差矩阵不是单位矩阵的倍数，则通常最小二乘估计是一致的，因此如果您拟合具有特定协方差结构的模型，首先需要测试它是否合适。

您可以在 gls 输出上使用 predict 。请参阅 ?predict.gls。您还可以通过相关结构中的“形式”项指定观察的顺序。例如：
corr=corAR1(form=~1)表示数据的顺序是它们在表中的顺序。 corr=corAR1(form=~Year)表示顺序是因子Year..最后的“0.5”值corr=corAR1(0.5,form=~1)?通常设置为估计的参数值以表示方差结构（phi，在AR的情况下，theta在MA的情况下...... .)。正如 Rob Hyndman 所说，设置它并用于优化是可选的。