'homework 2-1a'
## [1] "homework 2-1a"
library(fUnitRoots)
## Loading required package: timeDate
## Loading required package: timeSeries
## Loading required package: fBasics
da=read.table("m-unrate-4811.txt",header=T)
values = seq(from = as.Date("1948-01-01"), to = as.Date("2011-11-03"), by = 'month')
rownames(da) <- values
str(da)
## 'data.frame': 767 obs. of 1 variable:
## $ rate: num 3.4 3.8 4 3.9 3.5 3.6 3.6 3.9 3.8 3.7 ...
gdp=log(da[,1])
m1=ar(diff(gdp),method='mle')
m1$order
## [1] 12
adfTest(gdp,lags=10,type=c("c"))
##
## Title:
## Augmented Dickey-Fuller Test
##
## Test Results:
## PARAMETER:
## Lag Order: 10
## STATISTIC:
## Dickey-Fuller: -3.3319
## P VALUE:
## 0.01523
##
## Description:
## Thu Oct 17 18:17:10 2019 by user: grimr
"WHY? with Parameter 10 ADF-test statistic is -3.3319 with p-value 0.01523 Thus indicating that Unit-root hypothesis cannot be rejected"
## [1] "WHY? with Parameter 10 ADF-test statistic is -3.3319 with p-value 0.01523 Thus indicating that Unit-root hypothesis cannot be rejected"
"2-1B"
## [1] "2-1B"
eps=log(da$rate)
koeps=ts(eps,frequency=4,start=c(1983,1))
plot(koeps,type='l')
points(koeps,pch=0,cex=0.6)
par(mfcol=c(2,2))
koeps=log(da$rate)
deps=diff(koeps)
sdeps=diff(koeps,4)
ddeps=diff(sdeps)
acf(koeps,lag=20)
acf(deps,lag=20)
acf(sdeps,lag=20)
acf(ddeps,lag=20)
# Obtain time plots
c1=c("2","3","4","1")
c2=c("1","2","3","4")
par(mfcol=c(3,1))
plot(deps,xlab='year',ylab='diff',type='l')
points(deps,pch=c1,cex=0.7)
plot(sdeps,xlab='year',ylab='sea-diff',type='l')
points(sdeps,pch=c2,cex=0.7)
plot(ddeps,xlab='year',ylab='dd',type='l')
points(ddeps,pch=c1,cex=0.7)
# Estimation
m1=arima(koeps,order=c(0,1,1),seasonal=list(order=c(0,1,1),period=4))
m1
##
## Call:
## arima(x = koeps, order = c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 4))
##
## Coefficients:
## ma1 sma1
## 0.0955 -1.0000
## s.e. 0.0306 0.0074
##
## sigma^2 estimated as 0.001514: log likelihood = 1382.06, aic = -2758.12
tsdiag(m1,gof=20) # model checking
Box.test(m1$residuals,lag=12,type='Ljung')
##
## Box-Ljung test
##
## data: m1$residuals
## X-squared = 143.2, df = 12, p-value < 2.2e-16
pp=1-pchisq(13.30,10)
pp
## [1] 0.2073788
koeps=log(da$rate)
length(koeps)
## [1] 767
y=koeps[1:100]
m1=arima(y,order=c(0,1,1),seasonal=list(order=c(0,1,1),period=4))
m1
##
## Call:
## arima(x = y, order = c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 4))
##
## Coefficients:
## ma1 sma1
## 0.1522 -0.9996
## s.e. 0.0923 0.2234
##
## sigma^2 estimated as 0.004837: log likelihood = 112.04, aic = -218.08
pm1=predict(m1,7)
names(pm1)
## [1] "pred" "se"
pred=pm1$pred
se=pm1$se
ko=da$rate
fore=exp(pred+se^2/2)
v1=exp(2*pred+se^2)*(exp(se^2)-1)
s1=sqrt(v1)
eps=ko[80:107]
length(eps)
## [1] 28
tdx=(c(1:28)+3)/4+2005
upp=c(ko[100],fore+2*s1)
low=c(ko[100],fore-2*s1)
min(low,eps)
## [1] 2.291424
max(upp,eps)
## [1] 6.1
plot(tdx,eps,xlab='year',ylab='earnings',type='l')
points(tdx[22:28],fore,pch='*')
lines(tdx[21:28],upp,lty=2)
lines(tdx[21:28],low,lty=2)
points(tdx[22:28],ko[101:107],pch='o',cex=0.7)
"2-7"
## [1] "2-7"
data=read.table("q-jnj-earns-9211.txt",header=T)
data$date <- as.Date(with(data, paste(year, mon, day,sep="-")), "%Y-%m-%d")
data$day <- NULL
data$mon <- NULL
data$year <- NULL
rownames(data) <- data$date
data$date <- NULL
plot.ts( data, xlab='year', ylab='earnings')
