Chapter 8: Time Series Analysis and Forecasting

TS_OUTLIER_IDENTIFICATION

TextonSpedometerHPClogo.gif



Detects and determines outliers and simultaneously estimates the model parameters in a time series whose underlying outlier free series follows a general seasonal or nonseasonal ARMA model.

Required Arguments

MODEL —   Array of length 4 containing the order of the ARIMA model the outlier free series is following. Specifically, MODEL(1)= p, MODEL(2) = q, MODEL(3) = s, MODEL(4) = d. (Input)

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

X —   Array of length NOBS containing the outlier free series.   (Output)

Optional Arguments

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

DELTA — The dynamic dampening effect parameter used in the detection of a Temporary Change Outlier (TC), 0.0 < DELTA < 1.0.   (Input)
Default: DELTA = 0.7 .

CRITICAL —  Critical value used as a threshold for outlier detection, CRITICAL > 0.   (Input)
Default: CRITICAL = 3.0.

EPSILON — Positive tolerance value controlling the accuracy of parameter estimates during outlier detection.   (Input)
Default: EPSILON = 0.001.

RELERR — Stopping criterion for use in the nonlinear equation solver used by NSPE, see routine NSPE for more details.   (Input)
Default: RELERR = 1.0e-10.

TOLSS — Tolerance level used to determine convergence of the nonlinear least-squares algorithm used by NSLSE, see routine NSLSE for more details.   (Input)
TOLSS must be greater than zero and less than one.
Default: TOLSS = 0.9 Χ AMACH(4).

RESIDUAL — Array of length NOBS containing the residuals for the outlier free series.   (Output)

RESSIGMA — Residual standard error of the outlier free series.   (Output)

NOUTLIERS — The number of outliers detected.   (Output)

IOUTLIERSTATS — Pointer to an array of size NOUTLIERS by 2 containing outlier statistics. 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 no outliers are detected, an array of size 0 is returned.   (Output)

TAUSTAT — Pointer to an array of length NOUTLIERS containing the t value for each detected outlier.   (Output)

OMEGA — Pointer to an array of length NOUTLIERS containing the computed weights for the detected outliers.   (Output)

PARAMS — Array of length 1+p+q containing the estimated constant, AR and MA parameters, respectively.   (Output)

AIC — Akaike's Information Criterion (AIC) for the fitted model.   (Output)

AICC— Akaike's Corrected Information Criterion (AICC) for the fitted model.   (Output)

BIC — Bayesian Information Criterion (BIC) for the fitted model.   (Output)

FORTRAN 90 Interface

Generic:                              CALL TS_OUTLIER_IDENTIFICATION (MODEL, W, X [,…])

Specific:                             The specific interface names are S_TS_OUTLIER_IDENTIFICATION and D_TS_OUTLIER_IDENTIFICATION.

Description

Consider a univariate time series that can be described by the following multiplicative seasonal ARIMA model of order :

Here, ,  .  is the lag operator, ,  is a white noise process, and  denotes the mean of the series .

In general,  is not directly observable due to the influence of outliers. Chen and Liu (1993) distinguish between four types of outliers: innovational outliers (IO), additive outliers (AO), temporary changes (TC)  and level shifts (LS). If an outlier occurs as the last observation of the series, then Chen and Liu’s algorithm is unable to determine the outlier’s classification. In TS_OUTLIER_IDENTIFICATION, such an outlier is called a UI (unable to identify) and is treated as an innovational outlier.

In order to take the effects of multiple outliers occurring at time points  into account, Chen and Liu consider the following model:

Here,  is the observed outlier contaminated series, and  and  denote the magnitude and dynamic pattern of outlier , respectively.   is an indicator function that determines the temporal course of the outlier effect, ,  otherwise. Note that  operates on  via .

The last formula shows that the outlier free series  can be obtained from the original series  by removing all occurring outlier effects:

.

The different types of outliers are characterized by different values for :

1.   for an innovational outlier,

2.   for an additive outlier,

3.   for a level shift outlier and

4.   for a temporary change outlier.

TS_OUTLIER_IDENTIFICATION is an implementation of Chen and Liu’s algorithm. It determines the coefficients in  and the outlier effects in the model for the observed series jointly in three stages. The magnitude of the outlier effects is determined by least squares estimates. Outlier detection itself is realized by examination of the maximum value of the standardized statistics of the outlier effects. For a detailed description, see Chen and Liu’s original paper (1993).

