機能（シンボリック法優先）

Question

私は私のデータから多項式の係数を見つけた：

R <- c(0.256,0.512,0.768,1.024,1.28,1.437,1.594,1.72,1.846,1.972,2.098,2.4029) 
Ic <- c(1.78,1.71,1.57,1.44,1.25,1.02,0.87,0.68,0.54,0.38,0.26,0.17) 
NN <- 3 
ft <- lm(Ic ~ poly(R, NN, raw = TRUE)) 
pc <- coef(ft) 

だから私は多項式関数を作成することができます

f1 <- function(x) pc[1] + pc[2] * x + pc[3] * x^2 + pc[4] * x^3 

そして、たとえば、派生物を取る：

g1 <- Deriv(f1) 

どのようにしてeに書き直す必要がないように汎用関数を作成するか非常に新しい多項式度NN？

Answer 1

私の元の答えは、あなたが本当に欲しいものではないかもしれません。数値的ではなく象徴的であったからです。ここに象徴的な解決策があります。

## use `"x"` as variable name 
## taking polynomial coefficient vector `pc` 
## can return a string, or an expression by further parsing (mandatory for `D`) 
f <- function (pc, expr = TRUE) { 
    stringexpr <- paste("x", seq_along(pc) - 1, sep = "^") 
    stringexpr <- paste(stringexpr, pc, sep = " * ") 
    stringexpr <- paste(stringexpr, collapse = " + ") 
    if (expr) return(parse(text = stringexpr)) 
    else return(stringexpr) 
    } 

## an example cubic polynomial with coefficients 0.1, 0.2, 0.3, 0.4 
cubic <- f(pc = 1:4/10, TRUE) 

## using R base's `D` (requiring expression) 
dcubic <- D(cubic, name = "x") 
# 0.2 + 2 * x * 0.3 + 3 * x^2 * 0.4 

## using `Deriv::Deriv` 
library(Deriv) 

dcubic <- Deriv(cubic, x = "x", nderiv = 1L) 
# expression(0.2 + x * (0.6 + 1.2 * x)) 

Deriv(f(1:4/10, FALSE), x = "x", nderiv = 1L) ## use string, get string 
# [1] "0.2 + x * (0.6 + 1.2 * x)" 

もちろん、Derivは、より高次の派生語を得るのを容易にします。単にnderivと設定することができます。しかし、Dについては、再帰を使用する必要があります（例：?Dを参照）。私たちが式を使用する場合

Deriv(cubic, x = "x", nderiv = 2L) 
# expression(0.6 + 2.4 * x) 

Deriv(cubic, x = "x", nderiv = 3L) 
# expression(2.4) 

Deriv(cubic, x = "x", nderiv = 4L) 
# expression(0) 

、我々は後に、結果を評価することができるようになります。たとえば、

eval(cubic, envir = list(x = 1:4)) ## cubic polynomial 
# [1] 1.0 4.9 14.2 31.3 

eval(dcubic, envir = list(x = 1:4)) ## its first derivative 
# [1] 2.0 6.2 12.8 21.8 

上記のことは、関数の式をまとめることができることを意味します。関数の使用にはいくつかの利点があります.1つは、curveまたはplot.functionを使用してプロットできる点です。

fun <- function(x, expr) eval.parent(expr, n = 0L) 

注意、funの成功は、シンボルxの観点で表現するexprが必要です。例えばexprがyで定義された場合、funをfunction (y, expr)と定義する必要があります。今度は、範囲0 < x < 5に、cubicとdcubicをプロットするcurveを使用してみましょう：

curve(fun(x, cubic), from = 0, to = 5) ## colour "black" 
curve(fun(x, dcubic), add = TRUE, col = 2) ## colour "red" 

最も便利な方法は、むしろf + fun組み合わせを行うよりも、単一の関数FUNを定義することはもちろんです。このように、fとfunで使用される変数名の一貫性についても心配する必要はありません。

FUN <- function (x, pc, nderiv = 0L) { 
    ## check missing arguments 
    if (missing(x) || missing(pc)) stop ("arguments missing with no default!") 
    ## expression of polynomial 
    stringexpr <- paste("x", seq_along(pc) - 1, sep = "^") 
    stringexpr <- paste(stringexpr, pc, sep = " * ") 
    stringexpr <- paste(stringexpr, collapse = " + ") 
    expr <- parse(text = stringexpr) 
    ## taking derivatives 
    dexpr <- Deriv::Deriv(expr, x = "x", nderiv = nderiv) 
    ## evaluation 
    val <- eval.parent(dexpr, n = 0L) 
    ## note, if we take to many derivatives so that `dexpr` becomes constant 
    ## `val` is free of `x` so it will only be of length 1 
    ## we need to repeat this constant to match `length(x)` 
    if (length(val) == 1L) val <- rep.int(val, length(x)) 
    ## now we return 
    val 
    } 

我々は、係数pc <- c(0.1, 0.2, 0.3, 0.4)とx <- seq(0, 1, 0.2)上のその誘導体との三次多項式を評価したいと、私たちは、単に行うことができます。

FUN(x, pc) 
# [1] 0.1000 0.1552 0.2536 0.4144 0.6568 1.0000 

