I'm trying to fit a Boltzmann sigmoid 1/(1+exp((x-p1)/p2)) to this small experimental dataset:
xdata <- c(-60,-50,-40,-30,-20,-10,-0,10)
ydata <- c(0.04, 0.09, 0.38, 0.63, 0.79, 1, 0.83, 0.56)
I know that it is pretty simple to do it. For example, using nls:
fit <-nls(ydata ~ 1/(1+exp((xdata-p1)/p2)),start=list(p1=mean(xdata),p2=-5))
I get the following results:
Formula: ydata ~ 1/(1 + exp((xdata - p1)/p2))
Parameters:
Estimate Std. Error t value Pr(>|t|)
p1 -33.671 4.755 -7.081 0.000398 ***
p2 -10.336 4.312 -2.397 0.053490 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.1904 on 6 degrees of freedom
Number of iterations to convergence: 13
Achieved convergence tolerance: 7.079e-06
However, I need (due to theoretical reasons) the fitted curve to pass precisely through the point (-70, 0). Although the value of the fitted expression showed above passes near zero at x = -70, it is not exactly zero, which is not what I want.
So, the question is: Is there a way to tell nls (or some other function) to fit the same expression but forcing it to pass through a specified point?
Update:
As it has been mentioned in the comments, it is mathematically impossible to force the fit to go through the point (-70,0) using the function I provided (the Boltzmann sigmoid). On the other hand, @Cleb and @BenBolker have explained how to force the fit to go through any other point, for instance (-50, 0.09).
Building on @Cleb's answer, here's a way to pick a specified point the function must pass through and solve the resulting equation for one of the parameters:
dd <- data.frame(x=c(-60,-50,-40,-30,-20,-10,-0,10),
y=c(0.04, 0.09, 0.38, 0.63, 0.79, 1, 0.83, 0.56))
Initial fit (using plogis() rather than 1/(1+exp(-...)) for convenience):
fit <- nls(y ~ plogis(-(x-p1)/p2),
data=dd,
start=list(p1=mean(dd$x),p2=-5))
Now plug in (x3,y3) and solve for p2:
y3 = 1/(1+exp((x-p1)/p2))
logit(x) = qlogis(-x) = log(x/(1-x))
e.g. plogis(2)=0.88 -> qlogis(0.88)=2
qlogis(y3) = -(x-p1)/p2
p2 = -(x3-p1)/qlogis(y3)
Set up a function and plug it in for p2:
p2 <- function(p1,x,y) {
-(x-p1)/qlogis(y)
}
fit2 <- nls(y ~ plogis(-(x-p1)/p2(p1,dd$x[3],dd$y[3])),
data=dd,
start=list(p1=mean(dd$x)))
Plot the results:
plot(y~x,data=dd,ylim=c(0,1.1))
xr <- data.frame(x = seq(min(dd$x),max(dd$x),len=200))
lines(xr$x,predict(fit,newdata=xr))
lines(xr$x,predict(fit2,newdata=xr),col=2)