Chapter 8: Time Series Analysis and Forecasting

TS_OUTLIER_FORECAST

Computes forecasts, associated probability limits and  weights for an outlier contaminated time series.

Required Arguments

W —  Array of length NOBS containing the outlier free time series.  (Input)

RESIDUAL — Array of length NOBS containing the residuals of the outlier free time series determined from routine TS_OUTLIER_IDENTIFICATION. (Input)

NOUTLIERS — The number of outliers in W, determined from routine TS_OUTLIER_IDENTIFICATION.  (Input)

IOUTLIERSTATS Array of size NOUTLIERS by 2 containing outlier statistics from routine TS_OUTLIER_IDENTIFICATION.  (Input).
The first column contains the time at which the outlier was observed () and the second column contains an identifier indicating the type of outlier observed. Outlier types fall into one of five categories:

IOUTLIERSTATS

Category

0

Innovational Outliers (IO)

1

Additive Outliers  (AO)

2

Level Shift Outliers (LS)

3

Temporary Change Outliers (TC)

4

Unable to Identify (UI)

            If NOUTLIERS = 0 this array is ignored.

OMEGA — Array of length NOUTLIERS containing the omega weights for the outliers determined through routine TS_OUTLIER_IDENTIFICATION.  (Input)

DELTA — The dynamic dampening effect parameter used in the outlier detection,
 0.0 < DELTA < 1.0.  (Input)

MODEL Array of length four containing estimates for p, q, s and d in MODEL(1), MODEL(2), MODEL(3) and MODEL(4), respectively.  (Input)

PARAMS Array of length 1+p+q containing the estimated constant, AR and MA parameters as output from routine TS_OUTLIER_IDENTIFICATION.  (Input)

MXLEAD — Maximum lead time for forecasts.  (Input)
The forecasts are taken at origin , the time point of the last observed value in the series, for lead times 1, 2,…, MXLEAD.
MXLEAD must be greater than zero.

FCST An array of size MXLEAD by 3 containing the forecasted values for the outlier contaminated  series in the first column. The second column contains the deviations from each forecast that give the 100(1-ALPHA)% probability limits, and the third column contains the  weights of the infinite order moving average form of the model.   (Output)

Optional Arguments

NOBS —  Number of observations in the time series W.  (Input)
Default: NOBS = size (W).

ALPHA — Value in the exclusive interval (0,1) used to specify the 100(1-ALPHA)% probability limits of the forecast.  (Input)
Default: ALPHA = 0.05

OUTFREEFCST — An MXLEAD by 3 array containing the forecasted values for the original outlier free series in the first column. The second column contains standard errors for these forecasts, and the third column contains the  weights of the infinite order moving average form of the model.  (Output)

FORTRAN 90 Interface

Generic:          CALL TS_OUTLIER_FORECAST (W, RESIDUAL, NOUTLIERS, IOUTLIERSTATS, OMEGA, DELTA, MODEL, PARAMS, MXLEAD, FCST [,…])

Specific:         The specific interface names are S_TS_OUTLIER_FORECAST and D_TS_OUTLIER_FORECAST.

Description

Consider the following model for a given outlier contaminated univariate time series :

 

For an explanation of the notation, see the “Description” section for TS_OUTLIER_IDENTIFICATION.  It follows from the formula above that the Box-Jenkins forecast at origin  for lead time , , can be computed as:

Therefore, computation of the forecasts for  is done in two steps:

1.       Computation of  the forecasts for the outlier free series .

2.     Computation of  the forecasts for the original series  by adding the multiple outlier effects to the forecasts for  .

Step 1 above:

Since

where

the Box-Jenkins forecast at origin  for lead time , , can be computed recursively as:

Here,

 

and

 

Step 2 above:        

The formulas for  for the different types of outliers are as follows:

Innovational outliers (IO)

Additive outliers (AO)

Level shifts (LS)

Temporary changes (TC)

Assuming the outlier occurs at time point , the outlier impact is therefore:

Innovational outliers (IO)

Additive outliers (AO)

Level shifts (LS)

Temporary changes (TC)

 

From these formulas, the forecasts  can be computed easily.

The  percent probability limits for  and   are given by

where  is the  percentile of the standard normal distribution,  is an estimate of the variance  of the random shocks (returned from routine TS_OUTLIER_IDENTIFICATION), and the weights are the coefficients in

For a detailed explanation of these concepts, see Chapter 5: “Forecasting”, Box, Jenkins and Reinsel (1994).

Example

This example is a realization of an ARMA(2,1) process described by the model , a Gaussian white noise process.

Outliers were artificially added to the outlier free series  at time points (level shift, ) and (additive outlier, ), resulting in the outlier contaminated series .  For both series, forecasts were determined for time points  and compared with the actual values of the series.

 

USE TS_OUTLIER_IDENTIFICATION_INT

USE TS_OUTLIER_FORECAST_INT

USE WRRRL_INT

USE UMACH_INT

 

