Package plm.logit.fruehwirth.multi

Examples of plm.logit.fruehwirth.multi.FruehwirthMultiParticle


          public int compare(Double o1, Double o2) {
            return o1 < o2 ? 1 : -1;
          }
        });
    for (Entry<FruehwirthMultiParticle, ? extends Number> particleEntry : target.asMap().entrySet()) {
      final FruehwirthMultiParticle particle = particleEntry.getKey();

      /*
       * Fruewirth-Schnatter's method for sampling...
       * TODO do we sample new y's for each particle?
       * TODO when using FS aug sampling, we should probably go all the way
       * and replicate their beta sampling.
       * That would require a change here...
       */
      final double U = this.random.nextDouble();
      double lambdaSum = 0d;
      Vector sampledAugResponse = VectorFactory.getDefault().createVector(this.numCategories);
      for (int i = 0; i < this.numCategories; i++) {
  //        final Vector betaSample = particle.getLinearState().sample(this.random);
        final MultivariateGaussian predictivePrior = particle.getLinearState().clone();
        KalmanFilter kf = particle.getRegressionFilter(i);
        final Matrix G = kf.getModel().getA();
        predictivePrior.setMean(G.times(predictivePrior.getMean()));
        predictivePrior.setCovariance(
            G.times(predictivePrior.getCovariance()).times(G.transpose())
              .plus(kf.getModelCovariance()));
 
        // X * beta
        final double lambda = Math.exp(data.getObservedData().times(
            predictivePrior.getMean()).getElement(0));
        lambdaSum += lambda;
        final double dSampledAugResponse = -Math.log(
            -Math.log(U)/(1d+lambdaSum)
            - (data.getObservedValue().getElement(0) > 0d
                ? 0d : Math.log(this.random.nextDouble())/lambda));
        sampledAugResponse.setElement(i, dSampledAugResponse);
      }

      /*
       * Expand particle set over mixture components
       */
      for (int j = 0; j < 10; j++) {
        double categoriesTotalLogLikelihood = 0d;

        final UnivariateGaussian componentDist =
            this.evDistribution.getDistributions().get(j);
 
        for (int k = 0; k < this.numCategories; k++) {
          /*
           * TODO could avoid cloning if we didn't change the measurement covariance,
           * but instead used the componentDist explicitly.
           */
          final FruehwirthMultiParticle predictiveParticle = particle.clone();
          predictiveParticle.setAugResponseSample(sampledAugResponse);
 
          predictiveParticle.setEVcomponent(componentDist);
         
          /*
           * Update the observed data for the regression component.
           */
          predictiveParticle.getRegressionFilter(k).getModel().setC(data.getObservedData());
 
          final Matrix compVar = MatrixFactory.getDefault().copyArray(
              new double[][] {{componentDist.getVariance()}});
          predictiveParticle.getRegressionFilter(k).setMeasurementCovariance(compVar);
         
          final double logLikelihood = this.updater.computeLogLikelihood(predictiveParticle, data)
              + Math.log(this.evDistribution.getPriorWeights()[j])
              + (particleEntry.getValue().doubleValue() - prevTotalLogLikelihood);
 
          categoriesTotalLogLikelihood = LogMath.subtract(categoriesTotalLogLikelihood, logLikelihood);
          particleTotalLogLikelihood = LogMath.add(particleTotalLogLikelihood, logLikelihood);

          particleTree.put(logLikelihood, predictiveParticle);
        }
      }
    }

//    final WFCountedDataDistribution<FruehwirthMultiParticle> resampledParticles =
//        ExtSamplingUtils.waterFillingResample(logLikelihoods, particleTotalLogLikelihood,
//            particleSupport, random, this.numParticles);

    final WFCountedDataDistribution<FruehwirthMultiParticle> resampledParticles =
        ExtSamplingUtils.waterFillingResample(
            Doubles.toArray(particleTree.keySet()),
            particleTotalLogLikelihood,
            Lists.newArrayList(particleTree.values()),
            random, this.numParticles);

    /*
     * Propagate
     */
    target.clear();
    for (final Entry<FruehwirthMultiParticle, MutableDouble> particleEntry : resampledParticles.asMap().entrySet()) {
      final FruehwirthMultiParticle updatedParticle = sufficientStatUpdate(particleEntry.getKey(), data);
      target.set(updatedParticle, particleEntry.getValue().value);
    }
    Preconditions.checkState(target.getDomainSize() == this.numParticles);

  }
View Full Code Here


  }

  private FruehwirthMultiParticle sufficientStatUpdate(
      FruehwirthMultiParticle priorParticle, ObservedValue<Vector, Matrix> data) {
    final FruehwirthMultiParticle updatedParticle = priorParticle.clone();
   
    final Vector sampledAugResponse = priorParticle.getAugResponseSample();
    final KalmanFilter filter = updatedParticle.getRegressionFilter(
        priorParticle.getCategoryId());
    final UnivariateGaussian evComponent = updatedParticle.EVcomponent;
    // TODO we should've already set this, so it might be redundant.
    filter.setMeasurementCovariance(
        MatrixFactory.getDefault().copyArray(new double[][] {{
          evComponent.getVariance()}}));

    final MultivariateGaussian posteriorState = updatedParticle.getLinearState();
    filter.update(posteriorState, sampledAugResponse);
   
    return updatedParticle;
  }
View Full Code Here

       
        /*
         * Without an observation, start with any category...
         */
        final int categoryId = this.rng.nextInt(this.numCategories);
        final FruehwirthMultiParticle particleMvgDPDist =
            new FruehwirthMultiParticle(null,
                kf, initialPriorState, null,
                categoryId);

        initialParticles.increment(particleMvgDPDist);
      }
View Full Code Here

TOP

Related Classes of plm.logit.fruehwirth.multi.FruehwirthMultiParticle

Copyright © 2018 www.massapicom. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.