ARIMA differences in R and SAS

[Cpitson]

Hi everyone, 

I am working to understand some of the differences between SAS and R (using the fable, forecasting, or basic stats packages) to be able to run the same model in SAS and R. After being unable to replicate a more complex model, I tried to run a basic ARIMA (1,1,1) and work my way up. I know that this is not the right model for my data, the point of this is not to fit a model, rather to get the same forecasted values using SAS and R.

However, I am not understanding the SAS outcome - although the data has seasonality, I have not indicated that to SAS. But the forecasted values make me think that something is happening in SAS that I don't know about. Interestingly, when I add an intercept to my R model, I get similar values as the SAS model where I specify no intercept. 

I would appreciate any insight anyone has into this! I have read through a lot of R and SAS documentation, but have not found anything that sheds light on why this would happen. I have attached a screenshot of the forecasted values in SAS, R without an intercept, and R with an intercept. Thanks!

SAS Code: 

proc arima data=cld;

identify var=lgs1(1) nlag=13 minic scan esacf

stationarity=(adf=5)

crosscorr=(sspben)noprint;

estimate p=(1)q=(1) input=(sspben)

noint method=ml plot;

forecast lead=61 interval=month id=month alpha=0.05 out=s1out;

run;

R (fable) code without an intercept (with an intercept, replace ~0 with ~1):

arima111<-ow_input%>%

  model(arima = ARIMA(lgs1 ~ 0+ sspben +pdq(1,1,1)

                        +PDQ(0,0,0), method="ML"))

s1_forecast<-forecast(arima111, fc_data)%>%
rename(forecast=.mean)%>%
mutate(forecast=exp(forecast))%>%
as.data.frame()%>%
select(month, lgs1, forecast)

[sbxkoenk]

Hello,

in fact you fit a SARIMAX model (Seasonal ARIMA with regressors).

Can you try the statespace representation of your SARIMAX model?
You need PROC SSM (StateSpace Modelling) for that.

See here for an example :

SAS® 9.4 and SAS® Viya® 3.5 Programming Documentation
SAS/ETS 15.3 User's Guide
The SSM Procedure
Example 34.6 Model with Multiple ARIMA Components
https://go.documentation.sas.com/doc/en/pgmsascdc/9.4_3.5/etsug/etsug_ssm_examples06.htm

BR, Koen

[Cpitson]

Hi @sbxkoenk  thanks for your response! 

Interesting... so the proc arima command will identify seasonality even if you don't specify it? I thought that it would only incorporate seasonality if I were to add e.g. (12)., to indicate the seasonality. 

I am not understanding how to incorporate the arima into the SSM when I the input variable is in my dataset. I have tried a few versions to run the SSM (not too familar with them / SAS) and get different errors. For example, I get the error message "The name in the TREND statement must not be an analysis variable, such as a variable from the input data set."

proc ssm data=cld;
id date interval=month;
trend sspben(arma(p=1 d=1 q=1));
model lgs1= sspben;
run;

In the linked example, it looks like their trend statements are introducing new variables with arima models?

Thanks!

[sbxkoenk]

---

@Cpitson wrote:

I thought that it would only incorporate seasonality if I were to add e.g. (12)., to indicate the seasonality. 

That's how it works indeed.
Why do you think PROC ARIMA (SAS/ETS) is modelling seasonality when not asked to do so?

For a seasonal model on monthly data, you will indeed have to put that (12) explicitly.
Like here :

/*-- Seasonal Model for the Airline Series --*/
proc arima data=seriesg;
   identify var=xlog(1,12);
   estimate q=(1)(12) noint method=ml;
   forecast id=date interval=month printall out=b;
run;

With respect to statespace representation of SARIMAX using PROC SSM, I will check about including regressors (input variables).

BR,

Koen

[rselukar]

Hello,

ARIMAX model class specifications can be tedious and matching specs and results in two different languages is not easy.  Even if matching specs are created, the numbers may not match because the iterative model fitting phase may lead to different parameter estimates and then, in turn, different forecasts.  Therefore, for this illustration I am choosing a data set, Box & Jenkins's Series R, and a slightly simplified version of the Box & Tiao model (Example 4 in the ARIMA doc SAS Help Center: An Intervention Model for Ozone Data).  Please see the attached files, ozone.txt gives the R code, and ozone.sas gives the SAS code.  Essentially, the R spec:

fit = dt %>% model(arima=ARIMA(ozone~ 0 + x1 + pdq(0,0,1) + PDQ(0,1,1), method="ML"))

corresponds to

 identify var=ozone(12) crosscorr=(x1(12)) noprint;

estimate q=(1)(12) input=(x1) noint method=ml;

The parameter estimates and forecasts match well.

R:

Series: ozone
Model: LM w/ ARIMA(0,0,1)(0,1,1)[12] errors

Coefficients:
ma1 sma1 x1
0.3043 -0.6932 -1.2404
s.e. 0.0612 0.0662 0.2219

sigma^2 estimated as 0.6328: log likelihood=-245.25
AIC=498.5 AICc=498.7 BIC=511.77
Warning messages:
1: package ‘fable’ was built under R version 4.3.1
2: package ‘fabletools’ was built under R version 4.3.1
> fit %>% forecast(new_data=extra)
# A fable: 12 x 5 [1M]
# Key: .model [1]
.model date ozone .mean x1
<chr> <mth> <dist> <dbl> <dbl>
1 arima 1973 Jan N(1.6, 0.63) 1.55 1
2 arima 1973 Feb N(2.1, 0.69) 2.07 1
3 arima 1973 Mar N(2.7, 0.69) 2.72 1
4 arima 1973 Apr N(3, 0.69) 3.05 1
5 arima 1973 May N(3.4, 0.69) 3.40 1
6 arima 1973 Jun N(3.4, 0.69) 3.44 1
7 arima 1973 Jul N(4, 0.69) 4.01 1
8 arima 1973 Aug N(4.2, 0.69) 4.18 1
9 arima 1973 Sep N(3.6, 0.69) 3.58 1
10 arima 1973 Oct N(2.9, 0.69) 2.90 1
11 arima 1973 Nov N(2, 0.69) 1.97 1
12 arima 1973 Dec N(1.5, 0.69) 1.45 1

SAS:

Maximum Likelihood Estimation
Parameter	Estimate	Standard Error	t Value	Approx Pr > |t|	Lag	Variable	Shift
MA1,1	-0.30436	0.06618	-4.60	<.0001	1	ozone	0
MA2,1	0.69325	0.06273	11.05	<.0001	12	ozone	0
NUM1	-1.24034	0.21442	-5.78	<.0001	0	x1	0

Forecasts for variable ozone
Obs	Forecast	Std Error	95% Confidence Limits
217	1.5514	0.7955	-0.0076	3.1105
218	2.0651	0.8315	0.4354	3.6948
219	2.7236	0.8315	1.0939	4.3533
220	3.0485	0.8315	1.4188	4.6781
221	3.3960	0.8315	1.7663	5.0257
222	3.4374	0.8315	1.8077	5.0671
223	4.0125	0.8315	2.3828	5.6422
224	4.1797	0.8315	2.5500	5.8094
225	3.5754	0.8315	1.9457	5.2051
226	2.8963	0.8315	1.2666	4.5260
227	1.9678	0.8315	0.3381	3.5975
228	1.4521	0.8315	-0.1776	3.0818

[rselukar]

Some comments:

1. The sign convention for SAS and R ARMA parameters seems opposite.

2. In SAS/ETS, in addition to PROC ARIMA, you can do similar modeling with PROC UCM and PROC SSM.  Each has some advantages and disadvantages.