3
$\begingroup$

I'm trying to fit a 4 parameter boltzmann sigmoid and get an error: "Error in nls(y ~ a0 + (a1 - a0)/(1 + exp((a2 - x)/a3)), start = list(a0 = max(y), : singular gradient"

I have figured out that the code runs if the a3 parameter is set to .9 or greater, or -.9 or less. Does anyone have the reason this is? I want to provide a starting parameter for a3 as the slope according to the description on this website: http://www.originlab.com/doc/Origin-Help/Boltzmann-FitFunc . That is why I have the linear fit coefficient a3.s, but the result is < .9 and I get the error. Is there a way to estimate a3.s prior to use as starting parameter for nls? I am simply using a linear fit of the midpoint of the sigmoid +/- 10 x units - is that the correct interpretation of the a3 parameter?

Here is my code:

#fit boltzman sigmoid
a0.s=max(y); a1.s=min(y); a2.i=which.min( abs(((a0.s+a1.s)/2) - y) ); a2.s=x[a2.i]
lin.x.i=x<a2.s+10 & x>a2.s-10
a3.s=unname(coef(lm(y[lin.x.i]~so[lin.x.i]))[2])
fit <- nls(y ~ a0 + (a1-a0)/(1+exp((a2-x)/a3)), 
start=list(a0=max(y), a1=min(y), a2=a2.s,a3=.9) , trace=TRUE)
params=coef(fit)
curve(params[1]+(params[2]-params[1])/(1+exp((params[3]-x)/params[4])), 1,100,col='black',add=T,type='l')

Here is the data:

x=c( 75,  40,  90,  55, 15, 100,  10,  70,  90,  50,  15,   5,   5,  70, 100,  20,  60,  65,  20,  50,  30,  85,  60,  80,  55,  40,  45,  95,  10,  55, 60,  10,  35,  80,  75,  25,  30,   5,  35,  50, 100,  40,  30,  80,  20,  45,  25,  25,  95,  95,  65,  35,  90,  85,  70,  15,  75,  45,  85,  65);

y=c(4.673686, 0.034781, 5.014355, 0.843847, 0.013337, 4.214557, 0.015299, 5.017280, 4.327815, 0.041139, 0.008704, 0.007437, 0.005125, 4.725786, 3.869776, 0.018725, 4.514051, 3.232932, 0.012979, 0.257651, 0.028170, 4.723512, 2.676991, 5.018232, 0.633399, 0.040133, 0.051864, 5.019395, 0.006505, 0.642376, 2.752317, 0.010827, 0.029303, 4.050711, 3.698887, 0.018385, 0.029491, 0.013894, 0.032034, 0.053761, 5.029349, 0.038272, 0.032619, 5.030450, 0.022356, 0.053421, 0.025370, 0.024763, 4.948973, 3.254528, 1.149153, 0.038530, 4.612227, 4.048692, 4.809153, 0.016246, 5.014711, 0.062841, 5.026961, 2.951881)

Related to this question: the formula on the linked to website has a slight variation in the equation, where the the a2 parameter is used in the form "x-a2", while the equation I provided, and got from my data acquisition software's curve fitting function is the one I provided in the code with "a2-x". Which form of the Boltzmann is correct? Does the difference matter?

$\endgroup$

2 Answers 2

3
$\begingroup$

Why are you reinventing the wheel? Use the native function SSfpl for your model.

fitr <- nls(y ~ SSfpl(x, a1, a0, a2, ma3))

The parameter ma3 is -a3 in your notation, but otherwise the parametrization is identical, and you get slightly better convergence.

You should probably be using weighted least squares, since ordinary least squares assumes the variability of $y-E(y)$ does not depend on the value of $x$; which is clearly violated in your data.

$\endgroup$
-1
$\begingroup$

I guess this would also give an idea:

df <- data.frame(x,y)

# using drm() in "aomisc" package should give a solution
model.df <- drm(y ~ x, fct = L.4(), data = df)
summary(model.df)
plot(model.df, log="", main = "df")
$\endgroup$
1
  • $\begingroup$ Although implementation is often mixed with substantive content in questions, we are supposed to be a site for providing information about statistics, machine learning, etc., not code. It can be good to provide code as well, but please elaborate your substantive answer in text for people who don't read this language well enough to recognize & extract the answer from the code. $\endgroup$
    – Sycorax
    Oct 28, 2021 at 14:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.