使用 sapply() 了解 R 卷积代码

Question

我正在尝试分解 this post 中的 R 代码：

x <- c(0.17,0.46,0.62,0.08,0.40,0.76,0.03,0.47,0.53,0.32,0.21,0.85,0.31,0.38,0.69)


convolve.binomial <- function(p) {
  # p is a vector of probabilities of Bernoulli distributions.
  # The convolution of these distributions is returned as a vector
  # `z` where z[i] is the probability of i-1, i=1, 2, ..., length(p)+1.
  n <- length(p) + 1
  z <- c(1, rep(0, n-1))
  sapply(p, function(q) {z <<- (1 - q) * z + q * (c(0, z[-n])); q})
  z
}
convolve.binomial(x)
 [1] 5.826141e-05 1.068804e-03 8.233357e-03 3.565983e-02 9.775029e-02
 [6] 1.804516e-01 2.323855e-01 2.127628e-01 1.394564e-01 6.519699e-02
[11] 2.141555e-02 4.799630e-03 6.979119e-04 6.038947e-05 2.647052e-06
[16] 4.091095e-08

我在 RStudio 中尝试 debugging，但它仍然不透明。

问题在于行：sapply(p, function(q) {z <<- (1 - q) * z + q * (c(0, z[-n])); q})。

我猜想在调用的上下文中 convolve.binomial(x) p = q = x。至少如果我将线拉到函数之外我会得到相同的结果运行 sapply(x, function(x) {z <<- (1 - x) * z + x * (c(0, z[-n])); x}):

x <- c(0.17,0.46,0.62,0.08,0.40,0.76,0.03,0.47,0.53,0.32,0.21,0.85,0.31,0.38,0.69)
n <- length(x) + 1
z <- c(1, rep(0, n-1))
#  [1] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
sapply(x, function(x) {z <<- (1 - x) * z + x * (c(0, z[-n])); x})
z # Is extracted by calling it and contains the correct result

我的问题是：

;q} 以 sapply() 结尾的目的是什么？
它与 <<- 符号有何关系，旨在使 z 可在 "implicit" 循环之外访问，即 sapply()？

下面你可以看到我的问题"hacking"这行代码：

(x_complem = 1 - x)
sapply(x, function(x) {z <<- x_complem * z + x * (c(0, z[-n])); x})
z # Returns 16 values and warnings

z_offset = c(0, z[-n])
#  [1] 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
sapply(x, function(x) {z <<- (1 - x) * z + x * z_offset; x})
z # Returns different values.

Answer 1

;q} 以 sapply() 结尾的目的是什么？

sapplyreturnq里面的功能，其实不需要。以下功能将同样起作用。

convolve.binomial <- function(p) {
  n <- length(p) + 1
  z <- c(1, rep(0, n-1))
  sapply(p, function(q) {z <<- (1 - q) * z + q * (c(0, z[-n]))})
  z
}

它与 <<- 符号有何关系，旨在使 z 可在 "implicit" 循环之外访问，即 sapply()？

在 R 中，如果您使用 ?'<<-' 搜索 <<- 运算符的文档，它会说

The operators <<- and >>- are normally only used in function, and cause a search to be made through parent environments for an existing definition of the variable to be assigned. If such as variable is found (and its binding is not locked) then its value is redefined, otherwise assignment takes place in the global environment.

在函数 convolve.binomial 中，值 z 是在函数局部定义的。所以z <<-实际上是在convolve.binomial函数中重新定义了z。

总而言之，sapply 调用中的 z <<- 更改了 convolve.binomial 中已经定义的 z 变量，我们最终 return 这个 z。 sapply().

中不需要 ;q} 结尾

Answer 2

如果您想在函数执行过程中查看 z 的中间值，请在以下代码中插入 cat 或 print 命令：

sapply(x, function(x) {z <<- (1 - x) * z + x * (c(0, z[-n])); cat(z,"\n"); x})
#--------
0.83 0.17 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
0.4482 0.4736 0.0782 0 0 0 0 0 0 0 0 0 0 0 0 0 
0.170316 0.457852 0.323348 0.048484 0 0 0 0 0 0 0 0 0 0 0 0 
0.1566907 0.4348491 0.3341083 0.07047312 0.00387872 0 0 0 0 0 0 0 0 0 0 0 
0.09401443 0.3235858 0.3744046 0.1759272 0.03051648 0.001551488 0 0 0 0 0 0 0 0 0 0 
0.02256346 0.1491116 0.3357823 0.3267701 0.1410286 0.02356488 0.001179131 0 0 0 0 0 0 0 0 0 
snipped rest of output

我认为这更清楚地表明，每个中间步骤都代表一系列事件的一组概率。每行总和为 1.0，代表可能存在较少数量的二项式参数时个体计数存活的概率。最终结果显示组装完整序列后特定计数总和的概率。

另一个有趣的特征是这个结果在 x 中的概率随机重新排序下是不变的（对于原始问题应该是这样）。检查

的中间结果

plot(x)
lines(seq(length(z)), z)
z2 <- convolve.binomial(sample(x) )
lines(seq(length(z)), z2, col="red" )
z3 <- convolve.binomial(sample(x) )
lines(seq(length(z)), z3, col="blue" )

使用 sapply() 了解 R 卷积代码

Understanding R convolution code with sapply()

r

convolution

sapply