## EXAMPLES OF 2-D GAUSSIAN VECTORS

N<-5000
z.1 <- rnorm(N,0,1)
z.2 <- rnorm(N,0,1)
plot(z.1,z.2)

A <- matrix(nrow=2,ncol=2)
A[1,1] <- sqrt(2)
A[2,1] <- sqrt(2)
A[1,2] <- -1/sqrt(2)
A[2,2] <- 1/sqrt(2)
z <- matrix(nrow=2,ncol=N)
z[1,] <- z.1
z[2,] <- z.2
x <- A%*%z
plot(x[1,],x[2,])

cov(t(x))

Q<-A%*%t(A)
Q.inv <- solve(Q)
Q.inv

--------------------------------------

## 2-D DATA GAUSSIAN FIT

n <- 20
x <- rnorm(n,1,1)
y <- rnorm(n,1,1)

plot(x,y)

p.hat <- sum(((sign(x)+1)/2)*((sign(y)+1)/2))/n
p.hat

Q.12 <- cov(x,y)
Q.11 <- var(x)
Q.22 <- var(y)
mu.1 <- mean(x)
mu.2 <-mean(y)

Data <- matrix(nrow=20,ncol=2)
Data[,1] <- x
Data[,2] <- y
Q <- cov(Data)

require(mgcv)
A<-mroot(Q)
A
t(A)
A%*%t(A)-Q

N<-10000
z.1 <- rnorm(N,0,1)
z.2 <- rnorm(N,0,1)
xx <- A[1,1]*z.1 + A[1,2]*z.2 +mu.1
yy <- A[2,1]*z.1 + A[2,2]*z.2 +mu.2
plot(xx,yy)

p <- sum(((sign(xx)+1)/2)*((sign(yy)+1)/2))/N
p

------------------
## COMPARISON OF STABILITY BETWEEN ORIGINAL SAMPLE AND MODEL

p.hat <- 1:50
p <- 1:50
for (i in 1:50) {
n <- 20
x <- rnorm(n,1,1)
y <- rnorm(n,1,1)
mu.1 <- mean(x)
mu.2 <-mean(y)
Data <- matrix(nrow=n,ncol=2)
Data[,1] <- x
Data[,2] <- y
Q <- cov(Data)
A<-mroot(Q)
N<-1000
z.1 <- rnorm(N,0,1)
z.2 <- rnorm(N,0,1)
xx <- A[1,1]*z.1 + A[1,2]*z.2 +mu.1
yy <- A[2,1]*z.1 + A[2,2]*z.2 +mu.2
p.hat[i]<-sum(((sign(x)+1)/2)*((sign(y)+1)/2))/n
p[i]<-sum(((sign(xx)+1)/2)*((sign(yy)+1)/2))/N
}
sd(p.hat)
sd(p)


--------------------------------
--------------------------------