Getting Started

Minimum AIC Model Selection

The time series model is automatically selected by using the AIC. The TSUNIMAR call estimates the univariate autoregressive model and computes the AIC. You need to specify the maximum lag or order of the AR process with the MAXLAG= option or put the maximum lag as the sixth argument of the TSUNIMAR call. Here is an example:

   proc iml;
   y = { 2.430 2.506 2.767 2.940 3.169 3.450 3.594 3.774 3.695 3.411
         2.718 1.991 2.265 2.446 2.612 3.359 3.429 3.533 3.261 2.612
         2.179 1.653 1.832 2.328 2.737 3.014 3.328 3.404 2.981 2.557
         2.576 2.352 2.556 2.864 3.214 3.435 3.458 3.326 2.835 2.476
         2.373 2.389 2.742 3.210 3.520 3.828 3.628 2.837 2.406 2.675
         2.554 2.894 3.202 3.224 3.352 3.154 2.878 2.476 2.303 2.360
         2.671 2.867 3.310 3.449 3.646 3.400 2.590 1.863 1.581 1.690
         1.771 2.274 2.576 3.111 3.605 3.543 2.769 2.021 2.185 2.588
         2.880 3.115 3.540 3.845 3.800 3.579 3.264 2.538 2.582 2.907
         3.142 3.433 3.580 3.490 3.475 3.579 2.829 1.909 1.903 2.033
         2.360 2.601 3.054 3.386 3.553 3.468 3.187 2.723 2.686 2.821
         3.000 3.201 3.424 3.531 };
   call tsunimar(arcoef,ev,nar,aic) data=y opt={-1 1} 
                                    print=1 maxlag=20;

You can also invoke the TSUNIMAR subroutine as follows:

   call tsunimar(arcoef,ev,nar,aic,y,20,{-1 1},,1);

The optional arguments can be omitted. In this example, the argument MISSING is omitted, and thus the default option (MISSING=0) is used. The summary table of the minimum AIC method is displayed in Figure 13.4 and Figure 13.5. The final estimates are given in Figure 13.6.

Figure 13.4: Minimum AIC Table - I

SAS Output
line
 
         ORDER  INNOVATION VARIANCE                                       
            M           V(M)              AIC(M)                          
            0          0.31607294   -108.26753229                         
            1          0.11481982   -201.45277331                         
            2          0.04847420   -280.51201122                         
            3          0.04828185   -278.88576251                         
            4          0.04656506   -280.28905616                         
            5          0.04615922   -279.11190502                         
            6          0.04511943   -279.25356641                         
            7          0.04312403   -281.50543541                         
            8          0.04201118   -281.96304075                         
            9          0.04128036   -281.61262868                         
           10          0.03829179   -286.67686828                         
           11          0.03318558   -298.13013264                         
           12          0.03255171   -297.94298716                         
           13          0.03247784   -296.15655602                         
           14          0.03237083   -294.46677874                         
           15          0.03234790   -292.53337704                         
           16          0.03187416   -291.92021487                         
           17          0.03183282   -290.04220196                         
           18          0.03126946   -289.72064823                         
           19          0.03087893   -288.90203735                         
           20          0.02998019   -289.67854830                         
 


Figure 13.5: Minimum AIC Table - II

SAS Output
line
                                     AIC(M)-AICMIN (truncated at 40.0)    
                                 0        10        20        30        40
            M   AIC(M)-AICMIN    +---------+---------+---------+---------+
            0    189.862600      |                                       .
            1     96.677359      |                                       .
            2     17.618121      |                 *                     |
            3     19.244370      |                  *                    |
            4     17.841076      |                 *                     |
            5     19.018228      |                  *                    |
            6     18.876566      |                  *                    |
            7     16.624697      |                *                      |
            8     16.167092      |               *                       |
            9     16.517504      |                *                      |
           10     11.453264      |          *                            |
           11             0      *                                       |
           12      0.187145      *                                       |
           13      1.973577      | *                                     |
           14      3.663354      |   *                                   |
           15      5.596756      |     *                                 |
           16      6.209918      |     *                                 |
           17      8.087931      |       *                               |
           18      8.409484      |       *                               |
           19      9.228095      |        *                              |
           20      8.451584      |       *                               |
                                 +---------+---------+---------+---------+
 
        *****  MINIMUM AIC = -298.130133  ATTAINED AT M = 11 *****        


The minimum AIC order is selected as 11. Then the coefficients are estimated as in Figure 13.6. Note that the first 20 observations are used as presample values.

Figure 13.6: Minimum AIC Estimation

SAS Output
line
        ..........................M  A  I  C  E.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.181322                                           .  
        .     2    -0.551571                                           .  
        .     3     0.231372                                           .  
        .     4    -0.178040                                           .  
        .     5     0.019874                                           .  
        .     6    -0.062573                                           .  
        .     7     0.028569                                           .  
        .     8    -0.050710                                           .  
        .     9     0.199896                                           .  
        .    10     0.161819                                           .  
        .    11    -0.339086                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =   -298.1301326                                    .  
        .      Innovation Variance =    0.033186                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =    21  FINISH =   114        .  
        ................................................................  
 
 


You can estimate the AR(11) model directly by specifying OPT=$\{ -1~  0\} $ and using the first 11 observations as presample values. The AR(11) estimates shown in Figure 13.7 are different from the minimum AIC estimates in Figure 13.6 because the samples are slightly different. Here is the code:

   call tsunimar(arcoef,ev,nar,aic,y,11,{-1 0},,1);

Figure 13.7: AR(11) Estimation

