<p>You can do it by awk.</p>
<pre><code>awk -v OFS=&quot;\t&quot; '{ if ( (<span class="math-container">$3-$</span>2) &gt;= 500 &amp;&amp; (<span class="math-container">$3-$</span>2) &lt;= 1000) { print } }' in.bed &gt; out.bed
</code></pre>
<p>Assuming there is no header.</p>

CC BY-SA 4.0
2020-07-22T06:43:10.943Z
2020-07-22T06:43:10.943Z
geek_y
<p>You can use <code>rtracklayer</code> from <code>R bioconductor</code> to import it, filter and do much more stuff, i use the code from this <a href="https://charlesjb.github.io/How_to_import_narrowPeak/" rel="nofollow noreferrer">blog</a> for importing the narrow bed file. Below I use an example file since you did not provide yours:</p>
<pre><code>library(rtracklayer)
fl = &quot;zas1_default_peaks.narrowPeak&quot;

extraCols_narrowPeak &lt;- c(signalValue = &quot;numeric&quot;, pValue = &quot;numeric&quot;,
                          qValue = &quot;numeric&quot;, peak = &quot;integer&quot;)
gr &lt;- import(fl, format = &quot;BED&quot;,
                        extraCols = extraCols_narrowPeak)

head(gr)

GRanges object with 6 ranges and 6 metadata columns:
      seqnames        ranges strand |                name     score signalValue
         &lt;Rle&gt;     &lt;IRanges&gt;  &lt;Rle&gt; |         &lt;character&gt; &lt;numeric&gt;   &lt;numeric&gt;
  [1] AB325691     8982-9219      * | zas1_default_peak_1        18     1.31697
  [2] AB325691   16685-18050      * | zas1_default_peak_2        59     1.51283
  [3]        I 101971-102188      * | zas1_default_peak_3        17     1.23717
  [4]        I 102390-103538      * | zas1_default_peak_4        19     1.24967
  [5]        I 152286-152496      * | zas1_default_peak_5        23     1.37631
  [6]        I 175025-175577      * | zas1_default_peak_6        69     1.62063
         pValue    qValue      peak
      &lt;numeric&gt; &lt;numeric&gt; &lt;integer&gt;
  [1]   3.37358   1.89101       146
  [2]   8.12841   5.92232       419
  [3]   3.26012   1.79545        64
  [4]   3.47875   1.97673       708
  [5]   3.89921   2.32948       136
  [6]   9.34429    6.9608       374
