Error in Rcpp_get_top_K_or_more_matches_while_building_gamma

[jungleY078]

Hi
I am trying use QUILT to impute my data with a reference panel. But I get in some troubles.

This is my code:

for i in cat list2
do
hap_file="/data/Yujz/quilt-2021/ref/chr${i}.hap.gz"
legend_file="/data/Yujz/quilt-2021/ref/chr${i}.legend.gz"
sample_file="/data/Yujz/quilt-2021/ref/chr${i}.samples"

mkdir impute-chr1-1

R -e "args <- commandArgs(T); library('QUILT'); QUILT(chr=args[1], bamlist='bamlist-test.txt', outputdir='impute-chr1-1', tempdir = '/data/Yujz/quilt-2021/temp' ,reference_haplotype_file='${hap_file}', reference_legend_file='${legend_file}', nGen=120, nCores=1)" --args $i 

done

I met an error like:

[2024-03-03 12:44:40] Imputing sample: 1
[2024-03-03 12:45:52] downsample sample NKX06X40X21947 - 24 of 605786 reads removed
[2024-03-03 12:46:04] The average depth of this sample is:2.56574392710952
[2024-03-03 12:46:04] There are 605762 reads under consideration
[2024-03-03 12:46:08] i_gibbs=1, i_it = 1 small gibbs
[2024-03-03 12:49:39] i_gibbs=1, i_it = 1 full
[2024-03-03 12:50:20] i_gibbs=1, i_it = 2 small gibbs
[2024-03-03 12:53:48] i_gibbs=1, i_it = 2 full
[2024-03-03 12:54:30] i_gibbs=1, i_it = 3 small gibbs
[2024-03-03 12:57:58] i_gibbs=1, i_it = 3 full
[2024-03-03 12:58:57] i_gibbs=2, i_it = 1 small gibbs
[2024-03-03 13:02:22] i_gibbs=2, i_it = 1 full
problem with Rcpp_get_top_K_or_more_matches_while_building_gamma
iGrid = 27067
gamma_t_col(0) = -nan
alphaHat_t(0, iGrid) = -nan
betaHat_t_col(0) = 1
Error: upper value must be greater than lower value
Execution halted

This is an output of one single sample imputation, the depth of this sample is 2.5X. It seems like an error related samples, which means some samples would get in such error, but some else would not. I do not know how to deal with this problem, I have tracked in this error for a long time. Could you give me some advise？

Thanks！

[rwdavies]

Hey, this is an error on QUILT's part, that something has gone wrong with the math / the HMM. I don't see it much anymore, I thought I had sorted out all of those bugs. It's not immediately obvious to me what it is, though one thing, how big is the region you're imputing? 600,000 reads is a LOT of reads for low coverage sequencing! If this is an entire chromosome, can you try imputing in smaller chunks, say 5 Mbp, this will also lead to more accurate imputation, given the underlying assumptions of the model

[jungleY078]

Thanks for your timely reply!
I tried to separate my reference panel into smaller chunks (Chr1 was separated into 9 chunks as each chunks about 10Mbp). Then I ran quilt to impute my samples with these chunks, I got the same error as before. The reads numbers were less than 100000. There were 9 chunks, and the error happened on the second chunk as other chunks went well.

here is the log file：

