[R] solving x in a polynomial function

Peter Ehlers ehlers at ucalgary.ca
Fri Mar 1 22:30:25 CET 2013


On 2013-03-01 13:06, Mike Rennie wrote:
> Hi guys,
>
> Perhaps on the right track? However, not sure if it's correct. I fixed
> the bug that A.K. indicated above (== vs =), but the values don't seem
> right. From the original data, an a of 3 should give b of 2.5.
>
>> realroots <- function(model, b){
> +         is.zero <- function(x, tol = .Machine$double.eps^0.5) abs(x) < tol
> +         if(names(model)[1] == "(Intercept)")
> +                 r <- polyroot(c(coef(model)[1] - b, coef(model)[-1]))
> +         else
> +                 r <- polyroot(c(-b, coef(model)))
> +         Re(r[is.zero(Im(r))])
> + }
>>
>> r <- realroots(po.lm, 3)
>> predict(po.lm, newdata = data.frame(b = r)) # confirm
>     1
> 1.69
>
> So I think there's a calculation error somehwere.

You need to replace the following line

   if(names(model)[1] == "(Intercept)")

with

   if(names(coef(model))[1] == "(Intercept)")


Peter Ehlers


>
>
> On 3/1/13, arun <smartpink111 at yahoo.com> wrote:
>> Hi Rui,
>>
>> Looks like a bug:
>> ###if(names(model)[1] = "(Intercept)")
>> Should this be:
>>
>> if(names(model)[1] == "(Intercept)")
>>
>> A.K.
>>
>>
>>
>> ----- Original Message -----
>> From: Rui Barradas <ruipbarradas at sapo.pt>
>> To: Mike Rennie <mikerennie.r at gmail.com>
>> Cc: r-help Mailing List <r-help at r-project.org>
>> Sent: Friday, March 1, 2013 3:18 PM
>> Subject: Re: [R] solving x in a polynomial function
>>
>> Hello,
>>
>> Try the following.
>>
>>
>> a <- 1:10
>> b <- c(1, 2, 2.5, 3, 3.5, 4, 6, 7, 7.5, 8)
>>
>> dat <- data.frame(a = a, b = b)  # for lm(), it's better to use a df
>> po.lm <- lm(a~b+I(b^2)+I(b^3)+I(b^4), data = dat); summary(po.lm)
>>
>>
>> realroots <- function(model, b){
>>      is.zero <- function(x, tol = .Machine$double.eps^0.5) abs(x) < tol
>>      if(names(model)[1] = "(Intercept)")
>>          r <- polyroot(c(coef(model)[1] - b, coef(model)[-1]))
>>      else
>>          r <- polyroot(c(-b, coef(model)))
>>      Re(r[is.zero(Im(r))])
>> }
>>
>> r <- realroots(po.lm, 5.5)
>> predict(po.lm, newdata = data.frame(b = r))  # confirm
>>
>>
>>
>> Hope this helps,
>>
>> Rui Barradas
>>
>> Em 01-03-2013 18:47, Mike Rennie escreveu:
>>> Hi there,
>>>
>>> Does anyone know how I solve for x from a given y in a polynomial
>>> function? Here's some example code:
>>>
>>> ##example file
>>>
>>> a<-1:10
>>>
>>> b<-c(1,2,2.5,3,3.5,4,6,7,7.5,8)
>>>
>>> po.lm<-lm(a~b+I(b^2)+I(b^3)+I(b^4)); summary(po.lm)
>>>
>>> (please ignore that the model is severely overfit- that's not the point).
>>>
>>> Let's say I want to solve for the value b where a = 5.5.
>>>
>>> Any thoughts? I did come across the polynom package, but I don't think
>>> that does it- I suspect the answer is simpler than I am making it out
>>> to be. Any help would be welcome.
>>>
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
>



More information about the R-help mailing list