SAS Output
line
        ..........................M  A  I  C  E.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.149416                                           .  
        .     2    -0.533719                                           .  
        .     3     0.276312                                           .  
        .     4    -0.326420                                           .  
        .     5     0.169336                                           .  
        .     6    -0.164108                                           .  
        .     7     0.073123                                           .  
        .     8    -0.030428                                           .  
        .     9     0.151227                                           .  
        .    10     0.192808                                           .  
        .    11    -0.340200                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =   -318.7984105                                    .  
        .      Innovation Variance =    0.036563                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =    12  FINISH =   114        .  
        ................................................................  
 
 


The minimum AIC procedure can also be applied to the vector autoregressive (VAR) model by using the TSMULMAR subroutine. See the section Multivariate Time Series Analysis for details. Three variables are used as input. The maximum lag is specified as 10. Here is the code:

   data one;
      input invest income consum @@;
   datalines;
      . . . data lines omitted . . .
      ;
   proc iml;
   use one;
   read all into y var{invest income consum};
   mdel  = 1;
   maice = 2;
   misw  = 0;
   opt   = mdel ||  maice || misw;
   maxlag = 10;
   miss   = 0;
   print  = 1;
   call tsmulmar(arcoef,ev,nar,aic,y,maxlag,opt,miss,print);

The VAR(3) model minimizes the AIC and was selected as an appropriate model (see Figure 13.8). However, AICs of the VAR(4) and VAR(5) models show little difference from VAR(3). You can also choose VAR(4) or VAR(5) as an appropriate model in the context of minimum AIC since this AIC difference is much less than 1.

Figure 13.8: VAR Model Selection

SAS Output
line
 
         ORDER  INNOVATION VARIANCE                                       
            M         LOG(|V(M)|)         AIC(M)                          
            0         25.98001095   2136.36089828                         
            1         15.70406486   1311.73331883                         
            2         15.48896746   1312.09533158                         
            3         15.18567834   1305.22562428                         
            4         14.96865183   1305.42944974                         
            5         14.74838535   1305.36759889                         
            6         14.60269347   1311.42086432                         
            7         14.54981887   1325.08514729                         
            8         14.38596333   1329.64899297                         
            9         14.16383772   1329.43469312                         
           10         13.85377849   1322.00983656                         
 
                                     AIC(M)-AICMIN (truncated at 40.0)    
                                 0        10        20        30        40
            M   AIC(M)-AICMIN    +---------+---------+---------+---------+
            0    831.135274      |                                       .
            1      6.507695      |      *                                |
            2      6.869707      |      *                                |
            3             0      *                                       |
            4      0.203825      *                                       |
            5      0.141975      *                                       |
            6      6.195240      |     *                                 |
            7     19.859523      |                   *                   |
            8     24.423369      |                       *               |
            9     24.209069      |                       *               |
           10     16.784212      |                *                      |
                                 +---------+---------+---------+---------+
 


The TSMULMAR subroutine estimates the instantaneous response model with diagonal error variance. See the section Multivariate Time Series Analysis for details on the instantaneous response model. Therefore, it is possible to select the minimum AIC model independently for each equation. The best model is selected by specifying MAXLAG=5, as in the following code:

   call tsmulmar(arcoef,ev,nar,aic) data=y maxlag=5
        opt={1 1 0} print=1;

Figure 13.9: Model Selection via Instantaneous Response Model

variance
256.64375 29.803549 76.846777
29.803549 228.97341 119.60387
76.846777 119.60387 134.21764

arcoef
13.312109 1.5459098 15.963897
0.8257397 0.2514803 0
0.0958916 1.0057088 0
0.0320985 0.3544346 0.4698934
0.044719 -0.201035 0
0.0051931 -0.023346 0
0.1169858 -0.060196 0.0483318
0.1867829 0 0
0.0216907 0 0
-0.117786 0 0.3500366
0.1541108 0 0
0.0178966 0 0
0.0461454 0 -0.191437
-0.389644 0 0
-0.045249 0 0
-0.116671 0 0

aic
1347.6198


You can print the intermediate results of the minimum AIC procedure by using the PRINT=2 option.

Note that the AIC value depends on the MAXLAG=lag option and the number of parameters estimated. The minimum AIC VAR estimation procedure (MAICE=2) uses the following AIC formula:

\[  (T-lag) \log (|\hat{\Sigma }|) + 2(p \times n^2 + n \times \mi {intercept})  \]

In this formula, $p$ is the order of the $n$-variate VAR process, and intercept=1 if the intercept is specified; otherwise, intercept=0. When you use the MAICE=1 or MAICE=0 option, AIC is computed as the sum of AIC for each response equation. Therefore, there is an AIC difference of $n(n-1)$ since the instantaneous response model contains the additional $n(n-1)/2$ response variables as regressors.

The following code estimates the instantaneous response model. The results are shown in Figure 13.10.

   call tsmulmar(arcoef,ev,nar,aic) data=y 
                           maxlag=3 opt={1 0 0};
   print aic nar;
   print arcoef;

Figure 13.10: AIC from Instantaneous Response Model

aic nar
1403.0762 3

arcoef
4.8245814 5.3559216 17.066894
0.8855926 0.3401741 -0.014398
0.1684523 1.0502619 0.107064
0.0891034 0.4591573 0.4473672
-0.059195 -0.298777 0.1629818
0.1128625 -0.044039 -0.088186
0.1684932 -0.025847 -0.025671
0.0637227 -0.196504 0.0695746
-0.226559 0.0532467 -0.099808
-0.303697 -0.139022 0.2576405


The following code estimates the VAR model. The results are shown in Figure 13.11.

   call tsmulmar(arcoef,ev,nar,aic) data=y maxlag=3
        opt={1 2 0};
   print aic nar;
   print arcoef;

Figure 13.11: AIC from VAR Model

aic nar
1397.0762 3

