時間序列(數值模擬)

HW1

To generate the empirical distribution for the sum of the roll of two dice and compare it with the theoretical one.

我們可以發現實際值與理論值會有誤差，但擲骰500次比擲骰100次更接近理論值。

set.seed(1)
hw1A<-sample(x = 1:6,size = 100,replace = T,
             prob = c(1/6,1/6,1/6,1/6,1/6,1/6))
hw1B<-sample(x = 1:6,size = 100,replace = T,
             prob = c(1/6,1/6,1/6,1/6,1/6,1/6))
hw1C<-hw1A+hw1B
hw1<-as.data.frame(table(hw1C))
hw1$pro<-hw1$Freq/100
hw1$ther<-c(1,2,3,4,5,6,5,4,3,2,1)/36
names(hw1)<-c("加總值","次數","機率","理論值")

hw1A.2<-sample(x = 1:6,size = 500,replace = T,
               prob = c(1/6,1/6,1/6,1/6,1/6,1/6))
hw1B.2<-sample(x = 1:6,size = 500,replace = T,
               prob = c(1/6,1/6,1/6,1/6,1/6,1/6))
hw1C.2<-hw1A.2+hw1B.2
hw1.2<-as.data.frame(table(hw1C.2))
hw1.2$pro<-hw1.2$Freq/500
hw1.2$ther<-c(1,2,3,4,5,6,5,4,3,2,1)/36
names(hw1.2)<-c("加總值","次數","機率","理論值")
hw1

##    加總值 次數 機率     理論值
## 1       2    6 0.06 0.02777778
## 2       3    2 0.02 0.05555556
## 3       4    6 0.06 0.08333333
## 4       5    7 0.07 0.11111111
## 5       6   10 0.10 0.13888889
## 6       7   19 0.19 0.16666667
## 7       8   19 0.19 0.13888889
## 8       9   15 0.15 0.11111111
## 9      10    6 0.06 0.08333333
## 10     11    7 0.07 0.05555556
## 11     12    3 0.03 0.02777778

hw1.2

##    加總值 次數  機率     理論值
## 1       2   13 0.026 0.02777778
## 2       3   26 0.052 0.05555556
## 3       4   50 0.100 0.08333333
## 4       5   56 0.112 0.11111111
## 5       6   81 0.162 0.13888889
## 6       7   82 0.164 0.16666667
## 7       8   62 0.124 0.13888889
## 8       9   48 0.096 0.11111111
## 9      10   43 0.086 0.08333333
## 10     11   16 0.032 0.05555556
## 11     12   23 0.046 0.02777778

HW2

To generate the DF distribution and find out the critical values (at 1%,5% and 10% signicance levels) for the t-statistics for the null hypothesis = 0.

下圖為透過模擬Dickey-Fuller Distribution。以及critical values (at 1%,5% and 10% signicance levels)

hw2<-data.frame(1:10000)
names(hw2)<-"t.value"
for(j in 1:10000){
x<-rnorm(n = 50)
y0<-0
for(i in 1:50){
y0<-y0+x[i]
}

y<-data.frame(1:100)
names(y)<-"y"
e<-rnorm(n = 100)
for(i in 1:100){
  y0<-y0+e[i]
  y[i,1]<-y0
}
y$y.diff<-c(NA,diff(y$y,differences = 1))
df.t<-lm(y$y.diff[2:100]~y$y[1:99])
summary(df.t)
tstats <- coef(df.t) / sqrt(diag(vcov(df.t)))
hw2[j,1]<-tstats[2]
}
library(ggplot2)
ggplot(data = hw2)+
  geom_density(aes(x = hw2$t.value))

quantile(hw2$t.value,probs = 0.01)

##        1% 
## -3.509088

quantile(hw2$t.value,probs = 0.05)

##        5% 
## -2.921439

quantile(hw2$t.value,probs = 0.10)

##       10% 
## -2.598103