Skip to content

Mixlen refactoring: Remove dependency on ModelFactoryMixlen and cppoptlib, add support for the -ft option - #182

Open
StefanFlaumberg wants to merge 4 commits into
iqtree:masterfrom
StefanFlaumberg:fixmixlen
Open

Mixlen refactoring: Remove dependency on ModelFactoryMixlen and cppoptlib, add support for the -ft option#182
StefanFlaumberg wants to merge 4 commits into
iqtree:masterfrom
StefanFlaumberg:fixmixlen

Conversation

@StefanFlaumberg

Copy link
Copy Markdown
Contributor

This pull request simplifies the GHOST rate model implementation and allows using mixture models with GHOST rates for PMSF profile estimation on a guide tree. Specifically, it introduces the following changes:

  • Merge the ModelFactoryMixlen class into the ModelFactory class. The former class was just adding unnecessary complexity, its functionality is now seamlessly included into the ModelFactory class.
  • Remove dependency of the PhyloTreeMixlen class on the cppoptlib library and simplify initialization of PhyloTreeMixlen objects. The library has never been compiled and used, its version has long been deprecated, and if needed someday (very unlikely), it can be reimplemented without inheriting from PhyloTreeMixlen. The library required the number of heterotachy categories being set by the PhyloTreeMixlen() constructor (and used the ill-defined -mixlen option for this), this largely contradicted the general logic that implies that the number of categories should be passed within the model string (as -m "...+Hn").
    Now the logic is respected: the PhyloTreeMixlen() constructor takes only the alignment object as an argument, the number of heterotachy categories is set exclusively by the PhyloTreeMixlen::setModelFactory() function as given in the model string, and if a tree is read from a file with not enough or with extra sets of branch lengths, it is corrected by the PhyloTreeMixlen::initializeMixBranches() function as before.
  • With the new PhyloTreeMixlen class initialization rules (see above), mixlen is always 1 for a newly constructed object and is always higher than 1 after setting a RateHeterotachy site rate with setModelFactory(), it can never be 0. Also ASSERTs are added to ensure that a model factory provides a RateHeterotachy site rate (or a nullptr) and that a RateHeterotachy site rate uses a PhyloTreeMixlen tree.
  • Fix PhyloTree::getCurMixture() to return -1 instead of 0. The old output likely didn't affect anything for now, but clearly was a latent bug.
  • Fix and simplify switching to PhyloTreeMixlen in AliSimulator::initializeIQTreeFromTreeFile(). The fix comes from declaring IQTree *&current_tree as a reference to a pointer, not just a pointer in the supertree case -- this way the subtrees are reset, not just deleted in the supertree vector.
  • Use the IQTree or PhyloTreeMixlen object instead of the PhyloTree object for the guide tree in the computeSiteFrequencyModel() function of phyloanalysis.cpp to extend the -ft option to the mixture models that use heterotachy rates. This makes the fix for the Segmentation Fault when using PMSF with GHOST models iqtree2#459 issue complete.
    Note params.min_branch_length being temporarily set in the function: this is needed in the general case, but the need is especially pronounced for the heterotachy rate model, which can be using very short branch lengths (without the minimal length, this results in numerical underflows).
  • Minor: in the runStandardBootstrap() function of phyloanalysis.cpp, change the indices of the .bootlh file lines from always zero to using sample. This is the output of the -wtl option; it looks pretty useless, for the likelihoods are not phylogenetic, but polynomial, but let it at least have the correct indexing now.

These changes are intended not to affect tree and model optimization. The optimized likelihoods should be the same as in the upstream version wherever they can be compared.

@StefanFlaumberg

Copy link
Copy Markdown
Contributor Author

The remaining problems:
There are many. The PhyloTreeMixlen::retoreCheckpoint() function doesn't work at the beginning of the analysis, when mixlen = 1, due to checkpoint key mismatch. The IQTree -> PhyloTreeMixlen switching block, now placed in startTreeReconstruction(), is enough for ModelFinder, but not for MixtureFinder, which runs on tree copies, so the heterotachy rate model remains incompatible with MixtureFinder. The treeLengths() function accumulates values across successive calls, but I haven't touched this, as the fix affects analysis results. The logic behind the readTreeString() function setting relative_treelen while relative_treelen can also be initialized from a temp free rate model is quite complicated and a bit circular.
I opted not to fix these problems here to preserve the original algorithm and likelihoods. I will try to address them in the following PRs.

@StefanFlaumberg

Copy link
Copy Markdown
Contributor Author

On the new changes:

  • Rename getMixlen() to getNMixlen().
  • Merge the save/restoreBranchLengths() functions of PhyloSuperTreeUnlinked into the respective base functions of PhyloSuperTree. For both PhyloTree and PhyloSuperTree classes, make lenvec checking/resizing more robust in these functions (too small vectors now get properly resized).
  • Fix PhyloTree::getNBranchParameters() to properly account for mixlen and rooted trees. Previously, the return of the function was augmented in ModelFactory::getNParameters(), but this gave wrong results for rooted trees and branch scale optimization.

@StefanFlaumberg

Copy link
Copy Markdown
Contributor Author

@bqminh,
The conflict arises due to one newline having been modified in the last PR. If you are not going to propose any changes here, please try to resolve the conflict right within the merge commit if possible (I'd rather not add a dummy newline-removing commit).

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant