SNStatComp/rtrim

Error in update_alpha(method) in model = 2 with all changepoints set and serialcor = TRUE

Closed this issue · 2 comments

I tried to run model = 2 with all changepoints and serialcor = TRUE (the data variable is defined at the end of this post):

library(rtrim)

# set all possible changepoints
changepoints <- head(min(data$year):max(data$year) - min(data$year) + 1, -1) 
# run the model
c1 <- trim(count ~ site + year, data, model = 2, overdisp = TRUE, serialcor = TRUE, max_iter = 1000, 
    changepoints = changepoints, autodelete = TRUE, stepwise = FALSE)

but I got an error:

Error in update_alpha(method) : non-finite alpha problem 1a
In addition: Warning message:
In log(z_t %*% exp(B_i %*% beta - log(wt_i))) : NaNs produced

Setting serialcor = FALSE removes the error.

Now it is not clear if this is a bug or if I ran the function trim() with wrong parameter combinations. In that case I would expect a comprehensible error that would tell me so. The manual doesn't seem to mention it either.

Here are the data:

data <- structure(list(site = c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 
11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 
3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 
17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 
10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 
2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 
16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 
9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 
1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 
15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 
8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 
21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 
14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 
7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 
20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 
13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 
5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 
19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 
12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 
4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 
18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 
11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 
3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 
17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 
10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 1L, 
2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 13L, 14L, 15L, 
16L, 17L, 18L, 19L, 20L, 21L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 
9L, 10L, 11L, 12L, 13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L
), year = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 
2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 
3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 
4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 
5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 
6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 
7, 7, 7, 7, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 
9, 9, 9, 9, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 
10, 10, 10, 10, 10, 10, 10, 10, 11, 11, 11, 11, 11, 11, 11, 11, 
11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 11, 12, 12, 12, 
12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 
12, 12, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 13, 
13, 13, 13, 13, 13, 13, 13, 14, 14, 14, 14, 14, 14, 14, 14, 14, 
14, 14, 14, 14, 14, 14, 14, 14, 14, 14, 14, 14, 15, 15, 15, 15, 
15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 
15, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 
16, 16, 16, 16, 16, 16, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 
17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 17, 18, 18, 18, 18, 18, 
18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18
), count = c(575, 529, 982, 54, 54, 455, 302, 54, 1010, 50, 164, 
99, 104, 81, 92, 130, 1424, 11, 19, 873, 108, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, 10, 266, 99, 53, 53, 529, 198, 53, 1026, 
52, 45, 100, 94, 49, 102, 97, 231, 0, 0, NA, NA, 5, 148, 165, 
54, 0, 468, 65, 4, 19, 52, 29, 94, 106, 46, 99, 80, 61, 0, 0, 
2, NA, 11, 65, 153, 0, 0, 544, 62, 0, 16, 21, 28, 0, 109, 51, 
0, 5, 48, 0, 0, 10, 5, 4, 92, 0, 0, 0, 539, 62, 0, 7, 0, 21, 
53, 71, 0, 0, 0, 18, 0, 0, 9, NA, 4, 237, 5, 0, 0, 489, 19, 0, 
0, 1, 19, 49, 61, 0, 0, 0, 34, 0, 0, NA, NA, 0, 143, 0, 0, 0, 
262, 21, 0, 0, 6, 21, 51, 69, 0, 0, 0, 11, 0, 0, NA, NA, NA, 
95, NA, NA, 0, 217, 20, 0, 0, 6, 20, 47, 76, 0, 0, 0, 23, NA, 
0, 33, NA, 0, 29, 1, 0, 0, 156, 19, 0, 0, 7, 19, 55, NA, 0, 0, 
0, 179, 0, 0, 43, 3, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 0, 57, 0, 0, 0, 147, 
20, 0, 0, 5, 10, 0, 60, 0, 0, 0, 6, 0, 0, 0, 1, 0, 105, 0, 0, 
0, 183, 21, 0, 0, 5, 19, 46, 59, 0, 0, 0, 11, 0, 0, 37, 2, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, 0, 267, 0, 0, 0, 321, 48, 0, 0, 5, 45, 0, 54, 
0, 0, 0, 239, NA, 0, NA, NA, 0, 257, 0, 0, 0, 277, 50, 0, 0, 
15, 55, 0, 57, 0, 0, 0, 251, 0, 0, 0, 2, 0, 241, 0, 0, 0, 324, 
48, 0, 0, 19, 45, 0, 70, 0, 0, 0, 297, 0, 0, 0, 5)), row.names = c(NA, 
-378L), class = "data.frame")

Interestingly enough, if I add a covariate, the error disappears and the calculation goes through, even with serialcor = TRUE! Not something I would expect!

data$covar <- as.factor(ifelse(data$site <= 10, 1, 2))
c2 <- trim(count ~ site + year + covar, data, model = 2, overdisp = TRUE, serialcor = TRUE, max_iter = 5000, 
    changepoints = changepoints, autodelete = TRUE, stepwise = FALSE)

However, the model didn't converge properly (summary says "No convergence reached after 5000 iterations").

This data set seems to be an edge case, but rtrim should accept it.
I could not fix the original updating of the alpha parameters, but now the algorithm falls back to a second, more robust, method if this error occurs. In that case it wll now issue a warning, suggesting o use that alternative method throughout the run.

Implemented in the new version 2.3.0 (published today). Thanks for pointing this out.