if (FALSE) {
set.seed(1977)
exam2 <- function(n, mu, sigma, rho, a, b, trunc)
{
Sigma <- matrix(c(sigma[1]^2, rho*sigma[1]*sigma[2], rho*sigma[1]*sigma[2], sigma[2]^2), nrow=2)
x <- rTMVN(n, mean=mu, Sigma=Sigma, a=a, b=b, trunc=trunc)
x1.gr <- seq(mu[1]-3.5*sigma[1], mu[1]+3.5*sigma[1], length=100)
x2.gr <- seq(mu[2]-3.5*sigma[2], mu[2]+3.5*sigma[2], length=100)
z <- cbind(rep(x1.gr, 100), rep(x2.gr, each=100))
dens.z <- matrix(dMVN(z, mean=mu, Sigma=Sigma), ncol=100)
MEAN <- round(apply(x, 2, mean), 3)
SIGMA <- var(x)
SD <- sqrt(diag(SIGMA))
RHO <- round(SIGMA[1,2]/(SD[1]*SD[2]), 3)
SD <- round(SD, 3)
layout(matrix(c(0,1,1,0, 2,2,3,3), nrow=2, byrow=TRUE))
contour(x1.gr, x2.gr, dens.z, col="darkblue", xlab="x[1]", ylab="x[2]")
points(x[,1], x[,2], col="red")
title(sub=paste("Sample mean = ", MEAN[1], ", ", MEAN[2], ", Sample SD = ", SD[1], ", ", SD[2],
", Sample rho = ", RHO, sep=""))
plot(1:n, x[,1], type="l", xlab="Iteration", ylab="x[1]", col="darkgreen")
plot(1:n, x[,2], type="l", xlab="Iteration", ylab="x[2]", col="darkgreen")
return(x)
}
x1 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0, a=c(-6, -9), b=c(4, 11), trunc=c(3, 3))
x2 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-6, -9), b=c(4, 11), trunc=c(3, 3))
x3 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-100, -100), b=c(100, 100),
trunc=c(3, 3))
x4 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=-0.7, a=c(-6, -9), b=c(4, 11),
trunc=c(3, 3))
x5 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=-0.9, a=c(-6, -9), b=c(4, 11),
trunc=c(3, 3))
x6 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(0, 0), trunc=c(0, 2))
x7 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-1, 1), trunc=c(0, 2))
x8 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-1, 1), trunc=c(1, 2))
x9 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-1.5, 0.5), b=c(-0.5, 1.5),
trunc=c(3, 3))
x10 <- exam2(1000, mu=c(-1, 1), sigma=c(1, sqrt(2)), rho=0.7, a=c(-1.5, 0.5), b=c(-0.5, 1.5),
trunc=c(4, 3))
}
Run the code above in your browser using DataLab