Hi,

I have bam and bai files generated from my custom analysis pipeline, and I wanted to perform fragmentomics analysis on it using cfDNAPro. I am bit confused between 'readBam' vs 'read_bam_insert_metrics' commands. As mentioned in the documentation, both these commands input files is in 'bam' format. The 'readBam' output is 'curated GRanges object', while 'read_bam_insert_metrics' output is 'tibble'.

The confusion is:

readBam: "Read bam file into a curated GRanges object"

read_bam_insert_metrics: "Calculate insert sizes from a curated GRanges object"

Is the output of readBam is the input for read_bam_insert_metrics ?

Thanks

Rakesh

Depends on what you mean by input. If you mean to ask if you use readBam to generate something to use as an argument to read_bam_insert_metrics, then no. It's easy enough for you to answer that question yourself by reading the help files for both functions (which if you are planning to use R for much, you need to get into the habit of doing). Both functions take a path to a BAM as their first argument.

If you mean to ask if read_bam_insert_metrics uses readBam internally, then yes. You can tell by reading the code:

read_bam_insert_metrics <- function(bamfile = NULL,
                                    fragment_obj = NULL,
                                    chromosome_to_keep =paste0("chr", 1:22),
                                    strand_mode = 1,
                                    genome_label = "hg19",
                                    outdir = FALSE,
                                    isize_min = 1L,
                                    isize_max = 1000L,
                                    ...) {

  insert_size <- NULL

  # users are required to supply only one of 'bamfile' or 'fragment_obj" from
  # readBam() function
  switch(
    rlang::check_exclusive(bamfile, fragment_obj),
    bamfile = message("`bamfile` was supplied."),
    fragment_obj = message("`fragment_obj` was supplied.")
  )

  if(!is.null(bamfile)) {
    frag <- readBam(bamfile = bamfile, 
                  chromosome_to_keep = chromosome_to_keep,
                  strand_mode = strand_mode,
                  genome_label = genome_label,
                  outdir = outdir)
<snip>