r - pick 10 random numbers from standard normal distribution whose sum equals 5

I am having coding the following in R I want to pick 10 random numbers from a standard normal distribution whose sum equals 5. I have the following code so far (below), but this returns "numeric(0)" when the random numbers don't satisfy the condition. What I want is to choose 10 random numbers that do satisfy this condition. Is there a way to re-scale these numbers once they have been picked, or can I somehow insert a condition for this into "rnorm" ? Help very much appreciated!

a <- rnorm(10, 0, 1)

ones <- matrix(1, nrow=10, ncol=1 )

A <- a[t(a) %*% ones == 5]

The way I interpret this question, we are asking to generate from the distribution $(X_1, \ldots, X_{10} \mid \sum_{i = 1} ^ {10} X_i = 5)$ where $X_1, \ldots, X_n \stackrel {iid} \sim \mathcal N(0,1)$. This is easy enough - the joint distribution of $(X_1, \ldots, X_9, \sum_i X_i)$ is multivariate normal with mean vector $\mu = (0,\ldots,0)$ and covariance matrix

$$\Sigma = \pmatrix{1 & 0 & \ldots & 0 & 1 \\ 0 & 1 & \ldots & 0 & 1 \\ \vdots & \vdots & \cdots & \vdots & \vdots \\ 0 & 0 & \ldots & 1 & 1 \\ 1 & 1 & \ldots & 1 & 10}.$$

From this, it can be shown that $(X_1, \ldots, X_9 \mid \sum_i X_i)$ is, again, a normal distribution with mean $\mu^\star = (\bar X, \bar X, \ldots, \bar X)$ and covariance matrix $\Sigma^\star = \mathbf I - \frac 1 {10} \mathbf J$ where $\mathbf J$ is a matrix with all entries equal to $1$ - see, for example, here. So, to do this in R,

library(MASS)
Sigma  <- diag(9) - .1
X      <- numeric(10)
X[1:9] <- mvrnorm(1, rep(5 / 10, 9), Sigma)
X[10]  <- 5 - sum(X[1:9])


And now $(X_1, \ldots, X_{10})$ is drawn from the distribution $(X_1, \ldots, X_{10} \mid \sum_i X_i = 5)$.

• +1 for interpreting the question in a constructive manner and providing a good approach.
– whuber
May 12 '14 at 14:23

You can't actually do this. By adding the constraint of sum=5, you're reducing the degrees of freedom by 1 and therefore eliminating true randomness. You can either 1) pick 9 truly random numbers from a distribution, followed by one fixed number to get to a sum of 5, or 2) approximate sum=5 by doing a random sample from a distribution with a fixed mean and small variance, as Peter Flom has already explained.

You can generate 10 identically distributed normal random variables that sum to 5. But they won't be standard normal and they won't be independent.

Generate $X_1 \sim N({1 \over 2},1).$ Let $X_2 = 1-X_1.$ Generate $X_3 \sim N({1 \over 2},1).$ Let $X_4 = 1-X_3.$ And so on, generating 5 pairs of dependent random variables, whose sum is 5.

10 random numbers from a standard normal will have mean 0 and thus the sum will be about 0, not 5. If you ignore the "standard" part, you will still not be able to pick 10 random numbers from a normal distribution whose sum is exactly 5. You can get 10 such numbers with sum very close to 5 by e.g.

set.seed(123)
vars <- rnorm(10, 0.5, 0.001)
sum(vars)
vars


and if you make the sd (0.001) even smaller then the sum will tend to be even closer to 5.

• Oops, my bad. I will correct May 11 '14 at 11:09