MC模拟直方图
Histogram of MC simulation
我正在尝试比较马尔可夫链 (MC) 模拟和实际数据的直方图。我尝试使用下面的代码 运行 进行模拟,但我并不完全理解。 R 似乎已经接受了代码,但我不知道如何 运行 直方图...作为背景,数据是美国经济的扩张和收缩(可在此处找到:http://www.nber.org/cycles.html)。我已将这两个状态之间的转换矩阵设置为 P
,其中列总和为 1,并且状态之间的转换计算为 "transitions / months in each state"。我认为 n
对应于此处的转换,但我可能是错的...
P <- matrix(c(0.74961, 0.57291, 0.25039, 0.42709),2,2)
P <- t(P)
colSums(P)
n <- 33
MC.sim <- function(n,P) {
sim<-c()
m <- ncol(P)
sim[1] <- sample(1:m,1)
for(i in 2:n){
newstate <- sample(1:m,1,prob=P[,sim[i-1]])
sim[i] <- newstate
}
sim
}
如评论中所述,您可以使用 hist()
:
P <- matrix(c(0.74961, 0.57291, 0.25039, 0.42709),2,2)
P <- t(P)
colSums(P)
n <- 33
MC.sim <- function(n,P) {
sim<-c()
m <- ncol(P)
sim[1] <- sample(1:m,1)
for(i in 2:n){
newstate <- sample(1:m,1,prob=P[,sim[i-1]])
sim[i] <- newstate
}
return(sim)
}
# Save results of simulation
temp <- MC.sim(n,P)
#plot basic histogram
hist(temp)
我没有检查你的模拟。但它只生成 1 或 2 的值。
您也可以使用 ggplot()
:
# Save results of simulation
temp <- MC.sim(n,P)
# Make data frame for ggplot
t <- as.data.frame(temp)
names(t) <- "Sim"
p <- ggplot(t, aes(x = Sim))
p + geom_histogram(binwidth = 0.5)
我正在尝试比较马尔可夫链 (MC) 模拟和实际数据的直方图。我尝试使用下面的代码 运行 进行模拟,但我并不完全理解。 R 似乎已经接受了代码,但我不知道如何 运行 直方图...作为背景,数据是美国经济的扩张和收缩(可在此处找到:http://www.nber.org/cycles.html)。我已将这两个状态之间的转换矩阵设置为 P
,其中列总和为 1,并且状态之间的转换计算为 "transitions / months in each state"。我认为 n
对应于此处的转换,但我可能是错的...
P <- matrix(c(0.74961, 0.57291, 0.25039, 0.42709),2,2)
P <- t(P)
colSums(P)
n <- 33
MC.sim <- function(n,P) {
sim<-c()
m <- ncol(P)
sim[1] <- sample(1:m,1)
for(i in 2:n){
newstate <- sample(1:m,1,prob=P[,sim[i-1]])
sim[i] <- newstate
}
sim
}
如评论中所述,您可以使用 hist()
:
P <- matrix(c(0.74961, 0.57291, 0.25039, 0.42709),2,2)
P <- t(P)
colSums(P)
n <- 33
MC.sim <- function(n,P) {
sim<-c()
m <- ncol(P)
sim[1] <- sample(1:m,1)
for(i in 2:n){
newstate <- sample(1:m,1,prob=P[,sim[i-1]])
sim[i] <- newstate
}
return(sim)
}
# Save results of simulation
temp <- MC.sim(n,P)
#plot basic histogram
hist(temp)
我没有检查你的模拟。但它只生成 1 或 2 的值。
您也可以使用 ggplot()
:
# Save results of simulation
temp <- MC.sim(n,P)
# Make data frame for ggplot
t <- as.data.frame(temp)
names(t) <- "Sim"
p <- ggplot(t, aes(x = Sim))
p + geom_histogram(binwidth = 0.5)