[R] how to compute maximum of fitted polynomial?
David Winsemius
dwinsemius at comcast.net
Wed Jun 5 13:32:03 CEST 2013
On Jun 4, 2013, at 10:15 PM, Hans W Borchers wrote:
> Bert Gunter <gunter.berton <at> gene.com> writes:
>>
>> 1. This looks like a homework question. We should not do homework
>> here.
>> 2. optim() will only approximate the max.
>> 3. optim() is not the right numerical tool for this anyway.
>> optimize() is.
>> 4. There is never a guarantee numerical methods will find the max.
>> 5. This can (and should?) be done exactly using elementary math
>> rather
>> than numerical methods.
>>
>> Cheers,
>> Bert
>
> In the case of polynomials, "elementary math ... methods" can
> actually be
> executed with R:
>
> library(polynomial) # -6 + 11*x - 6*x^2 + x^3
> p0 <- polynomial(c(-6, 11, -6, 1)) # has zeros at 1, 2, and 3
> p1 <- deriv(p0); p2 <- deriv(p1) # first and second derivative
> xm <- solve(p1) # maxima and minima of p0
> xmax = xm[predict(p2, xm) < 0] # select the maxima
> xmax # [1] 1.42265
>
> Obviously, the same procedure will work for polynomials p0 of higher
> orders.
These look like the functions present in the 'polynom' package
authored by Bill Venables [aut] (S original), Kurt Hornik [aut, cre]
(R port), Martin Maechler. I wasn't able to find a 'polynomial'
package on CRAN. The 'mpoly' package by David Kahle offers
multivariate symbolic operations as well.
--
David.
>
> Hans Werner
>
>
>>> Em 04-06-2013 21:32, Joseph Clark escreveu:
>>>>
>>>> My script fits a third-order polynomial to my data with something
>>>> like
>>>> this:
>>>>
>>>> model <- lm( y ~ poly(x, 3) )
>>>>
>>>> What I'd like to do is find the theoretical maximum of the
>>>> polynomial
>>>> (i.e. the x at which "model" predicts the highest y).
>>>> Specifically, I'd
>>>> like to predict the maximum between 0 <= x <= 1.
>>>>
>>>> What's the best way to accomplish that in R?
>>>>
>>>> Bonus question: can R give me the derivative or 2nd derivative of
>>>> the
>>>> polynomial? I'd like to be able to compute these at that maximum
>>>> point.
>>>>
>>>> Thanks in advance!
--
David Winsemius, MD
Alameda, CA, USA
More information about the R-help
mailing list