arcoef
4.8245814 5.3559216 17.066894
0.8855926 0.3401741 -0.014398
0.1684523 1.0502619 0.107064
0.0891034 0.4591573 0.4473672
-0.059195 -0.298777 0.1629818
0.1128625 -0.044039 -0.088186
0.1684932 -0.025847 -0.025671
0.0637227 -0.196504 0.0695746
-0.226559 0.0532467 -0.099808
-0.303697 -0.139022 0.2576405


The AIC computed from the instantaneous response model is greater than that obtained from the VAR model estimation by 6. There is a discrepancy between Figure 13.11 and Figure 13.8 because different observations are used for estimation.

Nonstationary Data Analysis

The following example shows how to manage nonstationary data by using TIMSAC calls. In practice, time series are considered to be stationary when the expected values of first and second moments of the series do not change over time. This weak or covariance stationarity can be modeled by using the TSMLOCAR, TSMLOMAR, TSDECOMP, and TSTVCAR subroutines.

First, the locally stationary model is estimated. The whole series (1000 observations) is divided into three blocks of size 300 and one block of size 90, and the minimum AIC procedure is applied to each block of the data set. See the section Nonstationary Time Series for more details. Here is the code:

   data one;
      input y @@;
   datalines;
      . . . data lines omitted . . .
      ;

   proc iml;
   use one;
   read all var{y};

   mdel = -1;
   lspan = 300; /* local span of data */
   maice = 1;
   opt = mdel || lspan || maice;
   call tsmlocar(arcoef,ev,nar,aic,first,last)
                   data=y maxlag=10 opt=opt print=2;

Estimation results are displayed with the graphs of power spectrum $(\log 10(f_{YY}(g)))$, where $f_{YY}(g)$ is a rational spectral density function. See the section Spectral Analysis. The estimates for the first block and third block are shown in Figure 13.12 and Figure 13.15, respectively. As the first block and the second block do not have any sizable difference, the pooled model (AIC=45.892) is selected instead of the moving model (AIC=46.957) in Figure 13.13. However, you can notice a slight change in the shape of the spectrum of the third block of the data (observations 611 through 910). See Figure 13.14 and Figure 13.16 for comparison. The moving model is selected since the AIC (106.830) of the moving model is smaller than that of the pooled model (108.867).

Figure 13.12: Locally Stationary Model for First Block

SAS Output
line
        INITIAL LOCAL MODEL: N_CURR =   300                               
                             NAR_CURR =   8  AIC =   37.583203            
 
 
 
        ..........................CURRENT MODEL.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.605717                                           .  
        .     2    -1.245350                                           .  
        .     3     1.014847                                           .  
        .     4    -0.931554                                           .  
        .     5     0.394230                                           .  
        .     6    -0.004344                                           .  
        .     7     0.111608                                           .  
        .     8    -0.124992                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =     37.5832030                                    .  
        .      Innovation Variance =    1.067455                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =    11  FINISH =   310        .  
        ................................................................  


Figure 13.13: Locally Stationary Model Comparison

SAS Output
line
        --- THE FOLLOWING TWO MODELS ARE COMPARED ---                     
 
        MOVING MODEL: (N_PREV =  300, N_CURR = 300)                       
                       NAR_CURR   =  7  AIC =   46.957398                 
        CONSTANT MODEL: N_POOLED =   600                                  
                       NAR_POOLED =  8  AIC =   45.892350                 
 
        *****  CONSTANT MODEL ADOPTED  *****                              
 
 
 
        ..........................CURRENT MODEL.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.593890                                           .  
        .     2    -1.262379                                           .  
        .     3     1.013733                                           .  
        .     4    -0.926052                                           .  
        .     5     0.314480                                           .  
        .     6     0.193973                                           .  
        .     7    -0.058043                                           .  
        .     8    -0.078508                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =     45.8923501                                    .  
        .      Innovation Variance =    1.047585                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =    11  FINISH =   610        .  
        ................................................................  


Figure 13.14: Power Spectrum for First and Second Blocks

SAS Output
line
        POWER SPECTRAL DENSITY                                            
   20.00+                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |  XXXX                                                           
        XXX    XX    XXX                                                  
        |        XXXX                                                     
        |               X                                                 
        |                                                                 
   10.00+                                                                 
        |                X                                                
        |                                                                 
        |                 X                                               
        |                                                                 
        |                  X               XX                             
        |                                    X                            
        |                   X             X                               
        |                                                                 
        |                    X           X    X                           
       0+                     X                                           
        |                      X        X      X                          
        |                       XX    XX                                  
        |                         XXXX          X                         
        |                                                                 
        |                                        X                        
        |                                         X                       
        |                                                                 
        |                                          X                      
        |                                           X                     
   -10.0+                                            X                    
        |                                             XX                  
        |                                               XX                
        |                                                 XX              
        |                                                   XXX           
        |                                                      XXXXXX     
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
   -20.0+-----------+-----------+-----------+-----------+-----------+     
       0.0         0.1         0.2         0.3         0.4         0.5    
 
                                        FREQUENCY                         


Figure 13.15: Locally Stationary Model for Third Block

SAS Output
line
        --- THE FOLLOWING TWO MODELS ARE COMPARED ---                     
 
        MOVING MODEL: (N_PREV =  600, N_CURR = 300)                       
                       NAR_CURR   =  7  AIC =  106.829869                 
        CONSTANT MODEL: N_POOLED =   900                                  
                       NAR_POOLED =  8  AIC =  108.867091                 
 
        *************************************                             
        *****                           *****                             
        *****     NEW MODEL ADOPTED     *****                             
        *****                           *****                             
        *************************************                             
 
 
 
        ..........................CURRENT MODEL.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.648544                                           .  
        .     2    -1.201812                                           .  
        .     3     0.674933                                           .  
        .     4    -0.567576                                           .  
        .     5    -0.018924                                           .  
        .     6     0.516627                                           .  
        .     7    -0.283410                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =     60.9375188                                    .  
        .      Innovation Variance =    1.161592                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =   611  FINISH =   910        .  
        ................................................................  


