##############################################################################################
##This is ComboProject8.txt, a Maple package to generate and investigate  the US Electoral   #
#College   system                                                                            #
####It is the Maple package created by Team 8 in Dr. Z.'s Combinatorics Class at             #
#Rutgers University, Fall 2020.                                                              #
# Save this file as `ComboProject8.txt`, to use it                                           #
#Type, in a Maple session                                                                    #
#read `ComboProject8.txt`):                                                                  #
#and then to get a list of the functions type                                                #
#Help():                                                                                     #
#For Help with any of the functions, type                                                    # 
#Help(FunctionName)                                                                          #
#Team Leader:                                                                                #
#Team members:                                                                               #
##############################################################################################


with(combinat):

print(`This is ComboProject8.txt, a Maple package that is part of Project 8 in Dr. Z.'s Combinatorics Class at Rutgers University`):
print(`To study and simulate vote counting`):
print(``):
print(`Team Leader: tbd `):
print(``):
print(`Other Team members: tbd `):
print(``):
print(`For a list of all the functions type: Help(); `):
print(`For Help with any of the functions, type  Help(FunctionName):`):                                                                        
with(combinat):
Help:=proc()
if nargs=0 then
print(`The available procedures are: GFv, GFvp, IsConsistent, LC, SimuCount, SimuCount1, StatAnal, `):
print(`StatAnalS, USEC`):
print(` `): 


elif nargs=1 and args[1]=GFv then
print(`GFv(L,x): Given a list L of the numbers of electors in a country with nops(L) units (in the case of the US`):
print(` 50 states + DC) outputs the polynomial whose coefficient of x^i is the number of ways that `):
print(` one of the two candidates gotten i votes.  `):
print(`Try: `):
print(`GFv(USEC(),x);`):

elif nargs=1 and args[1]=GFvp then
print(`GFvp(L,p,x): Given a list L of the numbers of electors in a country with nops(L) units (in the case of the US`):
print(`50 states + DC) outputs the polynomial whose coefficient of x^i is the number of ways that`):
print(`one of the two candidates gotten i votes. We assume that the probability of the first candidate is p for each state`):
print(`indpendently. Try:  `):
print(`GFvp(USEC(),1/2,x);`):

elif nargs=1 and args[1]=LC then
print(`LC(p): Outputs 1 with probability p and 0 with prob. 1-p. Try: `):
print(`LC(2/3);`):

elif nargs=1 and args[1]=IsConsistent then
print(`IsConsistent(H): Given a list of partial scores, H, outputs true if and only if the ultimate winner was always ahead if it`):
print(`ended with a tie we declare the first candiate to be the winner. Try: `):
print(`IsConsistent([[0,0],[1,0],[2,0],[2,1]]); `);



elif nargs=1 and args[1]=SimuCount then
print(`SimuCount(L,p,N,K): Given a list L representing the  number of electors in nops(L) states, assuming that`):
print(`the probability of Trump winning is p independently. Simulates N random counts`):
print(`where every unpicked state is equally likely to be counted next. The output is the`):
print(`Statistical analysis for these N random runs up to K moments, followed by the ratio of conisistent histories.Try:`):
print(`SimuCount(USEC(),1/2,1000,4);`):

elif nargs=1 and args[1]=SimuCount1 then
print(`SimuCount1(L,p): Given a list L representing the  number of electors in nops(L) states, assuming that`):
print(`the probability of Trump winning is p independently. Simulates ONE random count`):
print(`where every unpicked state is equally likely to be counted next. The output is the`):
print(`list of partial scores Try:`):
print(`SimuCount1(USEC(),1/2);`):

elif nargs=1 and args[1]=StatAnal then
print(`StatAnal(f,x): Given a generating function f in x`):
print(`finds (i) the average, (ii) the  standard-deviation (iii) the third-throgh-K SCALED moments about the mean. `):
print(` (The r-th scaled moment about the mean of r.v. is the r-th moment divided by the r-th power of the standard-deviation) `):
print(` Try: `):
print(` StatAnal(GFv(USEC(),x), x,10); `):
print(` StatAnal(GFvp(USEC(),11/20,x), x,10); `):

