p=c(0,1,1)
L=rep(0,100)

f<-function(p){
for(i in 1:100){

L[i]=log(dnorm(Y[i],p[1]+p[2]*X[i],p[3]))
}
Lfinal=sum(L)

;-Lfinal
}
f(p)

m=nlminb(p,f,lower=c(-Inf,-Inf,0),upper=c(Inf,Inf,Inf))
m

library("nlme")
h=fdHess(m$par,f)
ih=solve(h$Hessian)
se=sqrt(diag(ih))
MLE=m$par
cbind(p,MLE,se)
se=sqrt(abs(diag(ih)))

Y=rbinom(100,1,.5)
X=rnorm(100,0,1)
Data=data.frame(Y,X)
Data

p=c(0,1)
L=rep(0,100)

f<-function(a){
for(i in 1:100){
p[i]=exp(a[1]+a[2]*X[i])/(1+exp(a[1]+a[2]*X[i]))
L[i]=log(dbinom(Y[i],1,p[i]))            
}
  
Lfinal=sum(L)
             
             
;-Lfinal
}

f(p)

library("nlme")
m=nlminb(a,f,lower=c(-Inf,-Inf),upper=c(Inf,Inf))

h=fdHess(m$par,f)
ih=solve(h$Hessian)
se=sqrt(diag(ih))
MLE=m$par
cbind(a,MLE,se)

x=c(1,-1,2,0,1,3)
M=matrix(x,ncol=1,nrow=6)
A=M%*%t(M)+diag(6)
det(A)