bfs counts always 0

[rpauly]

Hi!
All the counts from the bfs are NaN when normalized and 0 when not normalized. Below is the code I ran:

`> bfs <- c("./<path_to_bam>/1.bam") %>%
+     map(bin_fragment_size, bin_size = 5, normalized = TRUE) %>%
+     purrr::reduce(inner_join, by = "Breaks")

> tail (bfs, n=12)
    Breaks UNAFFECTED
69 341_346        NaN
70 346_351        NaN
71 351_356        NaN
72 356_361        NaN
73 361_366        NaN
74 366_371        NaN
75 371_376        NaN
76 376_381        NaN
77 381_386        NaN
78 386_391        NaN
79 391_396        NaN
80 396_400        NaN`

However, I know there are fragments as seen below

`fs1 <- get_fragment_size(bam = "/<path_to_bam>/1.bam", 
+                          mapqFilter = 30, 
+                          isProperPair = NA, 
+                          min_size = 1, 
+                          max_size = 400, 
+                          ignore_trimmed = FALSE, 
+                          simple_cigar = FALSE, 
+                          different_strands = TRUE)
> head(fs1)`
`      Sample                                                                 ID  chr    start      end size
1 UNAFFECTED   UNAFFECTED_A00450:173:HMNH2DSX2:3:2441:6560:29543:GAAGTG+AAGTCCA chr7 55242170 55242488  319
2 UNAFFECTED UNAFFECTED_A00450:173:HMNH2DSX2:4:1264:16206:35759:TAGCGTA+AGTGGTA chr7 55242186 55242365  180
3 UNAFFECTED  UNAFFECTED_A00450:173:HMNH2DSX2:3:2240:13901:14763:GTAACA+TAATGCG chr7 55242190 55242539  350
4 UNAFFECTED  UNAFFECTED_A00450:173:HMNH2DSX2:3:1204:4372:13823:GAACCTC+GCTGTCA chr7 55242197 55242511  315
5 UNAFFECTED  UNAFFECTED_A00450:173:HMNH2DSX2:3:1153:14389:34710:CGGTGTA+GTAACA chr7 55242203 55242561  359
6 UNAFFECTED   UNAFFECTED_A00450:173:HMNH2DSX2:3:2523:4453:13777:CCTCATA+AATGAG chr7 55242203 55242566  364`

Is it something I am not doing?
Thanks!
Rini

[alkodsi]

Hi @rpauly
The normalized counts are expected to be NaN if the not normalized counts are all 0.
As to why the counts are 0, I can't tell without debugging with some of your data.
You can try to run this internal function (https://github.com/alkodsi/ctDNAtools/blob/master/R/get_hist_bins.R) on the output of get_fragment_size() and see what happens.

[rpauly]

Thanks! Please find attached the the output of get_fragment_size().

FS1_get_fragment_size.txt