How to represent polynomials with numeric vectors in R

Question

In R how would one represent polynomial expressions and do polynomial math with the numeric vector objects? For example:

x1 <- c(2,1)  # 2 + x
x2 <- c(-1,3)  # -1 + 3*x

And want:

x1 * x2 # to return -2 + 5*x + 3*x^2 

Note: I answered a question this morning and then the poster apparently deleted it (making me wonder if it were homework.) So I am re-posting the question from memory.

Luisaluise

Try the mpoly package? – Bondie

That package does look more full-featured than 'polynom' in that it handles multivariate operations. I'm not intending to checkmark my own answer, so I invite better answers than the one I provided. – Luisaluise

Answer 1

One could multiply the coefficients directly using outer and then aggregate the results

x1 <- c(2,1)  # 2 + x
x2 <- c(-1,3)  # -1 + 3*x
tmp <- outer(x1, x2)
tapply(tmp, row(tmp) + col(tmp) - 1, sum)
# 1  2  3 
#-2  5  3

x1 <- c(2, 1) # 2 + x
x2 <- c(-1, 3, 2) # -1 + 3*x + 2*x^2
tmp <- outer(x1, x2)
tapply(tmp, row(tmp) + col(tmp) - 1, sum) # should give -2 + 5*x + 7*x^2 + 2*x^3
# 1  2  3  4 
#-2  5  7  2

as discussed in the comments the '-1' in the code isn't necessary. When coming up with the solution that helped me because it allowed me to map each location in the output of outer to where it would end up in the final vector. If we did a '-2' instead then it would map to the exponent on x in the resulting polynomial. But we really don't need it so something like the following would work just as well:

tmp <- outer(x1, x2)
tapply(tmp, row(tmp) + col(tmp), sum)

Peeve

That's kewl. Kind of a binomial expansion triangle tilted 45 degrees counterclockwise. – Luisaluise

This is very similar to the code in polynom:::Ops.polynomial, although they use tapply(m, row(m) + col(m), sum)... Why do both give the same result? – Ataraxia

I actually don't need the `-1 in the code but it helped me mentally think about the problem at first. I should probably just get rid of it. Those numbers are just being used to group together results so we could add/subtract any amount and it would stay the same. – Peeve

If you subtracted 2 instead of 1 you could use the sum as the exponent for x. – Luisaluise

Yeah - I was using it as the index of the resulting vector (thus only subtracting 1 instead of 2) - but either way it's unnecessary. – Peeve

Answer 2

Use the polynom package:

 require(polynom)
# Loading required package: polynom
# From the example for as.polynomial
 p <- as.polynomial(c(1,0,3,0))
 p
# 1 + 3*x^2 

 x1 <- c(2,1)
 x2 <- c(-1,3)
 px1 <- as.polynomial(x1)
 px2 <- as.polynomial(x2)

 px1*px2
# -2 + 5*x + 3*x^2 
 prod.p <- .Last.value
 str(prod.p)
# Class 'polynomial'  num [1:3] -2 5 3
 unclass(prod.p)
# [1] -2  5  3

Luisaluise

May I add that this does not allow symbolic manipulation of polynomial expressions, since the coefficients are stored as regular double values and precision is lost during arithmetic operations. If one needs to perform symbolic computations on such expressions, I'd suggest Mathematica instead. – Ataraxia

There is also an Ryacas package. – Luisaluise

If R cannot do what is wanted, I think Maxima does polynomials and is free. Although www.wolframalpha.com is free too and I think is based on Mathematica. – Clarindaclarine

R definitely can do what David asked: and the answer is using the polynom package, as long as you do not need exact arithmetic. – Salsify

For exact (or arbitrary precise) arithmetic, I'd also like to mention the gmp R package (based on the GNU MP ("GMP") library), for exact integer and rational arithmetic, and (my) package Rmpfr, based on the (GNU) MPFR library for arbitrary precision arithmetic. What really would be neat is a marriage of 'polynom' with gmp and Rmpfr. I'd be happy to collaborate on that. – Salsify

Answer 3

Polynomial multiplication is the convolution of the coefficients

convolve(c(2,1),rev(c(-1,3)),type="open")
#[1] -2  5  3

Ceilidh