[2024-03-06 14:04:39] Imputing sample: 1
[2024-03-06 14:04:39] Imputing sample: 2
[2024-03-06 14:04:39] Imputing sample: 3
[2024-03-06 14:05:13] The average depth of this sample is:2.12912785609929
[2024-03-06 14:05:13] There are 64105 reads under consideration
[2024-03-06 14:05:13] i_gibbs=1, i_it = 1 small gibbs
[2024-03-06 14:05:17] The average depth of this sample is:2.5822202322542
[2024-03-06 14:05:17] There are 76233 reads under consideration
[2024-03-06 14:05:17] i_gibbs=1, i_it = 1 small gibbs
[2024-03-06 14:05:27] The average depth of this sample is:3.14347761415989
[2024-03-06 14:05:27] There are 96249 reads under consideration
[2024-03-06 14:05:27] i_gibbs=1, i_it = 1 small gibbs
[2024-03-06 14:05:37] i_gibbs=1, i_it = 1 full
[2024-03-06 14:05:41] i_gibbs=1, i_it = 2 small gibbs
[2024-03-06 14:05:41] i_gibbs=1, i_it = 1 full
[2024-03-06 14:05:47] i_gibbs=1, i_it = 2 small gibbs
[2024-03-06 14:05:57] i_gibbs=1, i_it = 1 full
[2024-03-06 14:06:03] i_gibbs=1, i_it = 2 small gibbs
[2024-03-06 14:06:03] i_gibbs=1, i_it = 2 full
[2024-03-06 14:06:08] i_gibbs=1, i_it = 3 small gibbs
[2024-03-06 14:06:11] i_gibbs=1, i_it = 2 full
[2024-03-06 14:06:15] i_gibbs=1, i_it = 3 small gibbs
[2024-03-06 14:06:29] i_gibbs=1, i_it = 3 full
[2024-03-06 14:06:35] i_gibbs=2, i_it = 1 small gibbs
[2024-03-06 14:06:35] i_gibbs=1, i_it = 2 full
[2024-03-06 14:06:40] i_gibbs=1, i_it = 3 full
[2024-03-06 14:06:40] i_gibbs=1, i_it = 3 small gibbs
[2024-03-06 14:06:46] i_gibbs=2, i_it = 1 small gibbs
[2024-03-06 14:06:57] i_gibbs=2, i_it = 1 full
[2024-03-06 14:07:01] i_gibbs=2, i_it = 2 small gibbs
[2024-03-06 14:07:09] i_gibbs=1, i_it = 3 full
[2024-03-06 14:07:11] i_gibbs=2, i_it = 1 full
problem with Rcpp_get_top_K_or_more_matches_while_building_gamma
iGrid = 3124
gamma_t_col(0) = -nan
alphaHat_t(0, iGrid) = -nan
betaHat_t_col(0) = 1
[2024-03-06 14:07:16] i_gibbs=2, i_it = 1 small gibbs
[2024-03-06 14:07:22] i_gibbs=2, i_it = 2 full
[2024-03-06 14:07:27] i_gibbs=2, i_it = 3 small gibbs
[2024-03-06 14:07:44] i_gibbs=2, i_it = 1 full
[2024-03-06 14:07:48] i_gibbs=2, i_it = 3 full
problem with Rcpp_get_top_K_or_more_matches_while_building_gamma
iGrid = 3124
gamma_t_col(0) = -nan
alphaHat_t(0, iGrid) = -nan
betaHat_t_col(0) = 1
[2024-03-06 14:07:54] i_gibbs=3, i_it = 1 small gibbs
[2024-03-06 14:08:15] i_gibbs=3, i_it = 1 full
[2024-03-06 14:08:20] i_gibbs=3, i_it = 2 small gibbs
[2024-03-06 14:08:41] i_gibbs=3, i_it = 2 full
[2024-03-06 14:08:45] i_gibbs=3, i_it = 3 small gibbs
[2024-03-06 14:09:06] i_gibbs=3, i_it = 3 full
[2024-03-06 14:09:12] i_gibbs=4, i_it = 1 small gibbs
[2024-03-06 14:09:34] i_gibbs=4, i_it = 1 full
[2024-03-06 14:09:38] i_gibbs=4, i_it = 2 small gibbs
[2024-03-06 14:09:58] i_gibbs=4, i_it = 2 full
[2024-03-06 14:10:02] i_gibbs=4, i_it = 3 small gibbs
[2024-03-06 14:10:23] i_gibbs=4, i_it = 3 full
problem with Rcpp_get_top_K_or_more_matches_while_building_gamma
iGrid = 3124
gamma_t_col(0) = -nan
alphaHat_t(0, iGrid) = -nan
betaHat_t_col(0) = 1
[2024-03-06 14:10:27] Error : upper value must be greater than lower value