Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 5 additions & 15 deletions src/align/score.rs
Original file line number Diff line number Diff line change
Expand Up @@ -126,19 +126,12 @@ impl AlignmentScorer {
((genomic_span as f64).log2() * self.score_genomic_length_log2_scale - 0.5).ceil() as i32
}

/// Apply annotation bonus to junction score
///
/// # Arguments
/// * `base_score` - Base score from motif
/// * `annotated` - Whether the junction is annotated in GTF
///
/// # Returns
/// Adjusted score with annotation bonus applied
pub fn score_annotated_junction(&self, base_score: i32, annotated: bool) -> i32 {
/// Annotated junctions score `sjdb_score`; unannotated junctions score `motif_score`.
pub fn score_annotated_junction(&self, motif_score: i32, annotated: bool) -> i32 {
if annotated {
base_score + self.sjdb_score
self.sjdb_score
} else {
base_score
motif_score
}
}

Expand Down Expand Up @@ -959,17 +952,14 @@ mod tests {
out_filter_score_min_over_lread: 0.66,
};

// Annotated junction should get bonus
let annotated_score = scorer.score_annotated_junction(0, true);
assert_eq!(annotated_score, 2);

// Novel junction should not get bonus
let novel_score = scorer.score_annotated_junction(0, false);
assert_eq!(novel_score, 0);

// Bonus applies to any base score
let annotated_noncanon = scorer.score_annotated_junction(-8, true);
assert_eq!(annotated_noncanon, -6); // -8 + 2
assert_eq!(annotated_noncanon, 2);
}

#[test]
Expand Down
38 changes: 37 additions & 1 deletion src/align/stitch.rs
Original file line number Diff line number Diff line change
Expand Up @@ -1313,9 +1313,10 @@ fn stitch_align_to_transcript(
|| db.is_annotated(cluster.chr_idx, donor_fwd, acceptor_fwd, 2)
});

d_score += motif_score;
if is_annotated {
d_score += scorer.sjdb_score;
} else {
d_score += motif_score;
}

new_wt.n_junction += 1;
Expand Down Expand Up @@ -3430,4 +3431,39 @@ mod tests {
bin_start_edge = bin_start_edge.saturating_sub(win_flank_nbins);
assert_eq!(bin_start_edge, 0, "Flanking should saturate at 0");
}

#[test]
fn test_junction_score_annotated_uses_sjdb_not_motif() {
use crate::align::score::AlignmentScorer;

let scorer = AlignmentScorer::from_params_minimal();

for motif_score in [0_i32, scorer.score_gap_gcag, scorer.score_gap_atac] {
let baseline: i32 = 100;

let mut d_score_annot = baseline;
let is_annotated = true;
if is_annotated {
d_score_annot += scorer.sjdb_score;
} else {
d_score_annot += motif_score;
}
assert_eq!(d_score_annot, baseline + scorer.sjdb_score);

let mut d_score_novel = baseline;
let is_annotated = false;
if is_annotated {
d_score_novel += scorer.sjdb_score;
} else {
d_score_novel += motif_score;
}
assert_eq!(d_score_novel, baseline + motif_score);
}

let baseline: i32 = 100;
let motif_score = scorer.score_gap_gcag;
let additive_buggy = baseline + motif_score + scorer.sjdb_score;
let replacement_correct = baseline + scorer.sjdb_score;
assert_ne!(additive_buggy, replacement_correct);
}
}
Loading