Search in sources :

Example 1 with AbstractPeriodPriorDistribution

use of dr.evomodel.epidemiology.casetocase.periodpriors.AbstractPeriodPriorDistribution in project beast-mcmc by beast-dev.

the class CaseToCaseTransmissionLikelihood method getLogLikelihood.

public double getLogLikelihood() {
    if (!likelihoodKnown) {
        if (!treeProbKnown) {
            treeLikelihood.prepareTimings();
        }
        if (!transProbKnown) {
            try {
                transLogProb = 0;
                if (sortedTreeEvents == null) {
                    sortEvents();
                }
                double rate = transmissionRate.getParameterValue(0);
                ArrayList<AbstractCase> previouslyInfectious = new ArrayList<AbstractCase>();
                double currentEventTime;
                boolean first = true;
                for (TreeEvent event : sortedTreeEvents) {
                    currentEventTime = event.getTime();
                    AbstractCase thisCase = event.getCase();
                    if (event.getType() == EventType.INFECTION) {
                        if (first) {
                            if (indexCasePrior != null) {
                                transLogProb += Math.log(indexCasePrior.get(thisCase));
                            }
                            if (initialInfectionTimePrior != null) {
                                transLogProb += initialInfectionTimePrior.logPdf(currentEventTime);
                            }
                            if (!hasLatentPeriods) {
                                previouslyInfectious.add(thisCase);
                            }
                            first = false;
                        } else {
                            AbstractCase infector = event.getInfector();
                            if (thisCase.wasEverInfected()) {
                                if (previouslyInfectious.contains(thisCase)) {
                                    throw new BadPartitionException(thisCase.caseID + " infected after it was infectious");
                                }
                                if (event.getTime() > thisCase.endOfInfectiousTime) {
                                    throw new BadPartitionException(thisCase.caseID + " ceased to be infected before it was infected");
                                }
                                if (infector.endOfInfectiousTime < event.getTime()) {
                                    throw new BadPartitionException(thisCase.caseID + " infected by " + infector.caseID + " after the latter ceased to be infectious");
                                }
                                if (treeLikelihood.getInfectiousTime(infector) > event.getTime()) {
                                    throw new BadPartitionException(thisCase.caseID + " infected by " + infector.caseID + " before the latter became infectious");
                                }
                                if (!previouslyInfectious.contains(infector)) {
                                    throw new RuntimeException("Infector not previously infected");
                                }
                            }
                            for (AbstractCase nonInfector : previouslyInfectious) {
                                double timeDuringWhichNoInfection;
                                if (nonInfector.endOfInfectiousTime < event.getTime()) {
                                    timeDuringWhichNoInfection = nonInfector.endOfInfectiousTime - treeLikelihood.getInfectiousTime(nonInfector);
                                } else {
                                    timeDuringWhichNoInfection = event.getTime() - treeLikelihood.getInfectiousTime(nonInfector);
                                }
                                if (timeDuringWhichNoInfection < 0) {
                                    throw new RuntimeException("negative time");
                                }
                                double transRate = rate;
                                if (hasGeography) {
                                    transRate *= outbreak.getKernelValue(thisCase, nonInfector, spatialKernel);
                                }
                                transLogProb += -transRate * timeDuringWhichNoInfection;
                            }
                            if (thisCase.wasEverInfected()) {
                                double transRate = rate;
                                if (hasGeography) {
                                    transRate *= outbreak.getKernelValue(thisCase, infector, spatialKernel);
                                }
                                transLogProb += Math.log(transRate);
                            }
                            if (!hasLatentPeriods) {
                                previouslyInfectious.add(thisCase);
                            }
                        }
                    } else if (event.getType() == EventType.INFECTIOUSNESS) {
                        if (event.getTime() < Double.POSITIVE_INFINITY) {
                            if (event.getTime() > event.getCase().endOfInfectiousTime) {
                                throw new BadPartitionException(event.getCase().caseID + " noninfectious before" + "infectious");
                            }
                            if (first) {
                                throw new RuntimeException("First event is not an infection");
                            }
                            previouslyInfectious.add(thisCase);
                        }
                    }
                }
                transProbKnown = true;
            } catch (BadPartitionException e) {
                transLogProb = Double.NEGATIVE_INFINITY;
                transProbKnown = true;
                logLikelihood = Double.NEGATIVE_INFINITY;
                likelihoodKnown = true;
                return logLikelihood;
            }
        }
        if (!periodsProbKnown) {
            periodsLogProb = 0;
            HashMap<String, ArrayList<Double>> infectiousPeriodsByCategory = new HashMap<String, ArrayList<Double>>();
            for (AbstractCase aCase : outbreak.getCases()) {
                if (aCase.wasEverInfected()) {
                    String category = (outbreak).getInfectiousCategory(aCase);
                    if (!infectiousPeriodsByCategory.keySet().contains(category)) {
                        infectiousPeriodsByCategory.put(category, new ArrayList<Double>());
                    }
                    ArrayList<Double> correspondingList = infectiousPeriodsByCategory.get(category);
                    correspondingList.add(treeLikelihood.getInfectiousPeriod(aCase));
                }
            }
            for (String category : outbreak.getInfectiousCategories()) {
                Double[] infPeriodsInThisCategory = infectiousPeriodsByCategory.get(category).toArray(new Double[infectiousPeriodsByCategory.get(category).size()]);
                AbstractPeriodPriorDistribution hyperprior = outbreak.getInfectiousCategoryPrior(category);
                double[] values = new double[infPeriodsInThisCategory.length];
                for (int i = 0; i < infPeriodsInThisCategory.length; i++) {
                    values[i] = infPeriodsInThisCategory[i];
                }
                periodsLogProb += hyperprior.getLogLikelihood(values);
            }
            periodsProbKnown = true;
        }
        if (!treeProbKnown) {
            treeLogProb = treeLikelihood.getLogLikelihood();
            treeProbKnown = true;
        }
        if (transLogProb == Double.POSITIVE_INFINITY) {
            System.out.println("TransLogProb +INF");
            return Double.NEGATIVE_INFINITY;
        }
        if (periodsLogProb == Double.POSITIVE_INFINITY) {
            System.out.println("PeriodsLogProb +INF");
            return Double.NEGATIVE_INFINITY;
        }
        if (treeLogProb == Double.POSITIVE_INFINITY) {
            System.out.println("TreeLogProb +INF");
            return Double.NEGATIVE_INFINITY;
        }
        logLikelihood = treeLogProb + periodsLogProb + transLogProb;
        likelihoodKnown = true;
    }
    return logLikelihood;
}
Also used : AbstractPeriodPriorDistribution(dr.evomodel.epidemiology.casetocase.periodpriors.AbstractPeriodPriorDistribution)

