Java DistributionLikelihood类代码示例

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public static void main(String[] arg) throws IOException, TraceException {

        // constructing random variable representing mean of normal distribution
        Variable.D mean = new Variable.D("mean", 1.0);
        // give mean a uniform prior [-1000, 1000]
        mean.addBounds(new Parameter.DefaultBounds(1000, -1000, 1));

        // constructing random variable representing stdev of normal distribution
        Variable.D stdev = new Variable.D("stdev", 1.0);
        // give stdev a uniform prior [0, 1000]
        stdev.addBounds(new Parameter.DefaultBounds(1000, 0, 1));

        // construct normal distribution model
        NormalDistributionModel normal = new NormalDistributionModel(mean, stdev);

        // construct a likelihood for normal distribution
        DistributionLikelihood likelihood = new DistributionLikelihood(normal);

        // construct data
        Attribute.Default<double[]> d = new Attribute.Default<double[]>(
                "x", new double[]{1, 2, 3, 4, 5, 6, 7, 8, 9});

        // add data (representing 9 independent observations) to likelihood
        likelihood.addData(d);

        // construct two "operators" to be used as the proposal distribution
        MCMCOperator meanMove = new ScaleOperator(mean, 0.75);
        MCMCOperator stdevMove = new ScaleOperator(stdev, 0.75);

        // construct a logger to log progress of MCMC run to stdout (screen)
        MCLogger logger1 = new MCLogger(100);
        logger1.add(mean);
        logger1.add(stdev);

        // construct a logger to log to a log file for later analysis
        MCLogger logger2 = new MCLogger("tutorial1.log", 100, false, 0);
        logger2.add(mean);
        logger2.add(stdev);

        // construct MCMC object
        MCMC mcmc = new MCMC("tutorial1:normal");

        // initialize MCMC with chain length, likelihood, operators and loggers
        mcmc.init(100000, likelihood, new MCMCOperator[]{meanMove, stdevMove}, new Logger[]{logger1, logger2});

        // run the mcmc
        mcmc.chain();

        // perform post-analysis
        TraceAnalysis.report("tutorial1.log");
    }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public static void main(String[] arg) {

        // Define normal model
        Parameter meanParameter = new Parameter.Default(1.0); // Starting value
        Variable<Double> stdev = new Variable.D(1.0, 1); // Fixed value
        ParametricDistributionModel densityModel = new NormalDistributionModel(meanParameter, stdev);
        DistributionLikelihood likelihood = new DistributionLikelihood(densityModel);

        // Define prior
        DistributionLikelihood prior = new DistributionLikelihood(new NormalDistribution(0.0, 1.0)); // Hyper-priors
        prior.addData(meanParameter);

        // Define data
        likelihood.addData(new Attribute.Default<double[]>("Data", new double[] {0.0, 1.0, 2.0}));

        List<Likelihood> list = new ArrayList<Likelihood>();
        list.add(likelihood);
        list.add(prior);
        CompoundLikelihood posterior = new CompoundLikelihood(0, list);
        SliceOperator sliceSampler = new SliceOperator(meanParameter);

        final int length = 10000;
        double mean = 0;
        double variance = 0;

        for(int i = 0; i < length; i++) {
            sliceSampler.doOperation(posterior);
            double x = meanParameter.getValue(0);
            mean += x;
            variance += x*x;
        }
        mean /= length;
        variance /= length;
        variance -= mean*mean;
        System.out.println("E(x) = "+mean);
        System.out.println("V(x) = "+variance);
    }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {

            double mean = xo.getDoubleAttribute(MEAN);
            double offset = xo.getDoubleAttribute(OFFSET);

            DistributionLikelihood likelihood = new DistributionLikelihood(new PoissonDistribution(mean), offset);
            for (int j = 0; j < xo.getChildCount(); j++) {
                if (xo.getChild(j) instanceof Statistic) {
                    likelihood.addData((Statistic) xo.getChild(j));
                } else {
                    throw new XMLParseException("illegal element in " + xo.getName() + " element");
                }
            }

            return likelihood;
        }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
