I am interested in comparing three models
- a linear regression model: mod_linear <- lm (dv ~ iv)
- a polynomial regression model: mod_polynomial <- lm (dv ~ iv + I(iv^2))
- a one-breakpoint piecewise/segmented linear model: mod_segmented <- segmented(lm(dv~iv), seg = ~iv, psi=NA)
iv=independent variable, dv= dependent variable
Note for the segmented model: It estimates a breakpoint/changepoint/shift/knot in the slope, i.e., two lines and one inflection point in which the slope bends. This inflection point is itself a parameter (psi); it is not pre-specified, but is estimated via a robust fitting algorithm (built into the segmented package in R).
My question is: How do I compare the models to each other?
A. Is there a formal test that I can use? I know that nested models can be compared using the anova() function for a chi-square test in R, but I am not sure if it is appropriate to use anova() to compare the polynomial and piecewise models to each other.
B. Is it appropriate to use the R-squared values, residual standard error, AIC, and BIC, to compare the models?
Note: I am basing my reasoning for (A) on a paper where they compare piecewise and polynomial fits, and they state that "there is no formal test available that allows the direct comparison of a piecewise model to a power[polynomial] model" (pg 10, Ryan & Porth, 2007)
Ryan, S. E. & Porth, L. S. (2007). A tutorial on the piecewise regression approach applied to bedload transport data. U.S. Department of Agriculture & Forest Service Rocky Mountain Research Station, Fort Collins, CO, United States. https://www.arlis.org/docs/vol1/173845244.pdf
Below are fitted lines from each model polynomial=blue, segmented/piecewise=red, data points=gray
Reproducible code and sample data
## load packages required
library(segmented)
library(lme4)
## Sample data
iv <-c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,97,98,99,100,101,102,103,104,105,106,107,108,109,110,111,112,113,114,115,116,117,118,119,120,121,122,123,124,125,126,127,128,129,130,131,132,133,134,135,136,137,138,139,140,141,142,143,144,145,146,147,148,149,150,151,152,153,154,155,156,157,158,159,160,161,162,163,164,165,166,167,168,169,170,171,172,173,174,175,176,177,178,179,180,181,182,183,184,185,186,187,188,189,190,191,192,193,194,195,196,197,198,199,200,201,202,203,204,205,206,207,208,209,210,211,212,213,214,215,216,217,218,219,220,221,222,223,224,225,226,227,228,229,230,231,232,233,234,235,236,237,238,239,240,241,242,243,244,245,246,247,248,249,250,251,252,253,254,255,256,257,258,259,260,261,262,263,264,265,266,267,268,269,270,271,272,273,274,275,276,277,278,279,280,281,282,283,284,285,286,287,288,289,290,291,292,293,294,295,296,297,298,299,300,301,302,303,304,305,306,307,308,309,310,311,312)
dv <-c(0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8532,0.8916,0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.9372,0.9,0.9,0.9492,0.9492,0.9492,0.9492,0.9492,0.9984,0.9984,0.984,0.9984,0.984,0.9984,0.9984,0.9984,1.05,1.05,1.05,1.05,1.0812,1.0812,1.0812,1.0812,1.05,1.05,1.1016,1.1016,1.1016,1.1016,1.1016,1.1316,1.1316,1.1832,1.1832,1.1832,1.1832,1.1832,1.1544,1.1544,1.1544,1.2072,1.2072,1.2072,1.2072,1.236,1.2888,1.2888,1.2888,1.2888,1.236,1.2612,1.2612,1.2612,1.2612,1.2612,1.3164,1.3164,1.3416,1.3416,1.3416,1.3416,1.3956,1.3416,1.3704,1.3704,1.3704,1.3704,1.3704,1.4256,1.4256,1.4256,1.4256,1.4496,1.4496,1.4496,1.4496,1.482,1.482,1.482,1.482,1.482,1.482,1.482,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5384,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.5948,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.6512,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7076,1.7652,1.728,1.728,1.8408,1.7844,1.7844,1.7844,1.7652,1.8216,1.8216,1.8216,1.8792,1.8792,1.8792,1.8792,1.8792,1.9212,1.9368,1.9368,1.9788,1.9944,1.9788,1.9944,1.9944,2.052,2.052,2.1108,2.1108,2.1108,2.1684,2.1684,2.1684,2.1684,2.226,2.226,2.226,2.2848,2.2848,2.2848,2.2848,2.3436,2.3436,2.3436,2.4012,2.4012,2.4012,2.46,2.46,2.46,2.46,2.5188,2.5188,2.5188,2.5776,2.5776,2.5776,2.6244,2.6364,2.6832,2.6832,2.6832,2.6952,2.742,2.742,2.802,2.802,2.802,2.8608,2.8608,2.9196,2.9196,2.9196,2.9784,2.9784,2.9784,3.0384,3.0384,3.0972,3.0972,3.0972,3.0972,3.1572,3.1572,3.1572,3.216,3.216,3.2748,3.2676,3.2748,3.2748,3.3348,3.3264,3.3864,3.3864,3.4452,3.4452,3.5052,3.5052,3.5052,3.5652,3.624,3.624,3.624,3.684,3.6912,3.6912,3.744,3.8028,3.8028,3.87,3.87,3.8628,3.9228,3.9228,3.9828,3.9888,4.0416,4.056,4.1016,4.1076,4.1616,4.1676,4.2204,4.2276,4.2864,4.2936,4.3464,4.3464,4.4388,4.4064,4.4664,4.4664,4.4724,4.5324,4.5672,4.5996,4.626,4.7352,4.7352,4.7184,4.7784,4.794,4.8456,4.8372,4.8972,4.9728,5.016,5.0232,5.0832,5.142,5.0916,5.1348,5.1948,5.202,5.2692,5.2548,5.3292,5.3136,5.3736,5.3736,5.4336,5.4336,5.4996,5.4336,5.5068,5.4924,5.5752,5.6856,5.694,5.6184,5.7384,5.7972,5.7972,5.7972,5.8572,5.8716,5.9316,5.9172)
## Fit models
mod_linear <- lm (dv ~ iv)
mod_polynomial <- lm (dv ~ iv + I(iv^2))
mod_segmented <- segmented(mod_linear, seg.Z = ~iv, psi=NA,
it=20, control=seg.control(n.boot=50, size.boot=1000, K=1))
## Model comparison
# Chi-square test
anova(mod_linear, mod_polynomial, mod_segmented)
# Parsimony criteria: AIC and BIC
AIC(mod_linear, mod_polynomial, mod_segmented)
BIC (mod_linear, mod_polynomial, mod_segmented)
## Plot fitted models
plot(dv ~ iv , col="grey")
lines(predict(mod_polynomial), col="blue", lwd=3)
plot(mod_segmented, add=TRUE, col="red", lwd=3)
Results from model comparison
> # Chi-square test
> anova(mod_linear, mod_polynomial, mod_segmented)
Analysis of Variance Table
Model 1: dv ~ iv
Model 2: dv ~ iv + I(iv^2)
Model 3: dv ~ iv + U1.iv + psi1.iv
Res.Df RSS Df Sum of Sq F Pr(>F)
1 310 57.377
2 309 2.307 1 55.070 8338.356 < 2.2e-16 ***
3 308 2.034 1 0.273 41.334 4.883e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
> # Parsimony criteria: AIC and BIC
> AIC(mod_linear, mod_polynomial, mod_segmented)
df AIC
mod_linear 3 363.0917
mod_polynomial 4 -637.5626
mod_segmented 5 -674.8518
> BIC (mod_linear, mod_polynomial, mod_segmented)
df BIC
mod_linear 3 374.3207
mod_polynomial 4 -622.5906
mod_segmented 5 -656.1368
```