non-numeric argument to binary operator error in nls

Question

I am trying to run the following code but it gives me:

Error in qr.default(.swts * attr(rhs, "gradient")) : NA/NaN/Inf in foreign function call (arg 1) In addition: Warning message: In log(.expr4) : NaNs produced

Can you please help me about? Thanks!

model <- deriv( ~ c*(1+b*(q-1)*t)^(1/(1-q)),  c("c", "b", "q"), function (t, c, b,q){})

nls(Frequency ~ model(t, c, b, q), data=DF,start=list(c = 1, b = 1.5, q =0.5))

Following you can see a part of data, for which I am trying to fit a q-exponential distribution functin as I explained above. I am uiung nls function in R to obtain estimates (q-exponential) for the given data.

t Frequency 0 195746 1 93938 2 53181 3 31853 4 19856 5 12182 6 7847 7 5459 8 4325 9 3203 10 2750

Answer 1

tl;dr you need to work harder to think about/find reasonable starting values.

Set up data and model:

 dd <- data.frame(t=0:10,
            Frequency=c(195746,93938,53181,31853,19856,12182,
                        7847,5459,4325,3203,2750))
model <- deriv( ~ c*(1+b*(q-1)*t)^(1/(1-q)),
               c("c", "b", "q"), function (t, c, b,q){})

If you just evaluate model() at the initial parameters you can see that some of the derivatives with respect to q are NaN (it is always worthwhile to try evaluating the objective function and gradient at the starting values to make sure they make sense). In addition to this, if you look at the values of the objective function you can see that they're nowhere close to the data (they increase from 1 to 42, whereas the data decrease from 200000 to 3000 ...)

par0 <- list(c = 1, b = 1.5, q =0.5)
with(par0,model(dd$t,c,b,q))
##  [1]  1.0000  0.0625  0.2500  1.5625  4.0000  7.5625 12.2500 18.0625 25.0000
## [10] 33.0625 42.2500
## attr(,"gradient")
##             c     b         q
##  [1,]  1.0000  0.00 0.0000000
##  [2,]  0.0625 -0.25 0.4034264
##  [3,]  0.2500  1.00       NaN
##  [4,]  1.5625  3.75       NaN
...

The problem is that raising a negative value to a fractional power gives NaN. When q is <1, the base will go negative when t is sufficiently large ...

with(par1,model(dd$t,c,b,q))

I tried setting q=0 (so that the base value of the exponent is 1), but that still encounters the same problem.

What if we start from values where q-1 is >0 so that the base will be positive?

par2 <- list(c = 1, b = 1.5, q =2)
nls(Frequency ~ model(t, c, b, q), data=dd,start=par2)

Error in nls(Frequency ~ model(t, c, b, q), data = dd, start = par2) : step factor 0.000488281 reduced below 'minFactor' of 0.000976562

OK, that's getting a little better but still failing. What if we try to start on a more reasonable scale, i.e. set c to 100000?

par3 <- list(c = 1e5, b = 1.5, q =2)
fit3 <- nls(Frequency ~ model(t, c, b, q), data=dd,start=par3)

Plot results:

par(las=1,bty="l")
plot(Frequency/1000~t,data=dd)
newdat <- data.frame(t=seq(0,10,length=51))
lines(newdat$t,predict(fit3,newdata=newdat)/1000)