对于多项式，获取其所有极值并通过突出显示所有单调部分来绘制它 [英] For a polynomial, get all its extrema and plot it by highlighting all monotonic pieces

问题描述

有人问了我这个有趣的问题，我认为值得在这里发布，因为在 Stack Overflow 上没有任何相关主题.

假设我在长度为 n 的向量 pc 中有多项式系数，其中变量 的多项式为 n - 1x 可以以其原始形式表示:

pc[1] + pc[2] * x + pc[3] * x ^ 2 + ... + pc[n] * x ^ (n - 1)

R 核心函数 polyroot 可以在 complex 域中找到该多项式的所有根.但通常我们也对极值感兴趣，对于单变量函数，局部极小值和极大值交替出现，将函数分解成单调的部分.

我的问题是:

1. 如何获得多项式的实域中的所有极值(实际上是所有鞍点)?
2. 如何用 2 种配色方案绘制这个多项式:红色代表上升部分，绿色代表下降部分?

最好把它写成一个函数，这样我们就可以轻松探索/可视化多项式.

例如，考虑一个 5 次多项式:

pc <- c(1, -2.2, -13.4, -5.1, 1.9, 0.52)

解决方案

获得多项式的所有鞍点

事实上，可以通过对多项式的一阶导数使用polyroot来找到鞍点.这是一个执行此操作的函数.

SaddlePoly <- 函数 (pc) {##多项式至少需要二次才能有鞍点如果(长度(pc)<3L){message("多项式至少需要二次才能有鞍点！")返回(数字(0))}## 一阶导数的多项式系数pc1 <- pc[-1] * seq_len(length(pc) - 1)## 复杂域的根croots <- polyroot(pc1)## 保留真实域中的根## 测试 0 是否为浮点数时要小心rroots <- Re(croots)[abs(Im(croots)) <1e-14]## 请注意，`polyroot` 返回多个具有乘数的根## 返回唯一的实数根(升序)排序(唯一(rroots))}xs <- SaddlePoly(pc)#[1] -3.77435640 -1.20748286 -0.08654384 2.14530617

<小时>

计算多项式

我们需要能够评估多项式才能绘制它.

Someone asked me this interesting question and I think it worthwhile posting it here, as there has not been any relevant thread on Stack Overflow.

Suppose I have polynomial coefficients in a length-n vector pc, where a polynomial of degree n - 1 for variable x can be expressed in its raw form:

pc[1] + pc[2] * x + pc[3] * x ^ 2 + ... + pc[n] * x ^ (n - 1)

R core function polyroot can find all roots of this polynomial in complex domain. But often we are also interested in extrema, as for a univariate function, local minima and maxima turn up alternately, breaking the function into monotonic pieces.

My questions are:

1. How to obtain all extrema (actually all saddle points) in real domain of a polynomial?
2. How to sketch this polynomial with 2-colour scheme: red for ascending pieces and green for descending pieces?

It would be good to write this up as a function so that we can easily explore / visualize a polynomial.

As an example, consider a polynomial of degree 5:

pc <- c(1, -2.2, -13.4, -5.1, 1.9, 0.52)

解决方案

obtain all saddle points of a polynomial

In fact, saddle points can be found by using polyroot on the 1st derivative of the polynomial. Here is a function doing it.

SaddlePoly <- function (pc) {
  ## a polynomial needs be at least quadratic to have saddle points
  if (length(pc) < 3L) {
    message("A polynomial needs be at least quadratic to have saddle points!")
    return(numeric(0))
    }
  ## polynomial coefficient of the 1st derivative
  pc1 <- pc[-1] * seq_len(length(pc) - 1)
  ## roots in complex domain
  croots <- polyroot(pc1)
  ## retain roots in real domain
  ## be careful when testing 0 for floating point numbers
  rroots <- Re(croots)[abs(Im(croots)) < 1e-14]
  ## note that `polyroot` returns multiple root with multiplicies
  ## return unique real roots (in ascending order)
  sort(unique(rroots))
  }

xs <- SaddlePoly(pc)
#[1] -3.77435640 -1.20748286 -0.08654384  2.14530617

---

evaluate a polynomial

We need be able to evaluate a polynomial in order to plot it. My this answer has defined a function g that can evaluate a polynomial and its arbitrary derivatives. Here I copy this function in and rename it to PolyVal.

PolyVal <- function (x, pc, nderiv = 0L) {
  ## check missing aruments
  if (missing(x) || missing(pc)) stop ("arguments missing with no default!")
  ## polynomial order p
  p <- length(pc) - 1L
  ## number of derivatives
  n <- nderiv
  ## earlier return?
  if (n > p) return(rep.int(0, length(x)))
  ## polynomial basis from degree 0 to degree `(p - n)`
  X <- outer(x, 0:(p - n), FUN = "^")
  ## initial coefficients
  ## the additional `+ 1L` is because R vector starts from index 1 not 0
  beta <- pc[n:p + 1L]
  ## factorial multiplier
  beta <- beta * factorial(n:p) / factorial(0:(p - n))
  ## matrix vector multiplication
  base::c(X %*% beta)
  }

For example, we can evaluate the polynomial at all its saddle points:

PolyVal(xs, pc)
#[1]  79.912753  -4.197986   1.093443 -51.871351

---

sketch a polynomial with a 2-colour scheme for monotonic pieces

Here is a function to view / explore a polynomial.

ViewPoly <- function (pc, extend = 0.1) {
  ## get saddle points
  xs <- SaddlePoly(pc)
  ## number of saddle points (if 0 the whole polynomial is monotonic)
  n_saddles <- length(xs)
  if (n_saddles == 0L) {
    message("the polynomial is monotonic; program exits!")
    return(NULL)
    }
  ## set a reasonable xlim to include all saddle points
  if (n_saddles == 1L) xlim <- c(xs - 1, xs + 1)
  else xlim <- extendrange(xs, range(xs), extend)
  x <- c(xlim[1], xs, xlim[2])
  ## number of monotonic pieces
  k <- length(xs) + 1L 
  ## monotonicity (positive for ascending and negative for descending)
  y <- PolyVal(x, pc)
  mono <- diff(y)
  ylim <- range(y)
  ## colour setting (red for ascending and green for descending)
  colour <- rep.int(3, k)
  colour[mono > 0] <- 2
  ## loop through pieces and plot the polynomial
  plot(x, y, type = "n", xlim = xlim, ylim = ylim)
  i <- 1L
  while (i <= k) {
    ## an evaluation grid between x[i] and x[i + 1]
    xg <- seq.int(x[i], x[i + 1L], length.out = 20)
    yg <- PolyVal(xg, pc)
    lines(xg, yg, col = colour[i])
    i <- i + 1L
    }
  ## add saddle points
  points(xs, y[2:k], pch = 19)
  ## return (x, y)
  list(x = x, y = y)
  }

We can visualize the example polynomial in the question by:

ViewPoly(pc)
#$x
#[1] -4.07033952 -3.77435640 -1.20748286 -0.08654384  2.14530617  2.44128930
#
#$y
#[1]  72.424185  79.912753  -4.197986   1.093443 -51.871351 -45.856876

这篇关于对于多项式，获取其所有极值并通过突出显示所有单调部分来绘制它的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！