[R] Polynomial fitting
jm540 at york.ac.uk
Thu Aug 16 12:40:01 CEST 2007
Remember that polynomials of the form
y = b1*x + b2*x^2 + ... + bm*x^m
fit the linear regression equation form
Y = beta_1*x_1 + beta_2*x_2 + ... + beta_m*x_m
If one sets (from the 1st to the 2nd equation)
x -> x_1
x^2 -> x_2
x^3 -> x_3
In R this is easy, just use the identity operator I() when specifying the
e.g. for a 3rd order polynomial:
model <- lm(Y ~ x + I(x^2) + I(x^3) + I(x^4))
I'm looking some way to do in R a polynomial fit, say like polyfit
function of Octave/MATLAB.
For who don't know, c = polyfit(x,y,m) finds the coefficients of a
polynomial p(x) of degree m that fits the data, p(x[i]) to y[i], in a
least squares sense. The result c is a vector of length m+1 containing
the polynomial coefficients in descending powers:
p(x) = c*x^n + c*x^(n-1) + ... + c[n]*x + c[n+1]
For prediction, one can then use function polyval like the following:
y0 = polyval( polyfit( x, y, degree ), x0 )
y0 are the prediction values at points x0 using the given polynomial.
In R, we know there is lm for 1-degree polynomial:
lm( y ~ x ) == polyfit( x, y, 1 )
and for prediction I can just create a function like:
lsqfit <- function( model, xx ) return( xx * coefficients(model) +
and then: y0 <- lsqfit(x0)
(I've tried with predict.lm( model, newdata=x0 ) but obtain a bad result)
For a degree greater than 1, say m, what can I use.??
I've tried with
lm( y ~ poly(x, degree=m) )
I've also looked at glm, nlm, approx, ... but with these I can't
specify the polynomial degree.
Thank you so much!
Checked by AVG Free Edition.
More information about the R-help