Hi all, I'm having trouble calculating the likelihood of a time series with AR(1) errors. I am generating my covariance matrix according to page 2 of ( http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-timeseries-regression.pdf), using the library mvtnorm and the multivariate normal density function dmvnorm(). Here's some example code:
library(mvtnorm) # Generate a basic time series with AR(1) Errors: t <- 1:100 error <- as.numeric(arima.sim(n = length(t), list(ar = c(0.8897)), sd = 10)) series <- 5*t + error # Fit the series using a basic linear model assuming errors are IID Normal naive.model <- lm(series ~ t -1) # Examine and model the residuals residuals <- series - t*coef(naive.model) residual.model <- arima(residuals, c(1,0,0), include.mean=F) # Construct the covariance matrix, assuming the process variance (10^2) is known sigma <- diag(length(t)) sigma[(abs(row(sigma)-col(sigma)) == 1)] = as.numeric(coef(residual.model)) sigma <- sigma*10^2 # Calculate the MVN density... dmvnorm(series, t*coef(naive.model) ,sigma, log=T) Without fail, I get the following error message: Warning message: In log(eigen(sigma, symmetric = TRUE, only.values = TRUE)$values) : NaNs produced. It's worth noting that the matrix from the following ( https://stat.ethz.ch/pipermail/r-help/2007-May/131728.html) "works", but I think is actually for an MA(1) process rather than an AR(1) process. I gather the message means the proposed covariance matrix may not be invertible. This said I'm stuck on how to proceed and would be extremely appreciative of any thoughts. Thank you very much, Robert -- Robert Kindman Harvard College, Class of 2014 rkind...@college.harvard.edu +1.919.599.1921 [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
