auto_arima

Chapter 8: Time Series and Forecasting

auto_arima

Automatically identifies time series outliers, determines parameters of a multiplicative seasonal ARIMA model and produces forecasts that incorporate the effects of outliers whose effects persist beyond the end of the series.

Synopsis

#include <imsls.h>

float *imsls_f_auto_arima (int n_obs, int tpoints[], float x[],...,0)

The type double function is imsls_d_auto_arima.

Required Arguments

int n_obs (Input)
Number of observations in the original time series. Assuming that the series is defined at time points , the actual length of the series, including missing observations is .

int tpoints[] (Input)
A vector of length n_obs containing the time points the time series was observed. It is required that are in strictly ascending order.

float x[] (Input)
A vector of length n_obs containing the observed time series values . This series can contain outliers and missing observations. Outliers are identified by this routine and missing values are identified by the time values in vector tpoints. If the time interval between two consecutive time points is greater than one, i.e. , then missing values are assumed to exist between and at times . Therefore, the gap free series is assumed to be defined for equidistant time points . Missing values are automatically estimated prior to identifying outliers and producing forecasts. Forecasts are generated for both missing and observed values.

Return Value

Pointer to an array of length 1 + p + q with the estimated constant, AR and MA parameters used to fit the outlier-free series using an ARIMA model. Upon completion, if d=model[3]=0, then an ARMA(p, q) model or AR(p) model is fitted to the outlier-free version of the observed series . If d=model[3]>0, these parameters are computed for an ARMA(p,q) representation of the seasonally adjusted series , where and s=model[2]³1.
If an error occurred, NULL is returned.

Synopsis with Optional Arguments

#include <imsls.h>

float *imsls_f_auto_arima (int n_obs, int tpoints[], float x[],
IMSLS_RETURN_USER, float parameters[],
IMSLS_METHOD, int method,
IMSLS_MAX_LAG, int maxlag,
IMSLS_MODEL, int model[],
IMSLS_DELTA, float delta,
IMSLS_CRITICAL, float critical,
IMSLS_EPSILON, float epsilon,
IMSLS_RESIDUAL, float **residual,
IMSLS_RESIDUAL_USER, float residual[],
IMSLS_RESIDUAL_SIGMA, float *res_sigma,
IMSLS_NUM_OUTLIERS, int *num_outliers,
IMSLS_P_INITIAL, int n_p_initial, int p_initial[],
IMSLS_Q_INITIAL, int n_q_initial, int q_initial[],
IMSLS_S_INITIAL, int n_s_initial, int s_initial[],
IMSLS_D_INITIAL, int n_d_initial, int d_initial[],
IMSLS_OUTLIER_STATISTICS, int **outlier_stat,
IMSLS_OUTLIER_STATISTICS_USER, int outlier_stat[],
IMSLS_AIC, float *aic,
IMSLS_AICC, float *aicc,
IMSLS_BIC, float *bic,
IMSLS_MODEL_SELECTION_CRITERION, int criterion,
IMSLS_OUT_FREE_SERIES, float **outfree_series,
IMSLS_OUT_FREE_SERIES_USER, float outfree_series[],
IMSLS_CONFIDENCE, float confidence,
IMSLS_NUM_PREDICT, int n_predict,
IMSLS_OUT_FREE_FORECAST, float **outfree_forecast,
IMSLS_OUT_FREE_FORECAST_USER, float outfree_forecast[],
IMSLS_OUTLIER_FORECAST, float **outlier_forecast,
IMSLS_OUTLIER_FORECAST_USER, float outlier_forecast[],
0)

Optional Arguments

IMSLS_METHOD, int method (Input)
The method used in model selection:
1 — Automatic ARIMA selection
2 — Grid search
Requires arguments IMSLS_P_INITIAL and IMSLS_Q_INITIAL.
3 — Specified ARIMA model
Requires argument IMSLS_MODEL.
Default: method = 1
For more information, see the Description section.

IMSLS_MAX_LAG, int maxlag (Input)
The maximum lag allowed when fitting an AR(p) model.
Default: maxlag = 10