Intermediate and final estimates for the coefficients in  and  are computed by routines NSPE and NSLSE.  If the roots of or  lie on or within the unit circle, then the algorithm stops with an appropriate error message. In this case, different values for p and q should be tried.

Examples

Example 1

This example is based on estimates of the Canadian lynx population. TS_OUTLIER_IDENTIFICATION is used to fit an ARIMA(2,2,0) model of the form, , Gaussian White noise, to the given series. TS_OUTLIER_IDENTIFICATION computes parameters  and  and identifies a LS outlier at time point .

 

      use umach_int

      use ts_outlier_identification_int

 

      implicit none

 

      integer :: i, nout

      integer :: noutliers

      integer, dimension(4) :: model

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

      real(kind(1e0)) :: ressigma, aic

      real(kind(1e0)), dimension(3) :: parameters

      real(kind(1e0)), dimension(114) :: w, x

 

      w = (/ 0.24300E01,0.25060E01,0.27670E01,0.29400E01,0.31690E01,&

             0.34500E01,0.35940E01,0.37740E01,0.36950E01,0.34110E01,&

             0.27180E01,0.19910E01,0.22650E01,0.24460E01,0.26120E01,&

             0.33590E01,0.34290E01,0.35330E01,0.32610E01,0.26120E01,&

             0.21790E01,0.16530E01,0.18320E01,0.23280E01,0.27370E01,&

             0.30140E01,0.33280E01,0.34040E01,0.29810E01,0.25570E01,&

             0.25760E01,0.23520E01,0.25560E01,0.28640E01,0.32140E01,&

             0.34350E01,0.34580E01,0.33260E01,0.28350E01,0.24760E01,&

             0.23730E01,0.23890E01,0.27420E01,0.32100E01,0.35200E01,&

             0.38280E01,0.36280E01,0.28370E01,0.24060E01,0.26750E01,&

             0.25540E01,0.28940E01,0.32020E01,0.32240E01,0.33520E01,&

             0.31540E01,0.28780E01,0.24760E01,0.23030E01,0.23600E01,&

             0.26710E01,0.28670E01,0.33100E01,0.34490E01,0.36460E01,&

             0.34000E01,0.25900E01,0.18630E01,0.15810E01,0.16900E01,&

             0.17710E01,0.22740E01,0.25760E01,0.31110E01,0.36050E01,&

             0.35430E01,0.27690E01,0.20210E01,0.21850E01,0.25880E01,&

             0.28800E01,0.31150E01,0.35400E01,0.38450E01,0.38000E01,&

             0.35790E01,0.32640E01,0.25380E01,0.25820E01,0.29070E01,&

             0.31420E01,0.34330E01,0.35800E01,0.34900E01,0.34750E01,&

             0.35790E01,0.28290E01,0.19090E01,0.19030E01,0.20330E01,&

             0.23600E01,0.26010E01,0.30540E01,0.33860E01,0.35530E01,&

             0.34680E01,0.31870E01,0.27230E01,0.26860E01,0.28210E01,&

             0.30000E01,0.32010E01,0.34240E01,0.35310E01 /)

 

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

 

      call ts_outlier_identification(model, w, x, CRITICAL=3.5, &

                 RESSIGMA=ressigma, NOUTLIERS=noutliers,&

                 IOUTLIERSTATS=outlierstat, PARAMS=parameters,&

                 AIC=aic)

 

      call umach(2, nout)

      write(nout,FMT="(T2,A,I2)") 'Number of outliers:', noutliers

      write(nout,FMT="(T2,A)") 'Outlier statistics:'

      write(nout,FMT="(T4,A,TR10,A)") 'Time point','Outlier type'

      write(nout,FMT="(I8,TR20,I2)") (outlierstat(i,:),i=1,noutliers)

      write(nout,FMT="(/,T2,A)") 'ARMA parameters:'

      write(nout,FMT="(T3,I2,TR5,f10.6)") (i,parameters(i),i=1,3)

      write(nout,FMT="(/,T2,A,f10.6)") 'RSE: ', ressigma

      write(nout,FMT="(T2,A,f10.6)") 'AIC: ', aic

      write(nout,FMT="(/,T2,A)") 'Extract from the series:'

      write(nout,FMT="(T2,A,TR6,A,TR6,A)") 'time point',&

                           'original series', 'outlier free series'

      do i=1,36

        write(nout, FMT="(T5,I2,T15,F15.6,T35,F19.6)") i, w(i), x(i)

      end do

      end

