SKATO: rvtest frozen/stalled (3 variants) (minimal example provided)

[lindenb]

Hi, find in attachement where rvtest seems stalled/frozen (skato test ?) while there is a small number of variant.

$ unzip -l test.zip 
Archive:  test.zip
  Length      Date    Time    Name
---------  ---------- -----   ----
    34131  2022-01-24 13:59   test.ped
       38  2022-01-25 13:22   test.set
      218  2022-01-25 20:46   test.sh
     4150  2022-01-25 20:43   test.vcf.gz
      283  2022-01-25 20:43   test.vcf.gz.tbi
---------                     -------
    38820                     5 files

the test.sh contains:

#/!bin/sh
set -x
set -e
rvtest --noweb --inVcf test.vcf.gz \
        --pheno  test.ped --out chr3 \
        --setFile test.set \
        --burden cmc,zeggini,mb,fp,exactCMC,cmcWald,rarecover,cmat \
        --vt price,analytic \
        --kernel 'skato'

output:

[INFO]  Program version: 20190205      
[INFO]  Analysis started at: Tue Jan 25 20:52:45 2022                                                           
[INFO]  Loaded [ 1421 ] samples from genotype files                                                             
[INFO]  Loaded [ 1421 ] sample phenotypes                                                                       
[INFO]  Loaded 0 male, 0 female and 1421 sex-unknown samples from test.ped                                      
[INFO]  Loaded 351 cases, 1070 controls, and 0 missing phenotypes           
[WARN]  -- Enabling binary phenotype mode --                                                                    
[INFO]  Analysis begins with [ 1421 ] samples...                                                                
[INFO]  MadsonBrowning test significance will be evaluated using 10000 permutations at alpha = 0.05
[INFO]  Rare cover test significance will be evaluated using 10000 permutations at alpha = 0.05
[INFO]  cmat test significance will be evaluated using 10000 permutations at alpha = 0.05                       
[INFO]  Price's VT test significance will be evaluated using 10000 permutations at alpha = 0.05
[INFO]  SKAT-O test significance will be evaluated using weight = Beta[beta1 = 1.00, beta2 = 25.00]             
[INFO]  Loaded [ 1 ] set to tests.                                                                              
[INFO]  Impute missing genotype to mean (by default)                                                            
[INFO]  Analysis started                                                                                        
Matrix() called 1 times                                                                                         
[WARN]  Analytic VT test does not support binary outcomes. Results will be all NAs.                             
<===================FROZEN#########################

would it be a bug ? how can I prevent rvtest from freezing ? thanks !

test.zip

[lindenb]

as far as I can see the error is created by a NAN in MixtureChiSquare.cpp:getLiuPvalue

cdfchn(&which, &p, &q, &x, &l, &ncp, &status, &bound);

with X is NAN

the NAN value is created by a sqrt with a negative value sqrt(VarQ) in computeIntegrandDavies

modifying if (kappa > lambda.sum() * 10000) { with if (kappa > lambda.sum() * 10000 || VarQ <=0.0) { seems to fix the error but I don't know if it is the correct solution. Nevertheles, I'll submit a PR.