Brandon:LabNotes/Project1/2012-7-27

From ZhangLabWiki
Revision as of 02:26, 28 July 2012 by >Bsos (Created page with "==Peak calling and data analysis, continued from [http://genome-tech.ucsd.edu/LabNotes/index.php/Brandon:LabNotes/Project1/2012-7-20 IVT method accessibilty comp]== *seems to...")
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

Peak calling and data analysis, continued from IVT method accessibilty comp

  • seems to be working. comparing to GM12878 data to see effectiveness
  • macs seems to be calling more significant peaks when analyzing from extracted bed files, but less total bps in covered areas.
  • analyzed previous data wrong since data in bed files did not overlap.
    • thus if 2 reads overlapped, it would be 3 counts in the bed file.
    • if 3 did it would count as 5 reads
    • if 6 did, it would be read as ~10 reads.



bam, bed issues, peak calling discrepancies

  • issues found in bam files, bed files, using macs14, peak callind


peak calling issues BED files


  • when using bedtools to call genome coverage with the coverage.bed function, it results in formation of bed entries that to not overlap. Thus if 3 reads partially overlap, 5 entries will be made based on the points where reads start/begin
genomeCoverageBed -bg -ibam $uniqueBam -g /home/bsos/software/BEDTools-Version-2.12.0/genomes/human.hg19.genome > $bed
  • This results in artificially inflated read counts for for overlapping reads, which is why many more peak calls are made.
  • Though this is useful for viewing on UCSC genome browser which requires smaller files as inputs.


  • all analyses are being fixed by not using BEDs generated from using genomeCoverageBed. (use SAM files, or extracted BED files from SAM files)


  • 6 reads reported from BAM file (IVT_extracted_BED tract) are reported as 10 reads in the genomeCoverageBed generated file (IVT_1000_cell tract).

File:BAMorExtractedBed versus coverageBED.png

  • 3 reads reported from BAM file (IVT_extracted_BED tract) are reported as 5 reads in the genomeCoverageBed generated file (IVT_1000_cell tract).

File:BAMorExtractedBed versus coverageBED2.png



bam/extracted bed file calling discrepancies


  • macs14 seems to be better at calling peaks directly from the bam file. THus BED files should probabaly be only reserved for viewing reads etc on the UCSC genome browswer. As can be seen, more peaks are called in the accessibility assay on cells versus pure DNA, which is what we should be seeing.
    • first column is extracted BED file from BAM.
samtools view s_5_1_Indx49.unique.bam | awk '{print $3"\t"$4"\t"$4+60}'| sort -u > s_5_1_Indx49.unique.bed
    • second columns is unique bam file
    • 3rd column is bed file from using genomeCoverageBed
' bed_extracted bam bed_coverage '
peaks_run2 peaks_run2 peaks2 % difference
10000 348 656 3143 1.885057471
1000 2,549 3,343 11247 1.311494704
100 131 360 651 2.748091603
6ng 642 564 6529 0.878504673
600pg 60 20 447 0.333333333
60pg 24 24 25 1
all most all
  • check UCSC data and see how that compares. make graphs?
  • ran replicates to ensure peak calling was similar between runs and was reproducible.

differences in peaking calling nbtn extracted bed/bam

comparision between GM12878 datas