a = sum(beta*c(1,gammai[i,])) return(exp(a)/(1+exp(a))) }) U = runif(length(pydata)) ydata = as.numeric(U