Random proababilty of rolling down to 1 from 1000-CodePudding

If I roll a random number between 1 and a 1000, then use this random number as the new range from 1 and new number for example if it rolled 500 the new range would be 1 to 500. What is the mean attempts to get down to 1?

CodePudding user response：

The means satisfy the following recurrence: if f(n) denotes the mean number of steps to get to 1 from n, then f(1) = 0 and for n>1,

f(n) = 1   (f(1)   ...   f(n-1)   f(n))/n -->
f(n) = n/(n-1)   [f(1)   ...   f(n-1)]/(n-1).

We can find f(1000) with the following script:

import numpy as np

means = [0]
for n in range(2,1001):
    means.append(n/(n-1) np.average(means))
print(means[-1])

The result is 8.484470860550345. I'm not sure if there is any nice closed form for f(n).

CodePudding user response：

I get 8.48 for an answer. The logic is as follows: let f(n) denote the mean number of rows to get down to 1. Then, f satisfies the following recurrence:

f(1) = 1
f(n) = 1/n   (n-1)/n * (1   sum(f(i) for i in range(2, n 1)) / (n - 1))

For general n, there is 1/n chance of getting n one this turn, and (n-1)/n chance of one after this turn. In the latter case, the expected number of turns is one plus a weighted average. Rearranging the terms gives the following formula:

from functools import lru_cache

@lru_cache(None)
def f(n):
  if n == 1:
    return 1
  return n / (n - 1)   1 / (n - 1) * sum(f(i) for i in range(2, n))

You could also use fractions module to compute the answer exactly, although the resulting fraction is quite large.

After some fiddling, it seems that Mathematica can find a "closed form" solution in terms of polygamma function with

f[n_] := 1   EulerGamma   PolyGamma[0, n]

One way to derive this is to re-express the recurrence in terms of y(n) = sum(f(k) for k in range(1, n)). This also shows that f(n) is O(log(n)) asymptotically.

CodePudding user response：

First 100 digits:

8.484470860550344912656518204333900176521679169708803665773626749957699349165202440959934437411845081

The exact fraction:

60484649691685213005704263395601539057400162783408267244512899237619438652990632145483065609746993406733894830753023779634909455921438397706516749459866848443713860300689848494770432824401693665503983189312845928714471628150308917599333411703457439105339322525040288258228358123705825959020674166771023657949375130672985808604458424810309577532641272807764258202248144159846532965261433911670727905231051075533636840438428488803438997
/
7128865274665093053166384155714272920668358861885893040452001991154324087581111499476444151913871586911717817019575256512980264067621009251465871004305131072686268143200196609974862745937188343705015434452523739745298963145674982128236956232823794011068809262317708861979540791247754558049326475737829923352751796735248042463638051137034331214781746850878453485678021888075373249921995672056932029099390891687487672697950931603520000

Came up with the same formula as Ben, but used the fractions module to get the exact result:

from fractions import Fraction

attempts = [None, 0]
for hi in range(2, 1001):
    attempts.append(Fraction(hi, hi-1)   Fraction(sum(attempts[1:]), hi-1))

# Print the first 100 digits
a = attempts[-1]
print(a.numerator * 10**99 // a.denominator)

# Print the fraction
print(a)

CodePudding user response：

Tthe mean attempts to get “the number” down to one? So like EXP[len((1,1000), …, (1,1))]?

Possibly something like

from random import randint
from statistics import mean

trials = list()
while len(trials) < 30:
    value, attempts = 1000, 0
    while value != 1:
        value = randint(1, value)
        attempts  = 1
    trials.append(attempts)
print(mean(trials))

CodePudding user response：

I found ~8.5 with this code:

import random

def f():
    i = 0
    n = 0
    r = 1000
    while n != 1:
        n = random.randint(1, r)
        i  = 1
        r = n
    return i

N = 10000
l = [f() for _ in range(N)]
print(sum(l) / len(l))

Output:

8.5206

CodePudding user response：

we figured it out, I forgot to add a 1 to the number in my randrange in the while loop, I now get the value of 8.48....... with the code below, Thanks everyone who helped answer this question :).

import random
import math
random.seed(13)
import statistics

tries = []
for i in range(1000000):
  number = random.randrange(1, 1001)
  attempts = 1
  while number != 1:
    number = random.randrange(1, number 1)
    attempts  = 1
  tries.append(attempts)

print(statistics.mean(tries))