debruijn/colors.rs
1use std::{collections::HashMap, fmt::Display, iter::Sum, marker::PhantomData};
2
3use log::debug;
4
5use crate::{graph::DebruijnGraph, summarizer::{SummaryConfig, SummaryData, ID, Marker}, Kmer};
6use std::fmt::Debug;
7
8/// mode for coloring nodes in dot files - check compatibility with node data
9#[derive(Debug, Clone, PartialEq)]
10pub enum ColorMode<'a> {
11 /// only compatible with [`IDSumData`]
12 IDGroups {id_group_ids: &'a HashMap<ID, ID>, n_id_groups: usize},
13 /// only compatible with [`IDSumData`]
14 IDS {n_ids: usize},
15 /// compatible with all [`SummaryData`] containing `Tags`
16 SampleGroups,
17 /// compatible with all [`SummaryData`] containing `TagsCounts`
18 FoldChange
19}
20
21impl ColorMode<'_> {
22 /// get the hash map connecting IDs to their group IDs
23 pub fn id_group_ids(&self) -> Option<&HashMap<ID, ID>> {
24 match self {
25 Self::IDGroups { id_group_ids, n_id_groups: _ } => Some(id_group_ids),
26 _ => None
27 }
28 }
29}
30
31
32/// contains the hues, the markers signifying which tag belongs to which group,
33/// the maximun kmer count and the average kmer count
34#[derive(Clone, Debug, PartialEq)]
35pub struct Colors<'a, SD: SummaryData<DI>, DI> {
36 color_mode: ColorMode<'a>,
37 // 2-bit encoded group associations of labels
38 marker0: Marker,
39 marker1: Marker,
40 // factor (slope) for log2(fold change) to hue transformation
41 log2_fc_factor: Option<f32>,
42 // slope (m) and y intercept (b) for n obs to value transformation
43 _log2_nobs_mb: Option<(f32, f32)>,
44 // slope (m) and y intercept (b) for log10(p-value)) to saturation transformation
45 log10_p_mb: Option<(f32, f32)>,
46 // slope (m) and y intercept (b) for log10(edge multiplicity) to pen width transformation
47 log10_em_mb: Option<(f32, f32)>,
48 // phantom data who
49 phantom_data_sd: PhantomData<SD>,
50 phantom_data_di: PhantomData<DI>,
51}
52
53// add `'b, 'a: 'b` in case of lifetime error
54impl<'a, SD: SummaryData<DI> + Debug, DI> Colors<'a, SD, DI> {
55 const HUE_RED: f32 = 0.;
56 const HUE_YELLOW: f32 = 60. / 360.;
57 const HUE_GREEN: f32 = 120. / 360.;
58 const HUE_PURPLE: f32 = 289. / 360.;
59
60 const SAT_MIN: f32 = 0.5;
61 const SAT_MAX: f32 = 1.;
62 const SAT_DEF: f32 = 1.;
63
64 const VAL_MIN: f32 = 0.7;
65 const VAL_MAX: f32 = 1.;
66 const VAL_DEF: f32 = 1.;
67
68 const FC_MAX: f32 = 5.;
69 const FC_MIN: f32 = -5.;
70 const P_MAX: f32 = -4.;
71 const SIGN_P: f32 = 0.05;
72
73 const EDGE_WIDTH_MAX: f32 = 20.;
74 const EDGE_WIDTH_MIN: f32 = 3.;
75 const EDGE_WIDTH_DEF: f32 = 8.;
76
77
78 /// Creates a new [`Colors<SD>`].
79 pub fn new<'b: 'a, K: Kmer>(graph: &DebruijnGraph<K, SD>, summary_config: &SummaryConfig, color_mode: ColorMode<'b>) -> Self {
80 let (log2_fc_factor,
81 _log2_nobs_mb,
82 log10_p_mb,
83 log10_em_mb) = match color_mode {
84 ColorMode::FoldChange => {
85 // fold change factor: fold change of node will be multiplied by this factor to get hue
86 // factor is the slope of a linear function
87 let log2_fc_factor = match graph.get_node(0).data().p_value(summary_config) {
88 Some(_) => {
89 let (min_max, _, _) = get_min_max(
90 &graph,
91 &|&graph| Box::new(graph
92 .iter_nodes()
93 .map(|node| node.data().fold_change(summary_config).expect("error getting fold change"))
94 )
95 );
96 debug!("log2_fc min max: {:?}", min_max);
97
98 match min_max {
99 Some((min_fc, max_fc)) => {
100 // make symmetrical
101 let m_fc = if min_fc.abs() > max_fc.abs() {
102 min_fc.abs()
103 } else {
104 max_fc.abs()
105 };
106
107 // if too large, replace with max fold change
108 let val_fc = if m_fc > Self::FC_MAX {
109 Self::FC_MAX
110 } else {
111 m_fc
112 };
113
114 // yellow should be where log2(fc) = 0
115 // division by zero should not happen bc max has to be larger than min
116 // and larger absolute of the two is used
117 Some(Self::HUE_YELLOW / val_fc)
118 },
119 None => None
120 }
121 },
122 None => None
123 };
124
125 debug!("log2_fc factor: {:?}", log2_fc_factor);
126
127 // number of observations (nobs)
128 let (nobs, _, _) = get_min_max(
129 &graph,
130 &|&graph| Box::new(graph
131 .iter_nodes()
132 .map(|node| node.data().sum().unwrap_or(1) as f32)
133 )
134 );
135 debug!("nobs min max: {:?}", nobs);
136
137 // calculate m and b for value = m * log10(nobs) + b
138 let log2_nobs_mb = match nobs {
139 Some((min_nobs, max_nobs)) => {
140 let min = min_nobs.log2();
141 let max = max_nobs.log2();
142
143 let m = (Self::VAL_MAX - Self::VAL_MIN) / (max - min);
144 let b = Self::VAL_MIN - m * min;
145
146 Some((m, b))
147 }
148 None => None
149 };
150
151 debug!("log2_nobs m, b: {:?}", log2_nobs_mb);
152
153 // calculate m and b for saturation = m * log10(p-value) + b
154 let log10_p_mb = match graph.get_node(0).data().p_value(summary_config) {
155 Some(_) => {
156 // get min and max p-value
157 let (min_max, _, _) = get_min_max(
158 &graph,
159 &|&graph| Box::new(graph
160 .iter_nodes()
161 .map(|node| node.data().p_value(summary_config).expect("error getting p-value"))
162 )
163 );
164 debug!("p min max: {:?}", min_max);
165 match min_max {
166 Some((min_p, max_p)) => {
167 let min = max_p.log10();
168 let max = min_p.log10();
169
170 // replace "max" of p with P_MAX if too small
171 let max = if max < Self::P_MAX { Self::P_MAX } else { max };
172
173 // b should be 0 -> p = 1 -> log10(1) = 0 -> VAL_MIN=0
174 let m = (Self::SAT_MAX - Self::SAT_MIN) / (max - min);
175 let b = Self::SAT_MIN - m * min;
176
177 Some((m, b))
178 }
179 None => None
180 }
181 },
182 None => None
183 };
184
185 debug!("log10_p m, b: {:?}", log10_p_mb);
186
187 // calculate m and b for pen width = m * log10(edge mults) + b
188 let log10_em_mb = match graph.get_node(0).data().edge_mults() {
189 Some(_) => {
190 // get min and max
191 let (min_max, _, _) = get_min_max(&graph, &|graph| Box::new(graph
192 .iter_nodes()
193 .flat_map(|node| node.data().edge_mults().expect("error getting edge mult").edge_mults())
194 .filter(|element| *element != 0)
195 ));
196
197 debug!("edge mults min max: {:?}", min_max);
198
199 min_max.map(|(min, max)|{
200 let min = (min as f32).log10();
201 let max = (max as f32).log10();
202
203 // b should be 1 -> em = 1 -> log10(1) = 0 -> EDGE_WITH_MIN=1
204 let m = (Self::EDGE_WIDTH_MAX - Self::EDGE_WIDTH_MIN) / (max - min);
205 let b = Self::EDGE_WIDTH_MIN - m * min;
206
207 (m, b)
208 })
209 },
210 None => None
211 };
212
213 debug!("log10_em m, b: {:?}", log10_em_mb);
214
215 (log2_fc_factor, log2_nobs_mb, log10_p_mb, log10_em_mb)
216 },
217 _ => (None, None, None, None)
218 };
219
220 let (marker0, marker1) = summary_config.get_markers();
221
222 Colors {
223 color_mode,
224 marker0,
225 marker1,
226 log2_fc_factor,
227 _log2_nobs_mb,
228 log10_p_mb,
229 log10_em_mb,
230 phantom_data_sd: PhantomData,
231 phantom_data_di: PhantomData
232 }
233 }
234
235 /// calculate the saturation
236 fn saturation(&self, data: &SD, summary_config: &SummaryConfig) -> f32 {
237 match self.log10_p_mb {
238 Some((_m, _b)) => if data.p_value(summary_config).expect("error getting p-value") < Self::SIGN_P { Self::SAT_MAX } else { Self::SAT_MIN },
239 None => Self::SAT_DEF
240 }
241 }
242
243 /// get the color(s) for a nore in hsv format for dot
244 /// -> single "hue saturation value" or list
245 /// "hue saturation value:hue saturation value:hue saturation value"
246 /// for multiple colors
247 fn hsv_dot(&self, data: &SD, summary_config: &SummaryConfig, saturation: f32, value: f32) -> String {
248 match self.color_mode {
249 ColorMode::FoldChange => {
250 let hue = match self.log2_fc_factor {
251 // if fold change available calculate hue based on log2(fc)
252 Some(fc_factor) => {
253 match data.fold_change(summary_config).unwrap() {
254 Self::FC_MAX..=f32::INFINITY => Self::HUE_GREEN,
255 f32::NEG_INFINITY..=Self::FC_MIN => Self::HUE_RED,
256 fc => fc * fc_factor + Self::HUE_YELLOW
257 }
258 },
259 None => Self::HUE_PURPLE
260 };
261 format!("{hue} {saturation} {value}")
262 },
263 ColorMode::SampleGroups => {
264 let hue = match data.tags() {
265 Some(tag) => {
266 if tag.bit_and(self.marker0) & !tag.bit_and(self.marker1) {
267 // tags are only in marker0 group
268 Self::HUE_GREEN
269 } else if !tag.bit_and(self.marker0) & tag.bit_and(self.marker1) {
270 // tag is only in marker1 group
271 Self::HUE_RED
272 } else if tag.bit_and(self.marker0) & tag.bit_and(self.marker1) {
273 // tag is in both groups
274 Self::HUE_YELLOW
275 } else {
276 // tag is in neither group
277 // should only happen if the tags started with more than three distinct characters
278 // (overflow purple)
279 Self::HUE_PURPLE
280 }
281 },
282 None => Self::HUE_PURPLE
283 };
284 format!("{hue} {saturation} {value}")
285 },
286 ColorMode::IDGroups { id_group_ids, n_id_groups } => {
287 match data.ids() {
288 Some(ids) => {
289 format!("{:?}",
290 ids.iter().map(|id| format!("{} {saturation} {value}",
291 *(id_group_ids.get(id).expect("id was not in HM")) as f32 / n_id_groups as f32
292 )).collect::<Vec<_>>()
293 ).replace("\"", "").replace("[", "").replace("]", "").replace(", ", ":")
294 // turns '["hue saturation value", "hue saturation value", hue saturation value"]' into 'hue saturation value:hue saturation value:hue saturation value' -> DOT format
295 },
296 // if we don't have IDs, check if we have mapped IDs
297 None => match data.mapped_ids() {
298 Some(ids) => {
299 if ids.is_empty() {
300 // if we have a map summarizer, but no mapped ids, the vector is empty
301 // usually this means that the node is "false" -> use default purple and lower sv
302 format!("{} {} {}", Self::HUE_PURPLE, saturation/2., value/2.)
303 } else {
304 // we have mapped IDs, use the same procedure as IDs
305 format!("{:?}",
306 ids.iter().map(|id| format!("{} {saturation} {value}",
307 *(id_group_ids.get(id).expect("id was not in HM")) as f32 / n_id_groups as f32
308 )).collect::<Vec<_>>()
309 ).replace("\"", "").replace("[", "").replace("]", "").replace(", ", ":")
310 }
311 }
312 None => format!("{} {saturation} {value}", Self::HUE_PURPLE)
313 }
314
315
316 }
317 }
318 ColorMode::IDS { n_ids } => {
319 match data.ids() {
320 Some(ids) => {
321 format!("{:?}",
322 ids.iter().map(|id| format!("{} {saturation} {value}",
323 *id as f32 / n_ids as f32
324 )).collect::<Vec<_>>()
325 ).replace("\"", "").replace("[", "").replace("]", "").replace(", ", ":")
326 // turns '["hue saturation value", "hue saturation value", hue saturation value"]' into 'hue saturation value:hue saturation value:hue saturation value' -> DOT format
327 }
328 // if we don't have IDs, check if we have mapped IDs
329 None => match data.mapped_ids() {
330 Some(ids) => {
331 if ids.is_empty() {
332 // if we have a map summarizer, but no mapped ids, the vector is empty
333 // usually this means that the node is "false" -> use default purple and lower sv (probaby gray)
334 format!("{} {} {}", Self::HUE_PURPLE, saturation/2., value/2.)
335 } else {
336 // we have mapped IDs, use the same procedure as IDs
337 format!("{:?}",
338 ids.iter().map(|id| format!("{} {saturation} {value}",
339 *id as f32 / n_ids as f32
340 )).collect::<Vec<_>>()
341 ).replace("\"", "").replace("[", "").replace("]", "").replace(", ", ":")
342 }
343 }
344 None => format!("{} {saturation} {value}", Self::HUE_PURPLE)
345 }
346 }
347 }
348 }
349 }
350
351 /// get the color for a node in a HSV format
352 pub fn node_color_dot(&self, data: &SD, summary_config: &SummaryConfig, outline: bool) -> String {
353 // calculate saturation
354 let saturation = self.saturation(data, summary_config);
355
356 // set value as default value
357 let value = Self::VAL_DEF;
358
359 // get hue and color
360 let hsv_colors = self.hsv_dot(data, summary_config, saturation, value);
361
362 // adapt font color to value (currently always black)
363 let font_color= if value <= 0.5 { "white" } else { "black" };
364
365 // set outline (eg if it is in a path)
366 let mut prefix = if outline {
367 "black, penwidth=10, fillcolor=".to_string()
368 } else {
369 "".to_string()
370 };
371
372 // set shape and style
373 let shape_style = match self.color_mode {
374 ColorMode::FoldChange | ColorMode::SampleGroups => "style=filled", // default to oval shape
375 ColorMode::IDGroups { id_group_ids: _, n_id_groups: _ } | ColorMode::IDS { n_ids: _ } => "shape=rectangle, style=striped" // striped only works in rectangle
376 };
377
378 // overwrite generated path with mapped paths
379 if let Some(t_ids) = data.mapped_ids() {
380 if !t_ids.is_empty() {
381 match self.color_mode {
382 ColorMode::IDS { n_ids } => {
383 let outline_hue = t_ids.iter().map(|id| *id as f32 / (n_ids * t_ids.len()) as f32).sum::<f32>();
384 prefix = format!("\"{outline_hue} 1 0.6\", penwidth=30, fillcolor=");
385 }
386 ColorMode::IDGroups { id_group_ids, n_id_groups } => {
387 let outline_hue = t_ids.iter().map(|id| *(id_group_ids.get(id).expect("id was not in HM")) as f32 / (n_id_groups * t_ids.len()) as f32).sum::<f32>();
388 prefix = format!("\"{outline_hue} 1 0.6\", penwidth=30, fillcolor=");
389 }
390 _ => ()
391 }
392 }
393 }
394
395 // return formatted string for color, fillcolor, and fontcolor
396 format!("{shape_style}, color={prefix}\"{hsv_colors}\", fontcolor={font_color}")
397 }
398
399 /// get the hue for hsl color in json
400 pub fn hue_json(&self, data: &SD, summary_config: &SummaryConfig) -> i32 {
401 let hue = match self.color_mode {
402 ColorMode::FoldChange => {
403 match self.log2_fc_factor {
404 // if fold change available calculate hue based on log2(fc)
405 Some(fc_factor) => {
406 match data.fold_change(summary_config).unwrap() {
407 Self::FC_MAX..=f32::INFINITY => Self::HUE_GREEN,
408 f32::NEG_INFINITY..=Self::FC_MIN => Self::HUE_RED,
409 fc => fc * fc_factor + Self::HUE_YELLOW
410 }
411 },
412 None => Self::HUE_PURPLE
413 }
414 },
415 ColorMode::SampleGroups => {
416 match data.tags() {
417 Some(tag) => {
418 if tag.bit_and(self.marker0) & !tag.bit_and(self.marker1) {
419 // tags are only in marker0 group
420 Self::HUE_GREEN
421 } else if !tag.bit_and(self.marker0) & tag.bit_and(self.marker1) {
422 // tag is only in marker1 group
423 Self::HUE_RED
424 } else if tag.bit_and(self.marker0) & tag.bit_and(self.marker1) {
425 // tag is in both groups
426 Self::HUE_YELLOW
427 } else {
428 // tag is in neither group
429 // should only happen if the tags started with more than three distinct characters
430 // (overflow purple)
431 Self::HUE_PURPLE
432 }
433 },
434 None => Self::HUE_PURPLE
435 }
436 },
437 ColorMode::IDGroups { id_group_ids, n_id_groups } => {
438 match data.ids() {
439 Some(ids) => {
440 // calculate separate hues of group IDs
441 let id_hues = ids
442 .iter()
443 .map(|id| *(id_group_ids
444 .get(id)
445 .expect("id was not in HM")) as f32 / n_id_groups as f32)
446 .collect::<Vec<_>>();
447
448 // calculate average hue
449 id_hues.iter().sum::<f32>() / id_hues.len() as f32
450 },
451 None => Self::HUE_PURPLE
452 }
453 }
454 ColorMode::IDS { n_ids } => {
455 match data.ids() {
456 Some(ids) => {
457 // calculate separate hues of group IDs
458 let id_hues = ids.iter().map(|id| *id as f32 / n_ids as f32).collect::<Vec<_>>();
459
460 // calculate the average hue
461 id_hues.iter().sum::<f32>() / id_hues.len() as f32
462 }
463 None => Self::HUE_PURPLE
464 }
465 }
466 };
467
468 // transform into degrees
469 (hue * 365.) as i32
470 }
471
472 /// get the edge width based on the edge multiplicity
473 pub fn edge_width(&self, edge_mult: u32) -> f32 {
474 match self.log10_em_mb {
475 Some((m, b)) => (edge_mult as f32).log10() * m + b,
476 None => Self::EDGE_WIDTH_DEF
477 }
478 }
479
480 pub fn id_group_ids(&self) -> Option<&HashMap<ID, ID>> {
481 self.color_mode.id_group_ids()
482 }
483
484}
485
486/// get the minimum and maximum value of an iterator, and if they are three times as
487/// small/large as the next smallest/largest item, also return the latter
488/// filters for -inf and inf with floats
489pub fn get_min_max<I: Debug, F, II, N>(iter_struct: &I, iter_value: &F) -> (Option<(N, N)>, Option<N>, Option<N>)
490where
491 F: Fn(&I) -> Box<II>,
492 II: Iterator<Item = N>,
493 N: CFilter
494{
495 let min_o = iter_value(iter_struct)
496 .filter(|value| value.filter())
497 .min_by(|a, b| a.partial_cmp(b).unwrap());
498
499 // check if there are actually elements to process
500 if let Some(min) = min_o {
501 let max = iter_value(iter_struct)
502 .filter(|value| value.filter())
503 .max_by(|a, b| a.partial_cmp(b).unwrap())
504 .expect("error: empty iterator");
505
506 if max > min {
507 let snd_max = iter_value(iter_struct)
508 .filter(|value| *value != max && value.filter())
509 .max_by(|a, b| a.partial_cmp(b).unwrap())
510 .expect("error: empty iterator");
511 let snd_min = iter_value(iter_struct)
512 .filter(|value| *value != min && value.filter())
513 .min_by(|a, b| a.partial_cmp(b).unwrap())
514 .expect("error: empty iterator");
515
516 // if max is OUTL times bigger than snd_max, it is an outlier
517 const OUTL: f64 = 3.;
518 let outlier_max = if (max.to_f64() / snd_max.to_f64()) > OUTL && snd_max > snd_min { Some(snd_max) } else { None };
519 let outlier_min = if (min.to_f64() / snd_min.to_f64()) > OUTL && snd_max > snd_min { Some(snd_min) } else { None };
520
521 return (Some((min, max)), outlier_min, outlier_max)
522 }
523 }
524
525 (None, None, None)
526}
527
528/// trait for filtering the values in [`get_min_max`]
529pub trait CFilter: Sum + PartialOrd + Copy + Display {
530 fn filter(self) -> bool;
531 fn to_f64(self) -> f64;
532}
533
534impl CFilter for usize {
535 fn filter(self) -> bool {
536 self != usize::MAX
537 }
538 fn to_f64(self) -> f64 {
539 self as f64
540 }
541}
542
543impl CFilter for u32 {
544 fn filter(self) -> bool {
545 self != u32::MAX
546 }
547 fn to_f64(self) -> f64 {
548 self as f64
549 }
550}
551
552impl CFilter for f32 {
553 fn filter(self) -> bool {
554 self.is_finite()
555 }
556 fn to_f64(self) -> f64 {
557 self as f64
558 }
559}
560