Find Fragment Size of a bamfilee

Question

[GK1610]

I have a test.bam file which has paired end reads data

I want to find the fragment size of each pair to compare the results from macs2 fragment size determination.

[GenoMax]

Starting with your original data you can calculate the insert sizes by a couple of different ways using BBMap. BBMap is available here.

Answer 1

[Devon Ryan]

You can get the distribution with bamPEFragmentSize from deepTools.

Answer 2

[apnri]

Fragment size calculation is usually determined by cross-correlation analysis, wherein you compute correlation between the positive and negative strand densities from a particular position. Then you shift the position of both strands by one position - you should expect to see maximum correlation between the two strands at fragment length.

You can use spp (used by ENCODE consortium) package to do this in R.

Snippet of code:

library(spp)    
chip.data <- read.bam.tags(bamfile)
binding.characteristics <- get.binding.characteristics(chip.data,srange=c(50,500),bin=5, accept.all.tags = T)
print(paste("binding peak separation distance =",binding.characteristics$peak$x))

[Devon Ryan]

This is only done for single-end data (otherwise, the fragment length distribution can simply be measured) and the ENCODE cross-correlation stuff is...incredibly questionable (the amount of filter and mucking with the data that's done to get the results is, frankly, absurd).