Skip to content

predict.py fails with segmentation fault due to hicstraw incompatibility with URL hic file #258

Description

@jaavedm

Running the example K562 chr22 example fails with segmentation fault.
Specifically, the version of hicstraw 1.3.1 that I have installed does not work with URLs. Instead, the example hic file should be a local file instead inside of the config.tsv file. I'm creating this ticket so ABC devs can investigate if to fix this behavior in code, or modify docs so that others know how to workaround this issue.

Here is the error I received:

[Sun Jun 15 10:47:31 2025]
rule create_predictions:
    input: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/EnhancerList.txt, /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/GeneList.txt
    output: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutative.tsv.gz, /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutativeNonExpressedGenes.tsv.gz
    jobid: 1
    reason: Missing output files: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutative.tsv.gz, /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutativeNonExpressedGenes.tsv.gz
    wildcards: biosample=K562_chr22
    resources: tmpdir=/tmp, mem_mb=20000, mem_mib=19074

reading genes
reading enhancers
Making predictions for chromosome: chr22
Making putative predictions table...
Done. There are 2391900 putative enhancers for chromosome chr22
Elapsed time: 3.187191963195801
Begin HiC
[Sun Jun 15 10:47:35 2025]
Error in rule create_predictions:
    jobid: 1
    input: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/EnhancerList.txt, /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/GeneList.txt
    output: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutative.tsv.gz, /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions/EnhancerPredictionsAllPutativeNonExpressedGenes.tsv.gz
    conda-env: /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/.snakemake/conda/1dcda64c06900c4f07d893f5fbaf1a3b_
    shell:

                python /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/workflow/scripts/predict.py                         --enhancers /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/EnhancerList.txt                  --outdir /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions                        --score_column ABC.Score                    --chrom_sizes /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/reference/hg38/GRCh38_EBV.no_alt.chrom.sizes.tsv                         --accessibility_feature DHS                     --cellType K562_chr22                   --genes /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/GeneList.txt                  --hic_gamma 1.024238616787792                       --hic_scale 5.9594510043736655                  --hic_pseudocount_distance 5000                         --hic_file https://www.encodeproject.org/files/ENCFF621AIY/@@download/ENCFF621AIY.hic --hic_type hic --hic_resolution 5000                  --scale_hic_using_powerlaw

        (one of the commands exited with non-zero exit code; note that snakemake uses bash strict mode!)

Shutting down, this might take some time.
Exiting because a job execution failed. Look above for error message
Complete log: .snakemake/log/2025-06-15T104730.509317.snakemake.log

Snakemake log did not have anything about the segfault, so I ran it outside of snakemake:

python /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/workflow/scripts/predict.py \
	--enhancers /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/EnhancerList.txt \
	--outdir /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Predictions \
	--score_column ABC.Score \
	--chrom_sizes /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/reference/hg38/GRCh38_EBV.no_alt.chrom.sizes.tsv \
	--accessibility_feature DHS \
	--cellType K562_chr22 \
	--genes /storage/jaavedm/data1/packages/ABC-Enhancer-Gene-Prediction/results/K562_chr22/Neighborhoods/GeneList.txt \
	--hic_gamma 1.024238616787792 \
	--hic_scale 5.9594510043736655 \
	--hic_pseudocount_distance 5000 \
	--hic_file https://www.encodeproject.org/files/ENCFF621AIY/@@download/ENCFF621AIY.hic \
	--hic_type hic \
	--hic_resolution 5000 \
	--scale_hic_using_powerlaw
reading genes
reading enhancers
Making predictions for chromosome: chr22
Making putative predictions table...
Done. There are 2391900 putative enhancers for chromosome chr22
Elapsed time: 3.2994375228881836
Begin HiC
HiC file is : https://www.encodeproject.org/files/ENCFF621AIY/@@download/ENCFF621AIY.hic
Segmentation fault (core dumped)

Here is how I debugged it:

  1. Activate the abc conda environment
    conda activate abc-env

  2. Change to the ABC-Enhancer-Gene-Prediction git cloned directory

  3. Download the example .hic file.
    curl -L -o ENCFF621AIY.hic https://www.encodeproject.org/files/ENCFF621AIY/@@download/ENCFF621AIY.hic

  4. Check if hicstraw (the offending line inside of predictor.py) causes the segfault.
    open python and run these lines:

>>> import hicstraw
>>>print(hicstraw.__version__)
1.3.1
>>> hic = hicstraw.HiCFile("https://www.encodeproject.org/files/ENCFF621AIY/@@download/ENCFF621AIY.hic")
Segmentation fault (core dumped)

It does

After downloading the .hic file, run this code in python

>>> import hicstraw
>>> hic = hicstraw.HiCFile("ENCFF621AIY.hic")
>>> print(hic.getChromosomes()) 

IT WORKS

  1. Now modify the config/config_biosamples_chr22.tsv to look like this instead.
biosample       DHS     ATAC    H3K27ac default_accessibility_feature   HiC_file        HiC_type        HiC_resolution  alt_TSS alt_genes
K562_chr22      example_chr/chr22/ENCFF860XAE.chr22.sorted.se.bam               example_chr/chr22/ENCFF790GFL.chr22.sorted.se.bam       DHSENCFF621AIY.hic  hic     5000    example_chr/chr22/RefSeqCurated.170308.bed.CollapsedGeneBounds.chr22.hg38.TSS500bp.bed  example_chr/chr22/RefSeqCurated.170308.bed.CollapsedGeneBounds.chr22.hg38.bed

Metadata

Metadata

Assignees

No one assigned

    Labels

    bugSomething isn't working

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions