Examples with BeagleDataLikelihoodDelegate dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate used on opensource projects

Example 1 with BeagleDataLikelihoodDelegate

use of dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate in project beast-mcmc by beast-dev.

the class TreeDataLikelihoodParser method createTreeDataLikelihood.

protected Likelihood createTreeDataLikelihood(List<PatternList> patternLists, List<BranchModel> branchModels, List<SiteRateModel> siteRateModels, TreeModel treeModel, BranchRateModel branchRateModel, TipStatesModel tipStatesModel, boolean useAmbiguities, PartialsRescalingScheme scalingScheme, boolean delayRescalingUntilUnderflow) throws XMLParseException {
    if (tipStatesModel != null) {
        throw new XMLParseException("Tip State Error models are not supported yet with TreeDataLikelihood");
    }
    List<Taxon> treeTaxa = treeModel.asList();
    List<Taxon> patternTaxa = patternLists.get(0).asList();
    if (!patternTaxa.containsAll(treeTaxa)) {
        throw new XMLParseException("TreeModel contains more taxa than the partition pattern list.");
    }
    if (!treeTaxa.containsAll(patternTaxa)) {
        throw new XMLParseException("TreeModel contains fewer taxa than the partition pattern list.");
    }
    //        DataLikelihoodDelegate dataLikelihoodDelegate = new BeagleDataLikelihoodDelegate(
    //                treeModel,
    //                patternLists.get(0),
    //                branchModel,
    //                siteRateModel,
    //                useAmbiguities,
    //                scalingScheme,
    //                delayRescalingUntilUnderflow);
    boolean useBeagle3 = Boolean.parseBoolean(System.getProperty("USE_BEAGLE3", "true"));
    boolean useJava = Boolean.parseBoolean(System.getProperty("java.only", "false"));
    if (useBeagle3 && MultiPartitionDataLikelihoodDelegate.IS_MULTI_PARTITION_COMPATIBLE() && !useJava) {
        DataLikelihoodDelegate dataLikelihoodDelegate = new MultiPartitionDataLikelihoodDelegate(treeModel, patternLists, branchModels, siteRateModels, useAmbiguities, scalingScheme, delayRescalingUntilUnderflow);
        return new TreeDataLikelihood(dataLikelihoodDelegate, treeModel, branchRateModel);
    } else {
        // The multipartition data likelihood isn't available so make a set of single partition data likelihoods
        List<Likelihood> treeDataLikelihoods = new ArrayList<Likelihood>();
        for (int i = 0; i < patternLists.size(); i++) {
            DataLikelihoodDelegate dataLikelihoodDelegate = new BeagleDataLikelihoodDelegate(treeModel, patternLists.get(i), branchModels.get(i), siteRateModels.get(i), useAmbiguities, scalingScheme, delayRescalingUntilUnderflow);
            treeDataLikelihoods.add(new TreeDataLikelihood(dataLikelihoodDelegate, treeModel, branchRateModel));
        }
        if (treeDataLikelihoods.size() == 1) {
            return treeDataLikelihoods.get(0);
        }
        return new CompoundLikelihood(treeDataLikelihoods);
    }
}

Example 2 with BeagleDataLikelihoodDelegate

use of dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate in project beast-mcmc by beast-dev.

the class CheckPointTreeModifier method incorporateAdditionalTaxa.

/**
     * Add the remaining taxa, which can be identified through the TreeDataLikelihood XML elements.
     */
