/* nag_tsa_multi_filter_arima (g13bac) Example Program.
*
* Copyright 2023 Numerical Algorithms Group.
*
* Mark 29.2, 2023.
*/
#include <nag.h>
#include <stdio.h>
int main(void) {
/* Scalars */
double a1, a2, cx, cy;
Integer i, idd, ii, ij, iqxd, j, k, n, nb, ni, npar, nparx, nx, ny, nser,
npara, tdxxy, tdmrx, ldparx, tdparx;
Integer exit_status = 0;
/* Arrays */
double *b = 0, *fsd = 0, *fva = 0, *par = 0, *parx = 0, *x = 0, *y = 0,
*rms = 0, *parxx = 0;
Integer mr[14], mrx[7], *mrxx = 0;
Nag_ArimaOrder arimaj, arimaf, arimav;
Nag_TransfOrder transfj;
Nag_G13_Opt options;
NagError fail;
INIT_FAIL(fail);
exit_status = 0;
/* Initialize the options structure used by nag_tsa_multi_inputmod_forecast
* (g13bjc) */
/* nag_tsa_options_init (g13bxc).
* Initialization function for option setting
*/
nag_tsa_options_init(&options);
printf("nag_tsa_multi_filter_arima (g13bac) Example Program Results\n");
/* Skip heading in data file */
scanf("%*[^\n] ");
scanf("%" NAG_IFMT "%*[^\n] ", &nx);
if (nx > 0) {
/* Allocate array x */
if (!(x = NAG_ALLOC(nx + 2, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 1; i <= nx; ++i)
scanf("%lf", &x[i - 1]);
scanf("%*[^\n] ");
/* Read univariate Arima for series */
for (i = 1; i <= 7; ++i)
scanf("%" NAG_IFMT "", &mrx[i - 1]);
scanf("%*[^\n] ");
scanf("%lf%*[^\n] ", &cx);
nparx = mrx[0] + mrx[2] + mrx[3] + mrx[5];
arimaj.p = mrx[0];
arimaj.d = mrx[1];
arimaj.q = mrx[2];
arimaj.bigp = mrx[3];
arimaj.bigd = mrx[4];
arimaj.bigq = mrx[5];
arimaj.s = mrx[6];
nser = 1;
if (nparx > 0) {
/* Allocate array parx */
if (!(parx = NAG_ALLOC(nparx + nser, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 1; i <= nparx; ++i)
scanf("%lf", &parx[i - 1]);
scanf("%*[^\n] ");
/* Read model by which to filter series */
for (i = 1; i <= 7; ++i)
scanf("%" NAG_IFMT "", &mr[i - 1]);
scanf("%*[^\n] ");
arimaf.p = mr[0];
arimaf.d = mr[1];
arimaf.q = mr[2];
arimaf.bigp = mr[3];
arimaf.bigd = mr[4];
arimaf.bigq = mr[5];
arimaf.s = mr[6];
npar = mr[0] + mr[2] + mr[3] + mr[5];
if (npar > 0) {
/* Allocate array par */
if (!(par = NAG_ALLOC(npar + nparx, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
for (i = 1; i <= npar; ++i)
scanf("%lf", &par[i - 1]);
scanf("%*[^\n] ");
/* Initially backforecast QY values */
/* (1) Reverse series in situ */
n = nx / 2;
ni = nx;
for (i = 1; i <= n; ++i) {
a1 = x[i - 1];
a2 = x[ni - 1];
x[i - 1] = a2;
x[ni - 1] = a1;
--ni;
}
idd = mrx[1] + mrx[4];
/* (2) Possible sign reversal for ARIMA constant */
if (idd % 2 != 0)
cx = -cx;
/* (3) Calculate number of backforecasts required */
iqxd = mrx[2] + mrx[5] * mrx[6];
/* Calculate series length */
ny = nx + iqxd;
/* Allocate array y */
if (!(y = NAG_ALLOC(ny, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
if (iqxd != 0) {
/* Allocate arrays fsd, fva and st. */
if (!(fsd = NAG_ALLOC(iqxd, double)) ||
!(fva = NAG_ALLOC(iqxd, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
/* (4) Set up parameter list for call to forecast
* routine g13bjc
*/
npara = nparx + nser;
parx[npara - 1] = cx;
tdxxy = nser;
tdmrx = nser - 1;
ldparx = nser - 1;
tdparx = nser - 1;
if (!(rms = NAG_ALLOC(nser, double)) ||
!(parxx = NAG_ALLOC(nser, double)) ||
!(mrxx = NAG_ALLOC(7 * nser, Integer))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
/* nag_tsa_transf_orders (g13byc).
* Allocates memory to transfer function model orders
*/
nag_tsa_transf_orders(nser, &transfj, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_tsa_transf_orders (g13byc)"
".\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
rms[0] = 0;
transfj.nag_b = 0;
transfj.nag_q = 0;
transfj.nag_p = 0;
transfj.nag_r = 1;
for (i = 1; i <= 7; ++i)
mrxx[i - 1] = 0;
parxx[0] = 0;
/* Tell nag_tsa_multi_inputmod_forecast (g13bjc) not to
* print parameters on entry */
options.list = Nag_FALSE;
/* nag_tsa_multi_inputmod_forecast (g13bjc).
* Forecasting function
*/
nag_tsa_multi_inputmod_forecast(
&arimaj, nser, &transfj, parx, npara, nx, iqxd, x, tdxxy, rms,
mrxx, tdmrx, parxx, ldparx, tdparx, fva, fsd, &options, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_tsa_multi_inputmod_forecast "
"(g13bjc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
j = iqxd;
for (i = 1; i <= iqxd; ++i) {
y[i - 1] = fva[j - 1];
--j;
}
/* Move series into y */
j = iqxd + 1;
k = nx;
for (i = 1; i <= nx; ++i) {
if (j > 305)
goto END;
y[j - 1] = x[k - 1];
++j;
--k;
}
}
/* Move ARIMA for series into mr */
for (i = 1; i <= 7; ++i)
mr[i + 6] = mrx[i - 1];
arimav.p = mr[7];
arimav.d = mr[8];
arimav.q = mr[9];
arimav.bigp = mr[10];
arimav.bigd = mr[11];
arimav.bigq = mr[12];
arimav.s = mr[13];
/* Move parameters of ARIMA for y into par */
for (i = 1; i <= nparx; ++i)
par[npar + i - 1] = parx[i - 1];
npar += nparx;
/* Move constant and reset sign reversal */
cy = cx;
if (idd % 2 != 0)
cy = -cy;
/* Set parameters for call to filter routine
* nag_tsa_multi_filter_arima (g13bac) */
nb = ny + MAX(mr[2] + mr[5] * mr[6],
mr[0] + mr[1] + (mr[3] + mr[4]) * mr[6]);
/* Allocate array b */
if (!(b = NAG_ALLOC(nb, double))) {
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
/* Filter series by call to nag_tsa_multi_filter_arima (g13bac) */
/* nag_tsa_multi_filter_arima (g13bac).
* Multivariate time series, filtering (pre-whitening) by an
* ARIMA model
*/
nag_tsa_multi_filter_arima(y, ny, &arimaf, &arimav, par, npar, cy, b,
nb, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_tsa_multi_filter_arima (g13bac).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf("\n");
printf(" Original Filtered\n");
printf("Backforecasts y-series series\n");
if (iqxd != 0) {
ij = -iqxd;
for (i = 1; i <= iqxd; ++i) {
printf("%8" NAG_IFMT "%17.4f%15.4f\n", ij, y[i - 1], b[i - 1]);
++ij;
}
}
printf("\n");
printf(" Filtered Filtered "
"Filtered Filtered\n");
printf(" series series "
" series series\n");
for (i = iqxd + 1; i <= ny; i += 4) {
for (ii = i; ii <= MIN(ny, i + 3); ++ii) {
printf("%5" NAG_IFMT "", ii - iqxd);
printf("%9.4f ", b[ii - 1]);
}
printf("\n");
}
}
}
}
END:
/* Free the options structure used by nag_tsa_multi_inputmod_forecast
* (g13bjc) */
/* nag_tsa_free (g13xzc).
* Freeing function for use with g13 option setting
*/
nag_tsa_free(&options);
NAG_FREE(b);
NAG_FREE(fsd);
NAG_FREE(fva);
NAG_FREE(par);
NAG_FREE(parx);
NAG_FREE(x);
NAG_FREE(y);
NAG_FREE(rms);
NAG_FREE(parxx);
NAG_FREE(mrxx);
return exit_status;
}