From e58bb6a1cb10400b82d9731e7c2342746cc01df8 Mon Sep 17 00:00:00 2001 From: Sara <sara.schulte@uni-duesseldorf.de> Date: Sat, 7 Jan 2023 23:03:16 +0100 Subject: [PATCH] Add changes to snakemake pipeline --- Pipeline/scripts/create_vac_seq_from_sol.py | 2 +- Pipeline/scripts/vaccine_ilp_hog.py | 9 +++++---- 2 files changed, 6 insertions(+), 5 deletions(-) diff --git a/Pipeline/scripts/create_vac_seq_from_sol.py b/Pipeline/scripts/create_vac_seq_from_sol.py index 023fd9a..87da9d1 100644 --- a/Pipeline/scripts/create_vac_seq_from_sol.py +++ b/Pipeline/scripts/create_vac_seq_from_sol.py @@ -31,7 +31,7 @@ for edge in hog.edges(): edge_dict[key] = edge # create subgraph from edges chosen by ILP in order to find eulerian path to construct vaccine sequence -sub_hog = nx.DiGraph() +sub_hog = nx.MultiDiGraph() sub_hog.add_edges_from([edge_dict[e] for e in edge_dict if e in edges]) print(sub_hog) print(nx.is_eulerian(sub_hog)) diff --git a/Pipeline/scripts/vaccine_ilp_hog.py b/Pipeline/scripts/vaccine_ilp_hog.py index 6892456..ebb19fc 100644 --- a/Pipeline/scripts/vaccine_ilp_hog.py +++ b/Pipeline/scripts/vaccine_ilp_hog.py @@ -61,8 +61,8 @@ def solve_msks_hog(k, alleles, freq_vector, B_matrix, leaves, pep_count, graph, # create new model m = gp.Model('ivp_on_hog') - # set time limit to 10 hours - m.setParam('TimeLimit', 36000) + # set time limit to 48 hours + m.setParam('TimeLimit', 172800) m.setParam('Seed', 42) # m.setParam('PoolSearchMode', 2) # m.setParam('PoolSolutions', 5) @@ -78,7 +78,7 @@ def solve_msks_hog(k, alleles, freq_vector, B_matrix, leaves, pep_count, graph, log('Create edge variables') x_edges = gp.tupledict() for node_a, node_b in sorted(graph.edges): - x_edges[node_a, node_b] = m.addVar(vtype=GRB.BINARY, + x_edges[node_a, node_b] = m.addVar(vtype=GRB.INTEGER, name=graph.nodes[node_a]['string'] + '_' + graph.nodes[node_b]['string']) # the objective function is to maximize population coverage @@ -130,7 +130,8 @@ def solve_msks_hog(k, alleles, freq_vector, B_matrix, leaves, pep_count, graph, sol = m.getAttr('X', x_edges) for node_a, node_b in x_edges: # print('%s -> %s: %g' % (node_a, node_b, sol[node_a, node_b])) - if sol[node_a, node_b] >= 0.5: + traversals = round(sol[node_a, node_b]) + for i in range(traversals): if node_a in leaves.values(): chosen_pep.append(graph.nodes[node_a]['string']) if node_b in leaves.values(): -- GitLab