Error message "States for individual 0 are zero. Error during selection" but imputing a whole chromosome for 300 samples?

[robertwhbaldwin]

I included my log file below. I have 300 samples at around 1x coverage and am trying to impute chromsome1. Based on what I've read the error message I'm receiving should not occur as there will be plenty of variants. Any idea why I'm seeing this error message? Thank You - Robert

./glimpse.sh

[GLIMPSE2] Phase and impute low coverage sequencing data

• Authors : Simone RUBINACCI & Olivier DELANEAU, University of Lausanne
• Contact : simone.rubinacci@unil.ch & olivier.delaneau@unil.ch
• Version : GLIMPSE2_phase v2.0.0 / commit = 3bed6d9 / release = 2022-12-07
• Citation : BiorXiv, (2022). DOI: https://doi.org/10.1101/2022.11.28.518213
• : Nature Genetics 53, 120–126 (2021). DOI: https://doi.org/10.1038/s41588-020-00756-0
• Run date : 22/01/2023 - 10:15:00

Files:

• List BAM/CRAM : [../bam_list]
• Reference binary : [panel_1_1_95772355.bin]
• Output file : [/home/robert/IMP/Glimpse/imputed_chr1_100.bcf]
• Output format : [BCF format | ZLIB compression]

GLIMPSE_phase parameters:

• Imputation model : [Reference panel imputation]
• Input region : [Given by binary reference panel]
• Output region : [Given by binary reference panel]
• Sparse MAF : [Given by binary reference panel]
• Recombination rates : [Given by binary reference panel]
• Ploidy : [Only diploid samples in region]
• Keep monom. ref sites: [NO]

Model parameters:

• Burnin iterations : [5]

• Main iterations : [15]

• Ne [eff. pop. size] : [100000]
• Phase error rate : [0.0001]
• Imputation error rate: [1e-12]
• Min value for hap GLs: [1e-10]

Selection parameters:

• K init : [1000]
• K pbwt : [2000]
• PBWT depth : [12]
• PBWT modulo (cM) : [0.1]
• State list : [No list provided]

Genotype calling:

• Calling model : [standard]
• Perform calling]

BAM/CRAM filters and options:

• Min mapping quality : [10]
• Min base quality : [10]
• Max depth : [40]
• Keep failed QC : [NO]
• Keep orphan reads : [NO]
• Keep duplicate reads : [NO]
• Keep supp. alignment : [NO]
• Check pairing : [NO]
• Ignore orientation : [NO]
• Illumina-1.3+ : [NO]

Other parameters

• Seed : [15052011]

• Threads : [8]

Initialisation:

• Binary reference panel parsing [done] (0.38s)
• Binary reference panel [Nrh=358] [L=754558] [Lrare= 0 (0.0%) - Lcommon= 754558 (100.0%)]
• Input region : [1:1-95772355]
• Output region : [1:1-95772355]
• Sparse MAF : [0.001]
• Reading BAM files [W::hts_idx_load3] The index file is older than the data file: /archive/1xBAMS/TLUN-2018-913-010.bam.bai [W::hts_idx_load3] The index file is older than the data file: /archive/1xBAMS/TLUN-2018-913-053.bam.bai [W::hts_idx_load3] The index file is older than the data file: /archive/1xBAMS/TLUN-2018-913-189.bam.bai [W::hts_idx_load3] The index file is older than the data file: /archive/1xBAMS/TLUN-2018-913-243.bam.bai [W::hts_idx_load3] The index file is older than the data file: /archive/1xBAMS/TLUN-2018-913-326.bam.bai
• Reading BAM files done (75.62s)
• Coverage statistics reported in: [/home/robert/IMP/Glimpse/imputed_chr1_100_stats_coverage.txt.gz]
• Avg. sequencing coverage: 0.947329
• RARE transpose Ref (0.00s)
• Region spans 95698211 bp and 95.70 cM
• No PBWT allocated (Kpbwt=0 or Kpbwt >= n_ref_haps)

Initializing iteration

