Skip to content

Commit

Permalink
fixes
Browse files Browse the repository at this point in the history
  • Loading branch information
@fcunial
fcunial committed Aug 15, 2023
1 parent df0e8b6 commit 07fe605
Showing 1 changed file with 102 additions and 20 deletions.
122 changes: 102 additions & 20 deletions src/de/mpi_cbg/revant/apps/BuildAssemblyGraph.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,17 +49,20 @@ public static void main(String[] args) throws IOException {
final int IDENTITY_THRESHOLD = IO.quantum;
final int MIN_COMPONENT_SIZE = 10; // Arbitrary

boolean keep;
boolean keep, hasPrefix, hasSuffix, hasPrefixPrime, hasSuffixPrime, hasLoop;
int i, j, k;
int type, nAlignments, idGenerator, top, node, neighbor, nComponents, size, nextFullyUnique;
int type, nAlignments, idGenerator, top, node, neighbor, nNeighbors, nComponents, size, nextFullyUnique;
int prefixOrSuffixA, prefixOrSuffixB, max, nSingletons, nTips, nLoops;
double errorRate;
String str1, str2, str3;
RepeatAlphabet.Character tmpChar;
Edge tmpEdge;
BufferedReader br1, br2, br3;
BufferedWriter bw;
boolean[] containsUnique;
int[] component, componentSize, stack;
BufferedWriter[] bws;
Edge[] edges;

// Building the graph
System.err.println("Building the graph...");
Expand Down Expand Up @@ -89,6 +92,12 @@ public static void main(String[] args) throws IOException {
continue;
}
keep=str2.equalsIgnoreCase("1")&&str3.equalsIgnoreCase("1");
if (type==0) {
prefixOrSuffixA=Alignments.startA<=IDENTITY_THRESHOLD?0:1;
prefixOrSuffixB=Alignments.startB<=IDENTITY_THRESHOLD?0:1;
}
else { prefixOrSuffixA=2; prefixOrSuffixB=2; }




Expand All @@ -103,22 +112,67 @@ public static void main(String[] args) throws IOException {



addEdge(Alignments.readA-1,Alignments.readB-1,keep);
addEdge(Alignments.readB-1,Alignments.readA-1,keep);
addEdge(Alignments.readA-1,Alignments.readB-1,keep,prefixOrSuffixA);
addEdge(Alignments.readB-1,Alignments.readA-1,keep,prefixOrSuffixB);
str1=br1.readLine(); str2=br2.readLine(); str3=br3.readLine();
if (nAlignments%10000==0) System.err.println("Processed "+nAlignments+" alignments");
}
br1.close(); br2.close(); br3.close();

// Removing duplicated edges
System.err.println("Removing duplicated edges...");
max=0;
for (i=0; i<N_READS; i++) {
if (lastNeighbor[i]>max) max=lastNeighbor[i];
}
max=(max+1)>>1;
edges = new Edge[max];
for (i=0; i<max; i++) edges[i] = new Edge();
for (i=0; i<N_READS; i++) {
if (lastNeighbor[i]<=0) continue;
Arrays.sort(neighbors[i],0,lastNeighbor[i]+1);
if (lastNeighbor[i]<=1) continue;
nNeighbors=(lastNeighbor[i]+1)>>1;
for (j=0; j<=lastNeighbor[i]; j+=2) {
k=j>>1;
edges[k].neighbor=neighbors[i][j];
edges[k].prefixOrSuffix=neighbors[i][j+1];
}
Arrays.sort(edges,0,nNeighbors);
k=0;
for (j=1; j<=lastNeighbor[i]; j++) {
if (neighbors[i][j]!=neighbors[i][k]) neighbors[i][++k]=neighbors[i][j];
for (j=1; j<nNeighbors; j++) {
if (edges[j].neighbor!=edges[k].neighbor || edges[j].prefixOrSuffix!=edges[k].prefixOrSuffix) {
k++; tmpEdge=edges[k]; edges[k]=edges[j]; edges[j]=tmpEdge;
}
}
lastNeighbor[i]=k;
lastNeighbor[i]=-1;
for (j=0; j<=k; j++) {
neighbors[i][++lastNeighbor[i]]=edges[j].neighbor;
neighbors[i][++lastNeighbor[i]]=edges[j].prefixOrSuffix;
}
}

// Computing number of nodes with only-prefix or only-suffix alignments
nSingletons=0; nTips=0; nLoops=0;
for (i=0; i<N_READS; i++) {
if (lastNeighbor[i]==-1) { nSingletons++; continue; }
hasPrefix=false; hasSuffix=false; hasLoop=false;
hasPrefixPrime=false; hasSuffixPrime=false;
for (j=0; j<=lastNeighbor[i]; j+=2) {
if (neighbors[i][j+1]==0) hasPrefix=true;
else if (neighbors[i][j+1]==1) hasSuffix=true;
if (j>0 && neighbors[i][j]!=neighbors[i][j-2]) {
if (hasPrefixPrime && hasSuffixPrime) hasLoop=true;
hasPrefixPrime=false; hasSuffixPrime=false;
}
if (neighbors[i][j+1]==0) hasPrefixPrime=true;
else if (neighbors[i][j+1]==1) hasSuffixPrime=true;
}
if (hasPrefixPrime && hasSuffixPrime) hasLoop=true;
if (hasPrefix!=hasSuffix) nTips++;
if (hasLoop) nLoops++;
}
System.err.println("Reads with no alignment: "+nSingletons+" ("+((100.0*nSingletons)/N_READS)+"%)");
System.err.println("Reads with only prefix or only suffix alignments: "+nTips+" ("+((100.0*nTips)/N_READS)+"%)");
System.err.println("Reads with suffix and prefix alignments to the same other read: "+nLoops+" ("+((100.0*nLoops)/N_READS)+"%)");

// Building $containsUnique$.
RepeatAlphabet.deserializeAlphabet(ALPHABET_FILE,2);
Expand Down Expand Up @@ -148,7 +202,7 @@ public static void main(String[] args) throws IOException {
top=0; stack[0]=i; size=1;
while (top>=0) {
node=stack[top--];
for (j=0; j<=lastNeighbor[node]; j++) {
for (j=0; j<=lastNeighbor[node]; j+=2) {
neighbor=neighbors[node][j];
if (neighbor<0) neighbor=-1-neighbor;
if (component[neighbor]!=-1) continue;
Expand Down Expand Up @@ -190,7 +244,7 @@ public static void main(String[] args) throws IOException {
nextFullyUnique=str1!=null?Integer.parseInt(str1):Math.POSITIVE_INFINITY;
}
else bws[k].write(i+" [uniqueStatus=\""+(containsUnique[i]?1:0)+"\"];\n");
for (j=0; j<=lastNeighbor[i]; j++) {
for (j=0; j<=lastNeighbor[i]; j+=2) {
neighbor=neighbors[i][j]>=0?neighbors[i][j]:-1-neighbors[i][j];
if (i<neighbor) bws[k].write(i+" -- "+neighbor+";\n");
}
Expand All @@ -199,7 +253,7 @@ public static void main(String[] args) throws IOException {
for (i=0; i<nComponents; i++) { bws[i].write("}"); bws[i].close(); }
br1.close();

// Outputting filtered connected components
// Printing filtered connected components
System.err.println("Computing filtered connected components...");
Math.set(component,N_READS-1,-1);
idGenerator=-1;
Expand All @@ -210,7 +264,7 @@ public static void main(String[] args) throws IOException {
top=0; stack[0]=i; size=1;
while (top>=0) {
node=stack[top--];
for (j=0; j<=lastNeighbor[node]; j++) {
for (j=0; j<=lastNeighbor[node]; j+=2) {
neighbor=neighbors[node][j];
if (neighbor<0 || component[neighbor]!=-1) continue;
component[neighbor]=idGenerator;
Expand Down Expand Up @@ -251,7 +305,7 @@ public static void main(String[] args) throws IOException {
nextFullyUnique=str1!=null?Integer.parseInt(str1):Math.POSITIVE_INFINITY;
}
else bws[k].write(i+" [uniqueStatus=\""+(containsUnique[i]?1:0)+"\"];\n");
for (j=0; j<=lastNeighbor[i]; j++) {
for (j=0; j<=lastNeighbor[i]; j+=2) {
if (neighbors[i][j]>=0 && i<neighbors[i][j]) bws[k].write(i+" -- "+neighbors[i][j]+";\n");
}
}
Expand All @@ -263,18 +317,46 @@ public static void main(String[] args) throws IOException {

/**
* Adds $to$ to $from$.
*
* @param prefixOrSuffixFrom the alignment uses a prefix (0), suffix (1), or neither
* (2), of $from$.
*/
private static final void addEdge(int from, int to, boolean keep) {
final int CAPACITY = 2; // Arbitrary
private static final void addEdge(int from, int to, boolean keep, int prefixOrSuffixFrom) {
final int CAPACITY = 4; // Arbitrary

lastNeighbor[from]++;
if (neighbors[from]==null || neighbors[from].length==0) neighbors[from] = new int[CAPACITY];
else if (lastNeighbor[from]==neighbors[from].length) {
else if (lastNeighbor[from]+2>=neighbors[from].length) {
int[] newArray = new int[neighbors[from].length<<1];
System.arraycopy(neighbors[from],0,newArray,0,neighbors[from].length);
System.arraycopy(neighbors[from],0,newArray,0,lastNeighbor[from]+1);
neighbors[from]=newArray;
}
neighbors[from][lastNeighbor[from]]=keep?to:-1-to;
neighbors[from][++lastNeighbor[from]]=keep?to:-1-to;
neighbors[from][++lastNeighbor[from]]=prefixOrSuffixFrom;
}


private static class Edge implements Comparable {
public int neighbor, prefixOrSuffix;

public Edge() { }

public void set(int n, int p) {
this.neighbor=n; this.prefixOrSuffix=p;
}

public int compareTo(Object other) {
Edge otherEdge = (Edge)other;
if (neighbor<otherEdge.neighbor) return -1;
else if (neighbor>otherEdge.neighbor) return 1;
if (prefixOrSuffix<otherEdge.prefixOrSuffix) return -1;
else if (prefixOrSuffix>otherEdge.prefixOrSuffix) return 1;
return 0;
}

public boolean equals(Object other) {
Edge otherEdge = (Edge)other;
return neighbor==otherEdge.neighbor && prefixOrSuffix==otherEdge.prefixOrSuffix;
}
}

}

0 comments on commit 07fe605

Please sign in to comment.