0
$\begingroup$

It is known that X1, X2,...,Xn are random samples from the population X, where the probability of X taking 0 or 1 is equal, that is, $P(X=0)=P(X=1)=\frac{1}{2}$.

I want to find an approximation of $P\left(\sum_{i=1}^{100} X_{i} \leq 55\right)$ by means of the central limit theorem.

I use the following code to calculate this problem, but MMA keeps running and can't get the result:

NProbability[Sum[x[i], {i, 1, 100}] <= 55, 
 Table[x[i] \[Distributed] 
   EmpiricalDistribution[{0.5, 0.5} -> {0, 1.}], {i, 1, 100}]]

Is there any way to quickly find the approximate value of this problem?

Supplementary mathematical analysis process:

$$E\left(\sum_{i=1}^{100} X_{i}\right) X=100 E X=50 . \quad D\left(\sum_{i=1}^{100} X_{i}\right)=100 D X=25$$

$$\begin{aligned} &\text { According to the central limit theorem } \sum_{i=1}^{100} X_{i} \sim N(50,25)\\ &\therefore P\left\{\sum_{i=1}^{100} X_{i} \leq 55\right\}=P\left\{\frac{\sum_{i=1}^{100} X_{i}-55}{5} \leq \frac{55-50}{5}\right\}=\Phi(1) \end{aligned} $$

Count[Table[Total[RandomChoice[{0.5, 0.5} -> {0, 1}, 100]], 100000], 
  u_ /; u <= 55]/100000.
(*0.8655*)
Erf[1.]
(*0.84270079295*)
Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ This is basically like modelling the number of heads in a sequence of coin flips. It's a BinomialDistribution[100,1/2] and you're counting the number of successful trials: Probability[x <= 55, x \[Distributed] BinomialDistribution[100, 1/2]] which gives precisely 68482723177360620218041365161 / 79228162514264337593543950336 or roughly 0.864373 . If you want to demonstrate convergence to the normal distribution, then N@Probability[x <= 55 + .5, x \[Distributed] NormalDistribution[50, 5]] - the .5 is a correction, and stddev of binomial is Sqrt[n (1 - p) p] and mean is np for n=100,p=1/2 $\endgroup$
    – flinty
    Commented Aug 3, 2020 at 1:01
  • $\begingroup$ @flinty If I change the title to read: NProbability[Sum[x[i], {i, 1, 100}] <= 55, Table[x[i] \[Distributed] EmpiricalDistribution[{0.25, 0.5, 0.25} -> {0, 1., 2.}], {i, 1, 100}]], what should you do? $\endgroup$
    – A little mouse on the pampas
    Commented Aug 3, 2020 at 1:05
  • 1
    $\begingroup$ NProbability[Total[Array[x, 100]] <= 55, # \[Distributed] dist & /@ Array[x, 100]] will not complete for 100, but at that level I would just assume an almost-normal distribution anyway and use the limiting distribution: distapprox=NormalDistribution[Mean[dist]*100,StandardDeviation[dist]*Sqrt[100]] then Probability[x<=55+0.5,x \[Distributed] distapprox]] which gives 1.55443*10^-10 - and that low probability makes sense because Total@RandomVariate[dist, 100] is normally very close to 100. $\endgroup$
    – flinty
    Commented Aug 3, 2020 at 1:19
  • 1
    $\begingroup$ You can also do d2 = TransformedDistribution[ Total[Array[x, 100]], # \[Distributed] dist & /@ Array[x, 100]]; and get mu = Mean[d2] and sd = StandardDeviation[d2] then just plug mu*100 and sd*Sqrt[100] into NormalDistribution later. $\endgroup$
    – flinty
    Commented Aug 3, 2020 at 1:21
  • 3
    $\begingroup$ Yeah that will never work - don't use NProbability with a sum this big because it's like doing very high dimensional intergration - use the normal distribution approximation instead. $\endgroup$
    – flinty
    Commented Aug 3, 2020 at 1:36

3 Answers 3

Reset to default
2
$\begingroup$

Given a discrete distribution $D$ of finite variance, you only need the mean and standard deviation to apply the central limit theorem here:

dist = BinomialDistribution[100,1/2]
mu = Mean[dist];
sd = StandardDeviation[dist];

nclt = NormalDistribution[mu, sd];

It is important to apply a continuity correction, since we've gone from a discrete to continuous distribution, which is why I add 0.5 below:

Probability[x <= 55 + 0.5, x \[Distributed] nclt]
(* 0.864334 *)

Probability[x <= 55, x \[Distributed] BinomialDistribution[100, 1/2]]
(* 68482723177360620218041365161 / 79228162514264337593543950336 *)
(* approx 0.864373 *)

You can apply this reasoning for other distributions, but you must remember to scale the mean by a factor of $n$ and the standard deviation by a factor of $\sqrt{n}$:

dist2 = EmpiricalDistribution[{0.25, 0.5, 0.25} -> {0, 1, 2}];
mu2 = Mean[dist2];
sd2 = StandardDeviation[dist2];

nclt2 = NormalDistribution[mu2 * 100, sd2 * Sqrt[100]];
Probability[x <= 55 + 0.5, x \[Distributed] nclt2]
(* 1.55443*10^-10 *)
Improve this answer
$\endgroup$
4
$\begingroup$

I wonder if you went to the wrong solution because you used an inefficient approach to calculate the sum. Either of the following work almost instantaneously:

AbsoluteTiming[Sum[PDF[BinomialDistribution[100, 1/2], x], {x, 0, 55}]]
(* {0.0065887, 68482723177360620218041365161/79228162514264337593543950336} *)

AbsoluteTiming[CDF[BinomialDistribution[100, 1/2], 55]
(* {0.0001947, 68482723177360620218041365161/79228162514264337593543950336} *)
Improve this answer
$\endgroup$
1
$\begingroup$

Direct method to do what the OP was trying to do, but faster:

NProbability[
  Total[Array[x, 100]] <= 55,
  Array[x, 100] \[Distributed] ProductDistribution[{BernoulliDistribution[1/2], 100}]
]

0.864673

This is still not a very scaleable method, though.

Improve this answer
$\endgroup$

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.