public ArrayList<NodeRef> incorporateAdditionalTaxa(CheckPointUpdaterApp.UpdateChoice choice, BranchRates rateModel) {
    System.out.println("Tree before adding taxa:\n" + treeModel.toString() + "\n");
    ArrayList<NodeRef> newTaxaNodes = new ArrayList<NodeRef>();
    for (String str : newTaxaNames) {
        for (int i = 0; i < treeModel.getExternalNodeCount(); i++) {
            if (treeModel.getNodeTaxon(treeModel.getExternalNode(i)).getId().equals(str)) {
                newTaxaNodes.add(treeModel.getExternalNode(i));
                //always take into account Taxon dates vs. dates set through a TreeModel
                System.out.println(treeModel.getNodeTaxon(treeModel.getExternalNode(i)).getId() + " with height " + treeModel.getNodeHeight(treeModel.getExternalNode(i)) + " or " + treeModel.getNodeTaxon(treeModel.getExternalNode(i)).getHeight());
            }
        }
    }
    System.out.println("newTaxaNodes length = " + newTaxaNodes.size());
    ArrayList<Taxon> currentTaxa = new ArrayList<Taxon>();
    for (int i = 0; i < treeModel.getExternalNodeCount(); i++) {
        boolean taxonFound = false;
        for (String str : newTaxaNames) {
            if (str.equals((treeModel.getNodeTaxon(treeModel.getExternalNode(i))).getId())) {
                taxonFound = true;
            }
        }
        if (!taxonFound) {
            System.out.println("Adding " + treeModel.getNodeTaxon(treeModel.getExternalNode(i)).getId() + " to list of current taxa");
            currentTaxa.add(treeModel.getNodeTaxon(treeModel.getExternalNode(i)));
        }
    }
    System.out.println("Current taxa count = " + currentTaxa.size());
    //iterate over both current taxa and to be added taxa
    boolean originTaxon = true;
    for (Taxon taxon : currentTaxa) {
        if (taxon.getHeight() == 0.0) {
            originTaxon = false;
            System.out.println("Current taxon " + taxon.getId() + " has node height 0.0");
        }
    }
    for (NodeRef newTaxon : newTaxaNodes) {
        if (treeModel.getNodeTaxon(newTaxon).getHeight() == 0.0) {
            originTaxon = false;
            System.out.println("New taxon " + treeModel.getNodeTaxon(newTaxon).getId() + " has node height 0.0");
        }
    }
    //check the Tree(Data)Likelihoods in the connected set of likelihoods
    //focus on TreeDataLikelihood, which has getTree() to get the tree for each likelihood
    //also get the DataLikelihoodDelegate from TreeDataLikelihood
    ArrayList<TreeDataLikelihood> likelihoods = new ArrayList<TreeDataLikelihood>();
    ArrayList<Tree> trees = new ArrayList<Tree>();
    ArrayList<DataLikelihoodDelegate> delegates = new ArrayList<DataLikelihoodDelegate>();
    for (Likelihood likelihood : Likelihood.CONNECTED_LIKELIHOOD_SET) {
        if (likelihood instanceof TreeDataLikelihood) {
            likelihoods.add((TreeDataLikelihood) likelihood);
            trees.add(((TreeDataLikelihood) likelihood).getTree());
            delegates.add(((TreeDataLikelihood) likelihood).getDataLikelihoodDelegate());
        }
    }
    //suggested to go through TreeDataLikelihoodParser and give it an extra option to create a HashMap
    //keyed by the tree; am currently not overly fond of this approach
    ArrayList<PatternList> patternLists = new ArrayList<PatternList>();
    for (DataLikelihoodDelegate del : delegates) {
        if (del instanceof BeagleDataLikelihoodDelegate) {
            patternLists.add(((BeagleDataLikelihoodDelegate) del).getPatternList());
        } else if (del instanceof MultiPartitionDataLikelihoodDelegate) {
            MultiPartitionDataLikelihoodDelegate mpdld = (MultiPartitionDataLikelihoodDelegate) del;
            List<PatternList> list = mpdld.getPatternLists();
            for (PatternList pList : list) {
                patternLists.add(pList);
            }
        }
    }
    if (patternLists.size() == 0) {
        throw new RuntimeException("No patterns detected. Please make sure the XML file is BEAST 1.9 compatible.");
    }
    //aggregate all patterns to create distance matrix
    //TODO What about different trees for different partitions?
    Patterns patterns = new Patterns(patternLists.get(0));
    if (patternLists.size() > 1) {
        for (int i = 1; i < patternLists.size(); i++) {
            patterns.addPatterns(patternLists.get(i));
        }
    }
    //set the patterns for the distance matrix computations
    choice.setPatterns(patterns);
    //add new taxa one at a time
    System.out.println("Adding " + newTaxaNodes.size() + " taxa ...");
    for (NodeRef newTaxon : newTaxaNodes) {
        treeModel.setNodeHeight(newTaxon, treeModel.getNodeTaxon(newTaxon).getHeight());
        System.out.println("\nadding Taxon: " + newTaxon + " (height = " + treeModel.getNodeHeight(newTaxon) + ")");
        //check if this taxon has a more recent sampling date than all other nodes in the current TreeModel
        double offset = checkCurrentTreeNodes(newTaxon, treeModel.getRoot());
        System.out.println("Sampling date offset when adding " + newTaxon + " = " + offset);
        //AND set its current node height to 0.0 IF no originTaxon has been found
        if (offset < 0.0) {
            if (!originTaxon) {
                System.out.println("Updating all node heights with offset " + Math.abs(offset));
                updateAllTreeNodes(Math.abs(offset), treeModel.getRoot());
                treeModel.setNodeHeight(newTaxon, 0.0);
            }
        } else if (offset == 0.0) {
            if (!originTaxon) {
                treeModel.setNodeHeight(newTaxon, 0.0);
            }
        }
        //get the closest Taxon to the Taxon that needs to be added
        //take into account which taxa can currently be chosen
        Taxon closest = choice.getClosestTaxon(treeModel.getNodeTaxon(newTaxon), currentTaxa);
        System.out.println("\nclosest Taxon: " + closest + " with original height: " + closest.getHeight());
        //get the distance between these two taxa
        double distance = choice.getDistance(treeModel.getNodeTaxon(newTaxon), closest);
        System.out.println("at distance: " + distance);
        //TODO what if distance == 0.0 ??? how to choose closest taxon then (in absence of geo info)?
        //find the NodeRef for the closest Taxon (do not rely on node numbering)
        NodeRef closestRef = null;
        //careful with node numbering and subtract number of new taxa
        for (int i = 0; i < treeModel.getExternalNodeCount(); i++) {
            if (treeModel.getNodeTaxon(treeModel.getExternalNode(i)) == closest) {
                closestRef = treeModel.getExternalNode(i);
            }
        }
        System.out.println(closestRef + " with height " + treeModel.getNodeHeight(closestRef));
        //System.out.println("trying to set node height: " + closestRef + " from " + treeModel.getNodeHeight(closestRef) + " to " + closest.getHeight());
        //treeModel.setNodeHeight(closestRef, closest.getHeight());
        double timeForDistance = distance / rateModel.getBranchRate(treeModel, closestRef);
        System.out.println("timeForDistance = " + timeForDistance);
        //get parent node of branch that will be split
        NodeRef parent = treeModel.getParent(closestRef);
        //determine height of new node
        double insertHeight;
        if (treeModel.getNodeHeight(closestRef) == treeModel.getNodeHeight(newTaxon)) {
            insertHeight = treeModel.getNodeHeight(closestRef) + timeForDistance / 2.0;
            System.out.println("treeModel.getNodeHeight(closestRef) == treeModel.getNodeHeight(newTaxon): " + insertHeight);
            if (insertHeight >= treeModel.getNodeHeight(parent)) {
                insertHeight = treeModel.getNodeHeight(closestRef) + EPSILON * (treeModel.getNodeHeight(parent) - treeModel.getNodeHeight(closestRef));
            }
        } else {
            double remainder = (timeForDistance - Math.abs(treeModel.getNodeHeight(closestRef) - treeModel.getNodeHeight(newTaxon))) / 2.0;
            if (remainder > 0) {
                insertHeight = Math.max(treeModel.getNodeHeight(closestRef), treeModel.getNodeHeight(newTaxon)) + remainder;
                System.out.println("remainder > 0: " + insertHeight);
                if (insertHeight >= treeModel.getNodeHeight(parent)) {
                    insertHeight = treeModel.getNodeHeight(closestRef) + EPSILON * (treeModel.getNodeHeight(parent) - treeModel.getNodeHeight(closestRef));
                }
            } else {
                insertHeight = EPSILON * (treeModel.getNodeHeight(parent) - Math.max(treeModel.getNodeHeight(closestRef), treeModel.getNodeHeight(newTaxon)));
                insertHeight += Math.max(treeModel.getNodeHeight(closestRef), treeModel.getNodeHeight(newTaxon));
                System.out.println("remainder <= 0: " + insertHeight);
            }
        }
        System.out.println("insert at height: " + insertHeight);
        //pass on all the necessary variables to a method that adds the new taxon to the tree
        addTaxonAlongBranch(newTaxon, parent, closestRef, insertHeight);
        //option to print tree after each taxon addition
        System.out.println("\nTree after adding taxon " + newTaxon + ":\n" + treeModel.toString());
        //add newly added Taxon to list of current taxa
        currentTaxa.add(treeModel.getNodeTaxon(newTaxon));
    }
    return newTaxaNodes;
}

