CODEWITHSUNDEEP
rdplyrlinear-regressionsegmentpiecewise
Daniel Valencia C.
Daniel Valencia C.

Reputation: 2279

How to use segmented package when working with data frames with dplyr package to perform piecewise linear regression?

My data frame is separated by groups. I want to perform piecewise linear regression on each group and for that I intend to use the segmented package.

First I created the linear models for each group using the dplyr package. The next step is to segment these models, however this is where I'm stuck. Any tips or other way to do this? The ultimate goal is to use these segments to make a graph.

library(dplyr)
library(segmented)

Group <- c("A", "B")
x <- 0:10
y <- c(0, 0.4, 0.6, 0.8, 0.9, 0.9, 0.95, 0.97, 0.98, 0.99, 1,
       0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1)

df <- expand.grid(x = x,
                  Group = Group)

df$y <- y

Segment <- df %>%
  group_by(Group) %>%
  do(my.lm = lm(data = .,
                formula = y ~ x))

Unsuccessful attempts:

Segment <- df %>%
  group_by(Group) %>%
  do(my.lm = lm(data = .,
                formula = y ~ x),
     my.seg = segmented(my.lm,
                        seg.Z = x))


Segment <- df %>%
  group_by(Group) %>%
  do(my.lm = lm(data = .,
                formula = y ~ x)) %>%
  do(my.seg = segmented(my.lm,
                        seg.Z = x))


Segment <- df %>%
  group_by(Group) %>%
  do(my.lm = lm(data = .,
                formula = y ~ x)) %>%
  mutate(my.seg = segmented(my.lm,
                        seg.Z = x))

Upvotes: 1

Views: 932

Answers (1)

akrun
akrun

Reputation: 887891

One option is to wrap with tryCatch and return a NA for possible errors

library(dplyr)
out <- df %>% 
    nest_by(Group) %>%
    mutate(my.lm = list(lm(y ~ x, data = data)),
        my.seg = list(tryCatch(segmented(my.lm, seg.Z = ~ x),
         error = function(e) list(NA))))

-output

> out
# A tibble: 2 x 4
# Rowwise:  Group
  Group               data my.lm  my.seg    
  <fct> <list<tibble[,2]>> <list> <list>    
1 A               [11 × 2] <lm>   <segmentd>
2 B               [11 × 2] <lm>   <list [1]>
> out$my.seg
[[1]]
Call: segmented.lm(obj = my.lm, seg.Z = ~x)

Meaningful coefficients of the linear terms:
(Intercept)            x         U1.x  
    0.03333      0.30000     -0.27488  

Estimated Break-Point(s):
psi1.x  
 2.691  

[[2]]
[[2]][[1]]
[1] NA

Upvotes: 1

Related Questions