Dateien nach "ILP" hochladen

Now ILP solver picks optimal combination of two different similarities and similarity can also be calcuated using the mean of similarity instead of the comparitive value

ILP/caffeinesynthesis.py

@@ -33,6 +33,10 @@ VERTICES1 = {
     8: ([0, 0], ['Caffeine']),
 }
 
+#25014CX3
+
+#03
+
 VERTICES2 = {
     1: ([0], ['Xanthine']),
     2: ([0], ['p_{0,0}']),
@@ -47,69 +51,107 @@ VERTICES2 = {
 VERTICES = ['Xanthine', 'p_{0,0}', 'p_{0,1}', 'p_{0,2}', 'p_{0,3}', 'p_{0,4}', 'p_{0,5}',  'Caffeine']
 
 NMRLIKELYHOODS = [0.0, 0.2, 0.7, 0.1, 0.1, 0.2, 0.7, 0.0]
-NMRLIKELYHOODS1 = [0.75, 0.59, 0.36, 0.89, 0.59, 0.86, 0.78, 0.65]
+NMRLIKELYHOODS1 = [0.54, 0.52, 0.43, 0.62, 0.52, 0.56, 0.48, 0.5]
-NMRLIKELYHOODS2 = [0.56, 0.75, 0.45, 0.79, 0.75, 0.76, 0.84, 0.74]
+NMRLIKELYHOODS2 = [0.52, 0.52, 0.52, 0.61, 0.52, 0.56, 0.61, 0.59]
-NMRLIKELYHOODS3 = [0.54, 0.53, 0.7, 0.77, 0.53, 0.85, 0.94, 0.79]
+NMRLIKELYHOODS3 = [0.5, 0.53, 0.59, 0.59, 0.53, 0.56, 0.68, 0.62]
 
 FIXED_FLOWS = {
     1: 1,
     14: 1,
 }
 
-def build_model(name, hyperedges, vertices, nmrlikelihoods, excluded_support=None):
+def build_model(name, hyperedges, vertices, nmrlikelihoods1, nmrlikelihoods2, nmrlikelihoods3, excluded_support=None):
     model = Model(name)
 
-    x = {e_id: model.addVar(vtype=GRB.INTEGER, lb=0, name=f"x_{e_id}") for e_id in hyperedges}
+    x = {e_id: model.addVar(vtype=GRB.INTEGER, lb = 0, name = f"x_{e_id}") for e_id in hyperedges}
-    b = {e_id: model.addVar(vtype=GRB.BINARY, name=f"b_{e_id}") for e_id in hyperedges}
+    b = {e_id: model.addVar(vtype=GRB.BINARY, name = f"b_{e_id}") for e_id in hyperedges}
-    n = model.addVars(vertices, vtype=GRB.CONTINUOUS, lb=0.0, ub=1.0, name="nmr")
+    n1 = model.addVars(vertices, vtype=GRB.CONTINUOUS, lb = 0.0, ub = 1.0, name = "nmr1")
+    n2 = model.addVars(vertices, vtype=GRB.CONTINUOUS, lb = 0.0, ub = 1.0, name = "nmr2")
+    n3 = model.addVars(vertices, vtype=GRB.CONTINUOUS, lb = 0.0, ub = 1.0, name = "nmr3")
+    c1 = model.addVars(vertices, vtype = GRB.BINARY, name = "n1_choice")
+    c2 = model.addVars(vertices, vtype = GRB.BINARY, name = "n2_choice")
+    c3 = model.addVars(vertices, vtype = GRB.BINARY, name = "n3_choice")
+    
 
-    for v, nmr in zip(vertices, nmrlikelihoods):
+    for v, nmr in zip(vertices, nmrlikelihoods1):
-        n[v] = nmr
+        n1[v] = nmr
-        print(f'Vertice: {v}, Similarity: {n[v]}')
+        #print(f'Vertice: {v}, Similarity: {n1[v]}')
+    for v, nmr in zip(vertices, nmrlikelihoods2):
+        n2[v] = nmr
+
+    for v, nmr in zip(vertices, nmrlikelihoods3):
+        n3[v] = nmr
 
     vertices = set(v for tails, heads in hyperedges.values() for v in tails + heads)
 
     for v in vertices:
         inflow = quicksum(x[e_id] for e_id, (_, heads) in hyperedges.items() if v in heads)
         outflow = quicksum(x[e_id] for e_id, (tails, _) in hyperedges.items() if v in tails)
-        model.addConstr(inflow == outflow, name=f"flow_conservation_{v}")
+        model.addConstr(inflow == outflow, name = f"flow_conservation_{v}")
 
     for e_id, value in FIXED_FLOWS.items():
-        model.addConstr(x[e_id] == value, name=f"fixed_flow_{e_id}")
+        model.addConstr(x[e_id] == value, name = f"fixed_flow_{e_id}")
 
     for e_id in hyperedges:
-        model.addGenConstrIndicator(b[e_id], 0, x[e_id] == 0, name=f"unused_implies_zero_{e_id}")
+        model.addGenConstrIndicator(b[e_id], 0, x[e_id] == 0, name = f"unused_implies_zero_{e_id}")
-        model.addConstr(x[e_id] >= b[e_id], name=f"used_implies_positive_flow_{e_id}")
+        model.addConstr(x[e_id] >= b[e_id], name = f"used_implies_positive_flow_{e_id}")
 
     reaction_path = {}
 
-    #for v_id in vertices:
-        #model.addConstr()
 
     if excluded_support:
-        model.addConstr(quicksum(b[e_id] for e_id in excluded_support) <= len(excluded_support) - 1, name="different_hyperedges",)
+        model.addConstr(quicksum(b[e_id] for e_id in excluded_support) <= len(excluded_support) - 1, name = "different_hyperedges",)
     
     #Multiplizier den node Wert mit infow + outflow
-
-    #model.setObjective(quicksum(n[v_id] for v_id in vertices), GRB.MINIMIZE) #Maximize Similarity of Nodes
     model.ModelSense = GRB.MAXIMIZE
+    #Multiply with a list of binarys where only one value is one and the rest zero and not both the same position
+    '''
     model.setObjectiveN(
-        quicksum(n[t_id[0]] * x[e_id] for e_id, (_, t_id) in hyperedges.items() if t_id != []),
+        quicksum(n1[v_id] * c1[v_id] for v_id in vertices if v_id != []) 
-        index=0,
+        + quicksum(n3[v_id] * c3[v_id] for v_id in vertices if v_id != []),
-        priority=2,
+        index = 0,
-        name="maximize_nmr_similarity",
+        priority = 2,
+        name = "maximize_nmr_similarity",
     )   
+    '''
 
+    #Multiply node value with infow or outflow
+    '''
+    model.setObjectiveN(
+        quicksum(n1[t_id[0]] * c1[t_id[0]] * x[e_id] for e_id, (_, t_id) in hyperedges.items() if t_id != [])
+        +quicksum(n3[t_id[0]] * c3[t_id[0]] * x[e_id] for e_id, (_, t_id) in hyperedges.items() if t_id != []),
+        index = 0,
+        priority = 2,
+        name = "maximize_nmr_similarity",
+    ) 
+    '''
+
+    #Hybrid of both above
+    model.setObjective(
+        quicksum(n1[t_id[0]] * c1[t_id[0]] * x[e_id] for e_id, (_, t_id) in hyperedges.items() if t_id != [])
+        +quicksum(n3[t_id[0]] * c3[t_id[0]] * x[e_id] for e_id, (_, t_id) in hyperedges.items() if t_id != []),
+        sense = GRB.MAXIMIZE,
+    ) 
+    
+
+    model.addConstr(gp.quicksum(c1[v_id] for v_id in vertices) == 1)
+    model.addConstr(gp.quicksum(c3[v_id] for v_id in vertices) == 1)
+    for v_id in vertices:
+        model.addConstr(c1[v_id] + c3[v_id] <= 1)
+
+    model.addConstr(gp.quicksum(x[e_id] for e_id, (_, t_id) in hyperedges.items()) == 5)
+    ''''
     model.setObjectiveN(
         quicksum(-1 *b[e_id] for e_id in hyperedges),
         index=1,
         priority=1,
         name="minimize_used_hyperedges",
     )
+    '''
 
     return model, x, b
 
-def positive_entries(variable_dict, threshold=0.5):
+def positive_entries(variable_dict, threshold = 0.5):
     return {e_id: var.X for e_id, var in variable_dict.items() if var.X > threshold}
 
 def print_solution(title, flow_solution, binary_solution, hyperedges):
@@ -126,7 +168,7 @@ def print_solution(title, flow_solution, binary_solution, hyperedges):
     print(f"Number of used hyperedges: {len(binary_solution)}")
 
 def main():
-    model, x, b = build_model("HypergraphFlow", HYPEREDGES, VERTICES, NMRLIKELYHOODS)
+    model, x, b = build_model("HypergraphFlow", HYPEREDGES, VERTICES, NMRLIKELYHOODS1, NMRLIKELYHOODS2, NMRLIKELYHOODS3)
     model.optimize()
     if model.status != GRB.Status.OPTIMAL:
         print("No optimal solution found for the first model.")

ILP/nmrSimilarity.py

@@ -232,10 +232,10 @@ CNMR2 = {
 
 #Combination of methyl group and base purine rings from two papers
 CNMR3 = {
-    1: ([153.1], [1]),
+    1: ([154.9], [1]),
     2: ([150.0], [1]),
     3: ([108.1], [1]),
-    4: ([154.9], [1]),
+    4: ([153.1], [1]),
     5: ([142.6], [1]),
     6: ([29.3], [1]),
     7: ([33.9], [1]),
@@ -268,7 +268,7 @@ def define_border_values(spectraref, spectranew, bin_width):
         shifts.append(shift[0])
     highest_ppm = math.ceil(max(shifts))
     lowest_ppm = math.floor(min(shifts))
-    lowest_ppm = min(shifts) - bin_width/2 #Worse result. None of the previously wrong (except 0.6) become right
+    #lowest_ppm = min(shifts) - bin_width/2 #Worse result. None of the previously wrong (except 0.6) become right
     return (lowest_ppm, highest_ppm)
 
 def similarity_nmr(spectraref, spectranew, bin_width):
@@ -294,21 +294,25 @@ def correction(spectra, corretionppm):
     return newspectra
 
 def main():
-    spectrumref = CNMR3
+    spectrumref = CCAFFEINE2
-    spectrafalse = [CXANTHINE, C1XANTHINE, C3XANTHINE, C7XANTHINE, C1XANTHINE, C17XANTHINE, C37XANTHINE, C137XANTHINE]
+    spectra = [CXANTHINE, C1XANTHINE, C3XANTHINE, C7XANTHINE, C1XANTHINE, C17XANTHINE, C37XANTHINE, C137XANTHINE]
+    spectranames = ["CXANTHINE", "C1XANTHINE", "C3XANTHINE", "C7XANTHINE", "C1XANTHINE", "C17XANTHINE", "C37XANTHINE", "C137XANTHINE"]
     likelihood = []
-    for spectrumtrue in spectrafalse:
+    for spectrumtrue in spectra:
         #errorlist = {}
-        errorlist = []
+        #errorlist = []
-        for correctionvalue in range(8, 12): #also tested 8.4, 8.37, 11
+        similaritybycorrection = []
+        correctionvalues = range(8, 16) #also tested 8.4, 8.37, 11, 9.4 (for CNMR3)
+        for correctionvalue in correctionvalues: 
             spectrumrefcorrected = correction(spectrumref, correctionvalue) #CCAFFEINE 11 (klappt hier sehr gut) CCAFFEINE2 12 CPARAXANTHINE 10 CNMR1 9, 10 o 11 (sehr gut) CNMR2 10 o 11
-            error = 0
+            #Likelihood by number of higher similarity than all others.
+            '''error = 0
             total = 0
-            for spectrumfalse in spectrafalse:
+            for spectrumfalse in spectra:
                 positive = 0
                 negative = 0
                 bad_binwidth = []
-                for i in np.arange(0.1, 3.9, 0.1): #Irgendwie hat 0.99999999999 einen Fehler den 1.0 nicht hat 
+                for i in np.arange(0.1, 2.6, 0.1): #successfull at max 3.9, but max 1.7 is lowest where nmr3 correctly classified, 1.6 increases likelihood of 3,7 over 1,7 even with 8.4 correction
                     truesimilarity = similarity_nmr(spectrumtrue, spectrumrefcorrected, i)
                     falsesimilarity = similarity_nmr(spectrumfalse, spectrumrefcorrected, i)
                     #print(truesimilarity)
@@ -324,8 +328,19 @@ def main():
                 error += negative
             #errorlist[correctionvalue] = error
             errorlist.append(error)
-        #print(min(range(len(errorlist)), key=errorlist.__getitem__))
+        print(min(range(len(errorlist)), key=errorlist.__getitem__))
-        likelihood.append(round((total - min(errorlist))/total, 2))
+        likelihood.append(round((total - min(errorlist))/total, 2))'''
+        #Likelihood by mean similarity
+        #This method demonstrates the same problems as the other likelihood method
+            similaritylist = []
+            for i in np.arange(0.1, 3.9, 0.1):
+                similaritylist.append(similarity_nmr(spectrumtrue, spectrumrefcorrected, i))
+            similaritymean = sum(similaritylist) / len(similaritylist)
+            similaritybycorrection.append(similaritymean)
+        name = spectranames[spectra.index(spectrumtrue)]
+        correctionindex = max(range(len(similaritybycorrection)), key=similaritybycorrection.__getitem__)
+        print(f'{name}: {correctionindex} = {correctionvalues[correctionindex]}')
+        likelihood.append(round(max(similaritybycorrection), 2))
     print(likelihood)
     '''for i in np.arange(0.01, 0.07, 0.01):
         print(f'Increment i: {i}')