# 8: n=36 with data set X  beta=(-7,0,1,1,1,2) - Probit
# 9: n=100 with new X  beta=(-6,0.02,0,1,0) - Probit

# source("gen.sparse.corr.model.r")

scenario<-9  #8 or 9

if(scenario==8){
  #just load X and beta0 to generate the Y
  load("~/R/ERGM/Marginal/Psi_fixed_beta_degree(2)_n36_M30000.RData")
  # beta=(-7,0,1,1,1,2)
  n<-36
  sigma<-sqrt(0.5)
}

if(scenario==9){
  #just load X and beta0 to generate the Y
  load("~/R/ERGM/Marginal/Psi_fixed_beta_degree(2)_n100_M30000.RData")
  # beta=(-6,0.02,0,1,0)
  n <- 100
  sigma<- 1 #sqrt(0.5)
}

N<-n*(n-1)/2


# generate 30,000 data sets

M<-3e4
# 30000 times 630 (sim times N)
Y.MCMC <- matrix(NA,M,N)
eta<-as.numeric(X%*%beta.old)
ind.lt<-lower.tri(diag(n))

for(i in 1:M){
  
  ei <- rnorm(n,sd=sigma)
  ej <- rnorm(n,sd=sigma)
  
  hilf.matrix <- matrix(ei,n,n,byrow=F)+ matrix(ej,n,n,byrow=T)
  
  
  Y.star<-  1*((eta + hilf.matrix[ind.lt])>0)
  
  Y.MCMC[i,]<-Y.star
  
}#end for(i in 1:M){


save(M,Y.obs,X,Y.MCMC,beta.old,file=paste("Probit_n",n,"_M",M,".RData",sep=""))


cat("FINISHED\n")