Search in sources :

Example 41 with RealParameter

use of beast.core.parameter.RealParameter in project bacter by tgvaughan.

the class SkylinePopulationFunctionTest method test2.

@Test
public void test2() throws Exception {
    SkylinePopulationFunction skyline = new SkylinePopulationFunction();
    skyline.initByName("acg", acg, "popSizes", new RealParameter("5.0 1.0 5.0 1.0"), "groupSizes", new IntegerParameter("0 0 0 0"));
    for (CFEventList.Event cfEvent : acg.getCFEvents()) {
        double t = cfEvent.getHeight();
        assertTrue(Math.abs(t - skyline.getInverseIntensity(skyline.getIntensity(t))) < 1e-14);
    }
}
Also used : IntegerParameter(beast.core.parameter.IntegerParameter) CFEventList(bacter.CFEventList) RealParameter(beast.core.parameter.RealParameter) Test(org.junit.Test)

Example 42 with RealParameter

use of beast.core.parameter.RealParameter in project bacter by tgvaughan.

the class SkylinePopulationFunctionTest method test1.

@Test
public void test1() throws Exception {
    SkylinePopulationFunction skyline = new SkylinePopulationFunction();
    skyline.initByName("acg", acg, "popSizes", new RealParameter("5.0 1.0 5.0 1.0"), "groupSizes", new IntegerParameter("0 0 0 0"));
    for (double t = 0.0; t < 10; t += 0.01) assertTrue(Math.abs(t - skyline.getInverseIntensity(skyline.getIntensity(t))) < 1e-14);
}
Also used : IntegerParameter(beast.core.parameter.IntegerParameter) RealParameter(beast.core.parameter.RealParameter) Test(org.junit.Test)

Example 43 with RealParameter

use of beast.core.parameter.RealParameter in project bacter by tgvaughan.

the class AddRemoveConversionTest method testHR.

/**
 * Tests that probability density of forward move calculated
 * by drawNewConversion() matches probability density of backward
 * move calculated by getConversionProb().
 *
 * @throws Exception
 */
@Test
public void testHR() throws Exception {
    ConstantPopulation popFunc = new ConstantPopulation();
    popFunc.initByName("popSize", new RealParameter("1.0"));
    Locus locus = new Locus("locus", 10000);
    TaxonSet taxonSet = getTaxonSet(10);
    SimulatedACG acg = new SimulatedACG();
    acg.initByName("rho", 1.0 / locus.getSiteCount(), "delta", 50.0, "locus", locus, "taxonset", taxonSet, "populationModel", popFunc);
    AddRemoveConversion operator = new AddRemoveConversion();
    // Loop until a valid proposal is made
    double logP1;
    List<Conversion> oldConversions;
    do {
        operator.initByName("weight", 1.0, "acg", acg, "delta", new RealParameter("50.0"), "populationModel", popFunc);
        oldConversions = Lists.newArrayList(acg.getConversions(locus));
        logP1 = operator.drawNewConversion();
    } while (Double.isInfinite(logP1));
    System.out.println("logP1 = " + logP1);
    // Identify new recomination
    Conversion newRecomb = null;
    for (Conversion recomb : acg.getConversions(locus)) {
        if (!oldConversions.contains(recomb))
            newRecomb = recomb;
    }
    assertNotNull(newRecomb);
    double logP2 = operator.getConversionProb(newRecomb);
    System.out.println("logP2 = " + logP2);
    assertTrue(Math.abs(logP1 - logP2) < 1e-10);
}
Also used : ConstantPopulation(beast.evolution.tree.coalescent.ConstantPopulation) SimulatedACG(bacter.model.SimulatedACG) RealParameter(beast.core.parameter.RealParameter) Locus(bacter.Locus) TaxonSet(beast.evolution.alignment.TaxonSet) Conversion(bacter.Conversion) Test(org.junit.Test)

Example 44 with RealParameter

use of beast.core.parameter.RealParameter in project bacter by tgvaughan.

the class ACGScaler method proposal.

