Compare commits
27 Commits
| Author | SHA1 | Date | |
|---|---|---|---|
| b4cc240048 | |||
| ff72c9b359 | |||
| 88eb8aca50 | |||
| 98bf452891 | |||
| c2db4f87c1 | |||
| 8935407ade | |||
| 9fcc20343d | |||
| e4d094d796 | |||
| f385ebc31f | |||
| 8745550e11 | |||
| 41805135b3 | |||
| 373a5e02f9 | |||
| 7f18311054 | |||
| bcb816c3e6 | |||
| dad0fd35fd | |||
| 35d580cfcf | |||
| ab8d98ed81 | |||
| 3d9890e16a | |||
| dd64ac2731 | |||
| a5238624f1 | |||
| d8ba42b801 | |||
| 8edd89d784 | |||
| 2829b88689 | |||
| 108b0ec13f | |||
| a8b58d3f79 | |||
| bf64d57731 | |||
| c068c3db3c |
108
readme.md
108
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 [weightd bipartite graph](https://en.wikipedia.org/wiki/Bipartite_graph) representing the sample plate.
|
||||
BiGpairSEQ creates a [weighted 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.)
|
||||
@@ -29,17 +29,13 @@ Unfortunately, it's a fairly new algorithm, and not yet implemented by the graph
|
||||
So this program instead uses the Fibonacci heap-based algorithm of Fredman and Tarjan (1987), which has a worst-case
|
||||
runtime of **O(n (n log(n) + m))**. The algorithm is implemented as described in Melhorn and Näher (1999).
|
||||
|
||||
The current version of the program uses a pairing heap instead of a Fibonacci heap for its priority queue,
|
||||
which has lower theoretical efficiency but also lower complexity overhead, and is often equivalently performant
|
||||
in practice.
|
||||
|
||||
## USAGE
|
||||
|
||||
### 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/)
|
||||
BiGpairSEQ_Sim is an executable .jar file. Requires Java 14 or higher. [OpenJDK 17](https://jdk.java.net/17/)
|
||||
recommended.
|
||||
|
||||
Run with the command:
|
||||
@@ -58,17 +54,30 @@ main menu looks like this:
|
||||
```
|
||||
--------BiGPairSEQ SIMULATOR--------
|
||||
ALPHA/BETA T CELL RECEPTOR MATCHING
|
||||
USING WEIGHTED BIPARTITE GRAPHS
|
||||
USING WEIGHTED BIPARTITE GRAPHS
|
||||
------------------------------------
|
||||
Please select an option:
|
||||
1) Generate a population of distinct cells
|
||||
2) Generate a sample plate of T cells
|
||||
3) Generate CDR3 alpha/beta occupancy data and overlap graph
|
||||
4) Simulate bipartite graph CDR3 alpha/beta matching (BiGpairSEQ)
|
||||
8) Options
|
||||
9) About/Acknowledgments
|
||||
0) Exit
|
||||
```
|
||||
|
||||
By default, the Options menu looks like this:
|
||||
```
|
||||
--------------OPTIONS---------------
|
||||
1) Turn on cell sample file caching
|
||||
2) Turn on plate file caching
|
||||
3) Turn on graph/data file caching
|
||||
4) Turn off serialized binary graph output
|
||||
5) Turn on GraphML graph output
|
||||
6) Maximum weight matching algorithm options
|
||||
0) Return to main menu
|
||||
```
|
||||
|
||||
### INPUT/OUTPUT
|
||||
|
||||
To run the simulation, the program reads and writes 4 kinds of files:
|
||||
@@ -77,19 +86,26 @@ To run the simulation, the program reads and writes 4 kinds of files:
|
||||
* 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.
|
||||
These files are often generated in sequence. 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.
|
||||
|
||||
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.
|
||||
To save file I/O time, the most recent instance of each of these four
|
||||
files either generated or read from disk can be cached in program memory. When caching is active, 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,
|
||||
or caching is turned off for that file type. 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.
|
||||
(Note that cached Graph/Data files must be transformed back into their original state after a matching experiment, which
|
||||
may take some time. Whether file I/O or graph transformation takes longer for graph/data files is likely to be
|
||||
device-specific.)
|
||||
|
||||
The program's caching behavior can be controlled in the Options menu. By default, all caching is OFF.
|
||||
|
||||
The program can optionally output Graph/Data files in .GraphML format (.graphml) for data portability. This can be
|
||||
turned on in the Options menu. By default, GraphML output is OFF.
|
||||
|
||||
---
|
||||
#### Cell Sample Files
|
||||
Cell Sample files consist of any number of distinct "T cells." Every cell contains
|
||||
four sequences: Alpha CDR3, Beta CDR3, Alpha CDR1, Beta CDR1. The sequences are represented by
|
||||
@@ -107,7 +123,6 @@ Comments are preceded by `#`
|
||||
|
||||
Structure:
|
||||
|
||||
---
|
||||
# Sample contains 1 unique CDR1 for every 4 unique CDR3s.
|
||||
| Alpha CDR3 | Beta CDR3 | Alpha CDR1 | Beta CDR1 |
|
||||
|---|---|---|---|
|
||||
@@ -131,11 +146,14 @@ 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 approximately exponential with a lambda ~0.6. (Howie, et al. 2015))*
|
||||
* Total number of wells on the plate
|
||||
* Number of sections on plate
|
||||
* Number of T cells per well
|
||||
* per section, if more than one section
|
||||
* Well populations random or fixed
|
||||
* If random, minimum and maximum population sizes
|
||||
* If fixed
|
||||
* Number of sections on plate
|
||||
* Number of T cells per well
|
||||
* per section, if more than one section
|
||||
* Dropout rate
|
||||
|
||||
Files are in CSV format. There are no header labels. Every row represents a well.
|
||||
@@ -149,7 +167,6 @@ Dropout sequences are replaced with the value `-1`. Comments are preceded by `#`
|
||||
|
||||
Structure:
|
||||
|
||||
---
|
||||
```
|
||||
# Cell source file name:
|
||||
# Each row represents one well on the plate
|
||||
@@ -178,14 +195,19 @@ 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.)
|
||||
|
||||
These files do not have a human-readable structure, and are not portable to other programs. (Export of graphs in a
|
||||
portable data format may be implemented in the future. The tricky part is encoding the necessary metadata.)
|
||||
These files do not have a human-readable structure, and are not portable to other programs.
|
||||
|
||||
(For portability to other software, turn on GraphML output in the Options menu. This will produce a .graphml file
|
||||
for the weighted graph, with vertex attributes sequence, type, and occupancy data.)
|
||||
|
||||
---
|
||||
|
||||
#### Matching Results Files
|
||||
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.
|
||||
Matching results files consist of the results of a BiGpairSEQ matching simulation. Making them requires a serialized
|
||||
binary Graph/Data file (.ser). (Because .graphML files are larger than .ser files, BiGpairSEQ_Sim supports .graphML
|
||||
output only. Graph/data input must use a serialized binary.)
|
||||
|
||||
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:
|
||||
@@ -200,7 +222,6 @@ Options when running a BiGpairSEQ simulation of CDR3 alpha/beta matching:
|
||||
|
||||
Example output:
|
||||
|
||||
---
|
||||
```
|
||||
# Source Sample Plate file: 4MilCellsPlate.csv
|
||||
# Source Graph and Data file: 4MilCellsPlateGraph.ser
|
||||
@@ -251,27 +272,30 @@ 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~~ ~~UNABANDONED~~ DONE
|
||||
* ~~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.
|
||||
* 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.
|
||||
* It is possible, though the modifications to the graph incur their own performance penalties. Need testing to see which option is best.
|
||||
* ~~Test whether pairing heap (currently used) or Fibonacci heap is more efficient for priority queue in current matching algorithm~~ DONE
|
||||
* ~~in theory Fibonacci heap should be more efficient, but complexity overhead may eliminate theoretical advantage~~
|
||||
* ~~Add controllable heap-type parameter?~~
|
||||
* Parameter implemented. Fibonacci heap the current default.
|
||||
* ~~Implement sample plates with random numbers of T cells per well.~~ DONE
|
||||
* Possible BiGpairSEQ advantage over pairSEQ: BiGpairSEQ is resilient to variations in well population sizes on a sample plate; pairSEQ is not.
|
||||
* preliminary data suggests that BiGpairSEQ behaves roughly as though the whole plate had whatever the *average* well concentration is, but that's still speculative.
|
||||
* 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~~
|
||||
* _Got this working, but at the cost of a profoundly strange bug in graph occupancy filtering. Have reverted the repo until I can figure out what caused that. Given how easily Thingiverse transposes CSV matrices in R, might not even be worth fixing._
|
||||
* _Got this working, but at the cost of a profoundly strange bug in graph occupancy filtering. Have reverted the repo until I can figure out what caused that. Given how easily Thingiverse transposes CSV matrices in R, might not even be worth fixing.
|
||||
* ~~Enable GraphML output in addition to serialized object binaries, for data portability~~ DONE
|
||||
* ~~Custom vertex type with attribute for sequence occupancy?~~ ABANDONED
|
||||
* Have a branch where this is implemented, but there's a bug that broke matching. Don't currently have time to fix.
|
||||
* Re-implement command line arguments, to enable scripting and statistical simulation studies
|
||||
* Implement sample plates with random numbers of T cells per well.
|
||||
* Possible BiGpairSEQ advantage over pairSEQ: BiGpairSEQ is resilient to variations in well population sizes on a sample plate; pairSEQ is not.
|
||||
* preliminary data suggests that BiGpairSEQ behaves roughly as though the whole plate had whatever the *average* well concentration is, but that's still speculative.
|
||||
* Enable GraphML output in addition to serialized object binaries, for data portability
|
||||
* Custom vertex type with attribute for sequence occupancy?
|
||||
* Re-implement CDR1 matching method
|
||||
* Implement Duan and Su's maximum weight matching algorithm
|
||||
* Add controllable algorithm-type parameter?
|
||||
* Test whether pairing heap (currently used) or Fibonacci heap is more efficient for priority queue in current matching algorithm
|
||||
* in theory Fibonacci heap should be more efficient, but complexity overhead may eliminate theoretical advantage
|
||||
* Add controllable heap-type parameter?
|
||||
|
||||
|
||||
* Add controllable algorithm-type parameter?
|
||||
* This would be fun and valuable, but probably take more time than I have for a hobby project.
|
||||
|
||||
|
||||
## CITATIONS
|
||||
* Howie, B., Sherwood, A. M., et al. ["High-throughput pairing of T cell receptor alpha and beta sequences."](https://pubmed.ncbi.nlm.nih.gov/26290413/) Sci. Transl. Med. 7, 301ra131 (2015)
|
||||
|
||||
@@ -1,6 +1,6 @@
|
||||
import java.util.Random;
|
||||
|
||||
//main class. For choosing interface type and caching file data
|
||||
//main class. For choosing interface type and holding settings
|
||||
public class BiGpairSEQ {
|
||||
|
||||
private static final Random rand = new Random();
|
||||
@@ -10,6 +10,12 @@ public class BiGpairSEQ {
|
||||
private static String plateFilename = null;
|
||||
private static GraphWithMapData graphInMemory = null;
|
||||
private static String graphFilename = null;
|
||||
private static boolean cacheCells = false;
|
||||
private static boolean cachePlate = false;
|
||||
private static boolean cacheGraph = false;
|
||||
private static String priorityQueueHeapType = "FIBONACCI";
|
||||
private static boolean outputBinary = true;
|
||||
private static boolean outputGraphML = false;
|
||||
|
||||
public static void main(String[] args) {
|
||||
if (args.length == 0) {
|
||||
@@ -30,66 +36,141 @@ public class BiGpairSEQ {
|
||||
return cellSampleInMemory;
|
||||
}
|
||||
|
||||
public static void setCellSampleInMemory(CellSample cellSampleInMemory) {
|
||||
BiGpairSEQ.cellSampleInMemory = cellSampleInMemory;
|
||||
public static void setCellSampleInMemory(CellSample cellSample, String filename) {
|
||||
if(cellSampleInMemory != null) {
|
||||
clearCellSampleInMemory();
|
||||
}
|
||||
cellSampleInMemory = cellSample;
|
||||
cellFilename = filename;
|
||||
System.out.println("Cell sample file " + filename + " cached.");
|
||||
}
|
||||
|
||||
public static void clearCellSampleInMemory() {
|
||||
cellSampleInMemory = null;
|
||||
cellFilename = null;
|
||||
System.gc();
|
||||
System.out.println("Cell sample file cache cleared.");
|
||||
|
||||
}
|
||||
|
||||
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 setPlateInMemory(Plate plate, String filename) {
|
||||
if(plateInMemory != null) {
|
||||
clearPlateInMemory();
|
||||
}
|
||||
plateInMemory = plate;
|
||||
plateFilename = filename;
|
||||
System.out.println("Sample plate file " + filename + " cached.");
|
||||
}
|
||||
|
||||
public static void clearPlateInMemory() {
|
||||
plateInMemory = null;
|
||||
plateFilename = null;
|
||||
System.gc();
|
||||
System.out.println("Sample plate file cache cleared.");
|
||||
|
||||
}
|
||||
|
||||
public static String getPlateFilename() {
|
||||
return plateFilename;
|
||||
}
|
||||
|
||||
public static void setPlateFilename(String plateFilename) {
|
||||
BiGpairSEQ.plateFilename = plateFilename;
|
||||
|
||||
public static GraphWithMapData getGraphInMemory() {return graphInMemory;
|
||||
}
|
||||
|
||||
public static GraphWithMapData getGraphInMemory() {
|
||||
return graphInMemory;
|
||||
}
|
||||
|
||||
public static void setGraphInMemory(GraphWithMapData g) {
|
||||
public static void setGraphInMemory(GraphWithMapData g, String filename) {
|
||||
if (graphInMemory != null) {
|
||||
clearGraphInMemory();
|
||||
}
|
||||
graphInMemory = g;
|
||||
graphFilename = filename;
|
||||
System.out.println("Graph and data file " + filename + " cached.");
|
||||
}
|
||||
|
||||
public static void clearGraphInMemory() {
|
||||
graphInMemory = null;
|
||||
graphFilename = null;
|
||||
System.gc();
|
||||
System.out.println("Graph and data file cache cleared.");
|
||||
}
|
||||
|
||||
public static String getGraphFilename() {
|
||||
return graphFilename;
|
||||
}
|
||||
|
||||
public static void setGraphFilename(String filename) {
|
||||
graphFilename = filename;
|
||||
|
||||
public static boolean cacheCells() {
|
||||
return cacheCells;
|
||||
}
|
||||
|
||||
public static void setCacheCells(boolean cacheCells) {
|
||||
//if not caching, clear the memory
|
||||
if(!cacheCells){
|
||||
BiGpairSEQ.clearCellSampleInMemory();
|
||||
System.out.println("Cell sample file caching: OFF.");
|
||||
}
|
||||
else {
|
||||
System.out.println("Cell sample file caching: ON.");
|
||||
}
|
||||
BiGpairSEQ.cacheCells = cacheCells;
|
||||
}
|
||||
|
||||
public static boolean cachePlate() {
|
||||
return cachePlate;
|
||||
}
|
||||
|
||||
public static void setCachePlate(boolean cachePlate) {
|
||||
//if not caching, clear the memory
|
||||
if(!cachePlate) {
|
||||
BiGpairSEQ.clearPlateInMemory();
|
||||
System.out.println("Sample plate file caching: OFF.");
|
||||
}
|
||||
else {
|
||||
System.out.println("Sample plate file caching: ON.");
|
||||
}
|
||||
BiGpairSEQ.cachePlate = cachePlate;
|
||||
}
|
||||
|
||||
public static boolean cacheGraph() {
|
||||
return cacheGraph;
|
||||
}
|
||||
|
||||
public static void setCacheGraph(boolean cacheGraph) {
|
||||
//if not caching, clear the memory
|
||||
if(!cacheGraph) {
|
||||
BiGpairSEQ.clearGraphInMemory();
|
||||
System.out.println("Graph/data file caching: OFF.");
|
||||
}
|
||||
else {
|
||||
System.out.println("Graph/data file caching: ON.");
|
||||
}
|
||||
BiGpairSEQ.cacheGraph = cacheGraph;
|
||||
}
|
||||
|
||||
public static String getPriorityQueueHeapType() {
|
||||
return priorityQueueHeapType;
|
||||
}
|
||||
|
||||
public static void setPairingHeap() {
|
||||
priorityQueueHeapType = "PAIRING";
|
||||
}
|
||||
|
||||
public static void setFibonacciHeap() {
|
||||
priorityQueueHeapType = "FIBONACCI";
|
||||
}
|
||||
|
||||
public static boolean outputBinary() {return outputBinary;}
|
||||
public static void setOutputBinary(boolean b) {outputBinary = b;}
|
||||
|
||||
public static boolean outputGraphML() {return outputGraphML;}
|
||||
public static void setOutputGraphML(boolean b) {outputGraphML = b;}
|
||||
|
||||
}
|
||||
|
||||
@@ -1,10 +1,37 @@
|
||||
import java.util.ArrayList;
|
||||
import java.util.Collections;
|
||||
import java.util.List;
|
||||
import java.util.stream.IntStream;
|
||||
|
||||
public class CellSample {
|
||||
|
||||
private List<Integer[]> cells;
|
||||
private Integer cdr1Freq;
|
||||
|
||||
public CellSample(Integer numDistinctCells, Integer cdr1Freq){
|
||||
this.cdr1Freq = cdr1Freq;
|
||||
List<Integer> numbersCDR3 = new ArrayList<>();
|
||||
List<Integer> numbersCDR1 = new ArrayList<>();
|
||||
Integer numDistCDR3s = 2 * numDistinctCells + 1;
|
||||
IntStream.range(1, numDistCDR3s + 1).forEach(i -> numbersCDR3.add(i));
|
||||
IntStream.range(numDistCDR3s + 1, numDistCDR3s + 1 + (numDistCDR3s / cdr1Freq) + 1).forEach(i -> numbersCDR1.add(i));
|
||||
Collections.shuffle(numbersCDR3);
|
||||
Collections.shuffle(numbersCDR1);
|
||||
|
||||
//Each cell represented by 4 values
|
||||
//two CDR3s, and two CDR1s. First two values are CDR3s (alpha, beta), second two are CDR1s (alpha, beta)
|
||||
List<Integer[]> distinctCells = new ArrayList<>();
|
||||
for(int i = 0; i < numbersCDR3.size() - 1; i = i + 2){
|
||||
Integer tmpCDR3a = numbersCDR3.get(i);
|
||||
Integer tmpCDR3b = numbersCDR3.get(i+1);
|
||||
Integer tmpCDR1a = numbersCDR1.get(i % numbersCDR1.size());
|
||||
Integer tmpCDR1b = numbersCDR1.get((i+1) % numbersCDR1.size());
|
||||
Integer[] tmp = {tmpCDR3a, tmpCDR3b, tmpCDR1a, tmpCDR1b};
|
||||
distinctCells.add(tmp);
|
||||
}
|
||||
this.cells = distinctCells;
|
||||
}
|
||||
|
||||
public CellSample(List<Integer[]> cells, Integer cdr1Freq){
|
||||
this.cells = cells;
|
||||
this.cdr1Freq = cdr1Freq;
|
||||
|
||||
@@ -288,7 +288,7 @@ public class CommandLineInterface {
|
||||
|
||||
//for calling from command line
|
||||
public static void makeCells(String filename, Integer numCells, Integer cdr1Freq){
|
||||
CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
|
||||
CellSample sample = new CellSample(numCells, cdr1Freq);
|
||||
CellFileWriter writer = new CellFileWriter(filename, sample);
|
||||
writer.writeCellsToFile();
|
||||
}
|
||||
|
||||
@@ -1,35 +0,0 @@
|
||||
import org.jgrapht.graph.SimpleWeightedGraph;
|
||||
import org.jgrapht.nio.graphml.GraphMLImporter;
|
||||
|
||||
import java.io.BufferedReader;
|
||||
import java.io.IOException;
|
||||
import java.nio.file.Files;
|
||||
import java.nio.file.Path;
|
||||
|
||||
public class GraphMLFileReader {
|
||||
|
||||
private String filename;
|
||||
private SimpleWeightedGraph graph;
|
||||
|
||||
public GraphMLFileReader(String filename, SimpleWeightedGraph graph) {
|
||||
if(!filename.matches(".*\\.graphml")){
|
||||
filename = filename + ".graphml";
|
||||
}
|
||||
this.filename = filename;
|
||||
this.graph = graph;
|
||||
|
||||
try(//don't need to close reader bc of try-with-resources auto-closing
|
||||
BufferedReader reader = Files.newBufferedReader(Path.of(filename));
|
||||
){
|
||||
GraphMLImporter<SimpleWeightedGraph, BufferedReader> importer = new GraphMLImporter<>();
|
||||
importer.importGraph(graph, reader);
|
||||
}
|
||||
catch (IOException ex) {
|
||||
System.out.println("Graph file " + filename + " not found.");
|
||||
System.err.println(ex);
|
||||
}
|
||||
}
|
||||
|
||||
public SimpleWeightedGraph getGraph() { return graph; }
|
||||
|
||||
}
|
||||
@@ -1,4 +1,8 @@
|
||||
import org.jgrapht.graph.DefaultWeightedEdge;
|
||||
import org.jgrapht.graph.SimpleWeightedGraph;
|
||||
import org.jgrapht.nio.Attribute;
|
||||
import org.jgrapht.nio.AttributeType;
|
||||
import org.jgrapht.nio.DefaultAttribute;
|
||||
import org.jgrapht.nio.dot.DOTExporter;
|
||||
import org.jgrapht.nio.graphml.GraphMLExporter;
|
||||
|
||||
@@ -7,25 +11,69 @@ import java.io.IOException;
|
||||
import java.nio.file.Files;
|
||||
import java.nio.file.Path;
|
||||
import java.nio.file.StandardOpenOption;
|
||||
import java.util.HashMap;
|
||||
import java.util.LinkedHashMap;
|
||||
import java.util.Map;
|
||||
|
||||
public class GraphMLFileWriter {
|
||||
|
||||
String filename;
|
||||
SimpleWeightedGraph graph;
|
||||
GraphWithMapData data;
|
||||
|
||||
|
||||
public GraphMLFileWriter(String filename, SimpleWeightedGraph graph) {
|
||||
public GraphMLFileWriter(String filename, GraphWithMapData data) {
|
||||
if(!filename.matches(".*\\.graphml")){
|
||||
filename = filename + ".graphml";
|
||||
}
|
||||
this.filename = filename;
|
||||
this.graph = graph;
|
||||
this.data = data;
|
||||
}
|
||||
|
||||
// public void writeGraphToFile() {
|
||||
// try(BufferedWriter writer = Files.newBufferedWriter(Path.of(filename), StandardOpenOption.CREATE_NEW);
|
||||
// ){
|
||||
// GraphMLExporter<SimpleWeightedGraph, BufferedWriter> exporter = new GraphMLExporter<>();
|
||||
// exporter.exportGraph(graph, writer);
|
||||
// } catch(IOException ex){
|
||||
// System.out.println("Could not make new file named "+filename);
|
||||
// System.err.println(ex);
|
||||
// }
|
||||
// }
|
||||
|
||||
public void writeGraphToFile() {
|
||||
SimpleWeightedGraph graph = data.getGraph();
|
||||
Map<Integer, Integer> vertexToAlphaMap = data.getPlateVtoAMap();
|
||||
Map<Integer, Integer> vertexToBetaMap = data.getPlateVtoBMap();
|
||||
Map<Integer, Integer> alphaOccs = data.getAlphaWellCounts();
|
||||
Map<Integer, Integer> betaOccs = data.getBetaWellCounts();
|
||||
try(BufferedWriter writer = Files.newBufferedWriter(Path.of(filename), StandardOpenOption.CREATE_NEW);
|
||||
){
|
||||
GraphMLExporter<SimpleWeightedGraph, BufferedWriter> exporter = new GraphMLExporter<>();
|
||||
//create exporter. Let the vertex labels be the unique ids for the vertices
|
||||
GraphMLExporter<Integer, SimpleWeightedGraph<Vertex, DefaultWeightedEdge>> exporter = new GraphMLExporter<>(v -> v.toString());
|
||||
//set to export weights
|
||||
exporter.setExportEdgeWeights(true);
|
||||
//set type, sequence, and occupancy attributes for each vertex
|
||||
exporter.setVertexAttributeProvider( v -> {
|
||||
Map<String, Attribute> attributes = new HashMap<>();
|
||||
if(vertexToAlphaMap.containsKey(v)) {
|
||||
attributes.put("type", DefaultAttribute.createAttribute("CDR3 Alpha"));
|
||||
attributes.put("sequence", DefaultAttribute.createAttribute(vertexToAlphaMap.get(v)));
|
||||
attributes.put("occupancy", DefaultAttribute.createAttribute(
|
||||
alphaOccs.get(vertexToAlphaMap.get(v))));
|
||||
}
|
||||
else if(vertexToBetaMap.containsKey(v)) {
|
||||
attributes.put("type", DefaultAttribute.createAttribute("CDR3 Beta"));
|
||||
attributes.put("sequence", DefaultAttribute.createAttribute(vertexToBetaMap.get(v)));
|
||||
attributes.put("occupancy", DefaultAttribute.createAttribute(
|
||||
betaOccs.get(vertexToBetaMap.get(v))));
|
||||
}
|
||||
return attributes;
|
||||
});
|
||||
//register the attributes
|
||||
exporter.registerAttribute("type", GraphMLExporter.AttributeCategory.NODE, AttributeType.STRING);
|
||||
exporter.registerAttribute("sequence", GraphMLExporter.AttributeCategory.NODE, AttributeType.STRING);
|
||||
exporter.registerAttribute("occupancy", GraphMLExporter.AttributeCategory.NODE, AttributeType.STRING);
|
||||
//export the graph
|
||||
exporter.exportGraph(graph, writer);
|
||||
} catch(IOException ex){
|
||||
System.out.println("Could not make new file named "+filename);
|
||||
@@ -33,3 +81,4 @@ public class GraphMLFileWriter {
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@@ -4,61 +4,75 @@ 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 {
|
||||
public interface GraphModificationFunctions {
|
||||
|
||||
//remove over- and under-weight edges
|
||||
public static List<Integer[]> filterByOverlapThresholds(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
|
||||
int low, int high) {
|
||||
static List<Integer[]> filterByOverlapThresholds(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
|
||||
int low, int high, boolean saveEdges) {
|
||||
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 (DefaultWeightedEdge e : graph.edgeSet()) {
|
||||
if ((graph.getEdgeWeight(e) > high) || (graph.getEdgeWeight(e) < low)) {
|
||||
if(saveEdges) {
|
||||
Integer source = graph.getEdgeSource(e);
|
||||
Integer target = graph.getEdgeTarget(e);
|
||||
Integer weight = (int) graph.getEdgeWeight(e);
|
||||
Integer[] edge = {source, target, weight};
|
||||
removedEdges.add(edge);
|
||||
}
|
||||
else {
|
||||
graph.setEdgeWeight(e, 0.0);
|
||||
}
|
||||
}
|
||||
}
|
||||
for (Integer[] edge : removedEdges) {
|
||||
graph.removeEdge(edge[0], edge[1]);
|
||||
if(saveEdges) {
|
||||
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,
|
||||
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) {
|
||||
Integer maxOccupancyDifference, boolean saveEdges) {
|
||||
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);
|
||||
if (saveEdges) {
|
||||
Integer source = graph.getEdgeSource(e);
|
||||
Integer target = graph.getEdgeTarget(e);
|
||||
Integer weight = (int) graph.getEdgeWeight(e);
|
||||
Integer[] edge = {source, target, weight};
|
||||
removedEdges.add(edge);
|
||||
}
|
||||
else {
|
||||
graph.setEdgeWeight(e, 0.0);
|
||||
}
|
||||
}
|
||||
}
|
||||
for (Integer[] edge : removedEdges) {
|
||||
graph.removeEdge(edge[0], edge[1]);
|
||||
if(saveEdges) {
|
||||
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,
|
||||
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) {
|
||||
Integer minOverlapPercent,
|
||||
boolean saveEdges) {
|
||||
List<Integer[]> removedEdges = new ArrayList<>();
|
||||
for (DefaultWeightedEdge e : graph.edgeSet()) {
|
||||
Integer alphaOcc = alphaWellCounts.get(plateVtoAMap.get(graph.getEdgeSource(e)));
|
||||
@@ -66,20 +80,27 @@ public abstract class GraphModificationFunctions {
|
||||
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);
|
||||
if(saveEdges) {
|
||||
Integer source = graph.getEdgeSource(e);
|
||||
Integer target = graph.getEdgeTarget(e);
|
||||
Integer intWeight = (int) graph.getEdgeWeight(e);
|
||||
Integer[] edge = {source, target, intWeight};
|
||||
removedEdges.add(edge);
|
||||
}
|
||||
else {
|
||||
graph.setEdgeWeight(e, 0.0);
|
||||
}
|
||||
}
|
||||
}
|
||||
for (Integer[] edge : removedEdges) {
|
||||
graph.removeEdge(edge[0], edge[1]);
|
||||
if(saveEdges) {
|
||||
for (Integer[] edge : removedEdges) {
|
||||
graph.removeEdge(edge[0], edge[1]);
|
||||
}
|
||||
}
|
||||
return removedEdges;
|
||||
}
|
||||
|
||||
public static void addRemovedEdges(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
|
||||
static void addRemovedEdges(SimpleWeightedGraph<Integer, DefaultWeightedEdge> graph,
|
||||
List<Integer[]> removedEdges) {
|
||||
for (Integer[] edge : removedEdges) {
|
||||
DefaultWeightedEdge e = graph.addEdge(edge[0], edge[1]);
|
||||
|
||||
@@ -27,6 +27,7 @@ public class InteractiveInterface {
|
||||
//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("8) Options");
|
||||
System.out.println("9) About/Acknowledgments");
|
||||
System.out.println("0) Exit");
|
||||
try {
|
||||
@@ -37,9 +38,10 @@ public class InteractiveInterface {
|
||||
case 3 -> makeCDR3Graph();
|
||||
case 4 -> matchCDR3s();
|
||||
//case 6 -> matchCellsCDR1();
|
||||
case 8 -> mainOptions();
|
||||
case 9 -> acknowledge();
|
||||
case 0 -> quit = true;
|
||||
default -> throw new InputMismatchException("Invalid input.");
|
||||
default -> System.out.println("Invalid input.");
|
||||
}
|
||||
} catch (InputMismatchException | IOException ex) {
|
||||
System.out.println(ex);
|
||||
@@ -72,17 +74,15 @@ public class InteractiveInterface {
|
||||
System.out.println(ex);
|
||||
sc.next();
|
||||
}
|
||||
CellSample sample = Simulator.generateCellSample(numCells, cdr1Freq);
|
||||
CellSample sample = new CellSample(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();
|
||||
if(BiGpairSEQ.cacheCells()) {
|
||||
BiGpairSEQ.setCellSampleInMemory(sample, filename);
|
||||
}
|
||||
BiGpairSEQ.setCellSampleInMemory(sample);
|
||||
BiGpairSEQ.setCellFilename(filename);
|
||||
}
|
||||
|
||||
//Output a CSV of sample plate
|
||||
@@ -219,9 +219,9 @@ public class InteractiveInterface {
|
||||
System.out.println("Reading Cell Sample file: " + cellFile);
|
||||
CellFileReader cellReader = new CellFileReader(cellFile);
|
||||
cells = cellReader.getCellSample();
|
||||
BiGpairSEQ.clearCellSampleInMemory();
|
||||
BiGpairSEQ.setCellSampleInMemory(cells);
|
||||
BiGpairSEQ.setCellFilename(cellFile);
|
||||
if(BiGpairSEQ.cacheCells()) {
|
||||
BiGpairSEQ.setCellSampleInMemory(cells, cellFile);
|
||||
}
|
||||
}
|
||||
assert filename != null;
|
||||
Plate samplePlate;
|
||||
@@ -242,8 +242,9 @@ public class InteractiveInterface {
|
||||
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);
|
||||
if(BiGpairSEQ.cachePlate()) {
|
||||
BiGpairSEQ.setPlateInMemory(samplePlate, filename);
|
||||
}
|
||||
}
|
||||
|
||||
//Output serialized binary of GraphAndMapData object
|
||||
@@ -251,7 +252,6 @@ public class InteractiveInterface {
|
||||
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.");
|
||||
@@ -271,16 +271,16 @@ public class InteractiveInterface {
|
||||
assert cellFile != null;
|
||||
CellSample cellSample;
|
||||
//check if cells are already in memory
|
||||
if(cellFile.equals(BiGpairSEQ.getCellFilename())) {
|
||||
if(cellFile.equals(BiGpairSEQ.getCellFilename()) && BiGpairSEQ.getCellSampleInMemory() != null) {
|
||||
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);
|
||||
if(BiGpairSEQ.cacheCells()) {
|
||||
BiGpairSEQ.setCellSampleInMemory(cellSample, cellFile);
|
||||
}
|
||||
}
|
||||
|
||||
assert plateFile != null;
|
||||
@@ -290,12 +290,12 @@ public class InteractiveInterface {
|
||||
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(BiGpairSEQ.cachePlate()) {
|
||||
BiGpairSEQ.setPlateInMemory(plate, plateFile);
|
||||
}
|
||||
}
|
||||
if (cellSample.getCells().size() == 0){
|
||||
System.out.println("No cell sample found.");
|
||||
@@ -309,12 +309,20 @@ public class InteractiveInterface {
|
||||
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.");
|
||||
if(BiGpairSEQ.outputBinary()) {
|
||||
GraphDataObjectWriter dataWriter = new GraphDataObjectWriter(filename, data);
|
||||
dataWriter.writeDataToFile();
|
||||
System.out.println("Serialized binary graph/data file written to: " + filename);
|
||||
}
|
||||
if(BiGpairSEQ.outputGraphML()) {
|
||||
GraphMLFileWriter graphMLWriter = new GraphMLFileWriter(filename, data);
|
||||
graphMLWriter.writeGraphToFile();
|
||||
System.out.println("GraphML file written to: " + filename);
|
||||
}
|
||||
if(BiGpairSEQ.cacheGraph()) {
|
||||
BiGpairSEQ.setGraphInMemory(data, filename);
|
||||
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
@@ -366,17 +374,15 @@ public class InteractiveInterface {
|
||||
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);
|
||||
if(graphFilename.equals(BiGpairSEQ.getGraphFilename())) {
|
||||
data = BiGpairSEQ.getGraphInMemory();
|
||||
}
|
||||
else {
|
||||
data = BiGpairSEQ.getGraphInMemory();
|
||||
GraphDataObjectReader dataReader = new GraphDataObjectReader(graphFilename);
|
||||
data = dataReader.getData();
|
||||
if(BiGpairSEQ.cacheGraph()) {
|
||||
BiGpairSEQ.setGraphInMemory(data, graphFilename);
|
||||
}
|
||||
}
|
||||
//simulate matching
|
||||
MatchingResult results = Simulator.matchCDR3s(data, graphFilename, lowThreshold, highThreshold, maxOccupancyDiff,
|
||||
@@ -493,6 +499,79 @@ public class InteractiveInterface {
|
||||
// }
|
||||
// }
|
||||
|
||||
private static void mainOptions(){
|
||||
boolean backToMain = false;
|
||||
while(!backToMain) {
|
||||
System.out.println("\n--------------OPTIONS---------------");
|
||||
System.out.println("1) Turn " + getOnOff(!BiGpairSEQ.cacheCells()) + " cell sample file caching");
|
||||
System.out.println("2) Turn " + getOnOff(!BiGpairSEQ.cachePlate()) + " plate file caching");
|
||||
System.out.println("3) Turn " + getOnOff(!BiGpairSEQ.cacheGraph()) + " graph/data file caching");
|
||||
System.out.println("4) Turn " + getOnOff(!BiGpairSEQ.outputBinary()) + " serialized binary graph output");
|
||||
System.out.println("5) Turn " + getOnOff(!BiGpairSEQ.outputGraphML()) + " GraphML graph output");
|
||||
System.out.println("6) Maximum weight matching algorithm options");
|
||||
System.out.println("0) Return to main menu");
|
||||
try {
|
||||
input = sc.nextInt();
|
||||
switch (input) {
|
||||
case 1 -> BiGpairSEQ.setCacheCells(!BiGpairSEQ.cacheCells());
|
||||
case 2 -> BiGpairSEQ.setCachePlate(!BiGpairSEQ.cachePlate());
|
||||
case 3 -> BiGpairSEQ.setCacheGraph(!BiGpairSEQ.cacheGraph());
|
||||
case 4 -> BiGpairSEQ.setOutputBinary(!BiGpairSEQ.outputBinary());
|
||||
case 5 -> BiGpairSEQ.setOutputGraphML(!BiGpairSEQ.outputGraphML());
|
||||
case 6 -> algorithmOptions();
|
||||
case 0 -> backToMain = true;
|
||||
default -> System.out.println("Invalid input");
|
||||
}
|
||||
} catch (InputMismatchException ex) {
|
||||
System.out.println(ex);
|
||||
sc.next();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Helper function for printing menu items in mainOptions(). Returns a string based on the value of parameter.
|
||||
*
|
||||
* @param b - a boolean value
|
||||
* @return String "on" if b is true, "off" if b is false
|
||||
*/
|
||||
private static String getOnOff(boolean b) {
|
||||
if (b) { return "on";}
|
||||
else { return "off"; }
|
||||
}
|
||||
|
||||
private static void algorithmOptions(){
|
||||
boolean backToOptions = false;
|
||||
while(!backToOptions) {
|
||||
System.out.println("\n---------ALGORITHM OPTIONS----------");
|
||||
System.out.println("1) Use scaling algorithm by Duan and Su.");
|
||||
System.out.println("2) Use LEDA book algorithm with Fibonacci heap priority queue");
|
||||
System.out.println("3) Use LEDA book algorithm with pairing heap priority queue");
|
||||
System.out.println("0) Return to Options menu");
|
||||
try {
|
||||
input = sc.nextInt();
|
||||
switch (input) {
|
||||
case 1 -> System.out.println("This option is not yet implemented. Choose another.");
|
||||
case 2 -> {
|
||||
BiGpairSEQ.setFibonacciHeap();
|
||||
System.out.println("MWM algorithm set to LEDA with Fibonacci heap");
|
||||
backToOptions = true;
|
||||
}
|
||||
case 3 -> {
|
||||
BiGpairSEQ.setPairingHeap();
|
||||
System.out.println("MWM algorithm set to LEDA with pairing heap");
|
||||
backToOptions = true;
|
||||
}
|
||||
case 0 -> backToOptions = true;
|
||||
default -> System.out.println("Invalid input");
|
||||
}
|
||||
} catch (InputMismatchException ex) {
|
||||
System.out.println(ex);
|
||||
sc.next();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
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.");
|
||||
|
||||
@@ -3,6 +3,7 @@ import org.jgrapht.alg.matching.MaximumWeightBipartiteMatching;
|
||||
import org.jgrapht.generate.SimpleWeightedBipartiteGraphMatrixGenerator;
|
||||
import org.jgrapht.graph.DefaultWeightedEdge;
|
||||
import org.jgrapht.graph.SimpleWeightedGraph;
|
||||
import org.jheaps.tree.FibonacciHeap;
|
||||
import org.jheaps.tree.PairingHeap;
|
||||
|
||||
import java.math.BigDecimal;
|
||||
@@ -16,36 +17,12 @@ 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 {
|
||||
public class Simulator implements GraphModificationFunctions {
|
||||
private static final int cdr3AlphaIndex = 0;
|
||||
private static final int cdr3BetaIndex = 1;
|
||||
private static final int cdr1AlphaIndex = 2;
|
||||
private static final int cdr1BetaIndex = 3;
|
||||
|
||||
public static CellSample generateCellSample(Integer numDistinctCells, Integer cdr1Freq) {
|
||||
//In real T cells, CDR1s have about one third the diversity of CDR3s
|
||||
List<Integer> numbersCDR3 = new ArrayList<>();
|
||||
List<Integer> numbersCDR1 = new ArrayList<>();
|
||||
Integer numDistCDR3s = 2 * numDistinctCells + 1;
|
||||
IntStream.range(1, numDistCDR3s + 1).forEach(i -> numbersCDR3.add(i));
|
||||
IntStream.range(numDistCDR3s + 1, numDistCDR3s + 1 + (numDistCDR3s / cdr1Freq) + 1).forEach(i -> numbersCDR1.add(i));
|
||||
Collections.shuffle(numbersCDR3);
|
||||
Collections.shuffle(numbersCDR1);
|
||||
|
||||
//Each cell represented by 4 values
|
||||
//two CDR3s, and two CDR1s. First two values are CDR3s (alpha, beta), second two are CDR1s (alpha, beta)
|
||||
List<Integer[]> distinctCells = new ArrayList<>();
|
||||
for(int i = 0; i < numbersCDR3.size() - 1; i = i + 2){
|
||||
Integer tmpCDR3a = numbersCDR3.get(i);
|
||||
Integer tmpCDR3b = numbersCDR3.get(i+1);
|
||||
Integer tmpCDR1a = numbersCDR1.get(i % numbersCDR1.size());
|
||||
Integer tmpCDR1b = numbersCDR1.get((i+1) % numbersCDR1.size());
|
||||
Integer[] tmp = {tmpCDR3a, tmpCDR3b, tmpCDR1a, tmpCDR1b};
|
||||
distinctCells.add(tmp);
|
||||
}
|
||||
return new CellSample(distinctCells, cdr1Freq);
|
||||
}
|
||||
|
||||
//Make the graph needed for matching CDR3s
|
||||
public static GraphWithMapData makeGraph(List<Integer[]> distinctCells, Plate samplePlate, boolean verbose) {
|
||||
Instant start = Instant.now();
|
||||
@@ -146,8 +123,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<>();
|
||||
boolean saveEdges = BiGpairSEQ.cacheGraph();
|
||||
int numWells = data.getNumWells();
|
||||
Integer alphaCount = data.getAlphaCount();
|
||||
Integer betaCount = data.getBetaCount();
|
||||
@@ -160,33 +137,50 @@ public class Simulator {
|
||||
|
||||
//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");}
|
||||
removedEdges.addAll(GraphModificationFunctions.filterByOverlapThresholds(graph, lowThreshold, highThreshold));
|
||||
removedEdges.addAll(GraphModificationFunctions.filterByOverlapThresholds(graph, lowThreshold, highThreshold, saveEdges));
|
||||
if(verbose){System.out.println("Over- and under-weight edges removed");}
|
||||
|
||||
//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.");}
|
||||
removedEdges.addAll(GraphModificationFunctions.filterByOverlapPercent(graph, alphaWellCounts, betaWellCounts,
|
||||
plateVtoAMap, plateVtoBMap, minOverlapPercent));
|
||||
plateVtoAMap, plateVtoBMap, minOverlapPercent, saveEdges));
|
||||
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);}
|
||||
removedEdges.addAll(GraphModificationFunctions.filterByRelativeOccupancy(graph, alphaWellCounts, betaWellCounts,
|
||||
plateVtoAMap, plateVtoBMap, maxOccupancyDifference));
|
||||
plateVtoAMap, plateVtoBMap, maxOccupancyDifference, saveEdges));
|
||||
if(verbose){System.out.println("Edges between vertices of with excessively different occupancy values " +
|
||||
"removed");}
|
||||
|
||||
//Find Maximum Weighted Matching
|
||||
//using jheaps library class PairingHeap for improved efficiency
|
||||
if(verbose){System.out.println("Finding maximum weighted matching");}
|
||||
//Attempting to use addressable heap to improve performance
|
||||
MaximumWeightBipartiteMatching maxWeightMatching =
|
||||
new MaximumWeightBipartiteMatching(graph,
|
||||
MaximumWeightBipartiteMatching maxWeightMatching;
|
||||
//Use correct heap type for priority queue
|
||||
String heapType = BiGpairSEQ.getPriorityQueueHeapType();
|
||||
switch (heapType) {
|
||||
case "PAIRING" -> {
|
||||
maxWeightMatching = new MaximumWeightBipartiteMatching(graph,
|
||||
plateVtoAMap.keySet(),
|
||||
plateVtoBMap.keySet(),
|
||||
i -> new PairingHeap(Comparator.naturalOrder()));
|
||||
}
|
||||
case "FIBONACCI" -> {
|
||||
maxWeightMatching = new MaximumWeightBipartiteMatching(graph,
|
||||
plateVtoAMap.keySet(),
|
||||
plateVtoBMap.keySet(),
|
||||
i -> new FibonacciHeap(Comparator.naturalOrder()));
|
||||
}
|
||||
default -> {
|
||||
maxWeightMatching = new MaximumWeightBipartiteMatching(graph,
|
||||
plateVtoAMap.keySet(),
|
||||
plateVtoBMap.keySet());
|
||||
}
|
||||
}
|
||||
//get the matching
|
||||
MatchingAlgorithm.Matching<String, DefaultWeightedEdge> graphMatching = maxWeightMatching.getMatching();
|
||||
if(verbose){System.out.println("Matching completed");}
|
||||
Instant stop = Instant.now();
|
||||
@@ -242,6 +236,7 @@ public class Simulator {
|
||||
}
|
||||
|
||||
//Metadata comments for CSV file
|
||||
String algoType = "LEDA book with heap: " + heapType;
|
||||
int min = Math.min(alphaCount, betaCount);
|
||||
//rate of attempted matching
|
||||
double attemptRate = (double) (trueCount + falseCount) / min;
|
||||
@@ -272,6 +267,7 @@ public class Simulator {
|
||||
Map<String, String> metadata = new LinkedHashMap<>();
|
||||
metadata.put("sample plate filename", data.getSourceFilename());
|
||||
metadata.put("graph filename", dataFilename);
|
||||
metadata.put("algorithm type", algoType);
|
||||
metadata.put("well populations", wellPopulationsString);
|
||||
metadata.put("total alphas found", alphaCount.toString());
|
||||
metadata.put("total betas found", betaCount.toString());
|
||||
@@ -292,10 +288,11 @@ public class Simulator {
|
||||
}
|
||||
}
|
||||
|
||||
//put the removed edges back on the graph
|
||||
System.out.println("Restoring removed edges to graph.");
|
||||
GraphModificationFunctions.addRemovedEdges(graph, removedEdges);
|
||||
|
||||
if(saveEdges) {
|
||||
//put the removed edges back on the graph
|
||||
System.out.println("Restoring removed edges to graph.");
|
||||
GraphModificationFunctions.addRemovedEdges(graph, removedEdges);
|
||||
}
|
||||
//return MatchingResult object
|
||||
return output;
|
||||
}
|
||||
@@ -671,7 +668,7 @@ public class Simulator {
|
||||
|
||||
private static Map<Integer, Integer> makeVertexToSequenceMap(Map<Integer, Integer> sequences, Integer startValue) {
|
||||
Map<Integer, Integer> map = new LinkedHashMap<>(); //LinkedHashMap to preserve order of entry
|
||||
Integer index = startValue;
|
||||
Integer index = startValue; //is this necessary? I don't think I use this.
|
||||
for (Integer k: sequences.keySet()) {
|
||||
map.put(index, k);
|
||||
index++;
|
||||
|
||||
@@ -1,14 +1,20 @@
|
||||
|
||||
|
||||
public class Vertex {
|
||||
private final Integer peptide;
|
||||
private final Integer vertexLabel;
|
||||
private final Integer sequence;
|
||||
private final Integer occupancy;
|
||||
|
||||
public Vertex(Integer peptide, Integer occupancy) {
|
||||
this.peptide = peptide;
|
||||
public Vertex(Integer vertexLabel, Integer sequence, Integer occupancy) {
|
||||
this.vertexLabel = vertexLabel;
|
||||
this.sequence = sequence;
|
||||
this.occupancy = occupancy;
|
||||
}
|
||||
|
||||
public Integer getPeptide() {
|
||||
return peptide;
|
||||
public Integer getVertexLabel() { return vertexLabel; }
|
||||
|
||||
public Integer getSequence() {
|
||||
return sequence;
|
||||
}
|
||||
|
||||
public Integer getOccupancy() {
|
||||
|
||||
Reference in New Issue
Block a user