Cache graph data on creation

readme.md

@@ -159,8 +159,10 @@ Structure:
 Graph and Data files are serialized binaries of a Java object containing the weigthed bipartite graph representation of a
 Sample Plate, along with the necessary metadata for matching and results output. Making them requires a Cell Sample file 
 (to construct a list of correct sequence pairs for checking the accuracy of BiGpairSEQ simulations) and a 
-Sample Plate file (to construct the associated occupancy graph). These files can be several gigabytes in size.
+Sample Plate file (to construct the associated occupancy graph).
-Writing them to a file lets us generate a graph and its metadata once, then use it for multiple different BiGpairSEQ simulations.
+
+These files can be several gigabytes in size. Writing them to a file lets us generate a graph and its metadata once,
+then use it for multiple different BiGpairSEQ simulations.
 
 Options for creating a Graph and Data file:
 * The Cell Sample file to use
@@ -172,7 +174,11 @@ portable data format may be implemented in the future. The tricky part is encodi
 ---
 
 #### Matching Results Files 
-Matching results files consist of the results of a BiGpairSEQ matching simulation.
+Matching results files consist of the results of a BiGpairSEQ matching simulation. Making them requires a Graph and 
+Data file. To save file I/O time, the data from the most recent Graph and Data file read or generated is cached
+by the simulator. Subsequent BiGpairSEQ simulations run with the same input filename will use the cached version
+rather than reading in again from disk.
+
 Files are in CSV format. Rows are sequence pairings with extra relevant data. Columns are pairing-specific details.
 Metadata about the matching simulation is included as comments. Comments are preceded by `#`.
 
@@ -239,8 +245,9 @@ slightly less time than the simulation itself. Real elapsed time from start to f
 ## TODO
 
 * ~~Try invoking GC at end of workloads to reduce paging to disk~~ DONE
-* ~~Hold graph data in memory until another graph is read-in?~~ ABANDONED
+* Hold graph data in memory until another graph is read-in? ~~ABANDONED~~ ~~UNABANDONED~~ DONE
-  * *No, this won't work, because BiGpairSEQ simulations alter the underlying graph based on filtering constraints. Changes would cascade with multiple experiments.*
+  * ~~*No, this won't work, because BiGpairSEQ simulations alter the underlying graph based on filtering constraints. Changes would cascade with multiple experiments.*~~
+  * Might have figured out a way to do it, by taking edges out and then putting them back into the graph. This may actually be possible. If so, awesome.
 * See if there's a reasonable way to reformat Sample Plate files so that wells are columns instead of rows. 
   * ~~Problem is variable number of cells in a well~~
   * ~~Apache Commons CSV library writes entries a row at a time~~ 

src/main/java/BiGpairSEQ.java

