# Does autocorrelation cause bias in the regression parameters in piecewise regression?

In simple linear regression problems, autocorrelated residuals are supposed not to result in biased estimates for the regression parameters. Can the same be said for piecewise regression?

Suppose I want to fit a continuous, piecewise linear function of a single variable. Let's say for example we have data on shipping cost and weight of shipment. The function is piecewise because as the weight increases, at some point an additional rail car is required. We want to find the breakpoints and the slopes of the individual pieces. The model is fit, and for whatever reason, the residuals are found to be serially correlated in time. Could the regression parameters be biased?

Suppose it is known that there are two breakpoints at (unknown) points x1 and x2. We want to fit the data to a model f(x) given by:

x < x1:      f(x) = a + m1*x
x1 < x < x2: f(x) = a + m1*x1 + m2*(x - x1)
x > x2:      f(x) = a + m1*x1 + m2*(x2 - x1) + m3*(x - x2)


I use the nlm function in R to find the unknown parameters x1, x2, m1, m2 and m3:

sqerr <- function(prm,y,x) {
a <- prm[1]
x1 <- prm[2]
x2 <- prm[3]
m1 <- prm[4]
m2 <- prm[5]
m3 <- prm[6]
sqerr <- sum((y-(a+ifelse(x<x1,m1*x,
m1*x1+ifelse(x<x2,m2*(x-x1),
m2*(x2-x1)+m3*(x-x2)))))^2)
}
ai <- 0.4; x1i <- 0.4; x2i <- 0.7; m1i <- 0.0; m2i <- 0.8; m3i <- 3
prm <- c(ai,x1i,x2i,m1i,m2i,m3i)
uu <- nlm(sqerr,prm,data$Y,data$X)


Then I plot the residuals vs. the lag-1 residuals:

y <- data$Y x <- data$X
a <- uu$est[1] x1 <- uu$est[2]
x2 <- uu$est[3] m1 <- uu$est[4]
m2 <- uu$est[5] m3 <- uu$est[6]
resid <- (y-(a+ifelse(x<x1,m1*x,m1*x1+ifelse(x<x2,m2*(x-x1),m2*(x2-x1)+m3*(x-x2)))))
plot(resid[1:149]~resid[2:150])


There is clearly some sequential correlation. So my question is, are the regression parameters biased because of this? I have an old paper by Kadiyala (A Transformation Used to Circumvent the Problem of Autocorrelation, Econometrica Vol. 36, No. 1, Jan. 1968) that states:

"It is well known (see Watson [7] and Watson and Hannan [8]) that simple least squares estimators, though unbiased (when the independent variables are "fixed variates"),are, in general inefficient in the presence of autocorrelation among the disturbances."

It seems that by "simple least squares" he means linear equations of the form y = a + bx (that is the example used in the paper). But I have seen other papers that seem to imply that the estimators (i.e., regression parameters) are unbiased no matter what type of model you have. I don't think it's true in general.

• What precisely do you mean by "stepwise regression"? Is this a synonym for a model w/ a piecewise covariate, or are you referring to stepwise selection methods to determine the covariates that belong in the final model? – gung - Reinstate Monica May 26 '12 at 14:08
• Yes, I meant piecewise -- I edited the title. – Ringold May 26 '12 at 14:41
• Please post your data as we have seen this kind of problem resolved with causal input series and empirically identified level shifts and/or local time trends. A cost accountanting firm had retained us to analyze such data in order to estimate "fixed costs" and "variable costs". I will try and help. – IrishStat May 26 '12 at 16:40
• See my edits above to the original question – Ringold May 26 '12 at 22:17

A regression parameter that is often forgotten is the variance of the residuals. This one will be biased if residuals are correlated. This means that p-values of whatever test you are performing have to be handled with great care.

Otherwise, if you fit a single line through something that is not linear (your case), you should observe auto-correlation of the residuals, but through the X variable, not through time. In that case the parameters are not biased, they are just wrong.

However you specifically mention that your residuals are auto-correlated in time, so you could perhaps add time as a variable in your model and check whether this decorrelates the residuals.

