VARIMA ignores specials

[robjhyndman]

library(fable)

Dados <- data.frame(
  time = seq(240),
  Admissões = c(
    230L, 213L, 235L, 222L, 310L, 234L, 238L, 333L, 240L, 242L, 197L, 217L,
    104L, 213L, 255L, 220L, 251L, 233L, 181L, 168L, 129L, 188L, 149L,
    186L, 112L, 163L, 220L, 199L, 186L, 201L, 163L, 229L, 222L, 192L,
    172L, 198L, 109L, 120L, 190L, 149L, 181L, 168L, 182L, 168L, 227L,
    202L, 193L, 226L, 165L, 190L, 168L, 126L, 163L, 155L, 159L, 176L,
    159L, 153L, 198L, 189L, 86L, 148L, 123L, 160L, 156L, 121L, 209L,
    149L, 155L, 142L, 149L, 246L, 140L, 131L, 148L, 157L, 182L, 137L,
    128L, 137L, 164L, 105L, 155L, 125L, 220L, 263L, 248L, 289L, 454L,
    275L, 289L, 317L, 347L, 359L, 410L, 329L, 206L, 212L, 278L, 262L,
    337L, 283L, 303L, 268L, 285L, 353L, 271L, 327L, 170L, 217L, 254L,
    262L, 350L, 270L, 303L, 350L, 322L, 305L, 264L, 278L, 178L, 257L,
    262L, 215L, 265L, 321L, 285L, 259L, 215L, 247L, 277L, 245L, 123L,
    198L, 217L, 296L, 299L, 342L, 316L, 204L, 234L, 270L, 194L, 233L,
    171L, 217L, 255L, 300L, 293L, 242L, 257L, 307L, 354L, 243L, 259L,
    263L, 153L, 215L, 286L, 276L, 327L, 322L, 262L, 278L, 285L, 269L,
    279L, 358L, 180L, 301L, 505L, 317L, 537L, 441L, 464L, 368L, 366L,
    388L, 248L, 279L, 173L, 246L, 230L, 203L, 297L, 200L, 251L, 267L,
    207L, 217L, 239L, 302L, 394L, 381L, 257L, 70L, 140L, 186L, 189L,
    221L, 299L, 460L, 328L, 298L, 397L, 441L, 353L, 321L, 455L, 450L,
    419L, 444L, 494L, 526L, 628L, 476L, 462L, 611L, 562L, 504L, 588L,
    513L, 440L, 557L, 487L, 475L, 407L, 296L, 407L, 456L, 432L, 345L,
    479L, 386L, 354L, 411L, 467L, 397L, 430L, 280L
  ),
  Demissões = c(
    188L, 218L, 210L, 207L, 219L, 191L, 285L, 197L, 222L, 232L, 218L, 180L,
    177L, 219L, 204L, 178L, 169L, 120L, 145L, 150L, 114L, 204L, 178L,
    208L, 188L, 164L, 195L, 127L, 127L, 185L, 144L, 142L, 135L, 159L,
    154L, 195L, 116L, 159L, 142L, 152L, 144L, 165L, 207L, 201L, 170L,
    169L, 163L, 175L, 148L, 134L, 128L, 138L, 139L, 100L, 150L, 141L,
    121L, 151L, 115L, 155L, 175L, 149L, 135L, 135L, 126L, 127L, 132L,
    140L, 131L, 246L, 133L, 148L, 151L, 99L, 107L, 166L, 94L, 112L,
    93L, 71L, 131L, 103L, 132L, 101L, 300L, 348L, 374L, 322L, 308L,
    300L, 303L, 317L, 282L, 311L, 298L, 269L, 236L, 286L, 250L, 232L,
    245L, 248L, 276L, 258L, 283L, 316L, 279L, 277L, 266L, 260L, 270L,
    274L, 257L, 240L, 288L, 269L, 220L, 283L, 242L, 278L, 212L, 271L,
    247L, 217L, 245L, 269L, 277L, 197L, 220L, 263L, 303L, 226L, 286L,
    227L, 242L, 279L, 317L, 250L, 294L, 225L, 249L, 258L, 152L, 199L,
    322L, 237L, 287L, 238L, 239L, 293L, 219L, 274L, 314L, 288L, 256L,
    303L, 261L, 226L, 267L, 244L, 262L, 265L, 277L, 263L, 252L, 249L,
    241L, 420L, 306L, 696L, 355L, 349L, 404L, 307L, 357L, 320L, 271L,
    251L, 206L, 273L, 289L, 214L, 203L, 224L, 269L, 202L, 248L, 217L,
    215L, 332L, 188L, 215L, 256L, 271L, 318L, 226L, 204L, 180L, 151L,
    217L, 223L, 275L, 241L, 266L, 303L, 266L, 346L, 278L, 289L, 283L,
    319L, 385L, 379L, 408L, 382L, 468L, 501L, 530L, 516L, 408L, 431L,
    644L, 411L, 545L, 430L, 426L, 373L, 388L, 394L, 446L, 432L, 376L,
    550L, 498L, 398L, 484L, 407L, 414L, 682L, 399L
  )
) |>
  tsibble::as_tsibble(index = time)

Dados |>
  model(varima = VARIMA(vars(Admissões, Demissões) ~ pdq(0, 1, 0)))
#> # A mable: 1 x 1
#>                    varima
#>                   <model>
#> 1 <VARIMA(1,1,1) w/ mean>

^{Created on 2024-11-28 with reprex v2.1.1}

[mitchelloharawild]

A few issues,

1. The default identification method is/was "kronecker_indices", which passes on automatic selection to MTS::Kronid(). I've now changed the default to be "none" when p and q are fully specified, although this isn't ideal either for identification issues within the given lag orders. @GeorgeAthana: Is there an approach to constraining/identifying a VARMA model for given lag orders p and q?
2. The MTS package (which is our current implementation) doesn't play well with 0 lag orders. For instance, when finding kronecker indices if the maximum search order plag < 1 it instead uses floor(log(T)) + 1 (where T is the series length).
   https://github.com/d-/MTS/blob/5a274c3f9bc5d290be02e85c971ac72ce2f1dc88/R/MTS.R#L4629-L4631
   Similarly even when using no identification in this case (the updated default), it estimates VARMA(1,0) when p=0 and q=0 is specified
   https://github.com/d-/MTS/blob/5a274c3f9bc5d290be02e85c971ac72ce2f1dc88/R/MTS.R#L2047

Changes / Improvements

• 5337ff7 Changed the default identification method to "none" when p and q are fully specified.
• 8d67360 I've patched an issue with using integers instead of doubles, again this is problem arises from the MTS package where it assumes all numerics are integers in the cpp implementation.
• aacfb69 Added detailed documentation about identification options and defaults

[robjhyndman]

It sounds like we need our own implementation.

[mitchelloharawild]

It sounds like we need our own implementation.

Yes, that would be a nice long term solution.