library(Matrix)
library(sandwich)
mydata <- read.csv(http://www.ats.ucla.edu/stat/data/binary.csv)
mylogit <- glm(admit ~ gre + gpa, data = mydata, family = "binomial")
X <- as.matrix(cbind(1, mydata[,c('gre','gpa')]))
n <- nrow(X)
pi<-mylogit$fit
w<-pi*(1-pi)
v<-Diagonal(n, x = w)
var_b<-solve(t(X)%*%v%*%X)
var_b
vcov(mylogit)