@Override
public double proposal() {
    // Keep track of number of positively scaled elements minus
    // negatively scaled elements.
    int count = 0;
    // Choose scaling factor:
    double f = scaleParam + Randomizer.nextDouble() * (1.0 / scaleParam - scaleParam);
    // Scale clonal frame:
    if (rootOnly) {
        acg.getRoot().setHeight(acg.getRoot().getHeight() * f);
        count += 1;
    } else {
        for (Node node : acg.getInternalNodes()) {
            node.setHeight(node.getHeight() * f);
            count += 1;
        }
    }
    // Scale recombinant edges:
    for (Locus locus : acg.getLoci()) {
        for (Conversion recomb : acg.getConversions(locus)) {
            if (!rootOnly || recomb.getNode1().getParent().isRoot()) {
                recomb.setHeight1(recomb.getHeight1() * f);
                count += 1;
            }
            if (!rootOnly || recomb.getNode2().isRoot() || recomb.getNode2().getParent().isRoot()) {
                recomb.setHeight2(recomb.getHeight2() * f);
                count += 1;
            }
            if (recomb.getHeight1() < recomb.getNode1().getHeight() || recomb.getHeight2() < recomb.getNode2().getHeight() || recomb.getHeight1() > recomb.getHeight2())
                return Double.NEGATIVE_INFINITY;
        }
    }
    // Check for illegal node heights:
    if (rootOnly) {
        for (Node node : acg.getRoot().getChildren()) {
            if (node.getHeight() > node.getParent().getHeight())
                return Double.NEGATIVE_INFINITY;
        }
    } else {
        for (Node node : acg.getExternalNodes()) {
            if (node.getHeight() > node.getParent().getHeight())
                return Double.NEGATIVE_INFINITY;
        }
    }
    for (RealParameter param : parametersInput.get()) {
        try {
            param.startEditing(null);
            param.scale(f);
        } catch (Exception e) {
            // bounds.  Needs to change!
            return Double.NEGATIVE_INFINITY;
        }
        count += param.getDimension();
    }
    for (RealParameter paramInv : parametersInverseInput.get()) {
        try {
            paramInv.startEditing(null);
            paramInv.scale(1.0 / f);
        } catch (Exception e) {
            // bounds.  Needs to change!
            return Double.NEGATIVE_INFINITY;
        }
        count -= paramInv.getDimension();
    }
    // Return log of Hastings ratio:
    return (count - 2) * Math.log(f);
}
Also used : Node(beast.evolution.tree.Node) RealParameter(beast.core.parameter.RealParameter) Locus(bacter.Locus) Conversion(bacter.Conversion)

Example 45 with RealParameter

use of beast.core.parameter.RealParameter in project bacter by tgvaughan.

the class AddRemoveConversion method main.

public static void main(String[] args) throws Exception {
    ConversionGraph acg = new ConversionGraph();
    ConstantPopulation popFunc = new ConstantPopulation();
    AddRemoveConversion operator = new AddRemoveConversion();
    operator.initByName("weight", 1.0, "acg", acg, "populationModel", popFunc, "rho", new RealParameter(Double.toString(1.0 / 10000.0)), "delta", new RealParameter("50.0"));
    popFunc.initByName("popSize", new RealParameter("1.0"));
    TaxonSet taxonSet = new TaxonSet();
    taxonSet.taxonsetInput.setValue(new Taxon("t1"), taxonSet);
    taxonSet.taxonsetInput.setValue(new Taxon("t2"), taxonSet);
    Locus locus = new Locus("locus", 10000);
    try (PrintStream ps = new PrintStream("out.txt")) {
        for (int i = 0; i < 100000; i++) {
            acg.initByName("locus", locus, "taxonset", taxonSet, "fromString", "(0:1.0,1:1.0)2:0.0;");
            operator.drawNewConversion();
            ps.println(acg.getConversions(locus).get(0).getStartSite() + " " + acg.getConversions(locus).get(0).getEndSite());
        }
    }
}
Also used : PrintStream(java.io.PrintStream) ConstantPopulation(beast.evolution.tree.coalescent.ConstantPopulation) Taxon(beast.evolution.alignment.Taxon) RealParameter(beast.core.parameter.RealParameter) TaxonSet(beast.evolution.alignment.TaxonSet) Locus(bacter.Locus) ConversionGraph(bacter.ConversionGraph)

Aggregations

RealParameter (beast.core.parameter.RealParameter)97 Test (org.junit.Test)39 IntegerParameter (beast.core.parameter.IntegerParameter)23 Tree (beast.evolution.tree.Tree)16 State (beast.core.State)14 SCMigrationModel (beast.evolution.tree.SCMigrationModel)13 Alignment (beast.evolution.alignment.Alignment)11 TypeSet (beast.evolution.tree.TypeSet)11 MCMC (beast.core.MCMC)10 SiteModel (beast.evolution.sitemodel.SiteModel)10 Frequencies (beast.evolution.substitutionmodel.Frequencies)10 ConstantPopulation (beast.evolution.tree.coalescent.ConstantPopulation)10 StructuredCoalescentTreeDensity (multitypetree.distributions.StructuredCoalescentTreeDensity)10 TaxonSet (beast.evolution.alignment.TaxonSet)9 MultiTypeTreeStatLogger (multitypetree.util.MultiTypeTreeStatLogger)9 MultiTypeTreeFromNewick (beast.evolution.tree.MultiTypeTreeFromNewick)8 Node (beast.evolution.tree.Node)8 Operator (beast.core.Operator)7 RandomTree (beast.evolution.tree.RandomTree)7 Locus (bacter.Locus)6