R – Calculating Standard Error of Slopes in Piecewise Linear Regression with Known Breakpoints

piecewise linearrregressionstandard error

The situation

I have a dataset with one dependent $y$ and one independent variable $x$. I want to fit a continuous piecewise linear regression with $k$ known/fixed breakpoints occurring at $(a_{1}, a_{2}, \ldots, a_{k})$. The breakpoins are known without uncertainty, so I don't want to estimate them. Then I fit a regression (OLS) of the form
$$
y_{i} = \beta_{0} + \beta_{1}x_{i} + \beta_{2}\operatorname{max}(x_{i}-a_{1},0) + \beta_{3}\operatorname{max}(x_{i}-a_{2},0) +\ldots+ \beta_{k+1}\operatorname{max}(x_{i}-a_{k},0) +\epsilon_{i}
$$
Here is an example in R

set.seed(123)
x <- c(1:10, 13:22)
y <- numeric(20)
y[1:10] <- 20:11 + rnorm(10, 0, 1.5)
y[11:20] <- seq(11, 15, len=10) + rnorm(10, 0, 2)

Let's assume that the breakpoint $k_1$ occurs at $9.6$:

mod <- lm(y~x+I(pmax(x-9.6, 0)))
summary(mod)

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)          21.7057     1.1726  18.511 1.06e-12 ***
x                    -1.1003     0.1788  -6.155 1.06e-05 ***
I(pmax(x - 9.6, 0))   1.3760     0.2688   5.120 8.54e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

The intercept and slope of the two segments are: $21.7$ and $-1.1$ for the first and $8.5$ and $0.27$ for the second, respectively.

Breakpoint

Questions

  1. How to easily calculate the intercept and slope of each segment? Can the model be reparemetrized to do this in one calculation?
  2. How to calculate the standard error of each slope of each segment?
  3. How to test whether two adjacent slopes have the same slopes (i.e. whether the breakpoint can be omitted)?

Best Answer

  1. How to easily calculate the intercept and slope of each segment?

The slope of each segment is calculated by simply adding all the coefficients up to the current position. So the slope estimate at $x=15$ is $-1.1003 + 1.3760 = 0.2757\,$.

The intercept is a little harder, but it's a linear combination of coefficients (involving the knots).

In your example, the second line meets the first at $x=9.6$, so the red point is on the first line at $21.7057 -1.1003 \times 9.6 = 11.1428$. Since the second line passes through the point $(9.6, 11.428)$ with slope $0.2757$, its intercept is $11.1428 - 0.2757 \times 9.6 = 8.496$. Of course, you can put those steps together and it simplifies right down to the intercept for the second segment = $\beta_0 - \beta_2 k_1 = 21.7057 - 1.3760 \times 9.6$.

Can the model be reparameterized to do this in one calculation?

Well, yes, but it's probably easier in general to just compute it from the model.

2. How to calculate the standard error of each slope of each segment?

Since the estimate is a linear combination of regression coefficients $a^\top\hat\beta$, where $a$ consists of 1's and 0s, the variance is $a^\top\text{Var}(\hat\beta)a$. The standard error is the square root of that sum of variance and covariance terms.

e.g. in your example, the standard error of the slope of the second segment is:

Sb <- vcov(mod)[2:3,2:3]
sqrt(sum(Sb))

alternatively in matrix form:

Sb <- vcov(mod)
a <- matrix(c(0,1,1),nr=3)
sqrt(t(a) %*% Sb %*% a)

3. How to test whether two adjacent slopes have the same slopes (i.e. whether the breakpoint can be omitted)?

This is tested by looking at the coefficient in the table of that segment. See this line:

I(pmax(x - 9.6, 0))   1.3760     0.2688   5.120 8.54e-05 ***

That's the change in slope at 9.6. If that change is different from 0, the two slopes aren't the same. So the p-value for a test that the second segment has the same slope as the first is right at the end of that line.

Related Question