Search in sources :

Example 96 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class LogCombiner method rescaleTree.

private void rescaleTree(Tree tree, double scale) {
    if (tree instanceof MutableTree) {
        MutableTree mutableTree = (MutableTree) tree;
        for (int i = 0; i < tree.getNodeCount(); i++) {
            NodeRef node = tree.getNode(i);
            if (node != tree.getRoot()) {
                double length = tree.getBranchLength(node);
                mutableTree.setBranchLength(node, length * scale);
            }
        }
    } else {
        throw new IllegalArgumentException("Tree not mutable");
    }
}
Also used : MutableTree(dr.evolution.tree.MutableTree) NodeRef(dr.evolution.tree.NodeRef)

Example 97 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class DensityMap method addTree.

public void addTree(Tree tree, double sampleTime, String attributeName1, String attributeName2) {
    checkCalibration();
    double[][] variance = null;
    Object[] obj = (Object[]) tree.getAttribute(MultivariateDiffusionModel.PRECISION_TREE_ATTRIBUTE);
    if (obj != null) {
        variance = new Matrix(MatrixParameter.parseFromSymmetricDoubleArray(obj).getParameterAsMatrix()).inverse().toComponents();
    }
    for (int i = 0; i < tree.getNodeCount(); i++) {
        NodeRef node = tree.getNode(i);
        if (node != tree.getRoot()) {
            NodeRef parent = tree.getParent(node);
            double t1 = tree.getNodeHeight(node);
            double t2 = tree.getNodeHeight(parent);
            if (t1 <= sampleTime && t2 >= sampleTime) {
                Double valueX1 = transform((Double) tree.getNodeAttribute(node, attributeName1));
                Double valueY1 = transform((Double) tree.getNodeAttribute(node, attributeName2));
                Double valueX2 = transform((Double) tree.getNodeAttribute(parent, attributeName1));
                Double valueY2 = transform((Double) tree.getNodeAttribute(parent, attributeName2));
                if (valueX1 != null && valueY1 != null && valueX2 != null && valueY2 != null) {
                    addPoint(sampleTime, t1, t2, valueX1, valueY1, valueX2, valueY2, variance);
                }
            }
        }
    }
}
Also used : NodeRef(dr.evolution.tree.NodeRef) Matrix(dr.math.matrixAlgebra.Matrix)

Example 98 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class DensityMap method calibrate.

public void calibrate(Tree tree, String attributeName1, String attributeName2) {
    boolean foundAttribute1 = false;
    boolean foundAttribute2 = false;
    jointDensity = true;
    if (isCalibrated) {
        throw new RuntimeException("Already calibrated");
    }
    //		}
    for (int i = 0; i < tree.getNodeCount(); i++) {
        NodeRef node = tree.getNode(i);
        if (node != tree.getRoot()) {
            Double value = (Double) tree.getNodeAttribute(node, attributeName1);
            if (value != null) {
                value = transform(value);
                if (value < minX)
                    minX = value;
                if (value > maxX)
                    maxX = value;
                foundAttribute1 = true;
            }
            value = (Double) tree.getNodeAttribute(node, attributeName2);
            if (value != null) {
                value = transform(value);
                if (value < minY)
                    minY = value;
                if (value > maxY)
                    maxY = value;
                foundAttribute2 = true;
            }
        }
    }
    if (!foundAttribute1) {
        throw new RuntimeException("Can't find any attributes, " + attributeName1 + ", in tree " + tree.getId());
    }
    if (!foundAttribute2) {
        throw new RuntimeException("Can't find any attributes, " + attributeName2 + ", in tree " + tree.getId());
    }
//		System.err.printf("Calibrated: minY = %3.2f, maxY = %3.2f, minX = %3.2f, maxX = %3.2f\n",minY,maxY,minX,maxX);
//		System.exit(-1);
}
Also used : NodeRef(dr.evolution.tree.NodeRef)

Example 99 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class TaxaOriginTrait method getTraitName.

private String getTraitName() {
    FlexibleTree tree = trees[0];
    String traitName = null;
    for (int i = 0; i < tree.getExternalNodeCount(); i++) {
        NodeRef node = tree.getExternalNode(i);
        for (String taxaName : taxaNames) {
            if (tree.getNodeTaxon(node).getId().equals(taxaName)) {
                String attributeValue = (String) tree.getNodeAttribute(node, attributeName);
                if (traitName != null && !traitName.equals(attributeValue)) {
                    throw new RuntimeException("Not all taxa given have the same trait value");
                } else {
                    traitName = attributeValue;
                }
            }
        }
    }
    return traitName;
}
Also used : NodeRef(dr.evolution.tree.NodeRef) FlexibleTree(dr.evolution.tree.FlexibleTree)

