######################################################################
##CountChickens.txt: Save this file as  CountChickens.txt            #
## To use it, stay in the                                            #
##same directory, get into Maple (by typing: maple <Enter> )         #
##and then type:  CountChickens.txt<Enter>                           #
##Then follow the instructions given there                           #
##                                                                   #
##Written by Doron Zeilberger, Rutgers University ,                  #
#DoronZeil at gmail dot com                                          #
######################################################################
 
#Created: July 2019
 
print(`Created:  July , 2019`):
print(` This is CountChickens.txt `):
print(`It is the packages that accompanies the article `):
print(`Using Symbolic Computation to analyze Some Childish Board Games`):
print(`by Shalosh B. Ekhad  and Doron Zeilberger`):
print(`published in the Personal J. of Ekhad and Zeilberger and arxiv.org`):
print(``):
print(``):
print(`Please report bugs to DoronZeil at gmail dot com `):
print(``):
 print(`The most current version of this  package and paper`):
 print(` are  available from`):
 print(`http://sites.math.rutgers.edu/~zeilberg/  .`):


print(`---------------------------------------`):
 print(`For a list of the Supporting procedures type ezra1();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the Simulation rocedures`):
 print(` type ezraSi();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the Floating point procedures, that were supposed to make it faster (but they don't)`):
 print(` so don't use them`):
 print(` type ezraFloat();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the STORY procedures type ezraSt();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

print(`---------------------------------------`):
 print(`For a list of the MAIN procedures type ezra();, for help with`):
 print(`a specific procedure, type ezra(procedure_name);   .`):
 print(``):
print(`---------------------------------------`):

with(combinat):
Digits:=20:

ezraSi:=proc()

if args=NULL then
 print(` The Simulation  procedures are:  ManyPlays, OnePlay, OnePlayV `):
 print(``):

else
ezra(args):
fi:

end:

ezraFloat:=proc()

if args=NULL then
 print(` The Floating  procedures are:  ChSerFl `):
 print(``):

else
ezra(args):
fi:

end:

ezraSt:=proc()

if args=NULL then
 print(` The STORY procedures are:  ChickPaper `):
 print(``):

else
ezra(args):
fi:

end:

ezra1:=proc()

if args=NULL then
 print(` The supporting procedures are:  CCb, CCb1, FindNext, KillGadol, KillNeg, NeiState, Pashet, PreUmb, Umb `):
 print(``):

else
ezra(args):
fi:

end:

ezra:=proc()

if args=NULL then
 print(`The main procedures are: Alpha, ChSer, Info , RandGame `):

 print(` `):





elif nops([args])=1 and op(1,[args])=Alpha then
print(`Alpha(f,x,N): Given a probability generating function`):
print(`f(x) (a polynomial, rational function and possibly beyond)`):
print(`returns a list, of length N, whose `):
print(`(i) First entry is the average`):
print(` (ii): Second entry is the variance `):
print(` for i=3...N, the i-th entry is the so-called alpha-coefficients `):
print(` that is the i-th moment about the mean divided by the `):
print(` variance to the power i/2 (For example, i=3 is the skewness `):
print(` and i=4 is the Kurtosis) `):
print(`If f(1) is not 1, than it first normalizes it by dividing`):
print(`by f(1) (if it is not 0) .`):
print(` For example, try: `):
print(` Alpha(((1+x)/2)^100,x,4); `):

elif nops([args])=1 and op(1,[args])=CCb then
print(`CCb(C,D,P,S,T): The "Count Your Chickens board". It is a triple , the labels , followed by the board, followed by the set of blue squares`):
print(`An emptry square is indicated by 0. It is assumed that the last square (nops(gu[2])+1) has all of them. Try:`):
print(`CCb(C,D,P,S,T);`):

elif nops([args])=1 and op(1,[args])=CCb1 then
print(`CCb1(): A compact version of CCb(C,D,P,S,T). Try:`):
print(`CCb1();`):


elif nops([args])=1 and op(1,[args])=ChickPaper then
print(`ChickPaper(CB,K): Inputs a general "Count your Chikens Board", CB=[Spinner,Board,SetOfBlueSquares], and a positive`):
print(`integer K, outputs an article narrating (i) the probability of winning (ii) the expected number of chicks`):
print(`at the end, (iii) the expected number of moves to get to finish the game, higher moments of each `):
print(` and the correlation. All this under the assumption that it ended in <=K moves. It also tells you `):
print(` the probability that it will last longer than K moves. Please take K to be such that this last probability `):
print(` is really small (to make it a good approximation to the real thing (with K=infinity). `):
print(` Try: `):
print(`ChickPaper([[Sheep,Cow],[0,0,Sheep,Cow,0,Cow,0,Sheep], {3,6}] ,50);`):
print(` ChickPaper(CCb(Cow,Dog,Pig,Sheep,Tractor),60); `):

elif nops([args])=1 and op(1,[args])=ChSer then
print(`ChSer(CB,t,X,K): The grand-generating function according to the number of moves (given by X), number of chickens (given by t)`):
print(`in a Count Your Chicken Game with board CB.`):
print(`conditioned on it terminating in <=K moves. `):
print(`The output is a polynomial of degree K where the coeff. of X^i is the`):
print(`prob. generating function of ending after round i, (according to the number of chickens) (conditioned on it terminating  in <=K moves.`):
print(` It also returns`):
print(`the probability, in floating point of NOT terminating in <=K rounds. Try:`):
print(`ChSer(CCb1(),t,X,30);`):
print(`ChSer([[S,C],[0,0,S,C,0,C,0,S],{3,6}] ,t,X,30);`):

elif nops([args])=1 and op(1,[args])=ChSerFl then
print(`ChSerFl(CB,t,X,K): A floating point version of ChSer(CB,t,X,K) (q.v.) that was supposed to be faster, but it is not. Try: `):
print(`ChSerFl(CCb1(),t,X,30);`):
print(`ChSerFl([[S,C],[0,0,S,C,0,C,0,S],{3,6}] ,t,X,30);`):

elif nops([args])=1 and op(1,[args])=FindNext then
print(`FindNext(CB,i,j): inputs a "Count your Chicken board, and a location i, and an outcome of the spinner  from 1 to 6, outputs`):
print(` the next location of the chicken. Try: `):
print(` FindNext(CCb(C,D,P,S,T),4,4); `):


elif nops([args])=1 and op(1,[args])=Info then
print(`Info(CB,K): Inputs a Count Your Chickens! board CB, and  a positive integer K, for the number of allowed rounds,`):
print(`and outputs a list consisting of`):
print(`(i) The probability of winning, `):
print(`(ii) The stat. for the "number of chicks": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]`):
print(`(iii) The stat. for the "number of rounds": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]`):
print(`(vi) stat for the  "number of rounds if winning": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]`):
print(`(v) The correlation between number of chicks and number of rounds`):
print(` (vi) The prob. of taking more than K rounds . Try:`):
print(`Info(CCb1(),50);`):
print(`Info([[S,C],[0,0,S,C,0,C,0,S],{3,6}] ,50);`):

elif nops([args])=1 and op(1,[args])=KillGadol then
print(`KillGadol(f,t,N): replaces each power larger than N by t^N. Try:`):
print(`KillGadol(a*t+b*t^3+d*t^4,t,2);`):

elif nops([args])=1 and op(1,[args])=KillNeg then
print(`KillNeg(f,t): inputs a Laurent polynomial f in t, converts all negative powers to 1. Try`):
print(`KillNeg(3/t+1+5*t,t);`):

elif nops([args])=1 and op(1,[args])=ManyPlays then
print(`ManyPlays(G,K,M):Statistical analysis of M runs of the game G with max. moves K. Try:`):
print(`It returns the`):
print(`ratio of won games, expected number of chicks at the end, and expected number of turns, and and expected number of turns if winning. Try:`):
print(`ManyPlays(CCb1(),40,200);`):

elif nops([args])=1 and op(1,[args])=NeiState then
print(`NeiState(CB,st): Given a Count Your Chicken Board, CB, and a state a, where a is the location 1<=a<=nops(CB[2])+1, and b is`):
print(` outputs the six equally states that come from it together with the gain of chicks.`):
print(`The absorbing Losing state is indicated by 0.`):
print(`Try:`):
print(`NeiState(CCb1(),1);`):

elif nops([args])=1 and op(1,[args])=OnePlay then
print(`OnePlay(G,K): Simulating one play of the Count Your Chicken Style game G with a cap of K moves`):
print(`it returns the number of moves needed (K+1 if it does not end by K moves) followed by the`):
print(`number of chicks at the end (it is a win if there nops(G[2]) chicks. Try`):
print(`OnePlay([[S,C],[0,0,S,C,0,C,0,S],{3,6}], 30);`):
print(`OnePlay(CCb1(),40);`):

elif nops([args])=1 and op(1,[args])=OnePlayV then
print(`OnePlayV(G,K): verbose vesion of OnePlay(G,K) (q.v.) . Try: `):
print(`OnePlayV([[S,C],[0,0,S,C,0,C,0,S],{3,6}], 30);`):
print(`OnePlayV(CCb1(),40);`):
print(`OnePlayV(CCb(Cow,Dog,Pig,Sheep,Tractor),30); `):

elif nops([args])=1 and op(1,[args])=Pashet then
print(`Pashet(f,x,t,eps): inputs a polynomial f in the variables x and t and a cut-op eps, discards`):
print(`terms whose absolute value is less than eps. Try:`):
print(`Pashet(x*t+1/10^11*x-1/10^11,x,t,1/10^10);`):


elif nops([args])=1 and op(1,[args])=PreUmb then
print(`PreUmb(f,s,t,CB): Given a monomial f in s and t applies the pre-umbra operation implied by the Count Your Chickens Board CB.`):
print(` Try: `):
print(`PreUmb(s^5*t^11,s,t,CCb1());`):


elif nops([args])=1 and op(1,[args])=RandGame then
print(`RandGame(k,n,b): inputs positive integers k,n,b and outputs a Count Your Chicken Style game with `):
print(`a dial of size k, a board with n squares and a set of blue squares of size b, Try:`):
print(`RandGame(6,40,3);`):

elif nops([args])=1 and op(1,[args])=Umb then
print(`Umb(f,s,t,CB): Given a polynomial f in s and t applies the umbral operation implied by the Count Your Chickens Board CB.`):
print(` Try: `):
print(`Umb(s^5*t^11/2+s^8*t^16/2,s,t,CCb1());`):


else
print(`There is no ezra for`,args):
fi:
 
end:



##########From HISTABRUT

#AveAndMoms(f,x,N): Given a probability generating function
#f(x) (a polynomial, rational function and possibly beyond)
#returns a list whose first entry is the average 
#(alias expectation, alias mean)
#followed by the variance, the third moment (about the mean) and
#so on, until the N-th moment (about the mean).
#If f(1) is not 1, than it first normalizes it by dividing
#by f(1) (if it is not 0) .
#For example, try:
#AveAndMoms(((1+x)/2)^100,x,4);

AveAndMoms:=proc(f,x,N) local mu,F,memu1,gu,i:
mu:=simplify(subs(x=1,f)):

if mu=0 then
print(f, `is neither a prob. generating function nor can it be made so`):
RETURN(FAIL):
fi:

F:=f/mu:


memu1:=simplify(subs(x=1,diff(F,x))):

gu:=[memu1]:

F:=F/x^memu1:

F:=x*diff(F,x):
for i from 2 to N do
 F:=x*diff(F,x):
 gu:=[op(gu),simplify(subs(x=1,F))]:
od:

gu:

end:

#Alpha(f,x,N): Given a probability generating function
#f(x) (a polynomial, rational function and possibly beyond)
#returns a list, of length N, whose 
#(i) First entry is the average
#(ii): Second entry is the variance
#for i=3...N, the i-th entry is the so-called alpha-coefficients
#that is the i-th moment about the mean divided by the
#variance to the power i/2 (For example, i=3 is the skewness
#and i=4 is the Kurtosis)
#If f(1) is not 1, than it first normalizes it by dividing
#by f(1) (if it is not 0) .
#For example, try:
#Alpha(((1+x)/2)^100,x,4);
Alpha:=proc(f,x,N) local gu,i:
 gu:=AveAndMoms(f,x,N):
 if gu=FAIL then
  RETURN(gu):
 fi:

if gu[2]=0 then
print(`The variance is 0`):
  RETURN(FAIL):
fi:

[gu[1],gu[2],seq(gu[i]/gu[2]^(i/2),i=3..N)]:

end:

##########End From HISTABRUT




#Tsamek1(L): Given a  of list of destinations with their probability [j,Pr[i->j]]
#ony mentions those destination with non-zero prob. Try:
#Tsamek1([[2,1/2],[3,1/2],[4,0],[5,0]]));
Tsamek1:=proc(L) local gu, i:
gu:=[]:

for i from 1 to nops(L) do
 if L[i][2]<>0 then
  gu:=[op(gu), L[i]]:
 fi:
od:
gu:
end:



#Tsamek(L): Given a list of lists L where L[i] is the list of destinations with their probability [j,Pr[i->j]]
#shrinks it as to only mention  destination with non-zero probabilities. Try:
#Tsamek([[[2,1/3],[3,2/3],[4,0]]]);
Tsamek:=proc(L) local i:
[seq(Tsamek1(L[i]),i=1..nops(L))]:
end:


CCb1:=proc():
[[1,2,3,4,5],
[0, 0, 4, 3, 5, 1, 2, 3, 1, 2, 4, 5, 0, 1, 3, 0, 0, 0, 5, 0, 5, 2, 4, 1, 2, 3, 5, 0, 4, 1, 0, 0, 5, 3, 4, 2, 0, 4, 1, 3], 
{5, 9, 23, 36, 40}]
end:


#CCb(C,D,P,S,T): The "Count Your Chickens board". It is a triple , the labels , followed by the board, followed by the set of blue squares
#An emptry square is indicated by 0. It is assumed that the last square (nops(gu[2])+1) has all of them. Try:
#CCb();
CCb:=proc(C,D,P,S,T) local gu1,gu2,gu3:
gu1:=[C,D,P,S,T]:
gu2:=[0,0,S,P,T,C,D,P,C,D,
      S,T,0,C,P,0,0,0,T,0,
      T,D,S,C,D,P,T,0,S,C,
      0,0,T,P,S,D,0,S,C,P]:
gu3:={5,9,23,36,40}:
[gu1,gu2,gu3]:
end:

#FindNext(CB,i,j): inputs a "Count your Chicken board, and a location i, and an outcome of the spinner  from 1 to 6, outputs
#the next location of the chicken. Try:
#FindNext(CCb(C,D,P,S,T),4,4);

FindNext:=proc(CB,i,j) local mu,i1:


if j=6 then
 RETURN(i):
fi:

mu:=CB[1][j]:

for i1 from i+1 to nops(CB[2]) while CB[2][i1]<>mu do od:
i1:
end:


#NeiState(CB,a): Given a Count Your Chicken Board, CB, and a state a the location 1<=a<=nops(CB[2])+1,
#outputs the six equally likely states that come from it, together with the gain of chicks. 
#The absorbing Losing state is indicated by 0 and winning by 1
#Try:
#NeiState(CCb1(),1);
NeiState:=proc(CB,a) local j,lu,a1:

if a>=nops(CB[2])+1 then
 RETURN([a]):
fi:

lu:=[]:

for j from 1 to nops(CB[1]) do
 a1:=FindNext(CB,a,j):
 
 if member(a1,CB[3]) then
  lu:=[op(lu),[a1,a1-a+1] ]:
else
  lu:=[op(lu),[a1,a1-a]]:
 fi:

od:

 lu:=[op(lu),[a,-1]]:

lu:

end:


#PreUmb(f,s,t,CB): Given a monomial f in s and t applies the pre-umbra operation implied by the Count Your Chickens Board CB.
#Try:
#PreUmb(s^5*t^11,s,t,CCb1());
PreUmb:=proc(f,s,t,CB) local d1,gu,c,i,lu:

d1:=degree(f,s):
 c:=normal(f/(s^d1)):
if (normal(diff(c,s))<>0) then
 RETURN(FAIL):
fi:

gu:=NeiState(CB,d1):

lu:=c/nops(gu)*add(s^gu[i][1]*t^gu[i][2],i=1..nops(gu)):

lu:=expand(KillGadol(KillNeg(lu,t),t,nops(CB[2]))):

end:



#Umb(f,s,t,CB): Given a polynomial f in s and t applies the umbral operation implied by the Count Your Chickens Board CB.
#Try:
#Umb(s^5*t^11/2+ s^7*t^12/2,s,t,CCb1());
Umb:=proc(f,s,t,CB) local i:

if not type(f,`+`) then
 RETURN(PreUmb(f,s,t,CB)):
fi:

add(PreUmb(op(i,f),s,t,CB),i=1..nops(f)):
end:

#ChSer(CB,t,X,K): The grand-generating function according to the number of moves (given by X), number of chickens (given by t)
#in a Count Your Chicken Game with board CB.
#The output is a polynomial of degree K where the coeff. of X^i is the
#prob. generating function of ending after round i, (according to the number of chickens). It also returns
#the probability, in floating point of NOT terminating in <=K rounds. Try:
#ChSer(CCb1(),t,X,30);
ChSer:=proc(CB,t,X,K) local N,gu,mu,i,s:

N:=nops(CB[2])+1:


gu:=0:

mu:=s:


for i from 1 to K do


 mu:=Umb(mu,s,t,CB):
if mu=FAIL then
  RETURN(FAIL):
fi:

 gu:=expand(gu+coeff(mu,s,N)*X^i):
 mu:=expand(mu-coeff(mu,s,N)*s^N):
od:

if gu<>0 then
[normal(gu/subs({X=1,t=1},gu)),evalf(subs({s=1,t=1},mu))]:
else
[gu,evalf(subs({s=1,t=1},mu))]:
fi:

end:



#Pashet(f,x,t,eps): inputs a polynomial f in the variables x and t and a cut-op eps, discards
#terms whose absolute value is less than eps. Try:
#Pashet(x*t+1/10^11*x-1/10^11,x,t,1/10^10);
Pashet:=proc(f,x,t,eps) local i,j,gu,f1:
gu:=0:
for i from ldegree(f,x) to degree(f,x) do
 f1:=coeff(f,x,i):
 for j from ldegree(f1,t) to degree(f1,t) do
  if abs(coeff(f1,t,j))>eps then
    gu:=gu+coeff(f1,t,j)*x^i*t^j:
  fi:
 od:
od:
gu:
end:





#ChickPaper(CB,K): Inputs a general "Count your Chikens Board", CB=[Spinner,Board,SetOfBlueSquares], and a positive
#integer K, outputs an article narrating (i) the probability of winning (ii) the expected number of chicks
#at the end, (iii) the expected number of moves to get to finish the game, higher moments of each
#and the correlation. All this under the assumption that it ended in <=K moves. It also tells you
#the probability that it will last longer than K moves. Please take K to be such that this last probability
#is really small (to make it a good approximation to the real thing (with K=infinity).
#Try:
#ChickPaper(CCb1(),60);
ChickPaper:=proc(CB,K) local lu, gu,t,X,Fox,i,gu1,gu2,mu1,mu2,t0,m11,c,d:

t0:=time():


if gu=FAIL then
 RETURN(FAIL):
fi:

print(`Statistical Analysis of a Count Your Chickens! Board Game`):
print(``):
print(`By Shalosh B. Ekhad `):
print(``):

print(`Consider the "Cooperative" ( i.e. Solitaire) game where one spins a spinner with`, nops(CB[1])+1 , `possible outcomes labelled `):
 print([op(CB[1]),Fox]):
print(``):
print(` and there is board with `, nops(CB[2])+1, `squares that looks as follows` ):
print(``):
lu:=CB[2]:
#lu:=subs(0=`_`,lu):
lu:=subs(0=` `,lu):
lu:=[op(lu),{op(CB[1])}]:
print(lu):
print(``):
print(`where the last square contains all the animals `):
print(`There is also a set of locations marked blue. These are in the following locations`):
print(``):
print(CB[3]):

print(``):
print(`You start at the first square, and advance along the board, at the same time  gathering chicks, according to the rules. `):
print(`Your goal is to arrive at the coop located at the last square with a total of`, nops(CB[2]), `chicks. `):
print(``):
print(`The rules are as follows. If you get a Fox, then you remove one chick from the coop (unless it is emtpy, then you do nothing)`):
print(`and you stay at the same place. Otherwise you locate the next occurrence, on the board of the animal that you spun, and move`):
print(`the pawn (chicken) to that location. You add the number of squares thus advanced to the coop.`):
print(`If the new location is a blue square, then you add an additional chick to the coop. `):
print(`You keep going until you reach the last square (where the coop is)`):
print(`If the number of chicks in the coop is`, nops(CB[2]), `then you won, otherwise you lost. `):
print(`Here we will analyze this game limiting the game to <=`, K, `rounds. `):
print(``):
gu:=ChSer(CB,t,X,K):
print(`Note that the probability that the game will last more than`, K, `rounds is `, gu[2]):

if gu[2]<10^(-10) then
print(`Hence this is a very good approximation to the ideal case when there is no limit to the number of rounds.`):
fi:

gu:=gu[1]:


gu1:=subs(X=1,gu):
gu2:=subs(t=1,gu):

print(`The probabiliy of winning (under the above assumption) is`):
print(``):
print(evalf(coeff(gu1,t,nops(CB[2])))):
print(``):



mu1:=evalf(Alpha(gu1,t,6)):

mu2:=evalf(Alpha(gu2,X,6)):

print(`The expected  number of chicks in the coop at the end is`):
print(``):
print(mu1[1]):
print(``):

print(`The variance of the random variable "number of chicks at the end" is `):
print(``):
print(mu1[2]):
print(``):


print(`The skewness (aka scaled 3rd moment) of the random variable "number of chicks at the end" is `):
print(``):
print(mu1[3]):
print(``):


print(`The kurtosis  (aka scaled 4th moment) of the random variable "number of chicks at the end" is `):
print(``):
print(mu1[4]):
print(``):

for i from 5 to 6 do
print(``):
print(`The scaled `, i, `-th moment about the mean is this random variable is`):
print(``):
print(mu1[i]):
print(``):
od:
print(``):

print(`The expected  number of rounds until the end is`):
print(``):
print(mu2[1]):
print(``):
print(`The variance of the random variable "number of rounds" is ` ):
print(``):
print(mu2[2]):
print(``):

print(`The skewness (aka scaled 3rd moment) of the random variable "number of rounds " is `):
print(``):
print(mu2[3]):
print(``):


print(`The kurtosis  (aka scaled 4th moment) of the random variable "number of rounds" is `):
print(``):
print(mu2[4]):
print(``):


for i from 5 to 6 do
print(``):
print(`The scaled `, i, ` -th moment about the mean is this random variable is`):
print(``):
print(mu2[i]):
print(``):
od:

m11:=(evalf(subs({t=1,X=1},diff(diff(gu,X),t)))-mu1[1]*mu2[1])/sqrt(mu1[2]*mu2[2]):

print(`The correlation between the number of chicks and the number of rounds is`):
print(``):
print(m11):
print(``):

gu:=evalf(Pashet(gu,X,t,0.5/10^10),10):

print(`Finally, the bi-variate  probability generating function in the variables `, X, t, `whose coeficient of  `, X^c*t^d , `is the probability`):
print(`that the game ends after`, c, ` rounds and with `, d, `chicks , is `):
print(`We ignore terms smaller than `,  0.5/10^10 ):
gu:=add(sort(coeff(gu,t,degree(gu,t)-i),X)*t^(degree(gu,t)-i),i=0..degree(gu,t)):
print(``):
print(gu):
print(``):
print(`and in Maple notation `):
print(``):
lprint(gu):
print(``):

print(`------------------------------------`):
print(``):
print(`This ends this article that took`, time()-t0, `seconds to generate. `):
print(``):

end:

#Info(CB,K): Inputs a Count Your Chickens! board CB, and  a positive integer K, for the number of allowed rounds,
#and outputs a list consisting of
#(i) The probability of winning, 
#(ii) The stat. for the "number of chicks": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]
#(iii) The stat. for the "number of rounds": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]
#(vi) stat for the  "number of rounds if winning": [av,variance,skewness,kurtosis, scaled 5th moment about the mean, scaled 6th moment about the mean]
#(v) The correlation between number of chicks and number of rounds
#(vi) The prob. of taking more than K rounds
Info:=proc(CB,K) local t,X,mu,mu1,ProbW,A1,A2,A3,eps,m11:
mu:=ChSer(CB,t,X,K):
eps:=mu[2]:
mu:=mu[1]:

mu1:=coeff(mu,t,nops(CB[2])):

ProbW:=evalf(subs(X=1,mu1)):

A1:=evalf(Alpha(subs(X=1,mu),t,6)):

A2:=evalf(Alpha(subs(t=1,mu),X,6)):

A3:=evalf(Alpha(mu1/subs(X=1,mu1),X,6)):


m11:=(subs({t=1,X=1},diff(diff(mu,X),t))-A1[1]*A2[1])/sqrt(A1[2]*A2[2]):
[ProbW,A1,A2,A3,m11,eps]:
end:


#KillNeg(f,t): inputs a Laurent polynomial f in t, converts all negative powers to 1. Try
#KillNeg(3/t+1+5*t,t);
KillNeg:=proc(f,t) local i:

if ldegree(f,t)>=0 then
 RETURN(f):
else
add(coeff(f,t,i),i=ldegree(f,t)..-1)+add(coeff(f,t,i)*t^i,i=0..degree(f,t)):
fi:
end:

#KillGadol(f,t,N): replaces each power larger than N by t^N. Try:
#KillGadol(a*t+b*t^3+d*t^4,t,2);
KillGadol:=proc(f,t,N) local i:
if degree(f,t)<=N then
  RETURN(f):
else
add(coeff(f,t,i)*t^i,i=ldegree(f,t)..N)+add(coeff(f,t,i),i=N+1..degree(f,t))*t^N:
fi:
end:


#RandSet1(S,k): a random subset with k elements of the set S, try:
#RandSet1({1,2,3,4},2);
RandSet1:=proc(S,k) local gu,S1,i,a:
if k>nops(S) then
 RETURN(FAIL):
fi:

if k=0 then
 RETURN({}):
fi:

gu:={}:
S1:=S:

for i from 1 to k do
 a:=S1[rand(1..nops(S1))()]:
 gu:=gu union {a}:
 S1:=S1 minus {a}:
od:
gu:
end:


#RandGame(k,n,b): inputs positive integers k,n,b and outputs a Count Your Chicken Style game with 
#a dial of size k, a board with n squares and a set of blue squares of size b, Try:
#RandGame(6,40,3);
RandGame:=proc(k,n,b) local S,i,Bo,Bl,ra:
S:=[seq(i,i=1..k-1)]:
ra:=rand(1..k):
Bo:=[0,seq(ra(),i=1..n-1)]:
Bl:=RandSet1({seq(i,i=2..n-1)},b):
[S,Bo,Bl]:
end:



###Start Floating point version

#PreUmbFl(f,s,t,CB): Floating point version of PreUmb(f,s,t,CB) (q.v.). Try:
#PreUmbFl(s^5*t^11,s,t,CCb1());
PreUmbFl:=proc(f,s,t,CB) local d1,gu,c,i,lu:

d1:=degree(f,s):
 c:=evalf(normal(f/(s^d1))):
if (normal(diff(c,s))<>0) then
 RETURN(FAIL):
fi:

gu:=NeiState(CB,d1):

lu:=evalf(c/nops(gu)*add(s^gu[i][1]*t^gu[i][2],i=1..nops(gu))):

lu:=expand(KillGadol(KillNeg(lu,t),t,nops(CB[2]))):

end:



#UmbFl(f,s,t,CB): Floating point  version of Umb(f,s,t,CB) (q.v.). Try:
#UmbFl(s^5*t^11/2+ s^7*t^12/2,s,t,CCb1());
UmbFl:=proc(f,s,t,CB) local i:

if not type(f,`+`) then
 RETURN(PreUmbFl(f,s,t,CB)):
fi:

add(PreUmbFl(op(i,f),s,t,CB),i=1..nops(f)):
end:

#ChSerFl(CB,t,X,K):  Floating point version of ChSer(CB,t,X,K):  (q.v.) Try:
#ChSerFl(CCb1(),t,X,30);
ChSerFl:=proc(CB,t,X,K) local N,gu,mu,i,s:

N:=nops(CB[2])+1:


gu:=0:

mu:=s:


for i from 1 to K do


 mu:=UmbFl(mu,s,t,CB):
if mu=FAIL then
  RETURN(FAIL):
fi:

 gu:=expand(gu+coeff(mu,s,N)*X^i):
 mu:=expand(mu-coeff(mu,s,N)*s^N):
od:

if gu<>0 then
[normal(gu/subs({X=1,t=1},gu)),evalf(subs({s=1,t=1},mu))]:
else
[gu,evalf(subs({s=1,t=1},mu))]:
fi:

end:


####end floating point version

#OnePlayV(G,K): Simulating one play of the Count Your Chicken Style game G with a cap of K moves
#it returns the number of moves needed (K+1 if it does not end by K moves) followed by the
#number of chicks at the end (it is a win if there nops(G[2]) chicks. Try:
#OnePlayV(CCb1(),40);
OnePlayV:=proc(G,K) local st,ra,i,an,ch,st1:
print(`This is a simulation of one play of a game where there is a spinner with`, nops(G[1])+1, `equally likely animals`):
print(op(G[1])):
print(`together with a Fox`):
print(`The board consists of`, nops(G[2]), `squares, some of them are empty (indicated by 0), and other ones labelled by one of the animals`):
print(`The boad is`):
print(G[2]):
print(``):
print(`In addition the following squares`, op(G[3]), `are blue where you get a bonus if Mama Chicks lends on it`):
print(``):
print(`Mama Chicks proceeds along the board always going to the next square labelled by the animal that was spun`):
print(`and collects more chicks, according to the number of squares that it walked, and gets an extra chick if it lended on a blue square`):
print(`if the spinner spun a Fox, then it loses one chick (unless it has no chicks), and stays where it is`):
print(`You win if you reached the end location` , nops(G[2])+1, `with `, nops(G[2]) , `chicks, otherwise you lose. `):

st:=1:
ch:=0:
ra:=rand(1..nops(G[1])+1):

for i from 1 to K do
print(`This is move number `, i):
an:=ra():



if an=nops(G[1])+1 then

print(` it is a Fox`):
 ch:=max(ch-1,0):
print(`Now there are`, ch, `chicks `):
else
print(` the spinner gave`, G[1][an]):
 for st1 from st+1 to nops(G[2]) while G[2][st1]<>G[1][an] do od:
  print(`Now the location is`, st1):
   if st1=nops(G[2])+1 then
       ch:=ch+st1-st:
     ch:=min(ch,nops(G[2])):
      if ch<nops(G[2]) then
        print(`We finished with `, ch, `chicks that is less than`, nops(G[2]), ` so we lost. `):
        RETURN([i,ch]):
      else
      print(`We finished with `, ch, `chicks that is what we want!`, nops(G[2]), ` so we won. `):
          RETURN([i,ch]):
      fi:       
   else
    ch:=ch+st1-st:
    if member(st1,G[3]) then
    print(`Since location `, st1, `is blue, we get an extra chick`):
      ch:=ch+1:
      ch:=min(ch,nops(G[2])):
     print(`Now there are`, ch, `chicks `):
    fi:
     st:=st1:
  fi:
fi:
od:
print(`We did not finish after`, K, `moves, we are still in state`, st, `but with only`, ch, `chicks . `):
[K+1,ch]:
end:
 
#OnePlay(G,K): Simulating one play of the Count Your Chicken Style game G with a cap of K moves
#it returns the number of moves needed (K+1 if it does not end by K moves) followed by the
#number of chicks at the end (it is a win if there nops(G[2]) chicks. Try:
#OnePlay(CCb1(),40);
OnePlay:=proc(G,K) local st,ra,i,an,ch,st1:

st:=1:
ch:=0:
ra:=rand(1..nops(G[1])+1):

for i from 1 to K do
an:=ra():



if an=nops(G[1])+1 then

 ch:=max(ch-1,0):
else
 for st1 from st+1 to nops(G[2]) while G[2][st1]<>G[1][an] do od:
   if st1=nops(G[2])+1 then
       ch:=ch+st1-st:
     ch:=min(ch,nops(G[2])):
      if ch<nops(G[2]) then
             RETURN([i,ch]):
      else
               RETURN([i,ch]):
      fi:       
   else
    ch:=ch+st1-st:
    if member(st1,G[3]) then
      ch:=ch+1:
      ch:=min(ch,nops(G[2])):
    fi:
     st:=st1:
  fi:
fi:
od:
[K+1,ch]:
end:
  



#ManyPlays(G,K,M):Statistical analysis of M runs of the game G with max. moves K. Try:
#It returns the
#ratio of won games, expected number of chicks at the end, and expected number of turns, and expected number of turns if winning
#ManyPlay(CCb1(),40,200);

ManyPlays:=proc(G,K,M) local guL,guW,guN,X,t,i,P,gu:

guL:=0:
guW:=0:
guN:=0:

for i from 1 to M do
P:=OnePlay(G,K):
if P[1]=K+1 then
 guN:=guN+X^(K+1)*t^P[2]:
elif P[2]=nops(G[2]) then
   guW:=guW+X^P[1]*t^P[2]:
else
 guL:=guL+X^P[1]*t^P[2]:
fi:

gu:=guL+guW:
od:

if gu=0 then
  RETURN(FAIL):
fi:

gu:=gu/subs({X=1,t=1},gu):


[evalf(subs({X=1,t=1},guW)/M), evalf(subs({X=1,t=1},diff(gu,t))), 
evalf(subs({X=1,t=1},diff(gu,X))/subs({X=1,t=1},gu)),
evalf(subs({X=1,t=1},diff(guW,X))/subs({X=1,t=1},guW))]:

end: