3

I can use following MyFun (user defined) function to simulate observations. But could not figure out how to repeat this function with replicate function when the function returns a list of output.

MyFun <- function(nSim, Size, Prob) {
  M1 <- t(mapply(rbinom, prob = Prob, n = nSim, size = Size))
  dimnames(M1) <- list(Prob, paste0("V", 1:nSim))
  MeanM1 <- M1/Size
  Results1 <- list(M1, MeanM1)
  return(Results1)
}

MyFun(nSim=5, Size=4, Prob=c(0.2, 0.4))

[[1]]
    V1 V2 V3 V4 V5
0.2  2  2  1  2  1
0.4  2  3  0  3  1

[[2]]
     V1   V2   V3   V4   V5
0.2 0.5 0.50 0.25 0.50 0.25
0.4 0.5 0.75 0.00 0.75 0.25

But doesn't give the desired results with replicate function.

replicate(
    n=2
  , MyFun(nSim=5, Size=2, Prob=c(0.2, 0.4))
  , simplify = "array"
  )

    [,1]       [,2]      
[1,] Integer,10 Integer,10
[2,] Numeric,10 Numeric,10
G5W
  • 36,531
  • 10
  • 47
  • 80
MYaseen208
  • 22,666
  • 37
  • 165
  • 309
  • 1
    Maybe just `lapply(1:2, function(x) MyFun(nSim=5, Size=4, Prob=c(0.2, 0.4)))`? – David Arenburg Aug 02 '15 at 17:56
  • I answered but now I am realizing... how is it possible that the output is a function is a matrix in which every element is a matrix itself? – SabDeM Aug 02 '15 at 18:02

2 Answers2

3

Your function does work, it has only a "strange" output. Try to store the result into an object called aa.

 aa <- replicate(
+     n=2
+     , MyFun(nSim=5, Size=2, Prob=c(0.2, 0.4))
+     , simplify = T
+ )
 aa
     [,1]       [,2]      
[1,] Integer,10 Integer,10
[2,] Numeric,10 Numeric,10

but every element of that 2 x 2 matrix is itself a matrix. Try this:

 class(aa)
[1] "matrix"
 sapply(aa, class)
[1] "matrix" "matrix" "matrix" "matrix"

In fact if you extract the first "cell" of the previous matrix, you'll see that is a matrix that does contains the output of the MyFun:

 aa[1,1]
[[1]]
    V1 V2 V3 V4 V5
0.2  0  1  0  0  1
0.4  1  2  0  2  0
SabDeM
  • 7,050
  • 2
  • 25
  • 38
3

Try adding a simplify argument to your function, and pass it to mapply().

MyFun <- function(nSim, Size, Prob, simplify = "array") {
    M1 <- t(
        mapply(rbinom, prob = Prob, n = nSim, size = Size, SIMPLIFY = simplify)
    )
    dimnames(M1) <- list(Prob, paste0("V", 1:nSim))
    MeanM1 <- M1/Size
    Results1 <- list(M1, MeanM1)
    return(Results1)
}

Then call it like this, using simplify = FALSE in replicate(). This will produce a list of length 2 with a list of two arrays per element.

replicate(2, MyFun(nSim=5, Size=2, Prob=c(0.2, 0.4)), simplify = FALSE)

[[1]]
[[1]][[1]]
    V1 V2 V3 V4 V5
0.2  1  0  2  0  0
0.4  2  0  1  1  0

[[1]][[2]]
     V1 V2  V3  V4 V5
0.2 0.5  0 1.0 0.0  0
0.4 1.0  0 0.5 0.5  0


[[2]]
[[2]][[1]]
    V1 V2 V3 V4 V5
0.2  0  0  0  0  1
0.4  2  1  1  1  0

[[2]][[2]]
    V1  V2  V3  V4  V5
0.2  0 0.0 0.0 0.0 0.5
0.4  1 0.5 0.5 0.5 0.0
Rich Scriven
  • 97,041
  • 11
  • 181
  • 245
  • Excellent anwwer @Richard. Thanks a lot. Do you think my code is efficiently written for large number of simulations and replications? Thanks – MYaseen208 Aug 02 '15 at 17:56
  • I dunno, use it to write some benchmarks :). It seems okay to me, but the only thing I would be concerned about is `mapply()`. But it's not slow so you should be alright. Test it! – Rich Scriven Aug 02 '15 at 18:00