Example 100 with NodeRef

use of dr.evolution.tree.NodeRef in project beast-mcmc by beast-dev.

the class DiffusionRateStatistic method getStatisticValue.

public double getStatisticValue(int dim) {
    String traitName = traitLikelihoods.get(0).getTraitName();
    double treelength = 0;
    double treeDistance = 0;
    double maxDistanceFromRoot = 0;
    double maxDistanceOverTimeFromRoot = 0;
    //double[] rates =  null;
    List<Double> rates = new ArrayList<Double>();
    //double[] diffusionCoefficients =  null;
    List<Double> diffusionCoefficients = new ArrayList<Double>();
    double waDiffusionCoefficient = 0;
    double lowerHeight = heightLowers[dim];
    double upperHeight = Double.MAX_VALUE;
    if (heightLowers.length == 1) {
        upperHeight = heightUpper;
    } else {
        if (dim > 0) {
            if (!cumulative) {
                upperHeight = heightLowers[dim - 1];
            }
        }
    }
    for (AbstractMultivariateTraitLikelihood traitLikelihood : traitLikelihoods) {
        MultivariateTraitTree tree = traitLikelihood.getTreeModel();
        BranchRateModel branchRates = traitLikelihood.getBranchRateModel();
        for (int i = 0; i < tree.getNodeCount(); i++) {
            NodeRef node = tree.getNode(i);
            if (node != tree.getRoot()) {
                NodeRef parentNode = tree.getParent(node);
                if ((tree.getNodeHeight(parentNode) > lowerHeight) && (tree.getNodeHeight(node) < upperHeight)) {
                    double[] trait = traitLikelihood.getTraitForNode(tree, node, traitName);
                    double[] parentTrait = traitLikelihood.getTraitForNode(tree, parentNode, traitName);
                    double[] traitUp = parentTrait;
                    double[] traitLow = trait;
                    double timeUp = tree.getNodeHeight(parentNode);
                    double timeLow = tree.getNodeHeight(node);
                    double rate = (branchRates != null ? branchRates.getBranchRate(tree, node) : 1.0);
                    MultivariateDiffusionModel diffModel = traitLikelihood.diffusionModel;
                    double[] precision = diffModel.getPrecisionParameter().getParameterValues();
                    if (tree.getNodeHeight(parentNode) > upperHeight) {
                        timeUp = upperHeight;
                        //TODO: implement TrueNoise??
                        traitUp = imputeValue(trait, parentTrait, upperHeight, tree.getNodeHeight(node), tree.getNodeHeight(parentNode), precision, rate, false);
                    }
                    if (tree.getNodeHeight(node) < lowerHeight) {
                        timeLow = lowerHeight;
                        traitLow = imputeValue(trait, parentTrait, lowerHeight, tree.getNodeHeight(node), tree.getNodeHeight(parentNode), precision, rate, false);
                    }
                    double time = timeUp - timeLow;
                    treelength += time;
                    double[] rootTrait = traitLikelihood.getTraitForNode(tree, tree.getRoot(), traitName);
                    if (useGreatCircleDistances && (trait.length == 2)) {
                        // Great Circle distance
                        SphericalPolarCoordinates coord1 = new SphericalPolarCoordinates(traitLow[0], traitLow[1]);
                        SphericalPolarCoordinates coord2 = new SphericalPolarCoordinates(traitUp[0], traitUp[1]);
                        double distance = coord1.distance(coord2);
                        treeDistance += distance;
                        double dc = Math.pow(distance, 2) / (4 * time);
                        diffusionCoefficients.add(dc);
                        waDiffusionCoefficient += dc * time;
                        rates.add(distance / time);
                        SphericalPolarCoordinates rootCoord = new SphericalPolarCoordinates(rootTrait[0], rootTrait[1]);
                        double tempDistanceFromRoot = rootCoord.distance(coord2);
                        if (tempDistanceFromRoot > maxDistanceFromRoot) {
                            maxDistanceFromRoot = tempDistanceFromRoot;
                            maxDistanceOverTimeFromRoot = tempDistanceFromRoot / (tree.getNodeHeight(tree.getRoot()) - timeLow);
                            //distance between traitLow and traitUp for maxDistanceFromRoot
                            if (timeUp == upperHeight) {
                                maxDistanceFromRoot = distance;
                                maxDistanceOverTimeFromRoot = distance / time;
                            }
                        }
                    } else {
                        double distance = getNativeDistance(traitLow, traitUp);
                        treeDistance += distance;
                        double dc = Math.pow(distance, 2) / (4 * time);
                        diffusionCoefficients.add(dc);
                        waDiffusionCoefficient += dc * time;
                        rates.add(distance / time);
                        double tempDistanceFromRoot = getNativeDistance(traitLow, rootTrait);
                        if (tempDistanceFromRoot > maxDistanceFromRoot) {
                            maxDistanceFromRoot = tempDistanceFromRoot;
                            maxDistanceOverTimeFromRoot = tempDistanceFromRoot / (tree.getNodeHeight(tree.getRoot()) - timeLow);
                            //distance between traitLow and traitUp for maxDistanceFromRoot
                            if (timeUp == upperHeight) {
                                maxDistanceFromRoot = distance;
                                maxDistanceOverTimeFromRoot = distance / time;
                            }
                        }
                    }
                }
            }
        }
    }
    if (summaryStat == summaryStatistic.DIFFUSION_RATE) {
        if (summaryMode == Mode.AVERAGE) {
            return DiscreteStatistics.mean(toArray(rates));
        } else if (summaryMode == Mode.MEDIAN) {
            return DiscreteStatistics.median(toArray(rates));
        } else if (summaryMode == Mode.COEFFICIENT_OF_VARIATION) {
            // don't compute mean twice
            final double mean = DiscreteStatistics.mean(toArray(rates));
            return Math.sqrt(DiscreteStatistics.variance(toArray(rates), mean)) / mean;
        } else {
            return treeDistance / treelength;
        }
    } else if (summaryStat == summaryStatistic.DIFFUSION_COEFFICIENT) {
        if (summaryMode == Mode.AVERAGE) {
            return DiscreteStatistics.mean(toArray(diffusionCoefficients));
        } else if (summaryMode == Mode.MEDIAN) {
            return DiscreteStatistics.median(toArray(diffusionCoefficients));
        } else if (summaryMode == Mode.COEFFICIENT_OF_VARIATION) {
            // don't compute mean twice
            final double mean = DiscreteStatistics.mean(toArray(diffusionCoefficients));
            return Math.sqrt(DiscreteStatistics.variance(toArray(diffusionCoefficients), mean)) / mean;
        } else {
            return waDiffusionCoefficient / treelength;
        }
    } else if (summaryStat == summaryStatistic.WAVEFRONT_DISTANCE) {
        return maxDistanceFromRoot;
    } else {
        return maxDistanceOverTimeFromRoot;
    }
}
Also used : SphericalPolarCoordinates(dr.geo.math.SphericalPolarCoordinates) MultivariateTraitTree(dr.evolution.tree.MultivariateTraitTree) NodeRef(dr.evolution.tree.NodeRef) BranchRateModel(dr.evomodel.branchratemodel.BranchRateModel)

