pkg/tests/testthat/test-Moments.R

   1 test_that("both versions of Moments_Mi agree on various inputs",
   2 {
   3   for (n in c(20,200))
   4   {
   5     for (d in c(2,10,20))
   6     {
   7       E <- diag(d)
   8       Id <- as.double(E)
   9       X <- matrix( rnorm(n*d), nrow=n )
  10       Y <- rbinom(n, 1, .5)
  11       M2 <- as.double(.Moments_M2(X,Y))
  12       M2_R <- colMeans(Y * t( apply(X, 1, function(Xi) Xi %o% Xi - Id) ))
  13       expect_equal(max(abs(M2 - M2_R)), 0)
  14       M3 <- as.double(.Moments_M3(X,Y))
  15       M3_R <- colMeans(Y * t( apply(X, 1, function(Xi) {
  16         return (Xi %o% Xi %o% Xi -
  17           Reduce('+', lapply(1:d, function(j)
  18             as.double(Xi %o% E[j,] %o% E[j,])), rep(0, d*d*d)) -
  19           Reduce('+', lapply(1:d, function(j)
  20             as.double(E[j,] %o% Xi %o% E[j,])), rep(0, d*d*d)) -
  21           Reduce('+', lapply(1:d, function(j)
  22             as.double(E[j,] %o% E[j,] %o% Xi)), rep(0, d*d*d))) } ) ))
  23       expect_equal(max(abs(M3 - M3_R)), 0)
  24     }
  25   }
  26 })