The Purge Riddler

Another week, another Riddler. This weeks:

A town of 1,000 households has a strange law intended to prevent wealth-hoarding. On January 1 of every year, each household robs one other household, selected at random, moving all of that house’s money into their own house. The order in which the robberies take place is also random and is determined by a lottery. (Note that if House A robs House B first, and then C robs A, the houses of A and B would each be empty and C would have acquired the resources of both A and B.)

Two questions about this fateful day:

1. What is the probability that a house is not robbed over the course of the day
2. Suppose that every house has the same amount of cash to begin with — say $100. Which position in the lottery has the most expected cash at the end of the day, and what is that amount?

This one revolves around lots of conditional probability which I really enjoy.

Question One:

This one is the basis for main problem. Each player has a chance in robbing an individual. Therefore the chance to NOT get robbed is the probability that each player fails to rob you, or . This is approximately 0.3677.

Question Two:

This is the real problem here. For each player you need to know the chance that they get robbed before its their turn, then the amount that you are likely to steal from others and finally all multiplied by the chances your stuff gets stolen. This looks like this:

Expected Cash Returns in Relation to Position, n with 1000 players and 100$ starting sum

This, plotted, leaves you with the following graph:

To check the answer I plotted this against a simulated version of the game which was run a million times:

import random
import matplotlib.pyplot as plt

#Run a house robbing simulation
def house_rob(players):
  holdings = players*[1]
  for i in range(0,players):
    unselected = True
    to_take_from = 0
    while to_take_from == i or unselected:
      to_take_from = random.randint(0,players-1)
      unselected = False
    holdings[i] += holdings[to_take_from]
    holdings[to_take_from] = 0
    #print(holdings)
    unrobbed_people = holdings.count(1)
    return holdings,unrobbed_people

#Iteratively run trials
players = 1000
avg_holdings = players*[0]
num_trials = 300000
for trials in range(0,num_trials):
  results = house_rob(players)
  for i in range(len(results[0])):
    avg_holdings[i] += results[0][i]

#Make data ready to graph
player_titles = []
real = []
plottable_holdings = []
delta = []
for n in range(1,players):
  real.append((1*(998/999)**n+1+1*((1-(998/999)**n)*(1/999)))*(998/999)**(999-n))
  player_titles.append(n)
plottable_holdings.append(avg_holdings[n]/num_trials)
  delta.append(plottable_holdings[n-1]-real[n-1])

#And make those graphs
plt.scatter(player_titles,plottable_holdings)
plt.scatter(player_titles,real, color = "g")
plt.title("Expected Cash Returns in Relation to Position")
plt.xlabel("Player Position")
plt.ylabel("Expected Cash Returns, in hundreds")
plt.show()