require(mvtnorm)
mu <- c(1,2)
cov_mat <- matrix(c(1, 0.9, 0.9, 1), 2, 2, TRUE)
set.seed(1)
sample_chol <- rmvnorm(2000, mu, cov_mat)
x1<-sample_chol[,1]
y1<-sample_chol[,2]
sample_val <- function(x) {
delta_1 <- 0.6
delta_2 <- 0.4
return(2mu - x + c(runif(1, -delta_1, delta_1), runif(1, -delta_2, delta_2)))
}
results_4 <- matrix(NA, nrow = N, ncol = 2)
results_4[1,] <- x_start
i <- 2
for( i in 2:N ) {
new_val <- sample_val(results_4[i - 1, ])
alpha_x_y <- alpha(results_4[i - 1], new_val)
if(runif(1) < alpha_x_y){
results_4[i, ] <- new_val
}else{
results_4[i,] <- results_4[i-1, ]
}
}
results_4 <- results_4[round(0.1N):N,]
x<-results_4[,1]
y<-results_4[,2]
plot(
x,
y,
type=“p”,#线型
main=“cholesky”,W
sub=“MH???Random Walk”,
xlab=“x”,
ylab=“y”,
ylim=c(-6,6),
xlim=c(-6,6),
cex=0.3,#符号的大小
col="#cccccc",#必须要写全,六位的颜色标记;
pch=19,#设置标点的样式,19号表示实心;
asp=0.5)
#y/x的比例，y轴数值长度与x轴数值长度的比值
legend(“topleft”,inset=.02,c(“cholesky”,“MH-Pseudo dominating”),pch=c(20,20),col=c(“red”,“grey”))
points(x1,y1,pch=19, cex = 0.5, col= “red”)
输出结果: