update readme

readme.md

@@ -12,7 +12,7 @@ Unlike pairSEQ, which calculates p-values for every TCR alpha/beta overlap and c
 against a null distribution, BiGpairSEQ does not do any statistical calculations
 directly.
 
-BiGpairSEQ creates a [simple bipartite weighted graph](https://en.wikipedia.org/wiki/Bipartite_graph) representing the sample plate.
+BiGpairSEQ creates a [weightd bipartite graph](https://en.wikipedia.org/wiki/Bipartite_graph) representing the sample plate.
 The distinct TCRA and TCRB sequences form the two sets of vertices. Every TCRA/TCRB pair that share a well
 are connected by an edge, with the edge weight set to the number of wells in which both sequences appear.
 (Sequences present in *all* wells are filtered out prior to creating the graph, as there is no signal in their occupancy pattern.)
@@ -37,6 +37,8 @@ in practice.
 
 ### RUNNING THE PROGRAM
 
+[Download the current version of BiGpairSEQ_Sim.](https://gitea.ejsf.synology.me/efischer/BiGpairSEQ/releases)
+
 BiGpairSEQ_Sim is an executable .jar file. Requires Java 11 or higher. [OpenJDK 17](https://jdk.java.net/17/)
 recommended.
 
@@ -67,16 +69,26 @@ Please select an option:
 0) Exit
 ```
 
-### OUTPUT
+### INPUT/OUTPUT
 
 To run the simulation, the program reads and writes 4 kinds of files:
 * Cell Sample files in CSV format
 * Sample Plate files in CSV format
-* Graph and Data files in binary object serialization format
+* Graph/Data files in binary object serialization format
 * Matching Results files in CSV format
 
+These files are often generated in sequence. To save file I/O time, the most recent instance of each of these four
+files either generated or read from disk is cached in program memory. This is especially important for Graph/Data files,
+which can be several gigabytes in size. Since some simulations may require running multiple, 
+differntly-configured BiGpairSEQ matchings on the same graph, keeping the most recent graph cached drastically reduces 
+execution time. 
+
+Subsequent uses of the same data file won't need to be read in again until another file of that type is used or generated.
+The program checks whether it needs to update its cached data by comparing filenames as entered by the user. On 
+encountering a new filename, the program flushes its cache and reads in the new file.
+
 When entering filenames, it is not necessary to include the file extension (.csv or .ser). When reading or
-writing files, the program will automatically add the correct extension to any filename without one.
+writing files, the program will automatically add the correct extension to any filename without one. 
 
 #### Cell Sample Files
 Cell Sample files consist of any number of distinct "T cells." Every cell contains 
@@ -119,7 +131,7 @@ Options when making a Sample Plate file:
     * Standard deviation size 
   * Exponential
     * Lambda value
-      * (Based on the slope of the graph in Figure 4C of the pairSEQ paper, the distribution of the original experiment was exponential with a lambda of approximately 0.6. (Howie, et al. 2015))
+      * *(Based on the slope of the graph in Figure 4C of the pairSEQ paper, the distribution of the original experiment was exponential with a lambda of approximately 0.6. (Howie, et al. 2015))*
 * Total number of wells on the plate
 * Number of sections on plate
 * Number of T cells per well
@@ -127,7 +139,7 @@ Options when making a Sample Plate file:
 * Dropout rate
 
 Files are in CSV format. There are no header labels. Every row represents a well. 
-Every column represents an individual cell, containing four sequences, depicted as an array string:
+Every value represents an individual cell, containing four sequences, depicted as an array string:
 `[CDR3A, CDR3B, CDR1A, CDR1B]`. So a representative cell might look like this: 
 
 `[525902, 791533, -1, 866282]`
@@ -153,14 +165,16 @@ Structure:
 
 ---
 
-#### Graph and Data Files
-Graph and Data files are serialized binaries of a Java object containing the weigthed bipartite graph representation of a
+#### Graph/Data Files
+Graph/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.
-Writing them to a file lets us generate a graph and its metadata once, then use it for multiple different BiGpairSEQ simulations.
+Sample Plate file (to construct the associated occupancy graph).
 
-Options for creating a Graph and Data file:
+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/Data file:
 * The Cell Sample file to use
 * The Sample Plate file to use. (This must have been generated from the selected Cell Sample file.)
 
@@ -170,8 +184,8 @@ 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.
-Files are in CSV format. Rows are sequence pairings with extra relevant data. Columns are pairing-specific details.
+Matching results files consist of the results of a BiGpairSEQ matching simulation. Making them requires a Graph and 
+Data file. Matching results 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 `#`.
 
 Options when running a BiGpairSEQ simulation of CDR3 alpha/beta matching:
@@ -237,8 +251,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
-  * *No, this won't work, because BiGpairSEQ simulations alter the underlying graph based on filtering constraints. Changes would cascade with multiple experiments.*
+* 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.*~~
+  * 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

@@ -0,0 +1,95 @@
+import java.util.Random;
+
+//main class. For choosing interface type and caching file data
+public class BiGpairSEQ {
+
+    private static final Random rand = new Random();
+    private static CellSample cellSampleInMemory = null;
+    private static String cellFilename = null;
+    private static Plate plateInMemory = null;
+    private static String plateFilename = null;
+    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 re-implemented.
+            //CommandLineInterface.startCLI(args);
+            System.out.println("Command line arguments are still being re-implemented.");
+        }
+    }
+
+    public static Random getRand() {
+        return rand;
+    }
+
+    public static CellSample getCellSampleInMemory() {
+        return cellSampleInMemory;
+    }
+
+    public static void setCellSampleInMemory(CellSample cellSampleInMemory) {
+        BiGpairSEQ.cellSampleInMemory = cellSampleInMemory;
+    }
+
+    public static void clearCellSampleInMemory() {
+        cellSampleInMemory = null;
+        System.gc();
+    }
+
+    public static String getCellFilename() {
+        return cellFilename;
+    }
+
+    public static void setCellFilename(String cellFilename) {
+        BiGpairSEQ.cellFilename = cellFilename;
+    }
+
+    public static Plate getPlateInMemory() {
+        return plateInMemory;
+    }
+
+    public static void setPlateInMemory(Plate plateInMemory) {
+        BiGpairSEQ.plateInMemory = plateInMemory;
+    }
+
+    public static void clearPlateInMemory() {
+        plateInMemory = null;
+        System.gc();
+    }
+
+    public static String getPlateFilename() {
+        return plateFilename;
+    }
+
+    public static void setPlateFilename(String plateFilename) {
+        BiGpairSEQ.plateFilename = plateFilename;
+    }
+
+    public static GraphWithMapData getGraphInMemory() {
+        return graphInMemory;
+    }
+
+    public static void setGraphInMemory(GraphWithMapData g) {
+        if (graphInMemory != null) {
+            clearGraphInMemory();
+        }
+        graphInMemory = g;
+    }
+
+    public static void clearGraphInMemory() {
+        graphInMemory = null;
+        System.gc();
+    }
+
+    public static String getGraphFilename() {
+        return graphFilename;
+    }
+
+    public static void setGraphFilename(String filename) {
+        graphFilename = filename;
+    }
+
+}

src/main/java/CellFileReader.java

@@ -13,6 +13,7 @@ public class CellFileReader {
 
     private String filename;
     private List<Integer[]> distinctCells = new ArrayList<>();
+    private Integer cdr1Freq;
 
     public CellFileReader(String filename) {
         if(!filename.matches(".*\\.csv")){
@@ -38,19 +39,37 @@ public class CellFileReader {
                 cell[3] = Integer.valueOf(record.get("Beta CDR1"));
                 distinctCells.add(cell);
             }
+
+
         } catch(IOException ex){
             System.out.println("cell file " + filename + " not found.");
             System.err.println(ex);
         }
+
+        //get CDR1 frequency
+        ArrayList<Integer> cdr1Alphas = new ArrayList<>();
+        for (Integer[] cell : distinctCells) {
+            cdr1Alphas.add(cell[3]);
+        }
+        double count = cdr1Alphas.stream().distinct().count();
+        count = Math.ceil(distinctCells.size() / count);
+        cdr1Freq = (int) count;
+
+    }
+
+    public CellSample getCellSample() {
+        return new CellSample(distinctCells, cdr1Freq);
     }
 
     public String getFilename() { return filename;}
 
-    public List<Integer[]> getCells(){
+    //Refactor everything that uses this to have access to a Cell Sample and get the cells there instead.
+    public List<Integer[]> getListOfDistinctCellsDEPRECATED(){
         return distinctCells;
     }
 
-    public Integer getCellCount() {
+    public Integer getCellCountDEPRECATED() {
+        //Refactor everything that uses this to have access to a Cell Sample and get the count there instead.
         return distinctCells.size();
     }
 }

src/main/java/CellSample.java

@@ -18,7 +18,7 @@ public class CellSample {
         return cdr1Freq;
     }
 
-    public Integer population(){
+    public Integer getCellCount(){
         return cells.size();
     }
 

src/main/java/CommandLineInterface.java

@@ -0,0 +1,328 @@
+import org.apache.commons.cli.*;
+
+/*
+ * 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) {
+        //These command line options are a big mess
+        //Really, I don't think command line tools are expected to work in this many different modes
+        //making cells, making plates, and matching are the sort of thing that UNIX philosophy would say
+        //should be three separate programs.
+        //There might be a way to do it with option parameters?
+
+        //main options set
+        Options mainOptions = new Options();
+        Option makeCells = Option.builder("cells")
+                .longOpt("make-cells")
+                .desc("Makes a file of distinct cells")
+                .build();
+        Option makePlate = Option.builder("plates")
+                .longOpt("make-plates")
+                .desc("Makes a sample plate file")
+                .build();
+        Option makeGraph = Option.builder("graph")
+                .longOpt("make-graph")
+                .desc("Makes a graph and data file")
+                .build();
+        Option matchCDR3 = Option.builder("match")
+                .longOpt("match-cdr3")
+                .desc("Match CDR3s. Requires a cell sample file and any number of plate files.")
+                .build();
+        OptionGroup mainGroup = new OptionGroup();
+        mainGroup.addOption(makeCells);
+        mainGroup.addOption(makePlate);
+        mainGroup.addOption(makeGraph);
+        mainGroup.addOption(matchCDR3);
+        mainGroup.setRequired(true);
+        mainOptions.addOptionGroup(mainGroup);
+
+        //Reuse clones of this for other options groups, rather than making it lots of times
+        Option outputFile = Option.builder("o")
+                .longOpt("output-file")
+                .hasArg()
+                .argName("filename")
+                .desc("Name of output file")
+                .build();
+        mainOptions.addOption(outputFile);
+
+        //Options cellOptions = new Options();
+        Option numCells = Option.builder("nc")
+                .longOpt("num-cells")
+                .desc("The number of distinct cells to generate")
+                .hasArg()
+                .argName("number")
+                .build();
+        mainOptions.addOption(numCells);
+        Option cdr1Freq = Option.builder("d")
+                .longOpt("peptide-diversity-factor")
+                .hasArg()
+                .argName("number")
+                .desc("Number of distinct CDR3s for every CDR1")
+                .build();
+        mainOptions.addOption(cdr1Freq);
+        //Option cellOutput = (Option) outputFile.clone();
+        //cellOutput.setRequired(true);
+        //mainOptions.addOption(cellOutput);
+
+        //Options plateOptions = new Options();
+        Option inputCells = Option.builder("c")
+                .longOpt("cell-file")
+                .hasArg()
+                .argName("file")
+                .desc("The cell sample file used for filling wells")
+                .build();
+        mainOptions.addOption(inputCells);
+        Option numWells = Option.builder("w")
+                .longOpt("num-wells")
+                .hasArg()
+                .argName("number")
+                .desc("The number of wells on each plate")
+                .build();
+        mainOptions.addOption(numWells);
+        Option numPlates = Option.builder("np")
+                .longOpt("num-plates")
+                .hasArg()
+                .argName("number")
+                .desc("The number of plate files to output")
+                .build();
+        mainOptions.addOption(numPlates);
+        //Option plateOutput = (Option) outputFile.clone();
+        //plateOutput.setRequired(true);
+        //plateOutput.setDescription("Prefix for plate output filenames");
+        //mainOptions.addOption(plateOutput);
+        Option plateErr = Option.builder("err")
+                .longOpt("drop-out-rate")
+                .hasArg()
+                .argName("number")
+                .desc("Well drop-out rate. (Probability between 0 and 1)")
+                .build();
+        mainOptions.addOption(plateErr);
+        Option plateConcentrations = Option.builder("t")
+                .longOpt("t-cells-per-well")
+                .hasArgs()
+                .argName("number 1, number 2, ...")
+                .desc("Number of T cells per well for each plate section")
+                .build();
+        mainOptions.addOption(plateConcentrations);
+
+//different distributions, mutually exclusive
+        OptionGroup plateDistributions = new OptionGroup();
+        Option plateExp = Option.builder("exponential")
+                .desc("Sample from distinct cells with exponential frequency distribution")
+                .build();
+        plateDistributions.addOption(plateExp);
+        Option plateGaussian = Option.builder("gaussian")
+                .desc("Sample from distinct cells with gaussain frequency distribution")
+                .build();
+        plateDistributions.addOption(plateGaussian);
+        Option platePoisson = Option.builder("poisson")
+                .desc("Sample from distinct cells with poisson frequency distribution")
+                .build();
+        plateDistributions.addOption(platePoisson);
+        mainOptions.addOptionGroup(plateDistributions);
+
+        Option plateStdDev = Option.builder("stddev")
+                .desc("Standard deviation for gaussian distribution")
+                .hasArg()
+                .argName("number")
+                .build();
+        mainOptions.addOption(plateStdDev);
+
+        Option plateLambda = Option.builder("lambda")
+                .desc("Lambda for exponential distribution")
+                .hasArg()
+                .argName("number")
+                .build();
+        mainOptions.addOption(plateLambda);
+
+
+
+//
+//            String cellFile, String filename, Double stdDev,
+//                    Integer numWells, Integer numSections,
+//                    Integer[] concentrations, Double dropOutRate
+//
+
+        //Options matchOptions = new Options();
+        inputCells.setDescription("The cell sample file to be used for matching.");
+        mainOptions.addOption(inputCells);
+        Option lowThresh = Option.builder("low")
+                .longOpt("low-threshold")
+                .hasArg()
+                .argName("number")
+                .desc("Sets the minimum occupancy overlap to attempt matching")
+                .build();
+        mainOptions.addOption(lowThresh);
+        Option highThresh = Option.builder("high")
+                .longOpt("high-threshold")
+                .hasArg()
+                .argName("number")
+                .desc("Sets the maximum occupancy overlap to attempt matching")
+                .build();
+        mainOptions.addOption(highThresh);
+        Option occDiff = Option.builder("occdiff")
+                .longOpt("occupancy-difference")
+                .hasArg()
+                .argName("Number")
+                .desc("Maximum difference in alpha/beta occupancy to attempt matching")
+                .build();
+        mainOptions.addOption(occDiff);
+        Option overlapPer = Option.builder("ovper")
+                .longOpt("overlap-percent")
+                .hasArg()
+                .argName("Percent")
+                .desc("Minimum overlap percent to attempt matching (0 -100)")
+                .build();
+        mainOptions.addOption(overlapPer);
+        Option inputPlates = Option.builder("p")
+                .longOpt("plate-files")
+                .hasArgs()
+                .desc("Plate files to match")
+                .build();
+        mainOptions.addOption(inputPlates);
+
+
+
+        CommandLineParser parser = new DefaultParser();
+        try {
+            CommandLine line = parser.parse(mainOptions, args);
+            if(line.hasOption("match")){
+                //line = parser.parse(mainOptions, args);
+                //String cellFile = line.getOptionValue("c");
+                String graphFile = line.getOptionValue("g");
+                Integer lowThreshold = Integer.valueOf(line.getOptionValue(lowThresh));
+                Integer highThreshold = Integer.valueOf(line.getOptionValue(highThresh));
+                Integer occupancyDifference = Integer.valueOf(line.getOptionValue(occDiff));
+                Integer overlapPercent = Integer.valueOf(line.getOptionValue(overlapPer));
+                for(String plate: line.getOptionValues("p")) {
+                    matchCDR3s(graphFile, lowThreshold, highThreshold, occupancyDifference, overlapPercent);
+                }
+            }
+            else if(line.hasOption("cells")){
+                //line = parser.parse(mainOptions, args);
+                String filename = line.getOptionValue("o");
+                Integer numDistCells = Integer.valueOf(line.getOptionValue("nc"));
+                Integer freq = Integer.valueOf(line.getOptionValue("d"));
+                makeCells(filename, numDistCells, freq);
+            }
+            else if(line.hasOption("plates")){
+                //line = parser.parse(mainOptions, args);
+                String cellFile = line.getOptionValue("c");
+                String filenamePrefix = line.getOptionValue("o");
+                Integer numWellsOnPlate = Integer.valueOf(line.getOptionValue("w"));
+                Integer numPlatesToMake = Integer.valueOf(line.getOptionValue("np"));
+                String[] concentrationsToUseString = line.getOptionValues("t");
+                Integer numSections = concentrationsToUseString.length;
+
+                Integer[] concentrationsToUse = new Integer[numSections];
+                for(int i = 0; i <numSections; i++){
+                    concentrationsToUse[i] = Integer.valueOf(concentrationsToUseString[i]);
+                }
+                Double dropOutRate = Double.valueOf(line.getOptionValue("err"));
+                if(line.hasOption("exponential")){
+                    Double lambda = Double.valueOf(line.getOptionValue("lambda"));
+                    for(int i = 1; i <= numPlatesToMake; i++){
+                        makePlateExp(cellFile, filenamePrefix + i, lambda, numWellsOnPlate,
+                                concentrationsToUse,dropOutRate);
+                    }
+                }
+                else if(line.hasOption("gaussian")){
+                    Double stdDev = Double.valueOf(line.getOptionValue("std-dev"));
+                    for(int i = 1; i <= numPlatesToMake; i++){
+                        makePlate(cellFile, filenamePrefix + i, stdDev, numWellsOnPlate,
+                                concentrationsToUse,dropOutRate);
+                    }
+
+                }
+                else if(line.hasOption("poisson")){
+                    for(int i = 1; i <= numPlatesToMake; i++){
+                        makePlatePoisson(cellFile, filenamePrefix + i, numWellsOnPlate,
+                                concentrationsToUse,dropOutRate);
+                    }
+                }
+            }
+        }
+        catch (ParseException exp) {
+            System.err.println("Parsing failed.  Reason: " + exp.getMessage());
+        }
+    }
+
+    //for calling from command line
+    public static void makeCells(String filename, Integer numCells, Integer cdr1Freq){
+        CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
+        CellFileWriter writer = new CellFileWriter(filename, sample);
+        writer.writeCellsToFile();
+    }
+
+    public static void makePlateExp(String cellFile, String filename, Double lambda,
+                                    Integer numWells, Integer[] concentrations, Double dropOutRate){
+        CellFileReader cellReader = new CellFileReader(cellFile);
+        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
+        samplePlate.fillWellsExponential(cellReader.getFilename(), cellReader.getListOfDistinctCellsDEPRECATED(), lambda);
+        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
+        writer.writePlateFile();
+    }
+
+    private static void makePlatePoisson(String cellFile, String filename, Integer numWells,
+                                        Integer[] concentrations, Double dropOutRate){
+        CellFileReader cellReader = new CellFileReader(cellFile);
+        Double stdDev = Math.sqrt(cellReader.getCellCountDEPRECATED());
+        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
+        samplePlate.fillWells(cellReader.getFilename(), cellReader.getListOfDistinctCellsDEPRECATED(), stdDev);
+        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
+        writer.writePlateFile();
+    }
+
+    private static void makePlate(String cellFile, String filename, Double stdDev,
+                                 Integer numWells, Integer[] concentrations, Double dropOutRate){
+        CellFileReader cellReader = new CellFileReader(cellFile);
+        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
+        samplePlate.fillWells(cellReader.getFilename(), cellReader.getListOfDistinctCellsDEPRECATED(), stdDev);
+        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
+        writer.writePlateFile();
+    }
+
+    private static void matchCDR3s(String graphFile, Integer lowThreshold, Integer highThreshold,
+                                   Integer occupancyDifference, Integer overlapPercent) {
+
+    }
+}

src/main/java/Equations.java

@@ -4,10 +4,6 @@ import java.math.MathContext;
 
 public abstract class Equations {
 
-    public static int getRandomNumber(int min, int max) {
-        return (int) ((Math.random() * (max - min)) + min);
-    }
-
     //pValue calculation as described in original pairSEQ paper.
     //Included for comparison with original results.
     //Not used by BiGpairSEQ for matching.

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/GraphWithMapData.java

@@ -11,7 +11,7 @@ public class GraphWithMapData implements java.io.Serializable {
     private String sourceFilename;
     private final SimpleWeightedGraph graph;
     private Integer numWells;
-    private Integer[] wellConcentrations;
+    private Integer[] wellPopulations;
     private Integer alphaCount;
     private Integer betaCount;
     private final Map<Integer, Integer> distCellsMapAlphaKey;
@@ -31,7 +31,7 @@ public class GraphWithMapData implements java.io.Serializable {
                             Map<Integer, Integer> betaWellCounts, Duration time) {
         this.graph = graph;
         this.numWells = numWells;
-        this.wellConcentrations = wellConcentrations;
+        this.wellPopulations = wellConcentrations;
         this.alphaCount = alphaCount;
         this.betaCount = betaCount;
         this.distCellsMapAlphaKey = distCellsMapAlphaKey;
@@ -52,8 +52,8 @@ public class GraphWithMapData implements java.io.Serializable {
         return numWells;
     }
 
-    public Integer[] getWellConcentrations() {
-        return wellConcentrations;
+    public Integer[] getWellPopulations() {
+        return wellPopulations;
     }
 
     public Integer getAlphaCount() {

src/main/java/InteractiveInterface.java

@@ -0,0 +1,510 @@
+import java.io.IOException;
+import java.util.*;
+import java.util.regex.Matcher;
+import java.util.regex.Pattern;
+
+//
+public class InteractiveInterface {
+
+    private static final Random rand = BiGpairSEQ.getRand();
+    private static final Scanner sc = new Scanner(System.in);
+    private static int input;
+    private static boolean quit = false;
+
+    public static void startInteractive() {
+
+            while (!quit) {
+                System.out.println();
+                System.out.println("--------BiGPairSEQ SIMULATOR--------");
+                System.out.println("ALPHA/BETA T CELL RECEPTOR MATCHING");
+                System.out.println("  USING WEIGHTED BIPARTITE GRAPHS  ");
+                System.out.println("------------------------------------");
+                System.out.println("Please select an option:");
+                System.out.println("1) Generate a population of distinct cells");
+                System.out.println("2) Generate a sample plate of T cells");
+                System.out.println("3) Generate CDR3 alpha/beta occupancy data and overlap graph");
+                System.out.println("4) Simulate bipartite graph CDR3 alpha/beta matching (BiGpairSEQ)");
+                //Need to re-do the CDR3/CDR1 matching to correspond to new pattern
+                //System.out.println("5) Generate CDR3/CDR1 occupancy graph");
+                //System.out.println("6) Simulate CDR3/CDR1 T cell matching");
+                System.out.println("9) About/Acknowledgments");
+                System.out.println("0) Exit");
+                try {
+                    input = sc.nextInt();
+                    switch (input) {
+                        case 1 -> makeCells();
+                        case 2 -> makePlate();
+                        case 3 -> makeCDR3Graph();
+                        case 4 -> matchCDR3s();
+                        //case 6 -> matchCellsCDR1();
+                        case 9 -> acknowledge();
+                        case 0 -> quit = true;
+                        default -> throw new InputMismatchException("Invalid input.");
+                    }
+                } catch (InputMismatchException | IOException ex) {
+                    System.out.println(ex);
+                    sc.next();
+                }
+            }
+            sc.close();
+    }
+
+    private static void makeCells() {
+        String filename = null;
+        Integer numCells = 0;
+        Integer cdr1Freq = 1;
+        try {
+            System.out.println("\nSimulated T-Cells consist of integer values representing:\n" +
+                    "* a pair of alpha and beta CDR3 peptides (unique within simulated population)\n" +
+                    "* a pair of alpha and beta CDR1 peptides (not necessarily unique).");
+            System.out.println("\nThe cells will be written to a CSV file.");
+            System.out.print("Please enter a file name: ");
+            filename = sc.next();
+            System.out.println("\nCDR3 sequences are more diverse than CDR1 sequences.");
+            System.out.println("Please enter the factor by which distinct CDR3s outnumber CDR1s: ");
+            cdr1Freq = sc.nextInt();
+            System.out.print("\nPlease enter the number of T-cells to generate: ");
+            numCells = sc.nextInt();
+            if(numCells <= 0){
+                throw new InputMismatchException("Number of cells must be a positive integer.");
+            }
+        } catch (InputMismatchException ex) {
+            System.out.println(ex);
+            sc.next();
+        }
+        CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
+        assert filename != null;
+        System.out.println("Writing cells to file");
+        CellFileWriter writer = new CellFileWriter(filename, sample);
+        writer.writeCellsToFile();
+        System.out.println("Cell sample written to: " + filename);
+        if(BiGpairSEQ.getCellSampleInMemory() != null) {
+            BiGpairSEQ.clearCellSampleInMemory();
+        }
+        BiGpairSEQ.setCellSampleInMemory(sample);
+        BiGpairSEQ.setCellFilename(filename);
+    }
+
+    //Output a CSV of sample plate
+    private static void makePlate() {
+        String cellFile = null;
+        String filename = null;
+        Double stdDev = 0.0;
+        Integer numWells = 0;
+        Integer numSections;
+        Integer[] populations = {1};
+        Double dropOutRate = 0.0;
+        boolean poisson = false;
+        boolean exponential = false;
+        double lambda = 1.5;
+        try {
+            System.out.println("\nSimulated sample plates consist of:");
+            System.out.println("* a number of wells");
+            System.out.println("  * separated into one or more sections");
+            System.out.println("    * each of which has a set quantity of cells per well");
+            System.out.println("    * selected from a statistical distribution of distinct cells");
+            System.out.println("    * with a set dropout rate for individual sequences within a cell");
+            System.out.println("\nMaking a sample plate requires a population of distinct cells");
+            System.out.print("Please enter name of an existing cell sample file: ");
+            cellFile = sc.next();
+            System.out.println("\nThe sample plate will be written to a CSV file");
+            System.out.print("Please enter a name for the output file: ");
+            filename = sc.next();
+            System.out.println("\nSelect T-cell frequency distribution function");
+            System.out.println("1) Poisson");
+            System.out.println("2) Gaussian");
+            System.out.println("3) Exponential");
+            System.out.println("(Note: approximate distribution in original paper is exponential, lambda = 0.6)");
+            System.out.println("(lambda value approximated from slope of log-log graph in figure 4c)");
+            System.out.println("(Note: wider distributions are more memory intensive to match)");
+            System.out.print("Enter selection value: ");
+            input = sc.nextInt();
+            switch (input) {
+                case 1 -> poisson = true;
+                case 2 -> {
+                    System.out.println("How many distinct T-cells within one standard deviation of peak frequency?");
+                    System.out.println("(Note: wider distributions are more memory intensive to match)");
+                    stdDev = sc.nextDouble();
+                    if (stdDev <= 0.0) {
+                        throw new InputMismatchException("Value must be positive.");
+                    }
+                }
+                case 3 -> {
+                    exponential = true;
+                    System.out.print("Please enter lambda value for exponential distribution: ");
+                    lambda = sc.nextDouble();
+                    if (lambda <= 0.0) {
+                        lambda = 0.6;
+                        System.out.println("Value must be positive. Defaulting to 0.6.");
+                    }
+                }
+                default -> {
+                    System.out.println("Invalid input. Defaulting to exponential.");
+                    exponential = true;
+                }
+            }
+            System.out.print("\nNumber of wells on plate: ");
+            numWells = sc.nextInt();
+            if(numWells < 1){
+                throw new InputMismatchException("No wells on plate");
+            }
+            //choose whether to make T cell population/well random
+            boolean randomWellPopulations;
+            System.out.println("Randomize number of T cells in each well? (y/n)");
+            String ans = sc.next();
+            Pattern pattern = Pattern.compile("(?:yes|y)", Pattern.CASE_INSENSITIVE);
+            Matcher matcher = pattern.matcher(ans);
+            if(matcher.matches()){
+                randomWellPopulations = true;
+            }
+            else{
+                randomWellPopulations = false;
+            }
+            if(randomWellPopulations) { //if T cell population/well is random
+                numSections = numWells;
+                Integer minPop;
+                Integer maxPop;
+                System.out.print("Please enter minimum number of T cells in a well: ");
+                minPop = sc.nextInt();
+                if(minPop < 1) {
+                    throw new InputMismatchException("Minimum well population must be positive");
+                }
+                System.out.println("Please enter maximum number of T cells in a well: ");
+                maxPop = sc.nextInt();
+                if(maxPop < minPop) {
+                    throw new InputMismatchException("Max well population must be greater than min well population");
+                }
+                //maximum should be inclusive, so need to add one to max of randomly generated values
+                populations = rand.ints(minPop, maxPop + 1)
+                        .limit(numSections)
+                        .boxed()
+                        .toArray(Integer[]::new);
+                System.out.print("Populations: ");
+                System.out.println(Arrays.toString(populations));
+            }
+            else{ //if T cell population/well is not random
+                System.out.println("\nThe plate can be evenly sectioned to allow different numbers of T cells per well.");
+                System.out.println("How many sections would you like to make (minimum 1)?");
+                numSections = sc.nextInt();
+                if (numSections < 1) {
+                    throw new InputMismatchException("Too few sections.");
+                } else if (numSections > numWells) {
+                    throw new InputMismatchException("Cannot have more sections than wells.");
+                }
+                int i = 1;
+                populations = new Integer[numSections];
+                while (numSections > 0) {
+                    System.out.print("Enter number of T cells per well in section " + i + ": ");
+                    populations[i - 1] = sc.nextInt();
+                    i++;
+                    numSections--;
+                }
+            }
+            System.out.println("\nErrors in amplification can induce a well dropout rate for sequences");
+            System.out.print("Enter well dropout rate (0.0 to 1.0): ");
+            dropOutRate = sc.nextDouble();
+            if(dropOutRate < 0.0 || dropOutRate > 1.0) {
+                throw new InputMismatchException("The well dropout rate must be in the range [0.0, 1.0]");
+            }
+        }catch(InputMismatchException ex){
+            System.out.println(ex);
+            sc.next();
+        }
+        assert cellFile != null;
+        CellSample cells;
+        if (cellFile.equals(BiGpairSEQ.getCellFilename())){
+            cells = BiGpairSEQ.getCellSampleInMemory();
+        }
+        else {
+            System.out.println("Reading Cell Sample file: " + cellFile);
+            CellFileReader cellReader = new CellFileReader(cellFile);
+            cells = cellReader.getCellSample();
+            BiGpairSEQ.clearCellSampleInMemory();
+            BiGpairSEQ.setCellSampleInMemory(cells);
+            BiGpairSEQ.setCellFilename(cellFile);
+        }
+        assert filename != null;
+        Plate samplePlate;
+        PlateFileWriter writer;
+        if(exponential){
+            samplePlate = new Plate(numWells, dropOutRate, populations);
+            samplePlate.fillWellsExponential(cellFile, cells.getCells(), lambda);
+            writer = new PlateFileWriter(filename, samplePlate);
+        }
+        else {
+            if (poisson) {
+                stdDev = Math.sqrt(cells.getCellCount()); //gaussian with square root of elements approximates poisson
+            }
+            samplePlate = new Plate(numWells, dropOutRate, populations);
+            samplePlate.fillWells(cellFile, cells.getCells(), stdDev);
+            writer = new PlateFileWriter(filename, samplePlate);
+        }
+        System.out.println("Writing Sample Plate to file");
+        writer.writePlateFile();
+        System.out.println("Sample Plate written to file: " + filename);
+        BiGpairSEQ.setPlateInMemory(samplePlate);
+        BiGpairSEQ.setPlateFilename(filename);
+    }
+
+    //Output serialized binary of GraphAndMapData object
+    private static void makeCDR3Graph() {
+        String filename = null;
+        String cellFile = null;
+        String plateFile = null;
+
+        try {
+            String str = "\nGenerating bipartite weighted graph encoding occupancy overlap data ";
+            str = str.concat("\nrequires a cell sample file and a sample plate file.");
+            System.out.println(str);
+            System.out.print("\nPlease enter name of an existing cell sample file: ");
+            cellFile = sc.next();
+            System.out.print("\nPlease enter name of an existing sample plate file: ");
+            plateFile = sc.next();
+            System.out.println("\nThe graph and occupancy data will be written to a serialized binary file.");
+            System.out.print("Please enter a name for the output file: ");
+            filename = sc.next();
+        } catch (InputMismatchException ex) {
+            System.out.println(ex);
+            sc.next();
+        }
+
+        assert cellFile != null;
+        CellSample cellSample;
+        //check if cells are already in memory
+        if(cellFile.equals(BiGpairSEQ.getCellFilename())) {
+            cellSample = BiGpairSEQ.getCellSampleInMemory();
+        }
+        else {
+            BiGpairSEQ.clearCellSampleInMemory();
+            System.out.println("Reading Cell Sample file: " + cellFile);
+            CellFileReader cellReader = new CellFileReader(cellFile);
+            cellSample = cellReader.getCellSample();
+            BiGpairSEQ.setCellSampleInMemory(cellSample);
+            BiGpairSEQ.setCellFilename(cellFile);
+        }
+
+        assert plateFile != null;
+        Plate plate;
+        //check if plate is already in memory
+        if(plateFile.equals(BiGpairSEQ.getPlateFilename())){
+            plate = BiGpairSEQ.getPlateInMemory();
+        }
+        else {
+            BiGpairSEQ.clearPlateInMemory();
+            System.out.println("Reading Sample Plate file: " + plateFile);
+            PlateFileReader plateReader = new PlateFileReader(plateFile);
+            plate = new Plate(plateReader.getFilename(), plateReader.getWells());
+            BiGpairSEQ.setPlateInMemory(plate);
+            BiGpairSEQ.setPlateFilename(plateFile);
+        }
+        if (cellSample.getCells().size() == 0){
+            System.out.println("No cell sample found.");
+            System.out.println("Returning to main menu.");
+        }
+        else if(plate.getWells().size() == 0 || plate.getPopulations().length == 0){
+            System.out.println("No sample plate found.");
+            System.out.println("Returning to main menu.");
+        }
+        else{
+            List<Integer[]> cells = cellSample.getCells();
+            GraphWithMapData data = Simulator.makeGraph(cells, plate, true);
+            assert filename != null;
+            GraphDataObjectWriter dataWriter = new GraphDataObjectWriter(filename, data);
+            dataWriter.writeDataToFile();
+            System.out.println("Graph and Data file written to: " + filename);
+            BiGpairSEQ.setGraphInMemory(data);
+            BiGpairSEQ.setGraphFilename(filename);
+            System.out.println("Graph and Data file " + filename + " cached.");
+        }
+    }
+
+    //Simulate matching and output CSV file of results
+    private static void matchCDR3s() throws IOException {
+        String filename = null;
+        String graphFilename = null;
+        Integer lowThreshold = 0;
+        Integer highThreshold = Integer.MAX_VALUE;
+        Integer maxOccupancyDiff = Integer.MAX_VALUE;
+        Integer minOverlapPercent = 0;
+        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: ");
+            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();
+            System.out.println("\nWhat is the minimum number of CDR3 alpha/beta overlap wells to attempt matching?");
+            lowThreshold = sc.nextInt();
+            if(lowThreshold < 1){
+                lowThreshold = 1;
+                System.out.println("Value for low occupancy overlap threshold must be positive");
+                System.out.println("Value for low occupancy overlap threshold set to 1");
+            }
+            System.out.println("\nWhat is the maximum number of CDR3 alpha/beta overlap wells to attempt matching?");
+            highThreshold = sc.nextInt();
+            if(highThreshold < lowThreshold) {
+                highThreshold = lowThreshold;
+                System.out.println("Value for high occupancy overlap threshold must be >= low overlap threshold");
+                System.out.println("Value for high occupancy overlap threshold set to " + lowThreshold);
+            }
+            System.out.println("What is the minimum percentage of a sequence's wells in alpha/beta overlap to attempt matching? (0 - 100)");
+            minOverlapPercent = sc.nextInt();
+            if (minOverlapPercent < 0 || minOverlapPercent > 100) {
+                System.out.println("Value outside range. Minimum occupancy overlap percentage set to 0");
+            }
+            System.out.println("\nWhat is the maximum difference in alpha/beta occupancy to attempt matching?");
+            maxOccupancyDiff = sc.nextInt();
+            if (maxOccupancyDiff < 0) {
+                maxOccupancyDiff = 0;
+                System.out.println("Maximum allowable difference in alpha/beta occupancy must be nonnegative");
+                System.out.println("Maximum allowable difference in alpha/beta occupancy set to 0");
+            }
+        } catch (InputMismatchException ex) {
+            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.getGraphInMemory() == null) {
+            BiGpairSEQ.clearGraphInMemory();
+            //read object data from file
+            GraphDataObjectReader dataReader = new GraphDataObjectReader(graphFilename);
+            data = dataReader.getData();
+            //set new graph in memory and new filename
+            BiGpairSEQ.setGraphInMemory(data);
+            BiGpairSEQ.setGraphFilename(graphFilename);
+        }
+        else {
+            data = BiGpairSEQ.getGraphInMemory();
+        }
+        //simulate matching
+        MatchingResult results = Simulator.matchCDR3s(data, graphFilename, lowThreshold, highThreshold, maxOccupancyDiff,
+                minOverlapPercent, true);
+        //write results to file
+        assert filename != null;
+        MatchingFileWriter writer = new MatchingFileWriter(filename, results);
+        System.out.println("Writing results to file");
+        writer.writeResultsToFile();
+        System.out.println("Results written to file: " + filename);
+    }
+
+    ///////
+    //Rewrite this to fit new matchCDR3 method with file I/O
+    ///////
+//    public static void matchCellsCDR1(){
+//    /*
+//    The idea here is that we'll get the CDR3 alpha/beta matches first. Then we'll try to match CDR3s to CDR1s by
+//    looking at the top two matches for each CDR3. If CDR3s in the same cell simply swap CDR1s, we assume a correct
+//    match
+//     */
+//        String filename = null;
+//        String preliminaryResultsFilename = null;
+//        String cellFile = null;
+//        String plateFile = null;
+//        Integer lowThresholdCDR3 = 0;
+//        Integer highThresholdCDR3 = Integer.MAX_VALUE;
+//        Integer maxOccupancyDiffCDR3 = 96; //no filtering if max difference is all wells by default
+//        Integer minOverlapPercentCDR3 = 0; //no filtering if min percentage is zero by default
+//        Integer lowThresholdCDR1 = 0;
+//        Integer highThresholdCDR1 = Integer.MAX_VALUE;
+//        boolean outputCDR3Matches = false;
+//        try {
+//            System.out.println("\nSimulated experiment requires a cell sample file and a sample plate file.");
+//            System.out.print("Please enter name of an existing cell sample file: ");
+//            cellFile = sc.next();
+//            System.out.print("Please enter name of an existing sample plate file: ");
+//            plateFile = 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();
+//            System.out.println("What is the minimum number of CDR3 alpha/beta overlap wells to attempt matching?");
+//            lowThresholdCDR3 = sc.nextInt();
+//            if(lowThresholdCDR3 < 1){
+//                throw new InputMismatchException("Minimum value for low threshold is 1");
+//            }
+//            System.out.println("What is the maximum number of CDR3 alpha/beta overlap wells to attempt matching?");
+//            highThresholdCDR3 = sc.nextInt();
+//            System.out.println("What is the maximum difference in CDR3 alpha/beta occupancy to attempt matching?");
+//            maxOccupancyDiffCDR3 = sc.nextInt();
+//            System.out.println("What is the minimum CDR3 overlap percentage to attempt matching? (0 - 100)");
+//            minOverlapPercentCDR3 = sc.nextInt();
+//            if (minOverlapPercentCDR3 < 0 || minOverlapPercentCDR3 > 100) {
+//                throw new InputMismatchException("Value outside range. Minimum percent set to 0");
+//            }
+//            System.out.println("What is the minimum number of CDR3/CDR1 overlap wells to attempt matching?");
+//            lowThresholdCDR1 = sc.nextInt();
+//            if(lowThresholdCDR1 < 1){
+//                throw new InputMismatchException("Minimum value for low threshold is 1");
+//            }
+//            System.out.println("What is the maximum number of CDR3/CDR1 overlap wells to attempt matching?");
+//            highThresholdCDR1 = sc.nextInt();
+//            System.out.println("Matching CDR3s to CDR1s requires first matching CDR3 alpha/betas.");
+//            System.out.println("Output a file for CDR3 alpha/beta match results as well?");
+//            System.out.print("Please enter y/n: ");
+//            String ans = sc.next();
+//            Pattern pattern = Pattern.compile("(?:yes|y)", Pattern.CASE_INSENSITIVE);
+//            Matcher matcher = pattern.matcher(ans);
+//            if(matcher.matches()){
+//                outputCDR3Matches = true;
+//                System.out.println("Please enter filename for CDR3 alpha/beta match results");
+//                preliminaryResultsFilename = sc.next();
+//                System.out.println("CDR3 alpha/beta matches will be output to file");
+//            }
+//            else{
+//                System.out.println("CDR3 alpha/beta matches will not be output to file");
+//            }
+//        } catch (InputMismatchException ex) {
+//            System.out.println(ex);
+//            sc.next();
+//        }
+//        CellFileReader cellReader = new CellFileReader(cellFile);
+//        PlateFileReader plateReader = new PlateFileReader(plateFile);
+//        Plate plate = new Plate(plateReader.getFilename(), plateReader.getWells());
+//        if (cellReader.getCells().size() == 0){
+//            System.out.println("No cell sample found.");
+//            System.out.println("Returning to main menu.");
+//        }
+//        else if(plate.getWells().size() == 0){
+//            System.out.println("No sample plate found.");
+//            System.out.println("Returning to main menu.");
+//
+//        }
+//        else{
+//            if(highThresholdCDR3 >= plate.getSize()){
+//                highThresholdCDR3 = plate.getSize() - 1;
+//            }
+//            if(highThresholdCDR1 >= plate.getSize()){
+//                highThresholdCDR1 = plate.getSize() - 1;
+//            }
+//            List<Integer[]> cells = cellReader.getCells();
+//            MatchingResult preliminaryResults = Simulator.matchCDR3s(cells, plate, lowThresholdCDR3, highThresholdCDR3,
+//                    maxOccupancyDiffCDR3, minOverlapPercentCDR3, true);
+//            MatchingResult[] results = Simulator.matchCDR1s(cells, plate, lowThresholdCDR1,
+//                    highThresholdCDR1, preliminaryResults);
+//            MatchingFileWriter writer = new MatchingFileWriter(filename + "_FirstPass", results[0]);
+//            writer.writeResultsToFile();
+//            writer = new MatchingFileWriter(filename + "_SecondPass", results[1]);
+//            writer.writeResultsToFile();
+//            if(outputCDR3Matches){
+//                writer = new MatchingFileWriter(preliminaryResultsFilename, preliminaryResults);
+//                writer.writeResultsToFile();
+//            }
+//        }
+//    }
+
+    private static void acknowledge(){
+        System.out.println("This program simulates BiGpairSEQ, a graph theory based adaptation");
+        System.out.println("of the pairSEQ algorithm for pairing T cell receptor sequences.");
+        System.out.println();
+        System.out.println("For full documentation, view readme.md file distributed with this code");
+        System.out.println("or visit https://gitea.ejsf.synology.me/efischer/BiGpairSEQ.");
+        System.out.println();
+        System.out.println("pairSEQ citation:");
+        System.out.println("Howie, B., Sherwood, A. M., et. al.");
+        System.out.println("High-throughput pairing of T cell receptor alpha and beta sequences.");
+        System.out.println("Sci. Transl. Med. 7, 301ra131 (2015)");
+        System.out.println();
+        System.out.println("BiGpairSEQ_Sim by Eugene Fischer, 2021-2022");
+    }
+}

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

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

src/main/java/MatchingFileWriter.java

@@ -7,13 +7,10 @@ import java.nio.file.Files;
 import java.nio.file.Path;
 import java.nio.file.StandardOpenOption;
 import java.util.List;
-import java.util.regex.Pattern;
-
 
 public class MatchingFileWriter {
 
     private String filename;
-    private String sourceFileName;
     private List<String> comments;
     private List<String> headers;
     private List<List<String>> allResults;
@@ -23,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

@@ -1,18 +1,41 @@
 import java.time.Duration;
+import java.util.ArrayList;
 import java.util.List;
 import java.util.Map;
 
 public class MatchingResult {
-    private String sourceFile;
-    private List<String> comments;
-    private List<String> headers;
-    private List<List<String>> allResults;
-    private Map<Integer, Integer> matchMap;
-    private Duration time;
 
-    public MatchingResult(String sourceFileName, List<String> comments, List<String> headers, List<List<String>> allResults, Map<Integer, Integer>matchMap, Duration time){
-        this.sourceFile = sourceFileName;
-        this.comments = comments;
+    private final Map<String, String> metadata;
+    private final List<String> comments;
+    private final List<String> headers;
+    private final List<List<String>> allResults;
+    private final Map<Integer, Integer> matchMap;
+    private final Duration time;
+
+    public MatchingResult(Map<String, String> metadata, List<String> headers,
+                          List<List<String>> allResults, Map<Integer, Integer>matchMap, Duration time){
+        /*
+         * POSSIBLE KEYS FOR METADATA MAP ARE:
+         * sample plate filename *
+         * graph filename *
+         * well populations *
+         * total alphas found *
+         * total betas found *
+         * high overlap threshold
+         * low overlap threshold
+         * maximum occupancy difference
+         * minimum overlap percent
+         * pairing attempt rate
+         * correct pairing count
+         * incorrect pairing count
+         * pairing error rate
+         * simulation time
+         */
+        this.metadata = metadata;
+        this.comments = new ArrayList<>();
+        for (String key : metadata.keySet()) {
+            comments.add(key +": " + metadata.get(key));
+        }
         this.headers = headers;
         this.allResults = allResults;
         this.matchMap = matchMap;
@@ -20,6 +43,8 @@ public class MatchingResult {
 
     }
 
+    public Map<String, String> getMetadata() {return metadata;}
+
     public List<String> getComments() {
         return comments;
     }
@@ -40,7 +65,32 @@ public class MatchingResult {
         return time;
     }
 
-    public String getSourceFileName() {
-        return sourceFile;
+    public String getPlateFilename() {
+        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/Plate.java

@@ -1,26 +1,28 @@
-import java.util.*;
+
 
 /*
 TODO: Implement exponential distribution using inversion method - DONE
 TODO: Implement discrete frequency distributions using Vose's Alias Method
  */
 
+import java.util.*;
+
 public class Plate {
     private String sourceFile;
     private List<List<Integer[]>> wells;
-    private Random rand = new Random();
+    private final Random rand = BiGpairSEQ.getRand();
     private int size;
     private double error;
-    private Integer[] concentrations;
+    private Integer[] populations;
     private double stdDev;
     private double lambda;
     boolean exponential = false;
 
 
-    public Plate(int size, double error, Integer[] concentrations) {
+    public Plate(int size, double error, Integer[] populations) {
         this.size = size;
         this.error = error;
-        this.concentrations = concentrations;
+        this.populations = populations;
         wells = new ArrayList<>();
     }
 
@@ -35,9 +37,9 @@ public class Plate {
                 concentrations.add(w.size());
             }
         }
-        this.concentrations = new Integer[concentrations.size()];
-        for (int i = 0; i < this.concentrations.length; i++) {
-            this.concentrations[i] = concentrations.get(i);
+        this.populations = new Integer[concentrations.size()];
+        for (int i = 0; i < this.populations.length; i++) {
+            this.populations[i] = concentrations.get(i);
         }
     }
 
@@ -45,27 +47,19 @@ public class Plate {
         this.lambda = lambda;
         exponential = true;
         sourceFile = sourceFileName;
-        int numSections = concentrations.length;
+        int numSections = populations.length;
         int section = 0;
         double m;
         int n;
-        int test=0;
         while (section < numSections){
             for (int i = 0; i < (size / numSections); i++) {
                 List<Integer[]> well = new ArrayList<>();
-                for (int j = 0; j < concentrations[section]; j++) {
+                for (int j = 0; j < populations[section]; j++) {
                     do {
                         //inverse transform sampling: for random number u in [0,1), x = log(1-u) / (-lambda)
                         m = (Math.log10((1 - rand.nextDouble()))/(-lambda)) * Math.sqrt(cells.size());
                     } while (m >= cells.size() || m < 0);
                     n = (int) Math.floor(m);
-                    //n = Equations.getRandomNumber(0, cells.size());
-                    // was testing generating the cell sample file with exponential dist, then sampling flat here
-                    //that would be more realistic
-                    //But would mess up other things in the simulation with how I've coded it.
-                    if(n > test){
-                        test = n;
-                    }
                     Integer[] cellToAdd = cells.get(n).clone();
                     for(int k = 0; k < cellToAdd.length; k++){
                         if(Math.abs(rand.nextDouble()) < error){//error applied to each seqeunce
@@ -78,20 +72,19 @@ public class Plate {
             }
             section++;
         }
-        System.out.println("Highest index: " +test);
     }
 
     public void fillWells(String sourceFileName, List<Integer[]> cells, double stdDev) {
         this.stdDev = stdDev;
         sourceFile = sourceFileName;
-        int numSections = concentrations.length;
+        int numSections = populations.length;
         int section = 0;
         double m;
         int n;
         while (section < numSections){
             for (int i = 0; i < (size / numSections); i++) {
                 List<Integer[]> well = new ArrayList<>();
-                for (int j = 0; j < concentrations[section]; j++) {
+                for (int j = 0; j < populations[section]; j++) {
                     do {
                         m = (rand.nextGaussian() * stdDev) + (cells.size() / 2);
                     } while (m >= cells.size() || m < 0);
@@ -110,8 +103,8 @@ public class Plate {
         }
     }
 
-    public Integer[] getConcentrations(){
-        return concentrations;
+    public Integer[] getPopulations(){
+        return populations;
     }
 
     public int getSize(){

src/main/java/PlateFileWriter.java

@@ -7,7 +7,6 @@ import java.nio.file.Files;
 import java.nio.file.Path;
 import java.nio.file.StandardOpenOption;
 import java.util.*;
-import java.util.regex.Pattern;
 
 public class PlateFileWriter {
     private int size;
@@ -17,8 +16,7 @@ public class PlateFileWriter {
     private Double error;
     private String filename;
     private String sourceFileName;
-    private String[] headers;
-    private List<Integer> concentrations;
+    private Integer[] populations;
     private boolean isExponential = false;
 
     public PlateFileWriter(String filename, Plate plate) {
@@ -37,8 +35,8 @@ public class PlateFileWriter {
         }
         this.error = plate.getError();
         this.wells = plate.getWells();
-        this.concentrations = Arrays.asList(plate.getConcentrations());
-        concentrations.sort(Comparator.reverseOrder());
+        this.populations = plate.getPopulations();
+        Arrays.sort(populations);
     }
 
     public void writePlateFile(){
@@ -59,28 +57,32 @@ public class PlateFileWriter {
             }
         }
 
-        //this took forever
-        List<List<String>> rows = new ArrayList<>();
-        List<String> tmp = new ArrayList<>();
-        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++){
-            tmp = new ArrayList<>();
-            for(List<String> c: wellsAsStrings){
-                tmp.add(c.get(row));
-            }
-            rows.add(tmp);
-        }
-        //build string of well concentrations
-        StringBuilder concen = new StringBuilder();
-        for(Integer i: concentrations){
-            concen.append(i.toString());
-            concen.append(" ");
-        }
-        String concenString = concen.toString();
+//        //this took forever and I don't use it
+//        //if I wanted to use it, I'd replace printer.printRecords(wellsAsStrings) with printer.printRecords(rows)
+//        List<List<String>> rows = new ArrayList<>();
+//        List<String> tmp = new ArrayList<>();
+//        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++){
+//            tmp = new ArrayList<>();
+//            for(List<String> c: wellsAsStrings){
+//                tmp.add(c.get(row));
+//            }
+//            rows.add(tmp);
+//        }
 
+        //make string out of populations array
+        StringBuilder populationsStringBuilder = new StringBuilder();
+        populationsStringBuilder.append(populations[0].toString());
+        for(int i = 1; i < populations.length; i++){
+            populationsStringBuilder.append(", ");
+            populationsStringBuilder.append(populations[i].toString());
+        }
+        String wellPopulationsString = populationsStringBuilder.toString();
+
+        //set CSV format
         CSVFormat plateFileFormat = CSVFormat.Builder.create()
                 .setCommentMarker('#')
                 .build();
@@ -92,7 +94,7 @@ public class PlateFileWriter {
             printer.printComment("Each row represents one well on the plate.");
             printer.printComment("Plate size: " + size);
             printer.printComment("Error rate: " + error);
-            printer.printComment("Concentrations: " + concenString);
+            printer.printComment("Well populations: " + wellPopulationsString);
             if(isExponential){
                 printer.printComment("Lambda: " + lambda);
             }

src/main/java/Simulator.java

@@ -13,6 +13,8 @@ import java.time.Duration;
 import java.util.*;
 import java.util.stream.IntStream;
 
+import static java.lang.Float.*;
+
 //NOTE: "sequence" in method and variable names refers to a peptide sequence from a simulated T cell
 public class Simulator {
     private static final int cdr3AlphaIndex = 0;
@@ -49,6 +51,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 +66,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");}
-        filterByOccupancyThreshold(allAlphas, 1, numWells - 1);
-        filterByOccupancyThreshold(allBetas, 1, numWells - 1);
+        filterByOccupancyThresholds(allAlphas, 1, numWells - 1);
+        filterByOccupancyThresholds(allBetas, 1, numWells - 1);
         if(verbose){System.out.println("Sequences removed");}
         int pairableAlphaCount = allAlphas.size();
         if(verbose){System.out.println("Remaining alphas count: " + pairableAlphaCount);}
@@ -131,10 +130,15 @@ public class Simulator {
 
         Instant stop = Instant.now();
         Duration time = Duration.between(start, stop);
-        //return GraphWithMapData object
-        return new GraphWithMapData(graph, numWells, samplePlate.getConcentrations(), alphaCount, betaCount,
+
+        //create GraphWithMapData object
+        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;
     }
 
     //match CDR3s.
@@ -142,6 +146,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();
@@ -152,24 +158,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);
-        if(verbose){System.out.println("Over- and under-weight edges set to 0.0");}
+        removedEdges.addAll(GraphModificationFunctions.filterByOverlapThresholds(graph, lowThreshold, highThreshold));
+        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);
-        if(verbose){System.out.println("Edges with weights too far below vertex occupancy values set to 0.0");}
+        removedEdges.addAll(GraphModificationFunctions.filterByOverlapPercent(graph, alphaWellCounts, betaWellCounts,
+                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,
-                maxOccupancyDifference);
+        removedEdges.addAll(GraphModificationFunctions.filterByRelativeOccupancy(graph, alphaWellCounts, betaWellCounts,
+                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");}
@@ -233,351 +241,382 @@ public class Simulator {
             allResults.add(result);
         }
 
-        //Metadate comments for CSV file
+        //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> wellConcentrations = Arrays.asList(data.getWellConcentrations());
-        //make string out of concentrations list
-        StringBuilder concentrationStringBuilder = new StringBuilder();
-        for(Integer i: wellConcentrations){
-            concentrationStringBuilder.append(i.toString());
-            concentrationStringBuilder.append(" ");
+        BigDecimal pairingErrorRateTrunc;
+        if(pairingErrorRate == NaN || pairingErrorRate == POSITIVE_INFINITY || pairingErrorRate == NEGATIVE_INFINITY) {
+            pairingErrorRateTrunc = new BigDecimal(-1, mc);
         }
-        String concentrationString = concentrationStringBuilder.toString();
-
-        List<String> comments = new ArrayList<>();
-        comments.add("Source Sample Plate filename: " + data.getSourceFilename());
-        comments.add("Source Graph and Data filename: " + dataFilename);
-        comments.add("T cell counts in sample plate wells: " + concentrationString);
-        comments.add("Total alphas found: " + alphaCount);
-        comments.add("Total betas found: " + betaCount);
-        comments.add("High overlap threshold: " + highThreshold);
-        comments.add("Low overlap threshold: " + lowThreshold);
-        comments.add("Minimum overlap percent: " + minOverlapPercent);
-        comments.add("Maximum occupancy difference: " + maxOccupancyDifference);
-        comments.add("Pairing attempt rate: " + attemptRateTrunc);
-        comments.add("Correct pairings: " + trueCount);
-        comments.add("Incorrect pairings: " + falseCount);
-        comments.add("Pairing error rate: " + pairingErrorRateTrunc);
+        else{
+            pairingErrorRateTrunc = new BigDecimal(pairingErrorRate, mc);
+        }
+        //get list of well populations
+        Integer[] wellPopulations = data.getWellPopulations();
+        //make string out of populations list
+        StringBuilder populationsStringBuilder = new StringBuilder();
+        populationsStringBuilder.append(wellPopulations[0].toString());
+        for(int i = 1; i < wellPopulations.length; i++){
+            populationsStringBuilder.append(", ");
+            populationsStringBuilder.append(wellPopulations[i].toString());
+        }
+        String wellPopulationsString = populationsStringBuilder.toString();
+        //total simulation time
         Duration time = Duration.between(start, stop);
         time = time.plus(data.getTime());
-        comments.add("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
+
+        Map<String, String> metadata = new LinkedHashMap<>();
+        metadata.put("sample plate filename", data.getSourceFilename());
+        metadata.put("graph filename", dataFilename);
+        metadata.put("well populations", wellPopulationsString);
+        metadata.put("total alphas found", alphaCount.toString());
+        metadata.put("total betas found", betaCount.toString());
+        metadata.put("high overlap threshold", highThreshold.toString());
+        metadata.put("low overlap threshold", lowThreshold.toString());
+        metadata.put("minimum overlap percent", minOverlapPercent.toString());
+        metadata.put("maximum occupancy difference", maxOccupancyDifference.toString());
+        metadata.put("pairing attempt rate", attemptRateTrunc.toString());
+        metadata.put("correct pairing count", Integer.toString(trueCount));
+        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(metadata, header, allResults, matchMap, time);
         if(verbose){
-            for(String s: comments){
+            for(String s: output.getComments()){
                 System.out.println(s);
             }
         }
-        return new MatchingResult(data.getSourceFilename(), comments, header, allResults, matchMap, time);
-    }
 
-    //Simulated matching of CDR1s to CDR3s. Requires MatchingResult from prior run of matchCDR3s.
-    public static MatchingResult[] matchCDR1s(List<Integer[]> distinctCells,
-                                  Plate samplePlate, Integer lowThreshold,
-                                            Integer highThreshold, MatchingResult priorResult){
-        Instant start = Instant.now();
-        Duration previousTime = priorResult.getTime();
-        Map<Integer, Integer> previousMatches = priorResult.getMatchMap();
-        int numWells = samplePlate.getSize();
-        int[] cdr3Indices = {cdr3AlphaIndex, cdr3BetaIndex};
-        int[] cdr1Indices = {cdr1AlphaIndex, cdr1BetaIndex};
+        //put the removed edges back on the graph
+        System.out.println("Restoring removed edges to graph.");
+        GraphModificationFunctions.addRemovedEdges(graph, removedEdges);
 
-        System.out.println("Making previous match maps");
-        Map<Integer, Integer> cdr3AtoBMap = previousMatches;
-        Map<Integer, Integer> cdr3BtoAMap = invertVertexMap(cdr3AtoBMap);
-        System.out.println("Previous match maps made");
-
-        System.out.println("Making cell maps");
-        Map<Integer, Integer> alphaCDR3toCDR1Map = makeSequenceToSequenceMap(distinctCells, cdr3AlphaIndex, cdr1AlphaIndex);
-        Map<Integer, Integer> betaCDR3toCDR1Map = makeSequenceToSequenceMap(distinctCells, cdr3BetaIndex, cdr1BetaIndex);
-        System.out.println("Cell maps made");
-
-        System.out.println("Making well maps");
-        Map<Integer, Integer> allCDR3s = samplePlate.assayWellsSequenceS(cdr3Indices);
-        Map<Integer, Integer> allCDR1s = samplePlate.assayWellsSequenceS(cdr1Indices);
-        int CDR3Count = allCDR3s.size();
-        System.out.println("all CDR3s count: " + CDR3Count);
-        int CDR1Count = allCDR1s.size();
-        System.out.println("all CDR1s count: " + CDR1Count);
-        System.out.println("Well maps made");
-
-        System.out.println("Removing unpaired CDR3s from well maps");
-        List<Integer> unpairedCDR3s = new ArrayList<>();
-        for(Integer i: allCDR3s.keySet()){
-            if(!(cdr3AtoBMap.containsKey(i) || cdr3BtoAMap.containsKey(i))){
-                unpairedCDR3s.add(i);
-            }
-        }
-        for(Integer i: unpairedCDR3s){
-            allCDR3s.remove(i);
-        }
-        System.out.println("Unpaired CDR3s removed.");
-        System.out.println("Remaining CDR3 count: " + allCDR3s.size());
-
-        System.out.println("Removing below-minimum-overlap-threshold and saturating-occupancy CDR1s");
-        filterByOccupancyThreshold(allCDR1s, lowThreshold, numWells - 1);
-        System.out.println("CDR1s removed.");
-        System.out.println("Remaining CDR1 count: " + allCDR1s.size());
-
-        System.out.println("Making vertex maps");
-
-        //For the SimpleWeightedBipartiteGraphMatrixGenerator, all vertices must have
-        // distinct numbers associated with them. Since I'm using a 2D array, that means
-        // distinct indices between the rows and columns. vertexStartValue lets me track where I switch
-        // from numbering rows to columns, so I can assign unique numbers to every vertex, and then
-        // subtract the vertexStartValue from CDR1s to use their vertex labels as array indices
-        Integer vertexStartValue = 0;
-        //keys are sequential integer vertices, values are CDR3s
-        Map<Integer, Integer> plateVtoCDR3Map = makeVertexToSequenceMap(allCDR3s, vertexStartValue);
-        //New start value for vertex to CDR1 map should be one more than final vertex value in CDR3 map
-        vertexStartValue += plateVtoCDR3Map.size();
-        //keys are sequential integers vertices, values are CDR1s
-        Map<Integer, Integer> plateVtoCDR1Map = makeVertexToSequenceMap(allCDR1s, vertexStartValue);
-        //keys are CDR3s, values are sequential integer vertices from previous map
-        Map<Integer, Integer> plateCDR3toVMap = invertVertexMap(plateVtoCDR3Map);
-        //keys are CDR1s, values are sequential integer vertices from previous map
-        Map<Integer, Integer> plateCDR1toVMap = invertVertexMap(plateVtoCDR1Map);
-        System.out.println("Vertex maps made");
-
-        System.out.println("Creating adjacency matrix");
-        //Count how many wells each CDR3 appears in
-        Map<Integer, Integer> cdr3WellCounts = new HashMap<>();
-        //count how many wells each CDR1 appears in
-        Map<Integer, Integer> cdr1WellCounts = new HashMap<>();
-        //add edges, where weights are number of wells the peptides share in common.
-        //If this is too slow, can make a 2d array and use the SimpleWeightedGraphMatrixGenerator class
-        Map<Integer, Integer> wellNCDR3s = null;
-        Map<Integer, Integer> wellNCDR1s = null;
-        double[][] weights = new double[plateVtoCDR3Map.size()][plateVtoCDR1Map.size()];
-        countSequencesAndFillMatrix(samplePlate, allCDR3s, allCDR1s, plateCDR3toVMap, plateCDR1toVMap,
-                cdr3Indices, cdr1Indices, cdr3WellCounts, cdr1WellCounts, weights);
-        System.out.println("Matrix created");
-
-        System.out.println("Creating graph");
-        SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph =
-                new SimpleWeightedGraph<>(DefaultWeightedEdge.class);
-
-        SimpleWeightedBipartiteGraphMatrixGenerator graphGenerator = new SimpleWeightedBipartiteGraphMatrixGenerator();
-        List<Integer> cdr3Vertices = new ArrayList<>(plateVtoCDR3Map.keySet()); //This will work because LinkedHashMap preserves order of entry
-        graphGenerator.first(cdr3Vertices);
-        List<Integer> cdr1Vertices = new ArrayList<>(plateVtoCDR1Map.keySet());
-        graphGenerator.second(cdr1Vertices); //This will work because LinkedHashMap preserves order of entry
-        graphGenerator.weights(weights);
-        graphGenerator.generateGraph(graph);
-        System.out.println("Graph created");
-
-        System.out.println("Removing edges outside of weight thresholds");
-        filterByOccupancyThreshold(graph, lowThreshold, highThreshold);
-        System.out.println("Over- and under-weight edges set to 0.0");
-
-        System.out.println("Finding first maximum weighted matching");
-        MaximumWeightBipartiteMatching firstMaxWeightMatching =
-                new MaximumWeightBipartiteMatching(graph, plateVtoCDR3Map.keySet(), plateVtoCDR1Map.keySet());
-        MatchingAlgorithm.Matching<String, DefaultWeightedEdge> graphMatching = firstMaxWeightMatching.getMatching();
-        System.out.println("First maximum weighted matching found");
-
-
-        //first processing run
-        Map<Integer, Integer> firstMatchCDR3toCDR1Map = new HashMap<>();
-        Iterator<DefaultWeightedEdge> weightIter = graphMatching.iterator();
-        DefaultWeightedEdge e;
-        while(weightIter.hasNext()){
-            e = weightIter.next();
-//            if(graph.getEdgeWeight(e) < lowThreshold || graph.getEdgeWeight(e) > highThreshold) {
-//                continue;
-//            }
-            Integer source = graph.getEdgeSource(e);
-            Integer target = graph.getEdgeTarget(e);
-            firstMatchCDR3toCDR1Map.put(plateVtoCDR3Map.get(source), plateVtoCDR1Map.get(target));
-        }
-        System.out.println("First pass matches: " + firstMatchCDR3toCDR1Map.size());
-
-        System.out.println("Removing edges from first maximum weighted matching");
-        //zero out the edge weights in the matching
-        weightIter = graphMatching.iterator();
-        while(weightIter.hasNext()){
-            graph.removeEdge(weightIter.next());
-        }
-        System.out.println("Edges removed");
-
-        //Generate a new matching
-        System.out.println("Finding second maximum weighted matching");
-        MaximumWeightBipartiteMatching secondMaxWeightMatching =
-                new MaximumWeightBipartiteMatching(graph, plateVtoCDR3Map.keySet(), plateVtoCDR1Map.keySet());
-        graphMatching = secondMaxWeightMatching.getMatching();
-        System.out.println("Second maximum weighted matching found");
-
-
-        //second processing run
-        Map<Integer, Integer> secondMatchCDR3toCDR1Map = new HashMap<>();
-        weightIter = graphMatching.iterator();
-        while(weightIter.hasNext()){
-            e = weightIter.next();
-//            if(graph.getEdgeWeight(e) < lowThreshold || graph.getEdgeWeight(e) > highThreshold) {
-//                continue;
-//            }
-            Integer source = graph.getEdgeSource(e);
-//            if(!(CDR3AtoBMap.containsKey(source) || CDR3BtoAMap.containsKey(source))){
-//                continue;
-//            }
-            Integer target = graph.getEdgeTarget(e);
-            secondMatchCDR3toCDR1Map.put(plateVtoCDR3Map.get(source), plateVtoCDR1Map.get(target));
-        }
-        System.out.println("Second pass matches: " + secondMatchCDR3toCDR1Map.size());
-
-        System.out.println("Mapping first pass CDR3 alpha/beta pairs");
-        //get linked map for first matching attempt
-        Map<Integer, Integer> firstMatchesMap = new LinkedHashMap<>();
-        for(Integer alphaCDR3: cdr3AtoBMap.keySet()) {
-            if (!(firstMatchCDR3toCDR1Map.containsKey(alphaCDR3))) {
-                continue;
-            }
-            Integer betaCDR3 = cdr3AtoBMap.get(alphaCDR3);
-            if (!(firstMatchCDR3toCDR1Map.containsKey(betaCDR3))) {
-                continue;
-            }
-            firstMatchesMap.put(alphaCDR3, firstMatchCDR3toCDR1Map.get(alphaCDR3));
-            firstMatchesMap.put(betaCDR3, firstMatchCDR3toCDR1Map.get(betaCDR3));
-        }
-        System.out.println("First pass CDR3 alpha/beta pairs mapped");
-
-        System.out.println("Mapping second pass CDR3 alpha/beta pairs.");
-        System.out.println("Finding CDR3 pairs that swapped CDR1 matches between first pass and second pass.");
-        //Look for matches that simply swapped already-matched alpha and beta CDR3s
-        Map<Integer, Integer> dualMatchesMap = new LinkedHashMap<>();
-        for(Integer alphaCDR3: cdr3AtoBMap.keySet()) {
-            if (!(firstMatchCDR3toCDR1Map.containsKey(alphaCDR3) && secondMatchCDR3toCDR1Map.containsKey(alphaCDR3))) {
-                continue;
-            }
-            Integer betaCDR3 = cdr3AtoBMap.get(alphaCDR3);
-            if (!(firstMatchCDR3toCDR1Map.containsKey(betaCDR3) && secondMatchCDR3toCDR1Map.containsKey(betaCDR3))) {
-                continue;
-            }
-            if(firstMatchCDR3toCDR1Map.get(alphaCDR3).equals(secondMatchCDR3toCDR1Map.get(betaCDR3))){
-                if(firstMatchCDR3toCDR1Map.get(betaCDR3).equals(secondMatchCDR3toCDR1Map.get(alphaCDR3))){
-                    dualMatchesMap.put(alphaCDR3, firstMatchCDR3toCDR1Map.get(alphaCDR3));
-                    dualMatchesMap.put(betaCDR3, firstMatchCDR3toCDR1Map.get(betaCDR3));
-                }
-            }
-        }
-        System.out.println("Second pass mapping made. Dual CDR3/CDR1 pairings found.");
-
-        Instant stop = Instant.now();
-        //results for first map
-        System.out.println("RESULTS FOR FIRST PASS MATCHING");
-        List<List<String>> allResults = new ArrayList<>();
-        Integer trueCount = 0;
-        Iterator iter = firstMatchesMap.keySet().iterator();
-
-        while(iter.hasNext()){
-            Boolean proven = false;
-            List<String> tmp = new ArrayList<>();
-            tmp.add(iter.next().toString());
-            tmp.add(iter.next().toString());
-            tmp.add(firstMatchesMap.get(Integer.valueOf(tmp.get(0))).toString());
-            tmp.add(firstMatchesMap.get(Integer.valueOf(tmp.get(1))).toString());
-            if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(2)))){
-                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(3)))){
-                    proven = true;
-                }
-            }
-            else if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(3)))){
-                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(2)))){
-                    proven = true;
-                }
-            }
-            tmp.add(proven.toString());
-            allResults.add(tmp);
-            if(proven){
-                trueCount++;
-            }
-        }
-
-        List<String> comments = new ArrayList<>();
-        comments.add("Plate size: " + samplePlate.getSize() + " wells");
-        comments.add("Previous pairs found: " + previousMatches.size());
-        comments.add("CDR1 matches attempted: " + allResults.size());
-        double attemptRate = (double) allResults.size() / previousMatches.size();
-        comments.add("Matching attempt rate: " + attemptRate);
-        comments.add("Number of correct matches: " + trueCount);
-        double correctRate = (double) trueCount / allResults.size();
-        comments.add("Correct matching rate: " + correctRate);
-        NumberFormat nf = NumberFormat.getInstance(Locale.US);
-        Duration time = Duration.between(start, stop);
-        time = time.plus(previousTime);
-        comments.add("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
-        for(String s: comments){
-            System.out.println(s);
-        }
-
-
-
-        List<String> headers = new ArrayList<>();
-        headers.add("CDR3 alpha");
-        headers.add("CDR3 beta");
-        headers.add("first matched CDR1");
-        headers.add("second matched CDR1");
-        headers.add("Correct match?");
-
-        MatchingResult firstTest = new MatchingResult(samplePlate.getSourceFileName(),
-                comments, headers, allResults, dualMatchesMap, time);
-
-        //results for dual map
-        System.out.println("RESULTS FOR SECOND PASS MATCHING");
-        allResults = new ArrayList<>();
-        trueCount = 0;
-        iter = dualMatchesMap.keySet().iterator();
-        while(iter.hasNext()){
-            Boolean proven = false;
-            List<String> tmp = new ArrayList<>();
-            tmp.add(iter.next().toString());
-            tmp.add(iter.next().toString());
-            tmp.add(dualMatchesMap.get(Integer.valueOf(tmp.get(0))).toString());
-            tmp.add(dualMatchesMap.get(Integer.valueOf(tmp.get(1))).toString());
-            if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(2)))){
-                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(3)))){
-                    proven = true;
-                }
-            }
-            else if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(3)))){
-                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(2)))){
-                    proven = true;
-                }
-            }
-            tmp.add(proven.toString());
-            allResults.add(tmp);
-            if(proven){
-                trueCount++;
-            }
-        }
-
-        comments = new ArrayList<>();
-        comments.add("Plate size: " + samplePlate.getSize() + " wells");
-        comments.add("Previous pairs found: " + previousMatches.size());
-        comments.add("High overlap threshold: " + highThreshold);
-        comments.add("Low overlap threshold: " + lowThreshold);
-        comments.add("CDR1 matches attempted: " + allResults.size());
-        attemptRate = (double) allResults.size() / previousMatches.size();
-        comments.add("Matching attempt rate: " + attemptRate);
-        comments.add("Number of correct matches: " + trueCount);
-        correctRate = (double) trueCount / allResults.size();
-        comments.add("Correct matching rate: " + correctRate);
-        comments.add("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
-
-        for(String s: comments){
-            System.out.println(s);
-        }
-
-        System.out.println("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
-        MatchingResult dualTest = new MatchingResult(samplePlate.getSourceFileName(), comments, headers,
-                allResults, dualMatchesMap, time);
-        MatchingResult[] output = {firstTest, dualTest};
+        //return MatchingResult object
         return output;
     }
 
+    //Commented out CDR1 matching until it's time to re-implement it
+//    //Simulated matching of CDR1s to CDR3s. Requires MatchingResult from prior run of matchCDR3s.
+//    public static MatchingResult[] matchCDR1s(List<Integer[]> distinctCells,
+//                                  Plate samplePlate, Integer lowThreshold,
+//                                            Integer highThreshold, MatchingResult priorResult){
+//        Instant start = Instant.now();
+//        Duration previousTime = priorResult.getTime();
+//        Map<Integer, Integer> previousMatches = priorResult.getMatchMap();
+//        int numWells = samplePlate.getSize();
+//        int[] cdr3Indices = {cdr3AlphaIndex, cdr3BetaIndex};
+//        int[] cdr1Indices = {cdr1AlphaIndex, cdr1BetaIndex};
+//
+//        System.out.println("Making previous match maps");
+//        Map<Integer, Integer> cdr3AtoBMap = previousMatches;
+//        Map<Integer, Integer> cdr3BtoAMap = invertVertexMap(cdr3AtoBMap);
+//        System.out.println("Previous match maps made");
+//
+//        System.out.println("Making cell maps");
+//        Map<Integer, Integer> alphaCDR3toCDR1Map = makeSequenceToSequenceMap(distinctCells, cdr3AlphaIndex, cdr1AlphaIndex);
+//        Map<Integer, Integer> betaCDR3toCDR1Map = makeSequenceToSequenceMap(distinctCells, cdr3BetaIndex, cdr1BetaIndex);
+//        System.out.println("Cell maps made");
+//
+//        System.out.println("Making well maps");
+//        Map<Integer, Integer> allCDR3s = samplePlate.assayWellsSequenceS(cdr3Indices);
+//        Map<Integer, Integer> allCDR1s = samplePlate.assayWellsSequenceS(cdr1Indices);
+//        int CDR3Count = allCDR3s.size();
+//        System.out.println("all CDR3s count: " + CDR3Count);
+//        int CDR1Count = allCDR1s.size();
+//        System.out.println("all CDR1s count: " + CDR1Count);
+//        System.out.println("Well maps made");
+//
+//        System.out.println("Removing unpaired CDR3s from well maps");
+//        List<Integer> unpairedCDR3s = new ArrayList<>();
+//        for(Integer i: allCDR3s.keySet()){
+//            if(!(cdr3AtoBMap.containsKey(i) || cdr3BtoAMap.containsKey(i))){
+//                unpairedCDR3s.add(i);
+//            }
+//        }
+//        for(Integer i: unpairedCDR3s){
+//            allCDR3s.remove(i);
+//        }
+//        System.out.println("Unpaired CDR3s removed.");
+//        System.out.println("Remaining CDR3 count: " + allCDR3s.size());
+//
+//        System.out.println("Removing below-minimum-overlap-threshold and saturating-occupancy CDR1s");
+//        filterByOccupancyThreshold(allCDR1s, lowThreshold, numWells - 1);
+//        System.out.println("CDR1s removed.");
+//        System.out.println("Remaining CDR1 count: " + allCDR1s.size());
+//
+//        System.out.println("Making vertex maps");
+//
+//        //For the SimpleWeightedBipartiteGraphMatrixGenerator, all vertices must have
+//        // distinct numbers associated with them. Since I'm using a 2D array, that means
+//        // distinct indices between the rows and columns. vertexStartValue lets me track where I switch
+//        // from numbering rows to columns, so I can assign unique numbers to every vertex, and then
+//        // subtract the vertexStartValue from CDR1s to use their vertex labels as array indices
+//        Integer vertexStartValue = 0;
+//        //keys are sequential integer vertices, values are CDR3s
+//        Map<Integer, Integer> plateVtoCDR3Map = makeVertexToSequenceMap(allCDR3s, vertexStartValue);
+//        //New start value for vertex to CDR1 map should be one more than final vertex value in CDR3 map
+//        vertexStartValue += plateVtoCDR3Map.size();
+//        //keys are sequential integers vertices, values are CDR1s
+//        Map<Integer, Integer> plateVtoCDR1Map = makeVertexToSequenceMap(allCDR1s, vertexStartValue);
+//        //keys are CDR3s, values are sequential integer vertices from previous map
+//        Map<Integer, Integer> plateCDR3toVMap = invertVertexMap(plateVtoCDR3Map);
+//        //keys are CDR1s, values are sequential integer vertices from previous map
+//        Map<Integer, Integer> plateCDR1toVMap = invertVertexMap(plateVtoCDR1Map);
+//        System.out.println("Vertex maps made");
+//
+//        System.out.println("Creating adjacency matrix");
+//        //Count how many wells each CDR3 appears in
+//        Map<Integer, Integer> cdr3WellCounts = new HashMap<>();
+//        //count how many wells each CDR1 appears in
+//        Map<Integer, Integer> cdr1WellCounts = new HashMap<>();
+//        //add edges, where weights are number of wells the peptides share in common.
+//        //If this is too slow, can make a 2d array and use the SimpleWeightedGraphMatrixGenerator class
+//        Map<Integer, Integer> wellNCDR3s = null;
+//        Map<Integer, Integer> wellNCDR1s = null;
+//        double[][] weights = new double[plateVtoCDR3Map.size()][plateVtoCDR1Map.size()];
+//        countSequencesAndFillMatrix(samplePlate, allCDR3s, allCDR1s, plateCDR3toVMap, plateCDR1toVMap,
+//                cdr3Indices, cdr1Indices, cdr3WellCounts, cdr1WellCounts, weights);
+//        System.out.println("Matrix created");
+//
+//        System.out.println("Creating graph");
+//        SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph =
+//                new SimpleWeightedGraph<>(DefaultWeightedEdge.class);
+//
+//        SimpleWeightedBipartiteGraphMatrixGenerator graphGenerator = new SimpleWeightedBipartiteGraphMatrixGenerator();
+//        List<Integer> cdr3Vertices = new ArrayList<>(plateVtoCDR3Map.keySet()); //This will work because LinkedHashMap preserves order of entry
+//        graphGenerator.first(cdr3Vertices);
+//        List<Integer> cdr1Vertices = new ArrayList<>(plateVtoCDR1Map.keySet());
+//        graphGenerator.second(cdr1Vertices); //This will work because LinkedHashMap preserves order of entry
+//        graphGenerator.weights(weights);
+//        graphGenerator.generateGraph(graph);
+//        System.out.println("Graph created");
+//
+//        System.out.println("Removing edges outside of weight thresholds");
+//        filterByOccupancyThreshold(graph, lowThreshold, highThreshold);
+//        System.out.println("Over- and under-weight edges set to 0.0");
+//
+//        System.out.println("Finding first maximum weighted matching");
+//        MaximumWeightBipartiteMatching firstMaxWeightMatching =
+//                new MaximumWeightBipartiteMatching(graph, plateVtoCDR3Map.keySet(), plateVtoCDR1Map.keySet());
+//        MatchingAlgorithm.Matching<String, DefaultWeightedEdge> graphMatching = firstMaxWeightMatching.getMatching();
+//        System.out.println("First maximum weighted matching found");
+//
+//
+//        //first processing run
+//        Map<Integer, Integer> firstMatchCDR3toCDR1Map = new HashMap<>();
+//        Iterator<DefaultWeightedEdge> weightIter = graphMatching.iterator();
+//        DefaultWeightedEdge e;
+//        while(weightIter.hasNext()){
+//            e = weightIter.next();
+////            if(graph.getEdgeWeight(e) < lowThreshold || graph.getEdgeWeight(e) > highThreshold) {
+////                continue;
+////            }
+//            Integer source = graph.getEdgeSource(e);
+//            Integer target = graph.getEdgeTarget(e);
+//            firstMatchCDR3toCDR1Map.put(plateVtoCDR3Map.get(source), plateVtoCDR1Map.get(target));
+//        }
+//        System.out.println("First pass matches: " + firstMatchCDR3toCDR1Map.size());
+//
+//        System.out.println("Removing edges from first maximum weighted matching");
+//        //zero out the edge weights in the matching
+//        weightIter = graphMatching.iterator();
+//        while(weightIter.hasNext()){
+//            graph.removeEdge(weightIter.next());
+//        }
+//        System.out.println("Edges removed");
+//
+//        //Generate a new matching
+//        System.out.println("Finding second maximum weighted matching");
+//        MaximumWeightBipartiteMatching secondMaxWeightMatching =
+//                new MaximumWeightBipartiteMatching(graph, plateVtoCDR3Map.keySet(), plateVtoCDR1Map.keySet());
+//        graphMatching = secondMaxWeightMatching.getMatching();
+//        System.out.println("Second maximum weighted matching found");
+//
+//
+//        //second processing run
+//        Map<Integer, Integer> secondMatchCDR3toCDR1Map = new HashMap<>();
+//        weightIter = graphMatching.iterator();
+//        while(weightIter.hasNext()){
+//            e = weightIter.next();
+////            if(graph.getEdgeWeight(e) < lowThreshold || graph.getEdgeWeight(e) > highThreshold) {
+////                continue;
+////            }
+//            Integer source = graph.getEdgeSource(e);
+////            if(!(CDR3AtoBMap.containsKey(source) || CDR3BtoAMap.containsKey(source))){
+////                continue;
+////            }
+//            Integer target = graph.getEdgeTarget(e);
+//            secondMatchCDR3toCDR1Map.put(plateVtoCDR3Map.get(source), plateVtoCDR1Map.get(target));
+//        }
+//        System.out.println("Second pass matches: " + secondMatchCDR3toCDR1Map.size());
+//
+//        System.out.println("Mapping first pass CDR3 alpha/beta pairs");
+//        //get linked map for first matching attempt
+//        Map<Integer, Integer> firstMatchesMap = new LinkedHashMap<>();
+//        for(Integer alphaCDR3: cdr3AtoBMap.keySet()) {
+//            if (!(firstMatchCDR3toCDR1Map.containsKey(alphaCDR3))) {
+//                continue;
+//            }
+//            Integer betaCDR3 = cdr3AtoBMap.get(alphaCDR3);
+//            if (!(firstMatchCDR3toCDR1Map.containsKey(betaCDR3))) {
+//                continue;
+//            }
+//            firstMatchesMap.put(alphaCDR3, firstMatchCDR3toCDR1Map.get(alphaCDR3));
+//            firstMatchesMap.put(betaCDR3, firstMatchCDR3toCDR1Map.get(betaCDR3));
+//        }
+//        System.out.println("First pass CDR3 alpha/beta pairs mapped");
+//
+//        System.out.println("Mapping second pass CDR3 alpha/beta pairs.");
+//        System.out.println("Finding CDR3 pairs that swapped CDR1 matches between first pass and second pass.");
+//        //Look for matches that simply swapped already-matched alpha and beta CDR3s
+//        Map<Integer, Integer> dualMatchesMap = new LinkedHashMap<>();
+//        for(Integer alphaCDR3: cdr3AtoBMap.keySet()) {
+//            if (!(firstMatchCDR3toCDR1Map.containsKey(alphaCDR3) && secondMatchCDR3toCDR1Map.containsKey(alphaCDR3))) {
+//                continue;
+//            }
+//            Integer betaCDR3 = cdr3AtoBMap.get(alphaCDR3);
+//            if (!(firstMatchCDR3toCDR1Map.containsKey(betaCDR3) && secondMatchCDR3toCDR1Map.containsKey(betaCDR3))) {
+//                continue;
+//            }
+//            if(firstMatchCDR3toCDR1Map.get(alphaCDR3).equals(secondMatchCDR3toCDR1Map.get(betaCDR3))){
+//                if(firstMatchCDR3toCDR1Map.get(betaCDR3).equals(secondMatchCDR3toCDR1Map.get(alphaCDR3))){
+//                    dualMatchesMap.put(alphaCDR3, firstMatchCDR3toCDR1Map.get(alphaCDR3));
+//                    dualMatchesMap.put(betaCDR3, firstMatchCDR3toCDR1Map.get(betaCDR3));
+//                }
+//            }
+//        }
+//        System.out.println("Second pass mapping made. Dual CDR3/CDR1 pairings found.");
+//
+//        Instant stop = Instant.now();
+//        //results for first map
+//        System.out.println("RESULTS FOR FIRST PASS MATCHING");
+//        List<List<String>> allResults = new ArrayList<>();
+//        Integer trueCount = 0;
+//        Iterator iter = firstMatchesMap.keySet().iterator();
+//
+//        while(iter.hasNext()){
+//            Boolean proven = false;
+//            List<String> tmp = new ArrayList<>();
+//            tmp.add(iter.next().toString());
+//            tmp.add(iter.next().toString());
+//            tmp.add(firstMatchesMap.get(Integer.valueOf(tmp.get(0))).toString());
+//            tmp.add(firstMatchesMap.get(Integer.valueOf(tmp.get(1))).toString());
+//            if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(2)))){
+//                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(3)))){
+//                    proven = true;
+//                }
+//            }
+//            else if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(3)))){
+//                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(2)))){
+//                    proven = true;
+//                }
+//            }
+//            tmp.add(proven.toString());
+//            allResults.add(tmp);
+//            if(proven){
+//                trueCount++;
+//            }
+//        }
+//
+//        List<String> comments = new ArrayList<>();
+//        comments.add("Plate size: " + samplePlate.getSize() + " wells");
+//        comments.add("Previous pairs found: " + previousMatches.size());
+//        comments.add("CDR1 matches attempted: " + allResults.size());
+//        double attemptRate = (double) allResults.size() / previousMatches.size();
+//        comments.add("Matching attempt rate: " + attemptRate);
+//        comments.add("Number of correct matches: " + trueCount);
+//        double correctRate = (double) trueCount / allResults.size();
+//        comments.add("Correct matching rate: " + correctRate);
+//        NumberFormat nf = NumberFormat.getInstance(Locale.US);
+//        Duration time = Duration.between(start, stop);
+//        time = time.plus(previousTime);
+//        comments.add("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
+//        for(String s: comments){
+//            System.out.println(s);
+//        }
+//
+//
+//
+//        List<String> headers = new ArrayList<>();
+//        headers.add("CDR3 alpha");
+//        headers.add("CDR3 beta");
+//        headers.add("first matched CDR1");
+//        headers.add("second matched CDR1");
+//        headers.add("Correct match?");
+//
+//        MatchingResult firstTest = new MatchingResult(samplePlate.getSourceFileName(),
+//                comments, headers, allResults, dualMatchesMap, time);
+//
+//        //results for dual map
+//        System.out.println("RESULTS FOR SECOND PASS MATCHING");
+//        allResults = new ArrayList<>();
+//        trueCount = 0;
+//        iter = dualMatchesMap.keySet().iterator();
+//        while(iter.hasNext()){
+//            Boolean proven = false;
+//            List<String> tmp = new ArrayList<>();
+//            tmp.add(iter.next().toString());
+//            tmp.add(iter.next().toString());
+//            tmp.add(dualMatchesMap.get(Integer.valueOf(tmp.get(0))).toString());
+//            tmp.add(dualMatchesMap.get(Integer.valueOf(tmp.get(1))).toString());
+//            if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(2)))){
+//                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(3)))){
+//                    proven = true;
+//                }
+//            }
+//            else if(alphaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(0))).equals(Integer.valueOf(tmp.get(3)))){
+//                if(betaCDR3toCDR1Map.get(Integer.valueOf(tmp.get(1))).equals(Integer.valueOf(tmp.get(2)))){
+//                    proven = true;
+//                }
+//            }
+//            tmp.add(proven.toString());
+//            allResults.add(tmp);
+//            if(proven){
+//                trueCount++;
+//            }
+//        }
+//
+//        comments = new ArrayList<>();
+//        comments.add("Plate size: " + samplePlate.getSize() + " wells");
+//        comments.add("Previous pairs found: " + previousMatches.size());
+//        comments.add("High overlap threshold: " + highThreshold);
+//        comments.add("Low overlap threshold: " + lowThreshold);
+//        comments.add("CDR1 matches attempted: " + allResults.size());
+//        attemptRate = (double) allResults.size() / previousMatches.size();
+//        comments.add("Matching attempt rate: " + attemptRate);
+//        comments.add("Number of correct matches: " + trueCount);
+//        correctRate = (double) trueCount / allResults.size();
+//        comments.add("Correct matching rate: " + correctRate);
+//        comments.add("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
+//
+//        for(String s: comments){
+//            System.out.println(s);
+//        }
+//
+//        System.out.println("Simulation time: " + nf.format(time.toSeconds()) + " seconds");
+//        MatchingResult dualTest = new MatchingResult(samplePlate.getSourceFileName(), comments, headers,
+//                allResults, dualMatchesMap, time);
+//        MatchingResult[] output = {firstTest, dualTest};
+//        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
@@ -621,62 +660,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<>();

src/main/java/UserInterface.java

@@ -1,694 +0,0 @@
-import org.apache.commons.cli.*;
-
-import java.io.IOException;
-import java.util.List;
-import java.util.Scanner;
-import java.util.InputMismatchException;
-import java.util.regex.Matcher;
-import java.util.regex.Pattern;
-
-//
-public class UserInterface {
-
-    final static Scanner sc = new Scanner(System.in);
-    static int input;
-    static boolean quit = false;
-
-    public static void main(String[] args) {
-        //for now, commenting out all the command line argument stuff.
-        // Refactoring to output files of graphs, so it would all need to change anyway.
-
-//        if(args.length != 0){
-//            //These command line options are a big mess
-//            //Really, I don't think command line tools are expected to work in this many different modes
-//            //making cells, making plates, and matching are the sort of thing that UNIX philosophy would say
-//            //should be three separate programs.
-//            //There might be a way to do it with option parameters?
-//
-//            Options mainOptions = new Options();
-//            Option makeCells = Option.builder("cells")
-//                    .longOpt("make-cells")
-//                    .desc("Makes a file of distinct cells")
-//                    .build();
-//            Option makePlate = Option.builder("plates")
-//                    .longOpt("make-plates")
-//                    .desc("Makes a sample plate file")
-//                    .build();
-//            Option matchCDR3 = Option.builder("match")
-//                    .longOpt("match-cdr3")
-//                    .desc("Match CDR3s. Requires a cell sample file and any number of plate files.")
-//                    .build();
-//            OptionGroup mainGroup = new OptionGroup();
-//            mainGroup.addOption(makeCells);
-//            mainGroup.addOption(makePlate);
-//            mainGroup.addOption(matchCDR3);
-//            mainGroup.setRequired(true);
-//            mainOptions.addOptionGroup(mainGroup);
-//
-//            //Reuse clones of this for other options groups, rather than making it lots of times
-//            Option outputFile = Option.builder("o")
-//                    .longOpt("output-file")
-//                    .hasArg()
-//                    .argName("filename")
-//                    .desc("Name of output file")
-//                    .build();
-//            mainOptions.addOption(outputFile);
-//
-//            //Options cellOptions = new Options();
-//            Option numCells = Option.builder("nc")
-//                    .longOpt("num-cells")
-//                    .desc("The number of distinct cells to generate")
-//                    .hasArg()
-//                    .argName("number")
-//                    .build();
-//            mainOptions.addOption(numCells);
-//            Option cdr1Freq = Option.builder("d")
-//                    .longOpt("peptide-diversity-factor")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("Number of distinct CDR3s for every CDR1")
-//                    .build();
-//            mainOptions.addOption(cdr1Freq);
-//            //Option cellOutput = (Option) outputFile.clone();
-//            //cellOutput.setRequired(true);
-//            //mainOptions.addOption(cellOutput);
-//
-//            //Options plateOptions = new Options();
-//            Option inputCells = Option.builder("c")
-//                    .longOpt("cell-file")
-//                    .hasArg()
-//                    .argName("file")
-//                    .desc("The cell sample file used for filling wells")
-//                    .build();
-//            mainOptions.addOption(inputCells);
-//            Option numWells = Option.builder("w")
-//                    .longOpt("num-wells")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("The number of wells on each plate")
-//                    .build();
-//            mainOptions.addOption(numWells);
-//            Option numPlates = Option.builder("np")
-//                    .longOpt("num-plates")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("The number of plate files to output")
-//                    .build();
-//            mainOptions.addOption(numPlates);
-//            //Option plateOutput = (Option) outputFile.clone();
-//            //plateOutput.setRequired(true);
-//            //plateOutput.setDescription("Prefix for plate output filenames");
-//            //mainOptions.addOption(plateOutput);
-//            Option plateErr = Option.builder("err")
-//                    .longOpt("drop-out-rate")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("Well drop-out rate. (Probability between 0 and 1)")
-//                    .build();
-//            mainOptions.addOption(plateErr);
-//            Option plateConcentrations = Option.builder("t")
-//                    .longOpt("t-cells-per-well")
-//                    .hasArgs()
-//                    .argName("number 1, number 2, ...")
-//                    .desc("Number of T cells per well for each plate section")
-//                    .build();
-//            mainOptions.addOption(plateConcentrations);
-//
-////different distributions, mutually exclusive
-//            OptionGroup plateDistributions = new OptionGroup();
-//            Option plateExp = Option.builder("exponential")
-//                    .desc("Sample from distinct cells with exponential frequency distribution")
-//                    .build();
-//            plateDistributions.addOption(plateExp);
-//            Option plateGaussian = Option.builder("gaussian")
-//                    .desc("Sample from distinct cells with gaussain frequency distribution")
-//                    .build();
-//            plateDistributions.addOption(plateGaussian);
-//            Option platePoisson = Option.builder("poisson")
-//                    .desc("Sample from distinct cells with poisson frequency distribution")
-//                    .build();
-//            plateDistributions.addOption(platePoisson);
-//            mainOptions.addOptionGroup(plateDistributions);
-//
-//            Option plateStdDev = Option.builder("stddev")
-//                    .desc("Standard deviation for gaussian distribution")
-//                    .hasArg()
-//                    .argName("number")
-//                    .build();
-//            mainOptions.addOption(plateStdDev);
-//
-//            Option plateLambda = Option.builder("lambda")
-//                    .desc("Lambda for exponential distribution")
-//                    .hasArg()
-//                    .argName("number")
-//                    .build();
-//            mainOptions.addOption(plateLambda);
-//
-//
-//
-////
-////            String cellFile, String filename, Double stdDev,
-////                    Integer numWells, Integer numSections,
-////                    Integer[] concentrations, Double dropOutRate
-////
-//
-//            //Options matchOptions = new Options();
-//            inputCells.setDescription("The cell sample file to be used for matching.");
-//            mainOptions.addOption(inputCells);
-//            Option lowThresh = Option.builder("low")
-//                    .longOpt("low-threshold")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("Sets the minimum occupancy overlap to attempt matching")
-//                    .build();
-//            mainOptions.addOption(lowThresh);
-//            Option highThresh = Option.builder("high")
-//                    .longOpt("high-threshold")
-//                    .hasArg()
-//                    .argName("number")
-//                    .desc("Sets the maximum occupancy overlap to attempt matching")
-//                    .build();
-//            mainOptions.addOption(highThresh);
-//            Option occDiff = Option.builder("occdiff")
-//                    .longOpt("occupancy-difference")
-//                    .hasArg()
-//                    .argName("Number")
-//                    .desc("Maximum difference in alpha/beta occupancy to attempt matching")
-//                    .build();
-//            mainOptions.addOption(occDiff);
-//            Option overlapPer = Option.builder("ovper")
-//                    .longOpt("overlap-percent")
-//                    .hasArg()
-//                    .argName("Percent")
-//                    .desc("Minimum overlap percent to attempt matching (0 -100)")
-//                    .build();
-//            mainOptions.addOption(overlapPer);
-//            Option inputPlates = Option.builder("p")
-//                    .longOpt("plate-files")
-//                    .hasArgs()
-//                    .desc("Plate files to match")
-//                    .build();
-//            mainOptions.addOption(inputPlates);
-//
-//
-//
-//            CommandLineParser parser = new DefaultParser();
-//            try {
-//                CommandLine line = parser.parse(mainOptions, args);
-//                if(line.hasOption("match")){
-//                    //line = parser.parse(mainOptions, args);
-//                    String cellFile = line.getOptionValue("c");
-//                    Integer lowThreshold = Integer.valueOf(line.getOptionValue(lowThresh));
-//                    Integer highThreshold = Integer.valueOf(line.getOptionValue(highThresh));
-//                    Integer occupancyDifference = Integer.valueOf(line.getOptionValue(occDiff));
-//                    Integer overlapPercent = Integer.valueOf(line.getOptionValue(overlapPer));
-//                    for(String plate: line.getOptionValues("p")) {
-//                        matchCDR3s(cellFile, plate, lowThreshold, highThreshold, occupancyDifference, overlapPercent);
-//                    }
-//                }
-//                else if(line.hasOption("cells")){
-//                    //line = parser.parse(mainOptions, args);
-//                    String filename = line.getOptionValue("o");
-//                    Integer numDistCells = Integer.valueOf(line.getOptionValue("nc"));
-//                    Integer freq = Integer.valueOf(line.getOptionValue("d"));
-//                    makeCells(filename, numDistCells, freq);
-//                }
-//                else if(line.hasOption("plates")){
-//                    //line = parser.parse(mainOptions, args);
-//                    String cellFile = line.getOptionValue("c");
-//                    String filenamePrefix = line.getOptionValue("o");
-//                    Integer numWellsOnPlate = Integer.valueOf(line.getOptionValue("w"));
-//                    Integer numPlatesToMake = Integer.valueOf(line.getOptionValue("np"));
-//                    String[] concentrationsToUseString = line.getOptionValues("t");
-//                    Integer numSections = concentrationsToUseString.length;
-//
-//                    Integer[] concentrationsToUse = new Integer[numSections];
-//                    for(int i = 0; i <numSections; i++){
-//                        concentrationsToUse[i] = Integer.valueOf(concentrationsToUseString[i]);
-//                    }
-//                    Double dropOutRate = Double.valueOf(line.getOptionValue("err"));
-//                    if(line.hasOption("exponential")){
-//                        Double lambda = Double.valueOf(line.getOptionValue("lambda"));
-//                        for(int i = 1; i <= numPlatesToMake; i++){
-//                            makePlateExp(cellFile, filenamePrefix + i, lambda, numWellsOnPlate,
-//                                     concentrationsToUse,dropOutRate);
-//                        }
-//                    }
-//                    else if(line.hasOption("gaussian")){
-//                        Double stdDev = Double.valueOf(line.getOptionValue("std-dev"));
-//                        for(int i = 1; i <= numPlatesToMake; i++){
-//                            makePlate(cellFile, filenamePrefix + i, stdDev, numWellsOnPlate,
-//                                     concentrationsToUse,dropOutRate);
-//                        }
-//
-//                    }
-//                    else if(line.hasOption("poisson")){
-//                        for(int i = 1; i <= numPlatesToMake; i++){
-//                            makePlatePoisson(cellFile, filenamePrefix + i, numWellsOnPlate,
-//                                     concentrationsToUse,dropOutRate);
-//                        }
-//                    }
-//                }
-//            }
-//            catch (ParseException exp) {
-//                System.err.println("Parsing failed.  Reason: " + exp.getMessage());
-//            }
-//        }
-//        else {
-            while (!quit) {
-                System.out.println();
-                System.out.println("--------BiGPairSEQ SIMULATOR--------");
-                System.out.println("ALPHA/BETA T CELL RECEPTOR MATCHING");
-                System.out.println("  USING WEIGHTED BIPARTITE GRAPHS  ");
-                System.out.println("------------------------------------");
-                System.out.println("Please select an option:");
-                System.out.println("1) Generate a population of distinct cells");
-                System.out.println("2) Generate a sample plate of T cells");
-                System.out.println("3) Generate CDR3 alpha/beta occupancy data and overlap graph");
-                System.out.println("4) Simulate bipartite graph CDR3 alpha/beta matching (BiGpairSEQ)");
-                //Need to re-do the CDR3/CDR1 matching to correspond to new pattern
-                //System.out.println("5) Generate CDR3/CDR1 occupancy graph");
-                //System.out.println("6) Simulate CDR3/CDR1 T cell matching");
-                System.out.println("9) About/Acknowledgments");
-                System.out.println("0) Exit");
-                try {
-                    input = sc.nextInt();
-                    switch (input) {
-                        case 1 -> makeCells();
-                        case 2 -> makePlate();
-                        case 3 -> makeCDR3Graph();
-                        case 4 -> matchCDR3s();
-                        //case 6 -> matchCellsCDR1();
-                        case 9 -> acknowledge();
-                        case 0 -> quit = true;
-                        default -> throw new InputMismatchException("Invalid input.");
-                    }
-                } catch (InputMismatchException | IOException ex) {
-                    System.out.println(ex);
-                    sc.next();
-                }
-            }
-            sc.close();
-//        }
-    }
-
-    private static void makeCells() {
-        String filename = null;
-        Integer numCells = 0;
-        Integer cdr1Freq = 1;
-        try {
-            System.out.println("\nSimulated T-Cells consist of integer values representing:\n" +
-                    "* a pair of alpha and beta CDR3 peptides (unique within simulated population)\n" +
-                    "* a pair of alpha and beta CDR1 peptides (not necessarily unique).");
-            System.out.println("\nThe cells will be written to a CSV file.");
-            System.out.print("Please enter a file name: ");
-            filename = sc.next();
-            System.out.println("\nCDR3 sequences are more diverse than CDR1 sequences.");
-            System.out.println("Please enter the factor by which distinct CDR3s outnumber CDR1s: ");
-            cdr1Freq = sc.nextInt();
-            System.out.print("\nPlease enter the number of T-cells to generate: ");
-            numCells = sc.nextInt();
-            if(numCells <= 0){
-                throw new InputMismatchException("Number of cells must be a positive integer.");
-            }
-        } catch (InputMismatchException ex) {
-            System.out.println(ex);
-            sc.next();
-        }
-        CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
-        assert filename != null;
-        CellFileWriter writer = new CellFileWriter(filename, sample);
-        writer.writeCellsToFile();
-        System.gc();
-    }
-
-//    //for calling from command line
-//    private static void makeCells(String filename, Integer numCells, Integer cdr1Freq){
-//        CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
-//        CellFileWriter writer = new CellFileWriter(filename, sample);
-//        writer.writeCellsToFile();
-//    }
-//
-//    private static void makePlateExp(String cellFile, String filename, Double lambda,
-//                                  Integer numWells, Integer[] concentrations, Double dropOutRate){
-//        CellFileReader cellReader = new CellFileReader(cellFile);
-//        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
-//        samplePlate.fillWellsExponential(cellReader.getFilename(), cellReader.getCells(), lambda);
-//        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
-//        writer.writePlateFile();
-//    }
-//
-//    private static void makePlatePoisson(String cellFile, String filename, Integer numWells,
-//                                         Integer[] concentrations, Double dropOutRate){
-//        CellFileReader cellReader = new CellFileReader(cellFile);
-//        Double stdDev = Math.sqrt(cellReader.getCellCount());
-//        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
-//        samplePlate.fillWells(cellReader.getFilename(), cellReader.getCells(), stdDev);
-//        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
-//        writer.writePlateFile();
-//    }
-//
-//    private static void makePlate(String cellFile, String filename, Double stdDev,
-//                                  Integer numWells, Integer[] concentrations, Double dropOutRate){
-//        CellFileReader cellReader = new CellFileReader(cellFile);
-//        Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
-//        samplePlate.fillWells(cellReader.getFilename(), cellReader.getCells(), stdDev);
-//        PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
-//        writer.writePlateFile();
-//    }
-
-    //Output a CSV of sample plate
-    private static void makePlate() {
-        String cellFile = null;
-        String filename = null;
-        Double stdDev = 0.0;
-        Integer numWells = 0;
-        Integer numSections;
-        Integer[] concentrations = {1};
-        Double dropOutRate = 0.0;
-        boolean poisson = false;
-        boolean exponential = false;
-        double lambda = 1.5;
-        try {
-            System.out.println("\nSimulated sample plates consist of:");
-            System.out.println("* a number of wells");
-            System.out.println("  * separated into one or more sections");
-            System.out.println("    * each of which has a set quantity of cells per well");
-            System.out.println("    * selected from a statistical distribution of distinct cells");
-            System.out.println("    * with a set dropout rate for individual sequences within a cell");
-            System.out.println("\nMaking a sample plate requires a population of distinct cells");
-            System.out.print("Please enter name of an existing cell sample file: ");
-            cellFile = sc.next();
-            System.out.println("\nThe sample plate will be written to a CSV file");
-            System.out.print("Please enter a name for the output file: ");
-            filename = sc.next();
-            System.out.println("\nSelect T-cell frequency distribution function");
-            System.out.println("1) Poisson");
-            System.out.println("2) Gaussian");
-            System.out.println("3) Exponential");
-            System.out.println("(Note: approximate distribution in original paper is exponential, lambda = 0.6)");
-            System.out.println("(lambda value approximated from slope of log-log graph in figure 4c)");
-            System.out.println("(Note: wider distributions are more memory intensive to match)");
-            System.out.print("Enter selection value: ");
-            input = sc.nextInt();
-            switch (input) {
-                case 1 -> poisson = true;
-                case 2 -> {
-                    System.out.println("How many distinct T-cells within one standard deviation of peak frequency?");
-                    System.out.println("(Note: wider distributions are more memory intensive to match)");
-                    stdDev = sc.nextDouble();
-                    if (stdDev <= 0.0) {
-                        throw new InputMismatchException("Value must be positive.");
-                    }
-                }
-                case 3 -> {
-                    exponential = true;
-                    System.out.println("Please enter lambda value for exponential distribution.");
-                    lambda = sc.nextDouble();
-                    if (lambda <= 0.0) {
-                        throw new InputMismatchException("Value must be positive.");
-                    }
-                }
-                default -> {
-                    System.out.println("Invalid input. Defaulting to exponential.");
-                    exponential = true;
-                }
-            }
-            System.out.print("\nNumber of wells on plate: ");
-            numWells = sc.nextInt();
-            if(numWells < 1){
-                throw new InputMismatchException("No wells on plate");
-            }
-            System.out.println("\nThe plate can be evenly sectioned to allow multiple concentrations of T-cells/well");
-            System.out.println("How many sections would you like to make (minimum 1)?");
-            numSections = sc.nextInt();
-            if(numSections < 1) {
-                throw new InputMismatchException("Too few sections.");
-            }
-            else if (numSections > numWells) {
-                throw new InputMismatchException("Cannot have more sections than wells.");
-            }
-            int i = 1;
-            concentrations = new Integer[numSections];
-            while(numSections > 0) {
-                System.out.print("Enter number of T-cells per well in section " + i +": ");
-                concentrations[i - 1] = sc.nextInt();
-                i++;
-                numSections--;
-            }
-            System.out.println("\nErrors in amplification can induce a well dropout rate for sequences");
-            System.out.print("Enter well dropout rate (0.0 to 1.0): ");
-            dropOutRate = sc.nextDouble();
-            if(dropOutRate < 0.0 || dropOutRate > 1.0) {
-                throw new InputMismatchException("The well dropout rate must be in the range [0.0, 1.0]");
-            }
-        }catch(InputMismatchException ex){
-            System.out.println(ex);
-            sc.next();
-        }
-        System.out.println("Reading Cell Sample file: " + cellFile);
-        assert cellFile != null;
-        CellFileReader cellReader = new CellFileReader(cellFile);
-        if(exponential){
-            Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
-            samplePlate.fillWellsExponential(cellReader.getFilename(), cellReader.getCells(), lambda);
-            PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
-            writer.writePlateFile();
-        }
-        else {
-            if (poisson) {
-                stdDev = Math.sqrt(cellReader.getCellCount()); //gaussian with square root of elements approximates poisson
-            }
-            Plate samplePlate = new Plate(numWells, dropOutRate, concentrations);
-            samplePlate.fillWells(cellReader.getFilename(), cellReader.getCells(), stdDev);
-            assert filename != null;
-            PlateFileWriter writer = new PlateFileWriter(filename, samplePlate);
-            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 makeCDR3Graph() {
-        String filename = null;
-        String cellFile = null;
-        String plateFile = null;
-
-        try {
-            String str = "\nGenerating bipartite weighted graph encoding occupancy overlap data ";
-            str = str.concat("\nrequires a cell sample file and a sample plate file.");
-            System.out.println(str);
-            System.out.print("\nPlease enter name of an existing cell sample file: ");
-            cellFile = sc.next();
-            System.out.print("\nPlease enter name of an existing sample plate file: ");
-            plateFile = sc.next();
-            System.out.println("\nThe graph and occupancy data will be written to a serialized binary file.");
-            System.out.print("Please enter a name for the output file: ");
-            filename = sc.next();
-        } catch (InputMismatchException ex) {
-            System.out.println(ex);
-            sc.next();
-        }
-        System.out.println("Reading Cell Sample file: " + cellFile);
-        assert cellFile != null;
-        CellFileReader cellReader = new CellFileReader(cellFile);
-        System.out.println("Reading Sample Plate file: " + plateFile);
-        assert plateFile != null;
-        PlateFileReader plateReader = new PlateFileReader(plateFile);
-        Plate plate = new Plate(plateReader.getFilename(), plateReader.getWells());
-        if (cellReader.getCells().size() == 0){
-            System.out.println("No cell sample found.");
-            System.out.println("Returning to main menu.");
-        }
-        else if(plate.getWells().size() == 0 || plate.getConcentrations().length == 0){
-            System.out.println("No sample plate found.");
-            System.out.println("Returning to main menu.");
-        }
-        else{
-            List<Integer[]> cells = cellReader.getCells();
-            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();
-        }
-    }
-
-    //Simulate matching and output CSV file of results
-    private static void matchCDR3s() throws IOException {
-        String filename = null;
-        String dataFilename = null;
-        Integer lowThreshold = 0;
-        Integer highThreshold = Integer.MAX_VALUE;
-        Integer maxOccupancyDiff = Integer.MAX_VALUE;
-        Integer minOverlapPercent = 0;
-        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();
-            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();
-            System.out.println("\nWhat is the minimum number of CDR3 alpha/beta overlap wells to attempt matching?");
-            lowThreshold = sc.nextInt();
-            if(lowThreshold < 1){
-                throw new InputMismatchException("Minimum value for low threshold set to 1");
-            }
-            System.out.println("\nWhat is the maximum number of CDR3 alpha/beta overlap wells to attempt matching?");
-            highThreshold = sc.nextInt();
-            System.out.println("\nWhat is the maximum difference in alpha/beta occupancy to attempt matching?");
-            maxOccupancyDiff = sc.nextInt();
-            System.out.println("\nWell overlap percentage = pair overlap / sequence occupancy");
-            System.out.println("What is the minimum well overlap percentage to attempt matching? (0 to 100)");
-            minOverlapPercent = sc.nextInt();
-            if (minOverlapPercent < 0 || minOverlapPercent > 100) {
-                throw new InputMismatchException("Value outside range. Minimum percent set to 0");
-            }
-        } catch (InputMismatchException ex) {
-            System.out.println(ex);
-            sc.next();
-        }
-        //read object data from file
-        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");
-        assert dataFilename != null;
-        GraphDataObjectReader dataReader = new GraphDataObjectReader(dataFilename);
-        GraphWithMapData data = dataReader.getData();
-        //set source file name
-        data.setSourceFilename(dataFilename);
-        //simulate matching
-        MatchingResult results = Simulator.matchCDR3s(data, dataFilename, lowThreshold, highThreshold, maxOccupancyDiff,
-                minOverlapPercent, true);
-        //write results to file
-        assert filename != null;
-        MatchingFileWriter writer = new MatchingFileWriter(filename, results);
-        System.out.println("Writing results to file");
-        writer.writeResultsToFile();
-        System.out.println("Results written to file: " + filename);
-        System.gc();
-    }
-
-    ///////
-    //Rewrite this to fit new matchCDR3 method with file I/O
-    ///////
-//    public static void matchCellsCDR1(){
-//    /*
-//    The idea here is that we'll get the CDR3 alpha/beta matches first. Then we'll try to match CDR3s to CDR1s by
-//    looking at the top two matches for each CDR3. If CDR3s in the same cell simply swap CDR1s, we assume a correct
-//    match
-//     */
-//        String filename = null;
-//        String preliminaryResultsFilename = null;
-//        String cellFile = null;
-//        String plateFile = null;
-//        Integer lowThresholdCDR3 = 0;
-//        Integer highThresholdCDR3 = Integer.MAX_VALUE;
-//        Integer maxOccupancyDiffCDR3 = 96; //no filtering if max difference is all wells by default
-//        Integer minOverlapPercentCDR3 = 0; //no filtering if min percentage is zero by default
-//        Integer lowThresholdCDR1 = 0;
-//        Integer highThresholdCDR1 = Integer.MAX_VALUE;
-//        boolean outputCDR3Matches = false;
-//        try {
-//            System.out.println("\nSimulated experiment requires a cell sample file and a sample plate file.");
-//            System.out.print("Please enter name of an existing cell sample file: ");
-//            cellFile = sc.next();
-//            System.out.print("Please enter name of an existing sample plate file: ");
-//            plateFile = 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();
-//            System.out.println("What is the minimum number of CDR3 alpha/beta overlap wells to attempt matching?");
-//            lowThresholdCDR3 = sc.nextInt();
-//            if(lowThresholdCDR3 < 1){
-//                throw new InputMismatchException("Minimum value for low threshold is 1");
-//            }
-//            System.out.println("What is the maximum number of CDR3 alpha/beta overlap wells to attempt matching?");
-//            highThresholdCDR3 = sc.nextInt();
-//            System.out.println("What is the maximum difference in CDR3 alpha/beta occupancy to attempt matching?");
-//            maxOccupancyDiffCDR3 = sc.nextInt();
-//            System.out.println("What is the minimum CDR3 overlap percentage to attempt matching? (0 - 100)");
-//            minOverlapPercentCDR3 = sc.nextInt();
-//            if (minOverlapPercentCDR3 < 0 || minOverlapPercentCDR3 > 100) {
-//                throw new InputMismatchException("Value outside range. Minimum percent set to 0");
-//            }
-//            System.out.println("What is the minimum number of CDR3/CDR1 overlap wells to attempt matching?");
-//            lowThresholdCDR1 = sc.nextInt();
-//            if(lowThresholdCDR1 < 1){
-//                throw new InputMismatchException("Minimum value for low threshold is 1");
-//            }
-//            System.out.println("What is the maximum number of CDR3/CDR1 overlap wells to attempt matching?");
-//            highThresholdCDR1 = sc.nextInt();
-//            System.out.println("Matching CDR3s to CDR1s requires first matching CDR3 alpha/betas.");
-//            System.out.println("Output a file for CDR3 alpha/beta match results as well?");
-//            System.out.print("Please enter y/n: ");
-//            String ans = sc.next();
-//            Pattern pattern = Pattern.compile("(?:yes|y)", Pattern.CASE_INSENSITIVE);
-//            Matcher matcher = pattern.matcher(ans);
-//            if(matcher.matches()){
-//                outputCDR3Matches = true;
-//                System.out.println("Please enter filename for CDR3 alpha/beta match results");
-//                preliminaryResultsFilename = sc.next();
-//                System.out.println("CDR3 alpha/beta matches will be output to file");
-//            }
-//            else{
-//                System.out.println("CDR3 alpha/beta matches will not be output to file");
-//            }
-//        } catch (InputMismatchException ex) {
-//            System.out.println(ex);
-//            sc.next();
-//        }
-//        CellFileReader cellReader = new CellFileReader(cellFile);
-//        PlateFileReader plateReader = new PlateFileReader(plateFile);
-//        Plate plate = new Plate(plateReader.getFilename(), plateReader.getWells());
-//        if (cellReader.getCells().size() == 0){
-//            System.out.println("No cell sample found.");
-//            System.out.println("Returning to main menu.");
-//        }
-//        else if(plate.getWells().size() == 0){
-//            System.out.println("No sample plate found.");
-//            System.out.println("Returning to main menu.");
-//
-//        }
-//        else{
-//            if(highThresholdCDR3 >= plate.getSize()){
-//                highThresholdCDR3 = plate.getSize() - 1;
-//            }
-//            if(highThresholdCDR1 >= plate.getSize()){
-//                highThresholdCDR1 = plate.getSize() - 1;
-//            }
-//            List<Integer[]> cells = cellReader.getCells();
-//            MatchingResult preliminaryResults = Simulator.matchCDR3s(cells, plate, lowThresholdCDR3, highThresholdCDR3,
-//                    maxOccupancyDiffCDR3, minOverlapPercentCDR3, true);
-//            MatchingResult[] results = Simulator.matchCDR1s(cells, plate, lowThresholdCDR1,
-//                    highThresholdCDR1, preliminaryResults);
-//            MatchingFileWriter writer = new MatchingFileWriter(filename + "_FirstPass", results[0]);
-//            writer.writeResultsToFile();
-//            writer = new MatchingFileWriter(filename + "_SecondPass", results[1]);
-//            writer.writeResultsToFile();
-//            if(outputCDR3Matches){
-//                writer = new MatchingFileWriter(preliminaryResultsFilename, preliminaryResults);
-//                writer.writeResultsToFile();
-//            }
-//        }
-//    }
-
-    private static void acknowledge(){
-        System.out.println("This program simulates BiGpairSEQ, a graph theory based adaptation");
-        System.out.println("of the pairSEQ algorithm for pairing T cell receptor sequences.");
-        System.out.println();
-        System.out.println("For full documentation, view readme.md file distributed with this code");
-        System.out.println("or visit https://gitea.ejsf.synology.me/efischer/BiGpairSEQ.");
-        System.out.println();
-        System.out.println("pairSEQ citation:");
-        System.out.println("Howie, B., Sherwood, A. M., et. al.");
-        System.out.println("High-throughput pairing of T cell receptor alpha and beta sequences.");
-        System.out.println("Sci. Transl. Med. 7, 301ra131 (2015)");
-        System.out.println();
-        System.out.println("BiGpairSEQ_Sim by Eugene Fischer, 2021-2022");
-    }
-}