private List<Callable<DrawResult>> createTasks() {

        List<Callable<DrawResult>> callers = new ArrayList<Callable<DrawResult>>();

        int offset = 0;
        int index = 0;
        for (GaussianProcessRandomGenerator gp : gpList) {
            final int copies = copyList.get(index);
            if (likelihoodList.get(index) instanceof DistributionLikelihood) { // Univariate
                callers.add(new DrawCaller(gp, copies, offset, true));
                offset += copies;
            } else {
                for (int i = 0; i < copies; ++i) {
                    callers.add(new DrawCaller(gp, 1, offset, false));
                    offset += gp.getDimension();
                }
            }
        }

        return callers;
    }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public LoadingsGibbsTruncatedOperator(LatentFactorModel LFM, Likelihood prior, double weight, boolean randomScan, MatrixParameterInterface loadings, DistributionLikelihood cutoffPrior) {
    setWeight(weight);

    this.loadings=loadings;
    this.prior = prior;
    this.LFM = LFM;
    if(prior instanceof MomentDistributionModel){
    priorPrecision = (((MomentDistributionModel) this.prior).getScaleMatrix()[0][0]);
    priorMeanPrecision = ((MomentDistributionModel) this.prior).getMean()[0] * priorPrecision;}
    else if (prior instanceof DistributionLikelihood){
        priorPrecision = 1 / ((DistributionLikelihood) this.prior).getDistribution().variance();
        priorMeanPrecision = ((DistributionLikelihood) this.prior).getDistribution().mean() * priorPrecision;
    }
    this.cutoffPrior = cutoffPrior;
    missingIndicator = LFM.getMissingIndicator();
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {

            double weight = xo.getDoubleAttribute(WEIGHT);

            DistributionLikelihood likelihood = (DistributionLikelihood) ((XMLObject) xo.getChild(LIKELIHOOD)).getChild(DistributionLikelihood.class);
            DistributionLikelihood prior = (DistributionLikelihood) ((XMLObject) xo.getChild(PRIOR)).getChild(DistributionLikelihood.class);

//            System.err.println("class: " + prior.getDistribution().getClass());

            if (!((prior.getDistribution() instanceof GammaDistribution) ||
                    (prior.getDistribution() instanceof GammaDistributionModel)
            ) ||
                    !((likelihood.getDistribution() instanceof NormalDistributionModel) ||
                            (likelihood.getDistribution() instanceof LogNormalDistributionModel)
                    ))
                throw new XMLParseException("Gibbs operator assumes normal-gamma model");

            return new NormalGammaPrecisionGibbsOperator(likelihood, prior.getDistribution(), weight);
        }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public NormalNormalMeanGibbsOperator(DistributionLikelihood inLikelihood, Distribution prior,
                                     double weight) {

    if (!(prior instanceof NormalDistribution || prior instanceof NormalDistributionModel))
        throw new RuntimeException("Mean prior must be Normal");

    this.likelihood = inLikelihood.getDistribution();
    this.dataList = inLikelihood.getDataList();
    if (likelihood instanceof NormalDistributionModel)
        this.meanParameter = (Parameter) ((NormalDistributionModel) likelihood).getMean();
    else if (likelihood instanceof LogNormalDistributionModel) {
        if (((LogNormalDistributionModel) likelihood).getParameterization() == LogNormalDistributionModel.Parameterization.MEAN_STDEV) {
            this.meanParameter = ((LogNormalDistributionModel) likelihood).getMeanParameter();
        } else {
            this.meanParameter = ((LogNormalDistributionModel) likelihood).getMuParameter();
        }
        isLog = true;
    } else
        throw new RuntimeException("Likelihood must be Normal or log Normal");

    this.prior = prior;
    setWeight(weight);
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {

            double weight = xo.getDoubleAttribute(WEIGHT);

            DistributionLikelihood likelihood = (DistributionLikelihood) ((XMLObject) xo.getChild(LIKELIHOOD)).getChild(DistributionLikelihood.class);
            DistributionLikelihood prior = (DistributionLikelihood) ((XMLObject) xo.getChild(PRIOR)).getChild(DistributionLikelihood.class);

//            System.err.println("class: " + prior.getDistribution().getClass());

            if (!((prior.getDistribution() instanceof NormalDistribution) ||
                    (prior.getDistribution() instanceof NormalDistributionModel)
            ) ||
                    !((likelihood.getDistribution() instanceof NormalDistributionModel) ||
                            (likelihood.getDistribution() instanceof LogNormalDistributionModel)
                    ))
                throw new XMLParseException("Gibbs operator assumes normal-normal model");

            return new NormalNormalMeanGibbsOperator(likelihood, prior.getDistribution(), weight);
        }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    DistributionLikelihood priorLikelihood = (DistributionLikelihood)xo.getElementFirstChild(PRIOR);
    DistributionLikelihood pseudoPriorLikelihood = (DistributionLikelihood)xo.getElementFirstChild(PSEUDO_PRIOR);
    Distribution prior = priorLikelihood.getDistribution();
    Distribution pseudoPrior = pseudoPriorLikelihood.getDistribution();
    Parameter modelIndicator = (Parameter)xo.getElementFirstChild(MODEL_INDICATOR);
    int[] models = xo.getIntegerArrayAttribute(MODELS);
    Parameter selectedVariable = (Parameter)xo.getElementFirstChild(SELECTED_VARIABLE);


    ModelSpecificPseudoPriorLikelihood likelihood =
            new ModelSpecificPseudoPriorLikelihood(
                    prior,
                    pseudoPrior,
                    modelIndicator,
                    models
            );
    likelihood.addData(selectedVariable);
    

    return likelihood;
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {
    DistributionLikelihood priorLikelihood = (DistributionLikelihood)xo.getElementFirstChild(PRIOR);
    DistributionLikelihood pseudoPriorLikelihood = (DistributionLikelihood)xo.getElementFirstChild(PSEUDO_PRIOR);
    Distribution prior = priorLikelihood.getDistribution();
    Distribution pseudoPrior = pseudoPriorLikelihood.getDistribution();
    Parameter bitVector = (Parameter)xo.getElementFirstChild(PARAMETER_VECTOR);
    int paramIndex = xo.getIntegerAttribute(PARAMETER_INDEX);
    Parameter selectedVariable = (Parameter)xo.getElementFirstChild(SELECTED_VARIABLE);


    TwoPartsDistributionLikelihood likelihood =
            new TwoPartsDistributionLikelihood(
                    prior,
                    pseudoPrior,
                    bitVector,
                    paramIndex
            );
    likelihood.addData(selectedVariable);

    return likelihood;
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
@Override
public Object parseXMLObject(XMLObject xo) throws XMLParseException {

    Object obj = xo.getChild(0);

    if (obj instanceof MultivariateDistributionLikelihood) {
        final MultivariateDistributionLikelihood mdl = (MultivariateDistributionLikelihood)
                xo.getChild(MultivariateDistributionLikelihood.class);

        if (!(mdl.getDistribution() instanceof GradientProvider)) {
            throw new XMLParseException("Not a gradient provider");
        }

        final GradientProvider provider = (GradientProvider) mdl.getDistribution();
        final Parameter parameter = mdl.getDataParameter();

        return new GradientWrtParameterProvider.ParameterWrapper(provider, parameter, mdl);
    } else {
        DistributionLikelihood dl = (DistributionLikelihood) obj;
        if (!(dl.getDistribution() instanceof GradientProvider)) {
            throw new XMLParseException("Not a gradient provider");
        }

        throw new RuntimeException("Not yet implemented");
    }
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public Object parseXMLObject(XMLObject xo) throws XMLParseException {

            double mean = xo.getDoubleAttribute(MEAN);
            double stdev = xo.getDoubleAttribute(STDEV);

            DistributionLikelihood likelihood = new DistributionLikelihood(new NormalDistribution(mean, stdev));
            for (int j = 0; j < xo.getChildCount(); j++) {
                if (xo.getChild(j) instanceof Statistic) {
                    likelihood.addData((Statistic) xo.getChild(j));
                } else if (!(xo.getChild(j) instanceof LatentTruncation)) {
                    throw new XMLParseException("Illegal element in " + xo.getName() + " element");
                }
            }

            LatentTruncation truncation = (LatentTruncation) xo.getChild(LatentTruncation.class);

            return new TruncatedWorkingDistribution(likelihood, truncation);
        }
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public NormalGammaPrecisionGibbsOperator(DistributionLikelihood inLikelihood, Distribution prior,
                                         double weight) {

    if (!(prior instanceof GammaDistribution || prior instanceof GammaDistributionModel))
        throw new RuntimeException("Precision prior must be Gamma");

    Distribution likelihood = inLikelihood.getDistribution();
    this.dataList = inLikelihood.getDataList();
    if (likelihood instanceof NormalDistributionModel) {
        this.precisionParameter = (Parameter) ((NormalDistributionModel) likelihood).getPrecision();
        this.meanParameter = (Parameter) ((NormalDistributionModel) likelihood).getMean();
    } else if (likelihood instanceof LogNormalDistributionModel) {
        this.precisionParameter = ((LogNormalDistributionModel) likelihood).getPrecisionParameter();
        this.meanParameter = ((LogNormalDistributionModel) likelihood).getMeanParameter();
        isLog = true;
    } else
        throw new RuntimeException("Likelihood must be Normal or log Normal");

    if (precisionParameter == null)
        throw new RuntimeException("Must characterize likelihood in terms of a precision parameter");

    this.prior = prior;
    setWeight(weight);
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public NormalNormalMeanGibbsOperator(DistributionLikelihood inLikelihood, Distribution prior,
                                     double weight) {

    if (!(prior instanceof NormalDistribution || prior instanceof NormalDistributionModel))
        throw new RuntimeException("Mean prior must be Normal");

    this.likelihood = inLikelihood.getDistribution();
    this.dataList = inLikelihood.getDataList();
    if (likelihood instanceof NormalDistributionModel)
        this.meanParameter = (Parameter) ((NormalDistributionModel) likelihood).getMean();
    else if (likelihood instanceof LogNormalDistributionModel) {
        this.meanParameter = ((LogNormalDistributionModel) likelihood).getMeanParameter();
        isLog = true;
    } else
        throw new RuntimeException("Likelihood must be Normal or log Normal");

    this.prior = prior;
    setWeight(weight);
}
 

import dr.inference.distribution.DistributionLikelihood; //导入依赖的package包/类
public void testLognormalPrior() {
//        ConstantPopulation constant = new ConstantPopulation(Units.Type.YEARS);
//        constant.setN0(popSize); // popSize
        Parameter popSize = new Parameter.Default(6.0);
        popSize.setId(ConstantPopulationModelParser.POPULATION_SIZE);
        ConstantPopulationModel demo = new ConstantPopulationModel(popSize, Units.Type.YEARS);

        //Likelihood
        Likelihood dummyLikelihood = new DummyLikelihood(demo);

        // Operators
        OperatorSchedule schedule = new SimpleOperatorSchedule();

        MCMCOperator operator = new ScaleOperator(popSize, 0.75);
        operator.setWeight(1.0);
        schedule.addOperator(operator);

        // Log
        ArrayLogFormatter formatter = new ArrayLogFormatter(false);

        MCLogger[] loggers = new MCLogger[2];
        loggers[0] = new MCLogger(formatter, 1000, false);
//        loggers[0].add(treeLikelihood);
        loggers[0].add(popSize);

        loggers[1] = new MCLogger(new TabDelimitedFormatter(System.out), 100000, false);
//        loggers[1].add(treeLikelihood);
        loggers[1].add(popSize);

        // MCMC
        MCMC mcmc = new MCMC("mcmc1");
        MCMCOptions options = new MCMCOptions();
        options.setChainLength(1000000);
        options.setUseCoercion(true); // autoOptimize = true
        options.setCoercionDelay(100);
        options.setTemperature(1.0);
        options.setFullEvaluationCount(2000);

        DistributionLikelihood logNormalLikelihood = new DistributionLikelihood(new LogNormalDistribution(1.0, 1.0), 0); // meanInRealSpace="false"
        logNormalLikelihood.addData(popSize);

        List<Likelihood> likelihoods = new ArrayList<Likelihood>();
        likelihoods.add(logNormalLikelihood);
        Likelihood prior = new CompoundLikelihood(0, likelihoods);

        likelihoods.clear();
        likelihoods.add(dummyLikelihood);
        Likelihood likelihood = new CompoundLikelihood(-1, likelihoods);

        likelihoods.clear();
        likelihoods.add(prior);
        likelihoods.add(likelihood);
        Likelihood posterior = new CompoundLikelihood(0, likelihoods);

        mcmc.setShowOperatorAnalysis(true);
        mcmc.init(options, posterior, schedule, loggers);
        mcmc.run();

        // time
        System.out.println(mcmc.getTimer().toString());

        // Tracer
        List<Trace> traces = formatter.getTraces();
        ArrayTraceList traceList = new ArrayTraceList("LognormalPriorTest", traces, 0);

        for (int i = 1; i < traces.size(); i++) {
            traceList.analyseTrace(i);
        }

//      <expectation name="param" value="4.48168907"/>

        TraceCorrelation popSizeStats = traceList.getCorrelationStatistics(traceList.getTraceIndex(ConstantPopulationModelParser.POPULATION_SIZE));

        System.out.println("Expectation of Log-Normal(1,1) is e^(M+S^2/2) = e^(1.5) = " + Math.exp(1.5));
        assertExpectation(ConstantPopulationModelParser.POPULATION_SIZE, popSizeStats, Math.exp(1.5));
    }