################# non-independent_generators.sage
print('probability calculus for non independent generators (4 in this case) modeled on unit n-cube the vertices are atoms of statement calculus')
print('example variable description for 3-cube with generators A,B,C')
print('y000 = Prob(1), y100=Prob(A), y010=Prob(B), y001=Prob(C), y110=Prob(A*B), y011=Prob(B*C), y101=Prob(C*A), y111=Prob(A*B*C)')
print('x000 = Prob((1+A)*(1+B)*(1+C)), x100=Prob(A*(1+B)*(1+C)), x010=Prob((1+A)*B*(1+C)), x001=Prob((1+A)*(1+B)*C), x110=Prob(A*B*(1+C)), x011=Prob((1+A)*B*C), x101=Prob(A*(1+B)*C), x111=Prob(A*B*C)')
print(' ')
num_gens = 4 # 4-cube
num_atoms = 2**num_gens
print('# of generators ' + str(num_gens) + ', # of atoms ' + str(num_atoms)); print(' ')
x = ['x' + Integer(ii+num_atoms).binary()[1:num_gens + 1] for ii in range(num_atoms)];print('x probability variables');print(x);print(' ')
y = ['y' + Integer(ii+num_atoms).binary()[1:num_gens + 1] for ii in range(num_atoms)];print('y probability variables');print(y);print(' ')
MPR = PolynomialRing(QQ,x + y);print(MPR);print(' ')
MPR.inject_variables(); print(' ')
xvars = [MPR.gens()[k] for k in [i for i, j in enumerate(MPR.gens()) if (str(j)[0] == 'x')]];print('xvars');print(xvars);print(' ') # atoms correspond to the vertices of the unit n-cube
yvars = [MPR.gens()[k] for k in [i for i, j in enumerate(MPR.gens()) if (str(j)[0] == 'y')]];print('yvars');print(yvars);print(' ') # these correspond to expressions formed by generators and *, y101 = a*c, y111=a*b*c, y000 = 1 = x000+x001+x010+... also y101 = Prob(A*C) etc
#AND = lambda A, B: A*B
#FA = lambda mylist: 1 if not mylist else Integer(reduce(AND, mylist))
#MM = matrix([[Integer(FA([str(xvarsx)[i] == '1' for i in iyvars])) for xvarsx in xvars] for iyvars in iy]);print('MM');print('is MM totally unimodular? going out on limb and saying yes');print(MM);print(' ')
# just a note: Prob(A*~B) = Prob(A*(1+B)) = Prob(A+A*B)=Prob(A)-Prob(A*B)
#MMinv = MM.inverse();print('MMinv');print(MMinv);print(' ') # well almost, is MM totally unimodular? Integer matrix with integer inverse det = +/- 1
# recursive way of defining MM and MMinv
MM = matrix([1])
for ii in range(1,num_gens + 1):
MM = block_matrix(2,2,[MM,MM,0,MM],subdivide=False)
print('Sierpinski triangle Matrix, M');print(MM);print(' ')
print('y = M*x, is M totally unimodular? going out on limb and saying yes');
MMinv = matrix([1])
for ii in range(1,num_gens+1):
MMinv = block_matrix(2,2,[MMinv,-MMinv,0,MMinv],subdivide=False)
print('M^-1');print(MMinv);print(' ')
x_of_y = list(MMinv*vector(yvars) - vector(xvars));print('x variables in terms of y variables');print(x_of_y);print(' ')
y_of_x = list(MM*vector(xvars) - vector(yvars));print('y variables in terms of x variables');print(y_of_x);print(' ')
prob_atoms = [QQ(normalvariate(0, 1)^2) for i in range(num_atoms)]
normalize = sum(prob_atoms)
prob_atoms = [pa/normalize for pa in prob_atoms]
RR_prob_atoms = [RR(pa) for pa in prob_atoms]
print('random pick atom i.e. x variable probabilities'); print(RR_prob_atoms); print(' ')
y_var_prob = list(MM*vector(RR_prob_atoms)); print('y variable probabilities'); print(y_var_prob);print(' ')
atom_prob_polys = list(vector(xvars) - vector(prob_atoms)) #just assign rational probabilities to the X variables
print('sum probabilities check');print(sum(prob_atoms));print(' ')
iy = [[i for i,j in enumerate(str(a)) if j == '1'] for a in yvars]; print('iy'); print(iy); print(' ')
leniy = [len(c) for c in iy]
ind_poly_list = []
for ii in range(len(yvars)):
if leniy[ii] > 1:
str_yvars = str(yvars[ii])
gens_list = []
for jj in range(1,len(str_yvars)):
if str_yvars[jj] == '1':
gens_list.append(yvars[2^(num_gens - jj)])
ind_poly_list.append(reduce(lambda a,b: a*b, gens_list) - yvars[ii])
print('define polynomials for independent generators, generators=3 then y100,y010,y001 are the generators')
print(ind_poly_list);print(' ')
#myideal=ideal(ind_poly_list)
#QR = MPR.quotient_ring(myideal)
#QR.inject_variables()
#bool_poly = y0101*y1100+y1111*y1000*y0110+2*y1110; print(bool_poly) # make sure coefficients are in ZZ
bool_poly = MPR.remove_var(*tuple(xvars)).change_ring(ZZ).random_element(); print('random integer polynomial in y variables'); print(bool_poly); print(' ')
bool_poly = vector([QQ(mod(a,2)) for a in bool_poly.coefficients()])*vector(bool_poly.monomials()); print('treated as boolean polynomial') ;print(bool_poly); print(' ') #reduce modulo 2 (boolean polynomial remember) then coerce new coefficients back to QQ
parse_bool_poly = [[b.strip() for b in c] for c in [a.split('*') for a in str(bool_poly).split('+')]]; print('parsing and removing leading and trailing spaces'); print(parse_bool_poly); print(' ')
parse_bool_poly = [[a.split('^')[0] for a in c] for c in parse_bool_poly]; print('removing powers y^m --> y'); print(parse_bool_poly); print(' ')
listA = []
for pbp in parse_bool_poly:
listB = []
for b in pbp:
listB.append('1' + b[1:])
k = int(listB[0],2)
for b in listB:
k = k | int(b,2) # or-ing the bits in y101 with y110 gives y111 (A*C)*(A*B)=A*B*C for instance ie y101*y110=y111
listA.append('y' + Integer(k).binary()[1:])
#creating the equivalent reduced boolean polynomial
#y vars as string being converted back to the generators in the polynomial ring and we previously injected the variable so we can write expressions in the generators directly
print('or-ing the bits in y101 with y110 gives y111 (A*C)*(A*B)=A*B*C for instance ie y101*y110=y111')
reduced_bool_poly = sum([MPR(y) for y in listA]); print('reduced to 1st order - linear - boolean polynomial in yvariables'); print(reduced_bool_poly); print(' ')
#now the procedure to obtain the probability expression
vctr = vector([reduced_bool_poly.coefficient(g) for g in yvars])
atom_bool_poly = vctr*MM*vector(xvars); print('expressing in the x variables - still boolean - using vector(y var coefficients)*M*vector(xvars)'); print(atom_bool_poly); print(' ')
atom_bool_poly = sum([QQ(mod(atom_bool_poly.coefficient(a),2))*a for a in xvars]); print('reducing the x var polynomial'); print(atom_bool_poly); print(' ')
prob_poly = vector([atom_bool_poly.coefficient(a) for a in xvars])*MMinv*vector(yvars); print('converting to probability using xvar_coeff_vector*(M^-1)*vector(yvars)');print(prob_poly); print(' ')
#for independent generators prob_poly needs to be linear in the y variables use this with the output from above (testpoly4)
prob_poly_monomials = prob_poly.monomials(); print('prob_poly_monomials'); print(prob_poly_monomials); print(' ')
prob_poly_coefficients = prob_poly.coefficients(); print('prob_poly_coefficients'); print(prob_poly_coefficients); print(' ')
monomial_gens_list = []
for a in prob_poly_monomials:
str_a = str(a)
gens_list = []
for jj in range(1,len(str_a)):
if str_a[jj] == '1':
gens_list.append(yvars[2^(num_gens - jj)])
if len(gens_list) <> 0:
monomial_gens_list.append(reduce(lambda a,b:a*b ,gens_list))
else:
monomial_gens_list.append(yvars[0]) #yvars[0]=y000 then gens_list = [] so we put yvars[0] back
prob_poly_ind = vector(prob_poly_coefficients)*vector(monomial_gens_list); print('if the generators are independent the probability poly becomes'); print(prob_poly_ind); print(' ')
#lagrange multiplier problem
print('approximation to maximum entropy using least squares subject to constraints:')
print('minimize sum of squares of atom probabilities to an equiprobable distribution subject to constraints sum(prob(atoms)) = 1.0 and prob(boolean_polynomial) = observed_value'); print(' ')
#atom_bool_poly2 = vector([prob_poly.coefficient(a) for a in yvars])*MM*vector(xvars) # should equal atom_bool_poly
#atom_bool_poly_coefficients = matrix([atom_bool_poly2.coefficient(a) for a in xvars])
atom_bool_poly_coefficients = matrix([atom_bool_poly.coefficient(a) for a in xvars]); print('coefficient vector of boolean polynomial in terms of atoms - xvars'); print(atom_bool_poly_coefficients); print(' ')
observed_prob_of_bool_poly = QQ(round(.05 + 0.9*random(), 3)); print('suppose observed probability of boolean polynomial = ' + str(observed_prob_of_bool_poly)); print(' ')
lgr = block_matrix(2,1,[atom_bool_poly_coefficients, matrix([1 for i in range(num_atoms)])], subdivide = False)
lagrange_matrix = block_matrix(QQ, 2, 2, [2*identity_matrix(num_atoms,num_atoms), lgr.transpose(), lgr, zero_matrix(2,2)],subdivide = False); print('lagrange multiplier matrix'); print(lagrange_matrix); print(' ')
lagrange_inverse = lagrange_matrix.inverse();print('inverse of lagrange matrix'); print(lagrange_inverse); print(' ')
lagrange_vector = vector([QQ(2/num_atoms) for i in range(num_atoms)] + [observed_prob_of_bool_poly, 1])
print('left hand side of vector([2*equiprobable, observed_prob, 1]) = lagrange_matrix*vector([prob(atoms), lambda, mu])'); print(lagrange_vector); print(' ')
atom_probabilities = (lagrange_inverse*lagrange_vector)[0:num_atoms]; print('atom probabilities that minimize sum of squares subject to observed probability of boolean polynomial'); print(atom_probabilities); print(' ')
# some output
probability calculus for non independent generators (4 in this case) modeled on unit n-cube the vertices are atoms of statement calculus example variable description for 3-cube with generators A,B,C y000 = Prob(1), y100=Prob(A), y010=Prob(B), y001=Prob(C), y110=Prob(A*B), y011=Prob(B*C), y101=Prob(C*A), y111=Prob(A*B*C) x000 = Prob((1+A)*(1+B)*(1+C)), x100=Prob(A*(1+B)*(1+C)), x010=Prob((1+A)*B*(1+C)), x001=Prob((1+A)*(1+B)*C), x110=Prob(A*B*(1+C)), x011=Prob((1+A)*B*C), x101=Prob(A*(1+B)*C), x111=Prob(A*B*C) # of generators 4, # of atoms 16 x probability variables ['x0000', 'x0001', 'x0010', 'x0011', 'x0100', 'x0101', 'x0110', 'x0111', 'x1000', 'x1001', 'x1010', 'x1011', 'x1100', 'x1101', 'x1110', 'x1111'] y probability variables ['y0000', 'y0001', 'y0010', 'y0011', 'y0100', 'y0101', 'y0110', 'y0111', 'y1000', 'y1001', 'y1010', 'y1011', 'y1100', 'y1101', 'y1110', 'y1111'] Multivariate Polynomial Ring in x0000, x0001, x0010, x0011, x0100, x0101, x0110, x0111, x1000, x1001, x1010, x1011, x1100, x1101, x1110, x1111, y0000, y0001, y0010, y0011, y0100, y0101, y0110, y0111, y1000, y1001, y1010, y1011, y1100, y1101, y1110, y1111 over Rational Field Defining x0000, x0001, x0010, x0011, x0100, x0101, x0110, x0111, x1000, x1001, x1010, x1011, x1100, x1101, x1110, x1111, y0000, y0001, y0010, y0011, y0100, y0101, y0110, y0111, y1000, y1001, y1010, y1011, y1100, y1101, y1110, y1111 xvars [x0000, x0001, x0010, x0011, x0100, x0101, x0110, x0111, x1000, x1001, x1010, x1011, x1100, x1101, x1110, x1111] yvars [y0000, y0001, y0010, y0011, y0100, y0101, y0110, y0111, y1000, y1001, y1010, y1011, y1100, y1101, y1110, y1111] Sierpinski triangle Matrix, M [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1] [0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1] [0 0 1 1 0 0 1 1 0 0 1 1 0 0 1 1] [0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 1] [0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1] [0 0 0 0 0 1 0 1 0 0 0 0 0 1 0 1] [0 0 0 0 0 0 1 1 0 0 0 0 0 0 1 1] [0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1] [0 0 0 0 0 0 0 0 0 1 0 1 0 1 0 1] [0 0 0 0 0 0 0 0 0 0 1 1 0 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 1] [0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1] y = M*x, is M totally unimodular? going out on limb and saying yes M^-1 [ 1 -1 -1 1 -1 1 1 -1 -1 1 1 -1 1 -1 -1 1] [ 0 1 0 -1 0 -1 0 1 0 -1 0 1 0 1 0 -1] [ 0 0 1 -1 0 0 -1 1 0 0 -1 1 0 0 1 -1] [ 0 0 0 1 0 0 0 -1 0 0 0 -1 0 0 0 1] [ 0 0 0 0 1 -1 -1 1 0 0 0 0 -1 1 1 -1] [ 0 0 0 0 0 1 0 -1 0 0 0 0 0 -1 0 1] [ 0 0 0 0 0 0 1 -1 0 0 0 0 0 0 -1 1] [ 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 -1] [ 0 0 0 0 0 0 0 0 1 -1 -1 1 -1 1 1 -1] [ 0 0 0 0 0 0 0 0 0 1 0 -1 0 -1 0 1] [ 0 0 0 0 0 0 0 0 0 0 1 -1 0 0 -1 1] [ 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 -1] [ 0 0 0 0 0 0 0 0 0 0 0 0 1 -1 -1 1] [ 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 -1] [ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 -1] [ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1] x variables in terms of y variables [-x0000 + y0000 - y0001 - y0010 + y0011 - y0100 + y0101 + y0110 - y0111 - y1000 + y1001 + y1010 - y1011 + y1100 - y1101 - y1110 + y1111, -x0001 + y0001 - y0011 - y0101 + y0111 - y1001 + y1011 + y1101 - y1111, -x0010 + y0010 - y0011 - y0110 + y0111 - y1010 + y1011 + y1110 - y1111, -x0011 + y0011 - y0111 - y1011 + y1111, -x0100 + y0100 - y0101 - y0110 + y0111 - y1100 + y1101 + y1110 - y1111, -x0101 + y0101 - y0111 - y1101 + y1111, -x0110 + y0110 - y0111 - y1110 + y1111, -x0111 + y0111 - y1111, -x1000 + y1000 - y1001 - y1010 + y1011 - y1100 + y1101 + y1110 - y1111, -x1001 + y1001 - y1011 - y1101 + y1111, -x1010 + y1010 - y1011 - y1110 + y1111, -x1011 + y1011 - y1111, -x1100 + y1100 - y1101 - y1110 + y1111, -x1101 + y1101 - y1111, -x1110 + y1110 - y1111, -x1111 + y1111] y variables in terms of x variables [x0000 + x0001 + x0010 + x0011 + x0100 + x0101 + x0110 + x0111 + x1000 + x1001 + x1010 + x1011 + x1100 + x1101 + x1110 + x1111 - y0000, x0001 + x0011 + x0101 + x0111 + x1001 + x1011 + x1101 + x1111 - y0001, x0010 + x0011 + x0110 + x0111 + x1010 + x1011 + x1110 + x1111 - y0010, x0011 + x0111 + x1011 + x1111 - y0011, x0100 + x0101 + x0110 + x0111 + x1100 + x1101 + x1110 + x1111 - y0100, x0101 + x0111 + x1101 + x1111 - y0101, x0110 + x0111 + x1110 + x1111 - y0110, x0111 + x1111 - y0111, x1000 + x1001 + x1010 + x1011 + x1100 + x1101 + x1110 + x1111 - y1000, x1001 + x1011 + x1101 + x1111 - y1001, x1010 + x1011 + x1110 + x1111 - y1010, x1011 + x1111 - y1011, x1100 + x1101 + x1110 + x1111 - y1100, x1101 + x1111 - y1101, x1110 + x1111 - y1110, x1111 - y1111] random pick atom i.e. x variable probabilities [0.254012949353045, 0.0319311459788169, 0.0571729369657601, 0.0786812600291195, 0.188655212482421, 0.00990680542001836, 0.0464720011232607, 0.123391437395270, 0.0107745913537673, 0.0932831908041969, 0.0140401085923424, 0.0125507437727733, 0.0276594050866677, 0.0495474821267312, 0.00126746343620038, 0.000653266079608370] y variable probabilities [1.00000000000000, 0.399945331606535, 0.334229217394335, 0.215276707276771, 0.447553073150179, 0.183498991021628, 0.171784168034340, 0.124044703474879, 0.209776251252288, 0.156034682783310, 0.0285115818809245, 0.0132040098523817, 0.0791276167292077, 0.0502007482063396, 0.00192072951580875, 0.000653266079608370] sum probabilities check 1 define polynomials for independent generators, generators=3 then y100,y010,y001 are the generators [y0001*y0010 - y0011, y0001*y0010 - y0101, y0001*y0010 - y0110, y0001*y0010 - y0111, y0001*y0010 - y1001, y0001*y0010 - y1010, y0001*y0010 - y1011, y0001*y0010 - y1100, y0001*y0010 - y1101, y0001*y0010 - y1110, y0001*y0010 - y1111] iy [[], [4], [3], [3, 4], [2], [2, 4], [2, 3], [2, 3, 4], [1], [1, 4], [1, 3], [1, 3, 4], [1, 2], [1, 2, 4], [1, 2, 3], [1, 2, 3, 4]] random integer polynomial in y variables y0101*y1001 - 5*y1010^2 + y0001*y1110 treated as boolean polynomial y0101*y1001 + y1010^2 + y0001*y1110 parsing and removing leading and trailing spaces [['y0101', 'y1001'], ['y1010^2'], ['y0001', 'y1110']] removing powers y^m --> y [['y0101', 'y1001'], ['y1010'], ['y0001', 'y1110']] or-ing the bits in y101 with y110 gives y111 (A*C)*(A*B)=A*B*C for instance ie y101*y110=y111 reduced to 1st order - linear - boolean polynomial in y variables y1010 + y1101 + y1111 expressing in the x variables - still boolean - using vector(y var coefficients)*M*vector(xvars) x1010 + x1011 + x1101 + x1110 + 3*x1111 reducing the x var polynomial x1010 + x1011 + x1101 + x1110 + x1111 converting to probability using vector(x var coefficients)*(M^-1)*vector(yvars) y1010 + y1101 - y1111 prob_poly_monomials [y1010, y1101, y1111] prob_poly_coefficients [1, 1, -1] if the generators are independent the probability poly becomes -y0001*y0010*y0100*y1000 + y0001*y0100*y1000 + y0010*y1000
approximation to maximum entropy using least squares subject to constraints: minimize sum of squares of atom probabilities to an equiprobable distribution subject to constraints sum(prob(atoms)) = 1.0 and prob(boolean_polynomial) = observed_value coefficient vector of boolean polynomial in terms of atoms - xvars [0 0 0 0 0 0 0 0 0 0 1 1 0 1 1 1] suppose observed probability of boolean polynomial = 299/1000 lagrange multiplier matrix [2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1] [0 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1] [0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1] [0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 0 1] [0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 1] [0 0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 0 2 0 0 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 0 0 2 0 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 0 0 0 2 0 0 0 0 0 0 0 1] [0 0 0 0 0 0 0 0 0 0 2 0 0 0 0 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 2 0 0 0 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 0 0 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 1 1] [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 1 1] [0 0 0 0 0 0 0 0 0 0 1 1 0 1 1 1 0 0] [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0] inverse of lagrange matrix [ 5/11 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 5/11 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 5/11 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 5/11 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 5/11 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 -1/22 5/11 -1/22 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 5/11 -1/22 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 5/11 -1/22 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 5/11 -1/22 0 0 -1/22 0 0 0 -1/11 1/11] [ -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 5/11 0 0 -1/22 0 0 0 -1/11 1/11] [ 0 0 0 0 0 0 0 0 0 0 2/5 -1/10 0 -1/10 -1/10 -1/10 1/5 0] [ 0 0 0 0 0 0 0 0 0 0 -1/10 2/5 0 -1/10 -1/10 -1/10 1/5 0] [ -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 -1/22 0 0 5/11 0 0 0 -1/11 1/11] [ 0 0 0 0 0 0 0 0 0 0 -1/10 -1/10 0 2/5 -1/10 -1/10 1/5 0] [ 0 0 0 0 0 0 0 0 0 0 -1/10 -1/10 0 -1/10 2/5 -1/10 1/5 0] [ 0 0 0 0 0 0 0 0 0 0 -1/10 -1/10 0 -1/10 -1/10 2/5 1/5 0] [ -1/11 -1/11 -1/11 -1/11 -1/11 -1/11 -1/11 -1/11 -1/11 -1/11 1/5 1/5 -1/11 1/5 1/5 1/5 -32/55 2/11] [ 1/11 1/11 1/11 1/11 1/11 1/11 1/11 1/11 1/11 1/11 0 0 1/11 0 0 0 2/11 -2/11] left hand side of vector([2*equiprobable, observed_prob, 1]) = lagrange_matrix*vector([prob(atoms), lambda, mu]) (1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 1/8, 299/1000, 1) atom probabilities that minimize sum of squares subject to observed probability of boolean polynomial (701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 701/11000, 299/5000, 299/5000, 701/11000, 299/5000, 299/5000, 299/5000)