I'm not entirely familiar with the results of the qs() function in R's package seasonal. But, by what I understand roughly, the results below show me that the data does not show any seasonality in the first place -- if I am wrong, please correct me. But could someone give me an idea of how the variable 'qsorievadj' is constructed and how these tests are performed in R within the seasonal package?

> qs(var)
                 qs   p-val
qsori       0.27402 0.87196
qsorievadj  0.00000 1.00000
qsrsd       0.00000 1.00000
qssadj      0.27402 0.87196
qssadjevadj 0.00000 1.00000
qsirr       0.00000 1.00000
qsirrevadj  0.00000 1.00000
qssrsd      0.00000 1.00000

These tests are not actually computed in R. The seasonal package is an interface to X-13-ARIMA-SEATS, which is where the computations are performed.

The QS statistics are based on a draft research report by Maravall (2012, Research Report, Banco de España). Details are also available in the documentation of X-13-ARIMA-SEATS.

These statistics are relatively straightforward to implement in R. Below I give you the general sketch to compute the statistic. (I give an example of the QS statistic computed for the original series).

# run the whole procedure for the "AirPassengers" series
m <- seas(AirPassengers)

The steps shown below follow the description given in the documentation of X-13-ARIMA-SEATS. First, the series for which the QS statistic is computed is differenced according to the chosen ARIMA model and the following rule:

$$ ndif = \max(1, \min(d + D, 2)) $$

where $ndif$ is the number of regular differences to be taken; $d$ and $D$ are respectively the number of regular and seasonal differences in the chosen ARIMA model. (If the QS statistic is computed for the series of residuals no differences are applied, $ndif=0$.)

x <- AirPassengers
S <- frequency(x)
ndif <- max(1, min(sum(arimamodel(m)[c(2,5)]), 2))
dx <- filter(x, polynomial(c(1,-1))^ndif, sides=1)
dx <- window(dx, start=time(x)[ndif+1])
# alternatively, we can do it without package "polynom" by simply doing,
# if ndiff=1:
#dx <- diff(x)
# if ndiff=2:
#dx <- diff(diff(x)) 

Next, the first two autocorrelations of seasonal order (e.g. 12 and 24 in monthly data) are obtained. If these autocorrelations are lower or equal to zero, then they are set to zero.

R <- acf(dx, lag.max=S*2, plot=FALSE)$acf[-1,,1][c(S, 2*S)]
if (R[1] <= 0)
  R[1] <- 0
if (R[2] <= 0)
  R[2] <- 0

The statistic is defined as follows ($R_s$ and $R_{2s}$ denote the autocorrelations obtained in the previous step):

$$ QS = n(n+2)\left( \frac{R_s^2}{n-s} + \frac{R_{2s}^2}{n-2s} \right) \,, $$

where $n$ is the number of observations in the differenced series and $s$ is the periodicity of the data (12 in this case with monthly data). According to simulations exercises, this statistic follows approximately the $\chi^2$ distribuion with 2 degrees of freedom.

Thus, the statistic and the corresponding p-value can be obtained as follows:

n <- length(dx)
QS <- n*(n+2)*(R[1]^2/(n-S) + R[2]^2/(n-2*S))
pvalue <- pchisq(q=QS, df=2, lower.tail=FALSE)
round(c(QS=QS, p.value=pvalue), 4)
#       QS  p.value 
# 167.6486   0.0000 

which agrees with the output returned by seasonal::qs:

#             qs p-val
# qsori 167.6486     0

In order to apply this test to other series, e.g. qsorievadj, you just need to take this series from the output files and apply the operations shown above to it.

Given the large value of the QS statistic (and the implied low p-value), we can conclude that there is seasonality in the series.

  • $\begingroup$ And as for the output that I've showed I can conclude there is no significant seasonality? $\endgroup$ – John Doe Apr 28 '15 at 11:57
  • $\begingroup$ @JohnDoe The QS statistics are low. In particular, taking a 5% significance level, the p-values (larger than 0.05) suggest that the there isn't significant seasonality in these series. $\endgroup$ – javlacalle Apr 28 '15 at 16:25
  • $\begingroup$ So if I do have a good fit of the ARIMA model and the qs statistics do not show me any significant seasonality in the original series, I should use the original data (if possible) for further analysis? Also, for some of the variables that I get this result, I still get some significant regressions of the original series with seasonal dummies (I reject that their joint significance is zero in the original data). I posted a question about this here: stats.stackexchange.com/questions/151657/… $\endgroup$ – John Doe May 10 '15 at 22:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.