Serge Kashlik
Serge Kashlik

Reputation: 413

R code for convolution of multiple random variables

I am trying to empirically prove the central limit theorem by adding a simple uniform distribution multiple times.

I have mostly worked with VBA before so the naturally I used a for loop, however I am curious if there is any way to simplify the code in R, as I would like to be able to run it up to at least 50 times. Code is shown below:

A <- seq(1,8,1)
B <- seq(1,8,1)
C <- seq(1,8,1)
D <- seq(1,8,1)


p.Z <- as.integer()

for(a in 1:8){
    for(b in 1:8){
        for(c in 1:8){
            for(d in 1:8){
                Z <- A[a] + B[b] + C[c] + D[d]
                p.Z <- c(p.Z, Z)
            }
        }
    }
}

Upvotes: 1

Views: 222

Answers (2)

akrun
akrun

Reputation: 887213

We can use

p1 <- Reduce(`+`, expand.grid(rep(list(1:8), 4)))
all.equal(p1, p.Z)
#[1] TRUE

Or with CJ from data.table

CJ(a = 1:8, b = 1:8, c = 1:8, d = 1:8)[, Reduce(`+`, .SD)]

Upvotes: 1

Ronak Shah
Ronak Shah

Reputation: 389012

I don't know about the statistics behind it but the code can definitely be simplified using expand.grid and rowSums. expand.grid creates all combinations of input vector.

p1 <- rowSums(expand.grid(1:8, 1:8, 1:8, 1:8))

Comparing it with p.Z from OP's for loop

identical(p1, p.Z)
#[1] TRUE

To do this dynamically, we can use do.call with replicate

n <- 4
rowSums(do.call(expand.grid, replicate(n, 1:8, simplify = FALSE)))

The same behavior can also be achieved using crossing from tidyr

rowSums(tidyr::crossing(a = 1:8, b = 1:8, c = 1:8, d = 1:8))

Upvotes: 1

Related Questions