Utility function to replicate td and td1coef, using xreg

Key Questions

Can we replicate X-13 trading days adjustment
Can we include the regressors in seasonal, with an example

Conclusion

Replication is possible. Adding the regressors to seasonal may be useful.

suppressPackageStartupMessages(library(tidyverse))
library(tsbox)
library(seasonal)
#> 
#> Attaching package: 'seasonal'
#> The following object is masked from 'package:tibble':
#> 
#>     view

Examples

Constructing weekday regressors

dates <- seq(as.Date("1931-01-01"), as.Date("2030-12-31"), by = "day")

first_of_month <- function(x) {
  as.Date(paste(
    data.table::year(dates),
    data.table::month(dates),
    1,
    sep = "-"
  ))
}

td1nolpyear
Include the weekday-weekend contrast variable (monthly and quarterly flow eries only): (no. of weekdays) −(5/2) (no. of Saturdays and Sundays).
tdnolpyear
Include the six day-of-week contrast variables (monthly and quarterly flow series only): (no. of Mondays) − (no. of Sundays), . . . , (no. of Saturdays) − (no. of Sundays).

td_m_tbl <-
  tibble(dates, wd = as.POSIXlt(dates)$wday, name = weekdays(dates)) %>%
  group_by(time = first_of_month(dates)) %>%
  summarize(
    td1 = sum(wd %in% 1:5) - 5 / 2 * sum(wd %in% c(6, 0)),
    mon = sum(wd == 1) - sum(wd == 0),
    tue = sum(wd == 2) - sum(wd == 0),
    wed = sum(wd == 3) - sum(wd == 0),
    thu = sum(wd == 4) - sum(wd == 0),
    fri = sum(wd == 5) - sum(wd == 0),
    sat = sum(wd == 6) - sum(wd == 0)
  )

td_m_tbl
#> # A tibble: 1,200 x 8
#>    time         td1   mon   tue   wed   thu   fri   sat
#>    <date>     <dbl> <int> <int> <int> <int> <int> <int>
#>  1 1931-01-01  -0.5     0     0     0     1     1     1
#>  2 1931-02-01   0       0     0     0     0     0     0
#>  3 1931-03-01  -0.5     0     0    -1    -1    -1    -1
#>  4 1931-04-01   2       0     0     1     1     0     0
#>  5 1931-05-01  -4      -1    -1    -1    -1     0     0
#>  6 1931-06-01   2       1     1     0     0     0     0
#>  7 1931-07-01   3       0     0     1     1     1     0
#>  8 1931-08-01  -4       0    -1    -1    -1    -1     0
#>  9 1931-09-01   2       0     1     1     0     0     0
#> 10 1931-10-01  -0.5     0     0     0     1     1     1
#> # … with 1,190 more rows

‘Trading day adjustment’ removes the effect of the weekdays, and but does not include holidays, such as Christmas or Easter. These are handled separately (Easter) or dealt with by standard seasonal adjustment (Christmas).

td1nolpyear <-
  td_m_tbl %>%
  select(time, value = td1) %>%
  ts_ts()

tdnolpyear <-
  td_m_tbl %>%
  select(-td1) %>%
  ts_long() %>%
  ts_ts()

Single coef

m1 <- seas(
  AirPassengers,
  xreg = td1nolpyear,
  regression.aictest = NULL,
  outlier = NULL,
  regression.usertype = "td"
)
summary(m1)
#> 
#> Call:
#> seas(x = AirPassengers, xreg = td1nolpyear, regression.aictest = NULL, 
#>     outlier = NULL, regression.usertype = "td")
#> 
#> Coefficients:
#>                     Estimate Std. Error z value Pr(>|z|)    
#> td1nolpyear       -0.0025474  0.0006732  -3.784 0.000154 ***
#> MA-Nonseasonal-01  0.3292278  0.0813633   4.046 5.20e-05 ***
#> MA-Seasonal-12     0.5695911  0.0739360   7.704 1.32e-14 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> SEATS adj.  ARIMA: (0 1 1)(0 1 1)  Obs.: 144  Transform: log
#> AICc: 976.6, BIC: 987.7  QS (no seasonality in final):    0  
#> Box-Ljung (no autocorr.): 24.86   Shapiro (normality): 0.9805 *

