Quantcast

Autocorrelation in R

classic Classic list List threaded Threaded
2 messages Options
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Autocorrelation in R

Iuri Gavronski
Hi,

I am trying to learn time series, and I am attending a colleague's
course on Econometrics. However, he uses e-views, and I use R. I am
trying to reproduce his examples in R, but I am having problems
specifying a AR(1) model. Would anyone help me with my code?

Thanks in advance!

Reproducible code follows:

download.file("https://sites.google.com/a/proxima.adm.br/main/ex_32.csv
--no-check-certificate", "ex_32.csv", method="wget")

ex32=read.csv("ex_32.csv")

lm_ex32=lm(gc ~ yd, data=ex32)

summary(lm_ex32)

# Durbin-Watson (slide 26)
library(lmtest)

dwtest(gc ~ yd, data=ex32)
# or
dwtest(lm_ex32)

# Breusch-Godfrey
bgtest(lm_ex32, order=2)

# AR(1)

# In e-views, the specification was:
# GC = YD AR(1)
# and the output was:

# Dependent Variable: GC
# Method: Least Squares
# Sample: 1970Q2 1995Q2
# Included observations: 101
# Convergence achieved after 6 interations
# =========================================================
# Variable    Coefficient   Std.Error t-Statistic Prob.
# C           -56.99706     19.84692  -2.871835   0.0050
# YD          0.937035      0.006520  143.7170    0.0000
# AR(1)       0.752407      0.066565  11.30338    0.0000
# =========================================================
# R-squared 0.999691 Mean dependent var 2345.867
# Adjusted R-squared 0.999685 S.D. dependent var 1284.675
# S.E. of regression 22.81029 Akaike info criterion 9.121554
# Sum squared resid 50990.32 Schwarz criterion 9.199231
# Log likelihood -457.6385 F-statistic 158548.1
# Durbin-Watson stat 2.350440 Prob(F-statistic) 0.000000

# following code based on
http://www.stat.pitt.edu/stoffer/tsa2/R_time_series_quick_fix.htm
# "And now for some regression with autocorrelated errors."

# I've tried to follow the example in Pinheiro & Bates (2004), p.
239-244, with no success.

gc_ts = ts(ex32[66:166,"gc"])
yd_ts = ts(ex32[66:166,"yd"])

library(nlme)
trend = time(gc_ts)

fit_lm = lm(gc_ts ~ trend + yd_ts)
acf(resid(fit_lm))
pacf(resid(fit_lm))



gls_ex32_ar1 = gls(gc_ts ~ trend + yd_ts, correlation = corAR1(form=
~yd_ts),method="ML")
summary(gls_ex32_ar1)


_____________________________________________
Dr. Iuri Gavronski
Assistant Professor
Programa de Pós-Graduação em Administração
Universidade do Vale do Rio dos Sinos – UNISINOS
Av. Unisinos, 950 – São Leopoldo – RS – Brasil
Sala (Room) 5A 406 D
93022-000
www.unisinos.br

TEL +55-51-3591-1122 ext. 1589
FAX +55-51-3590-8447
Email: [hidden email]

CV Lattes: http://lattes.cnpq.br/8843390959025944

______________________________________________
[hidden email] 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.
Reply | Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: Autocorrelation in R

Achim Zeileis-4
On Wed, 8 Jun 2011, Iuri Gavronski wrote:

> Hi,
>
> I am trying to learn time series, and I am attending a colleague's
> course on Econometrics. However, he uses e-views, and I use R. I am
> trying to reproduce his examples in R, but I am having problems
> specifying a AR(1) model. Would anyone help me with my code?
>
> Thanks in advance!
>
> Reproducible code follows:
>
> download.file("https://sites.google.com/a/proxima.adm.br/main/ex_32.csv
> --no-check-certificate", "ex_32.csv", method="wget")
>
> ex32=read.csv("ex_32.csv")
>
> lm_ex32=lm(gc ~ yd, data=ex32)
>
> summary(lm_ex32)
>
> # Durbin-Watson (slide 26)
> library(lmtest)
>
> dwtest(gc ~ yd, data=ex32)
> # or
> dwtest(lm_ex32)
>
> # Breusch-Godfrey
> bgtest(lm_ex32, order=2)
>
> # AR(1)
>
> # In e-views, the specification was:
> # GC = YD AR(1)
> # and the output was:
>
> # Dependent Variable: GC
> # Method: Least Squares
> # Sample: 1970Q2 1995Q2
> # Included observations: 101
> # Convergence achieved after 6 interations
> # =========================================================
> # Variable    Coefficient   Std.Error t-Statistic Prob.
> # C           -56.99706     19.84692  -2.871835   0.0050
> # YD          0.937035      0.006520  143.7170    0.0000
> # AR(1)       0.752407      0.066565  11.30338    0.0000
> # =========================================================
> # R-squared 0.999691 Mean dependent var 2345.867
> # Adjusted R-squared 0.999685 S.D. dependent var 1284.675
> # S.E. of regression 22.81029 Akaike info criterion 9.121554
> # Sum squared resid 50990.32 Schwarz criterion 9.199231
> # Log likelihood -457.6385 F-statistic 158548.1
> # Durbin-Watson stat 2.350440 Prob(F-statistic) 0.000000
I'm not sure what exactly E-Views does here, but an ARIMAX(1,0,0) model
estimated by least squares seems to come rather close.

## create a time series object of your data
ex32ts <- ts(ex32[,-1], start = c(1954, 1), freq = 4)

## select relevant subset
ex32ts1 <- window(ex32ts, start = c(1970, 2))

## fit ARIMAX(1,0,0) model
m <- arima(ex32ts1[,"gc"], order = c(1, 0, 0),
   xreg = ex32ts1[,"yd"], method = "CSS")

## print output, coefficient tests, etc.
m
coeftest(m)
logLik(m)

It seems to be slightly different but that can well be due to different
fitting algorithms...

hth,
Z

> # following code based on
> http://www.stat.pitt.edu/stoffer/tsa2/R_time_series_quick_fix.htm
> # "And now for some regression with autocorrelated errors."
>
> # I've tried to follow the example in Pinheiro & Bates (2004), p.
> 239-244, with no success.
>
> gc_ts = ts(ex32[66:166,"gc"])
> yd_ts = ts(ex32[66:166,"yd"])
>
> library(nlme)
> trend = time(gc_ts)
>
> fit_lm = lm(gc_ts ~ trend + yd_ts)
> acf(resid(fit_lm))
> pacf(resid(fit_lm))
>
>
>
> gls_ex32_ar1 = gls(gc_ts ~ trend + yd_ts, correlation = corAR1(form=
> ~yd_ts),method="ML")
> summary(gls_ex32_ar1)
>
>
> _____________________________________________
> Dr. Iuri Gavronski
> Assistant Professor
> Programa de Pós-Graduação em Administração
> Universidade do Vale do Rio dos Sinos ? UNISINOS
> Av. Unisinos, 950 ? São Leopoldo ? RS ? Brasil
> Sala (Room) 5A 406 D
> 93022-000
> www.unisinos.br
>
> TEL +55-51-3591-1122 ext. 1589
> FAX +55-51-3590-8447
> Email: [hidden email]
>
> CV Lattes: http://lattes.cnpq.br/8843390959025944
>
> ______________________________________________
> [hidden email] 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.
______________________________________________
[hidden email] 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.
Loading...