Figure 13.16: Power Spectrum for Third Block

SAS Output
line
        POWER SPECTRAL DENSITY                                            
   20.00+             X                                                   
        |            X                                                    
        |           X                                                     
        |              X                                                  
        |        XXX                                                      
        |   XXXXX                                                         
        | XX                                                              
        XX              X                                                 
        |                                                                 
        |                                                                 
   10.00+                X                                                
        |                                                                 
        |                                                                 
        |                 X                                               
        |                                                                 
        |                  X                                              
        |                                   X                             
        |                   X              X                              
        |                                    X                            
        |                    X            X                               
       0+                     X          X    X                           
        |                      X                                          
        |                       X      XX      X                          
        |                        XXXXXX                                   
        |                                       X                         
        |                                                                 
        |                                        X                        
        |                                                                 
        |                                         X                       
        |                                          X                      
   -10.0+                                           X                     
        |                                            XX                   
        |                                              XX       XXXXX     
        |                                                XXXXXXX          
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
   -20.0+-----------+-----------+-----------+-----------+-----------+     
       0.0         0.1         0.2         0.3         0.4         0.5    
 
                                        FREQUENCY                         


Finally, the moving model is selected since there is a structural change in the last block of data (observations 911 through 1000). The final estimates are stored in variables ARCOEF, EV, NAR, AIC, FIRST, and LAST. The final estimates and spectrum are given in Figure 13.17 and Figure 13.18, respectively. The power spectrum of the final model (Figure 13.18) is significantly different from that of the first and second blocks (see Figure 13.14).

Figure 13.17: Locally Stationary Model for Last Block

SAS Output
line
        --- THE FOLLOWING TWO MODELS ARE COMPARED ---                     
 
        MOVING MODEL: (N_PREV =  300, N_CURR =  90)                       
                       NAR_CURR   =  6  AIC =  139.579012                 
        CONSTANT MODEL: N_POOLED =   390                                  
                       NAR_POOLED =  9  AIC =  167.783711                 
 
        *************************************                             
        *****                           *****                             
        *****     NEW MODEL ADOPTED     *****                             
        *****                           *****                             
        *************************************                             
 
 
 
        ..........................CURRENT MODEL.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients: AR(M)                             .  
        .                                                              .  
        .     1     1.181022                                           .  
        .     2    -0.321178                                           .  
        .     3    -0.113001                                           .  
        .     4    -0.137846                                           .  
        .     5    -0.141799                                           .  
        .     6     0.260728                                           .  
        .                                                              .  
        .                                                              .  
        .      AIC =     78.6414932                                    .  
        .      Innovation Variance =    2.050818                       .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =   911  FINISH =  1000        .  
        ................................................................  


Figure 13.18: Power Spectrum for Last Block

SAS Output
line
        POWER SPECTRAL DENSITY                                            
   30.00+                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
        |           X                                                     
        |                                                                 
        |            X                                                    
        |                                                                 
        |                                                                 
   20.00+          X                                                      
        |                                                                 
        |                                                                 
        |         X   X                                                   
        |                                                                 
        |        X                                                        
        XXX     X                                                         
        |  XXXXX       X                                                  
        |                                                                 
        |                                                                 
   10.00+               X                                                 
        |                                                                 
        |                X                                                
        |                                                                 
        |                 X                                               
        |                                                                 
        |                  X                                              
        |                   X                                             
        |                    X                                            
        |                     XX                                          
       0+                       XX      XXXXXX                            
        |                         XXXXXX      XX                          
        |                                       XX                        
        |                                         XX               XX     
        |                                           XX         XXXX       
        |                                             XXXXXXXXX           
        |                                                                 
        |                                                                 
        |                                                                 
        |                                                                 
   -10.0+-----------+-----------+-----------+-----------+-----------+     
       0.0         0.1         0.2         0.3         0.4         0.5    
 
                                        FREQUENCY                         
 


The multivariate analysis for locally stationary data is a straightforward extension of the univariate analysis. The bivariate locally stationary VAR models are estimated. The selected model is the VAR(7) process with some zero coefficients over the last block of data. There seems to be a structural difference between observations from 11 to 610 and those from 611 to 896. Here is the code:

   proc iml;
   rudder = {. . . data lines omitted . . .};
   yawing = {. . . data lines omitted . . .};

   y = rudder` || yawing`;
   c = {0.01795 0.02419};
   *-- calibration of data --*/
   y = y # (c @ j(n,1,1));
   mdel = -1;
   lspan = 300; /* local span of data */
   maice = 1;
   call tsmlomar(arcoef,ev,nar,aic,first,last) data=y maxlag=10
           opt=(mdel || lspan || maice) print=1;

The results of the analysis are shown in Figure 13.19.

Figure 13.19: Locally Stationary VAR Model Analysis