IMSLS_MODEL, int model[] (Input/Output)
Array of length 4 containing the values for p, q, s, d. If method=3 is chosen, then the values for p and q must be defined. If IMSLS_S_INITIAL and IMSLS_D_INITIAL are not defined, then also s and d must be given. If method=1 or method=2, then model is ignored as an input array. On output, model contains the optimum values for p, q, s, d in model[0], model[1], model[2] and model[3], respectively.

IMSLS_DELTA, float delta (Input)
The dampening effect parameter used in the detection of a Temporary Change Outlier (TC), 0<delta<1.
Default: delta = 0.7

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

IMSLS_EPSILON, float epsilon (Input)
Positive tolerance value controlling the accuracy of parameter estimates during outlier detection.
Default: epsilon = 0.001

IMSLS_RESIDUAL, float **residual (Output)
Address of a pointer to an internally allocated array of length
, containing , the estimates of the white noise in the outlier free original series.

IMSLS_RESIDUAL_USER, float residual[] (Output)
Storage for array residual is provided by the user. See IMSLS_RESIDUAL.

IMSLS_RESIDUAL_SIGMA, float *res_sigma (Output)
Residual standard error (RSE) of the outlier free original series.

IMSLS_NUM_OUTLIERS, int *num_outliers (Output)
The number of outliers detected.

IMSLS_P_INITIAL, int n_p_initial, int p_initial[] (Input)
An array with n_p_initial elements containing the candidate values for p, from which the optimum is being selected. All candidate values in p_initial[] must be non-negative and n_p_initial ³ 1. If method=2, then IMSLS_P_INITIAL must be defined. Otherwise, n_p_initial and p_initial are ignored.

IMSLS_Q_INITIAL, int n_q_initial, int q_initial[] (Input)
An array with n_q_initial elements containing the candidate values for q, from which the optimum is being selected. All candidate values in q_initial[] must be non-negative and n_q_initial ³ 1. If method=2, then IMSLS_Q_INITIAL must be defined. Otherwise, n_q_initial and q_initial are ignored.

IMSLS_S_INITIAL, int n_s_initial, int s_initial[] (Input)
A vector of length n_s_initial containing the candidate values for s, from which the optimum is being selected. All candidate values in s_initial[] must be positive and n_s_initial ³ 1.
Default: n_s_initial=1, s_initial={1}

IMSLS_D_INITIAL, int n_d_initial, int d_initial[] (Input)
A vector of length n_d_initial containing the candidate values for d, from which the optimum is being selected. All candidate values in d_initial[] must be non-negative and n_d_initial ³ 1.
Default: n_d_initial=1, d_initial={0}

IMSLS_OUTLIER_STATISTICS, int **outlier_stat (Output)
Address of a pointer to an internally allocated array of length num_outliers 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:

0	Innovational Outliers (IO)
1	Additive Outliers (AO)
2	Level Shift Outliers (LS)
3	Temporary Change Outliers (TC)
4	Unable to Identify (UI).

If num_outliers=0, NULL is returned.

IMSLS_OUTLIER_STATISTICS_USER, int outlier_stat[] (Output)
A user allocated array of length n ´ 2 containing outlier statistics in its first num_outliers rows. Here, .
See IMSLS_OUTLIER_STATISTICS.
If num_outliers = 0, outlier_stat stays unchanged.

