From 1402e624c4228c2f478a6e1a1fb3d365c070d86f Mon Sep 17 00:00:00 2001 From: xji3 Date: Wed, 15 Jun 2022 14:09:56 -0500 Subject: [PATCH 1/2] Bug fix: Use LinkedHashSet in Model sets to reserve insertion order for CheckPoint save and load --- src/dr/app/checkpoint/BeastCheckpointer.java | 9 +++++++-- src/dr/inference/model/Model.java | 4 ++-- 2 files changed, 9 insertions(+), 4 deletions(-) diff --git a/src/dr/app/checkpoint/BeastCheckpointer.java b/src/dr/app/checkpoint/BeastCheckpointer.java index 335a10872f..41c441b58a 100644 --- a/src/dr/app/checkpoint/BeastCheckpointer.java +++ b/src/dr/app/checkpoint/BeastCheckpointer.java @@ -288,7 +288,9 @@ protected boolean writeStateToFile(File file, long state, double lnL, MarkovChai //check up front if there are any TreeParameterModel objects for (Model model : Model.CONNECTED_MODEL_SET) { if (model instanceof TreeParameterModel) { - //System.out.println("\nDetected TreeParameterModel: " + ((TreeParameterModel) model).toString()); + if (DEBUG) { + System.out.println("\nSave TreeParameterModel: " + model.getClass().getSimpleName()); + } traitModels.add((TreeParameterModel) model); } } @@ -505,7 +507,7 @@ protected long readStateFromFile(File file, MarkovChain markovChain, double[] ln // load the tree models last as we get the node heights from the tree (not the parameters which // which may not be associated with the right node - Set expectedTreeModelNames = new HashSet(); + Set expectedTreeModelNames = new LinkedHashSet<>(); //store list of TreeModels for debugging purposes ArrayList treeModelList = new ArrayList(); @@ -529,6 +531,9 @@ protected long readStateFromFile(File file, MarkovChain markovChain, double[] ln //first add all TreeParameterModels to a list if (model instanceof TreeParameterModel) { + if (DEBUG) { + System.out.println("\nLoad TreeParameterModel: " + model.getClass().getSimpleName()); + } traitModels.add((TreeParameterModel)model); } diff --git a/src/dr/inference/model/Model.java b/src/dr/inference/model/Model.java index b67ac53839..96f37d305c 100644 --- a/src/dr/inference/model/Model.java +++ b/src/dr/inference/model/Model.java @@ -170,8 +170,8 @@ public int getListenerCount() { // set to store all created models - final static Set FULL_MODEL_SET = new HashSet(); - final static Set CONNECTED_MODEL_SET = new HashSet(); + final static Set FULL_MODEL_SET = new LinkedHashSet(); + final static Set CONNECTED_MODEL_SET = new LinkedHashSet(); } From 5c7f22ebb88f375cb3c323703f3e4142c1d2d10d Mon Sep 17 00:00:00 2001 From: yucais Date: Wed, 15 Jun 2022 12:58:29 -0700 Subject: [PATCH 2/2] add zero sampling rate model --- ...CriticalBirthDeathSerialSamplingModel.java | 60 ++++++++++++++----- ...ewBirthDeathSerialSamplingModelParser.java | 5 ++ 2 files changed, 51 insertions(+), 14 deletions(-) diff --git a/src/dr/evomodel/speciation/CriticalBirthDeathSerialSamplingModel.java b/src/dr/evomodel/speciation/CriticalBirthDeathSerialSamplingModel.java index b14e883c07..82a2b09acb 100644 --- a/src/dr/evomodel/speciation/CriticalBirthDeathSerialSamplingModel.java +++ b/src/dr/evomodel/speciation/CriticalBirthDeathSerialSamplingModel.java @@ -27,26 +27,58 @@ import dr.inference.model.Parameter; +import static java.lang.Math.exp; +import static java.lang.Math.log; + public class CriticalBirthDeathSerialSamplingModel extends NewBirthDeathSerialSamplingModel { - public CriticalBirthDeathSerialSamplingModel(String modelName, - Parameter birthDeathRate, - Parameter serialSamplingRate, - Parameter treatmentProbability, - Parameter samplingFractionAtPresent, - Parameter originTime, - boolean condition, - Type units) { - - super(modelName, - birthDeathRate, birthDeathRate, serialSamplingRate, - treatmentProbability, samplingFractionAtPresent, originTime, - condition, units); + private int n_events; + + public CriticalBirthDeathSerialSamplingModel( + String modelName, + Parameter birthRate, + Parameter deathRate, + Parameter serialSamplingRate, + Parameter treatmentProbability, + Parameter samplingFractionAtPresent, + Parameter originTime, + boolean condition, + Type units) { + + super(modelName, birthRate, deathRate, serialSamplingRate, treatmentProbability, samplingFractionAtPresent, originTime, condition, units); + n_events = 0; } + @Override public double processInterval(int model, double tYoung, double tOld, int nLineages) { // TODO Do something different return super.processInterval(model, tYoung, tOld, nLineages); } -} \ No newline at end of file + + @Override + public double processOrigin(int model, double rootAge) { + double lambda = lambda(); + double rho = rho(); + double mu = mu(); + double v = exp(-(lambda - mu) * rootAge); + double p_n = log(lambda*rho + (lambda*(1-rho) - mu)* v) - log(1- v); + return -2*logq(rootAge) + (n_events-1)*p_n; + } + + @Override + public double processCoalescence(int model, double tOld) { + n_events += 1; + return 0; + } + + @Override + public double processSampling(int model, double tOld) { + return 0; + } + + @Override + public double logConditioningProbability() { + return -log(n_events); + } +} diff --git a/src/dr/evomodelxml/speciation/NewBirthDeathSerialSamplingModelParser.java b/src/dr/evomodelxml/speciation/NewBirthDeathSerialSamplingModelParser.java index d136f6dcf8..405050977e 100644 --- a/src/dr/evomodelxml/speciation/NewBirthDeathSerialSamplingModelParser.java +++ b/src/dr/evomodelxml/speciation/NewBirthDeathSerialSamplingModelParser.java @@ -26,6 +26,7 @@ package dr.evomodelxml.speciation; import dr.evolution.util.Units; +import dr.evomodel.speciation.CriticalBirthDeathSerialSamplingModel; import dr.evomodel.speciation.NewBirthDeathSerialSamplingModel; import dr.evoxml.util.XMLUnits; import dr.inference.model.Parameter; @@ -76,6 +77,10 @@ public Object parseXMLObject(XMLObject xo) throws XMLParseException { Logger.getLogger("dr.evomodel").info(citeThisModel); + if (psi.getParameterValue(0) < Double.MIN_VALUE){ + return new CriticalBirthDeathSerialSamplingModel(modelName, lambda, mu, psi, r, rho, origin, condition, units); + } + return new NewBirthDeathSerialSamplingModel(modelName, lambda, mu, psi, r, rho, origin, condition, units); }