R：如何获得两个分布的总和？

Question

我有一个简单的问题。我想对两个非参数分布求和。

这是一个例子。有两个城市有10间房子。我们知道每个房子的能源消耗。 (edited) 我想得到从每个城市中随机选择的房子的总和的概率分布。

A1 <- c(1,2,3,3,3,4,4,5,6,7) #10 houses' energy consumption for city A
B1 <- c(11,13,15,17,17,18,18,19,20,22) #10 houses' energy consumption for city B

我有A1和B1的概率分布，如何得到A1+B1的概率分布？如果我只是在 R 中使用 A1+B1，它会给出 12 15 18 20 20 22 22 24 26 29。但是，我认为这是不对的。因为家里没有秩序。

当我改变房子的顺序时，它给出了另一个结果。

# Original
A1 <- c(1,2,3,3,3,4,4,5,6,7)
B1 <- c(11,13,15,17,17,18,18,19,20,22)
#change order 1
A2 <- c(7,6,5,4,4,3,3,3,2,1) 
B2 <- c(22,20,19,18,18,17,17,15,13,11)
#change order 2
A3 <- c(3,3,3,4,4,5,6,7,1,2) 
B3 <- c(17,17,18,18,19,13,20,11,22,15)
sum1 <- A1+B1; sum1
sum2 <- A1+B2; sum2
sum3 <- A3+B3; sum3

红线是 sum1、sum2 和 sum3。我不确定我怎样才能得到两个 distributions.Please 之和的分布给我任何 ideas.Thanks!

（如果这些分布是正态分布或均匀分布，我可以很容易地得到分布的总和，但这些不是正态分布并且没有顺序）

Answer 1

你可能想要这样的东西：

rowSums(expand.grid(A1, B1))

使用 expand.grid 将为您提供 A1 和 B1 的所有组合的数据框，rowSums 将添加它们。

Answer 2

编辑：

既然我更好地理解了这个问题，并且看到了@jeremycg 的回答，我想我有一个不同的方法，我认为它会随着样本量的增加而更好地扩展。

与其依赖 A1 和 B1 中的值作为分布中的唯一值，我们还可以推断这些只是分布中的样本。为避免对分布施加特定形式，我将使用经验值 'equivalent'：样本密度。如果我们使用 density 函数，我们可以推断出从任一城镇抽取连续范围的家庭能源使用样本的相对概率。我们可以从 density()$x 值中随机抽取任意数量的能量（有替换），其中我们采用的 sample 是用 prob=density()$y 加权的……即峰值密度图位于 x 值处，应该更频繁地重新采样。

作为启发式，一个过于简单的陈述可以说 mean(A1) 是 3.8，而 mean(B1) 是 17，所以这两个城市的能源使用总和应该是，平均而言，~20.8 .将其用作 "does it make sense test"/ 启发式方法，我认为以下方法符合您想要的结果类型。

sample_sum <- function(A, B, n, ...){
    qss <- function(X, n, ...){
        r_X <- range(X)
        dens_X <- density(X, ...)
        sample(dens_X$x, size=n, prob=dens_X$y, replace=TRUE)
    }

    sample_A <- qss(A, n=n, ...)
    sample_B <- qss(B, n=n, ...)

    sample_A + sample_B
}

ss <- sample_sum(A1, B1, n=100, from=0)

png("~/Desktop/answer.png", width=5, height=5, units="in", res=150)
plot(density(ss))
dev.off()

请注意，我将密度图限制为 0，因为我假设您不想推断负能量。我看到合成密度的峰值正好在 20 以上，所以 'it makes sense'.

这里的潜在优势在于，您无需查看两个城市房屋的每一种可能的能源组合即可了解能源使用总和的分布。如果你能定义两者的分布，你就可以定义成对和的分布。

最后，计算时间是微不足道的，特别是与寻找所有组合的方法相比。例如，每个城市有 1000 万栋房屋，如果我尝试使用 expand.grid 方法，我会得到 Error: cannot allocate vector of size 372529.0 Gb 错误，而 sample_sum 方法需要 0.12 秒。

当然，如果答案对你没有帮助，那速度就一文不值了;)

Answer 3

理论上，两个随机变量的总和分布是它们PDF的卷积，details，如：

PDF(Z) = PDF(Y) * PDF(X)

所以，我认为这种情况可以通过convolution来计算。

# your data
A1 <- c(1,2,3,3,3,4,4,5,6,7) #10 houses' energy consumption for city A
B1 <- c(11,13,15,17,17,18,18,19,20,22) #10 houses' energy consumption for city B

# compute PDF/CDF
PDF_A1 <- table(A1)/length(A1)
CDF_A1 <- cumsum(PDF_A1)

PDF_B1 <- table(B1)/length(B1)
CDF_B1 <- cumsum(PDF_B1)

# compute the sum distribution 
PDF_C1 <- convolve(PDF_B1, PDF_A1, type = "open")

# plotting
plot(PDF_C1, type="l", axe=F, main="PDF of A1+B1")
box()
axis(2)
# FIXME: is my understand for X correct?
axis(1, at=seq(1:14), labels=(c(names(PDF_A1)[-1],names(PDF_B1))))

注：

CDF: cumulative distribution function

PDF: probability density function

## To make the x-values correspond to actually sums, consider
## compute PDF
## pad zeros in probability vectors to convolve
r <- range(c(A1, B1))
pdfA <- pdfB <- vector('numeric', diff(r)+1L)
PDF_A1 <- table(A1)/length(A1)                        # same as what you have done
PDF_B1 <- table(B1)/length(B1)
pdfA[as.numeric(names(PDF_A1))] <- as.vector(PDF_A1)  # fill the values
pdfB[as.numeric(names(PDF_B1))] <- as.vector(PDF_B1)

## compute the convolution and plot
res <- convolve(pdfA, rev(pdfB), type = "open")
plot(res, type="h", xlab='Sum', ylab='')

## In this simple case (with discrete distribution) you can compare
## to previous solution
tst <- rowSums(expand.grid(A1, B1))
plot(table(tst) / sum(as.vector(table(tst))), type='h')

Answer 4

在添加之前对分布进行排序不是解决了这个问题吗？

A1 <- c(1,2,3,3,3,4,4,5,6,7) #10 houses' energy consumption for city A
B1 <- c(11,13,15,17,17,18,18,19,20,22) #10 houses' energy consumption for city B
sort(A1)+sort(B1)

R：如何获得两个分布的总和？

R: How to get a sum of two distributions?

r

sum

distribution