SAS Output
line
        --- THE FOLLOWING TWO MODELS ARE COMPARED ---                     
 
        MOVING MODEL: (N_PREV =  600, N_CURR = 286)                       
                       NAR_CURR   =  7  AIC = -823.845234                 
        CONSTANT MODEL: N_POOLED =   886                                  
                       NAR_POOLED = 10  AIC = -716.818588                 
 
        *************************************                             
        *****                           *****                             
        *****     NEW MODEL ADOPTED     *****                             
        *****                           *****                             
        *************************************                             
 
 
 
        ..........................CURRENT MODEL.........................  
        .                                                              .  
        .                                                              .  
        .                                                              .  
        .     M     AR Coefficients                                    .  
        .                                                              .  
        .     1     0.932904    -0.130964                              .  
        .          -0.024401     0.599483                              .  
        .     2     0.163141     0.266876                              .  
        .          -0.135605     0.377923                              .  
        .     3    -0.322283     0.178194                              .  
        .           0.188603    -0.081245                              .  
        .     4     0.166094    -0.304755                              .  
        .          -0.084626    -0.180638                              .  
        .     5            0            0                              .  
        .                  0    -0.036958                              .  
        .     6            0            0                              .  
        .                  0     0.034578                              .  
        .     7            0            0                              .  
        .                  0     0.268414                              .  
        .                                                              .  
        .                                                              .  
        .      AIC =   -114.6911872                                    .  
        .                                                              .  
        .           Innovation Variance                                .  
        .                                                              .  
        .           1.069929     0.145558                              .  
        .           0.145558     0.563985                              .  
        .                                                              .  
        .                                                              .  
        .            INPUT DATA   START =   611  FINISH =   896        .  
        ................................................................  
 


Consider the time series decomposition

\[  y_ t = T_ t + S_ t + u_ t + \epsilon _ t  \]

where $T_ t$ and $S_ t$ are trend and seasonal components, respectively, and $u_ t$ is a stationary AR($p$) process. The annual real GNP series is analyzed under second difference stochastic constraints on the trend component and the stationary AR(2) process.

\begin{eqnarray*}  T_ t &  = &  2T_{t-1} - T_{t-2} + w_{1t} \\[0.05in] u_ t &  = &  \alpha _1 u_{t-1} + \alpha _2 u_{t-2} + w_{2t} \end{eqnarray*}

The seasonal component is ignored if you specify SORDER=0. Therefore, the following state space model is estimated:

\begin{eqnarray*}  y_ t &  = &  \textbf{H}\textbf{z}_ t + \epsilon _ t \\[0.05in] \textbf{z}_ t &  = &  \textbf{F}\textbf{z}_{t-1} + \textbf{w}_ t \end{eqnarray*}

where

\begin{eqnarray*}  \textbf{H} &  = &  \left[ \begin{array}{cccc} 1 &  0 &  1 &  0 \end{array} \right] \\[0.05in] \textbf{F} &  = &  \left[ \begin{array}{cccc} 2 &  -1 &  0 &  0 \\ 1 &  0 &  0 &  0 \\ 0 &  0 &  \alpha _1 &  \alpha _2 \\ 0 &  0 &  1 &  0 \end{array} \right] \\[0.05in] \textbf{z}_ t &  = &  (T_ t, T_{t-1}, u_ t, u_{t-1})^{\prime } \\[0.05in] \textbf{w}_ t &  = &  (w_{1t}, 0, w_{2t}, 0)^{\prime } \sim \left(0, \left[ \begin{array}{cccc} \sigma _1^2 &  0 &  0 &  0 \\ 0 &  0 &  0 &  0 \\ 0 &  0 &  \sigma _2^2 &  0 \\ 0 &  0 &  0 &  0 \end{array} \right] \right) \end{eqnarray*}

The parameters of this state space model are $\sigma ^2_1$, $\sigma ^2_2$, $\alpha _1$, and $\alpha _2$. Here is the code:

   proc iml;
   y = { 116.8 120.1 123.2 130.2 131.4 125.6 124.5 134.3
         135.2 151.8 146.4 139.0 127.8 147.0 165.9 165.5
         179.4 190.0 189.8 190.9 203.6 183.5 169.3 144.2
         141.5 154.3 169.5 193.0 203.2 192.9 209.4 227.2
         263.7 297.8 337.1 361.3 355.2 312.6 309.9 323.7
         324.1 355.3 383.4 395.1 412.8 406.0 438.0 446.1
         452.5 447.3 475.9 487.7 497.2 529.8 551.0 581.1
         617.8 658.1 675.2 706.6 724.7 };
   y = y`; /*-- convert to column vector --*/
   mdel  = 0;
   trade = 0;
   tvreg = 0;
   year  = 0;
   period= 0;
   log   = 0;
   maxit = 100;
   update = .; /* use default update method      */
   line   = .; /* use default line search method */
   sigmax = 0; /* no upper bound for variances   */
   back = 100;
   opt = mdel || trade || year || period || log || maxit ||
         update || line || sigmax || back;
   call tsdecomp(cmp,coef,aic) data=y order=2 sorder=0 nar=2
        npred=5 opt=opt icmp={1 3} print=1;
   y = y[52:61];
   cmp = cmp[52:66,];
   print y cmp;

The estimated parameters are printed when you specify the PRINT= option. In Figure 13.20, the estimated variances are printed under the title of TAU2(I), showing that $\hat{\sigma }_1^2=2.915$ and $\hat{\sigma }_2^2=113.9577$. The AR coefficient estimates are $\hat{\alpha }_1=1.397$ and $\hat{\alpha }_2=-0.595$. These estimates are also stored in the output matrix COEF.

Figure 13.20: Nonstationary Time Series and State Space Modeling

SAS Output
line
 
                       <<                          
 
                  ---  PARAMETER VECTOR  ---                              
 
                   1.607423E-01 6.283820E+00 8.761627E-01 -5.94879E-01    
 
                  ---  GRADIENT  ---                                      
 
                   3.352158E-04 5.237221E-06 2.907539E-04 -1.24376E-04    
 
 
                  LIKELIHOOD = -249.937193     SIG2 =       18.135085     
                  AIC        =  509.874385                                
 
          I   TAU2(I)       AR(I)      PARCOR(I)                          
          1     2.915075     1.397374     0.876163                        
          2   113.957607    -0.594879    -0.594879                        
 


The trend and stationary AR components are estimated by using the smoothing method, and out-of-sample forecasts are computed by using a Kalman filter prediction algorithm. The trend and AR components are stored in the matrix CMP since the ICMP={1 3} option is specified. The last 10 observations of the original series Y and the last 15 observations of two components are shown in Figure 13.21. Note that the first column of CMP is the trend component and the second column is the AR component. The last 5 observations of the CMP matrix are out-of-sample forecasts.

Figure 13.21: Smoothed and Predicted Values of Two Components

y cmp  
487.7 514.01141 -26.94342
497.2 532.62744 -32.48672
529.8 552.02402 -24.46593
551 571.90121 -20.15112
581.1 592.31944 -10.58646
617.8 613.21855 5.2504401
658.1 634.43665 20.799207
675.2 655.70431 22.161604
706.6 677.2125 27.927978
724.7 698.72364 25.957962
  720.23478 19.6592
  741.74593 12.029396
  763.25707 5.1147111
  784.76821 -0.008876
  806.27935 -3.05504


Seasonal Adjustment

Consider the simple time series decomposition

\[  y_ t = T_ t + S_ t + \epsilon _ t  \]

The TSBAYSEA subroutine computes seasonally adjusted series by estimating the seasonal component. The seasonally adjusted series is computed as $y_ t^* = y_ t - \hat{S}_ t$. The details of the adjustment procedure are given in the section Bayesian Seasonal Adjustment.

The monthly labor force series (1972–1978) are analyzed. You do not need to specify the options vector if you want to use the default options. However, you should change OPT[2] when the data frequency is not monthly (OPT[2]=12). The NPRED= option produces the multistep forecasts for the trend and seasonal components. The stochastic constraints are specified as ORDER=2 and SORDER=1.

\begin{eqnarray*}  T_ t &  = &  2T_{t-1} - T_{t-2} + w_{1t} \\ S_ t &  = &  -S_{t-1} - \cdots - S_{t-11} + w_{2t} \end{eqnarray*}

In Figure 13.22, the first column shows the trend components; the second column shows the seasonal components; the third column shows the forecasts; the fourth column shows the seasonally adjusted series; the last column shows the value of ABIC. The last 12 rows are the forecasts. The figure is generated by using the following statements:

   proc iml;
   y = { 5447 5412 5215 4697 4344 5426   
         5173 4857 4658 4470 4268 4116    
         4675 4845 4512 4174 3799 4847    
         4550 4208 4165 3763 4056 4058    
         5008 5140 4755 4301 4144 5380    
         5260 4885 5202 5044 5685 6106    
         8180 8309 8359 7820 7623 8569    
         8209 7696 7522 7244 7231 7195    
         8174 8033 7525 6890 6304 7655    
         7577 7322 7026 6833 7095 7022    
         7848 8109 7556 6568 6151 7453    
         6941 6757 6437 6221 6346 5880 }; 
   y = y`;

   call tsbaysea(trend,season,series,adj,abic)
        data=y order=2 sorder=1 npred=12 print=2;
   print trend season series adj abic;

