x.py

#!/usr/bin/env python3.1
#
#	x.py ... Classical Parrondo paradox.
#	
#	Copyright (C) 2010  KennyTM~
#	
#	This program is free software: you can redistribute it and/or modify
#	it under the terms of the GNU General Public License as published by
#	the Free Software Foundation, either version 3 of the License, or
#	(at your option) any later version.
#	
#	This program is distributed in the hope that it will be useful,
#	but WITHOUT ANY WARRANTY; without even the implied warranty of
#	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#	GNU General Public License for more details.
#	
#	You should have received a copy of the GNU General Public License
#	along with this program.  If not, see <http://www.gnu.org/licenses/>.
#

from collections import defaultdict, Iterable
from itertools import starmap, chain
from math import fsum
from operator import mul

class Parrondo(Iterable):
	"This class represents the payoff distribution in the classical Parrondo's game."

	def __init__(self, p, p1, p2):
		self.distr = {0: 1}
		self.p = p
		self.p1 = p1
		self.p2 = p2
	
	def __iter__(self):
		return sorted(self.distr.items())
	
	def mean(self):
		"Compute the mean of the current payoff distribution."
		total_prob = fsum(self.distr.values())
		total_payoff = fsum(starmap(mul, self.distr.items()))
		return total_payoff / total_prob
		
	def A(self, q):
		"Flip the 'A' coin, yield the new payoff distribution."
		if q:	
			pq = self.p * q
			pcq = q - pq
			for payoff, prob in self.distr.items():
				yield (payoff+1, prob * pq)
				yield (payoff-1, prob * pcq)
		
	def B(self, q):
		"Flip the 'B' coin."
		if q:
			p1q = self.p1 * q
			p2q = self.p2 * q
			for payoff, prob in self.distr.items():
				pq = p2q if payoff % 3 else p1q
				yield (payoff+1, prob * pq)
				yield (payoff-1, prob * (q - pq))


	def play(self, q):
		"Play A with probability q, play B with probability 1-q."
		res = {}
		for payoff, prob in chain(self.A(q), self.B(1-q)):
			if payoff in res:
				res[payoff] += prob
			else:
				res[payoff] = prob
		self.distr = res

	def chop(self):
		"Remove extremely low probability items."
		self.distr = {payoff: prob for payoff, prob in self.distr.items() if payoff == 0 or payoff*prob > 1e-5}

fibon = [1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987]

x = Parrondo(0.5 - 0.005,  0.1 - 0.005,  0.75 - 0.005)
for i in range(1000):
	print("{}\t{}".format(i, x.mean()))
	x.play( i not in fibon )