# Make some data for testing code:

n=40
x=1:n*2/n
a=20
b=-5
c=4
mu=a*exp(b+c*x)/(1+exp(b+c*x))
y=mu+rnorm(n)*4
plot(x,y)
lines(x,mu,lty=3)

# initial guess: 

beta=c(18,-5,6)
yhat = beta[1]*exp(beta[2]+beta[3]*x)/(1+exp(beta[2]+beta[3]*x))

# tolerance: 

sm=0.0001

# loop

d=1

while(d>sm){
	
# build the W matrix:

wmat=matrix(1:(3*n)*0,nrow=n)
for(i in 1:n){
	ce=beta[2]+beta[3]*x[i]
	wmat[i,1]= exp(ce)/(1+exp(ce))
	wmat[i,2]= beta[1]*exp(ce)*(1+exp(ce))^(-2)
	wmat[i,3]= beta[1]*exp(ce)*x[i]*(1+exp(ce))^(-2)
}	

# iterate!

beta1 = beta + solve(t(wmat)%*%wmat)%*%t(wmat)%*%(y-yhat)

d=sum((beta1-beta)^2)
beta=beta1
yhat = beta[1]*exp(beta[2]+beta[3]*x)/(1+exp(beta[2]+beta[3]*x))
lines(x,yhat,lty=2)
}