locked
Problems modeling multivariate Gaussian distribution RRS feed

  • Question

  • Hello,

    I'm trying to extend Bayesian PCA example to Bayesian mixture of FA, with full covariance matrix on the rows of matrix of factor loadings (vW in this case). But when I'm trying to define vW this way:

    vW[rC][rM] = Variable.Random(new VectorGaussian(Vector.Zero(vD),
                                                 PositiveDefiniteMatrix.IdentityScaledBy(vD, vAlpha[rC][rM])));

    compiler throws a bunch of errors like this: "cannot convert from 'MicrosoftResearch.Infer.Models.Variable<int>' to 'int' ". I understand what it means, but I can't figure out how to solve this. Could you please hint on how to do that?

    Here is the full code:

    public class BayesianMFA
        {
            // Inference engine
            public InferenceEngine engine = null;
            // Model variables
            public Variable<int> vN = null; // Number of data points
            public Variable<int> vD = null; // Dimensionality of data
            public Variable<int> vM = null; // Dimensionality of latent variables
            public Variable<int> vC = null; // Number of components
            public VariableArray2D<double> vData = null; // Data
    
            public VariableArray<VariableArray<Vector>, Vector[][]> vW = null;
            public VariableArray<VariableArray2D<double>, double[][,]> vZ = null;
            public VariableArray<VariableArray2D<double>, double[][,]> vT = null;
            public VariableArray<VariableArray2D<double>, double[][,]> vU = null;
    
            public VariableArray<VariableArray<double>, double[][]> vMu = null;
            public VariableArray<VariableArray<double>, double[][]> vPi = null;
            public VariableArray<VariableArray<double>, double[][]> vAlpha = null;
    
            // Priors - these are declared as distribution variables
            // so that we can set them at run-time.
            public VariableArray<Gamma> priorAlpha = null;
            public VariableArray<Gaussian> priorMu = null;
            public VariableArray<Gamma> priorPi = null;
            // Index of mixture component for each data point
            // and weights of the mixture
            public VariableArray<int> index = null;
            public Variable<Vector> weights = null;
            // Ranges of vN, vD, vM, vC
            public Range rN = null; 
            public Range rD = null; 
            public Range rM = null; 
            public Range rC = null; 
            /// <summary>
            /// Model constructor
            /// </summary>
            public BayesianMFA()
            {
                // The various dimensions will be set externally...
                vN = Variable.New<int>().Named("NumObs");
                vD = Variable.New<int>().Named("NumFeats");
                vM = Variable.New<int>().Named("MaxComponents");
                vC = Variable.New<int>().Named("NumOfMixComponents");
    
                rN = new Range(vN).Named("N");
                rD = new Range(vD).Named("D");
                rM = new Range(vM).Named("M");
                rC = new Range(vC).Named("C");
                
                // ... as will the data
                vData = Variable.Array<double>(rN, rD).Named("data");
                
                // ... and the priors
                priorAlpha = Variable.Array<Gamma>(rC).Named("PriorAlpha");
                priorMu = Variable.Array<Gaussian>(rC).Named("PriorMu");
                priorPi = Variable.Array<Gamma>(rC).Named("PriorPi");
    
                // Mixing matrix
                vAlpha = Variable.Array(Variable.Array<double>(rM), rC).Named("Alpha");
                vW = Variable.Array(Variable.Array<Vector>(rM), rC).Named("W");
                vAlpha[rC][rM] = Variable.Random<double, Gamma>(priorAlpha[rC]).ForEach(rM);
                vW[rC][rM] = Variable.Random(new VectorGaussian(Vector.Zero(vD),
                                                 PositiveDefiniteMatrix.IdentityScaledBy(vD, vAlpha[rC][rM])));
                
                // Latent variables are drawn from a standard Gaussian
                vZ = Variable.Array(Variable.Array<double>(rN,rM), rC).Named("Z");
                vZ[rC][rN,rM] = Variable.GaussianFromMeanAndPrecision(0.0, 1.0).ForEach(rN, rM, rC);
                
                // Multiply the latent variables with the mixing matrix...
                vT = Variable.Array(Variable.Array<double>(rN, rD), rC).Named("T");
                VariableArray<VariableArray<VariableArray<double>, double[][]>, double[][][]> vWMatrix;
                vWMatrix = Variable.Array(Variable.Array(Variable.Array<double>(rD),rM),rC).Named("U");
                vWMatrix[rC][rM] = Variable.ArrayFromVector(vW[rC][rM], rD);
                // Explicit Matrix multiply
                using (Variable.ForEach(rC))
                {
                    using (Variable.ForEach(rN))
                    {
                        using (Variable.ForEach(rD))
                        {
                            var vProds = Variable.Array<double>(rM);
                            vProds[rM] = vZ[rC][rN,rM] * vWMatrix[rC][rM][rD];
                            vT[rC][rN,rD] = Variable.Sum(vProds);
                        }
                    }
                }
                
                // ... add in a bias ...
                vMu = Variable.Array(Variable.Array<double>(rD), rC).Named("mu");
                vMu[rC][rD] = Variable.Random<double, Gaussian>(priorMu[rC]).ForEach(rD);
                vU = Variable.Array(Variable.Array<double>(rN, rD), rC).Named("U");
                vU[rC][rN, rD] = vT[rC][rN, rD] + vMu[rC][rD];
                
                // ... and add in some observation noise ...
                vPi = Variable.Array(Variable.Array<double>(rD), rC).Named("pi");
                vPi[rC][rD] = Variable.Random<double, Gamma>(priorPi[rC]).ForEach(rD);
                
                // Working with the mixture
                weights = Variable.New<Vector>().Named("weights");
                weights = Variable.DirichletSymmetric(rC, 1);
                index = Variable.Array<int>(rN).Named("index");
                using (Variable.ForEach(rN))
                {
                    index[rN] = Variable.Discrete(weights);
                    using (Variable.Switch(index[rN]))
                    {
                        vData[rN, rD] = Variable.GaussianFromMeanAndPrecision(vU[index[rN]][rN, rD], vPi[index[rN]][rD]);
                    }
                }
                // Inference engine
                engine = new InferenceEngine(new VariationalMessagePassing());
                return;
            }
        }

    Thanks a lot.

    Friday, April 27, 2012 1:16 PM