Aggregations

NodeRef (dr.evolution.tree.NodeRef)426 ArrayList (java.util.ArrayList)38 LinkedList (java.util.LinkedList)20 Taxon (dr.evolution.util.Taxon)18 Tree (dr.evolution.tree.Tree)14 TreeModel (dr.evomodel.tree.TreeModel)14 Parameter (dr.inference.model.Parameter)14 Clade (dr.evolution.tree.Clade)13 MutableTree (dr.evolution.tree.MutableTree)9 Node (dr.evomodel.arg.ARGModel.Node)9 MultivariateTraitTree (dr.evolution.tree.MultivariateTraitTree)8 BranchMapModel (dr.evomodel.epidemiology.casetocase.BranchMapModel)8 BitSet (java.util.BitSet)8 FixedBitSet (jebl.util.FixedBitSet)8 FlexibleTree (dr.evolution.tree.FlexibleTree)7 AbstractCase (dr.evomodel.epidemiology.casetocase.AbstractCase)7 Matrix (dr.math.matrixAlgebra.Matrix)6 CompoundParameter (dr.inference.model.CompoundParameter)5 TimeScale (dr.evolution.util.TimeScale)4 MicrosatelliteSamplerTreeModel (dr.evomodel.tree.MicrosatelliteSamplerTreeModel)4