Output

 

 Number of outliers: 1

 Outlier statistics:

   Time point          Outlier type

      16                     2

 

 ARMA parameters:

   1       0.000000

   2       0.124390

   3      -0.179959

 

 RSE:   0.319650

 AIC: 282.995575

 

 Extract from the series:

 time point      original series      outlier free series

     1               2.430000                2.430000

     2               2.506000                2.506000

     3               2.767000                2.767000

     4               2.940000                2.940000

     5               3.169000                3.169000

     6               3.450000                3.450000

     7               3.594000                3.594000

     8               3.774000                3.774000

     9               3.695000                3.695000

    10               3.411000                3.411000

    11               2.718000                2.718000

    12               1.991000                1.991000

    13               2.265000                2.265000

    14               2.446000                2.446000

    15               2.612000                2.612000

    16               3.359000                2.701984

    17               3.429000                2.771984

    18               3.533000                2.875984

    19               3.261000                2.603984

    20               2.612000                1.954984

    21               2.179000                1.521984

    22               1.653000                0.995984

    23               1.832000                1.174984

    24               2.328000                1.670984

    25               2.737000                2.079984

    26               3.014000                2.356984

    27               3.328000                2.670984

    28               3.404000                2.746984

    29               2.981000                2.323984

    30               2.557000                1.899984

    31               2.576000                1.918984

    32               2.352000                1.694984

    33               2.556000                1.898984

    34               2.864000                2.206984

    35               3.214000                2.556984

    36               3.435000                2.777984

 

Example 2

This example is an artificial realization of an ARMA(1,1) process via formula  Gaussian white noise, .

