### Simulate Bivariate Normal Data: BVN(mu.x, mu.y, sigmasq.x, sigmasq.y, rho)

### BVN(5,10,4,9,0.4)

n <- 1000

mu.x <- 5

mu.y <- 10

mu <- c(mu.x, mu.y)

sigmasq.x <- 4

sigmasq.y <- 9

sigma <- c(sigmasq.x, sigmasq.y)

rho <- 0.4

X <- rmvnorm(n, mean = mu, sd = sigma, rho = rho)

# Plot the data in a scatterplot

plot(X)

# Compute the mean

apply(X, 2, mean) # apply the mean finction to the columns (2)

# Compute the correlation coefficient

cor(X)

# Fit a regression line to the data

x <- X[,1]

y <- X[,2]

mvn.fit <- lm(y ~ x)

summary(mvn.fit)

y.fit <- mvn.fit$fit

plot(x,y) # plot y vs x

lines(x,y.fit) # plot y.fit vs x on the same plot.

# Exercise: BVN(200, 500, 625, 900, -0.6)