year = seq(1961 + 1/12, 1972 + 10/12, by = 1/12)
par(mar = c(4, 4, 2, 2) + 0.1, mfrow = c(2, 2))
for(ii in 1:4) {
plot(year, usgrain[, ii], main = names(usgrain)[ii],
type = "l", xlab = "", ylab = "")
points(year, usgrain[,ii], pch = "*")
}
apply(usgrain, 2, mean) # mu vector
cgrain = scale(usgrain, scale = FALSE) # Center the time series only
fit = vglm(cgrain ~ 1, rrar(Ranks = c(4, 1)), trace = TRUE)
summary(fit)
print(fit@misc$Ak1, dig = 2)
print(fit@misc$Cmatrices, dig = 3)
print(fit@misc$Dmatrices, dig = 3)
print(fit@misc$omegahat, dig = 3)
print(fit@misc$Phimatrices, dig = 2)
par(mar = c(4, 4, 2, 2) + 0.1, mfrow = c(4, 1))
for(ii in 1:4) {
plot(year, fit@misc$Z[,ii], main = paste("Z", ii, sep = ""),
type = "l", xlab = "", ylab = "")
points(year, fit@misc$Z[,ii], pch = "*")
}
Run the code above in your browser using DataLab