Example 2 with AbstractPeriodPriorDistribution

use of dr.evomodel.epidemiology.casetocase.periodpriors.AbstractPeriodPriorDistribution in project beast-mcmc by beast-dev.

the class CaseToCaseTransmissionLikelihood method getColumns.

// Not the most elegant solution, but you want two types of log out of this model, one for numerical parameters
// (which Tracer can read) and one for the transmission tree (which it cannot). This is set up so that C2CTransL
// is the numerical log and C2CTreeL the TT one.
public LogColumn[] getColumns() {
    ArrayList<LogColumn> columns = new ArrayList<LogColumn>();
    columns.add(new LogColumn.Abstract("trans_LL") {

        protected String getFormattedValue() {
            return String.valueOf(transLogProb);
        }
    });
    columns.add(new LogColumn.Abstract("period_LL") {

        protected String getFormattedValue() {
            return String.valueOf(periodsLogProb);
        }
    });
    columns.addAll(Arrays.asList(treeLikelihood.passColumns()));
    for (AbstractPeriodPriorDistribution hyperprior : (outbreak).getInfectiousMap().values()) {
        columns.addAll(Arrays.asList(hyperprior.getColumns()));
    }
    columns.add(new LogColumn.Abstract("FirstInfectionTime") {

        protected String getFormattedValue() {
            if (sortedTreeEvents == null) {
                sortEvents();
            }
            return String.valueOf(treeLikelihood.getInfectionTime(indexCase));
        }
    });
    columns.add(new LogColumn.Abstract("IndexCaseIndex") {

        protected String getFormattedValue() {
            return String.valueOf(treeLikelihood.getOutbreak().getCaseIndex(indexCase));
        }
    });
    return columns.toArray(new LogColumn[columns.size()]);
}
Also used : LogColumn(dr.inference.loggers.LogColumn) AbstractPeriodPriorDistribution(dr.evomodel.epidemiology.casetocase.periodpriors.AbstractPeriodPriorDistribution)

Aggregations

AbstractPeriodPriorDistribution (dr.evomodel.epidemiology.casetocase.periodpriors.AbstractPeriodPriorDistribution)2 LogColumn (dr.inference.loggers.LogColumn)1