FUN(x, pc, nderiv = 1L) 
# [1] 0.200 0.368 0.632 0.992 1.448 2.000 

FUN(x, pc, nderiv = 2L) 
# [1] 0.60 1.08 1.56 2.04 2.52 3.00 

FUN(x, pc, nderiv = 3L) 
# [1] 2.4 2.4 2.4 2.4 2.4 2.4 

FUN(x, pc, nderiv = 4L) 
# [1] 0 0 0 0 0 0 

は今プロットも簡単です：

curve(FUN(x, pc), from = 0, to = 5) 
curve(FUN(x, pc, 1), from = 0, to = 5, add = TRUE, col = 2) 
curve(FUN(x, pc, 2), from = 0, to = 5, add = TRUE, col = 3) 
curve(FUN(x, pc, 3), from = 0, to = 5, add = TRUE, col = 4) 

Answer 2

私の最終ソルス最終的には象徴的な派生を伴うイオンが長すぎるため、私は数値計算のために別のセッションを使用します。これを多項式の場合と同様に行うことができ、派生物が明示的に知られているので、それらをコード化することができます。注意してください、ここでR表現の使用はありません。すべてが関数を使用して直接行われます。

だから我々は最初の係数及び階乗乗数を掛けた後、程度p - nに度0から多項式基底を生成します。ここでよりpolyを使用する方が便利です。結果は、我々はシンボリック溶液中FUNを持っているものと一致している

x <- seq(0, 1, by = 0.2) 
pc <- 1:4/10 

g(x, pc, 0) 
# [1] 0.1000 0.1552 0.2536 0.4144 0.6568 1.0000 

g(x, pc, 1) 
# [1] 0.200 0.368 0.632 0.992 1.448 2.000 

g(x, pc, 2) 
# [1] 0.60 1.08 1.56 2.04 2.52 3.00 

g(x, pc, 3) 
# [1] 2.4 2.4 2.4 2.4 2.4 2.4 

g(x, pc, 4) 
# [1] 0 0 0 0 0 0 

：

## use `outer` 
g <- 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 
    drop(X %*% beta) 
    } 

我々はまだ例xとシンボリック溶液中で定義されたpcを使用しています。

同様に、我々はcurveを使用してgをプロットすることができます

curve(g(x, pc), from = 0, to = 5) 
curve(g(x, pc, 1), from = 0, to = 5, col = 2, add = TRUE) 
curve(g(x, pc, 2), from = 0, to = 5, col = 3, add = TRUE) 
curve(g(x, pc, 3), from = 0, to = 5, col = 4, add = TRUE) 

Answer 3

今、私たちはRパッケージpolynomを使用することを検討して、この質問自分自身をうまくできる方法を示すにはかなり多くの努力の後。小さなパッケージとして、単変量多項式の構築、派生、統合、算術、および根発見を実装することを目指しています。このパッケージはコンパイルされたコードなしでR言語で完全に書かれています。

## install.packages("polynom") 
library(polynom) 

以前に使用された3次多項式の例を考えます。

pc <- 1:4/10 

## step 1: making a "polynomial" object as preparation 
pcpoly <- polynomial(pc) 
#0.1 + 0.2*x + 0.3*x^2 + 0.4*x^3 

## step 2: compute derivative 
expr <- deriv(pcpoly) 

## step 3: convert to function 
g1 <- as.function(expr) 

#function (x) 
#{ 
# w <- 0 
# w <- 1.2 + x * w 
# w <- 0.6 + x * w 
# w <- 0.2 + x * w 
# w 
#} 
#<environment: 0x9f4867c> 

ステップバイステップの構成では、結果として得られる関数の中にはすべてのパラメータがあります。 xの値には単一の引数が必要です。対照的に、他の2つの答えの関数は、係数と派生順序も必須の引数として使用します。私たちは、他の2つの答えで見る同じグラフを生成するには

g1(seq(0, 1, 0.2)) 
# [1] 0.200 0.368 0.632 0.992 1.448 2.000 

この関数を呼び出すことができ、我々としても他の誘導体を得る：

g0 <- as.function(pcpoly) ## original polynomial 

## second derivative 
expr <- deriv(expr) 
g2 <- as.function(expr) 
#function (x) 
#{ 
# w <- 0 
# w <- 2.4 + x * w 
# w <- 0.6 + x * w 
# w 
#} 
#<environment: 0x9f07c68> 

## third derivative 
expr <- deriv(expr) 
g3 <- as.function(expr) 
#function (x) 
#{ 
# w <- 0 
# w <- 2.4 + x * w 
# w 
#} 
#<environment: 0x9efd740> 

はおそらく、あなたはすでに私がnderivを指定していないことに気づきました、再帰的に一度に1微分をとる。これはこのパッケージの欠点かもしれません。高次の導関数を容易にするものではありません。

今、私たちは、プロット

## As mentioned, `g0` to `g3` are parameter-free 
curve(g0(x), from = 0, to = 5) 
curve(g1(x), add = TRUE, col = 2) 
curve(g2(x), add = TRUE, col = 3) 
curve(g3(x), add = TRUE, col = 4)