<p>You can exclude the <code>which</code> statement and get the entirety of all of the chromosomes. Note however that R is a very bad platform to use for this, since you end up reading everything into memory (I hope the BAM file isn't huge).</p>
<p>Note that you're seriously reinventing the wheel. There are a number of packages that will directly produce the plot you want, such as <code>bamPEFragmentSize</code> from deepTools.</p>

CC BY-SA 4.0
2021-01-06T14:59:13.767Z
"bamPEFragmentSize" will see this
kcm
CC BY-SA 4.0
2021-01-06T15:51:54.760Z
" you're seriously reinventing the wheel." yes since I didn;t read about your tool "bamPEFragmentSize" I used it , Works absolutely fine
kcm
CC BY-SA 4.0
2021-01-09T15:32:57.093Z
is there any way by using this "bamPEFragmentSize" i can get smooth curve ?
kcm
CC BY-SA 4.0
2021-01-09T15:48:39.823Z
Save the values to a file and replot however you like.
Devon Ryan
CC BY-SA 4.0
2021-01-09T16:20:56.070Z
okay i will see the option to save the file ....right now i ran this "bamPEFragmentSize -hist fragmentSize.png --maxFragmentLength 100 -b HSC1.bam Mono1.bam LSC1.bam Blast1.bam --numberOfProcessors 10 --samplesLabel HSC Mono LSC Blast
"
kcm
CC BY-SA 4.0
2021-01-06T12:28:15.470Z
2021-01-06T12:28:15.470Z
Devon Ryan
<p>My tutorial <a href="https://rockefelleruniversity.github.io/RU_ATAC_Workshop.html" rel="nofollow noreferrer">link</a></p>
<pre><code>library(GenomicAlignments)

atacReads &lt;- readGAlignmentPairs(&quot;MACS2_peak_call/ALL_ATAC_data/HSC1.bam&quot;, param = ScanBamParam(mapqFilter = 1, 
            flag = scanBamFlag(isPaired = TRUE, isProperPair = TRUE), what = c(&quot;qname&quot;, 
        &quot;mapq&quot;, &quot;isize&quot;), which = GRanges(&quot;chr20&quot;, IRanges(1, 63025520))))
# length(atacReads)
atacReads


atacReads_read1 &lt;- GenomicAlignments::first(atacReads)
insertSizes &lt;- abs(elementMetadata(atacReads_read1)$isize)
head(insertSizes)


library(magrittr)
library(dplyr)
library(ggplot2)
fragLenPlot &lt;- table(insertSizes) %&gt;% data.frame %&gt;% dplyr::rename(InsertSize = insertSizes, 
    Count = Freq) %&gt;% mutate(InsertSize = as.numeric(as.vector(InsertSize)), 
 Count = as.numeric(as.vector(Count))) %&gt;% ggplot(aes(x = InsertSize, y = Count)) + 
  geom_line()

fragLenPlot + theme_bw()

fragLenPlot + scale_y_continuous(trans = &quot;log2&quot;) + theme_bw()
</code></pre>
<p>How to i specify the Grange to span the total chromosome. Now I use only &quot;chr20&quot; as shown in the example.</p>

CC BY-SA 4.0
2021-01-06T11:09:18.610Z
sorry don't quite get the question. is it about knowing the chromosome length from the bam file, or it's something to do with the plotting
StupidWolf
CC BY-SA 4.0
2021-01-06T11:24:13.500Z
Now Im too bit confused ,after reading paper where they don't mention which chromosome they taking and plotting.  It would be essentially this one " is it about knowing the chromosome length from the bam file" can i like  "chr1-chr20" and see the insert size distribution . I was not able to figure out how to define more that one chromosome in that Grange argument //
kcm
CC BY-SA 4.0
2021-01-06T11:26:32.773Z
@StupidWolf for example in this paper https://www.nature.com/articles/s41597-019-0071-0#Sec1  Fig2(b).
kcm
<atac-seq>
ATAC seq : Plotting insert sizes
# ATAC seq : Plotting insert sizes

(Asked by: username_0 on 2021-01-06)

My tutorial link
library(GenomicAlignments)

atacReads <- readGAlignmentPairs("MACS2_peak_call/ALL_ATAC_data/HSC1.bam", param = ScanBamParam(mapqFilter = 1, 
            flag = scanBamFlag(isPaired = TRUE, isProperPair = TRUE), what = c("qname", 
        "mapq", "isize"), which = GRanges("chr20", IRanges(1, 63025520))))
# length(atacReads)
atacReads


atacReads_read1 <- GenomicAlignments::first(atacReads)
insertSizes <- abs(elementMetadata(atacReads_read1)$isize)
head(insertSizes)


library(magrittr)
library(dplyr)
library(ggplot2)
fragLenPlot <- table(insertSizes) %>% data.frame %>% dplyr::rename(InsertSize = insertSizes, 
    Count = Freq) %>% mutate(InsertSize = as.numeric(as.vector(InsertSize)), 
 Count = as.numeric(as.vector(Count))) %>% ggplot(aes(x = InsertSize, y = Count)) + 
  geom_line()

fragLenPlot + theme_bw()

fragLenPlot + scale_y_continuous(trans = "log2") + theme_bw()

How to i specify the Grange to span the total chromosome. Now I use only "chr20" as shown in the example.

**Tags:** atac-seq
> **Comments:**
> **username_1**: sorry don't quite get the question. is it about knowing the chromosome length from the bam file, or it's something to do with the plotting  
> **username_0**: Now Im too bit confused ,after reading paper where they don't mention which chromosome they taking and plotting.  It would be essentially this one " is it about knowing the chromosome length from the bam file" can i like  "chr1-chr20" and see the insert size distribution . I was not able to figure out how to define more that one chromosome in that Grange argument //  
> **username_0**: @StupidWolf for example in this paper https://www.nature.com/articles/s41597-019-0071-0#Sec1  Fig2(b).  

---

## 1 answers

### Answer by username_2 on 2021-01-06. Score: 1 | Accepted answer

You can exclude the `which` statement and get the entirety of all of the chromosomes. Note however that R is a very bad platform to use for this, since you end up reading everything into memory (I hope the BAM file isn't huge).

Note that you're seriously reinventing the wheel. There are a number of packages that will directly produce the plot you want, such as `bamPEFragmentSize` from deepTools.

> **Comments:**
> **username_0**: "bamPEFragmentSize" will see this  
> **username_0**: " you're seriously reinventing the wheel." yes since I didn;t read about your tool "bamPEFragmentSize" I used it , Works absolutely fine  
> **username_0**: is there any way by using this "bamPEFragmentSize" i can get smooth curve ?  
> **username_2**: Save the values to a file and replot however you like.  
> **username_0**: okay i will see the option to save the file ....right now i ran this "bamPEFragmentSize -hist fragmentSize.png --maxFragmentLength 100 -b HSC1.bam Mono1.bam LSC1.bam Blast1.bam --numberOfProcessors 10 --samplesLabel HSC Mono LSC Blast "  

---