Example 3 with BeagleDataLikelihoodDelegate

use of dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate in project beast-mcmc by beast-dev.

the class LocationScaleGradientParser method parseTreeDataLikelihood.

private GradientWrtParameterProvider parseTreeDataLikelihood(XMLObject xo, TreeDataLikelihood treeDataLikelihood, String traitName, boolean useHessian) throws XMLParseException {
    BranchRateModel branchRateModel = treeDataLikelihood.getBranchRateModel();
    DataLikelihoodDelegate delegate = treeDataLikelihood.getDataLikelihoodDelegate();
    if (delegate instanceof ContinuousDataLikelihoodDelegate) {
        throw new XMLParseException("Not yet implemented! ");
    } else if (delegate instanceof BeagleDataLikelihoodDelegate) {
        BeagleDataLikelihoodDelegate beagleData = (BeagleDataLikelihoodDelegate) delegate;
        BranchModel branchModel = beagleData.getBranchModel();
        if (branchRateModel instanceof DefaultBranchRateModel || branchRateModel instanceof ArbitraryBranchRates) {
            if (xo.hasChildNamed(LOCATION)) {
                BranchSpecificFixedEffects location = parseLocation(xo);
                return new LocationGradient(traitName, treeDataLikelihood, beagleData, location, useHessian);
            } else if (xo.hasChildNamed(SCALE)) {
                Parameter scale = (Parameter) xo.getElementFirstChild(SCALE);
                return new ScaleGradient(traitName, treeDataLikelihood, beagleData, scale, useHessian);
            } else {
                throw new XMLParseException("Poorly formed");
            }
        } else if (branchModel instanceof ArbitrarySubstitutionParameterBranchModel) {
            BranchParameter branchParameter = (BranchParameter) xo.getChild(BranchParameter.class);
            if (xo.hasChildNamed(LOCATION)) {
                BranchSpecificFixedEffects location = parseLocation(xo);
                return new BranchSubstitutionParameterLocationGradient(traitName, treeDataLikelihood, beagleData, branchParameter, useHessian, location);
            } else if (xo.hasChildNamed(SCALE)) {
                Parameter scale = (Parameter) xo.getElementFirstChild(SCALE);
                return new BranchSubstitutionParameterScaleGradient(traitName, treeDataLikelihood, beagleData, branchParameter, scale, useHessian);
            } else {
                throw new XMLParseException("Not yet implemented.");
            }
        } else {
            throw new XMLParseException("Only implemented for an arbitrary rates model");
        }
    } else {
        throw new XMLParseException("Unknown likelihood delegate type");
    }
}

