避免在构建分数多项式函数时使用eval(parse())

4
我的目标是编写一个R函数,接受分数多项式(FP)的系数,并返回一个向量化的函数,该函数可以为给定的输入数字计算指定的FP。FP定义有两个重要规则:
  • x ^ 0被定义为log(x)
  • 幂可以有多个系数,其中幂p的第二个系数将一个log(x)因子添加到加法项(x ^ p * log(x)),第三个系数添加一个log(x)^ 2因子(x ^ p * log(x)^ 2),依此类推

我的当前解决方案是构建FP函数作为字符串,解析字符串并返回评估表达式的函数。我的问题是是否有更好/更快的方法,可以避免使用eval(parse()) - 可能使用一些substitute()魔法。

该函数必须处理在调用时未知每个幂的系数数量,但在调用时指定。最终的FP评估需要快速,因为它经常被调用。

希望不仅限于标准的幂次-2,-1,-0.5,0,0.5,1,2,3,这将是一个好的选择。理想情况下,所需函数可以同时执行两个步骤:接受FP系数和数字向量,并返回输入的FP值,同时仍然保持快速。

getFP <- function(p_2, p_1, p_0.5, p0, p0.5, p1, p2, p3, ...) {
    p <- as.list(match.call(expand.dots=TRUE)[-1])         # all args
    names(p) <- sub("^p", "", names(p))     # strip "p" from arg names
    names(p) <- sub("_", "-", names(p))     # replace _ by - in arg names

    ## for one power and the i-th coefficient: build string
    getCoefStr <- function(i, pow, coefs) {
        powBT  <- ifelse(as.numeric(pow), paste0("x^(", pow, ")"), "log(x)")
        logFac <- ifelse(i-1,             paste0("*log(x)^", i-1), "")
        paste0("(", coefs[i], ")*", powBT, logFac)
    }

    onePwrStr <- function(pow, p) { # for one power: build string for all coefs
        coefs  <- eval(p[[pow]])
        pwrStr <- sapply(seq(along=coefs), getCoefStr, pow, coefs)
        paste(pwrStr, collapse=" + ")
    }

    allPwrs <- sapply(names(p), onePwrStr, p)  # for each power: build string
    fpExpr  <- parse(text=paste(allPwrs, collapse=" + "))
    function(x) { eval(fpExpr) }
}

一个例子是-1.5*x^(-1) - 14*log(x) - 13*x^(0.5) + 6*x^0.5*log(x) + 1*x^3,其中指定了幂次 (-1, 0, 0.5, 0.5, 3),对应的系数为 (-1.5, -14, -13, 6, 1)。
> fp <- getFP(p_1=-1.5, p0=-14, p0.5=c(-13, 6), p3=1)
> fp(1:3)
[1] -13.50000000 -14.95728798   0.01988127

是的,挖掘变量名称以删除或替换它们看起来确实需要很多工作。那么要求输入为向量怎么样?例如 function(pwrs,coeffs),其中 pwrscoeffs 是幂和系数的等长向量? - Carl Witthoft
你只需要将多项式写成一个用expression包装的方程,然后使用substitute在适当的位置交换值,例如substitute(expression(p1*x + p2*x),list(p1=p1, p0=p0)) - Thomas
@CarlWitthoft 可能有这个可能性,但是我不确定如何处理当幂有多个系数时(不知道有多少系数)添加 log(x) 因子的规则。似乎需要识别重复的幂,并在它们的长度上进行 for() 循环。 - caracal
你可以使用与当前示例相同的条件逻辑,有条件地用 log(x) 来替换某些内容。您可以替换任何幂和系数...我将上面留下的作为注释,因为它远非完整的解决方案。 - Thomas
1
尝试简化您的问题,然后再添加这些复杂的细节可能会有所帮助;另外,在提问中增加更多关于您的函数在系数和幂方面需要多灵活的细节也是有帮助的。 - Thomas
显示剩余2条评论
1个回答

6
首先我们创建一个函数,用于生成序列中的单个项。
one <- function(p, c = 1, repeated = 1) {
  if (p == 0) {
    lhs <- substitute(c * log(x), list(c = c))
  } else {
    lhs <- substitute(c * x ^ p, list(c = c, p = p))
  }

  if (repeated == 1) return(lhs)
  substitute(lhs * log(x) ^ pow, list(lhs = lhs, pow = repeated - 1))
}
one(0)
# 1 * log(x)
one(2)
# 1 * x^2

one(2, 2)
# 2 * x^2

one(2, r = 2)
# 1 * x ^ 2 * log(x)^1
one(2, r = 3)
# 1 * x ^ 2 * log(x)^2

这里的重要工具是substitute(),可以在这里找到其解释。
接下来,我们编写一个将两个术语相加的函数。同样,这使用了substitute:
add_expr_1 <- function(x, y) {
  substitute(x + y, list(x = x, y = y))
}

add_expr_1(one(0, 1), one(2, 1))

我们可以使用这个方法来创建一个函数,将任意数量的项相加:
add_expr <- function(x) Reduce(add_expr_1, x)
add_expr(list(one(0, 1), one(1, 1), one(2, 3)))

有了这些部分,最终函数就很简单——我们计算出重复次数,然后使用Map()调用one()函数来处理每个powerscoefsreps的组合:

fp <- function(powers, coefs) {
  # Determine number of times each power is repeated. This is too
  # clever approach but I think it works
  reps <- ave(powers, powers, FUN = seq_along)

  # Now generate a list of expressions using one
  components <- Map(one, powers, coefs, reps)

  # And combine them together with plus
  add_expr(components)
}

powers <- c(-1, 0, 0.5, 0.5, 3)
coefs <-  c(-1.5, -14, -13, 6, 1)
fp(powers, coefs)
# -1.5 * x^-1 + -14 * log(x) + -13 * x^0.5 + 6 * x^0.5 * log(x)^1 + 
#   1 * x^3

整洁,非常整洁,谢谢!这正是我所期望的substitute()魔法。 - caracal
另外,我认为也可以使用 sequence(rle(powers)$lengths) 找到每个幂重复的次数。 - caracal
1
@caracal 如果重复不是连续的,我认为那样做不会起作用。 - hadley

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接