humanlongevity / HLA

xHLA: Fast and accurate HLA typing from short read sequence data
Other
101 stars 52 forks source link

Script typing.r breaks #38

Open abhisheknrl opened 5 years ago

abhisheknrl commented 5 years ago

Hi,

I managed to generate a new bam file using get-reads-alt-unmap.sh. When I use this new bam file as an input, the code breaks while running typing.r. When I follow the cat statements, I discovered that it breaks between line 310 and 350 in the section starting

more <- do.call(rbind, mclapply(solution, function(s){..............}

The error message reads like this.

/opt/bin/typer.sh: line 45:   822 Killed                  $BIN/typing.r $OUT/${ID}.tsv $OUT/${ID}.hla
Traceback (most recent call last):
  File "/opt/bin/run.py", line 64, in <module>
    check_call(bin_args)
  File "/usr/lib/python2.7/subprocess.py", line 540, in check_call
    raise CalledProcessError(retcode, cmd)
subprocess.CalledProcessError: Command '['/opt/bin/typer.sh', '/myd/data/bam/NA06985.aligned.sorted_xhla.bam', 'NA06985.aligned.sorted', 'full']' returned non-zero exit status 137

Can anyone replicate this? Thanks!

abhisheknrl commented 5 years ago

@tanghaibao Any thoughts on this?

Thanks!

cbrunoels commented 5 years ago

I'm having the same issue. I had a GRCh37 referenced BAM file of a WES sample that I converted to hg38 using CrossMap. I tried running the hg38 sorted BAM file but it returned the exact same error @abhisheknrl got, so I tried to fix it using get-reads-alt-unmap.sh and still received the same exact error. Here is the full output when I ran the "fixed" BAM file (after using get-reads-alt-unmap.sh):

[28/Mar/2019 23:16:37] INFO - Xie Chao's HLA typing algorithm
[28/Mar/2019 23:16:37] INFO - Sample_id: SAMPLE Input file: SAMPLE.hg38.fixed.bam
typer.sh parameters: DELETE=false FULL=false
Extracting reads from S3
Aligning reads to IMGT database
processing FASTQ file
    found 22993 reads
processing MSA file
    found 26381 HLA exons
processing FASTQ file
    matched to 25219 HLA exons
    4304 reads matched to HLA types not in the MSA file
translating matches to MSA
Typing
Read 1948430 rows and 16 (of 16) columns from 0.254 GB file in 00:00:03
[1] 1347262
[1] 1347262
[1] 7284
[1] 3389
[1] 11220
[1] 338793
[1] 289591
[1] 87944
[1] 201647
dcast
done dcast
weight
done weight
set to 1
done set to 1
dcast2
done dcast2
lpsolve
/opt/bin/typer.sh: line 45:   541 Killed                  $BIN/typing.r $OUT/${ID}.tsv $OUT/${ID}.hla
Traceback (most recent call last):
  File "/opt/bin/run.py", line 64, in <module>
    check_call(bin_args)
  File "/usr/lib/python2.7/subprocess.py", line 540, in check_call
    raise CalledProcessError(retcode, cmd)
subprocess.CalledProcessError: Command '['/opt/bin/typer.sh', 'SAMPLE.hg38.fixed.bam', 'SAMPLE']' returned non-zero exit status 137

For the record, the example executes perfectly, so I am sure that xHLA was successfully installed before I tried to run it. Here are some details about my OS if that might help:

macOS: 10.14.4-x86_64
CLT: 10.2.0.0.1.1552586384
Xcode: 10.2
XQuartz: 2.7.11 => /opt/X11
R: 3.4.0
abhisheknrl commented 5 years ago

Hi, I received an email with a suggestion in response to this post. Have not tried this but posting it here in case it helps @cbrunoels or others with the same problem.

Hi you need to install sambamba
http://lomereiter.github.io/sambamba/

if still not found sambamba command
you can edit get-reads-alt-unmap.sh
add absolute path before sambamba command
cbrunoels commented 5 years ago

Thanks @abhisheknrl. Unfortunately I already have sambamba up and running so this didn't help. For clarification I was able to run the get-reads-alt-unmap.sh script successfully, but when I try to actually impute the HLA alleles the script fails where I pointed out above.

Based on the error, it seems like the issue comes from typing.r line 220: lps <- lp('max', f.obj, f.con, f.dir, f.rhs, int.vec = which(f.type == 'i'), binary.vec = which(f.type == 'b')) The line before this in typing.r prints lpsolve and line 221 prints done lpsolve, so the issue must be localized to here. Is there any further troubleshooting I can do to find what's going on here?

abhisheknrl commented 5 years ago

No reply from the devs in 4 months. Probably means 'Forget about this, use something else'. At least that is how I have understood.