m2 <- seas(
  AirPassengers,
  regression.aictest = NULL,
  outlier = NULL,
  regression.variables = c("td1nolpyear", outlier = NULL)
)
summary(m2)
#> 
#> Call:
#> seas(x = AirPassengers, regression.aictest = NULL, outlier = NULL, 
#>     regression.variables = c("td1nolpyear", outlier = NULL))
#> 
#> Coefficients:
#>                     Estimate Std. Error z value Pr(>|z|)    
#> Weekday           -0.0025474  0.0006732  -3.784 0.000154 ***
#> MA-Nonseasonal-01  0.3292278  0.0813633   4.046 5.20e-05 ***
#> MA-Seasonal-12     0.5695911  0.0739360   7.704 1.32e-14 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> SEATS adj.  ARIMA: (0 1 1)(0 1 1)  Obs.: 144  Transform: log
#> AICc: 976.6, BIC: 987.7  QS (no seasonality in final):    0  
#> Box-Ljung (no autocorr.): 24.86   Shapiro (normality): 0.9805 *

All coefs

m3 <- seas(
  AirPassengers,
  xreg = tdnolpyear,
  regression.aictest = NULL,
  outlier = NULL,
  regression.usertype = "td"
)
summary(m3)
#> 
#> Call:
#> seas(x = AirPassengers, xreg = tdnolpyear, regression.aictest = NULL, 
#>     outlier = NULL, regression.usertype = "td")
#> 
#> Coefficients:
#>                    Estimate Std. Error z value Pr(>|z|)    
#> xreg1             -0.004982   0.004731  -1.053 0.292359    
#> xreg2             -0.004589   0.004762  -0.964 0.335172    
#> xreg3             -0.001612   0.004745  -0.340 0.734094    
#> xreg4             -0.003817   0.004680  -0.816 0.414652    
#> xreg5              0.003958   0.004706   0.841 0.400272    
#> xreg6              0.003164   0.004829   0.655 0.512342    
#> MA-Nonseasonal-01  0.298942   0.082193   3.637 0.000276 ***
#> MA-Seasonal-12     0.579965   0.073855   7.853 4.07e-15 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> SEATS adj.  ARIMA: (0 1 1)(0 1 1)  Obs.: 144  Transform: log
#> AICc: 983.9, BIC:  1008  QS (no seasonality in final):    0  
#> Box-Ljung (no autocorr.): 29.43   Shapiro (normality): 0.9781 *

m4 <- seas(
  AirPassengers,
  regression.aictest = NULL,
  outlier = NULL,
  regression.variables = c("tdnolpyear", outlier = NULL)
)
summary(m4)
#> 
#> Call:
#> seas(x = AirPassengers, regression.aictest = NULL, outlier = NULL, 
#>     regression.variables = c("tdnolpyear", outlier = NULL))
#> 
#> Coefficients:
#>                    Estimate Std. Error z value Pr(>|z|)    
#> Mon               -0.004982   0.004731  -1.053 0.292359    
#> Tue               -0.004589   0.004762  -0.964 0.335172    
#> Wed               -0.001612   0.004745  -0.340 0.734094    
#> Thu               -0.003817   0.004680  -0.816 0.414652    
#> Fri                0.003958   0.004706   0.841 0.400272    
#> Sat                0.003164   0.004829   0.655 0.512342    
#> MA-Nonseasonal-01  0.298942   0.082193   3.637 0.000276 ***
#> MA-Seasonal-12     0.579965   0.073855   7.853 4.07e-15 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> SEATS adj.  ARIMA: (0 1 1)(0 1 1)  Obs.: 144  Transform: log
#> AICc: 983.9, BIC:  1008  QS (no seasonality in final):    0  
#> Box-Ljung (no autocorr.): 29.43   Shapiro (normality): 0.9781 *

^{Created on 2019-12-28 by the reprex package (v0.3.0)}

christophsax / seasonal