Skip to content

Commit 240a63f

Browse files
steguastegua
committed
Update solution for xor
1 parent 3096c92 commit 240a63f

File tree

2 files changed

+197
-0
lines changed

2 files changed

+197
-0
lines changed

‎README.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,7 @@ This repository is maintained by the [Computational Optimization Research Group]
2020

2121
| Data | Solution | Link |
2222
|:-|:-|:-|
23+
|**[2024/05/03]**|*Training a BNN for the XOR logical function (non linearly separable*|[NN_and.py](https://github.com/mathcoding/opt4ds/blob/master/scripts/Xor_aula.py)|
2324
|**[2024/04/22]**|*Training a BNN for the AND logical function*|[NN_and.py](https://github.com/mathcoding/opt4ds/blob/master/scripts/NN_and.py)|
2425
|**[2024/04/15]**|*Optimal Color Transfer*|[colorTransfer.py](https://github.com/mathcoding/opt4ds/blob/master/scripts/colorTransfer.py)|
2526
|**[2024/04/12]**|*Linear regression Diabete dataset*|[regression_diabete.py](https://github.com/mathcoding/opt4ds/blob/master/scripts/regression_diabete.py)|

‎scripts/Xor_aula.py

Lines changed: 196 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,196 @@
1+
from gurobipy import Model, GRB, quicksum
2+
3+
def AddBias(X):
4+
return [x for x in X] + [1]
5+
6+
def Sign(x):
7+
return 1 if x >= 0.01 else -1
8+
9+
def MIPxor(Xs, Ys, Phi=lambda x: x, NH=2):
10+
Xs = [Phi(x) for x in Xs]
11+
12+
m = len(Xs[0])
13+
n = len(Xs)
14+
15+
nh = NH
16+
17+
# Create model ILP
18+
model = Model()
19+
model.setParam(GRB.Param.TimeLimit, 5) # In seconds
20+
model.setParam(GRB.Param.OutputFlag, 1) # 0: silent, 1: normal, 2: verbose
21+
22+
wp, wn = {}, {}
23+
for i in range(m):
24+
for h in range(nh):
25+
wp[i, h] = model.addVar(obj=0.001, vtype=GRB.BINARY)
26+
wn[i, h] = model.addVar(obj=0.001, vtype=GRB.BINARY)
27+
28+
up = [model.addVar(obj=0.001, vtype=GRB.BINARY) for h in range(nh)]
29+
un = [model.addVar(obj=0.001, vtype=GRB.BINARY) for h in range(nh)]
30+
31+
# Bias on output
32+
ubp = model.addVar(obj=0.001, vtype=GRB.BINARY)
33+
ubn = model.addVar(obj=0.001, vtype=GRB.BINARY)
34+
35+
z = {}
36+
vp, vn = {}, {}
37+
for k in range(n):
38+
for h in range(nh):
39+
z[k, h] = model.addVar(vtype=GRB.BINARY)
40+
# Variable to linearize the product of binary variables
41+
vp[k, h] = model.addVar(vtype=GRB.BINARY)
42+
vn[k, h] = model.addVar(vtype=GRB.BINARY)
43+
44+
y_hat = [model.addVar(vtype=GRB.BINARY) for k in range(n)]
45+
alpha = [model.addVar(lb=0, obj=1) for k in range(n)]
46+
47+
# First layer constraints
48+
M = 1000
49+
for k in range(n):
50+
for h in range(nh):
51+
model.addConstr( quicksum(Xs[k][i]*(wp[i,h]-wn[i,h]) for i in range(m)) >= 0.01 - M*(1 - z[k,h]) )
52+
model.addConstr( quicksum(Xs[k][i]*(wp[i,h]-wn[i,h]) for i in range(m)) <= 0.00 + M*z[k,h] )
53+
54+
# Second layer constraints
55+
for k in range(n):
56+
# (2z-1)*(up - un) = 2z*up - 2z*un - up + un
57+
model.addConstr( (ubp - ubn) + quicksum((2*vp[k,h] - 2*vn[k,h] - up[h] + un[h]) for h in range(nh)) >= 0.01 - M*(1 - y_hat[k]))
58+
model.addConstr( (ubp - ubn) + quicksum((2*vp[k,h] - 2*vn[k,h] - up[h] + un[h]) for h in range(nh)) <= 0.00 + M*y_hat[k])
59+
60+
for k in range(n):
61+
for h in range(nh):
62+
model.addConstr( vp[k,h] >= z[k,h] + up[h] - 1 )
63+
model.addConstr( vp[k,h] <= z[k,h] )
64+
model.addConstr( vp[k,h] <= up[h] )
65+
66+
model.addConstr( vn[k,h] >= z[k,h] + un[h] - 1 )
67+
model.addConstr( vn[k,h] <= z[k,h] )
68+
model.addConstr( vn[k,h] <= un[h] )
69+
70+
for k in range(n):
71+
model.addConstr( (2*y_hat[k]-1) - Ys[k] <= alpha[k])
72+
model.addConstr( Ys[k] - (2*y_hat[k]-1) <= alpha[k])
73+
74+
model.optimize()
75+
76+
if model.status != GRB.OPTIMAL and model.status != GRB.TIME_LIMIT:
77+
return None
78+
79+
wbar = {(i, h): wp[i, h].X - wn[i,h].X for i, h in wp}
80+
ubar = [up[h].X - un[h].X for h in range(nh)]
81+
ubias = ubp.X - ubn.X
82+
print('wbar nonzero', sum(wbar[i,h] != 0 for i, h in wbar))
83+
print('ubar nonzero', sum(u != 0 for u in ubar) + (ubias != 0))
84+
85+
def F(x):
86+
z = Phi(x)
87+
return Sign(sum(ubar[h]*Sign(sum(z[i]*wbar[i,h] for i in range(m))) for h in range(nh)))
88+
89+
return F
90+
91+
def MIPxor_noncnvex(Xs, Ys, Phi=lambda x: x, NH=2):
92+
Xs = [Phi(x) for x in Xs]
93+
94+
m = len(Xs[0])
95+
n = len(Xs)
96+
97+
nh = NH
98+
99+
# Create model ILP
100+
model = Model()
101+
model.setParam(GRB.Param.TimeLimit, 5) # In seconds
102+
model.setParam(GRB.Param.OutputFlag, 1) # 0: silent, 1: normal, 2: verbose
103+
104+
# For automatic linearization techniques
105+
model.setParam(GRB.Param.NonConvex, 2)
106+
107+
wp, wn = {}, {}
108+
for i in range(m):
109+
for h in range(nh):
110+
wp[i, h] = model.addVar(obj=0.001, vtype=GRB.BINARY)
111+
wn[i, h] = model.addVar(obj=0.001, vtype=GRB.BINARY)
112+
113+
up = [model.addVar(obj=0.001, vtype=GRB.BINARY) for h in range(nh)]
114+
un = [model.addVar(obj=0.001, vtype=GRB.BINARY) for h in range(nh)]
115+
116+
z = {}
117+
for k in range(n):
118+
for h in range(nh):
119+
z[k, h] = model.addVar(vtype=GRB.BINARY)
120+
121+
y_hat = [model.addVar(vtype=GRB.BINARY) for k in range(n)]
122+
alpha = [model.addVar(lb=0, obj=1) for k in range(n)]
123+
124+
# First layer constraints
125+
M = 1000
126+
for k in range(n):
127+
for h in range(nh):
128+
model.addConstr( quicksum(Xs[k][i]*(wp[i,h]-wn[i,h]) for i in range(m)) >= 0.01 - M*(1 - z[k,h]) )
129+
model.addConstr( quicksum(Xs[k][i]*(wp[i,h]-wn[i,h]) for i in range(m)) <= 0.00 + M*z[k,h] )
130+
131+
# Second layer constraints
132+
for k in range(n):
133+
model.addConstr( quicksum(z[k,h]*(up[h] - un[h]) for h in range(nh)) >= 0.01 - M*(1 - y_hat[k]))
134+
model.addConstr( quicksum(z[k,h]*(up[h] - un[h]) for h in range(nh)) <= 0.00 + M*y_hat[k])
135+
136+
for k in range(n):
137+
model.addConstr( (2*y_hat[k]-1) - Ys[k] <= alpha[k])
138+
model.addConstr( Ys[k] - (2*y_hat[k]-1) <= alpha[k])
139+
140+
model.optimize()
141+
142+
if model.status != GRB.OPTIMAL and model.status != GRB.TIME_LIMIT:
143+
return None
144+
145+
wbar = {(i, h): wp[i, h].X - wn[i,h].X for i, h in wp}
146+
ubar = [up[h].X - un[h].X for h in range(nh)]
147+
148+
print('wbar nonzero', sum(wbar[i,h] != 0 for i, h in wbar))
149+
print('ubar nonzero', sum(u != 0 for u in ubar))
150+
151+
def F(x):
152+
z = Phi(x)
153+
return Sign(sum(ubar[h]*Sign(sum(z[i]*wbar[i,h] for i in range(m))) for h in range(nh)))
154+
155+
return F
156+
157+
158+
from numpy.random import normal, seed
159+
seed(13)
160+
def AddNoise(X, mu=0.1):
161+
return list(map(lambda x: x + normal(0, mu), X))
162+
163+
Xor = [(-1, -1), (-1, 1), (1, -1), (1, 1)]
164+
Yor = [-1, 1, 1, -1]
165+
166+
# Sample points
167+
Xtrain, Ytrain = [], []
168+
for _ in range(250):
169+
Xtrain.extend([AddNoise(x) for x in Xor])
170+
Ytrain.extend(Yor)
171+
172+
Xtest, Ytest = [], []
173+
for _ in range(250):
174+
Xtest.extend([AddNoise(x) for x in Xor])
175+
Ytest.extend(Yor)
176+
177+
nh=2
178+
F1 = MIPxor(Xtrain, Ytrain, Phi=AddBias, NH=nh)
179+
F2 = MIPxor(Xtrain, Ytrain, Phi=AddBias, NH=nh)
180+
181+
acc1, acc2 = 0, 0
182+
for x, y in zip(Xtest, Ytest):
183+
acc1 = acc1 + (F1(x) == y)
184+
acc2 = acc2 + (F2(x) == y)
185+
# print(x, F(x), y)
186+
187+
print('Accuracy:', acc1/len(Xtest), len(Xtest), len(Xtrain), 'NH',nh)
188+
print('Accuracy Non convex:', acc2ù/len(Xtest), len(Xtest), len(Xtrain), 'NH',nh)
189+
190+
191+
192+
193+
194+
195+
196+

0 commit comments

Comments
 (0)