In segmented regression, how to choose so the second line cant have a negative slope?

Viewed 13

I am using the R package segmented to detrend my data and to obtain the residuals. I want to fit 2 linear lines with one breakpoint, for each of the groups of cyl, but I don't want the second line to have a negative slope (it should be positive or 0). How can I do that?

This is my code to create 2 segmented lines for each cyl:

library(segmented)

    mtcars%>% 
    nest_by(cyl) %>%
    mutate(mod = list(segmented(lm(hp~disp, data = data)))) %>% 
    summarize(augment(mod))
1 Answers

You could try setting the priors using the mcp package. In the summary tables, cp is the change point and hp_* are the slopes. Note that I set the priors so that the second slope is always positive.

library(tidyverse)
library(mcp)


prior = list(
  disp_2 = "dnorm(0, 10) T(0,)"
)

model = list(
  hp ~ disp,          
     ~ 0 + disp  
)

models_by_cyl <- mtcars |>
  nest(data = -cyl) |>
  mutate(fit = map(data, ~mcp(model,data = .x, prior = prior)))



models_by_cyl$fit
#> [[1]]
#> Family: gaussian(link = 'identity')
#> Iterations: 9000 from 3 chains.
#> Segments:
#>   1: hp ~ disp
#>   2: hp ~ 1 ~ 0 + disp
#> 
#> Population-level parameters:
#>     name  mean    lower  upper Rhat n.eff
#>     cp_1 237.5  1.8e+02 258.00    1   742
#>   disp_1  -0.1 -5.8e-01   0.36    1   176
#>   disp_2   3.2  4.9e-05  12.25    1  1889
#>    int_1 136.4  5.2e+01 226.81    1   181
#>  sigma_1  27.3  1.3e+01  43.89    1  1081
#> 
#> [[2]]
#> Family: gaussian(link = 'identity')
#> Iterations: 9000 from 3 chains.
#> Segments:
#>   1: hp ~ disp
#>   2: hp ~ 1 ~ 0 + disp
#> 
#> Population-level parameters:
#>     name   mean    lower upper Rhat n.eff
#>     cp_1 115.67  7.5e+01 146.7  1.0   991
#>   disp_1   0.35 -2.3e-01   1.1  1.1   156
#>   disp_2   1.55  1.8e-06   7.6  1.0  1649
#>    int_1  45.41 -1.8e+01 102.8  1.1   169
#>  sigma_1  21.39  1.2e+01  32.3  1.0  1937
#> 
#> [[3]]
#> Family: gaussian(link = 'identity')
#> Iterations: 9000 from 3 chains.
#> Segments:
#>   1: hp ~ disp
#>   2: hp ~ 1 ~ 0 + disp
#> 
#> Population-level parameters:
#>     name  mean    lower  upper Rhat n.eff
#>     cp_1 418.3  3.0e+02 472.00    1  1120
#>   disp_1   0.1 -2.7e-01   0.51    1   146
#>   disp_2   2.4  9.5e-05  10.24    1  2168
#>    int_1 168.5  3.2e+01 295.27    1   143
#>  sigma_1  55.6  3.6e+01  78.32    1  3593
Related