Figure 13.22: Trend and Seasonal Component Estimates and Forecasts

obs trend season series adj abic
1 4843.2502 576.86675 5420.1169 4870.1332 874.04585
2 4848.6664 612.79607 5461.4624 4799.2039  
3 4871.2876 324.02004 5195.3077 4890.98  
4 4896.6633 -198.7601 4697.9032 4895.7601  
5 4922.9458 -572.5562 4350.3896 4916.5562  
. . . . .  
71 6551.6017 -266.2162 6285.3855 6612.2162  
72 6388.9012 -440.3472 5948.5539 6320.3472  
73 6226.2006 650.7707 6876.9713    
74 6063.5001 800.93733 6864.4374    
75 5900.7995 396.19866 6296.9982    
76 5738.099 -340.2852 5397.8137    
77 5575.3984 -719.1146 4856.2838    
78 5412.6979 553.19764 5965.8955    
79 5249.9973 202.06582 5452.0631    
80 5087.2968 -54.44768 5032.8491    
81 4924.5962 -295.2747 4629.3215    
82 4761.8957 -487.6621 4274.2336    
83 4599.1951 -266.1917 4333.0034    
84 4436.4946 -440.3354 3996.1591    


The estimated spectral density function of the irregular series $\hat{\epsilon }_ t$ is shown in Figure 13.23 and Figure 13.24.

Figure 13.23: Spectrum of Irregular Component

SAS Output
line
    I   Rational   0.0       10.0      20.0      30.0      40.0      50.0      60.0                 
          Spectrum  +---------+---------+---------+---------+---------+---------+                   
    0 1.366798E+00  |*                                                           ===>X              
    1 1.571261E+00  |*                                                                              
    2 2.414836E+00  |  *                                                                            
    3 5.151906E+00  |      *                                                                        
    4 1.634887E+01  |           *                                                                   
    5 8.085674E+01  |                  *                                                            
    6 3.805530E+02  |                        *                                                      
    7 8.082536E+02  |                            *                                                  
    8 6.366350E+02  |                           *                                                   
    9 3.479435E+02  |                        *                                                      
   10 3.872650E+02  |                        *                                   ===>X              
   11 1.264805E+03  |                              *                                                
   12 1.726138E+04  |                                         *                                     
   13 1.559041E+03  |                              *                                                
   14 1.276516E+03  |                              *                                                
   15 3.861089E+03  |                                  *                                            
   16 9.593184E+03  |                                      *                                        
   17 3.662145E+03  |                                  *                                            
   18 5.499783E+03  |                                    *                                          
   19 4.443303E+03  |                                   *                                           
   20 1.238135E+03  |                             *                              ===>X              
   21 8.392131E+02  |                            *                                                  
   22 1.258933E+03  |                              *                                                
   23 2.932003E+03  |                                 *                                             
   24 1.857923E+03  |                               *                                               
   25 1.171437E+03  |                             *                                                 
   26 1.611958E+03  |                               *                                               
   27 4.822498E+03  |                                   *                                           
   28 4.464961E+03  |                                   *                                           
   29 1.951547E+03  |                               *                                               
   30 1.653182E+03  |                               *                            ===>X              
   31 2.308152E+03  |                                *                                              
   32 5.475758E+03  |                                    *                                          
   33 2.349584E+04  |                                          *                                    
   34 5.266969E+03  |                                    *                                          
   35 2.058667E+03  |                                *                                              
   36 2.215595E+03  |                                *                                              
   37 8.181540E+03  |                                      *                                        
   38 3.077329E+03  |                                 *                                             
   39 7.577961E+02  |                           *                                                   
   40 5.057636E+02  |                          *                                 ===>X              
   41 7.312090E+02  |                           *                                                   
   42 3.131377E+03  |                                 *                          ===>T              
   43 8.173276E+03  |                                      *                                        
   44 1.958359E+03  |                               *                                               
   45 2.216458E+03  |                                *                                              
   46 4.215465E+03  |                                   *                                           
   47 9.659340E+02  |                            *                                                  
   48 3.758466E+02  |                        *                                                      
   49 2.849326E+02  |                       *                                                       
   50 3.617848E+02  |                        *                                   ===>X              
   51 7.659839E+02  |                           *                                                   
   52 3.191969E+03  |                                  *                                            


