compute_fragments_length {FunChIP} | R Documentation

## Computing the length of the fragments in the .bam file.

### Description

Given a .bam file and GRanges object, it computes the positive and negative coverage for each GRanges element, estimates the distance between positive and negative peaks, and finally the fragment length `d`, i.e. the sum of the length of the reads and the distance between positive and negative peaks. See Details and the package vignette for the description of the method.

### Usage

    compute_fragments_length(object, bamf, min.d = 0, max.d = 200)

### Arguments

`object` | GRanges object of length N.
`bamf` | Path to the .bam file used to compute the coverage function. The associated .bam.bai index file must also be present.
`min.d` | integer. Minimum value for the distance between positive and negative peaks. Default is 0.
`max.d` | integer. Maximum value for the distance between positive and negative peaks. Default is 200.

### Details

Given a set of n= 1, … N regions, characterized by their positive and negative coverages, the function computes the distance between the positive peak f_{n+} and the negative peak, shifted by δ f_{n-}^{δ}:

D(f_{n+}, f_{n-}^{δ}) = \frac{\| f_{n+} - f_{n-}^{δ}\|^2_{L^2}}{\textrm{width}(\textrm{union}(f_{n+}, f_{n-}^{δ}))}

The function computes the d_{pn} minimizing the distance between postive and negative peaks

d_{pn} = \textrm{argmin}_{δ \in [\code{min.d}, \code{max.d}]} ∑_{n=1 }^N D(f_{n+}, f_{n-}^{δ})

The function returns both the plot of the global distance vs the fragment length d = d_{pn} + r, where r is the length of the reads, and the optimum value for d. r is also estimated from the .bam file as the average of the read lengths.

### Value

optimum value of the parameter d, to be used in the pileup_peak method.

### Author(s)

Alice Parodi, Marco J. Morelli, Laura M. Sangalli, Piercesare Secchi, Simone Vantini

### Examples

    # load the data
    # GRanges object

    data(GR100)

    # import the .bam file

    bamf <- system.file("extdata", "test.bam", package="FunChIP",
                          mustWork=TRUE)

    # compute the estimated fragment length

    d <- compute_fragments_length(GR[1:10], bamf, min.d = 0, max.d = 200)

[Package FunChIP version 1.13.0 Index]