Running 'doRUnit.R'
  Running 'forecastTests.R'
Running the tests in 'tests/forecastTests.R' failed.
Complete output:
  > require(glarma)
  Loading required package: glarma
  > require(zoo)
  Loading required package: zoo
  Warning message:
  In library(package, lib.loc = lib.loc, character.only = TRUE, logical.return = TRUE,  :
    there is no package called 'zoo'
  > ### Model number of deaths
  > data(DriverDeaths)
  > y <- DriverDeaths[, "Deaths"]
  > X <- as.matrix(DriverDeaths[, 2:5])
  > Population <- DriverDeaths[, "Population"]
  > 
  > ### Offset included
  > glarmamod <- glarma(y, X, offset = log(Population/100000),
  +                     phiLags = c(12),
  +                     thetaLags = c(1),
  +                     type = "Poi", method = "FS",
  +                     residuals = "Pearson", maxit = 100, grad = 1e-6)
  > print(summary(glarmamod))
  
  Call: glarma(y = y, X = X, offset = log(Population/1e+05), type = "Poi", 
      method = "FS", residuals = "Pearson", phiLags = c(12), thetaLags = c(1), 
      maxit = 100, grad = 1e-06)
  
  Pearson Residuals:
      Min       1Q   Median       3Q      Max  
  -1.7016  -0.7938   0.1189   0.5897   3.7199  
  
  GLARMA Coefficients:
            Estimate  Std.Error z-ratio Pr(>|z|)
  phi_12  -0.0338742  0.0899514  -0.377    0.706
  theta_1 -0.0005256  0.0827289  -0.006    0.995
  
  Linear Model Coefficients:
             Estimate Std.Error z-ratio Pr(>|z|)   
  Intercept  -2.13559   0.80711  -2.646  0.00815 **
  ReducedBAC -0.38582   0.15598  -2.474  0.01338 * 
  FriSat      0.08663   0.09091   0.953  0.34060   
  lnOMVDRate  0.55792   0.23779   2.346  0.01896 * 
  
      Null deviance: 95.499  on 71  degrees of freedom
  Residual deviance: 67.560  on 66  degrees of freedom
  AIC: 264.1927 
  
  Number of Fisher Scoring iterations: 17
  
  LRT and Wald Test:
  Alternative hypothesis: model is a GLARMA process
  Null hypothesis: model is a GLM with the same regression structure
            Statistic p-value
  LR Test       0.160   0.923
  Wald Test     0.142   0.932
  
  > 
  > XT1 <- matrix(X[72,], nrow = 1)
  > offsetT1 <- log(Population/100000)[72]
  > 
  > mu <- forecast(glarmamod, 1, XT1, offsetT1)$mu
  > print(mu)
  [1] 1.933074
  > 
  > 
  > ### Save some values
  > allX <- X
  > allFits <- fitted(glarmamod)
  > ally <- y
  > 
  > ### Look at a succession of forecasts
  > ### Using actual values in forecasts
  > forecasts <- numeric(72)
  > for (i in (62:71)){
  +     y <- DriverDeaths[1:i, "Deaths"]
  +     X <- as.matrix(DriverDeaths[1:i, 2:5])
  +     Population <- DriverDeaths[1:i, "Population"]
  + 
  +     ## Offset included
  +     glarmamod <- glarma(y, X, offset = log(Population/100000),
  +                         phiLags = c(12),
  +                         thetaLags = c(1),
  +                         type = "Poi", method = "FS",
  +                         residuals = "Pearson", maxit = 100, grad = 1e-6)
  +     XT1 <- matrix(allX[i + 1, ], nrow = 1)
  +     offsetT1 <- log(DriverDeaths$Population[i + 1]/100000)
  +     mu <- forecast(glarmamod, 1, XT1, offsetT1)$mu
  +     if (i == 62){
  +         forecasts[1:62] <- fitted(glarmamod)
  +     }
  +     forecasts[i+1] <- mu
  + }
  > par(mfrow = c(1,1))
  > forecasts <- ts(forecasts[63:72], start = c(1985, 10), deltat = 1/12)
  > fitted <- ts(allFits, start = c(1980, 8), deltat = 1/12)
  > obs <- ts(DriverDeaths$Deaths, start = c(1980, 8), deltat = 1/12)
  > plot(obs, ylab = "Driver Deaths", lty = 2,
  +      main = "Single Vehicle Nighttime Driver Deaths in Utah")
  > points(obs)
  > lines(fitted, lwd = 2)
  > lines(forecasts, col = "red")
  > par(xpd = NA)
  > graph.param <-
  +     legend("top",
  +            legend = c("observations",expression(estimated~mu[t]),
  +                       expression(predicted~mu[t])),
  +            ncol = 3,
  +            cex = 0.7,
  +            bty = "n", plot = FALSE)
  > legend(graph.param$rect$left,
  +        graph.param$rect$top + graph.param$rect$h,
  +        legend = c("observations", expression(estimated~mu[t]),
  +                   expression(predicted~mu[t])),
  +        col = c("black","black","red"),
  +        lwd = c(1,2,1), lty = c(2,1,1),
  +        pch = c(1, NA_integer_, NA_integer_),
  +        ncol = 3,
  +        cex = 0.7,
  +        bty = "n",
  +        text.font = 4)
  > par(xpd = FALSE)
  > 
  > ### Polio:
  > data(Polio)
  > y  <-  Polio[, 2]
  > X  <-  as.matrix(Polio[, 3:8])
  > glarmamod <- glarma(y, X, thetaLags = c(1,2,5),
  +                     type = "NegBin", method = "FS",
  +                     residuals = "Score", maxit = 100, grad = 1e-6)
  > 
  > summary(glarmamod)
  
  Call: glarma(y = y, X = X, type = "NegBin", method = "FS", residuals = "Score", 
      thetaLags = c(1, 2, 5), maxit = 100, grad = 1e-06)
  
  Score Residuals:
      Min       1Q   Median       3Q      Max  
  -1.0000  -1.0000  -0.2070   0.6889   7.8095  
  
  Negative Binomial Parameter:
        Estimate Std.Error z-ratio Pr(>|z|)   
  alpha    2.856     1.039   2.748    0.006 **
  
  GLARMA Coefficients:
          Estimate Std.Error z-ratio Pr(>|z|)    
  theta_1  0.27621   0.05708   4.840 1.30e-06 ***
  theta_2  0.23545   0.05812   4.051 5.09e-05 ***
  theta_5 -0.00880   0.06265  -0.140    0.888    
  
  Linear Model Coefficients:
                Estimate Std.Error z-ratio Pr(>|z|)    
  Intcpt         0.09423   0.13827   0.681 0.495572    
  Trend         -4.73820   2.65014  -1.788 0.073791 .  
  CosAnnual     -0.04591   0.15980  -0.287 0.773869    
  SinAnnual     -0.58596   0.17514  -3.346 0.000821 ***
  CosSemiAnnual  0.29616   0.12207   2.426 0.015263 *  
  SinSemiAnnual -0.26011   0.13196  -1.971 0.048705 *  
  
      Null deviance: 201.22  on 167  degrees of freedom
  Residual deviance: 255.41  on 158  degrees of freedom
  AIC: 506.134 
  
  Number of Fisher Scoring iterations: 26
  
  LRT and Wald Test:
  Alternative hypothesis: model is a GLARMA process
  Null hypothesis: model is a GLM with the same regression structure
            Statistic  p-value    
  LR Test       21.52 8.20e-05 ***
  Wald Test     31.07 8.21e-07 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
  > 
  > XT1 <- matrix(X[168, ], nrow = 1)
  > mu <- forecast(glarmamod, 1, XT1, offsetT1 = 0)$mu
  > 
  > 
  > plot(0:11, dnbinom(0:11, mu = mu, size =  coef(glarmamod)$NB))
  > 
  > ### Save some values
  > allX <- X
  > allFits <- fitted(glarmamod)
  > ally <- y
  > 
  > ### Look at a succession of forecasts
  > ### Using actual values in forecasts
  > forecasts <- numeric(168)
  > for (i in (156:167)){
  +     y  <-  Polio[1:i, 2]
  +     X  <-  as.matrix(Polio[1:i, 3:8])
  +     glarmamod <- glarma(y, X, thetaLags = c(1,2,5),
  +                         type = "NegBin", method = "FS",
  +                         residuals = "Score", maxit = 100, grad = 1e-6)
  +     XT1 <- matrix(c(1, (i + 1 - 72)/1000, cos(2*pi*(i + 1)/12),
  +              sin(2*pi*(i + 1)/12), cos(4*pi*(i + 1)/12),
  +              sin(4*pi*(i + 1)/12)), nrow = 1)
  +     offsetT1 <- 0
  +     mu <- forecast(glarmamod, 1, XT1, offsetT1)$mu
  +     if (i == 156){
  +         forecasts[1:156] <- fitted(glarmamod)
  +     }
  +     forecasts[i+1] <- mu
  + }
  > par(mfrow = c(1,1))
  > forecasts <- ts(forecasts[157:168], start = c(1983, 1), deltat = 1/12)
  > fitted <- ts(allFits, start = c(1970, 1), deltat = 1/12)
  > data <- ts(Polio$Cases, start = c(1970, 1), deltat = 1/12)
  > plot(data, ylab = "Polio Cases",
  +      main = "Monthly Numbers of Polio Cases in the USA")
  > points(data)
  > lines(fitted, lwd = 2)
  > lines(forecasts, col = "red")
  > par(xpd = NA)
  > graph.param <-
  +     legend("top",
  +            legend = c("observations",expression(estimated~mu[t]),
  +                       expression(predicted~mu[t])),
  +            lty = c(2,1,1),
  +            ncol = 3,
  +            cex = 0.7,
  +            bty = "n", plot = FALSE)
  > legend(graph.param$rect$left,
  +        graph.param$rect$top + graph.param$rect$h,
  +        legend = c("observations", expression(estimated~mu[t]),
  +                   expression(predicted~mu[t])),
  +        col = c("black","black","red"),
  +        lwd = c(1,2,1), lty = c(2,1,1),
  +        pch = c(1, NA_integer_, NA_integer_),
  +        ncol = 3,
  +        cex = 0.7,
  +        bty = "n",
  +        text.font = 4)
  > par(xpd = FALSE)
  > 
  > ### Boat Race:
  > data(OxBoatRace)
  > 
  > y1 <- OxBoatRace$Camwin
  > n1 <- rep(1, length(OxBoatRace$Year))
  > Y <- cbind(y1, n1 - y1)
  > X <- cbind(OxBoatRace$Intercept, OxBoatRace$Diff)
  > colnames(X) <- c("Intercept", "Weight Diff")
  > 
  > oxcamglm <- glm(Y ~ Diff + I(Diff^2),
  +                 data = OxBoatRace,
  +                 family = binomial(link = "logit"), x = TRUE)
  > 
  > X <- oxcamglm$x
  > 
  > glarmamod <- glarma(Y, X, thetaLags = c(1, 2), type = "Bin", method = "NR",
  +                     residuals = "Pearson", maxit = 100, grad = 1e-6)
  > 
  > summary(glarmamod)
  
  Call: glarma(y = Y, X = X, type = "Bin", method = "NR", residuals = "Pearson", 
      thetaLags = c(1, 2), maxit = 100, grad = 1e-06)
  
  Pearson Residuals:
      Min       1Q   Median       3Q      Max  
  -1.8865  -0.7968   0.4395   0.8436   2.9658  
  
  GLARMA Coefficients:
          Estimate Std.Error z-ratio Pr(>|z|)    
  theta_1   0.3396    0.1709   1.987 0.046975 *  
  theta_2   0.5552    0.1459   3.804 0.000142 ***
  
  Linear Model Coefficients:
               Estimate Std.Error z-ratio Pr(>|z|)   
  (Intercept)  0.349954  0.267249   1.309  0.19038   
  Diff         0.114755  0.038238   3.001  0.00269 **
  I(Diff^2)   -0.011333  0.004916  -2.305  0.02114 * 
  
      Null deviance: 216.16  on 155  degrees of freedom
  Residual deviance: 148.53  on 151  degrees of freedom
  AIC: 193.0553 
  
  Number of Newton Raphson iterations: 5
  
  LRT and Wald Test:
  Alternative hypothesis: model is a GLARMA process
  Null hypothesis: model is a GLM with the same regression structure
            Statistic  p-value    
  LR Test       15.56 0.000417 ***
  Wald Test     17.51 0.000158 ***
  ---
  Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  
  > 
  > par(mfrow=c(3,2))
  > plot(glarmamod)
  > 
  > 
  > ### Save some values
  > allX <- X
  > allFits <- fitted(glarmamod)
  > ally1 <- y1
  > 
  > 
  > ### Look at a succession of forecasts
  > ### Using actual values in forecasts
  > forecasts <- numeric(156)
  > for (i in (146:155)){
  +     y1 <- OxBoatRace$Camwin[1:i]
  +     n1 <- rep(1, i)
  +     Y <- cbind(y1, n1 - y1)
  +     X <- cbind(OxBoatRace$Intercept[1:i], OxBoatRace$Diff[1:i])
  +     colnames(X) <- c("Intercept", "Weight Diff")
  + 
  +     oxcamglm <- glm(Y ~ Diff + I(Diff^2),
  +                     data = OxBoatRace[1:i, ],
  +                     family = binomial(link = "logit"), x = TRUE)
  +     X <- oxcamglm$x
  +     glarmamod <- glarma(Y, X, thetaLags = c(1, 2), type = "Bin", method = "NR",
  +                     residuals = "Pearson", maxit = 100, grad = 1e-6)
  +     XT1 <- matrix(allX[i + 1, ], nrow = 1)
  +     offsetT1 <- 0
  +     mu <- forecast(glarmamod, 1, XT1, offsetT1 = 0)$mu
  +     if (i == 146){
  +         forecasts[1:146] <- fitted(glarmamod)
  +     }
  +     forecasts[i+1] <- mu
  + }
  > 
  > par(mfrow = c(1,1))
  > forecasts <- zoo(forecasts, order.by = OxBoatRace$Year)
  Error in zoo(forecasts, order.by = OxBoatRace$Year) : 
    could not find function "zoo"
  Execution halted