• Init rare Tar [#non-ref gls: 0] (0.55s)
• HMM imputation [#states=358.0 / %poly=100.0%] (131.09s)

Burn-in iteration [1/5]

• HAP update Tar (1.50s)
• RARE transpose Tar (0.00s)
• No PBWT selection (Kpbwt=0 or Kpbwt >= n_ref_haps)

ERROR: States for individual 0 are zero. Error during selection

ERROR: States for individual 1 are zero. Error during selection.

ERROR: States for individual 2 are zero. Error during selection.

E(base

[srubinacci]

Hi,

Indeed, something is off here. My guess is that you have a small reference panel and the default threshold does not work for you. Would it be possible for you to share the data of this chunk (by email)? I can look at it.

[srubinacci]

I should probably add that GLIMPSE2 is designed for VERY large reference panels, and maybe for your use case GLIMPSE1 is more appropriate.

[robertwhbaldwin]

thanks for the response! I'll have to get permission to share the data. But the reference panel was for rainbow trout Oncorhynchus mykiss and was downloaded from here http://mgb.qnlm.ac/ under downloads. It contained 179 samples, which is small, yes.

On Sun, Jan 29, 2023 at 5:29 PM Simone Rubinacci @.***> wrote:

I should probably add that GLIMPSE2 is designed for VERY large reference panels, and maybe for your use case GLIMPSE1 is more appropriate.

— Reply to this email directly, view it on GitHub https://github.com/odelaneau/GLIMPSE/issues/121#issuecomment-1407774666, or unsubscribe https://github.com/notifications/unsubscribe-auth/AITBRRSDCOV26DKIB3N5ZMDWU3ONHANCNFSM6AAAAAAUDAOUHU . You are receiving this because you authored the thread.Message ID: @.***>

[srubinacci]

I see - thank you for your response. I recommend to use GLIMPSE1 (release: https://github.com/odelaneau/GLIMPSE/releases/tag/v1.1.1, documentation: https://odelaneau.github.io/GLIMPSE/glimpse1/) then.

Thank you for getting in touch about this, will add flags/warnings in G2. Indeed I should make more clear when to use one or the other software. I apologise for the inconvenience.

Let me know if you get any issues with GLIMPSE1 (maybe in a new issue).

[robertwhbaldwin]

thanks for the advice I'll try G1. Does G1 take bam files as input? I'm confused about that part.

[srubinacci]

You're welcome. I should actually check your case more in details. It's likely that you trigger a silly G2 bug. The thing is that you only have 3I58 haplotypes and Kpbwt is set to 2000. So no selection is needed and the software should work on the 358 haps directly. I will likely look at this in details on Wednesday and will update you here.

Regarding reading directly BAM files, unfortunately no, G1 uses only genotype likelihoods file in VCF/BCF for now (it'll be part of future updates). For low-coverage data, the easiest and fast method to get genotype likelihoods is bcftools.

Just provide the bams to bcftools following these instructions: https://odelaneau.github.io/GLIMPSE/glimpse1/tutorial_b38.html#run_likelihoods

It's relatively simple: bcftools requires a BAM file and does calling at the positions you want: in your case reference panel positions. Only annoying thing is that bcftools requires these positions in two formats: TSV and VCF.gz (this explains 3.1). It also requires the reference genome (fasta) in the right build. It's a bit tedious, but now you should be able to understand step3.2 as well. Then you will need to run 3.3. as well. it's indeed a bit annoying compared to G2, but it should be pretty quick as bcftools is an amazing tool.

I can help if needed. Best,

Simone

[srubinacci]

It's likely a small bug, thanks for reporting it. I think GLIMPSE2 should work fine here. I will look at it very soon and report here. Sorry for the delays.

[srubinacci]

Hi @robertwhbaldwin @karinkumar Thank you for reporting this issue. Indeed, there was a bug in GLIMPSE when reference panels are small. I fixed and pushed the code (it's already online - feel free to pull the code). I will prepare a new release in the next few days containing few bugfixes if you can wait.

[srubinacci]

Prerelease v2.0.1 has been made available, solving this issue. A full release will be coming soon.