IMPLICIT NONE

 

              real(kind(1e0)), dimension(290) :: w = (/ &

                41.6699982,41.6699982,42.0752144,42.6123962,43.6161919,42.1932831, &

                43.1055450,44.3518715,45.3961258,45.0790215,41.8874397,40.2159805, &

                40.2447319,39.6208458,38.6873589,37.9272423,36.8718872,36.8310852, &

                37.4524879,37.3440933,37.9861374,40.3810501,41.3464622,42.6495285, &

                42.6096764,40.3134537,39.7971268,41.5401535,40.7160759,41.0363541, &

                41.8171883,42.4190292,43.0318832,43.9968109,44.0419617,44.3225212, &

                44.6082611,43.2199631,42.0419197,41.9679718,42.4926224,43.2091255, &

                43.2512283,41.2301674,40.1057358,40.4510574,41.5329170,41.5678177, &

                43.0090141,42.1592140,39.9234505,38.8394127,40.4319878,40.8679352, &

                41.4551926,41.9756317,43.9878922,46.5736389,45.5939293,42.4487762, &

                41.5325394,42.8830910,44.5771217,45.8541985,46.8249474,47.5686378, &

                46.6700745,45.4120026,43.2305107,42.7635345,43.7112923,42.0768661, &

                41.1835632,40.3352280,37.9761467,35.9550056,36.3212509,36.9925880, &

                37.2625008,37.0040665,38.5232544,39.4119797,41.8316803,43.7091446, &

                42.9381447,42.1066780,40.3771248,38.6518707,37.0550499,36.9447708, &

                38.1017685,39.4727097,39.8670387,39.3820763,38.2180786,37.7543488, &

                37.7265244,38.0290642,37.5531158,37.4685936,39.8233147,42.0480766, &

                42.4053535,43.0117416,44.1289330,45.0393829,45.1114540,45.0086479, &

                44.6560631,45.0278931,46.7830849,48.7649765,47.7991905,46.5339661, &

                43.3679199,41.6420822,41.2694893,41.5959740,43.5330009,43.3643608, &

                42.1471291,42.5552788,42.4521446,41.7629128,39.9476891,38.3217010, &

                40.5318718,42.8811569,44.4796944,44.6887932,43.1670265,41.2226143, &

                41.8330154,44.3721924,45.2697029,44.4174194,43.5068550,44.9793015, &

                45.0585403,43.2746620,40.3317070,40.3880501,40.2627106,39.6230278, &

                41.0305252,40.9262009,40.8326912,41.7084885,42.9038048,45.8650513, &

                46.5231590,47.9916115,47.8463135,46.5921936,45.8854408,45.9130440, &

                45.7450371,46.2964249,44.9394569,45.8141251,47.5284042,48.5527802, &

                48.3950577,47.8753052,45.8880005,45.7086983,44.6174774,43.5567932, &

                44.5891113,43.1778679,40.9405632,40.6206894,41.3330421,42.2759552, &

                42.4744949,43.0719833,44.2178459,43.8956337,44.1033440,45.6241455, &

                45.3724861,44.9167595,45.9180603,46.9077835,46.1666603,46.6013489, &

                46.6592331,46.7291603,47.1908340,45.9784355,45.1215782,45.6791115, &

                46.7379875,47.3036957,45.9968834,44.4669495,45.7734680,44.6315041, &

                42.9911766,46.3842583,43.7214432,43.5276833,41.3946495,39.7013168, &

                39.1033401,38.5292892,41.0096245,43.4535828,44.6525154,45.5725899, &

                46.2815285,45.2766647,45.3481712,45.5039482,45.6745682,44.0144806, &

                42.9305000,43.6785469,42.2500534,40.0007210,40.4477005,41.4432716, &

                42.0058670,42.9357758,45.6758842,46.8809929,46.8601494,47.0449791, &

                46.5420647,46.8939934,46.2963371,43.5479164,41.3864059,41.4046364, &

                42.3037987,43.6223717,45.8602371,47.3016396,46.8632469,45.4651413, &

                45.6275482,44.9968376,42.7558670,42.0218239,41.9883728,42.2571678, &

                44.3708687,45.7483635,44.8832512,44.7945862,44.8922577,44.7409401, &

                45.1726494,45.5686874,45.9946709,47.3151054,48.0654068,46.4817467, &

                42.8618279,42.4550323,42.5791168,43.4230957,44.7787971,43.8317108, &

                43.6481781,42.4183960,41.8426285,43.3475227,44.4749908,46.3498306, &

                47.8599319,46.2449913,43.6044006,42.4563484,41.2715340,39.8492508, &

                39.9997292,41.4410820,42.9388237,42.5687332,42.6384087,41.7088661, &

                43.9399033,45.4284401,44.4558411,45.1761856,45.3489113,45.1892662, &

                46.3754730,45.6082802 /)

   

              integer :: i, nout

integer, dimension(4) :: model

integer :: noutliers, nobs=280, mxlead=10

