On Tue, Feb 3, 2009 at 12:58 PM, Harsh <singhalblr at gmail.com>
wrote:> Hi,
> I am newbie user of ggplot and would like some assistance in
> implementing time series plots.
>
> I'd like to know how the tsdiag plot can be made in ggplot?
It's not particularly easy because the code for tsdiag interweaves
computing statistics and then plotting them. However, with a bit of
effort you can disentangle them and create the equivalent plots with
ggplot2.
# What does tsdiag look like?
fit <- arima(lh, c(1,0,0))
tsdiag(fit)
# Need to extract the calculations from the existing
# plotting code. I started by looking at the contents of
# stats:::tsdiag.Arima
# Standardised residuals
stdres <- as.numeric(resid(fit) / sqrt(fit$sigma2))
qplot(seq_along(stdres), stdres, yend = 0, xend = seq_along(stdres),
geom="segment") +
geom_hline(colour = "grey50")
# ACF
acf <- acf(resid(fit), plot = F)
acf_df <- data.frame(
acf = acf$acf,
lag = acf$lag
)
qplot(lag, acf, data = acf_df, yend = 0, xend = lag, geom="segment") +
geom_hline() + geom_point()
# To figure out the confidence bands, we need to dig into the code
# for stats:::plot.acf. Unfortunately this isn't available in its own
# method, so I'm not entirely sure I've extracted it correctly
ci <- 0.95
clim <- qnorm((1 + ci)/2) / sqrt(acf$n.used)
qplot(lag, acf, data = acf_df, yend = 0, xend = lag, geom="segment") +
geom_hline(colour = "grey50") +
geom_point() +
geom_hline(yintercept = c(-clim, clim), colour = "darkblue")
# Finally the code for the Ljung-Box statistic
ljungbox <- function(i, rs) {
Box.test(rs, i, type = "Ljung-Box")$p.value
}
acf_df$pvalue <- sapply(acf_df$lag, ljungbox, series = resid(fit))
qplot(lag, pvalue, data = acf_df) +
geom_hline(colour = "grey50") +
geom_hline(yintercept = 0.05, colour = "darkblue")
Regards,
Hadley
--
http://had.co.nz/