Expected value -- infinite power series sum?

Discussion in 'Physics & Math' started by Jennifer Murphy, Aug 22, 2014.

  1. Jennifer Murphy Registered Senior Member

    Messages:
    239
    In another thread about weighted averages, we have the odds that the four aces will be dealt in suits order either forward (SHDC) or backwards (CDHS) as 2/24 = 0.08333. Call that P.

    Code:
    P = 2/24
    
    Now then, the odds that the aces will be dealt in one of those two orders on the Nth deal and not before, is

    Code:
    Pn = (1-P)^(N-1) x P
    
    That is, it is the odds that it won't happen for N-1 deals and then it will happen on the Nth deal. Correct so far?

    Now I'd like to calculate the expected number of deals to get the aces in order. I believe the expected value, E, of the random variable Pn is the infinite sum

    Code:
    E = 1x((1-P))^0xP + 2x((1-P))^1xP + 3x((1-P))^2xP + 4x((1-P))^3xP + ...
    
    From this Wolfram page http://mathworld.wolfram.com/PowerSum.html, I found a power series sum that's close to what I need:

    Code:
    Sum [as k=0 to N] of (kx^k) = [P - (N+1)xP^(N+1) + NxP^(N+2)] / ((P-1)^2)
    
    What I think I need is

    Code:
    Sum [as k=0 to N] of ((k+1)x^k)
    
    Is there a way to convert the first one into the second one?

    Is this even the right way to go about this?

    Thanks
     
  2. Google AdSense Guest Advertisement



    to hide all adverts.
  3. przyk squishy Valued Senior Member

    Messages:
    3,203
    Yes, because

    \(\sum_{k = 0}^{N} (k \,+\, 1) x^{k} \,=\, \sum_{k = 0}^{N} k x^{k} \,+\, \sum_{k = 0}^{N} x^{k} \,,\)​

    and you already know how to solve the geometric series.

    Another way is to use that \((k \,+\, 1) x^{k} \,=\, \frac{\mathrm{d}}{\mathrm{d} x} x^{k+1}\) and pull the derivative out of the sum:

    \( \begin{eqnarray} E &=& P \sum_{k = 0}^{\infty} (k \,+\, 1) Q^{k} \\ &=& P \sum_{k = 0}^{\infty} \frac{\mathrm{d}}{\mathrm{d}Q} \, Q^{k + 1} \\ &=& P \, \frac{\mathrm{d}}{\mathrm{d}Q} \sum_{k = 0}^{\infty} Q^{k + 1} \\ &=& P \, \frac{\mathrm{d}}{\mathrm{d}Q} \, \frac{Q}{1 - Q} \\ &=& P \, \frac{1}{(1 - Q)^{2}} \\ &=& \frac{1}{P} \,, \end{eqnarray} \)​

    using \(P \,=\, 1 \,-\, Q\) like in your previous thread.
     
  4. Google AdSense Guest Advertisement



    to hide all adverts.
  5. Jennifer Murphy Registered Senior Member

    Messages:
    239
    Thanks. I'll work on that.

    Is that LaTeX you are using to generate those equations?
     
  6. Google AdSense Guest Advertisement



    to hide all adverts.
  7. Aqueous Id flat Earth skeptic Valued Senior Member

    Messages:
    6,152
    Hi there. In this case, the odds that any 4 particular cards will be drawn is (1/52) * (1/51) * (1/50) * (1/49) = 1.54 e-7.
     
  8. Jennifer Murphy Registered Senior Member

    Messages:
    239
    Yes, but I'm only interested in the aces and only in the order that they appear. I don't care if they are the first 4 cards, the last four cards, the 9th, 21st, 39th, and 40th, or any other positions.

    To simplify the problem, I could use a deck with just the 4 aces -- or just generate a random ordering of the numbers 1-2-3-4. There are 24 possible orders for these 4 cards: 1234 1243 1324 1342 etc.

    I am only interested in two orders: 1234 and 4321. There is a 2/24 or 0.08333 chance of one of these occurring on any one deal. I would like to calculate the expected number of deals to get one of those sequences.

    Is that clearer?
     
    Last edited: Aug 23, 2014
  9. przyk squishy Valued Senior Member

    Messages:
    3,203
    Yes. Some LaTeX math support was added to this forum several years ago. Anything you write between [noparse]\(...\)[/noparse] tags gets interpreted as LaTeX code. For instance, [noparse]\(P_{n} = (1-P)^{N-1} P\)[/noparse] renders as \(P_{n} = (1-P)^{N-1} P\).
     
  10. Jennifer Murphy Registered Senior Member

    Messages:
    239
    Thanks, once again, for pointing out what should have been obvious to me.

    Please Register or Log in to view the hidden image!

    I actually do know that \((k+1)x^k = kx^k + x^k\).

    For my application, I need
    \( E=\sum_{k = 1}^{\infty} kP(1-P)^{(k-1)} \)​

    On Wikipedia, http://en.wikipedia.org/wiki/List_of_mathematical_series, I found these infinite series for values of |z| < 1, which should work for for my application.
    \( \sum_{k=1}^{\infty} z^{k} = \frac z {1-z} \text{ ; |z|<1} \sum_{k=1}^{\infty} kz^{k} = \frac z {(1-z)^2} \text{ ; |z|<1} \)​

    Using your tip, I can adapt this to my situation with a little tweaking:
    \( \begin{eqnarray} E &=& \sum_{k = 1}^{\infty} kP(1-P)^{(k-1)} &=& \sum_{k = 1}^{\infty} kPQ^{(k-1)} &=& \sum_{k = 1}^{\infty} kP \frac {Q^k}{Q} &=& \frac P Q \sum_{k = 1}^{\infty} k Q^k &=& \frac P Q \frac Q {(1-Q)^2} &=& \frac P Q \frac Q {P^2} &=& \frac 1 P \end{eqnarray} \)​

    Is this correct?

    If it is, then E=12.

    Intuitively, such a simple formula resulting in a nice round number is a surprise, but I checked the math several times and can't find any errors.


    PS: Thanks for the little LaTeX tutorial. I think I starting to get the hang of it.

    Please Register or Log in to view the hidden image!

     
  11. rpenner Fully Wired Valued Senior Member

    Messages:
    4,833
    \( \sum_{k = 1}^{n} kP(1-P)^{(k-1)} = \frac{1 - (1 + n P) (1-P)^n}{P} ; \quad \quad \forall n \in \mathbb{n} \)

    Proof:
    Let n = 1, then
    \(\sum_{k = 1}^{n} kP(1-P)^{(k-1)} = P(1-P)^0 = P = \frac{P^2}{P} = \frac{1 - ( 1 - P^2) }{P} = \frac{1 - (1 + P) (1-P)}{P} = \frac{1 - (1 + n P) (1-P)^n}{P} ; \quad \quad n = 1\)​
    Assume the formula is true for all integers 1..m, then we can use that to prove the case for \(n=m+1\):
    \(\sum_{k = 1}^{n} kP(1-P)^{(k-1)} = \sum_{k = 1}^{m} kP(1-P)^{(k-1)} + (m+1)P(1-P)^{(m+1-1)} = \frac{1 - (1 + m P) (1-P)^m}{P} + \frac{(m+1)P^2(1-P)^{m}}{P} = \frac{ 1 - ( 1 + m P - (m+1)P^2 ) (1-P)^{m}}{P} = \frac{ 1 - ( 1 + m P + P ) (1 - P) (1-P)^{m}}{P} = \frac{ 1 - ( 1 + (m +1) P ) (1-P)^{m+1}}{P} = \frac{1 - (1 + n P) (1-P)^n}{P} ; \quad \quad n = m+1\)​
    Since we proved this for n=1 and proved that if it was true for n=m then it is also true for n=m+1, then it follow that it is true for all positive integers. QED.​

    \(E=\sum_{k = 1}^{\infty} kP(1-P)^{(k-1)} = \lim_{n\to\infty} \sum_{k = 1}^{n} kP(1-P)^{(k-1)} = \lim_{n\to\infty} \frac{1 - (1 + n P) (1-P)^n}{P} = \frac{1}{P} - \lim_{n\to\infty} \frac{1 + n P}{P} (1-P)^n = \frac{1}{P} ; \quad \quad 0 \lt P \lt 2\)
     
  12. Jennifer Murphy Registered Senior Member

    Messages:
    239
    Imagine that, I was actually right about something.

    I tried it your way, but when I got to the limit, I wasn't sure about a term like \({nP}P^n\). I thought that the \({P^n}\) factor, being a power, would go to zero faster than the \(nP\) factor, being just a product, would go to infinity, but that was just a hunch.

    Thanks for the confirmation.
     
  13. Aqueous Id flat Earth skeptic Valued Senior Member

    Messages:
    6,152
    My bad. I see I spoke too soon. I thought you were asking the odds of drawing four particular cards from the top of the deck.

    OK so if you have a deck of four cards, the probability of drawing any particular sequence is (1/4) * (1/3) * (1/2) = 1/24. Got it.
    Right.
    Ok.
    P(1234 or 4321) = P(1234) + P(4321) = 1/24 + 1/24 = 2/24. That's the odds of drawing either sequence from the deck of four aces.

    I think you meant to say 2/24. 2/14 is the odds of drawing either of two particular aces from the top of a full deck.

    From a full deck or from a deck of four aces? Let me do it for the full deck and you can tell me if I'm following you (since I missed the earlier thread). I'll go with this simpler way of numbering the ordered deck 1 through 52.

    The odds of drawing 1234 or 4321 on the first four cards drawn from the full shuffled deck is 2 * (1/52)(1/51)(1/50)(1/49). The odds of drawing them from the 2nd through 5th cards dealt is 2 * (1/51)(1/50)(1/49)(1/48). The odds of drawing one of the sequences beginning on the Nth card dealt is

    \({ p }_{ N }\quad =\quad \prod _{ k=0 }^{ 3 }{ \frac { 1 }{ 52-k-(N-1) } } ,\quad 1\le N\le 49.\)

    (BTW I use Google Chrome with the Daum LaTex equation editor. I prefer rendering expressions in another tab graphically then pasting the code to my post - I just have to enclose it in the [tex}[/tex] tags as przk mentioned.)

    You are asking about the expected value of N. That's the same as asking the average of the probabilities for each of the 49 possible positions in the deck in which the desired sequence can be dealt:

    \(E\{ N\} \quad =\quad \frac { 1 }{ 49 } \sum _{ i\quad =\quad 1 }^{ 49 }{ { p }_{ i } } =\quad \frac { 1 }{ 49 } \sum _{ i\quad =\quad 1 }^{ 49 }{ \prod _{ k=0 }^{ 3 }{ \frac { 1 }{ 53-k-i } } }\)

    I know you do spreadsheets. Here's a possible solution doing it that way:

    1. Cell A1 = 1
    2. Cell A2 = A1 + 1
    3. Copy (2) and paste to cells A3:A49. (i.e., write 1 though 49 in column A)
    4. Cell B1 = 1 / (53 - A1)
    5. Cell C1 = 1 / (52 - A1)
    6. Cell D1 = 1 / (51 - A1)
    7. Cell E1 = 1 / (50 - A1)
    8. Cell F1 = B1 * C1 * D1 * E1
    9. Copy B1:F1 and paste to B2:F49
    10 Cell F49 = 1. (The 52nd card dealt is known.)
    11. Cell F50 = SUM(F1:F49) / 49.

    I get 0.06597. That's the odds of getting one of the sequences (1234 or 4321). The odds of getting either of them is twice that much, or about 13.2%.

    Therefore the odds of never getting either sequence after dealing all the cards is (1-P) or 86.8%.

    That's one possible brute force approach (unless I goofed somewhere). Did I model the problem you had in mind correctly?

    I am only halfway following the discussion here. I'm of the impression that the preceding was a sidebar, and that the original question is still open.

    More importantly, when are we going to Vegas? If you'll liquidate everything and bring the playing money, and if przk & rpenner will bring their masterful skills to clean out all the dealers, I will certainly be glad to supervise the whole operation for . . . say . . . 86.8% of the take? That will still leave you with 13.2% of a trillion to split with przk & rpenner - some serious change, and just a taste of what you guys can clear after I book us at Monaco (after I deduct expenses from your share, as any diligent agent would do -- so, naturally, you could enjoy the tax deduction on the higher reported earnings. I would of course have you guys' best interests at heart as your agent. Oh you'll need to sign durable powers of attorney. But we can take that to another thread or hash it out on the plane).

    Please Register or Log in to view the hidden image!



    Hope I didn't derail this, if I simply have not gotten on page with you folks. If so, I at least had fun taking my spreadsheet and LaTex apps on a joyride and imagining how your collective resources might be applied to my retirement in Monaco. :mufc:
     
  14. Aqueous Id flat Earth skeptic Valued Senior Member

    Messages:
    6,152
    Wait: what did I just say? The expected number of deals is 0.06597? Ok I retract that. What did I just do? Oops.

    And 2/14 isn't the odds of drawing either of two particular aces from the top of a full deck! That's 2/52. Uh... well, whatever you were thinking. Maybe you just meant 2/24.

    Funny these thoughts didn't cross my mind BEFORE I pushed the Submit button. I could have just edited the post and said "I misunderstood you" but my mind is scattered today so I'll leave for anyone who wants to set me straight.
     
  15. Jennifer Murphy Registered Senior Member

    Messages:
    239
    We are still not on the same problem. One more try.

    1. A deck is a standard poker or bridge deck of 52 cards. It could also be just the 4 aces from such a deck.
    2. A deal is an event that entails shuffling the cards and then dealing them out until all 4 aces have been dealt. I don't care about the other cards, if any.
    3. Suit order is the standard bridge order of the suits: S(pades) H(earts) D(iamonds) C(lubs). It could be either foreward (SHDC) or backwards (CDHS).
    4. A run is a sequence of deals ending when the 4 aces appear in suit order.
    5. P = The probability that the 4 aces will be dealt in suit order on any one deal = 2/24 = 0.083333.
    6. Q = The probability that the 4 aces will not be dealt in suit order on any one deal = 1-P.
    7. Pn = The probability that the 4 aces will be dealt in suit order on Deal #N, and not before.

    Pn = the probability that the aces will be in order on Deal N (P) times the probability that they will not be dealt in order on the previous N-1 deals (\(Q^{N-1}\)).

    \( P_n = P Q^{N-1} = P (1-P)^{N-1) \)​

    I was trying to calculate E, the expected value of N, the average number of deals for the aces to be in order. It turns out to be 1/P = 12.


    And yes, I meant 2/24, not 2/14. I just fixed that.
     
  16. Aqueous Id flat Earth skeptic Valued Senior Member

    Messages:
    6,152
    Yes the expected number of tries should be 1/P. Apologies for losing my mind as I was walking through this.
     
  17. Jennifer Murphy Registered Senior Member

    Messages:
    239
    At least you found it, I'm still looking for mine...

    I ran a little simulation to test the results.

    Code:
       Hits    Deals     Hit%  Ave#Deals
        834   10,000   8.3400    11.9904
      4,327   50,000   8.6540    11.5554
      8,623  100,000   8.6230    11.5969
     17,086  200,000   8.5430    11.7055
     25,422  300,000   8.4740    11.8008
     33,810  400,000   8.4525    11.8308
     42,160  500,000   8.4320    11.8596
     50,539  600,000   8.4232    11.8720
    
    I'm a little surprised that the Hit% is not a lot closer to 8.3333 and the Ave#Deals is not a lot closer to 12.0000.

    Here's the code:

    Code:
    Hdr = right("Hits",7)||right("Deals",9)||right("Hit%",9)||right("Ave#Deals",11)
    Call LineOut DatFid, Hdr
    
    #deals = 0
    #hits  = 0
    LoopLimit = 10000
    Do i=1                   /* Display loop */
      Do j = 1 to LoopLimit    /* Tally loop */
        #deals=#deals+1          /* Count the deals */
        If RandInt(1,4)<=2 then  /* S or C? (2 chances in 4) */
          If RandInt(1,3)=1 then   /* SH or CD? (1 chance in 3) */
            If RandInt(1,2)=1 then   /* SHD or CDH? (1 chance in 2) */
              #hits=#hits+1            /* Aces in order, count it */
      End j
      Line = fmt(#hits,7)||fmt(#deals,9)||fmt(100*#hits/#deals,9,4)||fmt(#deals/#hits,11,4)
      Call LineOut DatFid, Line
    End i
    
     
  18. RJBeery Natural Philosopher Valued Senior Member

    Messages:
    4,222
    That is surprising; maybe a poor PRNG? Try making a simple loop of RandInt(1,4) for, say, 600,000 iterations and tallying the results?
     
  19. Jennifer Murphy Registered Senior Member

    Messages:
    239
    Yes, it was a poor PRNG, or, more accurately, a poor front end I wrote for one.

    I upgraded to a better PRNG and fixed my front end. Then I did as you suggested. Here is the distribution of the four random integers 1,2,3,4:

    Code:
         Trials      Tally 1s   Err1  Err1%      Tally 2s   Err2  Err2%      Tally 3s   Err3  Err3%      Tally 4s   Err4  Err4%
         10,000         2,495     -5 -0.20%         2,484    -16 -0.64%         2,512    +12 +0.48%         2,509     +9 +0.36%
         50,000        12,462    -38 -0.30%        12,517    +17 +0.14%        12,569    +69 +0.55%        12,452    -48 -0.38%
        100,000        25,022    +22 +0.09%        25,065    +65 +0.26%        25,189   +189 +0.76%        24,724   -276 -1.10%
        500,000       125,029    +29 +0.02%       125,233   +233 +0.19%       125,385   +385 +0.31%       124,353   -647 -0.52%
      1,000,000       250,002     +2 +0.00%       249,753   -247 -0.10%       250,831   +831 +0.33%       249,414   -586 -0.23%
      2,000,000       499,567   -433 -0.09%       499,458   -542 -0.11%       500,898   +898 +0.18%       500,077    +77 +0.02%
      3,000,000       750,416   +416 +0.06%       749,117   -883 -0.12%       750,782   +782 +0.10%       749,685   -315 -0.04%
      4,000,000     1,000,437   +437 +0.04%       998,671 -1,329 -0.13%     1,001,443 +1,443 +0.14%       999,449   -551 -0.06%
      5,000,000     1,250,646   +646 +0.05%     1,248,938 -1,062 -0.08%     1,251,115 +1,115 +0.09%     1,249,301   -699 -0.06%
      6,000,000     1,499,652   -348 -0.02%     1,499,668   -332 -0.02%     1,501,185 +1,185 +0.08%     1,499,495   -505 -0.03%
      7,000,000     1,749,926    -74 -0.00%     1,749,154   -846 -0.05%     1,751,349 +1,349 +0.08%     1,749,571   -429 -0.02%
      8,000,000     2,000,726   +726 +0.04%     1,997,988 -2,012 -0.10%     2,001,500 +1,500 +0.08%     1,999,786   -214 -0.01%
      9,000,000     2,251,233 +1,233 +0.05%     2,247,343 -2,657 -0.12%     2,251,940 +1,940 +0.09%     2,249,484   -516 -0.02%
     10,000,000     2,501,385 +1,385 +0.06%     2,497,044 -2,956 -0.12%     2,501,997 +1,997 +0.08%     2,499,574   -426 -0.02%
    
    Then I reran my previous simulation. Here are those results:

    Code:
         Hits      Deals     Hit%  Ave#Deals
        8,370    100,000   8.3700    11.9474
       16,662    200,000   8.3310    12.0034
       25,080    300,000   8.3600    11.9617
       33,359    400,000   8.3398    11.9908
       41,841    500,000   8.3682    11.9500
       83,648  1,000,000   8.3648    11.9549
      167,079  2,000,000   8.3540    11.9704
      250,664  3,000,000   8.3555    11.9682
      334,156  4,000,000   8.3539    11.9705
      417,190  5,000,000   8.3438    11.9849
      500,619  6,000,000   8.3437    11.9852
      584,251  7,000,000   8.3464    11.9812
      667,661  8,000,000   8.3458    11.9821
      751,380  9,000,000   8.3487    11.9780
      834,712 10,000,000   8.3471    11.9802
    1,669,170 20,000,000   8.3459    11.9820
    2,501,564 30,000,000   8.3385    11.9925
    
    That's a lot closer. Is it close enough?
     

Share This Page