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 t for lead time 
, 
, can be computed as:
Therefore, computation of the forecasts for 
is done in two steps:
-
Computation of the forecasts for the outlier free series
. -
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 Lj(B) for the different types of outliers are as follows:
|
|
|
|
Additive outliers (AO) |
|
|
Level shifts (LS) |
|
|
Temporary changes (TC) |
|
Assuming the outlier occurs at time point tj, the outlier impact is therefore:
|
|
|
|
|
|
|
|
|
|
|
|
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
