2
$\begingroup$

I know there are many arguments against intercept-free regression but we are very interested in the convex or concave shape of these lines. Because I got suspicious real results, I switched my thinking hat from EDA to SW testing and found clearly wrong results with straight lines from degenerated ns, bs, and for default poly even worse, but people tell me it is my misunderstanding of intercept and orthogonal and not a bug. (I know poly with raw=TRUE does not have this problem)

wrong bs (blue) and poly (red) regression lines

x <- c(0.2,1,2,3) # minor change in first point
y <- c(0.2,1,2,3)
fs <- lm(y ~ 0+splines::bs(x, degree=1))
fp <- lm(y ~ 0+poly(x, 1))
plot(x,y)
lines(x,fitted(fs), col="blue")
lines(x,fitted(fp), col="red")
summary(fs)
summary(fp)
sessionInfo()

This is the session info:

R version 4.1.1 (2021-08-10) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 19043)

Matrix products: default

locale: 1 LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C
[5] LC_TIME=English_United States.1252

attached base packages: 1 splines stats graphics grDevices utils datasets methods base

other attached packages: 1 onls_0.1-1 minpack.lm_1.2-1

loaded via a namespace (and not attached): 1 compiler_4.1.1 tools_4.1.1

$\endgroup$
4
  • $\begingroup$ Yes, it is your misunderstanding of orthogonal polynomials. The poly function constructs them with the assumption that an intercept is present in the model. It doesn't know that there is no intercept and you can't specify this. $\endgroup$
    – Roland
    Nov 1, 2021 at 9:01
  • $\begingroup$ I should have left out the poly example, but you say the default raw=F can not be used for an intercept-free regression. But please comment on the spline. $\endgroup$ Nov 1, 2021 at 11:22
  • $\begingroup$ I don't have a deep understanding of these, but it appears you must specify 0 as a boundary knot to ensure that the spline passes through zero: fs1 <- lm(y ~ 0+splines::bs(x, degree=1, Boundary.knots = c(0, 3))) $\endgroup$
    – Roland
    Nov 1, 2021 at 11:47
  • $\begingroup$ Yes, thank you very much, specifying boundary.knots indeed is the solution. And I also intuitively understand this. $\endgroup$ Nov 2, 2021 at 16:34

2 Answers 2

1
$\begingroup$

The poly function does not place the intercept at $x=0$.

It creates a first order polynomial $a + bx$ that is orthogonal to the zero order polynomial (the constant value $1$) and has a sum of squares equal to 1.

So you are fitting with that $a+bx$ instead of with $x$.

Check it out like this

x = c(0.2,1,2,3)

### sum equal to 0
> sum(poly(x, order = 1))
[1] 2.775558e-17
### sum of squares equal to 1
> sum(poly(x, order = 1)^2)
[1] 1

### So these are the values that you use to fit
> poly(x, order = 1)
              1
[1,] -0.6414044
[2,] -0.2613129
[3,]  0.2138015
[4,]  0.6889158  

See also https://en.wikipedia.org/wiki/Discrete_Chebyshev_polynomials

$\endgroup$
1
$\begingroup$

I admit it was my misunderstanding of orthogonal polynomials. Now I see the default, orthogonal, poly by definition can not be intercept-free so specifying 0+ results only in erroneously leaving out the necessary intercept term from the model.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.