MLE cumulative gaussian

library("ggplot2")

cumNorm <- function(x, p, m) (1/m) + (1 - 1/m) * pnorm(x, p[1], p[2])

likeL <- function(p, d, FUN, m) {
    pr <- FUN(d$x, p, m)
    -sum(d$nyes * log(pr) + d$nno * log(1 - pr))
}

curveFun <- function(pIni, d, FUN, m, xmin = min(d$x), xmax = max(d$x)) {
    p <- optim(pIni, likeL, d = d, FUN = FUN, m = m)$par
    xseq <- seq(xmin, xmax, len = 100)
    yseq <- FUN(xseq, p, m)
    data.frame(x = xseq, y = yseq)
}

x <- c(1, 2, 3, 4)
nyes <- c(70, 100, 100, 100)
nno <- c(30, 0, 0, 0)

dat <- data.frame(x, nyes, nno)
dat$p <- dat$nyes/(dat$nyes + dat$nno)

psycho <- curveFun(c(1, 1), dat, cumNorm, 2, xmin = 0)

p1 <- ggplot() + geom_point(data = dat, aes(x = x, y = p)) + geom_line(data = psycho, 
    aes(x = x, y = y))
p1