3
$\begingroup$

Inspired by the probability calculation here, I am trying to solve a little general one:

$$\mathbb{P} (\sum_{i=1}^{m-1} A_i + \sum_{i=1}^{m} S_i < L < \sum_{i=1}^{m} A_i + \sum_{i=1}^{m+1} S_i)$$

where,

$$A_i \sim \textrm{exp}(\lambda), \; S_i \sim \textrm{exp}(\mu), \; L \sim \textrm{exp}(\lambda), \textrm{ and } \lambda \neq \mu \textrm{ are two integers}.$$

All $A_i, S_i, \textrm{ and } L$ are mutually independent. $m$ is an integer parameter.


In this calculation, I fix $m = 4$. Because $A_i \sim \textrm{exp} (\lambda)$, we have $\sum_{i=1}^{i=3} A_i \sim \textrm{Erlang}(3, \lambda)$.

Probability[
  a3 + s4 < middle < a4 + s5, 
  {a3 \[Distributed] ErlangDistribution[3, λ], 
   s4 \[Distributed] ErlangDistribution[4, μ], 
   a4 \[Distributed] ErlangDistribution[4, λ], 
   s5 \[Distributed] ErlangDistribution[5, μ], 
   middle \[Distributed] ExponentialDistribution[λ]}]

After a few minutes, I got a long result from which I could not read anything useful:

$\frac{\mu ^4 \left(186624 \lambda ^{16}+3825792 \lambda ^{15} \mu +35971776 \lambda ^{14} \mu ^2+205589664 \lambda ^{13} \mu ^3+798109200 \lambda ^{12} \mu ^4+2217713192 \lambda ^{11} \mu ^5+4537481548 \lambda ^{10} \mu ^6+6954729890 \lambda ^9 \mu ^7+8071898695 \lambda ^8 \mu ^8+7133161040 \lambda ^7 \mu ^9+4799247376 \lambda ^6 \mu ^{10}+2441453824 \lambda ^5 \mu ^{11}+923808416 \lambda ^4 \mu ^{12}+252053248 \lambda ^3 \mu ^{13}+46855424 \lambda ^2 \mu ^{14}+5307904 \lambda \mu ^{15}+276224 \mu ^{16}\right)}{11664 (\lambda +\mu )^5 (2 \lambda +\mu )^7 (\lambda +2 \mu )^8}$

In addition, the command FullSimplify[Probability[]] does not help.


Edit (2014-11-17):

By paper-and-pencil, I have obtained a closed form $\frac{\mu + 2 \lambda}{\mu + \lambda} (\frac{1}{2} \frac{\mu}{\mu + \lambda})^{m}$. Notice that the sum $\sum_{m=1}^{\infty} \frac{\mu + 2 \lambda}{\mu + \lambda} (\frac{1}{2} \frac{\mu}{\mu + \lambda})^{m} = \frac{\mu}{\mu + \lambda}.$ The remaining $\frac{\lambda}{\mu + \lambda}$ for special case $m = 0$ has been omitted here.


Therefore, the question is

How to get a more compact form of the above probability? Specifically, is it consistent with my manual calculation (which may be wrong)?

$\endgroup$
3
  • 1
    $\begingroup$ I think your program is wrong. In your program, a3 and a4, s4 and s5 are independent, but actually they're not, a3 is part of a4, s4 is part of s5! $\endgroup$
    – xzczd
    Commented Nov 17, 2014 at 12:23
  • $\begingroup$ @xzczd Thanks. Even "error" is helpful to me. How should I program the non-independent? Is it right to first define each $A_i$ and $S_i$ separately and then sum them to get the Erlang distributions? $\endgroup$
    – hengxin
    Commented Nov 17, 2014 at 12:37
  • $\begingroup$ Network is a bit slow today. See my answer. $\endgroup$
    – xzczd
    Commented Nov 17, 2014 at 13:39

1 Answer 1

4
$\begingroup$

The approach that kguler suggest in your precedent question is completely suitable for the current one:

FindSequenceFunction[Table[With[{
       d1 = HypoexponentialDistribution[Flatten[Table[{λ, μ}, {i - 1}]]],
       d2 = HypoexponentialDistribution[{λ, μ}],
       d3 = ExponentialDistribution[μ], 
       d = ExponentialDistribution[λ]},
      Probability[l1 + l3 < l < l1 + l2 + l3, 
        {l1 \[Distributed] d1, l2 \[Distributed] d2, 
         l3 \[Distributed] d3, l \[Distributed] d}]], {i, 2, 5}]][m - 1] // Simplify //
      AbsoluteTiming

$\left\{108.573200,\frac{2^{-m} (2 \lambda +\mu ) \left(\frac{\mu }{\lambda +\mu }\right)^m}{\lambda +\mu }\right\}$

Of course ErlangDistribution can also be used, and it's faster than the above approach actually:

FindSequenceFunction[
    Table[Probability[
      a + s < middle < a + aadd + s + sadd, 
      {a \[Distributed] ErlangDistribution[i, λ], 
       s \[Distributed] ErlangDistribution[i + 1, μ], 
       aadd \[Distributed] ExponentialDistribution[λ], 
       sadd \[Distributed] ExponentialDistribution[μ], 
       middle \[Distributed] ExponentialDistribution[λ]}], 
       {i, 4}]][m - 1] // Simplify // AbsoluteTiming

$\left\{37.932000,\frac{2^{-m} (2 \lambda +\mu ) \left(\frac{\mu }{\lambda +\mu }\right)^m}{\lambda +\mu }\right\}$

You can use ParallelTable instead of Table to speed up even more.

BTW, there exists a more straight-forward but much slower approach:

With[{i = 2}, Probability[
   Total@Array[a, i - 1] + Total@Array[s, i] <l< Total@Array[a, i] + Total@Array[s, i + 1], 
   Join[{l \[Distributed] ExponentialDistribution[λ]}, 
    a[#] \[Distributed] ExponentialDistribution[λ]&/@Range[i], 
    s[#] \[Distributed] ExponentialDistribution[μ]&/@Range[i + 1]]]] // AbsoluteTiming

$\left\{96.188400,\frac{\mu ^2 (2 \lambda +\mu )}{4 (\lambda +\mu )^3}\right\}$

$\endgroup$
2
  • $\begingroup$ A long way for me to go to learn Mathematica. Thanks a lot. (I have an even worse network these days. I should have accepted your answer yesterday.) $\endgroup$
    – hengxin
    Commented Nov 18, 2014 at 2:13
  • $\begingroup$ @hengxin It's OK to wait 24 hours or even longer before accepting so your question may attract better answers :) . I added a ErlangDistribution approach, have a look. $\endgroup$
    – xzczd
    Commented Nov 20, 2014 at 4:24

Your Answer

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