• API
• FAQ
• Tools
• Archive
SHARE
TWEET # Untitled dechicom  Sep 6th, 2018 70 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. #! /usr/bin/python
2. # -*- coding: utf-8 -*-
3.
4. import copy as CO
5. import bisect as BI
6. import math as MA
7. import matplotlib.pyplot as PLT
8. '''
9. dice game
11. this one add matplotlib for visualizing the distri
12. '''
13.
14. class CountingValueError(ValueError):
15.     pass
16.
17. def arangements_ie( n, l ):
18.     """ arrangements_ie ( n,( k1, k2, k3,...) )
19.    Return the number of possible arrangements of n elements
20.    with k1, k2,.... elements being indistinguishable.
21.    """
22.     a= MA.factorial( n )
23.     K=0
24.     for k in l:
25.         K += k
26.         if K > n:
27.             raise CountingValueError ("Can't have k indistinguishable objects taken from n objects if k > n")
28.         a= a / MA.factorial( k )
29.     return a
30. def result(result):
31.     print "gain\t\tresult"
32.     for r in result:
33.         pass
34.
35.
36. def storegain(trialgain,trialproba):
37.     #print "### trial gain:",trialgain,"|proba:",trialproba*100,"%"
38.     global tp
39.     global gains
40.     gains.append([trialgain,trialproba])
41.     tp+=trialproba
42.     i = BI.bisect_right(totwin,trialgain) -1
43.     totproba [i] += trialproba
44.
45. def trial(t):
46.     global g_no
47.     global g_N
48.     trialproba = 0.0
49.     trialgain = 0.0
50.     for i  in range (g_no):
51.         k=t[i]
52.         if k!=0:
53.             p = gamewin[i]
54.             gain = gamewin[i]
55.             trialgain += gain * k
56.             proba=MA.pow(p,k)
57.             if trialproba == 0:
58.                 trialproba= proba
59.             else:
60.                trialproba*= proba
61.             #print "-- outcome:",i,"|k:",k,"|p:",p,"|totgain:", gain * k
62.     trialproba*=arangements_ie(g_N,t)
63.
64.     storegain(trialgain, trialproba)
65.
66.
67. def callback ( t):
68.
69.     #print t
70.     #print "."
71.     global cpt
72.     cpt+=1
73.     trial(t)
74.
75. def rec_gen ( co , t ,s ):
76.     '''
77.     N: number of experiment still  available
78.     co: current outcome (position in t )
79.     t : list of occurence number for each outcome
80.    '''
81.
82.     for i in xrange ( 0,g_N-s+1 ):
83.         t[ co ] = i
84.         ss=s+i
85.         if ss == g_N :
86.              callback ( t )
87.              continue
88.
89.         if  co > 0 :
90.             t2 = CO.copy ( t )
91.             rec_gen ( co-1 , t2,ss )
92.
93.     del(t)
94.     return
95.
96. def gen ( no , N ):
97.     '''
98.    no: number outcome
99.    N : Number of subsequent experiment (trial of N experiment)
100.    t : contain the number of occurence for each outcomes
101.    '''
102.     global g_no  # number of outcome
103.     g_no = no
104.     global g_N  # number of experiment
105.     g_N = N
106.     print "===", g_no,g_N
107.     t = [ 0 ] * no
108.     rec_gen ( no-1 , t,0 )
109.
110. def plot1(totwin,totproba):
111.
112.     PLT.title("distri")
113.     PLT.xlabel("win")
114.     PLT.ylabel("proba")
115.
116.     PLT.plot(totwin,totproba)
117.
118.     PLT.grid(True)
119.     PLT.show()
120.
121.
122.
123.
124.
125.
126. if __name__ == "__main__":
127.
128.     # number of expiriment (how many time I throw the dices)
129.     N = 100
130.
131.     #list of outcomes  (proba, gain)
132.     gamewin = [ [1.0/6, 200], [2.0/6, 10], [3.0/6, 0] ]
133.
134.     # game expectation E(X)
135.     gameexpect= 0
136.     for g in gamewin:
137.         gameexpect += g[ 0 ] * g[ 1 ]
138.     print "Game expectation = ", gameexpect
139.     # game variance var(X)=E(X²)-E²(X)
140.     gamevar = 0
141.     for g in gamewin:
142.         gamevar +=  MA.pow( g[ 1 ], 2 ) * g[ 0 ]
143.     gamevar -= MA.pow( gameexpect , 2 )
144.     print "Game variance = ", gamevar
145.     print "Game Standard Error =", MA.sqrt(gamevar)
146.
147.     # maxwin after N experiments
148.     maxwin = 200 * N
149.
150.     print "Maxwin after N throws:", maxwin
151.     print "Expectation after  N throws:", gameexpect * N
152.
153.     # to have the total gains  we get after N experiments  well distributed in
154.     #  the result  list [(gains,proba), .... ]
155.     totwin = range( 0, int(maxwin), int(maxwin/500))
156.     totproba = [0.0] * len( totwin )
157.
158.     cpt = 0
159.     tp=0
160.     gains=[]
161.     gen( len(gamewin),N )
162.     print "cpt",cpt,"tp",tp
163.     result = zip( totproba, totwin )
164.
165.     # verifying that the calculus are corrects:
166.     # total probability for all final gains = 1
167.     print "verify the sum of all probabilty = 1:",sum(totproba)
168.
169.     # game expectation E(Xn) after N throws
170.     totgameexpect= 0
171.     for g in result:
172.         totgameexpect += g[ 0 ] * g[ 1 ]
173.     print "Game expectation after N throws= ", totgameexpect
174.     # game variance after N throws var(Xn)=E(Xn²)-E²(Xn)
175.     totgamevar = 0
176.     for g in result:
177.         totgamevar +=  MA.pow( g[ 1 ], 2 ) * g[ 0 ]
178.     totgamevar -= MA.pow( totgameexpect , 2 )
179.     print "Game variance N throws = ", totgamevar
180.     print "Game Standard Error N throws =", MA.sqrt(totgamevar)
181.
182.     # have
183.     plot1(totwin,totproba)
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.

Top