Determining the power associated with the logistic model

Question

Although there has been some detailed discussions about power analysis on this website (for example here and here), the answer provided to this question has outlines the steps to simulating a power analysis, here.

Say we take some data (data was linked to a bootstrapping question)

We create a regression that will predict admit based on the two continous variables gpa and gre

• Now we have a n=400.
• We can then elect our power level, alpha = 0.5
• The effect size you would like to detect, e.g., odds ratios (we obtain this from our regression)

So in following the detailed method provided by @gung here, I want to run the simulation. Here is the code I have adjusted, but my output is not correct. Can someone outline what I have not understood

mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv")
head(mydata)

set.seed(1234)

my.mod <- glm(admit ~ gre + gpa , data = mydata, family = "binomial")


repetitions <- length(mydata$admit)

gre <- mydata$gre
gpa <- mydata$gpa


significant = matrix(nrow=repetitions, ncol=4)

for(i in 1:repetitions){
  responses          = mydata$admit
  #responses          = rbinom(n=N, size=1, prob=mydata$admit)      # we can interchange this comment
  model              = glm(responses ~ gre + gpa, family = binomial(link="logit"))
  significant[i,1:2] = (summary(model)$coefficients[2:3,4]<.05)
  significant[i,3]   = sum(significant[i,1:2])
  modelDev           = model$null.deviance-model$deviance
  significant[i,4]   = (1-pchisq(modelDev, 2))<.05
}



sum(significant[,1])/repetitions      # pre-specified effect power for gre

sum(significant[,2])/repetitions      # pre-specified effect power for gpa

sum(significant[,4])/repetitions  # power for likelihood ratio test of model

sum(significant[,3]==2)/repetitions   # all effects power

sum(significant[,3]>0)/repetitions    # any effect power