• adding "time" is a very bad practice as there is an implicit assumption of 1 and only 1 trend that starts at the beginning . If there is no "trend" for the first k periods , the correct series to use would be 0,0,0,0,,,,k,1,2,3,4,5..... additionally there may be level shifts (intercept changes) that shouldn't be confused with trends. – IrishStat May 26 '12 at 16:36
• (+1) Thanks @IrishStat, I was not aware of this. I am not sure I get the point though. Could you give some link or reference? It seems to me that this assumption is also implicitly made for every variable. – gui11aume May 26 '12 at 16:43
• You are right there is an assumption of linearity. What you are doing could be referred to as Model Specification Bias. Consider the following series 1,2,3,4,7,8,9,10 ... this series has a level shift (intercept change) in addition to a trend. Consider,1,2,3,4,5,14,16,18,20 which has both an intercept change and a slope change. Hope this helps – IrishStat May 26 '12 at 19:33

Thanks for sharing your data. It raises some interesting answers. To begin with a potentially useful model between y and x is which suggests a strong relationship between y and two previous y's and both a contemporaneous and lag 1 effect of X. The plot of actual/fit and forecast is and the cleansed ( outlier adjusted series) is does not suggest level shifts and/or local time trends but rather a few one-time anomalies. The ACF of the original series is while the ACF of the model residuals is . In conclusion there is no need for "local splines", local trends , level shifts in the presence of the x variable which carries the load for the visually suggestive non-conditional plot of y .Now if we ignore the x variable and any memory in y (the ARIMA structure ) and only focus on detecting and incorporating any needed pulses, level shifts, seasonal pulses and.or local time trends we get a totally different answer. Here is the equation showing two time trends and two level shifts and a few pulses reflecting unknown exceptional activity. The actual/fit/forecast is . The acf of the residuals suggests some omitted ( by design/specification of no ARIMA ) memory structure ! and a plot of the residuals . The fit/fore graphic tells the story of the equation in a visual way . As I commented to my friend Michael, the data suggests appropriate remedies. In summary in the absence of x and the memory of y there are a few "local splines" whose length and time range can be found analytically. If x is included and the past of y is considered then there is no additional need for these "local time trends' I hope this has been of help to the original poster and to others who have commented here. I notice that I misnamed the data ringwald rather than ringold. For transparency I am one of the developers of AUTOBOX , the software that I used for this analysis. The methods used are based upon the pioneering work of G.C.Tiao and others.

• +1 Nice example of an appropriate analysis with the OP's data. – whuber May 28 '12 at 19:50

I think piecewise regression means fitting several different lines at various cut points. It is not clear whether the number of cutoffs is prespecified and whether their locations are prespecified. Even if they are all prespecified it seems that each piece would be fit by ordinary regression and the problem of correlationed residuals and under or overestimated residual variance would be there in each line. There is also another assumption not stated. Are the residuals for each line assumed to have the same variance as for all the others? So the problem exists and may be worse in this more complex type of regression. Regarding IrishStat's comment. I think he is right for the piecewise model because the breakpoints could be time interventions that affect the stationary component and possibly also the nonstationary component of the model. But in ordinary harmonic regression the nonstationary seasonal component is modelled by sine and cosine functions of t and the residuals could be white noises or they could be model as a stationary process in time such as an AR(p) model.

• Yes, the number of cut points is specified (but not their location -- those are parameters), and it can be assumed that the residuals all have the same variance. – Ringold May 26 '12 at 17:35
• @Ringold that is why I wanted to look at your data to possibly identify the "cut points" as Michael alluded to. – IrishStat May 26 '12 at 18:33
• @Ringold I believe that IrishStat's time series software both identifies the number and the location of time points where a time series changes (an intervention occurs). If you give him your data he can show you how he would identify the points of intervantion and model the time series in the intervals between interventions. This would be a little different from your piecewise linear regerssion because he probably would fit an ARIMA model between interventions rather than a linear function of time. – Michael R. Chernick May 26 '12 at 19:04
• @Michael Not necessarily so . The data, properly analyzed might suggest the appropriate remedy . – IrishStat May 26 '12 at 19:37
• @IrishStat I was assuming that your software could identify the breakpoints but does not estimate linear functions of time. So you would solve part of the problem for him. Also and Arima fir in each interval could be better than the regression lines. That would depend on the data. Am I right? – Michael R. Chernick May 26 '12 at 20:18