elif nargs=1 and args[1]=StatAnalS then
print(`StatAnalS(f,x): Given a generating function f in x with SYMBOLIC coefficients`):
print(`finds (i) the average, (ii) the  standard-deviation (iii) the third-throgh-K SCALED moments about the mean. `):
print(` (The r-th scaled moment about the mean of r.v. is the r-th moment divided by the r-th power of the standard-deviation) `):
print(` Try: `):
print(` StatAnalS(GFvp(USEC(),p,x), x,2); `):

elif nargs=1 and args[1]=USEC then
print(`USEC(): The increasing list of the number of electors in the electoral state for each the 50 states (+DC)`):
print(`Try: `):
print(` USEC(); `):
print(``):

else
 print(`There is no Help for`):

print(args):
fi:

end:




#USEC(): The increasing list of the number of electors in the electoral state for each the 50 states
USEC:=proc():
[3$8,4$5,5$3,6$6,7$3,8$2,9$3,10$4,11$4,12,13,14,15,16$2,18,20$2,29$2,38,55];
end:




#GFv(L,p,x): Given a list L of the numbers of electors in a country with nops(L) units (in the case of the US
#50 states + DC) outputs the polynomial whose coefficient of x^i is the number of ways that
#one of the two candidates gotten i votes. 
#Try:
#GFv(USEC(),1/2,x);
GFv:=proc(L,x) local i:
sort(expand(mul(1+x^L[i],i=1..nops(L)))):
end:

#GFvp(L,p,x): Given a list L of the numbers of electors in a country with nops(L) units (in the case of the US
#50 states + DC) outputs the polynomial whose coefficient of x^i is the number of ways that
#one of the two candidates gotten i votes. We assume that the probability of the first candidate is p for each state
#indpendently
#Try:
#GFvp(USEC(),1/2,x);
GFvp:=proc(L,p,x) local i:
sort(expand(mul(1-p+p*x^L[i],i=1..nops(L)))):
end:

#StatAnal(f,x): Given a generating function f in x
#finds (i) the average, (ii) the  standard-deviation (iii) the third-throgh-K SCALED moments about the mean. 
#(The r-th scaled moment about the mean of r.v. is the r-th moment divided by the r-th power of the standard-deviation)
#Try:
#StatAnal(GFv(USEC(),x), x,10);
StatAnal:=proc(f,x,K) local f1,i,mu,sig,M,r:

if min(seq(coeff(f,x,i),i=0..degree(f,x)))<0  or subs(x=1,f)=0 then
 print(`Can't be made into a probability generating function `):
 RETURN(FAIL):
fi:

#turn it into a probability genearting function
f1:=expand(f/subs(x=1,f)):

#mu is the average
mu:=subs(x=1,diff(f1,x)):

#We redefine f1 to be the probability generating function for the CENTRALIZED prob. generating function
#where the average is 1
f1:=expand(f1/x^mu):

#We keep applying the operation f->x*diff(f,x) in order to computer higher moments

#Doing it once
f1:=x*diff(f1,x):

#Doing it again (for the second moment)
f1:=x*diff(f1,x):

#The standar deviation is the square-root of the variance (alias second-moment about the mean)
sig:=sqrt(subs(x=1,f1)):

#The first two entries of M are the most important part, the expectation (alias mean, alias average) and the standard-deviation

M:=[mu,sig]:



#We now find the higher moments about the mean
for r from 3 to K do
f1:=expand(x*diff(f1,x)):
M:=[op(M),subs(x=1,f1)/sig^r]:
od:

#The final output is the list of length K as desired
M:

end:

#LC(p): Outputs 1 with probability p and 0 with prob. 1-p. Try: 
#LC(2/3);
LC:=proc(p) local a,b,ra:
a:=numer(p):
b:=denom(p):
ra:=rand(1..b)():

if ra<=a then
 1:
else
 0:
fi:
end:


