MLE for the logistic regression: Programming problem I'm learning GLM models and found the following log-likelihood by hand and wrote the code below it as an inference exercise. However my result for the parameter is different from the glm() function's result. I'd appreciate if someone pointed out what is wrong, the calculations or my code. Thanks!
\begin{align}
\mathcal{L}(p_{i}, n_{i}, y) &= \prod_{i=1}^{N}{n_{i} \choose k_{i}} p_{i}^{k_{i}}(1-p_{i})^{n_{i}-k_{i}} \\
\ell(p_{i}, n_{i}, y) &= \sum_{i=1}^{N}\log {n_{i} \choose k_{i}} + k_{i} \log p_{i} + (n_{i} - k_{i})\log(1-p_{i}) \\
&= \sum_{i=1}^{N}\log {n_{i} \choose k_{i}} + k_{i}\log \frac{p_{i}}{1-p_{i}} + n_{i}\log(1-p_{i})\\
&= \sum_{i=1}^{N}\log {n_{i} \choose k_{i}} + k_{i} x'_{i}\beta - n_{i}\log(1 + \exp(x'_{i}\beta))\\
\frac{\partial\ell}{\partial \beta_{j}} &= \sum_{i=1}^{N}k_{i}x_{ij} - n_{i}x_{ij}\frac{\exp(x'_{i}\beta)}{1+ \exp(x'_{i}\beta)}
\end{align}
set.seed(100)
 dados <- rbinom(n = 40, size = 5, prob = 0.1)
 cov <- rep(1, times = 40)
 n <- 5
 soma <- 0

    ll <-  function(beta){
            inv_log <- exp(cov*beta)/(1 + exp(cov*beta))
            fn_dados <- dados*cov - n*inv_log*cov 
            soma <- sum(fn_dados)
              return(-soma)
    }

    optim(ll, par = c(0), method = 'BFGS')

    fail <- n - dados
    glm(formula = cbind(dados, fail) ~ 1, family= binomial())

EDIT: I found the solution by myself, me and all my disattention were optimizing the derivative, not the log-likelihood.
GitHub link as said in the comments:
https://github.com/sergioandrad/RegressaoLogistica-
 A: As you mentioned in your edit, you didn't have the correct objective function. I get the same answer as glm when i minimize the following negative log likelihood function:
neg_log_like <- function(beta){
  -sum(dados*cov*beta - n*log(1+exp(cov*beta)))
}

A: As said in my edit, I was optimizing the wrong function. My new code is the following:
set.seed(100)
dados <- rbinom(n = 40, size = 5, prob = 0.1)
cov <- rep(1, times = 40)
n <- 5
soma <- 0

ll <-  function(beta){
        fato <- log(factorial(n)/(factorial(n-dados)*factorial(dados)))
        fn_dados <- fato + dados*cov*beta - n*log(1+ exp(cov*beta))
        soma <- sum(fn_dados)
          return(-soma)
}

optim(ll, par = c(0), method = 'Brent', upper = 100, lower = -100)

fail <- n - dados

model <- glm(formula = cbind(dados, fail) ~ 1, family= binomial())
logLik(model)

You can see that altough the factorials do not change the parameter estimation, I left them there because most softwares release the full value, so this way we can compare likelihoods, deviances and other interesting measure of fit criteria.
