set.seed(961)
mu <- rep(0, 4)
Sigma <- matrix(c(1, -1, 0, 0,
                  -1, 1, 0, 0,
                  0, 0, 1, -1,
                  0, 0, -1, 1), nrow = 4, ncol = 4)
Sigma
#      [,1] [,2] [,3] [,4]
# [1,]    1   -1    0    0
# [2,]   -1    1    0    0
# [3,]    0    0    1   -1
# [4,]    0    0   -1    1

library(MASS)
rawvars <- mvrnorm(n = 1e4, mu = mu, Sigma = Sigma)

# cov(rawvars)

cor(rawvars)
#             [,1]         [,2]         [,3]         [,4]
# [1,]  1.000000000 -1.000000000 -0.006839596  0.006839597
# [2,] -1.000000000  1.000000000  0.006839597 -0.006839598
# [3,] -0.006839596  0.006839597  1.000000000 -1.000000000
# [4,]  0.006839597 -0.006839598 -1.000000000  1.000000000

pvars <- pnorm(rawvars)

binvars <- qbinom(pvars, 4, .3)

summary(as.factor(binvars))
#    0     1     2     3     4
# 9704 16283 10676  2994   343

hist(binvars)