Nonlinear Models

### Least Square Estimation
x <- c(0.02,0.02,0.06,0.06,0.11,0.11,0.22,0.22,0.56,0.56,1.1,1.1)
y <- c(76,47,97,107,123,139,159,152,191,201,207,200)

# mu(x) = p[1]*x / (p[2]+x)

fn <- function(p) sum((y-(p[1]*x)/(p[2]+x))^2)
out <- nlm(fn, p=c(200, 0.1), hessian=T)

xfit <-seq(.02, 1.1, .05)
yfit <- out$estimate[1]*xfit/(out$estimate[2]+xfit)

plot(x,y)
lines(spline(xfit,yfit))


### Maximum Likelihood Estimation

# N(p[1], p[2])

fn <- function(p) length(x)/2*log(p[2])+sum((x-p[1])^2)/(2*p[2])
gr <- function(p) c(sum(x-p[1])/p[2], 0.5*length(x)/p[2]-0.5*sum((x-p[1])^2)/p[2]^2)

# 100 random samples from N(20, 9)
x <- rnorm(100, 20, 3)

# MLEs
mean(x); var(x)*(length(x)-1)/length(x)

nlm(p=c(20,9), fn, gr)$estimate

optim(c(20,9), fn, gr)$par

# nls()
out<-nls(y ~ alpha*x/(beta+x), start=list(alpha=200, beta=0.1))

plot(y~x)
lines(x, predict(out), col='blue')