[R] X-12-ARIMA
Alistair Gray
alistair at statsresearch.co.nz
Thu Jun 10 04:37:50 CEST 2004
Dear All,
I've used the X-12-ARIMA or its earlier versions from S+ and R under both Unix
and Windows platforms for many years using the klugey approach of calling an
executable using in R the system function. I've found this serviceable for the
following reasons.
1) Paul Gilbert's hunch is correct that many of the subroutines have extensive
IO calls (especially the X-11 engine) and so it is not straightfoward to call
using .Fortran.
2) X-12-ARIMA has in the Unix final version 0.2.10 which I use 15 different
commands with their own parameters covering about six pages of description so
it's not clear to me having a long R function parameter list is a great advantage.
My experience with people who have ported X-12-ARIMA into other econometrics
software is that they port a limited range of the commands and options. (I'd be
interested to know from Dirk Eddelbuettel if this is the case for GNU-Gretl.)
3) X-12_ARIMA through its specification files (particularly the metafiles) is
set up to handle multiple runs on different time series with different parameter
settings. I'm not sure I'd want to re-invent it. Secondly I find it an OK way to
keep track of what I've done.
However, I'd be a grateful user of a less klugey approach.
Below is a clearer explanation of my klugey approach and functions which I've
used successfully. Feel free to use them.
# create data file for the series to be read by x12 fortran pgm
# blp is a object of class ts written to file blp.dat in the x12 datevalue
# format
writex12in(blp,"blp")
# create the specification file using your favourite editor say blp.spc to be
# read by x12 fortran pgm
# a vanilla specification file might have the following
# series{
# title="building consents"
# start=1973.01
# span=(1973.01, 2000.12)
# period=12
# file="blp.dat"
# format="datevalue"
# }
# x11{
# mode=mult
# sigmalim=(1.8 2.8)
# seasonalma=x11default
# trendma=13
# appendfcst=no
# save=(b1 c17 d10 d11 d12 d13)
# savelog=(m1 m2 m3 m4 m5 m6 m7 m8 m9 m10 m11 q q2 msr icr fb1 fd8 msf ids)
# }
# execute the x12 fortran pgm. executable x12a stored in /home/fred/x12a
system("/home/fred/x12a/x12a blp")
# read x12 fortran pgm output tables back into R
blp.x12 <- readx12out("blp", adtype="M", calendar=F, tblnames=NULL)
# basic versions of writex12in and readx12out
writex12in <- function(tso,file){
write.table(cbind(time(tso) %/% 1, cycle(tso), tso),
file=paste(file,".dat",sep=""),
sep=" ", quote=F, row.names=F, col.names=F)
}
readx12out <- function(file, adtype = "M", calendar = F, tblnames = NULL) {
notbls <- 6 + (calendar != F) + length(tblnames)
comp <- vector("list", notbls)
if(calendar == F) {
names(comp) <- c("original", "seasonal", "adjusted", "trend",
"irregular", "weights", tblnames)
tblnames <- c("b1", "d10", "d11", "d12", "d13", "c17", tblnames)
}
else {
names(comp) <- c("original", "seasonal", "adjusted", "trend",
"irregular", "calendar", "weights", tblnames)
if(calendar == "C")
tblnames <- c("b1", "d10", "d11", "d12", "d13", "xca", "c17",
tblnames)
else tblnames <- c("b1", "d10", "d11", "d12", "d13", "c16", "c17",
tblnames)
}
for(i in seq(along = comp)) {
series <- read.table(paste(file, ".", tblnames[i], sep = ""),
as.is = T, skip = 2)
times <- series[,1]
begin <- c(times[1]%/%100,times[1]%%100)
freq <- max(times%%100)
comp[[i]] <- ts(data = series[,2], start = begin, frequency = freq)
}
# for the log transform the adjusted series is the log of the original
# also x12 transforms trend, seasonal, etc back into original scale
# i.e. making a multiplicative decomposition in which case 1.0 needs to be
# subtracted from seasonal, irregular and calendar to be consistent with
# multiplicative. Also x12 makes correction for bias to trend in original scale
# so best to treat as multiplicative decomposition rather than decomposition in
# transformed scale
if(adtype == "L") {
comp$transformed <- comp$original
comp$seasonal <- comp$seasonal - 1
comp$irregular <- comp$irregular - 1
if(calendar != F)
comp$calendar <- comp$calendar - 1
comp$power <- "log additive presented as multiplicative. 1.0 has been
subtracted from seasonal, irregular & calendar (if present)"
}
else if(adtype == "M") {
comp$transformed <- comp$original
comp$seasonal <- comp$seasonal - 1
comp$irregular <- comp$irregular - 1
if(calendar != F)
comp$calendar <- comp$calendar - 1
comp$power <- "multiplicative. 1.0 has been subtracted from seasonal,
irregular & calendar (if present)"
}
else if(adtype == "A") {
comp$original <- rtso
comp$power <- "additive"
}
return(comp)
}
--
Alistair Gray Email: alistair at statsresearch.co.nz
Statistics Research Associates Ltd Web: www.statsresearch.co.nz
PO Box 12 649, Thorndon, Wellington Phone: +64 +4 972 6531
NEW ZEALAND Mobile: +64 +21 610 569
More information about the R-help
mailing list