veg / hyphy

HyPhy: Hypothesis testing using Phylogenies
http://www.hyphy.org
Other
205 stars 69 forks source link

GARD execution error #1739

Open liamxg opened 2 weeks ago

liamxg commented 2 weeks ago

Dear @spond,

When I run GARD for my dataset, it executes error as below:
Analysis Description
--------------------
  GARD : Genetic Algorithms for Recombination Detection. Implements a
heuristic approach to screening alignments of sequences for
recombination, by using the CHC genetic algorithm to search for
phylogenetic incongruence among different partitions of the data. The
number of partitions is determined using a step-up procedure, while the
placement of breakpoints is searched for with the GA. The best fitting
model (based on c-AIC) is returned; and additional post-hoc tests run to
distinguish topological incongruence from rate-variation. v0.2 adds and
spooling results to JSON after each breakpoint search conclusion

- __Requirements__: A sequence alignment.

- __Citation__: **Automated Phylogenetic Detection of Recombination Using a Genetic
Algorithm**, _Mol Biol Evol 23(10), 1891–1901

- __Written by__: Sergei L Kosakovsky Pond

- __Contact Information__: spond@temple.edu

- __Analysis Version__: 0.2

type: nucleotide
rv: GDD
>Maximum number of breakpoints to consider (permissible range = [1,100000], default value = 10000, integer): max-breakpoints: 10000
>How many site rate classes to use (permissible range = [2,10], default value = 4, integer): rate-classes: 4
mode: Normal
>Loaded a nucleotide multiple sequence alignment with **16** sequences, **864** sites (125 of which are variable) from `HIV-sets.nex`
>Minimum size of a partition is set to be 29 sites

### Fitting the baseline (single-partition; no breakpoints) model
* Log(L) = -2081.41, AIC-c =  4255.66 (44 estimated parameters)

### Performing an exhaustive single breakpoint analysis
[GARD] Breakpoint         60 of 124. Best cAIC =    4193.5157 [delta =      62.1480] with breakpoint at site        571Error:
  Internal error 
Internal error in _LikelihoodFunction::ConjugateGradientDescent. The function evaluated at current parameter values [-2027.698029793105] does not match the last recorded LF maximum [-2024.242028554538]

Function call stack
1 :  [namespace = HZNq_CHs] Optimize(mles, ^lf_id);

Keyword arguments:
  {
    "model":"GTR"
  }
-------
  2 :  [namespace = uGGLqrzj] res=estimators.FitExistingLF(&likelihoodFunction,modelObjects);

Keyword arguments:
  {
    "model":"GTR"
  }
-------
  3 :  [namespace = iMODMRzA] fit=gard.fitPartitionedModel(breakPoints,model,initialValues,null,FALSE);

Keyword arguments:
  {
    "model":"GTR"
  }
-------
  4 :  [namespace = mpi.KUvmxtbj] Call(result_callback,0,Eval(job+'('+Join(",",utility.Map(arguments,"_value_","utility.convertToArgumentString (_value_)"))+')'),arguments);

Keyword arguments:
  {
    "model":"GTR"
  }
-------
  5 :  [namespace = gard] mpi.QueueJob(queue,"gard.obtainModel_cAIC",{"0":{{siteIndex__}},"1":model,"2":baseLikelihoodInfo},"gard.storeSingleBreakPointModelResults");

Keyword arguments:
  {
    "model":"GTR"
  }
-------
  6 :  namespace 

Step 0.singleBreakPointBest_cAIC=^"math.Infinity";

Step 1.breakPointIndex=0;

Step 2.Branch under condition 'breakPointIndex<variableSites-1'
to
siteIndex=variableSiteMap[breakPointIndex];
else
  mpi.QueueComplete(queue);

Step 3.siteIndex=variableSiteMap[breakPointIndex];

Step 4.Branch under condition 'singleBreakPointBest_cAIC<baseline_cAIC'
to
io.ReportProgressBar("GARD","Breakpoint "+Format(1+breakPointIndex,10,0)+" of "+(variableSites-1)+". Best cAIC = "+Format(singleBreakPointBest_cAIC,12,4)+" [delta = "+Format(baseline_cAIC-singleBreakPointBest_cAIC,12,4)+"] with breakpoint at site "+Format(singleBreakPointBestLocation,10,0));
else
  io.ReportProgressBar("GARD","Breakpoint "+Format(1+breakPointIndex,10,0)+" of "+(variableSites-1)+". Best cAIC = "+Format(baseline_cAIC,12,4)+" with no breakpoints.");

Step 5.io.ReportProgressBar("GARD","Breakpoint "+Format(1+breakPointIndex,10,0)+" of "+(variableSites-1)+". Best cAIC = "+Format(singleBreakPointBest_cAIC,12,4)+" [delta = "+Format(baseline_cAIC-singleBreakPointBest_cAIC,12,4)+"] with breakpoint at site "+Format(singleBreakPointBestLocation,10,0));

Step 6.Go to step 8

Step 7.io.ReportProgressBar("GARD","Breakpoint "+Format(1+breakPointIndex,10,0)+" of "+(variableSites-1)+". Best cAIC = "+Format(baseline_cAIC,12,4)+" with no breakpoints.");

Step 8.Branch under condition 'gard.validatePartititon({{siteIndex}},minPartitionSize,numSites)'
to
mpi.QueueJob(queue,"gard.obtainModel_cAIC",{"0":{{siteIndex__}},"1":model,"2":baseLikelihoodInfo},"gard.storeSingleBreakPointModelResults");
else
  breakPointIndex+=1;

Step 9.mpi.QueueJob(queue,"gard.obtainModel_cAIC",{"0":{{siteIndex__}},"1":model,"2":baseLikelihoodInfo},"gard.storeSingleBreakPointModelResults");

Step 10.breakPointIndex+=1;

Step 11.Go to step 2

Step 12.mpi.QueueComplete(queue);

Step 13.io.ClearProgressBar();

Step 14.io.ReportProgressMessageMD('GARD','single-breakpoint','Done with single breakpoint analysis.');

Step 15.io.ReportProgressMessageMD('GARD','single-breakpoint',("   Best sinlge break point location: "+singleBreakPointBestLocation));

Step 16.io.ReportProgressMessageMD('GARD','single-breakpoint',("   c-AIC  = "+singleBreakPointBest_cAIC));;

Keyword arguments:
  {
    "model":"GTR"
  }
-------

  Check errors.log for execution error details.
spond commented 2 weeks ago

Dear @liamxg,

Add ENV="TOLERATE_NUMERICAL_ERRORS=1;" to the command line.

What's your HyPhy version (hyphy --version)?

Also, make sure you run the HYPHYMPI version for GARD. If you have openmpi installed (and it's free to install), you can use something like

mpirun -np N HYPHYMPI gard...

where N is the number of cores on your system.

Best, Sergei

liamxg commented 2 weeks ago

Dear @spond,

hyphy --version
HYPHY 2.5.48(MP) for Darwin on x86_64

and my MAC OS:

image
spond commented 2 weeks ago

Dear @liamxg,

Could you please update to 2.5.62, the latest HyPhy release?

You can run GARD with 8 threads on your CPU.

Best, Sergei

liamxg commented 2 weeks ago

Dear Prof. @spond,

My GARD can run now, thanks.

Could I know one of the outputs small.fas.best-gard.fit.bf can be used for? Thanks.

Best, Liam