Error using sequenza-utils with WES

Question

[Camilo Andres]

Hi! I'm trying to calculate the homologous recombination deficiency score of a cell line (MDA-MB-231) using whole exome sequencing data. To do this, I pretend to use the scarHRD package in R, but I first need a "*.seqz.gz" archive, which is made using sequenza-utils, but I have the following error message when I run this tool:

[mpileup] 1 samples in 1 input files
[mpileup] 1 samples in 1 input files
Traceback (most recent call last):
  File "/Users/camilo/micromamba/envs/bioinfo/bin/sequenza-utils", line 33, in <module>
    sys.exit(load_entry_point('sequenza-utils==2.1.9999b1', 'console_scripts', 'sequenza-utils')())
  File "/Users/camilo/micromamba/envs/bioinfo/lib/python3.8/site-packages/sequenza/commands.py", line 39, in main
    modules[args.module](subparsers, args.module, extra, log)
  File "/Users/camilo/micromamba/envs/bioinfo/lib/python3.8/site-packages/sequenza/programs/seqz_binning.py", line 35, in seqz_binning
    out.write(next(bins))
  File "/Users/camilo/micromamba/envs/bioinfo/lib/python3.8/site-packages/sequenza/seqz.py", line 208, in binned_seqz
    yield format_seqz_window(window_i) + replace_gc(
  File "/Users/camilo/micromamba/envs/bioinfo/lib/python3.8/site-packages/sequenza/seqz.py", line 157, in format_seqz_window
    int(round(window[2] / window[6], 0)),
ZeroDivisionError: division by zero

The code is as follows: ` $ sequenza-utils bam2seqz -gc hg38/hg38.wig.gz --fasta hg38/hg38.fa -n MCF10A/SRR3090727_recall_reads.bam -t MDA-MB-231/SRR1021654_recall_reads.bam | sequenza-utils seqz_binning -w 50 -s - | gzip > tumor_small.seqz.gz

I am quite new to this kind of analyses, so I am sorry if it is a very simple issue to solve.

Thanks in advance

[GenoMax]

A quick look at the sequenza-utils cookbook does not show examples of pipe you are using. Do you know if this piping works? You may want to try doing this in two separate steps. Instead of plain gzip you should probably use bgzip.

[Camilo Andres]

This pipeline is the one recommended by scar HRD, which is the package I pretend to use aftewards. I only removed the chromosome specifying part to see whether it worked that way, but it didn't, just got the same error. I will try using your suggestion. Thanks.

PS: I just tried using bgzip, but didn't work.

[GenoMax]

I only removed the chromosome specifying part to see whether it worked that way,

What does this mean?

Can you try doing this in 2 steps:

sequenza-utils bam2seqz -gc hg38/hg38.wig.gz --fasta hg38/hg38.fa -n MCF10A/SRR3090727_recall_reads.bam -t MDA-MB-231/SRR1021654_recall_reads.bam  --output sample.seqz.gz

followed by

sequenza-utils seqz_binning -w 50 -s --seqz sample.seqz.gz -o sample_bin50.seqz.gz

[Camilo Andres]

Thanks, this make me notice that is the first part of the pipeline that was giving me the error... "Segmentation fault: 11" So, I will look for solving this issue now. Thanks.