</code></pre>
<p>The narrowBed and peaks.xls contains the same information. If you really really need to read in the xls, you can do:</p>
<pre><code>library(GenomicRanges)
peaks = read.table(&quot;zas1_default_peaks.xls&quot;,comment=&quot;#&quot;,header=TRUE)
gr = makeGRangesFromDataFrame(peaks,keep.extra.columns=TRUE)
</code></pre>
<p>Once you have the genomic ranges object, filtering is the same:</p>
<pre><code> gr[width(gr)&gt;500 &amp; width(gr)&lt;1000]
GRanges object with 208 ranges and 7 metadata columns:
        seqnames          ranges strand |    length abs_summit    pileup
           &lt;Rle&gt;       &lt;IRanges&gt;  &lt;Rle&gt; | &lt;integer&gt;  &lt;integer&gt; &lt;numeric&gt;
    [1]        I   175025-175577      * |       553     175399    186.99
    [2]        I   186807-187489      * |       683     187009    208.66
    [3]        I   339070-339700      * |       631     339258    222.91
    [4]        I   431557-432088      * |       532     431734    195.66
</code></pre>

CC BY-SA 4.0
2020-07-22T22:06:49.100Z
Thank you so much for the answer. it works fine with narrow peaks and bed files but does not work with xls files. import function supports other formats like GFF,TwoBit, Wig, bedGRaph but not xls.
Mariam
CC BY-SA 4.0
2020-07-23T08:11:51.000Z
it contains the same information. I have edited my answer to show how u can do the import
StupidWolf
CC BY-SA 4.0
2020-07-23T18:46:35.420Z
Many Thanks, Is there any way to export the files as .xls to further analyze them using diffBind?
Mariam
CC BY-SA 4.0
2020-07-22T09:34:37.803Z
2020-07-23T08:11:30.537Z
2020-07-23T08:11:30.537Z
StupidWolf
<p>My data contain several chip-seq results. I have the peaks called by MACS2.I wanted to only look at those peaks that their size is e.g 500bp to 1000bp. How can I separate those peaks efficiently?</p>
<p>I know one way is to take the peaks.xls file sort them based on the length column and take out the peaks less than 1kb. But I was looking for a more efficient way to do that. Any suggestion?</p>
<p>Many Thanks</p>

CC BY-SA 4.0
2020-07-21T21:37:33.240Z
I would just write an R script or a Python to read the peaks and write those lines to new files according to their lengths.
Phoenix Mu
<sequence-analysis><differential-expression><chip-seq><peak-calling><macs2>
Separating peaks of chip-seq with specific length
# Separating peaks of chip-seq with specific length

(Asked by: username_0 on 2020-07-21)

My data contain several chip-seq results. I have the peaks called by MACS2.I wanted to only look at those peaks that their size is e.g 500bp to 1000bp. How can I separate those peaks efficiently?

I know one way is to take the peaks.xls file sort them based on the length column and take out the peaks less than 1kb. But I was looking for a more efficient way to do that. Any suggestion?

Many Thanks

**Tags:** sequence-analysis, differential-expression, chip-seq, peak-calling, macs2
> **Comments:**
> **username_1**: I would just write an R script or a Python to read the peaks and write those lines to new files according to their lengths.  

---

## 2 answers

### Answer by username_2 on 2020-07-22. Score: 1 | Accepted answer

You can do it by awk.
awk -v OFS="\t" '{ if ( ($3-$2) >= 500 && ($3-$2) <= 1000) { print } }' in.bed > out.bed

Assuming there is no header.


---

### Answer by username_3 on 2020-07-22. Score: 1

You can use rtracklayer from R bioconductor to import it, filter and do much more stuff, i use the code from this blog for importing the narrow bed file. Below I use an example file since you did not provide yours:
library(rtracklayer)
fl = "zas1_default_peaks.narrowPeak"

extraCols_narrowPeak <- c(signalValue = "numeric", pValue = "numeric",
                          qValue = "numeric", peak = "integer")
gr <- import(fl, format = "BED",
                        extraCols = extraCols_narrowPeak)

head(gr)

GRanges object with 6 ranges and 6 metadata columns:
      seqnames        ranges strand |                name     score signalValue
                 |             
  [1] AB325691     8982-9219      * | zas1_default_peak_1        18     1.31697
  [2] AB325691   16685-18050      * | zas1_default_peak_2        59     1.51283
  [3]        I 101971-102188      * | zas1_default_peak_3        17     1.23717
  [4]        I 102390-103538      * | zas1_default_peak_4        19     1.24967
  [5]        I 152286-152496      * | zas1_default_peak_5        23     1.37631
  [6]        I 175025-175577      * | zas1_default_peak_6        69     1.62063
         pValue    qValue      peak
        
  [1]   3.37358   1.89101       146
  [2]   8.12841   5.92232       419
  [3]   3.26012   1.79545        64
  [4]   3.47875   1.97673       708
  [5]   3.89921   2.32948       136
  [6]   9.34429    6.9608       374

The narrowBed and peaks.xls contains the same information. If you really really need to read in the xls, you can do:
library(GenomicRanges)
peaks = read.table("zas1_default_peaks.xls",comment="#",header=TRUE)
gr = makeGRangesFromDataFrame(peaks,keep.extra.columns=TRUE)

Once you have the genomic ranges object, filtering is the same:
 gr[width(gr)>500 & width(gr)<1000]
GRanges object with 208 ranges and 7 metadata columns:
        seqnames          ranges strand |    length abs_summit    pileup
                     |    
    [1]        I   175025-175577      * |       553     175399    186.99
    [2]        I   186807-187489      * |       683     187009    208.66
    [3]        I   339070-339700      * |       631     339258    222.91
    [4]        I   431557-432088      * |       532     431734    195.66

> **Comments:**
> **username_0**: Thank you so much for the answer. it works fine with narrow peaks and bed files but does not work with xls files. import function supports other formats like GFF,TwoBit, Wig, bedGRaph but not xls.  
> **username_3**: it contains the same information. I have edited my answer to show how u can do the import  
> **username_0**: Many Thanks, Is there any way to export the files as .xls to further analyze them using diffBind?  

---

