Suppose I have $n$ independent normal random variables

$$X_1 \sim \mathrm{N}(\mu_1, \sigma_1^2)\\X_2 \sim \mathrm{N}(\mu_2, \sigma_2^2)\\\vdots\\X_n \sim \mathrm{N}(\mu_n, \sigma_n^2)$$

and $Y=X_1+X_2+\dotsm+X_n$. How would I characterize the density of $Y$ if the distribution of each $X_i$ is each truncated to within $(\mu_i - 2\sigma_i, \mu_i + 2\sigma_i)$? In other words, I'm sampling from $n$ independent normal distributions, discarding samples not within $2\sigma_i$ of each mean, and summing them.

Right now, I'm doing this with the R code below:

x_mu <- c(12, 18, 7)
x_sd <- c(1.5, 2, 0.8)
a <- x_mu - 2 * x_sd
b <- x_mu + 2 * x_sd

samples <- sapply(1:3, function(i) {
  return(rtruncnorm(100000, a[i], b[i], x_mu[i], x_sd[i]))

y <- rowSums(samples)

Is there any method for generating the density of $Y$ directly?

  • 2
    $\begingroup$ Your question implies you know all the $\sigma_i$. Is that really the case or are you estimating them? There's a huge difference! Out of curiosity, why are you throwing away such data? Depending on your objectives, I suspect there exist (much) better procedures. $\endgroup$
    – whuber
    Jul 3, 2014 at 21:14
  • $\begingroup$ I do know all of the means and SDs for my data, yes. $\endgroup$
    – Devin
    Jul 3, 2014 at 21:41
  • 7
    $\begingroup$ I believe that you could characterize it as "a mess". This paper, jstor.org/stable/2236545 , examines the matter, with more scientific rigor. $\endgroup$ Jul 4, 2014 at 3:27
  • 2
    $\begingroup$ Outside approximation via CLT, this is relatively tricky. I guess if $n$ is small enough you could try numerical convolution. $\endgroup$
    – Glen_b
    Jul 5, 2014 at 7:04
  • 2
    $\begingroup$ @Silverfish Depending on implementation, platform and how fine a grid is tolerable, hundreds should be fine (perhaps more); besides speed, though, with enough terms you have to be much more careful about details of implementation or a number of numerical issues can start to crop up. $\endgroup$
    – Glen_b
    Jan 21, 2015 at 15:29

2 Answers 2


You could use approximation by saddlepoint methods, for the sum of truncated normals. I will not give the details now, you can look at my answer to General sum of Gamma distributions for hints. What we need is to find the moment-generating function for a truncated normal, which is easy. I will do it here for a standard normal truncated at $\pm 2$, which has density $$ f(x) =\begin{cases} \frac1{C} \phi(x), & |x| \le 2 \\ 0, & |x| > 2 \end{cases} $$ where $C=\Phi(2) - \Phi(-2)$ here $\phi(x), \Phi(x)$ are density and cdf for a standard normal, respectively.

The moment generating function can be calculated as $$ \DeclareMathOperator{\E}{\mathbb{E}} M(t) = \E e^{tX}=\frac1{C}\int_{-2}^2 e^{tx} \phi(x)\; dx=\frac1{C}e^{\frac12 t^2} [\Phi(2-t)-\Phi(-2-t) ] $$ and then we can use saddlepoint approximations.


I'm curious why, but yes, there is a simple way to generate the pdf of this sum of distributions:

## install.packages("truncnorm")
## install.packages("caTools")

x.mu <- c(12, 18, 7)
x.sd <- c(1.5, 2, 0.8)
x.a <- x.mu - 2*x.sd
x.b <- x.mu + 2*x.sd

dmulti <- function(x, a, b, mu, sd)
             dtruncnorm(x, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)
pmulti <- function(q, a, b, mu, sd)
             ptruncnorm(q, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)

pointrange <- range(c(x.a, x.b))
pointseq <- seq(pointrange[1], pointrange[2], length.out=100)
## Plot the probability density function
plot(pointseq, dmulti(pointseq, x.a, x.b, x.mu, x.sd),

## Plot the cumulative distribution function
plot(pointseq, pmulti(pointseq, x.a, x.b, x.mu, x.sd),
  • $\begingroup$ If I read this code correctly, you seem to be implementing something like a mixture rather than a summation. The plot this code produces is woefully incorrect. It's not even a valid probability density function! $\endgroup$
    – whuber
    Nov 3, 2014 at 18:57
  • $\begingroup$ @whuber, thanks for the catch. I normalized the pdf and added the cdf. $\endgroup$ Nov 8, 2014 at 12:38
  • 3
    $\begingroup$ Thank you. However, the basic error persists: you are computing a mixture distribution rather than the sum. $\endgroup$
    – whuber
    Nov 8, 2014 at 13:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.