I have some issue using nls()
to estimate parameters. I have a following set of functions to explain some data in hand:
funk1 <- function(a,x) { x^2*exp(-(l*(1-exp(-r*a))/r)) }
funk2 <- function(x) { sapply(x, function (s)
{ integrate(funk1, lower = 0, upper = s, x=s)$value }) }
I am trying to fit funk2
to y
:
y <- sort(runif(100, 0, 10^8))
When I use nls()
:
nls(y ~ funk2(z1$days.post.bmt), data= z1, start=list(l=0.02, r=0.002), trace=T)
it shows me following error:
Error in f(x, ...) : object 'l' not found
Isn't the whole point of nls()
to substitute different values for parameters l
and r
from parameter space to fit the function by minimizing SSR and give the parameter estimates? why it needs value of l
for it to work? I am definitely missing something big here. Please help!
Thanks in advance!
You must pass parameters l
and r
as function arguments of funk1
and funk2
.
funk1 <- function(a,x,l,r) {
x^2*exp(-(l*(1-exp(-r*a))/r))
}
funk2 <- function(x,l,r) {
sapply(x, function (s) {
integrate(funk1, lower = 0, upper = s, x=s, l=l, r=r)$value
})
}
I will generate some data to test:
z <- data.frame(days.post.bmt = 1:100,
y = funk2(1:100, l = 1, r = 1) + rpois(100, 1:100))
nls(y ~ funk2(days.post.bmt,l,r), data = z, start = list(l = 0.5, r = 0.5))
#Nonlinear regression model
# model: y ~ funk2(days.post.bmt, l, r)
# data: z
# l r
#0.9405 0.9400
# residual sum-of-squares: 6709
#Number of iterations to convergence: 5
#Achieved convergence tolerance: 2.354e-07
As a counter example, consider:
bad_funk1 <- function(a,x) {
x^2*exp(-(l*(1-exp(-r*a))/r))
}
bad_funk2 <- function(x) {
sapply(x, function (s) {
integrate(funk1, lower = 0, upper = s, x=s)$value
})
}
nls(y ~ bad_funk2(days.post.bmt), data = z, start = list(l = 0.5, r = 0.5))
# Error in f(x, ...) (from #2) : argument "l" is missing, with no default