Seqanswers Leaderboard Ad

Collapse

Announcement

Collapse
No announcement yet.
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • bedGraphToBigWig does not work with single base coverage

    Hi,

    I've identified a "bug" in UCSC's bedGraphToBigWig.
    It does not work with single base coverage.

    If I have a file where the start and end are always equal (i.e. single base coverage), the bigWig file generated by bedGraphToBigWig is empty.
    Does anyone have a workaround, or am I simply not respecting the format for bedgraph files? If I want to include the count just for base 1, how am I supposed to write it?
    1 1 count
    or
    1 2 count
    The 1st option appears correct but bedGraphToBigWig will not take it.
    It would be a big job to convert all my bedgraph files to the format bedGraphToBigWig appears to want.

    Here are 2 example bedgraph files: test1.bedgraph with the start and the end equal to 1, and test2.bedgraph where the start and end are 1 and 2, respectively. The bigWig files generated by bedGraphToBigWig is empty in the first case, and with 1 base covered in the second case.

    blancha@sw-2r12-n36 ~]$ more test1.bedgraph
    chr1 1 1 1
    [blancha@sw-2r12-n36 ~]$ more test2.bedgraph
    chr1 1 2 1
    [blancha@sw-2r12-n36 ~]$ bedGraphToBigWig test1.bedgraph UCSCChromInfo.txt test1.bw
    [blancha@sw-2r12-n36 ~]$ bedGraphToBigWig test2.bedgraph UCSCChromInfo.txt test2.bw
    [blancha@sw-2r12-n36 ~]$ bigWigInfo test1.bw
    version: 4
    isCompressed: yes
    isSwapped: 0
    primaryDataSize: 30
    primaryIndexSize: 6,204
    zoomLevels: 1
    chromCount: 1
    basesCovered: 0
    mean: -nan
    min: 1.000000
    max: 1.000000
    std: -nan
    [blancha@sw-2r12-n36 ~]$ bigWigInfo test2.bw
    version: 4
    isCompressed: yes
    isSwapped: 0
    primaryDataSize: 32
    primaryIndexSize: 6,204
    zoomLevels: 1
    chromCount: 1
    basesCovered: 1
    mean: 1.000000
    min: 1.000000
    max: 1.000000
    std: 0.000000

    Thank you for your help.
    Last edited by blancha; 05-28-2014, 01:51 PM.

  • #2
    I think your bedgraph should be the second way, 1 2 not 1 1, because bedgraph files are "half open".

    Comment


    • #3
      Great. Thanks BioComputer.
      That is exactly the formal definition I was looking for.

      I did waste half a day figuring out the "bug", though.
      Thankfully, you had a quick answer to my question, as I would never have figured out the reason behind the "bug" on my own.
      I guess it's more of a "feature".
      Now, I just need to decide if I should write a Python script to correct the coordinates, or whether if I run bismark_methylation_extractor again with the --zero_based parameter I will get the desired format.

      I used bismark_methylation_extractor to generate the BED files. The output was a bit confusing since, by default, the bedgraph file outputted has 0-based start coords and 1-based end coords, whereas the other BED files generated have 1-based genomic coords.

      To be fair, it is in the documentation, although I'm still unsure if adding the --zero-based parameter will result in 0-based start and end coordinates, or just 0-based start coordinates, and 1-based end coordinates. I imagine it is the latter.
      Bioinformatics is so much fun.

      Thanks again!

      Comment


      • #4
        So you just have to add 1 to each value in the end position column? I'm still learning this stuff so there's probably a better way to do this, but this command below worked when I tested it on a small file. It assumes tab-separated columns, should add 1 to each value in the 3rd column, and should print out any non-position lines (ie. any lines that don't start with "chr" like a file header) as is. Test first and use at your own risk

        Code:
        awk '/^chr/ {FS="\t";$3=$3+1;print} !/^chr/{print}' inputFile > outputFile
        Or instead of adding 1 to the end position do you need to substract 1 from each start position?

        Code:
        awk '/^chr/ {FS="\t";$2=$2-1;print} !/^chr/{print}' inputFile > outputFile
        Here's another version that should do the same thing, maybe it'll be faster?

        Code:
        awk 'BEGIN { OFS="\t"}{if (/^chr/) $3=$3+1; print}' inputFile > outputFile
        Last edited by biocomputer; 05-28-2014, 05:35 PM.

        Comment


        • #5
          If you are going to use UCSC browser to display the files consider this: http://genome.ucsc.edu/FAQ/FAQtracks.html#tracks1

          Comment

          Latest Articles

          Collapse

          • seqadmin
            Latest Developments in Precision Medicine
            by seqadmin



            Technological advances have led to drastic improvements in the field of precision medicine, enabling more personalized approaches to treatment. This article explores four leading groups that are overcoming many of the challenges of genomic profiling and precision medicine through their innovative platforms and technologies.

            Somatic Genomics
            “We have such a tremendous amount of genetic diversity that exists within each of us, and not just between us as individuals,”...
            05-24-2024, 01:16 PM
          • seqadmin
            Recent Advances in Sequencing Analysis Tools
            by seqadmin


            The sequencing world is rapidly changing due to declining costs, enhanced accuracies, and the advent of newer, cutting-edge instruments. Equally important to these developments are improvements in sequencing analysis, a process that converts vast amounts of raw data into a comprehensible and meaningful form. This complex task requires expertise and the right analysis tools. In this article, we highlight the progress and innovation in sequencing analysis by reviewing several of the...
            05-06-2024, 07:48 AM

          ad_right_rmr

          Collapse

          News

          Collapse

          Topics Statistics Last Post
          Started by seqadmin, 05-24-2024, 07:15 AM
          0 responses
          15 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 05-23-2024, 10:28 AM
          0 responses
          17 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 05-23-2024, 07:35 AM
          0 responses
          21 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 05-22-2024, 02:06 PM
          0 responses
          10 views
          0 likes
          Last Post seqadmin  
          Working...
          X