@@ -44,6 +44,7 @@ pub struct CellStats {
4444 pub fraction : f32 ,
4545}
4646
47+ #[ derive( Default ) ]
4748pub struct DisplayStrings {
4849 /// The display for the profile portion of the alignment
4950 pub profile_string : String ,
@@ -55,6 +56,40 @@ pub struct DisplayStrings {
5556 pub posterior_string : String ,
5657}
5758
59+ #[ derive( Default ) ]
60+ pub struct AlignmentStats {
61+ // The length of the alignment
62+ pub length : usize ,
63+ // The number of exact matches
64+ pub n_match : usize ,
65+ // The number of mismatches/substitutions
66+ pub n_mismatch : usize ,
67+ // The number of gap openings
68+ pub n_gap_open : usize ,
69+ // The total number of gaps
70+ pub n_gap : usize ,
71+ }
72+
73+ impl AlignmentStats {
74+ pub fn pid_with_gaps ( & self ) -> anyhow:: Result < f32 > {
75+ if self . length == 0 {
76+ bail ! ( "alignment length 0 in AlignmentStats::pid_with_gaps()" )
77+ }
78+
79+ Ok ( self . n_match as f32 / self . length as f32 )
80+ }
81+
82+ pub fn pid_no_gaps ( & self ) -> anyhow:: Result < f32 > {
83+ if self . length == 0 {
84+ bail ! ( "alignment length 0 in AlignmentStats::pid_no_gaps()" )
85+ } else if self . n_gap > self . length {
86+ bail ! ( "n_gap > length in AlignmentStats::pid_no_gaps()" )
87+ }
88+
89+ Ok ( self . n_match as f32 / ( self . length - self . n_gap ) as f32 )
90+ }
91+ }
92+
5893pub struct Alignment {
5994 /// The name of the profile/model
6095 pub profile_name : Option < String > ,
@@ -68,6 +103,8 @@ pub struct Alignment {
68103 pub cell_stats : Option < CellStats > ,
69104 /// The strings used for alignment display
70105 pub display_strings : Option < DisplayStrings > ,
106+ ///
107+ pub stats : Option < AlignmentStats > ,
71108}
72109
73110impl AsRef < Alignment > for & Alignment {
@@ -87,7 +124,7 @@ impl AsRef<Alignment> for &Alignment {
87124/// 0.85 - 0.95 -> "9"
88125///
89126/// 0.95 - 1.00 -> "*"
90- fn map_posterior_probability_to_bin_byte ( probability : f32 ) -> u8 {
127+ fn map_posterior_probability_to_bin_byte_utf8 ( probability : f32 ) -> u8 {
91128 let bin = ( probability * 10.0 ) . round ( ) ;
92129 UTF8_NUMERIC [ bin as usize ]
93130}
@@ -115,6 +152,8 @@ impl ScoreParams {
115152pub struct AlignmentBuilder < ' a > {
116153 target : Option < & ' a Sequence > ,
117154 profile : Option < & ' a Profile > ,
155+ profile_name : Option < String > ,
156+ target_name : Option < String > ,
118157 trace : Option < & ' a Trace > ,
119158 database_size : Option < usize > ,
120159 forward_score : Option < Bits > ,
@@ -133,11 +172,21 @@ impl<'a> AlignmentBuilder<'a> {
133172 self
134173 }
135174
175+ pub fn with_profile_name ( mut self , name : & str ) -> Self {
176+ self . profile_name = Some ( name. to_string ( ) ) ;
177+ self
178+ }
179+
136180 pub fn with_target ( mut self , target : & ' a Sequence ) -> Self {
137181 self . target = Some ( target) ;
138182 self
139183 }
140184
185+ pub fn with_target_name ( mut self , name : & str ) -> Self {
186+ self . target_name = Some ( name. to_string ( ) ) ;
187+ self
188+ }
189+
141190 pub fn with_database_size ( mut self , count : usize ) -> Self {
142191 self . database_size = Some ( count) ;
143192 self
@@ -225,83 +274,101 @@ impl<'a> AlignmentBuilder<'a> {
225274 None => None ,
226275 } ;
227276
228- let display_strings = match ( self . trace , self . profile , self . target ) {
277+ let ( display_strings, stats ) = match ( self . trace , self . profile , self . target ) {
229278 ( Some ( trace) , Some ( profile) , Some ( target) ) => {
230- let mut profile_bytes = vec ! [ ] ;
231- let mut target_bytes = vec ! [ ] ;
232- let mut middle_bytes = vec ! [ ] ;
233- let mut posteriors = vec ! [ ] ;
279+ let mut stats = AlignmentStats :: default ( ) ;
280+ stats . length = trace . core_len ( ) ;
281+
282+ let mut display = DisplayStrings :: default ( ) ;
234283
235- trace
284+ let mut prev_state = Trace :: INVALID_STATE ;
285+ for step in trace
236286 . iter ( )
237- . filter ( |s| {
238- s. state == Trace :: M_STATE
239- || s. state == Trace :: I_STATE
240- || s. state == Trace :: D_STATE
241- } )
242- . for_each ( |step| {
243- posteriors. push ( map_posterior_probability_to_bin_byte ( step. posterior ) ) ;
244- let profile_byte = profile. consensus_seq_bytes_utf8 [ step. profile_idx ] ;
245- let target_byte = target. utf8_bytes [ step. target_idx ] ;
246-
247- match step. state {
248- Trace :: I_STATE => {
249- profile_bytes. push ( Alignment :: PROFILE_GAP_BYTE ) ;
250- target_bytes. push ( target_byte) ;
251- middle_bytes. push ( UTF8_SPACE ) ;
252- }
253- Trace :: D_STATE => {
254- profile_bytes. push ( profile_byte) ;
255- target_bytes. push ( Alignment :: TARGET_GAP_BYTE ) ;
256- middle_bytes. push ( UTF8_SPACE ) ;
257- }
258- Trace :: M_STATE => {
259- let target_byte_digital = target. digital_bytes [ step. target_idx ] ;
260-
261- profile_bytes. push ( profile_byte) ;
262- target_bytes. push ( target_byte) ;
263-
264- if profile_byte. to_ascii_lowercase ( )
265- == target_byte. to_ascii_lowercase ( )
266- {
267- middle_bytes. push ( profile_byte) ;
268- } else if profile
269- . match_score ( target_byte_digital as usize , step. profile_idx )
270- > 0.0
271- {
272- middle_bytes. push ( UTF8_PLUS ) ;
287+ . skip_while ( |s| !s. is_core_state ( ) )
288+ . take_while ( |s| s. is_core_state ( ) )
289+ {
290+ display
291+ . posterior_string
292+ . push ( map_posterior_probability_to_bin_byte_utf8 ( step. posterior ) as char ) ;
293+
294+ let prf_char = profile. consensus_seq_bytes_utf8 [ step. profile_idx ] as char ;
295+ let seq_char = target. utf8_bytes [ step. target_idx ] as char ;
296+
297+ match step. state {
298+ Trace :: M_STATE => {
299+ display. profile_string . push ( prf_char) ;
300+ display. target_string . push ( seq_char) ;
301+
302+ let matches_consensus =
303+ prf_char. to_ascii_lowercase ( ) == seq_char. to_ascii_lowercase ( ) ;
304+
305+ let match_score = profile. match_score (
306+ target. digital_bytes [ step. target_idx ] as usize ,
307+ step. profile_idx ,
308+ ) ;
309+
310+ let mid_char = if matches_consensus {
311+ stats. n_match += 1 ;
312+ prf_char
313+ } else {
314+ stats. n_mismatch += 1 ;
315+ if match_score > 0.0 {
316+ UTF8_PLUS as char
273317 } else {
274- middle_bytes . push ( UTF8_SPACE ) ;
318+ UTF8_SPACE as char
275319 }
320+ } ;
321+
322+ display. middle_string . push ( mid_char) ;
323+ }
324+ Trace :: I_STATE => {
325+ stats. n_gap += 1 ;
326+ if step. state != prev_state {
327+ stats. n_gap_open += 1 ;
276328 }
277- _ => {
278- panic ! ( "invalid trace state in AlignmentBuilder: {}" , step. state)
329+
330+ // ---
331+
332+ display
333+ . profile_string
334+ . push ( Alignment :: PROFILE_GAP_BYTE as char ) ;
335+ display. target_string . push ( seq_char) ;
336+ display. middle_string . push ( UTF8_SPACE as char ) ;
337+ }
338+ Trace :: D_STATE => {
339+ stats. n_gap += 1 ;
340+ if step. state != prev_state {
341+ stats. n_gap_open += 1 ;
279342 }
343+
344+ // ---
345+
346+ display. profile_string . push ( prf_char) ;
347+ display
348+ . target_string
349+ . push ( Alignment :: TARGET_GAP_BYTE as char ) ;
350+ display. middle_string . push ( UTF8_SPACE as char ) ;
280351 }
281- } ) ;
282-
283- let profile_string = String :: from_utf8 ( profile_bytes) ?;
284- let target_string = String :: from_utf8 ( target_bytes) ?;
285- let middle_string = String :: from_utf8 ( middle_bytes) ?;
286- let posterior_string = String :: from_utf8 ( posteriors) ?;
287-
288- Some ( DisplayStrings {
289- profile_string,
290- target_string,
291- middle_string,
292- posterior_string,
293- } )
352+
353+ _ => bail ! ( "unexpected non-core state in filtered trace step iter" ) ,
354+ }
355+
356+ prev_state = step. state ;
357+ }
358+
359+ ( Some ( display) , Some ( stats) )
294360 }
295- _ => None ,
361+ _ => ( None , None ) ,
296362 } ;
297363
298364 Ok ( Alignment {
299- profile_name : self . profile . map ( |profile| profile . name . clone ( ) ) ,
300- target_name : self . target . map ( |target| target . name . clone ( ) ) ,
365+ profile_name : self . profile_name ,
366+ target_name : self . target_name ,
301367 boundaries,
302368 scores,
303369 cell_stats,
304370 display_strings,
371+ stats,
305372 } )
306373 }
307374}
0 commit comments