#SimuCount1(L,p): Given a list L representing the  number of electors in nops(L) `states', assuming that
#the probability of Trump winning is p independently. Simulates ONE random count
#where every unpicked state is equally likely to be counted next. The output is the
#list of partial scores Try:
#SimuCount1(USEC(),1/2);

SimuCount1:=proc(L,p) local H,cu,L1,r,v:

#At the start the score is [0,0] and the current history is [[]]:
cu:=[0,0]:
H:=[[0,0]]:

#L1 is the list of not-yet-counted states, it starts with L
L1:=L:

while L1<>[] do

#we pick a random state to count
r:=rand(1..nops(L1))():

#We toss a loaded coin with probability p
v:=LC(p):

#We update the current score
if v=1 then
 cu:=cu+[L1[r],0]:
else
 cu:=cu+[0,L1[r]]:
fi:

#We append the current score to the history list
H:=[op(H),cu]:

#We remove the current state from the still-to-count list

L1:=[op(1..r-1,L1),op(r+1..nops(L1),L1)]:
od:

#We return the  history. The last element is the final score
H:
end:

#IsConsistent(H): Given a list of partial scores, outputs true if and only if the ultimate winner was always ahead if it
#ended with a tie we declare the first candiate to be the winner. Try:

IsConsistent:=proc(H) local cu,i:
cu:=H[nops(H)]:

if cu[1]>=cu[2] then
 if min(seq(H[i][1]-H[i][2],i=1..nops(H)-1))>=0 then
  RETURN(true):
 else
  RETURN(false):
fi:

else

 if min(seq(H[i][2]-H[i][1],i=1..nops(H)-1))>=0 then
  RETURN(true):
 else
  RETURN(false):
fi:

fi:

end:

# SimuCount(L,p,N,K): simulates N random counting with the list L and probability of the first candiate  winning being p
#The output is the list consisting of the average, standard deviation, and first K moments, followed by the
#ratio of consistend histories. Try:
#SimuCount(USEC(),1/2,1000,4);
SimuCount:=proc(L,p,N,K) local f,x,co,H,i:

#The initial count is of the number of consistend votes was 0
co:=0:

#We compute at the same time the generating function according to the number of votes of the first candidate  at the end
#and keep track of those  where the vote was consistent

f:=0:
for i from 1 to N do
 H:=SimuCount1(L,p):

 f:=f+x^H[nops(H)][1]:

 if IsConsistent(H) then
   co:=co+1:
 fi:
od:
 
evalf(StatAnal(f,x,K)),evalf(co/N):

end:

#StatAnalS(f,x): Given a generating function f in x with SYMBOLIC coefficients
#finds (i) the average, (ii) the  standard-deviation (iii) the third-throgh-K SCALED moments about the mean. 
#(The r-th scaled moment about the mean of r.v. is the r-th moment divided by the r-th power of the standard-deviation)
#Try:
#StatAnalS(GFv(USEC(),x), x,10);
StatAnalS:=proc(f,x,K) local f1,i,mu,sig,M,r:



#turn it into a probability genearting function
f1:=expand(f/subs(x=1,f)):

#mu is the average
mu:=subs(x=1,diff(f1,x)):

#We redefine f1 to be the probability generating function for the CENTRALIZED prob. generating function
#where the average is 1
f1:=expand(f1/x^mu):

#We keep applying the operation f->x*diff(f,x) in order to computer higher moments

#Doing it once
f1:=x*diff(f1,x):

#Doing it again (for the second moment)
f1:=x*diff(f1,x):

#The standar deviation is the square-root of the variance (alias second-moment about the mean)
sig:=sqrt(subs(x=1,f1)):

#The first two entries of M are the most important part, the expectation (alias mean, alias average) and the standard-deviation

M:=[mu,sig]:



#We now find the higher moments about the mean
for r from 3 to K do
f1:=expand(x*diff(f1,x)):
M:=[op(M),subs(x=1,f1)/sig^r]:
od:

#The final output is the list of length K as desired
M:

end: