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)]
No need to mark, too long in a comment.
Something like this:
DT[, c("sigma", "a", "b", "c") :=
tryCatch({
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)))
}, error=function(e) NA_real_),
by=.(groups)]
output:
x y groups sigma a b c
1: 1 15.4 1 0.2986243 0.5265405 0.4565363 14.56728
2: 2 16.0 1 0.2986243 0.5265405 0.4565363 14.56728
3: 3 16.4 1 0.2986243 0.5265405 0.4565363 14.56728
4: 4 17.7 1 0.2986243 0.5265405 0.4565363 14.56728
5: 5 20.0 1 0.2986243 0.5265405 0.4565363 14.56728
6: 6 23.0 1 0.2986243 0.5265405 0.4565363 14.56728
7: 7 27.0 1 0.2986243 0.5265405 0.4565363 14.56728
8: 8 35.0 1 0.2986243 0.5265405 0.4565363 14.56728
9: 1 25.4 2 NA NA NA NA
10: 2 26.0 2 NA NA NA NA
11: 3 26.4 2 NA NA NA NA
12: 4 27.7 2 NA NA NA NA
13: 5 30.0 2 NA NA NA NA
14: 6 33.0 2 NA NA NA NA
15: 7 37.0 2 NA NA NA NA
16: 8 3.0 2 NA NA NA NA
17: 1 27.4 3 0.2986243 0.5265401 0.4565364 26.56728
18: 2 28.0 3 0.2986243 0.5265401 0.4565364 26.56728
19: 3 28.4 3 0.2986243 0.5265401 0.4565364 26.56728
20: 4 29.7 3 0.2986243 0.5265401 0.4565364 26.56728
21: 5 32.0 3 0.2986243 0.5265401 0.4565364 26.56728
22: 6 35.0 3 0.2986243 0.5265401 0.4565364 26.56728
23: 7 39.0 3 0.2986243 0.5265401 0.4565364 26.56728
24: 8 47.0 3 0.2986243 0.5265401 0.4565364 26.56728
x y groups sigma a b c