Home > machines > prt_machine_gpml.m

prt_machine_gpml

PURPOSE ^

Run Gaussian process model - wrapper for gpml toolbox

SYNOPSIS ^

function output = prt_machine_gpml(d,args)

DESCRIPTION ^

 Run Gaussian process model - wrapper for gpml toolbox
 FORMAT output = prt_machine_gpml(d,args)
 Inputs:
   d         - structure with data information, with mandatory fields:
     .train      - training data (cell array of matrices of row vectors,
                   each [Ntr x D]). each matrix contains one representation
                   of the data. This is useful for approaches such as
                   multiple kernel learning.
     .test       - testing data  (cell array of matrices row vectors, each
                   [Nte x D])
     .testcov    - testing covariance (cell array of matrices row vectors,
                   each [Nte x Nte])
     .tr_targets - training labels (for classification) or values (for
                   regression) (column vector, [Ntr x 1])
     .use_kernel - flag, is data in form of kernel matrices (true) or in 
                   form of features (false)
    args     - argument string, where
       -h         - optimise hyperparameters (otherwise don't)
       -f iter    - max # iterations for optimiser (ignored if -h not set)
       -l likfun  - likelihood function:
                       'likErf' - erf/probit likelihood (binary only)
       -c covfun  - covariance function:
                       'covLINkcell' - simple dot product
                       'covLINglm'   - construct a GLM
       -m meanfun - mean function:
                       'meanConstcell' - suitable for dot product
                       'meanConstglm'  - suitable for GLM
       -i inffun  - inference function:
                       'prt_infEP' - Expectation Propagation
    experimental args (use at your own risk):
       -p         - use priors for the hyperparameters. If specified, this
                    indicates that a maximum a posteriori (MAP) approach
                    will be used to set covariance function
                    hyperparameters. The priors are obtained by calling
                    prt_gp_priors('covFuncName')

       N.B.: for the arguments specifying functions, pass in a string, not
       a function handle. This script will generate a function handle
 
 Output:
    output  - output of machine (struct).
     * Mandatory fields:
      .predictions - predictions of classification or regression [Nte x D]
     * Optional fields:
      .type     - which type of machine this is (here, 'classifier')
      .func_val - predictive probabilties
      .mu       - test latent means
      .s2       - test latent variances
      .loghyper - log hyperparameters
      .nlml     - negative log marginal likelihood
      .alpha    - GP weighting coefficients
      .sW       - likelihood matrix (see Rasmussen & Williams, 2006)
      .L        - Cholesky factor
__________________________________________________________________________
 Copyright (C) 2011 Machine Learning & Neuroimaging Laboratory

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function output = prt_machine_gpml(d,args)
0002 % Run Gaussian process model - wrapper for gpml toolbox
0003 % FORMAT output = prt_machine_gpml(d,args)
0004 % Inputs:
0005 %   d         - structure with data information, with mandatory fields:
0006 %     .train      - training data (cell array of matrices of row vectors,
0007 %                   each [Ntr x D]). each matrix contains one representation
0008 %                   of the data. This is useful for approaches such as
0009 %                   multiple kernel learning.
0010 %     .test       - testing data  (cell array of matrices row vectors, each
0011 %                   [Nte x D])
0012 %     .testcov    - testing covariance (cell array of matrices row vectors,
0013 %                   each [Nte x Nte])
0014 %     .tr_targets - training labels (for classification) or values (for
0015 %                   regression) (column vector, [Ntr x 1])
0016 %     .use_kernel - flag, is data in form of kernel matrices (true) or in
0017 %                   form of features (false)
0018 %    args     - argument string, where
0019 %       -h         - optimise hyperparameters (otherwise don't)
0020 %       -f iter    - max # iterations for optimiser (ignored if -h not set)
0021 %       -l likfun  - likelihood function:
0022 %                       'likErf' - erf/probit likelihood (binary only)
0023 %       -c covfun  - covariance function:
0024 %                       'covLINkcell' - simple dot product
0025 %                       'covLINglm'   - construct a GLM
0026 %       -m meanfun - mean function:
0027 %                       'meanConstcell' - suitable for dot product
0028 %                       'meanConstglm'  - suitable for GLM
0029 %       -i inffun  - inference function:
0030 %                       'prt_infEP' - Expectation Propagation
0031 %    experimental args (use at your own risk):
0032 %       -p         - use priors for the hyperparameters. If specified, this
0033 %                    indicates that a maximum a posteriori (MAP) approach
0034 %                    will be used to set covariance function
0035 %                    hyperparameters. The priors are obtained by calling
0036 %                    prt_gp_priors('covFuncName')
0037 %
0038 %       N.B.: for the arguments specifying functions, pass in a string, not
0039 %       a function handle. This script will generate a function handle
0040 %
0041 % Output:
0042 %    output  - output of machine (struct).
0043 %     * Mandatory fields:
0044 %      .predictions - predictions of classification or regression [Nte x D]
0045 %     * Optional fields:
0046 %      .type     - which type of machine this is (here, 'classifier')
0047 %      .func_val - predictive probabilties
0048 %      .mu       - test latent means
0049 %      .s2       - test latent variances
0050 %      .loghyper - log hyperparameters
0051 %      .nlml     - negative log marginal likelihood
0052 %      .alpha    - GP weighting coefficients
0053 %      .sW       - likelihood matrix (see Rasmussen & Williams, 2006)
0054 %      .L        - Cholesky factor
0055 %__________________________________________________________________________
0056 % Copyright (C) 2011 Machine Learning & Neuroimaging Laboratory
0057 
0058 % Written by A Marquand
0059 % $Id$
0060 
0061 % configure default parameters for GP optimisation
0062 meanfunc  = @meanConstcell;
0063 covfunc   = @covLINkcell; 
0064 maxeval   = -20;
0065 if strcmp(d.pred_type,'classification')
0066     mode = 'classifier'; % it's good to be consistent!
0067     likfunc   = @likErf;
0068     inffunc   = @prt_infEP;
0069 else
0070     mode = 'regression';
0071     likfunc   = @likGauss;
0072     inffunc   = @prt_infExact;
0073     mtr       = mean(d.tr_targets);      % mean of the training data
0074 end
0075 
0076 % Error checks
0077 % -------------------------------------------------------------------------
0078 SANITYCHECK=true; % can turn off for "speed". Expert only.
0079 
0080 if SANITYCHECK==true
0081     % args should be a string (empty or otherwise)
0082     if ~ischar(args)
0083         error('prt_machine_gpml:argsNotString',['Error: gpml'...
0084             ' args should be a string. ' ...
0085             ' SOLUTION: Please do XXX']);
0086     end
0087     
0088     % check we can reach the binary library
0089     if ~exist('prt_gp','file')
0090         error('prt_machine_gpml:libNotFound',['Error:'...
0091             ' ''prt_gp'' function could not be found !' ...
0092             ' SOLUTION: Please check your path.']);
0093     end
0094     % check whether it is a two-class classification problem
0095     uTL=unique(d.tr_targets(:));
0096     k=numel(uTL); % number of classes
0097     if strcmp(mode,'classifier') && k > 2
0098         warning('prt_machine_gpml:classificationWithMoreThanTwoClasses',...
0099                ['Classification specified with > 2 classes. ',...
0100                 'Defaulting to multiclass Laplace approximation.']);
0101         output = prt_machine_gpclap(d,args);
0102         return;
0103     end
0104     % are we using a kernel ?
0105     if ~d.use_kernel
0106         error('prt_machine_gpml:useKernelIsFalse',['Error:'...
0107             ' This machine is currently only implemented for kernel data ' ...
0108             'SOLUTION: Please set use_kernel to true']);
0109     end
0110 end
0111 
0112 % parse input arguments
0113 % -------------------------------------------------------------------------
0114 % hyperparameters
0115 if ~isempty(regexp(args,'-h','once'))
0116     optimise_theta = true;
0117     eargs = regexp(args,'-f\s+[0-9]*','match');
0118     if ~isempty(eargs)
0119         eargs = regexp(cell2mat(eargs),'-f\s+','split');
0120         maxeval  = str2num(['-',cell2mat(eargs(2))]);
0121     end
0122 else
0123     optimise_theta = false;
0124 end
0125 % likelihood function
0126 largs = regexp(args,'-l\s+[a-zA-Z0-9_]*','match');
0127 if ~isempty(largs)
0128     largs = regexp(cell2mat(largs),'-l\s+','split');
0129     likfunc = str2func(cell2mat(largs(2)));
0130     if strcmpi(cell2mat(largs(2)),'Erf')
0131         likfunc   = @likErf;
0132     end
0133 end
0134 % covariance function
0135 cargs = regexp(args,'-c\s+[a-zA-Z0-9_]*','match');
0136 if ~isempty(cargs)
0137     cargs = regexp(cell2mat(cargs),'-c\s+','split');
0138     covfunc = str2func(cell2mat(cargs(2)));
0139 end
0140 % mean function
0141 margs = regexp(args,'-m\s+[a-zA-Z0-9_]*','match');
0142 if ~isempty(margs)
0143     margs = regexp(cell2mat(margs),'-m\s+','split');
0144     meanfunc = str2func(cell2mat(margs(2)));
0145 end
0146 % inference function
0147 iargs = regexp(args,'-i\s+[a-zA-Z0-9_]*','match');
0148 if ~isempty(iargs)
0149     iargs = regexp(cell2mat(iargs),'-i\s+','split');
0150     inffunc = str2func(cell2mat(iargs(2)));
0151 end
0152 % priors
0153 if ~isempty(regexp(args,'-p','once'))
0154     disp('Empirical priors specified. Using MAP for hyperparameters')
0155     priors = prt_gp_priors(func2str(covfunc));
0156     map = true;
0157 else
0158     map = false;
0159 end
0160 
0161 % Set default hyperparameters
0162 % -------------------------------------------------------------------------
0163 nhyp = str2num([feval(covfunc); feval(likfunc); feval(meanfunc)]);
0164 if nhyp(1) > 0
0165     hyp.cov = zeros(nhyp(1),1);
0166 end
0167 if nhyp(2) > 0
0168     hyp.lik = zeros(nhyp(2),1);
0169 end
0170 if nhyp(3) > 0 
0171     hyp.mean = zeros(nhyp(3),1);
0172 end
0173 
0174 % Assemble data matrices
0175 % -------------------------------------------------------------------------
0176 % handle the glm as a special case (for now)
0177 if strcmpi(func2str(covfunc),'covLINglm') || strcmpi(func2str(covfunc),'covLINglm_2class')
0178     % configure covariances
0179     K   = [d.train(:)'   {d.tr_param}];
0180     Ks  = [d.test(:)'    {d.te_param}];
0181     Kss = [d.testcov(:)' {d.te_param}];
0182     
0183     % get default hyperparamter values
0184     hyp.cov = log(prt_glm_design);
0185     
0186     [tmp1 tmp2 tmp3 tr_lbs] = prt_glm_design(hyp.cov, d.tr_param);
0187     [tmp1 tmp2 tmp3 te_lbs] = prt_glm_design(hyp.cov, d.te_param);   
0188 else
0189     % configure covariances
0190     K   = d.train;
0191     Ks  = d.test;
0192     Kss = d.testcov;
0193         
0194     tr_lbs = d.tr_targets;
0195     te_lbs = d.te_targets;
0196 end
0197 
0198 % configure targets
0199 if strcmp(mode,'classifier')
0200     % convert targets to +1/-1
0201     y = -1*(2 * tr_lbs - 3);
0202 else
0203     y = tr_lbs - mtr;
0204 end
0205     
0206 % Train and test GP model
0207 % -------------------------------------------------------------------------
0208 % train
0209 if optimise_theta
0210     if map
0211         [hyp,nlmls] = minimize(hyp, @prt_gp_map, maxeval, inffunc, meanfunc, covfunc, likfunc, K, y, priors);
0212     else
0213         [hyp nlmls] = minimize(hyp, @prt_gp, maxeval, inffunc, meanfunc, covfunc, likfunc, K, y);
0214     end
0215 else
0216     nlmls = prt_gp(hyp, inffunc, meanfunc, covfunc, likfunc, K, y);
0217 end
0218 
0219 % make predictions
0220 [ymu ys2 fmu fs2 lp post] = prt_gp(hyp, inffunc, meanfunc, covfunc, likfunc,K, y, Ks, zeros(size(Ks{1},1),1), Kss);
0221 
0222 % Outputs
0223 % -------------------------------------------------------------------------
0224 if strcmp(mode,'classifier')
0225     p = exp(lp);
0226     output.predictions = (1-real(p > 0.5)) + 1;
0227     output.func_val    = p;
0228 else % regression
0229     output.predictions = ymu + mtr;
0230     output.func_val    = output.predictions;
0231 end
0232 output.type        = mode;
0233 output.loghyper    = hyp;
0234 output.mu          = ymu;
0235 output.s2          = ys2;
0236 output.nlml        = min(nlmls);
0237 output.tr_targets  = tr_lbs;
0238 output.te_targets  = te_lbs;
0239 output.alpha       = post.alpha;
0240 %output.sW          = post.sW;
0241 %output.L           = post.L;
0242 
0243 end
0244

Generated on Tue 10-Feb-2015 18:16:33 by m2html © 2005