Figure 13.24: continued

SAS Output
line
   53 1.768107E+04  |                                         *                                     
   54 5.281385E+03  |                                    *                                          
   55 2.959704E+03  |                                 *                                             
   56 3.783522E+03  |                                  *                                            
   57 1.896625E+04  |                                         *                                     
   58 1.041753E+04  |                                       *                                       
   59 2.038940E+03  |                                *                                              
   60 1.347568E+03  |                              *                             ===>X              
 
X: If peaks (troughs) appear                                                                        
   at these frequencies,                                                                            
   try lower (higher) values                                                                        
   of rigid and watch ABIC                                                                          
T: If a peaks appears here                                                                          
   try trading-day adjustment                                                                       


Miscellaneous Time Series Analysis Tools

The forecast values of multivariate time series are computed by using the TSPRED call. In the following example, the multistep-ahead forecasts are produced from the VARMA(2,1) estimates. Since the VARMA model is estimated by using the mean deleted series, you should specify the CONSTANT=$-1$ option. You need to provide the original series instead of the mean deleted series to get the correct predictions. The forecast variance MSE and the impulse response function IMPULSE are also produced.

The VARMA($p,q$) model is written

\[  \textbf{y}_ t + \sum _{i=1}^ p \textbf{A}_ i\textbf{y}_{t-i} = \epsilon _ t + \sum _{i=1}^ q \textbf{M}_ i\epsilon _{t-i}  \]

Then the COEF matrix is constructed by stacking matrices $\bA _1,\ldots ,\bA _ p,$ $\bM _1,\ldots ,\bM _ q$. Here is the code:

   proc iml;
   c = { 264 235 239 239 275 277 274 334 334 306
         308 309 295 271 277 221 223 227 215 223
         241 250 270 303 311 307 322 335 335 334
         309 262 228 191 188 215 215 249 291 296 };
   f = { 690 690 688 690 694 702 702 702 700 702
         702 694 708 702 702 708 700 700 702 694
         698 694 700 702 700 702 708 708 710 704
         704 700 700 694 702 694 710 710 710 708 };
   t = { 1152 1288 1288 1288 1368 1456 1656 1496 1744 1464
         1560 1376 1336 1336 1296 1296 1280 1264 1280 1272
         1344 1328 1352 1480 1472 1600 1512 1456 1368 1280
         1224 1112 1112 1048 1176 1064 1168 1280 1336 1248 };
   p = { 254.14 253.12 251.85 250.41 249.09 249.19 249.52 250.19
         248.74 248.41 249.95 250.64 250.87 250.94 250.96 251.33
         251.18 251.05 251.00 250.99 250.79 250.44 250.12 250.19
         249.77 250.27 250.74 250.90 252.21 253.68 254.47 254.80
         254.92 254.96 254.96 254.96 254.96 254.54 253.21 252.08 };

   y = c` || f` || t` || p`;
   ar = {  .82028   -.97167    .079386   -5.4382,
          -.39983    .94448    .027938   -1.7477,
          -.42278  -2.3314    1.4682    -70.996,
           .031038  -.019231  -.0004904   1.3677,
          -.029811   .89262   -.047579    4.7873,
           .31476    .0061959 -.012221    1.4921,
           .3813    2.7182    -.52993    67.711,
          -.020818   .01764    .00037981  -.38154 };
   ma = {  .083035 -1.0509     .055898   -3.9778,
          -.40452    .36876    .026369    -.81146,
           .062379 -2.6506     .80784   -76.952,
           .03273   -.031555  -.00019776  -.025205 };
   coef = ar // ma;
   ev = { 188.55   6.8082    42.385   .042942,
          6.8082   32.169    37.995  -.062341,
          42.385   37.995    5138.8  -.10757,
          .042942  -.062341  -.10757  .34313 };

   nar = 2; nma = 1;
   call tspred(forecast,impulse,mse,y,coef,nar,nma,ev,
                  5,nrow(y),-1);

Figure 13.25: Multivariate ARMA Prediction