All replies

  • The issue is how to create a VectorGaussian variable with a diagonal precision matrix.  This is possible by first creating a uniform VectorGaussian and then applying constraints, as illustrated in the code below:
    public void VectorGaussianFromPrecisionDiagonalExample()
    {
    	Range d = new Range(3);
    	VariableArray<double> prec = Variable.Array<double>(d);
    	prec[d] = Variable.GammaFromShapeAndRate(1,1).ForEach(d);
    	Variable<Vector> x = VectorGaussianFromPrecisionDiagonal(prec);
    	InferenceEngine engine = new InferenceEngine(new VariationalMessagePassing());
    	Console.WriteLine(engine.Infer(x));
    }
    public Variable<Vector> VectorGaussianFromPrecisionDiagonal(VariableArray<double> precisionDiagonal)
    {
    	Range r = precisionDiagonal.Range;
    	int dimension = r.SizeAsInt;
    	Variable<Vector> x = Variable.VectorGaussianFromMeanAndPrecision(Vector.Zero(dimension), PositiveDefiniteMatrix.IdentityScaledBy(dimension, 1e-20));
    	using (var fb = Variable.ForEach(r)) {
    		Variable<double> xi = Variable.GetItem(x, fb.Index);
    		var y = Variable.GaussianFromMeanAndPrecision(xi, precisionDiagonal[r]);
    		y.ObservedValue = 0;
    	}
    	return x;
    }
    
    This technique allows essentially any structure to be imposed on the precision matrix.
    Friday, April 27, 2012 6:50 PM
    Owner