Segmented Package Errors: Checkpoint Confusion

Question

Using a segmented package to create piecewise linear regression, I see an error when I try to set my own breakpoints; It seems only when I try to install more than two.

(EDIT) Here is the code I'm using:

# data bullard <- structure(list(Rt = c(0, 4.0054, 25.1858, 27.9998, 35.7259, 39.0769, 45.1805, 45.6717, 48.3419, 51.5661, 64.1578, 66.828, 111.1613, 114.2518, 121.8681, 146.0591, 148.8134, 164.6219, 176.522, 177.9578, 180.8773, 187.1846, 210.5131, 211.483, 230.2598, 262.3549, 266.2318, 303.3181, 329.4067, 335.0262, 337.8323, 343.1142, 352.2322, 367.8386, 380.09, 388.5412, 390.4162, 395.6409), Tem = c(15.248, 15.4523, 16.0761, 16.2013, 16.5914, 16.8777, 17.3545, 17.3877, 17.5307, 17.7079, 18.4177, 18.575, 19.8261, 19.9731, 20.4074, 21.2622, 21.4117, 22.1776, 23.4835, 23.6738, 23.9973, 24.4976, 25.7585, 26.0231, 28.5495, 30.8602, 31.3067, 37.3183, 39.2858, 39.4731, 39.6756, 39.9271, 40.6634, 42.3641, 43.9158, 44.1891, 44.3563, 44.5837)), .Names = c("Rt", "Tem"), class = "data.frame", row.names = c(NA, -38L)) library(segmented) # create a linear model out.lm <- lm(Tem ~ Rt, data=bullard) o<-segmented(out.lm, seg.Z=~Rt, psi=list(Rt=c(200,300)), control=seg.control(display=FALSE)) 

Using the psi parameter, I tried the following:

 psi = list(x = c(150, 300)) -- OK psi = list(x = c(100, 200)) -- OK psi = list(x = c(200, 300)) -- OK psi = list(x = c(100, 300)) -- OK psi = list(x = c(120, 150, 300)) -- error 1 below psi = list(x = c(120, 300)) -- OK psi = list(x = c(120, 150)) -- OK psi = list(x = c(150, 300)) -- OK psi = list(x = c(100, 200, 300)) -- error 2 below 

(1) Error in segmented.lm(out.lm, seg.Z = ~Rt, psi = list(Rt = c(120, 150, : only 1 datum in an interval: breakpoint(s) at the boundary or too close

(2) Error in diag(Cov[id, id]) : subscript out of bounds

I have already listed my data on this subject , but as a guide, restrictions on x data are around 0-400.

The second question that relates to this is this: how can I fix breakpoints using this segmented package?

Answer 1

The problem here seems to be bad capture of errors in the segmented package. Looking at the code for segmented.lm , you can debug a little. For example, in the case of psi = list(x = c(100, 200, 300)) an extended linear model is installed, as shown below:

 lm(formula = Tem ~ Rt + U1.Rt + U2.Rt + U3.Rt + psi1.Rt + psi2.Rt + psi3.Rt, data = mf) Call: lm(formula = Tem ~ Rt + U1.Rt + U2.Rt + U3.Rt + psi1.Rt + psi2.Rt + psi3.Rt, data = mf) Coefficients: (Intercept) Rt U1.Rt U2.Rt U3.Rt psi1.Rt 15.34303 0.04149 0.04591 742.74186 -742.74499 1.02252 psi2.Rt psi3.Rt NA NA 

As you can see, the correspondence has NA values, which then lead to a degenerate variance-covariance matrix (called Cov in the code). The function does not check this and tries to pull the diagonal entries from the Cov and does not work with the error message. At least the first error, although perhaps not very useful, falls into the function of the function itself and assumes that the breakpoints are too close.

In the absence of better catching errors in a function, I believe that all you can do is use a trial and error approach (and avoid breakpoints that are too close). For example, psi = list(x = c(50, 200, 300)) works fine.

Answer 2

If you use while and tryCatch , you can tryCatch command until you decide that there is no error in the @jaySf model. I assume that this depends on the randomizer settings in the function, which can be seen in seg.control .

 lm.model <- lm(xdat ~ ydat, data = x) if.false <- F while(if.false == F){ tryCatch({ s <- segmented(lm.model, seg.Z =~ydata, psi = NA) if.false <- T }, error = function(e){ }, finally = {}) }