aryeelab / hichipper

A preprocessing and QC pipeline for HiChIP data
MIT License
33 stars 12 forks source link

Error using GRCh38.p12 human genome assembly #59

Closed bressin closed 5 years ago

bressin commented 5 years ago

Hi, Unfortunately, I get consistently an error when I run hichipper using the GRCh38 p12 human genome assembly. Apparently it works fine when I'm using the hg19.

Mon Nov 19 12:18:38 CET 2018: Starting hichipper pipeline v0.7.3 Mon Nov 19 12:18:38 CET 2018: Parsing user parameters Mon Nov 19 12:18:38 CET 2018: .yaml file detected Mon Nov 19 13:40:13 CET 2018: Performing hichipper-modified background peak calling Error in [[<-(*tmp*, name, value = c(366L, 677L, 705L, 723L, 729L, : 130901107 elements in value to replace 130901117 elements Calls: computeRatioEtc -> $<- -> $<- -> [[<- -> [[<- Execution halted Mon Nov 19 14:05:06 CET 2018: Modified background signal; performing peak calling Mon Nov 19 14:05:06 CET 2018: MACS2 TEXT OUTPUT INCOMING INFO @ Mon, 19 Nov 2018 14:05:06: Read and build treatment bedGraph... Traceback (most recent call last): File "/pkg/python-2.7.15-2/bin/macs2", line 622, in main() File "/pkg/python-2.7.15-2/bin/macs2", line 73, in main run( args ) File "/home/bressin/.local/lib/python2.7/site-packages/MACS2/bdgcmp_cmd.py", line 40, in run tbtrack = tbio.build_bdgtrack() File "MACS2/IO/BedGraphIO.pyx", line 97, in MACS2.IO.BedGraphIO.bedGraphIO.build_bdgtrack (MACS2/IO/BedGraphIO.c:1009) IOError: [Errno 2] No such file or directory: 'hichipper_all/adjustedTreatment.bdg.tmp'

Do you have any idea what happens? I use only canonical chromosomes chr1-chr22,chrX,chrY and chrM. Surprisingly it works when I add ChIP-seq data but I need also to results for the EACH,ALL mode.

Here is my *.yaml file

peaks:

  • EACH,ALL resfrags:
  • /project/owlmayer/Applications/05_HiChIP_Pipeline/data/annotation/GRCh38.p12.MboI.annotation.bed hicpro_output:
  • HiCProOutput

and the log file:

Mon Nov 19 12:18:38 CET 2018: Starting hichipper pipeline v0.7.3 Mon Nov 19 12:18:38 CET 2018: Executed from: /project/HiChIP_Pipeline/HiChIP_dTAG7_K562/GRCh38_p12/DMSO Mon Nov 19 12:18:38 CET 2018: Output folder: /project/HiChIP_Pipeline/HiChIP_dTAG7_K562/GRCh38_p12/DMSO/hichipper_all Mon Nov 19 12:18:38 CET 2018: Parsed manifest as follows: {'peaks': ['EACH,ALL'], 'hicpro_output': ['HiCProOutput'], 'resfrags': ['/project/owlmayer/Applications/05_HiChIP_Pipeline/data/annotation/GRCh38.p12.MboI.annotation.bed']} Mon Nov 19 12:18:38 CET 2018: Determined that the following samples are good to go: ['MRA14', 'MRA15', 'MRA13'] Mon Nov 19 12:18:38 CET 2018: User defined peaks specification: EACH,ALL Mon Nov 19 12:18:38 CET 2018: Calling peaks from all reads for each sample Mon Nov 19 12:34:07 CET 2018: macs2 command: macs2 callpeak -t hichipper_all/MRA14.all.Pairs.bed.tmp --keep-dup all -q 0.01 --extsize 147 --nomodel -g hs -B -f BED --verbose 0 -n hichipper_all/MRA14_temporary Mon Nov 19 14:21:39 CET 2018: macs2 command: macs2 callpeak -t hichipper_all/MRA15.all.Pairs.bed.tmp --keep-dup all -q 0.01 --extsize 147 --nomodel -g hs -B -f BED --verbose 0 -n hichipper_all/MRA15_temporary Mon Nov 19 16:23:33 CET 2018: macs2 command: macs2 callpeak -t hichipper_all/MRA13.all.Pairs.bed.tmp --keep-dup all -q 0.01 --extsize 147 --nomodel -g hs -B -f BED --verbose 0 -n hichipper_all/MRA13_temporary Mon Nov 19 18:12:23 CET 2018: Performing restriction fragment-aware padding Mon Nov 19 18:14:12 CET 2018: Processing MRA14 Mon Nov 19 18:14:12 CET 2018: Total_PETs=110890120 Mon Nov 19 18:14:32 CET 2018: Mapped_unique_quality_pairs=146526416 Mon Nov 19 18:14:39 CET 2018: Mapped_unique_quality_valid_pairs=59035095 Mon Nov 19 18:14:39 CET 2018: Intersecting PETs with anchors Mon Nov 19 18:14:39 CET 2018: Something went wrong in determining peaks for anchor inference; rerun with the flag to debug.

Thank you for your help. Annkatrin

bressin commented 5 years ago

Ok, I think I found my mistake. I went through the code and I see that the script is using the resfrags where it fails. I had Scaffolds in the resfrags file ... is this a potential reason?

bressin commented 5 years ago

Fixed it. Sorry.

caleblareau commented 5 years ago

Did it wind up being the scaffolds?

On Nov 20, 2018, at 1:36 AM, bressin notifications@github.com wrote:

Fixed it. Sorry.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub, or mute the thread.

bressin commented 5 years ago

Yes.

armenabnousi commented 4 years ago

I'm having the same problem when using GRCh38 (restriction fragments marked by HiC-Pro utils tool). I'm wondering how did you end up working around this. Thanks!

bressin commented 4 years ago

Remove the Scaffolds :)