Polynomial regression prediction works but the formula is wrong

Question

Why am I getting different results when I use stats::predict() versus using the fitted polynomial?

set.seed(999999)

library(tidyverse)

fit <- tibble(rand = rbeta(1e6, 5, 2.5)) %>%
  mutate(rand = cut(rand, 150)) %>%
  count(rand) %>%
  mutate(x = 1:150) %>%
  lm(n ~ poly(x, 10), .)

plot(predict(fit), main = "prediction")

coef(fit)

poly_fn <- polynom::polynomial(coef(fit)) # convert to function

f <- as.function(poly_fn)
plot(f(1:150))

See this question https://stackoverflow.com/questions/66614091/r-how-to-plot-custom-range-of-polynomial-produced-by-lm-poly-fit/66614696#66614696, you need to use the raw=TRUE opinion in the poly function. — Dave2e, Jun 25 '21 at 14:51

score 0 · Answer 1 · answered Jun 25 '21 at 14:51

You need to use poly(..., raw = TRUE) to get the coefficients in the form that polynom::polynomial uses:

set.seed(999999)

library(tidyverse)

fit <- tibble(rand = rbeta(1e6, 5, 2.5)) %>%
  mutate(rand = cut(rand, 150)) %>%
  count(rand) %>%
  mutate(x = 1:150) %>%
  lm(n ~ poly(x, 10, raw = TRUE), .)

plot(predict(fit), main = "prediction")

coef(fit)
#>               (Intercept)  poly(x, 10, raw = TRUE)1  poly(x, 10, raw = TRUE)2 
#>             -1.502706e+00              2.081423e+00             -3.652070e-01 
#>  poly(x, 10, raw = TRUE)3  poly(x, 10, raw = TRUE)4  poly(x, 10, raw = TRUE)5 
#>              5.063552e-02             -1.219081e-03              3.267963e-05 
#>  poly(x, 10, raw = TRUE)6  poly(x, 10, raw = TRUE)7  poly(x, 10, raw = TRUE)8 
#>             -4.348967e-07              2.143856e-09              3.615717e-12 
#>  poly(x, 10, raw = TRUE)9 poly(x, 10, raw = TRUE)10 
#>             -7.321806e-14              1.996748e-16

poly_fn <- polynom::polynomial(coef(fit)) # convert to function

f <- as.function(poly_fn)
plot(f(1:150))

^{Created on 2021-06-25 by the reprex package (v2.0.0)}

Polynomial regression prediction works but the formula is wrong

1 Answers1