Post Reply 
(11C) Probability of No Repetitions
01-11-2020, 02:15 PM (This post was last modified: 04-28-2020 11:50 PM by Albert Chan.)
Post: #3
RE: (11C) Probability of No Repetitions
(01-10-2020 12:42 PM)Gamo Wrote:  P = [1- (1/m)][1- (2/m)]....[1- (n-1)/m]

Perhaps it is more useful by adding a (1-0/m) term, we have P = perm(m,n) / mn

Using the same idea to approximate ln(perm(m,n)):

XCas> f := ln(m-x)
XCas> expand(subst(int(f) - f/2 + diff(f)/12, x = m-u))     → -1/(12*u) - ln(u)/2 - u*ln(u) + u
XCaS> g(u) := -1/(12*u) - ln(u)/2 - u*ln(u) + u
XCas> ln_perm(m,n) := g(m-n) - g(m)                                 // integral limit x = 0 to n, u=m-x

ln_perm() can be simplified and more accurate by putting tiny terms up front, and use logp1()

\(\large \ln(perm(m,n)) ≈ {-n\over 12m(m-n)} - [n + (m-n+½)\ln(1 - n/m)]+ n \ln(m)
\)

For ln(P), above last term cancelled out.

XCas> ln_P(m,n) := -n/(12*m*(m-n)) - (n + (m-n+0.5)*ln(1-n/m))
XCas> ln_perm(m,n) := ln_P(m,n) + n*ln(m)

XCas> e^ln_P(365, 4)       → 0.9836440875
XCas> e^ln_P(365,23)      → 0.4927027657
XCas> e^ln_P(365,48)      → 0.03940202712

Trivia: ln_perm(m,n) = ln_perm(m,a) + ln_perm(m-a, n-a), same as real permutation function.

if m ≈ n, ln_perm(m,n) is not very accurate.
With above trivia, splitting ln_perm(m,n) does not help.
However, we can replace the 2nd term with real ln(perm(m-a,n-a))

Example: ln(perm(69,68)) = ln(69!) ≈ 226.1905483
XCas> ln_perm(69,68)                   → 226.1882765 // slightly worse than log of Stirling's formula.
XCas> ln_perm(68,67) + ln(2)        → 226.1902224 // perm(69-67, 68-67) = perm(2,1) = 2
XCas> ln_perm(68,66) + ln(6)        → 226.1904485 // perm(69-66, 68-66) = perm(3,2) = 6

Edit: above ln_P formula is numerically very bad, use this instead
Find all posts by this user
Quote this message in a reply
Post Reply 


Messages In This Thread
(11C) Probability of No Repetitions - Gamo - 01-10-2020, 12:42 PM
RE: (11C) Probability of No Repetitions - Albert Chan - 01-11-2020 02:15 PM



User(s) browsing this thread: 1 Guest(s)