@@ -355,8 +355,6 @@ public void processFilesRemote(List<SequenceOutputFile> inputFiles, JobContext c
355355 {
356356 Integer gatkDepth = null ;
357357 Integer lofreqDepth = null ;
358- Double totalAltAf = 0.0 ;
359- int totalAltDepth = 0 ;
360358 Map <String , Double > alleleToAf = new HashMap <>();
361359 Map <String , Integer > alleleToDp = new HashMap <>();
362360
@@ -414,9 +412,6 @@ public void processFilesRemote(List<SequenceOutputFile> inputFiles, JobContext c
414412 String a = siteDef .getRenamedAllele (vc , vc .getAlternateAlleles ().get (0 ));
415413 if (a != null )
416414 {
417- totalAltAf += af ;
418- totalAltDepth += alleleDepth ;
419-
420415 double val = alleleToAf .getOrDefault (a , 0.0 );
421416 if (val == NO_DATA_VAL )
422417 {
@@ -441,7 +436,40 @@ public void processFilesRemote(List<SequenceOutputFile> inputFiles, JobContext c
441436 List <String > toWrite = new ArrayList <>(line );
442437 toWrite .add (String .valueOf (gatkDepth ));
443438 toWrite .add (String .valueOf (lofreqDepth ));
444- toWrite .add (String .valueOf (1 - totalAltAf ));
439+
440+ //NOTE: because an upstream deletion could be merged with this site (and their depths / AFs computed slightly differently), recalculate here:
441+
442+ //First establish true ALT depth, which could include upstream deletions:
443+ int totalAltDepth = 0 ;
444+ double totalAltAf = 0.0 ;
445+ for (String a : siteDef ._alternates )
446+ {
447+ double af = alleleToAf .getOrDefault (a , 0.0 );
448+ int dp = alleleToDp .getOrDefault (a , 0 );
449+
450+ if (dp == (int )NO_DATA_VAL )
451+ {
452+ continue ;
453+ }
454+
455+ totalAltDepth += dp ;
456+
457+ double adjAF = (double )dp / lofreqDepth ;
458+ if (Math .abs (adjAF - af ) > 0.01 )
459+ {
460+ ctx .getLogger ().error ("Significant AF adjustment: " + line .get (0 ) + "/" + line .get (4 ) + "/" + a + "/" + af + "/" + adjAF );
461+ }
462+
463+ totalAltAf += adjAF ;
464+ alleleToAf .put (a , adjAF );
465+ }
466+
467+ double refAF = 1 - totalAltAf ;
468+ toWrite .add (String .valueOf (refAF ));
469+ if (refAF < 0 )
470+ {
471+ ctx .getLogger ().error ("Negative REF AF: " + line .get (0 ) + "/" + line .get (4 ) + "/" + refAF );
472+ }
445473
446474 //Add AFs in order:
447475 List <Object > toAdd = new ArrayList <>();
0 commit comments