An additive outlier with  was added at time point , a temporary change outlier with  was added at time point .

 

      use umach_int

      use ts_outlier_identification_int

 

      implicit none

 

      integer :: i, nout

      integer :: noutliers

      integer, dimension(4) :: model

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

      real(kind(1e0)) :: ressigma, aic

      real(kind(1e0)), dimension(3) :: parameters

      real(kind(1e0)), dimension(300) :: w, x

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

 

      w = (/ 50.0000000,50.2728081,50.6242599,51.0373917,51.9317627,&

             50.3494759,51.6597252,52.7004929,53.5499802,53.1673279,&

             50.2373505,49.3373871,49.5516472,48.6692696,47.6606636,&

             46.8774185,45.7315445,45.6469727,45.9882355,45.5216560,&

             46.0479660,48.1958656,48.6387749,49.9055367,49.8077278,&

             47.7858467,47.9386749,49.7691956,48.5425873,49.1239853,&

             49.8518791,50.3320694,50.9146347,51.8772049,51.8745689,&

             52.3394470,52.7273712,51.4310036,50.6727448,50.8370399,&

             51.2843437,51.8162918,51.6933670,49.7038231,49.0189247,&

             49.455703,50.2718010,49.9605980,51.3775749,50.2285385,&

             48.2692299,47.6495590,49.2938499,49.1924858,49.6449242,&

             50.0446815,51.9972496,54.2576981,52.9835434,50.4193535,&

             50.3617897,51.8276901,53.1239929,54.0682144,54.9238319,&

             55.6877632,54.8896332,54.0701065,52.2754097,52.2522354,&

             53.1248703,51.1287193,50.5003815,49.6504173,47.2453079,&

             45.4555626,45.8449707,45.9765129,45.7682228,45.2343674,&

             46.6496811,47.0894432,49.3368340,50.8058052,49.9132500,&

             49.5893288,48.2470627,46.9779968,45.6760864,45.7070389,&

             46.6158409,47.5303612,47.5630417,47.0389214,46.0352287,&

             45.8161545,45.7974396,46.0015373,45.3796463,45.3461685,&

             47.6444016,49.3327446,49.3810692,50.2027817,51.4567032,&

             52.3986320,52.5819206,52.7721825,52.6919098,53.3274345,&

             55.1345940,56.8962631,55.7791634,55.0616989,52.3551178,&

             51.3264084,51.0968323,51.1980476,52.8001442,52.0545082,&

             50.8742943,51.5150337,51.2242050,50.5033989,48.7760124,&

             47.4179192,49.7319527,51.3320541,52.3918304,52.4140434,&

             51.0845947,49.6485748,50.6893463,52.9840813,53.3246994,&

             52.4568024,51.9196091,53.6683121,53.4555359,51.7755814,&

             49.2915611,49.8755112,49.4546776,48.6171913,49.9643021,&

             49.3766441,49.2551308,50.1021881,51.0769119,55.8328133,&

             52.0212708,53.4930801,53.2147255,52.2356453,51.9648819,&

             52.1816330,51.9898071,52.5623627,51.0717278,52.2431946,&

             53.6943054,54.3752098,54.1492615,53.8523254,52.1093712,&

             52.3982697,51.2405128,50.3018112,51.3819618,49.5479546,&

             47.5024452,47.4447708,47.8939056,48.4070015,48.2440681,&

             48.7389755,49.7309227,49.1998024,49.5798340,51.1196213,&

             50.6288414,50.3971405,51.6084099,52.4564743,51.6443901,&

             52.4080658,52.4643364,52.6257210,53.1604691,51.9309731,&

             51.4137230,52.1233368,52.9867249,53.3180733,51.9647636,&

             50.7947655,52.3815842,50.8353729,49.4136009,52.8355217,&

             52.2234840,51.1392517,48.5245132,46.8700218,46.1607285,&

             45.2324257,47.4157829,48.9989090,49.6230736,50.4352913,&

             51.1652985,50.2588654,50.7820129,51.0448799,51.2880516,&

             49.6898804,49.0288200,49.9338837,48.2214432,46.2103348,&

             46.9550171,47.5595894,47.7176018,48.4502945,50.9816895,&

             51.6950073,51.6973495,52.1941261,51.8988075,52.5617599,&

             52.0218391,49.5236053,47.9684906,48.2445183,48.8275146,&

             49.7176971,51.5649338,52.5627213,52.0182419,50.9688835,&

             51.5846901,50.9486771,48.8685837,48.5600624,48.4760094,&

             48.5348396,50.4187813,51.2542381,50.1872864,50.4407692,&

             50.6222687,50.4972000,51.0036087,51.3367500,51.7368202,&

             53.0463791,53.6261253,52.0728683,48.9740753,49.3280830,&

             49.2733917,49.8519020,50.8562126,49.5594254,49.6109200,&

             48.3785629,48.0026474,49.4874268,50.1596375,51.8059540,&

             53.0288620,51.3321075,49.3114815,48.7999306,47.7201881,&

             46.3433914,46.5303612,47.6294632,48.6012459,47.8567657,&

             48.0604057,47.1352806,49.5724792,50.5566483,49.4182968,&

             50.5578079,50.6883736,50.6333389,51.9766159,51.0595245,&

             49.3751640,46.9667702,47.1658173,47.4411278,47.5360374,&

             48.9914742,50.4747620,50.2728043,51.9117165,53.7627792 /)

 

 

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

 

      call ts_outlier_identification(model, w, x, CRITICAL=3.5,&

                        RESSIGMA=ressigma, NOUTLIERS=noutliers,&

                        IOUTLIERSTATS=outlierstat,OMEGA=omega,&

                        PARAMS=parameters, AIC=aic, RELERR=1.0e-05)

 

      call umach(2, nout)

      write(nout,FMT="(/,T2,A)") 'ARMA parameters:'

      write(nout,FMT="(T3,I2,TR5,f10.6)") (i,parameters(i),i=1,3)

      write(nout,FMT="(/,T2,A,I2)") 'Number of outliers:', noutliers

      write(nout,FMT="(/,T2,A)") 'Outlier statistics:'

      write(nout,FMT="(T4,A,TR10,A)") 'Time point','Outlier type'

      write(nout,FMT="(I8,TR20,I2)") (outlierstat(i,:),i=1,noutliers)

      write(nout,FMT="(/,T2,A)") 'Omega statistics:'

      write(nout,FMT="(T4,A,TR10,A)") 'Time point','Omega'

      write(nout,FMT="(I9,TR10,f10.6)") (outlierstat(i,1),omega(i),&

                                           i=1,noutliers)

      write(nout,FMT="(/,T2,A,f10.6)") 'RSE: ', ressigma

      write(nout,FMT="(T2,A,f12.6)") 'AIC: ', aic

      end

 

Output

 

ARMA parameters:

   1      10.700689

   2       0.787765

   3      -0.498039

 

 Number of outliers: 2

 

 Outlier statistics:

   Time point          Outlier type

     150                     1

     200                     3

 

 Omega statistics:

   Time point          Omega

      150            4.478198

      200            3.381984

 

 RSE:   1.007209

 AIC:  1417.036377



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