NSBJF
Computes Box‑Jenkins forecasts and their associated probability limits for a nonseasonal ARMA model.
Required Arguments
W — Vector of length NOBS containing the time series. (Input)
PAR — Vector of length NPAR containing the autoregressive parameters. (Input)
LAGAR — Vector of length NPAR containing the order of the autoregressive parameters. (Input)
The elements of LAGAR must be greater than zero.
PMA — Vector of length NPMA containing the moving average parameters. (Input)
LAGMA — Vector of length NPMA containing the order of the moving average parameters. (Input)
The elements of LAGMA must be greater than zero.
ICNST — Option for including the overall constant in the model. (Input)
ICNST | Action |
|---|
0 | No overall constant is included. |
1 | The overall constant is included. |
CNST — Estimate of the overall constant. (Input)
AVAR — Estimate of the random shock variance. (Input)
AVAR must be greater than 0.
ALPHA — Value in the exclusive interval (0, 1) used to specify the 100(1 ‑ ALPHA)% probability limits of the forecasts. (Input)
Typical choices for ALPHA are 0.10, 0.05, and 0.01.
MXBKOR — Maximum backward origin. (Input)
MXBKOR must be greater than or equal to zero and less than or equal to
NOBS ‑ max(MAXAR, MAXMA) where MAXAR = max(LAGAR(i)) and MAXMA = max(LAGMA(j)). Forecasts at origins NOBS ‑ MXBKOR through NOBS are generated.
MXLEAD — Maximum lead time for forecasts. (Input)
MXLEAD must be greater than zero.
FCST — MXLEAD by (MXBKOR + 3) matrix defined below. (Output)
Column | Content |
|---|
j | Forecasts for lead times l = 1, …, MXLEAD at origins NOBS ‑ MXBKOR ‑ 1 + j, j = 1, …, MXBKOR + 1. |
MXBKOR + 2 | Deviations from each forecast that give the 100(1 ‑ ALPHA)% probability limits. |
MXBKOR + 3 | Psi weights of the infinite order moving average form of the model. |
Optional Arguments
NOBS — Number of observations in the time series W. (Input)
NOBS must be greater than ICONST + max(LAGAR(i)) + max(LAGMA(j)).
Default: NOBS = size (W,1).
IPRINT — Printing option. (Input)
Default: IPRINT = 0.
IPRINT | Action |
|---|
0 | No printing is performed. |
1 | Prints the forecasts for lead times l = 1, …, MXLEAD at each origin t = (NOBS ‑ MXBKOR ), …, NOBS, the 100(1 ‑ ALPHA)% probability limit deviations, and the psi weights. |
NPAR — Number of autoregressive parameters. (Input)
NPAR must be greater than or equal to zero.
Default: NPAR = size (PAR,1).
NPMA — Number of moving average parameters. (Input)
NPMA must be greater than or equal to zero.
Default: NPMA = size (PMA,1).
LDFCST — Leading dimension of FCST exactly as specified in the dimension statement in the calling program. (Input)
LDFCST must be greater than or equal to MXLEAD.
Default: LDFCST = size (FCST,1).
FORTRAN 90 Interface
Generic: CALL NSBJF (W, PAR, LAGAR, PMA, LAGMA, ICNST, CNST, AVAR, ALPHA, MXBKOR, MXLEAD, FCST [, …])
Specific: The specific interface names are S_NSBJF and D_NSBJF.
FORTRAN 77 Interface
Single: CALL NSBJF (NOBS, W, IPRINT, NPAR, PAR, LAGAR, NPMA, PMA, LAGMA, ICNST, CNST, AVAR, ALPHA, MXBKOR, MXLEAD, FCST, LDFCST)
Double: The double precision name is DNSBJF.
Description
Routine NSBJF computes Box‑Jenkins forecasts and their associated probability limits for a nonseasonal ARMA model given a sample of n = NOBS observations {Wt} for t = 1, 2, …, n.
Suppose the time series {Wt} is generated by a nonseasonal ARMA model of the form
where B is the backward shift operator, θ0 = CONST,
p = NPAR and q = NPMA. Without loss of generality, we assume
so that the nonseasonal ARMA model is of order (pʹ, qʹ) where pʹ = lɸ (p) and qʹ = lθ (q). Note that the usual hierarchal model assumes
The Box‑Jenkins forecast at origin t for lead time l of Wt+l is defined in terms of the difference equation
where
The 100(1 ‑ α)% probability limits for Wt+l are given by
where za/2 is the 100(1 ‑ α/2) percentile of the standard normal distribution,
and { Ψj } are the parameters of the random shock form of the difference equation. Note that the forecasts are computed for lead times l = 1, 2, …, L at origins t = (n ‑ b), (n ‑ b + 1), …, n where L = MXLEAD and b = MXBKOR.
The Box‑Jenkins forecasts minimize the mean square error
Also, the forecasts may be easily updated according to the equation
This approach and others are given in Chapter 5 of Box and Jenkins (1976).
Comments
1. Workspace may be explicitly provided, if desired, by use of N2BJF/DN2BJF. The reference is:
CALL N2BJF (NOBS, W, IPRINT, NPAR, PAR, LAGAR, NPMA, PMA, LAGMA, ICNST, CNST, AVAR, ALPHA, MXBKOR, MXLEAD, FCST, LDFCST, PARH, PMAH, PSIH, PSI, LAGPSI)
The additional arguments are as follows:
PARH — Work vector of length equal to IARDEG + 1.
PMAH — Work vector of length equal to IMADEG + 1.
PSIH — Work vector of length equal to MXLEAD + 1.
PSI — Work vector of length equal to MXLEAD + 1.
LAGPSI — Work vector of length equal to MXLEAD + 1.
2. If the
W series has been transformed using a Box‑Cox transformation with parameters
POWER and
SHIFT, the forecasts and probability limits for the original series may be obtained by application of routine
BCTR with the same parameters and argument
IDIR set equal to one.
Example
Consider the Wölfer Sunspot Data (Anderson 1971, page 660) consisting of the number of sunspots observed each year from 1749 through 1924. The data set for this example consists of the number of sunspots observed from 1770 through 1869. Routine
NSBJF is used to computed forecasts and 95% probability limits for the forecasts for an
ARMA(2, 1) model fit using routine
NSPE. With
MXBKOR = 3, columns one through four of
FCST give forecasts given the data through 1866, 1867, 1868, and 1869, respectively. Column 5 gives the deviations from the forecast for computing probability limits, and column six gives the psi weights, which can be used to update forecasts once more data is available. For example, the forecast for the 102‑nd observation (year 1871) given the data through the 100‑th observation (year 1869) is 77.21, and 95% probability limits are given by 77.21
∓ 56.30. After observation 101 (
W101 for year 1870) is available, the forecast can be updated by using equation 7 with the psi weight (
Ψ1 = 1.37) and the one‑step‑ahead forecast error for observation 101 (
W101 ‑ 83.72) to give
77.21 + 1.37 (W101 − 83.72)
Since this updated forecast is one step ahead, the 95% probability limits are now given by the forecast ∓ 33.22.
USE GDATA_INT
USE NSPE_INT
USE NSBJF_INT
USE WRRRL_INT
IMPLICIT NONE
INTEGER LDFCST, MXBKOR, MXLEAD, NOBS, NPAR, NPMA
PARAMETER (MXBKOR=3, MXLEAD=12, NOBS=100, NPAR=2, NPMA=1, &
LDFCST=MXLEAD)
!
INTEGER ICNST, LAGAR(NPAR), LAGMA(NPMA), NCOL, NROW
REAL ALPHA, AVAR, CNST, FCST(LDFCST,MXBKOR+3), PAR(NPAR), &
PMA(NPMA), RDATA(176,2), W(NOBS), WMEAN
CHARACTER CLABEL(MXBKOR+4)*40, RLABEL(1)*6
!
EQUIVALENCE (W(1), RDATA(22,2))
!
DATA LAGAR(1), LAGAR(2)/1, 2/
DATA LAGMA(1)/1/
DATA RLABEL/'NUMBER'/, CLABEL/'%/Lead%/Time', &
'%/Forecast%/From 1866', '%/Forecast%/From 1867', &
'%/Forecast%/From 1868', '%/Forecast%/From 1869', &
' Deviation %/ for 95% %/Prob. Limits', '%/%/Psi'/
! Wolfer Sunspot Data for
! years 1770 through 1869
CALL GDATA (2, RDATA, NROW, NCOL)
! Compute preliminary parameter
! estimates for ARMA(2,1) model
CALL NSPE (W, CNST, PAR, PMA, AVAR)
!
! Include constant in forecast model
ICNST = 1
! Specify 95 percent probability
! limits for forecasts
ALPHA = 0.05
! Compute forecasts
CALL NSBJF (W, PAR, LAGAR, PMA, LAGMA, &
ICNST, CNST, AVAR, ALPHA, MXBKOR, MXLEAD, FCST)
! Print results
CALL WRRRL ('FCST', FCST, RLABEL, CLABEL, FMT='(5F9.2, F6.3)')
!
END
Output
FCST
Deviation
Lead Forecast Forecast Forecast Forecast for 95%
Time From 1866 From 1867 From 1868 From 1869 Prob. Limits Psi
1 18.28 16.62 55.19 83.72 33.22 1.368
2 28.92 32.02 62.76 77.21 56.30 1.127
3 41.01 45.83 61.89 63.46 67.62 0.616
4 49.94 54.15 56.46 50.10 70.64 0.118
5 54.09 56.56 50.19 41.38 70.75 -0.208
6 54.13 54.78 45.53 38.22 71.09 -0.326
7 51.78 51.17 43.32 39.30 71.91 -0.286
8 48.84 47.71 43.26 42.46 72.53 -0.169
9 46.53 45.47 44.46 45.77 72.75 -0.045
10 45.35 44.69 45.98 48.08 72.77 0.041
11 45.21 44.99 47.18 49.04 72.78 0.077
12 45.71 45.82 47.81 48.91 72.82 0.072