@@ -1,14 +1,42 @@
-//main class. Only job is to choose which interface to use.
+//main class. Only job is to choose which interface to use, and hold graph data in memory
 public class BiGpairSEQ {
 
-    private static void main(String[] args) {
+    private static GraphWithMapData graphInMemory = null;
+    private static String graphFilename = null;
+
+    public static void main(String[] args) {
         if (args.length == 0) {
             InteractiveInterface.startInteractive();
         }
         else {
-            //This will be uncommented when command line arguments are fixed.
+            //This will be uncommented when command line arguments are re-implemented.
             //CommandLineInterface.startCLI(args);
             System.out.println("Command line arguments are still being re-implemented.");
         }
     }
+
+    public static GraphWithMapData getGraph() {
+        return graphInMemory;
+    }
+
+    public static void setGraph(GraphWithMapData g) {
+        if (graphInMemory != null) {
+            clearGraph();
+        }
+        graphInMemory = g;
+    }
+
+    public static void clearGraph() {
+        graphInMemory = null;
+        System.gc();
+    }
+
+    public static String getGraphFilename() {
+        return graphFilename;
+    }
+
+    public static void setGraphFilename(String filename) {
+        graphFilename = filename;
+    }
+
 }

src/main/java/CommandLineInterface.java

@@ -1,6 +1,45 @@
 import org.apache.commons.cli.*;
 
-//Class for parsing options passed to program from command line
+/*
+ * Class for parsing options passed to program from command line
+ *
+ * Top-level flags:
+ * cells : to make a cell sample file
+ * plate : to make a sample plate file
+ * graph : to make a graph and data file
+ * match : to do a cdr3 matching (WITH OR WITHOUT MAKING A RESULTS FILE. May just want to print summary for piping.)
+ *
+ * Cell flags:
+ * count : number of cells to generate
+ * diversity factor : factor by which CDR3s are more diverse than CDR1s
+ * output : name of the output file
+ *
+ * Plate flags:
+ * cellfile : name of the cell sample file to use as input
+ * wells : the number of wells on the plate
+ * dist : the statistical distribution to use
+ *      (if exponential) lambda : the lambda value of the exponential distribution
+ *      (if gaussian) stddev : the standard deviation of the gaussian distribution
+ * rand : randomize well populations, take a minimum argument and a maximum argument
+ * populations : number of t cells per well per section (number of arguments determines number of sections)
+ * dropout : plate dropout rate, double from 0.0 to 1.0
+ * output : name of the output file
+ *
+ * Graph flags:
+ * cellfile : name of the cell sample file to use as input
+ * platefile : name of the sample plate file to use as input
+ * output : name of the output file
+ *
+ * Match flags:
+ * graphFile : name of graph and data file to use as input
+ * min : minimum number of overlap wells to attempt a matching
+ * max : the maximum number of overlap wells to attempt a matching
+ * maxdiff : (optional) the maximum difference in occupancy to attempt a matching
+ * minpercent : (optional) the minimum percent overlap to attempt a matching.
+ * writefile : (optional) the filename to write results to
+ * output : the values to print to System.out for piping
+ *
+ */
 public class CommandLineInterface {
 
     public static void startCLI(String[] args) {
@@ -20,7 +59,7 @@ public class CommandLineInterface {
                 .longOpt("make-plates")
                 .desc("Makes a sample plate file")
                 .build();
-        Option makeGraph = Option.builder("graoh")
+        Option makeGraph = Option.builder("graph")
                 .longOpt("make-graph")
                 .desc("Makes a graph and data file")
                 .build();

src/main/java/GraphDataObjectReader.java

@@ -13,6 +13,8 @@ public class GraphDataObjectReader {
             BufferedInputStream fileIn = new BufferedInputStream(new FileInputStream(filename));
                 ObjectInputStream in = new ObjectInputStream(fileIn))
         {
+            System.out.println("Reading graph data from file. This may take some time");
+            System.out.println("File I/O time is not included in results");
             data = (GraphWithMapData) in.readObject();
         } catch (FileNotFoundException | ClassNotFoundException ex) {
             ex.printStackTrace();

src/main/java/GraphDataObjectWriter.java

@@ -18,8 +18,11 @@ public class GraphDataObjectWriter {
 
     public void writeDataToFile() {
         try (BufferedOutputStream bufferedOut = new BufferedOutputStream(new FileOutputStream(filename));
+
              ObjectOutputStream out = new ObjectOutputStream(bufferedOut);
         ){
+            System.out.println("Writing graph and occupancy data to file. This may take some time.");
+            System.out.println("File I/O time is not included in results.");
             out.writeObject(data);
         } catch (IOException ex) {
             ex.printStackTrace();

src/main/java/GraphModificationFunctions.java

@@ -0,0 +1,90 @@
+import org.jgrapht.graph.DefaultWeightedEdge;
+import org.jgrapht.graph.SimpleWeightedGraph;
+
+import java.util.ArrayList;
+import java.util.List;
+import java.util.Map;
+import java.util.Set;
+
+public abstract class GraphModificationFunctions {
+
+    //remove over- and under-weight edges
+    public static List<Integer[]> filterByOverlapThresholds(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
+                                                            int low, int high) {
+        List<Integer[]> removedEdges = new ArrayList<>();
+        for(DefaultWeightedEdge e: graph.edgeSet()){
+            if ((graph.getEdgeWeight(e) > high) || (graph.getEdgeWeight(e) < low)){
+                Integer source = graph.getEdgeSource(e);
+                Integer target = graph.getEdgeTarget(e);
+                Integer weight = (int) graph.getEdgeWeight(e);
+                Integer[] edge = {source, target, weight};
+                removedEdges.add(edge);
+            }
+        }
+        for (Integer[] edge : removedEdges) {
+            graph.removeEdge(edge[0], edge[1]);
+        }
+        return removedEdges;
+    }
+
+    //Remove edges for pairs with large occupancy discrepancy
+    public static List<Integer[]> filterByRelativeOccupancy(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
+                                                  Map<Integer, Integer> alphaWellCounts,
+                                                  Map<Integer, Integer> betaWellCounts,
+                                                  Map<Integer, Integer> plateVtoAMap,
+                                                  Map<Integer, Integer> plateVtoBMap,
+                                                  Integer maxOccupancyDifference) {
+        List<Integer[]> removedEdges = new ArrayList<>();
+        for (DefaultWeightedEdge e : graph.edgeSet()) {
+            Integer alphaOcc = alphaWellCounts.get(plateVtoAMap.get(graph.getEdgeSource(e)));
+            Integer betaOcc = betaWellCounts.get(plateVtoBMap.get(graph.getEdgeTarget(e)));
+            if (Math.abs(alphaOcc - betaOcc) >= maxOccupancyDifference) {
+                Integer source = graph.getEdgeSource(e);
+                Integer target = graph.getEdgeTarget(e);
+                Integer weight = (int) graph.getEdgeWeight(e);
+                Integer[] edge = {source, target, weight};
+                removedEdges.add(edge);
+            }
+        }
+        for (Integer[] edge : removedEdges) {
+            graph.removeEdge(edge[0], edge[1]);
+        }
+        return removedEdges;
+    }
+
+    //Remove edges for pairs where overlap size is significantly lower than the well occupancy
+    public static List<Integer[]> filterByOverlapPercent(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
+                                                         Map<Integer, Integer> alphaWellCounts,
+                                                         Map<Integer, Integer> betaWellCounts,
+                                                         Map<Integer, Integer> plateVtoAMap,
+                                                         Map<Integer, Integer> plateVtoBMap,
+                                                         Integer minOverlapPercent) {
+        List<Integer[]> removedEdges = new ArrayList<>();
+        for (DefaultWeightedEdge e : graph.edgeSet()) {
+            Integer alphaOcc = alphaWellCounts.get(plateVtoAMap.get(graph.getEdgeSource(e)));
+            Integer betaOcc = betaWellCounts.get(plateVtoBMap.get(graph.getEdgeTarget(e)));
+            double weight = graph.getEdgeWeight(e);
+            double min = minOverlapPercent / 100.0;
+            if ((weight / alphaOcc < min) || (weight / betaOcc < min)) {
+                Integer source = graph.getEdgeSource(e);
+                Integer target = graph.getEdgeTarget(e);
+                Integer intWeight = (int) graph.getEdgeWeight(e);
+                Integer[] edge = {source, target, intWeight};
+                removedEdges.add(edge);
+            }
+        }
+        for (Integer[] edge : removedEdges) {
+            graph.removeEdge(edge[0], edge[1]);
+        }
+        return removedEdges;
+    }
+
+    public static void addRemovedEdges(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
+                                       List<Integer[]> removedEdges) {
+        for (Integer[] edge : removedEdges) {
+            DefaultWeightedEdge e = graph.addEdge(edge[0], edge[1]);
+            graph.setEdgeWeight(e, (double) edge[2]);
+        }
+    }
+
+}

src/main/java/InteractiveInterface.java

@@ -31,10 +31,10 @@ public class InteractiveInterface {
                 try {
                     input = sc.nextInt();
                     switch (input) {
-                        case 1 -> makeCellsInteractive();
+                        case 1 -> makeCells();
-                        case 2 -> makePlateInteractive();
+                        case 2 -> makePlate();
-                        case 3 -> makeCDR3GraphInteractive();
+                        case 3 -> makeCDR3Graph();
-                        case 4 -> matchCDR3sInteractive();
+                        case 4 -> matchCDR3s();
                         //case 6 -> matchCellsCDR1();
                         case 9 -> acknowledge();
                         case 0 -> quit = true;
@@ -48,7 +48,7 @@ public class InteractiveInterface {
             sc.close();
     }
 
-    private static void makeCellsInteractive() {
+    private static void makeCells() {
         String filename = null;
         Integer numCells = 0;
         Integer cdr1Freq = 1;
@@ -75,11 +75,10 @@ public class InteractiveInterface {
         assert filename != null;
         CellFileWriter writer = new CellFileWriter(filename, sample);
         writer.writeCellsToFile();
-        System.gc();
     }
 
     //Output a CSV of sample plate
-    private static void makePlateInteractive() {
+    private static void makePlate() {
         String cellFile = null;
         String filename = null;
         Double stdDev = 0.0;
@@ -187,12 +186,11 @@ public class InteractiveInterface {
             System.out.println("Writing Sample Plate to file");
             writer.writePlateFile();
             System.out.println("Sample Plate written to file: " + filename);
-            System.gc();
         }
     }
 
     //Output serialized binary of GraphAndMapData object
-    private static void makeCDR3GraphInteractive() {
+    private static void makeCDR3Graph() {
         String filename = null;
         String cellFile = null;
         String plateFile = null;
@@ -232,18 +230,18 @@ public class InteractiveInterface {
             GraphWithMapData data = Simulator.makeGraph(cells, plate, true);
             assert filename != null;
             GraphDataObjectWriter dataWriter = new GraphDataObjectWriter(filename, data);
-            System.out.println("Writing graph and occupancy data to file. This may take some time.");
-            System.out.println("File I/O time is not included in results.");
             dataWriter.writeDataToFile();
             System.out.println("Graph and Data file written to: " + filename);
-            System.gc();
+            BiGpairSEQ.setGraph(data);
+            BiGpairSEQ.setGraphFilename(filename);
+            System.out.println("Graph and Data file " + filename + " cached.");
         }
     }
 
     //Simulate matching and output CSV file of results
-    private static void matchCDR3sInteractive() throws IOException {
+    private static void matchCDR3s() throws IOException {
         String filename = null;
-        String dataFilename = null;
+        String graphFilename = null;
         Integer lowThreshold = 0;
         Integer highThreshold = Integer.MAX_VALUE;
         Integer maxOccupancyDiff = Integer.MAX_VALUE;
@@ -251,7 +249,7 @@ public class InteractiveInterface {
         try {
             System.out.println("\nBiGpairSEQ simulation requires an occupancy data and overlap graph file");
             System.out.println("Please enter name of an existing graph and occupancy data file: ");
-            dataFilename = sc.next();
+            graphFilename = sc.next();
             System.out.println("The matching results will be written to a file.");
             System.out.print("Please enter a name for the output file: ");
             filename = sc.next();
@@ -274,16 +272,23 @@ public class InteractiveInterface {
             System.out.println(ex);
             sc.next();
         }
+        assert graphFilename != null;
+        //check if this is the same graph we already have in memory.
+        GraphWithMapData data;
+        if(!(graphFilename.equals(BiGpairSEQ.getGraphFilename())) || BiGpairSEQ.getGraph() == null) {
+            BiGpairSEQ.clearGraph();
             //read object data from file
-        System.out.println("Reading graph data from file. This may take some time");
+            GraphDataObjectReader dataReader = new GraphDataObjectReader(graphFilename);
-        System.out.println("File I/O time is not included in results");
+            data = dataReader.getData();
-        assert dataFilename != null;
+            //set new graph in memory and new filename
-        GraphDataObjectReader dataReader = new GraphDataObjectReader(dataFilename);
+            BiGpairSEQ.setGraph(data);
-        GraphWithMapData data = dataReader.getData();
+            BiGpairSEQ.setGraphFilename(graphFilename);
-        //set source file name
+        }
-        data.setSourceFilename(dataFilename);
+        else {
+            data = BiGpairSEQ.getGraph();
+        }
         //simulate matching
-        MatchingResult results = Simulator.matchCDR3s(data, dataFilename, lowThreshold, highThreshold, maxOccupancyDiff,
+        MatchingResult results = Simulator.matchCDR3s(data, graphFilename, lowThreshold, highThreshold, maxOccupancyDiff,
                 minOverlapPercent, true);
         //write results to file
         assert filename != null;
@@ -291,7 +296,6 @@ public class InteractiveInterface {
         System.out.println("Writing results to file");
         writer.writeResultsToFile();
         System.out.println("Results written to file: " + filename);
-        System.gc();
     }
 
     ///////

src/main/java/META-INF/MANIFEST.MF

@@ -0,0 +1,3 @@
+Manifest-Version: 1.0
+Main-Class: BiGpairSEQ
+

src/main/java/MatchingFileWriter.java

@@ -11,7 +11,6 @@ import java.util.List;
 public class MatchingFileWriter {
 
     private String filename;
-    private String sourceFileName;
     private List<String> comments;
     private List<String> headers;
     private List<List<String>> allResults;
@@ -21,7 +20,6 @@ public class MatchingFileWriter {
             filename = filename + ".csv";
         }
         this.filename = filename;
-        this.sourceFileName = result.getSourceFileName();
         this.comments = result.getComments();
         this.headers = result.getHeaders();
         this.allResults = result.getAllResults();

src/main/java/MatchingResult.java

@@ -4,7 +4,7 @@ import java.util.List;
 import java.util.Map;
 
 public class MatchingResult {
-    private final String sourceFile;
+
     private final Map<String, String> metadata;
     private final List<String> comments;
     private final List<String> headers;
@@ -12,16 +12,15 @@ public class MatchingResult {
     private final Map<Integer, Integer> matchMap;
     private final Duration time;
 
-    public MatchingResult(String sourceFileName, Map<String, String> metadata, List<String> headers,
+    public MatchingResult(Map<String, String> metadata, List<String> headers,
                           List<List<String>> allResults, Map<Integer, Integer>matchMap, Duration time){
-        this.sourceFile = sourceFileName;
         /*
          * POSSIBLE KEYS FOR METADATA MAP ARE:
-         * sample plate filename
+         * sample plate filename *
-         * graph filename
+         * graph filename *
-         * well populations
+         * well populations *
-         * total alphas found
+         * total alphas found *
-         * total betas found
+         * total betas found *
          * high overlap threshold
          * low overlap threshold
          * maximum occupancy difference
@@ -66,7 +65,32 @@ public class MatchingResult {
         return time;
     }
 
-    public String getSourceFileName() {
+    public String getPlateFilename() {
-        return sourceFile;
+        return metadata.get("sample plate filename");
     }
+
+    public String getGraphFilename() {
+        return metadata.get("graph filename");
+    }
+
+    public Integer[] getWellPopulations() {
+        List<Integer> wellPopulations = new ArrayList<>();
+        String popString = metadata.get("well populations");
+        for (String p : popString.split(", ")) {
+            wellPopulations.add(Integer.parseInt(p));
+        }
+        Integer[] popArray = new Integer[wellPopulations.size()];
+        return wellPopulations.toArray(popArray);
+    }
+
+    public Integer getAlphaCount() {
+        return Integer.parseInt(metadata.get("total alpha count"));
+    }
+
+    public Integer getBetaCount() {
+        return Integer.parseInt(metadata.get("total beta count"));
+    }
+
+    //put in the rest of these methods following the same pattern
+
 }

src/main/java/PlateFileWriter.java

@@ -16,7 +16,6 @@ public class PlateFileWriter {
     private Double error;
     private String filename;
     private String sourceFileName;
-    private String[] headers;
     private Integer[] concentrations;
     private boolean isExponential = false;
 
@@ -58,20 +57,21 @@ public class PlateFileWriter {
             }
         }
 
-        //this took forever and I don't use it
+//        //this took forever and I don't use it
-        List<List<String>> rows = new ArrayList<>();
+//        //if I wanted to use it, I'd replace printer.printRecords(wellsAsStrings) with printer.printRecords(rows)
-        List<String> tmp = new ArrayList<>();
+//        List<List<String>> rows = new ArrayList<>();
-        for(int i = 0; i < wellsAsStrings.size(); i++){//List<Integer[]> w: wells){
+//        List<String> tmp = new ArrayList<>();
-            tmp.add("well " + (i+1));
+//        for(int i = 0; i < wellsAsStrings.size(); i++){//List<Integer[]> w: wells){
-        }
+//            tmp.add("well " + (i+1));
-        rows.add(tmp);
+//        }
-        for(int row = 0; row < maxLength; row++){
+//        rows.add(tmp);
-            tmp = new ArrayList<>();
+//        for(int row = 0; row < maxLength; row++){
-            for(List<String> c: wellsAsStrings){
+//            tmp = new ArrayList<>();
-                tmp.add(c.get(row));
+//            for(List<String> c: wellsAsStrings){
-            }
+//                tmp.add(c.get(row));
-            rows.add(tmp);
+//            }
-        }
+//            rows.add(tmp);
+//        }
 
         //get list of well populations
         List<Integer> wellPopulations = Arrays.asList(concentrations);

src/main/java/Simulator.java

@@ -1,3 +1,4 @@
+import org.jgrapht.Graph;
 import org.jgrapht.alg.interfaces.MatchingAlgorithm;
 import org.jgrapht.alg.matching.MaximumWeightBipartiteMatching;
 import org.jgrapht.generate.SimpleWeightedBipartiteGraphMatrixGenerator;
@@ -49,6 +50,7 @@ public class Simulator {
         Instant start = Instant.now();
         int[] alphaIndex = {cdr3AlphaIndex};
         int[] betaIndex = {cdr3BetaIndex};
+
         int numWells = samplePlate.getSize();
 
         if(verbose){System.out.println("Making cell maps");}
@@ -63,15 +65,11 @@ public class Simulator {
         if(verbose){System.out.println("All alphas count: " + alphaCount);}
         int betaCount = allBetas.size();
         if(verbose){System.out.println("All betas count: " + betaCount);}
-
         if(verbose){System.out.println("Well maps made");}
 
-        //Remove saturating-occupancy sequences because they have no signal value.
-        //Remove sequences with total occupancy below minimum pair overlap threshold
         if(verbose){System.out.println("Removing sequences present in all wells.");}
-        //if(verbose){System.out.println("Removing sequences with occupancy below the minimum overlap threshold");}
+        filterByOccupancyThresholds(allAlphas, 1, numWells - 1);
-        filterByOccupancyThreshold(allAlphas, 1, numWells - 1);
+        filterByOccupancyThresholds(allBetas, 1, numWells - 1);
-        filterByOccupancyThreshold(allBetas, 1, numWells - 1);
         if(verbose){System.out.println("Sequences removed");}
         int pairableAlphaCount = allAlphas.size();
         if(verbose){System.out.println("Remaining alphas count: " + pairableAlphaCount);}
@@ -136,6 +134,7 @@ public class Simulator {
         GraphWithMapData output = new GraphWithMapData(graph, numWells, samplePlate.getPopulations(), alphaCount, betaCount,
                 distCellsMapAlphaKey, plateVtoAMap, plateVtoBMap, plateAtoVMap,
                 plateBtoVMap, alphaWellCounts, betaWellCounts, time);
+        //Set source file name in graph to name of sample plate
         output.setSourceFilename(samplePlate.getSourceFileName());
         //return GraphWithMapData object
         return output;
@@ -146,6 +145,8 @@ public class Simulator {
                                             Integer highThreshold, Integer maxOccupancyDifference,
                                             Integer minOverlapPercent, boolean verbose) {
         Instant start = Instant.now();
+        //Integer arrays will contain TO VERTEX, FROM VERTEX, and WEIGHT (which I'll need to cast to double)
+        List<Integer[]> removedEdges = new ArrayList<>();
         int numWells = data.getNumWells();
         Integer alphaCount = data.getAlphaCount();
         Integer betaCount = data.getBetaCount();
@@ -156,24 +157,26 @@ public class Simulator {
         Map<Integer, Integer> betaWellCounts = data.getBetaWellCounts();
         SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph = data.getGraph();
 
-        //remove weights outside given overlap thresholds
+        //remove edges with weights outside given overlap thresholds, add those to removed edge list
         if(verbose){System.out.println("Eliminating edges with weights outside overlap threshold values");}
-        filterByOccupancyThreshold(graph, lowThreshold, highThreshold);
+        removedEdges.addAll(GraphModificationFunctions.filterByOverlapThresholds(graph, lowThreshold, highThreshold));
-        if(verbose){System.out.println("Over- and under-weight edges set to 0.0");}
+        if(verbose){System.out.println("Over- and under-weight edges removed");}
 
-        //Filter by overlap size
+        //remove edges between vertices with too small an overlap size, add those to removed edge list
         if(verbose){System.out.println("Eliminating edges with weights less than " + minOverlapPercent.toString() +
                 " percent of vertex occupancy value.");}
-        filterByOverlapSize(graph, alphaWellCounts, betaWellCounts, plateVtoAMap, plateVtoBMap, minOverlapPercent);
+        removedEdges.addAll(GraphModificationFunctions.filterByOverlapPercent(graph, alphaWellCounts, betaWellCounts,
-        if(verbose){System.out.println("Edges with weights too far below vertex occupancy values set to 0.0");}
+                plateVtoAMap, plateVtoBMap, minOverlapPercent));
+        if(verbose){System.out.println("Edges with weights too far below a vertex occupancy value removed");}
 
         //Filter by relative occupancy
         if(verbose){System.out.println("Eliminating edges between vertices with occupancy difference > "
                 + maxOccupancyDifference);}
-        filterByRelativeOccupancy(graph, alphaWellCounts, betaWellCounts, plateVtoAMap, plateVtoBMap,
+        removedEdges.addAll(GraphModificationFunctions.filterByRelativeOccupancy(graph, alphaWellCounts, betaWellCounts,
-                maxOccupancyDifference);
+                plateVtoAMap, plateVtoBMap, maxOccupancyDifference));
         if(verbose){System.out.println("Edges between vertices of with excessively different occupancy values " +
-                "set to 0.0");}
+                "removed");}
+
         //Find Maximum Weighted Matching
         //using jheaps library class PairingHeap for improved efficiency
         if(verbose){System.out.println("Finding maximum weighted matching");}
@@ -239,18 +242,20 @@ public class Simulator {
 
         //Metadata comments for CSV file
         int min = Math.min(alphaCount, betaCount);
+        //rate of attempted matching
         double attemptRate = (double) (trueCount + falseCount) / min;
         BigDecimal attemptRateTrunc = new BigDecimal(attemptRate, mc);
+        //rate of pairing error
         double pairingErrorRate = (double) falseCount / (trueCount + falseCount);
         BigDecimal pairingErrorRateTrunc = new BigDecimal(pairingErrorRate, mc);
         //get list of well concentrations
-        List<Integer> wellPopulations = Arrays.asList(data.getWellConcentrations());
+        Integer[] wellPopulations = data.getWellConcentrations();
         //make string out of concentrations list
         StringBuilder populationsStringBuilder = new StringBuilder();
-        populationsStringBuilder.append(wellPopulations.remove(0).toString());
+        populationsStringBuilder.append(wellPopulations[0].toString());
-        for(Integer i: wellPopulations){
+        for(int i = 1; i < wellPopulations.length; i++){
             populationsStringBuilder.append(", ");
-            populationsStringBuilder.append(i.toString());
+            populationsStringBuilder.append(wellPopulations[i].toString());
         }
         String wellPopulationsString = populationsStringBuilder.toString();
         //total simulation time
@@ -272,13 +277,19 @@ public class Simulator {
         metadata.put("incorrect pairing count", Integer.toString(falseCount));
         metadata.put("pairing error rate", pairingErrorRateTrunc.toString());
         metadata.put("simulation time", nf.format(time.toSeconds()));
-
+        //create MatchingResult object
-        MatchingResult output = new MatchingResult(data.getSourceFilename(), metadata, header, allResults, matchMap, time);
+        MatchingResult output = new MatchingResult(metadata, header, allResults, matchMap, time);
         if(verbose){
             for(String s: output.getComments()){
                 System.out.println(s);
             }
         }
+
+        //put the removed edges back on the graph
+        System.out.println("Restoring removed edges to graph.");
+        GraphModificationFunctions.addRemovedEdges(graph, removedEdges);
+
+        //return MatchingResult object
         return output;
     }
 
@@ -587,6 +598,18 @@ public class Simulator {
 //        return output;
 //    }
 
+    //Remove sequences based on occupancy
+    public static void filterByOccupancyThresholds(Map<Integer, Integer> wellMap, int low, int high){
+        List<Integer> noise = new ArrayList<>();
+        for(Integer k: wellMap.keySet()){
+            if((wellMap.get(k) > high) || (wellMap.get(k) < low)){
+                noise.add(k);
+            }
+        }
+        for(Integer k: noise) {
+            wellMap.remove(k);
+        }
+    }
 
     //Counts the well occupancy of the row peptides and column peptides into given maps, and
     //fills weights in the given 2D array
@@ -630,62 +653,6 @@ public class Simulator {
         }
     }
 
-
-    private static void filterByOccupancyThreshold(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
-                                                   int low, int high) {
-        for(DefaultWeightedEdge e: graph.edgeSet()){
-            if ((graph.getEdgeWeight(e) > high) || (graph.getEdgeWeight(e) < low)){
-                graph.setEdgeWeight(e, 0.0);
-            }
-        }
-    }
-    private static void filterByOccupancyThreshold(Map<Integer, Integer> wellMap, int low, int high){
-        List<Integer> noise = new ArrayList<>();
-        for(Integer k: wellMap.keySet()){
-            if((wellMap.get(k) > high) || (wellMap.get(k) < low)){
-                noise.add(k);
-            }
-        }
-        for(Integer k: noise) {
-            wellMap.remove(k);
-        }
-    }
-
-    //Remove edges for pairs with large occupancy discrepancy
-    private static void filterByRelativeOccupancy(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
-                                                  Map<Integer, Integer> alphaWellCounts,
-                                                  Map<Integer, Integer> betaWellCounts,
-                                                  Map<Integer, Integer> plateVtoAMap,
-                                                  Map<Integer, Integer> plateVtoBMap,
-                                                  Integer maxOccupancyDifference) {
-        for (DefaultWeightedEdge e : graph.edgeSet()) {
-            Integer alphaOcc = alphaWellCounts.get(plateVtoAMap.get(graph.getEdgeSource(e)));
-            Integer betaOcc = betaWellCounts.get(plateVtoBMap.get(graph.getEdgeTarget(e)));
-            //Adjust this to something cleverer later
-            if (Math.abs(alphaOcc - betaOcc) >= maxOccupancyDifference) {
-                graph.setEdgeWeight(e, 0.0);
-            }
-        }
-    }
-
-    //Remove edges for pairs where overlap size is significantly lower than the well occupancy
-    private static void filterByOverlapSize(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
-                                            Map<Integer, Integer> alphaWellCounts,
-                                            Map<Integer, Integer> betaWellCounts,
-                                            Map<Integer, Integer> plateVtoAMap,
-                                            Map<Integer, Integer> plateVtoBMap,
-                                            Integer minOverlapPercent) {
-        for (DefaultWeightedEdge e : graph.edgeSet()) {
-            Integer alphaOcc = alphaWellCounts.get(plateVtoAMap.get(graph.getEdgeSource(e)));
-            Integer betaOcc = betaWellCounts.get(plateVtoBMap.get(graph.getEdgeTarget(e)));
-            double weight = graph.getEdgeWeight(e);
-            double min = minOverlapPercent / 100.0;
-            if ((weight / alphaOcc < min) || (weight / betaOcc < min)) {
-                graph.setEdgeWeight(e, 0.0);
-            }
-        }
-    }
-
     private static Map<Integer, Integer> makeSequenceToSequenceMap(List<Integer[]> cells, int keySequenceIndex,
                                                                    int valueSequenceIndex){
         Map<Integer, Integer> keySequenceToValueSequenceMap = new HashMap<>();