Example 4 with BeagleDataLikelihoodDelegate

use of dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate in project beast-mcmc by beast-dev.

the class NodeHeightGradientParser method parseXMLObject.

public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    String traitName = xo.getAttribute(TRAIT_NAME, DEFAULT_TRAIT_NAME);
    final TreeDataLikelihood treeDataLikelihood = (TreeDataLikelihood) xo.getChild(TreeDataLikelihood.class);
    BranchRateModel branchRateModel = treeDataLikelihood.getBranchRateModel();
    if (branchRateModel instanceof DefaultBranchRateModel || branchRateModel instanceof ArbitraryBranchRates) {
        Parameter branchRates = null;
        if (branchRateModel instanceof ArbitraryBranchRates) {
            branchRates = ((ArbitraryBranchRates) branchRateModel).getRateParameter();
        }
        DataLikelihoodDelegate delegate = treeDataLikelihood.getDataLikelihoodDelegate();
        if (delegate instanceof ContinuousDataLikelihoodDelegate) {
            throw new XMLParseException("Not yet implemented! ");
        } else if (delegate instanceof BeagleDataLikelihoodDelegate) {
            BeagleDataLikelihoodDelegate beagleData = (BeagleDataLikelihoodDelegate) delegate;
            return new NodeHeightGradientForDiscreteTrait(traitName, treeDataLikelihood, beagleData, branchRates);
        } else {
            throw new XMLParseException("Unknown likelihood delegate type");
        }
    } else {
        throw new XMLParseException("Only implemented for an arbitrary rates model");
    }
}

Example 5 with BeagleDataLikelihoodDelegate

use of dr.evomodel.treedatalikelihood.BeagleDataLikelihoodDelegate in project beast-mcmc by beast-dev.

the class TreeDataLikelihoodParser method createTreeDataLikelihood.

