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$
  • 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 Nov 17 '14 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 Nov 17 '14 at 12:37
  • $\begingroup$ Network is a bit slow today. See my answer. $\endgroup$ – xzczd Nov 17 '14 at 13:39
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\}$

| improve this answer | |
$\endgroup$
  • $\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 Nov 18 '14 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 Nov 20 '14 at 4:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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