# Iterating distributions

I want to calculate the distribution that comes from summing up random numbers between -0.5 and 0.5.

A computationally intensive version of doing this is to simply do the following:

SmoothHistogram[
Table[Table[
Sum[RandomReal[] - 0.5, {i, 1,


This won't work however if the number of numbers to be summed is very large, like a 1000000 for example.

So I tried doing it by first calculating the distribution of one random number, then calculating the distribution of two random numbers from the first distribution. The first distribution is simply:

UnitBox[x]


The next distribution can be calculated by integrating the previous distribution from x-0.5 to x+0.5:

Integrate[UnitBox[a], {a, x - 1/2, x + 1/2}]


I can't figure out how to iterate this to, for example, draw a plot from the result.

I also tried the following:

NestList[PDF[
ProbabilityDistribution[
Integrate[# /. x -> a, {a, x - 1/2, x + 1/2},
Assumptions -> x \[Element] Reals], {x, -2, 2}], x] &,
UnitBox[x], 5]
Plot[%, {x, -2, 2}]


• Any reason not to use built-ins? – ciao Apr 19 '15 at 3:26
• No. What built-in would help? – user Apr 19 '15 at 3:27
• maybe UniformSumDistribution. E.g., dis[n_] := UniformSumDistribution[n, {-.5, .5}]; Plot[Evaluate[PDF[dis@#, x] & /@ Range[10]], {x, -2, 2}]? – kglr Apr 19 '15 at 3:30
• It works! Thanks a lot. – user Apr 19 '15 at 3:32
• @MathematicaUser39386, posted the comment as an answer. – kglr Apr 19 '15 at 4:03

• summing random numbers between -0.5 and 0.5.

• the number of numbers to be summed is very large, like a 1000000 for example.

The sum of $n$ identical Uniform random variables is known as a generalised Irwin-Hall distribution, implemented in Mathematica as the UniformSumDistribution [ see kguler's answer]. The latter takes an $n$-part piecewise form, so, for example, if $n = 10$, the pdf is:

PDF[UniformSumDistribution[10, {-1/2, 1/2}], x]


The above is exact and works beautifully for small $n$. However, the OP poses the problem of $n$ being very large -- such as $n$ = 1 million -- which involves creating a 1 million part piecewise structure which is certainly going to fail. In fact, for $n = 1000$, simply making a plot of the pdf takes about 400 seconds on my Mac Pro:

(AA = Plot[
PDF[UniformSumDistribution[1000, {-1/2, 1/2}], x], {x, -50, 50},
PlotRange -> All]) // AbsoluteTiming


390 seconds

• For $n=10,000$, the above plot is simply unworkable ... which leaves the OP's question: how to proceed for large $n$??

Large $n$: apply Central Limit Theorem

For large $n$, I would suggest applying the (Lindeberg-Lévy) Central Limit Theorem:

• If the random variables $(X_1, X_2, \dots)$ are independent and identically distributed, each with finite mean $\mu$ and finite variance $\sigma^2$, then the sample sum:

$$S_n \overset{a}{\sim} N\left(n \mu,n \sigma ^2\right)$$

In our case, $X \sim \text{Uniform}(-\frac12,\frac12)$ so $\mu = 0$ and $\sigma^2 = \frac{1}{12}$, so the asymptotic distribution of the sample sum is:

$$S_n \overset{a}{\sim} N\left(0, \frac{n}{12} \right)$$

with pdf $f(x)$:

       f = Exp[-6 (x^2/n)] / Sqrt[n Pi /6];
domain[f] = {x, -Infinity, Infinity} && {n > 0};


Here is a plot of pdf $f(x)$ when $n = 1000$:

BB = Plot[f /. n -> 1000, {x, -50, 50}, PlotRange -> All, PlotStyle -> Red]


This is, of course, instantaneous and works for arbitrarily large $n$. The following plot compares the exact solution AA to the asymptotic solution BB when $n = 1000$:

Show[BB, AA]


There is no discernible visual difference between the two plots here. Whereas the exact AA solution fails to evaluate for very large $n$, the asymptotic BB solution will always evaluate immediately, and with ever improving accuracy as $n$ increases.

The built-in function UniformSumDistribution may be useful:

usd[n_] := UniformSumDistribution[n, {-.5, .5}];
Plot[Evaluate[PDF[usd@#, x] & /@ Range[20]], {x, -2, 2},
PlotStyle -> (ColorData[{"Rainbow", {1, 20}}] /@ Range[20]),
Exclusions -> None, PlotRange -> All, ImageSize -> 500,
PlotLegends -> ("usd (" <> ToString[#] <> ")" & /@ Range[20])]


ClearAll[skd]
skd[n_] := SmoothKernelDistribution[RandomVariate[usd[n], 5000]];
Plot[Evaluate[PDF[#, x] & /@ {skd[5], usd[5], skd[3], usd[3]}], {x, -2, 2},
ImageSize -> 500, Exclusions -> None,
Filling -> {1 -> {{2}, {Orange, Yellow}}, 3 -> {{4}, {Blue, Green}}}]