[R] Polynomial fitting

Jon Minton 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
etc.

In R this is easy, just use the identity operator I() when specifying the
equation.
e.g. for a 3rd order polynomial:

model <- lm(Y ~ x + I(x^2) + I(x^3) + I(x^4))

hth, Jon

***********

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[1]*x^n + c[2]*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)[2] +
coefficients(model)[1] );
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!

Sincerely,

-- Marco




Checked by AVG Free Edition. 

16:55



More information about the R-help mailing list