IMSLS_AIC, float *aic (Output)
The AIC (Akaike's Information Criterion ) value for the optimum model. Uses an approximation of the maximum log-likelihood based on an estimate of the innovation variance of the series.

IMSLS_AICC, float *aicc (Output)
The AICC (corrected AIC) value for the optimum model. Uses an approximation of the maximum log-likelihood based on an estimate of the innovation variance of the series.

IMSLS_BIC, float *bic (Output)
The BIC (Bayesian Information Criterion) value for the optimum model. Uses an approximation of the maximum log-likelihood based on an estimate of the innovation variance of the series.

IMSLS_MODEL_SELECTION_CRITERION, int criterion (Input)
The information criterion used for optimum model selection.

criterion	selected information criterion
0	Akaike's Information Criterion (AIC)
1	Akaike's Corrected Information Criterion (AICC)
2	Bayesian Information Criterion (BIC)

Default: criterion = 0.

IMSLS_OUT_FREE_SERIES, float **outfree_series (Output)
Address of a pointer to an internally allocated array of length n by 2, where . The first column of outfree_series contains the n_obs observations from the original series, , plus estimated values for any time gaps. The second column contains the same values as the first column adjusted by removing any outlier effects. In effect, the second column contains estimates of the underlying outlier-free series, . If no outliers are detected then both columns will contain identical values.

IMSLS_OUT_FREE_SERIES_USER, float outfree_series[] (Output)
A user allocated array of length n by 2, where . For further details, see IMSLS_OUT_FREE_SERIES.

IMSLS_CONFIDENCE, float confidence (Input)
Confidence level for computing forecast confidence limits, taken from the exclusive interval (0, 100). Typical choices for confidence are 90.0, 95.0 and 99.0.
Default: confidence = 95.0

IMSLS_NUM_PREDICT, int n_predict (Input)
The number of forecasts requested. Forecasts are made at origin , i.e. from the last observed value of the series.
Default: n_predict = 0

IMSLS_OUT_FREE_FORECAST, float **outfree_forecast (Output)
Address of a pointer to an internally allocated array of length n_predict by 3. The first column contains the forecasted values for the original outlier free series for t=+1, + 2,..., + n_predict. The second column contains standard errors for these forecasts, and the third column contains the psi weights of the infinite order moving average form of the model.

IMSLS_OUT_FREE_FORECAST_USER, float outfree_forecast[] (Output)
A user allocated array of length n_predict by 3. For more information, see IMSLS_OUT_FREE_FORECAST.

IMSLS_OUTLIER_FORECAST, float **outlier_forecast (Output)
Address of a pointer to an internally allocated array of length n_predict by 3. The first column contains the forecasted values for the original series for t=+1, +2,..., +n_predict. 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.

IMSLS_OUTLIER_FORECAST_USER, float outlier_forecast[] (Output)
A user allocated array of length n_predict by 3. For more information, see IMSLS_OUTLIER_FORECAST.

IMSLS_RETURN_USER, float x[] (Output)
A user allocated array containing the estimated constant, AR and MA parameters in its first 1+p+q locations. The values p and q can be estimated by upper bounds: If method=1, an upper bound for p would be maxlag, and q= 0. If method=2, upper bounds for p and q would be the maximum values in arrays p_initial and q_initial, respectively. If method=3,
p= model[0] and q= model[1].

Description

Function imsls_f_auto_arima determines the parameters of a multiplicative seasonal ARIMA model, and then uses the fitted model to identify outliers and prepare forecasts. The order of this model can be specified or automatically determined.
The ARIMA model handled by imsls_f_auto_arima has the following form:

where

and

It is assumed that all roots of and lie outside the unit circle. Clearly, if this reduces to the traditional ARIMA(p, d, q) model.

is the unobserved, outlier-free time series with mean , and white noise . This model is referred to as the underlying, outlier-free model. Function imsls_f_auto_arima does not assume that this series is observable. It assumes that the observed values might be contaminated by one or more outliers, whose effects are added to the underlying outlier-free series:

Outlier identification uses the algorithm developed by Chen and Liu (1993). Outliers are classified into 1 of 5 types:

1. innovational

2. additive

3. level shift

4. temporary change and

5. unable to identify

Once outliers are identified, imsls_f_auto_arima estimates , the outlier-free series representation of the data, by removing the estimated outlier effects.

Using the information about the adjusted ARIMA model and the removed outliers, forecasts are then prepared for the outlier-free series. Outlier effects are added to these forecasts to produce a forecast for the observed series, . If there are no outliers, then the forecasts for the outlier-free series and the observed series will be identical.

Model Selection

Users have an option of either specifying specific values for p, q , s and d or have imsls_f_auto_arima automatically select best fit values. Model selection can be conducted in one of three methods listed below depending upon the value of variable method.

Method 1: Automatic ARIMA Selection

This method initially searches for the AR(p) representation with minimum AIC for the noisy data, where p =0,...,maxlag.

If IMSLS_D_INITIAL is defined then the values in s_initial and d_initial are included in the search to find an optimum ARIMA representation of the series. Here, every possible combination of values for p, s in s_initial and d in d_initial is examined. The best found ARIMA representation is then used as input for the outlier detection routine.

The optimum values for p, q, s and d are returned in model[0], model[1], model[2] and model[3], respectively.

Method 2: Grid Search

The second automatic method conducts a grid search for p and q using all possible combinations of candidate values in p_initial and q_initial. Therefore, for this method the definition of IMSLS_P_INITIAL and IMSLS_Q_INITIAL is required.

If IMSLS_D_INITIAL is defined, the grid search is extended to include the candidate values for s and d given in s_initial and d_initial, respectively.

If IMSLS_D_INITIAL is not defined, no seasonal adjustment is attempted, and the grid search is restricted to searching for optimum values of p and q only.

The optimum values of p, q, s and d are returned in model[0], model[1], model[2] and model[3], respectively.

Method 3: Specified ARIMAModel

In the third method, specific values for p, q, s and d are given. The values for p and q must be defined in model[0] and model[1], respectively. If IMSLS_S_INITIAL and IMSLS_D_INITIAL are not defined, then values and must be specified in model[2] and model[3]. If IMSLS_S_INITIAL and IMSLS_D_INITIAL are defined, then a grid search for the optimum values of s and d is conducted using all possible combinations of input values in s_initial and d_initial. The optimum values of s and d can be found in model[2] and model[3], respectively.

Outliers

The algorithm of Chen and Liu (1993) is used to identify outliers. The number of outliers identified is returned in num_outliers. Both the time and classification for these outliers are returned in outlier_stat[]. Outliers are classified into one of five categories based upon the standardized statistic for each outlier type. The time at which the outlier occurred is given in the first column of outlier_stat. The outlier identifier returned in the second column is according to the descriptions in the following table:

Outlier Identifier	Name	General Description
0	(IO) Innovational Outlier	Innovational outliers persist. That is, there is an initial impact at the time the outlier occurs. This effect continues in a lagged fashion with all future observations. The lag coefficients are determined by the coefficient of the underlying ARIMA model.
1	(AO) Additive Outlier	Additive outliers do not persist. As the name implies, an additive outlier effects only the observation at the time the outlier occurs. Hence additive outliers have no effect on future forecasts.
2	(LS) Level Shift	Level shift outliers persist. They have the effect of either raising or lowering the mean of the series starting at the time the outlier occurs. This shift in the mean is abrupt and permanent.
3	(TC) Temporary Change	Temporary change outliers persist and are similar to level shift outliers with one major exception. Like level shift outliers, there is an abrupt change in the mean of the series at the time this outlier occurs. However, unlike level shift outliers, this shift is not permanent. The TC outlier gradually decays, eventually bringing the mean of the series back to its original value. The rate of this decay is modeled using the parameter delta. The default of delta= 0.7 is the value recommended for general use by Chen and Liu (1993).
4	(UI) Unable to Identify	If an outlier is identified as the last observation, then the algorithm is unable to determine the outlier's classification. For forecasting, a UI outlier is treated as an IO outlier. That is, its effect is lagged into the forecasts.

Except for additive outliers (AO), the effect of an outlier persists to observations following that outlier. Forecasts produced by imsls_f_auto_arima take this into account.

Examples

Example 1

This example uses time series LNU03327709 from the US Department of Labor, Bureau of Labor Statistics. It contains the unadjusted special unemployment rate, taken monthly from January 1994 through September 2005. The values 01/2004 – 03/2005 are used by imsls_f_auto_arima for outlier detection and parameter estimation. In this example, Method 1 without seasonal adjustment is chosen to find an appropriate AR(p) model. A forecast is done for the following six months and compared with the actual values 04/2005 – 09/2005.

Link to example source

#include <imsls.h>

#include <stdlib.h>

#include <stdio.h>

int main()

{

float *parameters = NULL, *outlier_forecast = NULL;

int *outlier_stat = NULL;

int n_obs, n_predict, i, num_outliers;

float aic, res_sigma;

int model[4];

float forecast_table[24];

float x[141] = {

12.8,12.2,11.9,10.9,10.6,11.3,11.1,10.4,10.0,9.7,9.7,9.7,

11.1,10.5,10.3,9.8,9.8,10.4,10.4,10.0,9.7,9.3,9.6,9.7,

10.8,10.7,10.3,9.7,9.5,10.0,10.0,9.3,9.0,8.8,8.9,9.2,

10.4,10.0,9.6,9.0,8.5,9.2,9.0,8.6,8.3,7.9,8.0,8.2,

9.3,8.9,8.9,7.7,7.6,8.4,8.5,7.8,7.6,7.3,7.2,7.3,

8.5,8.2,7.9,7.4,7.1,7.9,7.7,7.2,7.0,6.7,6.8,6.9,

7.8,7.6,7.4,6.6,6.8,7.2,7.2,7.0,6.6,6.3,6.8,6.7,

8.1,7.9,7.6,7.1,7.2,8.2,8.1,8.1,8.2,8.7,9.0,9.3,

10.5,10.1,9.9,9.4,9.2,9.8,9.9,9.5,9.0,9.0,9.4,9.6,

11.0,10.8,10.4,9.8,9.7,10.6,10.5,10.0,9.8,9.5,9.7,9.6,

10.9,10.3,10.4,9.3,9.3,9.8,9.8,9.3,8.9,9.1,9.1,9.1,

10.2,9.9,9.4,8.7,8.6,9.3,9.1,8.8,8.5};

int times[141] = {

1,2,3,4,5,6,7,8,9,10,11,12,

13,14,15,16,17,18,19,20,21,22,23,24,

25,26,27,28,29,30,31,32,33,34,35,36,

37,38,39,40,41,42,43,44,45,46,47,48,

49,50,51,52,53,54,55,56,57,58,59,60,

61,62,63,64,65,66,67,68,69,70,71,72,

73,74,75,76,77,78,79,80,81,82,83,84,

85,86,87,88,89,90,91,92,93,94,95,96,

97,98,99,100,101,102,103,104,105,106,107,108,

109,110,111,112,113,114,115,116,117,118,119,120,

121,122,123,124,125,126,127,128,129,130,131,132,

133,134,135,136,137,138,139,140,141};

n_predict = 6;

n_obs = 135;

parameters = imsls_f_auto_arima(n_obs, times, x, IMSLS_MODEL, model,

IMSLS_AIC, &aic,

IMSLS_MAX_LAG, 5,

IMSLS_CRITICAL, 4.0,

IMSLS_NUM_OUTLIERS, &num_outliers,

IMSLS_OUTLIER_STATISTICS, &outlier_stat,

IMSLS_RESIDUAL_SIGMA, &res_sigma,

IMSLS_NUM_PREDICT, n_predict,

IMSLS_OUTLIER_FORECAST, &outlier_forecast,

0);

printf("\nMethod 1: Automatic ARIMA model selection,"

" no differencing\n");

printf("\nModel chosen: p=%d, q=%d, s=%d, d=%d\n", model[0],

model[1], model[2], model[3]);

printf("\nNumber of outliers: %d\n\n", num_outliers);

printf("Outlier statistics:\n\n");

printf("Time point\t\tOutlier type\n");

for (i=0; i<num_outliers; i++)

printf("%d\t\t%d\n", outlier_stat[2*i], outlier_stat[2*i+1]);

printf("\nAIC = %lf\n", aic);

printf("RSE = %lf\n\n", res_sigma);

printf("Parameters:\n");

for (i=0; i<=model[0]+model[1]; i++)

printf("parameters[%d]=%lf\n", i, parameters[i]);

for (i=0; i<n_predict; i++)

{

forecast_table[4*i] = x[n_obs+i];

forecast_table[4*i+1] = outlier_forecast[3*i];

forecast_table[4*i+2] = outlier_forecast[3*i+1];

forecast_table[4*i+3] = outlier_forecast[3*i+2];

}

imsls_f_write_matrix("\t* * * Forecast Table * * *"

"\nOrig. series\t forecast\tprob. limits\tpsi weights\n",

n_predict, 4, forecast_table, IMSLS_WRITE_FORMAT, "%11.4f", 0);

if (parameters)

{

imsls_free(parameters);

parameters = NULL;

}

if (outlier_forecast)

{

imsls_free(outlier_forecast);

outlier_forecast = NULL;

}

if (outlier_stat)

{

imsls_free(outlier_stat);

outlier_stat = NULL;

}

return;

}

Output

Method 1: Automatic ARIMA model selection, no differencing

Model chosen: p=5, q=0, s=1, d=0

Number of outliers: 2

Outlier statistics:

Time point Outlier type

97 0

109 0

AIC = 418.118835

RSE = 0.428036

Parameters:

parameters[0]=0.166748

parameters[1]=0.857952

parameters[2]=-0.048224

parameters[3]=-0.239799

parameters[4]=0.200456

parameters[5]=0.211686

* * * Forecast Table * * *

Orig. series forecast prob. limits psi weights

1 2 3 4

1 8.7000 9.0586 0.8389 0.8580

2 8.6000 9.0823 1.1054 0.6879

3 9.3000 9.4117 1.2469 0.3090

4 9.1000 9.6112 1.2736 0.2266

5 8.8000 9.5866 1.2877 0.3983

6 8.5000 9.4094 1.3304 0.5762

Example 2

This is the same as Example 1 , except now imsls_f_auto_arima uses Method 2 with a possible seasonal adjustment. As a result, the unadjusted model with is chosen as optimum.

Link to example source

#include <imsls.h>

#include <stdlib.h>

#include <stdio.h>

int main()

{

int n_obs, n_predict, i, num_outliers;

float aic, res_sigma;

int model[4];

int n_s_initial = 2;

int n_d_initial = 3;

int s_initial[2] = {1,2};

int d_initial[3] = {0,1,2};

int n_p_initial = 4, n_q_initial = 4;

int p_initial[4] = {0,1,2,3};

int q_initial[4] = {0,1,2,3};

float parameters_user[141];

float outfree_series_user[282];

int outlier_stat_user[282];

float outlier_forecast_user[24];

float forecast_table[24];

float x[141] = {

12.8,12.2,11.9,10.9,10.6,11.3,11.1,10.4,10.0,9.7,9.7,9.7,

11.1,10.5,10.3,9.8,9.8,10.4,10.4,10.0,9.7,9.3,9.6,9.7,

10.8,10.7,10.3,9.7,9.5,10.0,10.0,9.3,9.0,8.8,8.9,9.2,

10.4,10.0,9.6,9.0,8.5,9.2,9.0,8.6,8.3,7.9,8.0,8.2,

9.3,8.9,8.9,7.7,7.6,8.4,8.5,7.8,7.6,7.3,7.2,7.3,

8.5,8.2,7.9,7.4,7.1,7.9,7.7,7.2,7.0,6.7,6.8,6.9,

7.8,7.6,7.4,6.6,6.8,7.2,7.2,7.0,6.6,6.3,6.8,6.7,

8.1,7.9,7.6,7.1,7.2,8.2,8.1,8.1,8.2,8.7,9.0,9.3,

10.5,10.1,9.9,9.4,9.2,9.8,9.9,9.5,9.0,9.0,9.4,9.6,

11.0,10.8,10.4,9.8,9.7,10.6,10.5,10.0,9.8,9.5,9.7,9.6,

10.9,10.3,10.4,9.3,9.3,9.8,9.8,9.3,8.9,9.1,9.1,9.1,

10.2,9.9,9.4,8.7,8.6,9.3,9.1,8.8,8.5};

int times[141] = {

1,2,3,4,5,6,7,8,9,10,11,12,

13,14,15,16,17,18,19,20,21,22,23,24,

25,26,27,28,29,30,31,32,33,34,35,36,

37,38,39,40,41,42,43,44,45,46,47,48,

49,50,51,52,53,54,55,56,57,58,59,60,

61,62,63,64,65,66,67,68,69,70,71,72,

73,74,75,76,77,78,79,80,81,82,83,84,

85,86,87,88,89,90,91,92,93,94,95,96,

97,98,99,100,101,102,103,104,105,106,107,108,

109,110,111,112,113,114,115,116,117,118,119,120,

121,122,123,124,125,126,127,128,129,130,131,132,

133,134,135,136,137,138,139,140,141};

n_predict = 6;

n_obs = 135;

imsls_f_auto_arima(n_obs, times, x, IMSLS_MODEL, model,

IMSLS_AIC, &aic,

IMSLS_CRITICAL, 4.0,

IMSLS_MAX_LAG, 5,

IMSLS_METHOD, 2,

IMSLS_P_INITIAL, n_p_initial, p_initial,

IMSLS_Q_INITIAL, n_q_initial, q_initial,

IMSLS_S_INITIAL, n_s_initial, s_initial,

IMSLS_D_INITIAL, n_d_initial, d_initial,

IMSLS_NUM_OUTLIERS, &num_outliers,

IMSLS_OUTLIER_STATISTICS_USER, outlier_stat_user,

IMSLS_RESIDUAL_SIGMA, &res_sigma,

IMSLS_NUM_PREDICT, 6,

IMSLS_OUTLIER_FORECAST_USER, outlier_forecast_user,

IMSLS_RETURN_USER, parameters_user,

0);

for (i=0; i<n_predict; i++)

{

forecast_table[4*i] = x[n_obs+i];

forecast_table[4*i+1] = outlier_forecast_user[3*i];

forecast_table[4*i+2] = outlier_forecast_user[3*i+1];

forecast_table[4*i+3] = outlier_forecast_user[3*i+2];

}

printf("\nMethod 2: Grid search, differencing allowed\n");

printf("\nModel chosen: p=%d, q=%d, s=%d, d=%d\n", model[0],

model[1], model[2], model[3]);

printf("\nNumber of outliers: %d\n\n", num_outliers);

printf("Outlier statistics:\n\n");

printf("Time point\t\tOutlier type\n");

for (i=0; i<num_outliers; i++)

printf("%d\t\t%d\n", outlier_stat_user[2*i],

outlier_stat_user[2*i+1]);

printf("\nAIC = %lf\n", aic);

printf("RSE = %lf\n\n", res_sigma);

printf("Parameters:\n");

for (i=0; i<=model[0]+model[1]; i++)

printf("parameters[%d]=%lf\n", i, parameters_user[i]);

imsls_f_write_matrix("\n\t* * * Forecast Table * * *"

"\nOrig. series\t forecast\tprob. limits\tpsi weights\n",

n_predict, 4, forecast_table, IMSLS_WRITE_FORMAT, "%11.4f", 0);

return;

}

Output

Method 2: Grid search, differencing allowed

Model chosen: p=3, q=2, s=1, d=0

Number of outliers: 1

Outlier statistics:

Time point Outlier type

109 0

AIC = 408.076813

RSE = 0.412409

Parameters:

parameters[0]=0.509478

parameters[1]=1.944665

parameters[2]=-1.901104

parameters[3]=0.901657

parameters[4]=1.113017

parameters[5]=-0.914998

* * * Forecast Table * * *

Orig. series forecast prob. limits psi weights

1 2 3 4

1 8.7000 9.1109 0.8083 0.8316

2 8.6000 9.1811 1.0513 0.6312

3 9.3000 9.5185 1.1686 0.5480

4 9.1000 9.7804 1.2497 0.6157

5 8.8000 9.7117 1.3451 0.7245

6 8.5000 9.3842 1.4671 0.7326

Example 3

This example is the same as Example 2 but now Method 3 with the optimum model parameters from Example 2 are chosen for outlier detection and forecasting.

Link to example source

#include <imsls.h>

#include <stdlib.h>

#include <stdio.h>

int main()

{

float *parameters = NULL, *outlier_forecast = NULL;

int *outlier_stat = NULL;

int n_obs, n_predict, i, num_outliers;

float aic, res_sigma;

int model[4];

float forecast_table[24];

float x[141] = {

12.8,12.2,11.9,10.9,10.6,11.3,11.1,10.4,10.0,9.7,9.7,9.7,

11.1,10.5,10.3,9.8,9.8,10.4,10.4,10.0,9.7,9.3,9.6,9.7,

10.8,10.7,10.3,9.7,9.5,10.0,10.0,9.3,9.0,8.8,8.9,9.2,

10.4,10.0,9.6,9.0,8.5,9.2,9.0,8.6,8.3,7.9,8.0,8.2,

9.3,8.9,8.9,7.7,7.6,8.4,8.5,7.8,7.6,7.3,7.2,7.3,

8.5,8.2,7.9,7.4,7.1,7.9,7.7,7.2,7.0,6.7,6.8,6.9,

7.8,7.6,7.4,6.6,6.8,7.2,7.2,7.0,6.6,6.3,6.8,6.7,

8.1,7.9,7.6,7.1,7.2,8.2,8.1,8.1,8.2,8.7,9.0,9.3,

10.5,10.1,9.9,9.4,9.2,9.8,9.9,9.5,9.0,9.0,9.4,9.6,

11.0,10.8,10.4,9.8,9.7,10.6,10.5,10.0,9.8,9.5,9.7,9.6,

10.9,10.3,10.4,9.3,9.3,9.8,9.8,9.3,8.9,9.1,9.1,9.1,

10.2,9.9,9.4,8.7,8.6,9.3,9.1,8.8,8.5};

int times[141] = {

1,2,3,4,5,6,7,8,9,10,11,12,

13,14,15,16,17,18,19,20,21,22,23,24,

25,26,27,28,29,30,31,32,33,34,35,36,

37,38,39,40,41,42,43,44,45,46,47,48,

49,50,51,52,53,54,55,56,57,58,59,60,

61,62,63,64,65,66,67,68,69,70,71,72,

73,74,75,76,77,78,79,80,81,82,83,84,

85,86,87,88,89,90,91,92,93,94,95,96,

97,98,99,100,101,102,103,104,105,106,107,108,

109,110,111,112,113,114,115,116,117,118,119,120,

121,122,123,124,125,126,127,128,129,130,131,132,

133,134,135,136,137,138,139,140,141};

n_predict = 6;

n_obs = 135;

model[0] = 3;

model[1] = 2;

model[2] = 1;

model[3] = 0;

parameters = imsls_f_auto_arima(n_obs, times, x, IMSLS_MODEL, model,

IMSLS_AIC, &aic,

IMSLS_CRITICAL, 4.0,

IMSLS_METHOD, 3,

IMSLS_NUM_OUTLIERS, &num_outliers,

IMSLS_OUTLIER_STATISTICS, &outlier_stat,

IMSLS_RESIDUAL_SIGMA, &res_sigma,

IMSLS_NUM_PREDICT, 6,

IMSLS_OUTLIER_FORECAST, &outlier_forecast,

0);

printf("\nMethod 3: Specified ARIMA model\n");

printf("\nModel: p=%d, q=%d, s=%d, d=%d\n", model[0], model[1],

model[2], model[3]);

printf("\nNumber of outliers: %d\n\n", num_outliers);

printf("Outlier statistics:\n\n");

printf("Time point\t\tOutlier type\n");

for (i=0; i<num_outliers; i++)

printf("%d\t\t%d\n", outlier_stat[2*i], outlier_stat[2*i+1]);

printf("\nAIC = %lf\n", aic);

printf("RSE = %lf\n", res_sigma);

printf("\nParameters:\n");

for (i=0; i<=model[0]+model[1]; i++)

printf("parameters[%d]=%lf\n", i, parameters[i]);

for (i=0; i<n_predict; i++)

{

forecast_table[4*i] = x[n_obs+i];

forecast_table[4*i+1] = outlier_forecast[3*i];

forecast_table[4*i+2] = outlier_forecast[3*i+1];

forecast_table[4*i+3] = outlier_forecast[3*i+2];

}

imsls_f_write_matrix("\t* * * Forecast Table * * *"

"\nOrig. series\t forecast\tprob. limits\tpsi weights\n",

n_predict, 4, forecast_table, IMSLS_WRITE_FORMAT, "%11.4f", 0);

if (parameters)

{

imsls_free(parameters);

parameters = NULL;

}

if (outlier_forecast)

{

imsls_free(outlier_forecast);

outlier_forecast = NULL;

}

if (outlier_stat)

{

imsls_free(outlier_stat);

outlier_stat = NULL;

}

return;

}

Output

Method 3: Specified ARIMA model

Model: p=3, q=2, s=1, d=0

Number of outliers: 1

Outlier statistics:

Time point Outlier type

109 0

AIC = 408.076813

RSE = 0.412409

Parameters:

parameters[0]=0.509478

parameters[1]=1.944665

parameters[2]=-1.901104

parameters[3]=0.901657

parameters[4]=1.113017

parameters[5]=-0.914998

* * * Forecast Table * * *

Orig. series forecast prob. limits psi weights

1 2 3 4

1 8.7000 9.1109 0.8083 0.8316

2 8.6000 9.1811 1.0513 0.6312

3 9.3000 9.5185 1.1686 0.5480

4 9.1000 9.7804 1.2497 0.6157

5 8.8000 9.7117 1.3451 0.7245

6 8.5000 9.3842 1.4671 0.7326

Visual Numerics, Inc.
Visual Numerics - Developers of IMSL and PV-WAVE
http://www.vni.com/
PHONE: 713.784.3131
FAX:713.781.9260