integer, dimension(:,:), pointer :: outlierstat

real(kind(1e0)), dimension(:), pointer :: omega

real(kind(1e0)) :: delta = 0.7, res_sigma, aic

real(kind(1e0)), dimension(280) :: x, residual

real(kind(1e0)), dimension(10,4) :: forecast_table

real(kind(1e0)), dimension(10,3) :: fcst, outfreefcst

real(kind(1e0)), dimension(4) :: params

character (len = 10) :: fmt

character (len = 20), dimension(5) :: clabel

character (len = 4), dimension(10) :: rlabel

 

             fmt = '(F11.4)'

             clabel(1) = '   '

             clabel(2) = 'Orig. series'

             clabel(3) = 'forecast'

             clabel(4) = 'prob. limits'

             clabel(5) = 'psi weights'

 

             rlabel = (/ ' 1',' 2',' 3',' 4',' 5',' 6',' 7',' 8',' 9','10' /)

 

             model = (/ 2,1,1,0 /)

 

             CALL TS_OUTLIER_IDENTIFICATION (model, w, x, NOBS=nobs, DELTA=delta,  &

                        RELERR=1.0e-5, RESIDUAL=residual,  &

                        RESSIGMA=res_sigma, NOUTLIERS=noutliers, &

                        IOUTLIERSTATS=outlierstat, OMEGA=omega, &

                        PARAMS=params, AIC=aic)

 

             CALL UMACH(2, nout)

             WRITE (nout,*) 'ARMA parameters:'

             DO i=1, 1+model(1)+model(2)

               WRITE (nout,*) i,'    ', params(i)

             END DO

 

             WRITE (nout,*)

             WRITE (nout,*) 'Number of outliers: ', noutliers

             WRITE (nout,*)

             WRITE (nout,*) 'Outlier statistics: '

             WRITE (nout,FMT="(T2, A, T22, A)") 'Time point','Outlier type'

             DO i=1,noutliers

               WRITE (nout,FMT="(T2, I10, T22, I12)") outlierstat(i,1), outlierstat(i,2)

             END DO

 

             WRITE (nout,*)

             WRITE (nout,*) 'RSE: ', res_sigma

             WRITE (nout,*) 'AIC: ', aic

             WRITE (nout,*)

 

             CALL TS_OUTLIER_FORECAST (x, residual, noutliers, outlierstat, &

                            omega, delta, model, params, mxlead, fcst,  &

                            NOBS=nobs, OUTFREEFCST=outfreefcst)

 

             forecast_table(1:mxlead,1) = w(281:290)

             forecast_table(1:mxlead,2:4) = fcst(1:mxlead,1:3)

 

             CALL WRRRL ('* * * Forecast Table for outlier contaminated series * * *',&

                          forecast_table, rlabel, clabel, FMT=fmt)

 

             forecast_table(1:mxlead,1) = w(281:290) - 2.5

             forecast_table(1:mxlead,2:4) = outfreefcst(1:mxlead,1:3)

 

             WRITE (nout,*)

 

             clabel(2) = 'Outlier free series'

             CALL WRRRL ('* * * Forecast Table for outlier free series * * *', &

                          forecast_table, RLABEL, CLABEL, FMT=fmt)

             END

             

Output

 

ARMA parameters:

               1      8.837544

               2      0.9461826

               3      -0.1512835

               4      -0.5606939

 

              Number of outliers:  2

 

              Outlier statistics:

              Time point          Outlier type

                     150                     2

                     200                     1

 

              RSE:  1.0042976

              AIC:  1323.6127

 

             * * * Forecast Table for outlier contaminated series * * *

                  Orig. series     forecast  prob. limits  psi weights

               1       42.6384      42.3113        1.9684       1.5069

               2       41.7089      42.7868        3.5598       1.2745

               3       43.9399      43.2756        4.3550       0.9779

               4       45.4284      43.6662        4.7615       0.7325

               5       44.4558      43.9618        4.9750       0.5451

               6       45.1762      44.1825        5.0894       0.4050

               7       45.3489      44.3465        5.1514       0.3007

               8       45.1893      44.4683        5.1853       0.2233

               9       46.3755      44.5588        5.2039       0.1658

              10       45.6083      44.6259        5.2141       0.1231

 

             * * * Forecast Table for outlier free series * * *

                 Outlier free     forecast  prob. limits  psi weights

                       series

              1       40.1384      40.5805        1.9684       1.5069

              2       39.2089      41.0560        3.5598       1.2745

              3       41.4399      41.5449        4.3550       0.9779

              4       42.9284      41.9355        4.7615       0.7325

              5       41.9558      42.2311        4.9750       0.5451

              6       42.6762      42.4517        5.0894       0.4050

              7       42.8489      42.6158        5.1514       0.3007

              8       42.6893      42.7376        5.1853       0.2233

              9       43.8755      42.8281        5.2039       0.1658

             10       43.1083      42.8952        5.2141       0.1231



http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260