observed
Y1
Y2 predicted
P1
P2
264 690 269.950 700.750
235 690 256.764 691.925
239 688 239.996 693.467
239 690 242.320 690.951
275 694 247.169 693.214
277 702 279.024 696.157
274 702 284.041 700.449
334 702 286.890 701.580
334 700 321.798 699.851
306 702 330.355 702.383
308 702 297.239 700.421
309 694 302.651 701.928
295 708 294.570 696.261
271 702 283.254 703.936
277 702 269.600 703.110
221 708 270.349 701.557
223 700 231.523 705.438
227 700 233.856 701.785
215 702 234.883 700.185
223 694 229.156 701.837
241 698 235.054 697.060
250 694 249.288 698.181
270 700 257.644 696.665
303 702 272.549 699.281
311 700 301.947 701.667
307 702 306.422 700.708
322 708 304.120 701.204
335 708 311.590 704.654
335 710 320.570 706.389
334 704 315.127 706.439
309 704 311.083 703.735
262 700 288.159 702.801
228 700 251.352 700.805
191 694 226.749 700.247
188 702 199.775 696.570
215 694 202.305 700.242
215 710 222.951 696.451
249 710 226.553 704.483
291 710 259.927 707.610
296 708 291.446 707.861
    293.899 707.430
    293.477 706.933
    292.564 706.190
    290.313 705.384
    286.559 704.618


The first 40 forecasts in Figure 13.25 are one-step predictions. The last observation is the five-step forecast values of variables C and F. You can construct the confidence interval for these forecasts by using the mean square error matrix, MSE. See the section Multivariate Time Series Analysis for more details about impulse response functions and the mean square error matrix.

The TSROOT call computes the polynomial roots of the AR and MA equations. When the AR($p$) process is written

\[  y_ t = \sum _{i=1}^ p \alpha _ i y_{t-i} + \epsilon _ t  \]

you can specify the following polynomial equation:

\[  z^ p - \sum _{i=1}^ p \alpha _ i z^{p-i} = 0  \]

When all $p$ roots of the preceding equation are inside the unit circle, the AR($p$) process is stationary. The MA($q$) process is invertible if the following polynomial equation has all roots inside the unit circle:

\[  z^ q + \sum _{i=1}^ q \theta _ i z^{q-i} = 0  \]

where $\theta _ i$ are the MA coefficients. For example, the best AR model is selected and estimated by the TSUNIMAR subroutine (see Figure 13.26). You can obtain the roots of the preceding equation by calling the TSROOT subroutine. Since the TSROOT subroutine can handle the complex AR or MA coefficients, note that you should add zero imaginary coefficients for the second column of the MATIN matrix for real coefficients. Here is the code:

   proc iml;
   y = { 2.430 2.506 2.767 2.940 3.169 3.450 3.594 3.774 3.695 3.411
         2.718 1.991 2.265 2.446 2.612 3.359 3.429 3.533 3.261 2.612
         2.179 1.653 1.832 2.328 2.737 3.014 3.328 3.404 2.981 2.557
         2.576 2.352 2.556 2.864 3.214 3.435 3.458 3.326 2.835 2.476
         2.373 2.389 2.742 3.210 3.520 3.828 3.628 2.837 2.406 2.675
         2.554 2.894 3.202 3.224 3.352 3.154 2.878 2.476 2.303 2.360
         2.671 2.867 3.310 3.449 3.646 3.400 2.590 1.863 1.581 1.690
         1.771 2.274 2.576 3.111 3.605 3.543 2.769 2.021 2.185 2.588
         2.880 3.115 3.540 3.845 3.800 3.579 3.264 2.538 2.582 2.907
         3.142 3.433 3.580 3.490 3.475 3.579 2.829 1.909 1.903 2.033
         2.360 2.601 3.054 3.386 3.553 3.468 3.187 2.723 2.686 2.821
         3.000 3.201 3.424 3.531 };

   call tsunimar(ar,v,nar,aic) data=y maxlag=5
        opt=({-1 1}) print=1;
   /*-- set up complex coefficient matrix --*/
   ar_cx = ar || j(nrow(ar),1,0);
   call tsroot(root) matin=ar_cx nar=nar nma=0 print=1;

In Figure 13.27, the roots and their lengths from the origin are shown. The roots are also stored in the matrix ROOT. All roots are within the unit circle, while the MOD values of the fourth and fifth roots appear to be sizable (0.9194).

Figure 13.26: Minimum AIC AR Estimation

SAS Output
line
                                    lag   ar_coef                         
 
                                      1 1.3003068                         
                                      2  -0.72328                         
                                      3 0.2421928                         
                                      4 -0.378757                         
                                      5 0.1377273                         
                               aic innovation_varinace                    
 
                         -318.6138           0.0490554                    
 


Figure 13.27: Roots of AR Characteristic Polynomial Equation

SAS Output
line
 
                      Roots of AR Characteristic Polynomial               
 
      I           Real      Imaginary       MOD(Z)    ATAN(I/R)       Degr
 
      1       -0.29755        0.55991       0.6341       2.0593     117.98
      2       -0.29755       -0.55991       0.6341      -2.0593    -117.98
      3        0.40529              0       0.4053            0           
      4        0.74505        0.53866       0.9194       0.6260      35.86
      5        0.74505       -0.53866       0.9194      -0.6260     -35.86
 
            Z**5-AR(1)*Z**4-AR(2)*Z**3-AR(3)*Z**2-AR(4)*Z**1-AR(5)=0      
 


The TSROOT subroutine can also recover the polynomial coefficients if the roots are given as an input. You should specify the QCOEF=1 option when you want to compute the polynomial coefficients instead of polynomial roots. You can compare the result with the preceding output of the TSUNIMAR call. Here is the code:

   call tsroot(ar_cx) matin=root nar=nar qcoef=1
                       nma=0 print=1;

The results are shown in Figure 13.28.

Figure 13.28: Polynomial Coefficients

SAS Output
line
 
                             Polynomial Coefficents                       
 
                         I       AR(real)       AR(imag)                  
 
                         1        1.30031              0                  
                         2       -0.72328    1.11022E-16                  
                         3        0.24219    8.32667E-17                  
                         4       -0.37876    2.77556E-17                  
                         5        0.13773              0