protected Likelihood createTreeDataLikelihood(List<PatternList> patternLists, List<BranchModel> branchModels, List<SiteRateModel> siteRateModels, Tree treeModel, BranchRateModel branchRateModel, TipStatesModel tipStatesModel, boolean useAmbiguities, boolean preferGPU, PartialsRescalingScheme scalingScheme, boolean delayRescalingUntilUnderflow, PreOrderSettings settings) throws XMLParseException {
    if (tipStatesModel != null) {
        throw new XMLParseException("Tip State Error models are not supported yet with TreeDataLikelihood");
    }
    List<Taxon> treeTaxa = treeModel.asList();
    List<Taxon> patternTaxa = patternLists.get(0).asList();
    if (!patternTaxa.containsAll(treeTaxa)) {
        throw new XMLParseException("TreeModel " + treeModel.getId() + " contains more taxa (" + treeModel.getExternalNodeCount() + ") than the partition pattern list (" + patternTaxa.size() + ").");
    }
    if (!treeTaxa.containsAll(patternTaxa)) {
        throw new XMLParseException("TreeModel " + treeModel.getId() + " contains fewer taxa (" + treeModel.getExternalNodeCount() + ") than the partition pattern list (" + patternTaxa.size() + ").");
    }
    boolean useBeagle3MultiPartition = false;
    if (patternLists.size() > 1) {
        // will currently recommend true if using GPU, CUDA or OpenCL.
        useBeagle3MultiPartition = MultiPartitionDataLikelihoodDelegate.IS_MULTI_PARTITION_RECOMMENDED();
        if (System.getProperty("USE_BEAGLE3_EXTENSIONS") != null) {
            useBeagle3MultiPartition = Boolean.parseBoolean(System.getProperty("USE_BEAGLE3_EXTENSIONS"));
        }
        if (System.getProperty("beagle.multipartition.extensions") != null && !System.getProperty("beagle.multipartition.extensions").equals("auto")) {
            useBeagle3MultiPartition = Boolean.parseBoolean(System.getProperty("beagle.multipartition.extensions"));
        }
    }
    boolean useJava = Boolean.parseBoolean(System.getProperty("java.only", "false"));
    int threadCount = -1;
    int beagleThreadCount = -1;
    if (System.getProperty(BEAGLE_THREAD_COUNT) != null) {
        beagleThreadCount = Integer.parseInt(System.getProperty(BEAGLE_THREAD_COUNT));
    }
    if (beagleThreadCount == -1) {
        if (System.getProperty(THREAD_COUNT) != null) {
            threadCount = Integer.parseInt(System.getProperty(THREAD_COUNT));
        }
    }
    if (useBeagle3MultiPartition && !useJava) {
        if (beagleThreadCount == -1 && threadCount >= 0) {
            System.setProperty(BEAGLE_THREAD_COUNT, Integer.toString(threadCount));
        }
        try {
            DataLikelihoodDelegate dataLikelihoodDelegate = new MultiPartitionDataLikelihoodDelegate(treeModel, patternLists, branchModels, siteRateModels, useAmbiguities, scalingScheme, delayRescalingUntilUnderflow);
            return new TreeDataLikelihood(dataLikelihoodDelegate, treeModel, branchRateModel);
        } catch (DataLikelihoodDelegate.DelegateTypeException dte) {
            useBeagle3MultiPartition = false;
        }
    }
    // The multipartition data likelihood isn't available so make a set of single partition data likelihoods
    List<Likelihood> treeDataLikelihoods = new ArrayList<Likelihood>();
    // Todo: allow for different number of threads per beagle instance according to pattern counts
    if (beagleThreadCount == -1 && threadCount >= 0) {
        System.setProperty(BEAGLE_THREAD_COUNT, Integer.toString(threadCount / patternLists.size()));
    }
    for (int i = 0; i < patternLists.size(); i++) {
        DataLikelihoodDelegate dataLikelihoodDelegate = new BeagleDataLikelihoodDelegate(treeModel, patternLists.get(i), branchModels.get(i), siteRateModels.get(i), useAmbiguities, preferGPU, scalingScheme, delayRescalingUntilUnderflow, settings);
        treeDataLikelihoods.add(new TreeDataLikelihood(dataLikelihoodDelegate, treeModel, branchRateModel));
    }
    if (treeDataLikelihoods.size() == 1) {
        return treeDataLikelihoods.get(0);
    }
    return new CompoundLikelihood(treeDataLikelihoods);
}