This question builds off a previous one I asked here Exponential curve fitting with nls using data.table groups.
I am fitting exponential curves using nls
to a data table object by multiple groups. Not all of the data appears to fit an exponential model and nls
sometimes throws an error, stopping all further calculations for the remaining groups.
I've attached a MWE below with my attempt at using a tryCatch
to skip problematic groups but I end up with the error output for all of my new columns. How can I skip calculating nls
values for my problematic groups?
## Example data table
DT <- data.table(
x = c(1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8,
1,2,3,4,5,6,7,8),
y = c(15.4,16,16.4,17.7,20,23,27,35,
25.4,26,26.4,27.7,30,33,37,45,
27.4,28,28.4,29.7,32,35,39,47),
groups = c(1,1,1,1,1,1,1,1,
2,2,2,2,2,2,2,2,
3,3,3,3,3,3,3,3)
)
## Fit exponential curve using starting values a,b,c for each group
DT[, c("sigma", "a", "b", "c") := {
c.0 <- min(y) * 0.5
model.0 <- lm(log(y - c.0) ~ x, data=.SD)
start <- list(a=exp(coef(model.0)[1]), b=coef(model.0)[2], c=c.0)
model <- nls(y ~ a * exp(b * x) + c,
data=.SD,
start=start,
control=nls.control(maxiter=500))
c(.(sigma=summary(model)$sigma), as.list(coef(model)))
},
by=.(groups)]
## Modify data table to ruin nls model for group 2
set(DT, i=16L, j="y", value=3)
## Calculation works for group 1 but stops for group 2 and onwards
DT[, c("sigma", "a", "b", "c") := {
c.0 <- min(y) * 0.5
model.0 <- lm(log(y - c.0) ~ x, data=.SD)
start <- list(a=exp(coef(model.0)[1]), b=coef(model.0)[2], c=c.0)
model <- nls(y ~ a * exp(b * x) + c,
data=.SD,
start=start,
control=nls.control(maxiter=500))
c(.(sigma=summary(model)$sigma), as.list(coef(model)))
},
by=.(groups)]
## My poor attempt at using a tryCatch just gives NA to every column
DT[, c("sigma","a", "b", "c") := {
c.0 <- min(y) * 0.5
model.0 <- lm(log(y - c.0) ~ x, data=.SD)
start <- list(a=exp(coef(model.0)[1]), b=coef(model.0)[2], c=c.0)
model <- tryCatch(
{
nls(y ~ a * exp(b * x) + c,
data=.SD,
start=start,
control=nls.control(maxiter=500))
c(.(sigma=summary(model)$sigma), as.list(coef(model)))
},
error=function(err){
return(NA_real_)
}
)
},
by=.(groups)]