PLearn 0.1
IIDNoiseKernel.cc
Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 // IIDNoiseKernel.cc
00004 //
00005 // Copyright (C) 2006-2007 Nicolas Chapados
00006 //
00007 // Redistribution and use in source and binary forms, with or without
00008 // modification, are permitted provided that the following conditions are met:
00009 //
00010 //  1. Redistributions of source code must retain the above copyright
00011 //     notice, this list of conditions and the following disclaimer.
00012 //
00013 //  2. Redistributions in binary form must reproduce the above copyright
00014 //     notice, this list of conditions and the following disclaimer in the
00015 //     documentation and/or other materials provided with the distribution.
00016 //
00017 //  3. The name of the authors may not be used to endorse or promote
00018 //     products derived from this software without specific prior written
00019 //     permission.
00020 //
00021 // THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR
00022 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
00023 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN
00024 // NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
00025 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
00026 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00027 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00028 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00029 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00030 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00031 //
00032 // This file is part of the PLearn library. For more information on the PLearn
00033 // library, go to the PLearn Web site at www.plearn.org
00034 
00035 // Authors: Nicolas Chapados
00036 
00039 #include "IIDNoiseKernel.h"
00040 #include <plearn/base/lexical_cast.h>
00041 #include <plearn/math/TMat_maths.h>
00042 
00043 namespace PLearn {
00044 using namespace std;
00045 
00046 PLEARN_IMPLEMENT_OBJECT(
00047     IIDNoiseKernel,
00048     "Kernel representing independent and identically-distributed observation noise",
00049     "This Kernel is typically used as a base class for covariance functions used\n"
00050     "in gaussian processes (see GaussianProcessRegressor).  It represents simple\n"
00051     "i.i.d. additive noise that applies to 'identical training cases' i and j:\n"
00052     "\n"
00053     "  k(D_i,D_j) = delta_i,j * sn\n"
00054     "\n"
00055     "where D_i and D_j are elements from the current data set (established by\n"
00056     "the setDataForKernelMatrix function), delta_i,j is the Kronecker delta\n"
00057     "function, and sn is softplus(isp_noise_sigma), with softplus(x) =\n"
00058     "log(1+exp(x)).  Note that 'identity' is not equivalent to 'vector\n"
00059     "equality': in particular, at test-time, this noise is NEVER added.\n"
00060     "Currently, two vectors are considered identical if and only if they are the\n"
00061     "SAME ROW of the current data set, and hence the noise term is added only at\n"
00062     "TRAIN-TIME across the diagonal of the Gram matrix (when the\n"
00063     "computeGramMatrix() function is called).  This is why at test-time, no such\n"
00064     "noise term is added.  The idea (see the book \"Gaussian Processes for\n"
00065     "Machine Learning\" by Rasmussen and Williams for details) is that\n"
00066     "observation noise only applies when A SPECIFIC OBSERVATION is drawn from\n"
00067     "the GP distribution: if we sample a new point at the same x, we will get a\n"
00068     "different realization for the noise, and hence the correlation between the\n"
00069     "two noise realizations is zero.  This class can only be sure that two\n"
00070     "observations are \"identical\" when they are presented all at once through\n"
00071     "the data matrix.\n"
00072     "\n"
00073     "The Kronecker terms computed by the base class are ADDDED to the noise\n"
00074     "computed by this kernel (at test-time also).\n"
00075     );
00076 
00077 
00078 IIDNoiseKernel::IIDNoiseKernel()
00079     : m_isp_noise_sigma(-100.0), /* very close to zero... */
00080       m_isp_kronecker_sigma(-100.0)
00081 { }
00082 
00083 
00084 //#####  declareOptions  ######################################################
00085 
00086 void IIDNoiseKernel::declareOptions(OptionList& ol)
00087 {
00088     declareOption(
00089         ol, "isp_noise_sigma", &IIDNoiseKernel::m_isp_noise_sigma,
00090         OptionBase::buildoption,
00091         "Inverse softplus of the global noise variance.  Default value=-100.0\n"
00092         "(very close to zero after we take softplus).");
00093 
00094     declareOption(
00095         ol, "isp_kronecker_sigma", &IIDNoiseKernel::m_isp_kronecker_sigma,
00096         OptionBase::buildoption,
00097         "Inverse softplus of the noise variance term for the product of\n"
00098         "Kronecker deltas associated with kronecker_indexes, if specified.");
00099     
00100     // Now call the parent class' declareOptions
00101     inherited::declareOptions(ol);
00102 }
00103 
00104 
00105 //#####  build  ###############################################################
00106 
00107 void IIDNoiseKernel::build()
00108 {
00109     // ### Nothing to add here, simply calls build_
00110     inherited::build();
00111     build_();
00112 }
00113 
00114 
00115 //#####  build_  ##############################################################
00116 
00117 void IIDNoiseKernel::build_()
00118 { }
00119 
00120 
00121 //#####  evaluate  ############################################################
00122 
00123 real IIDNoiseKernel::evaluate(const Vec& x1, const Vec& x2) const
00124 {
00125     // Assume that if x1 and x2 are identical, they are actually the same
00126     // instance of a data point.  This should not be called to compare a train
00127     // point against a test point (use evaluate_i_x for this purpose).
00128     return (x1 == x2? softplus(m_isp_noise_sigma) : 0.0) +
00129         softplus(m_isp_kronecker_sigma) * inherited::evaluate(x1,x2);
00130 }
00131 
00132 
00133 //#####  evaluate_i_x  ########################################################
00134 
00135 real IIDNoiseKernel::evaluate_i_x(int i, const Vec& x, real) const
00136 {
00137     // Noise component is ZERO between a test and any train example.
00138     // Just compute the Kronecker part is not necessarily zero
00139     Vec* train_row = dataRow(i);
00140     PLASSERT( train_row );
00141     return softplus(m_isp_kronecker_sigma) * inherited::evaluate(*train_row, x);
00142 }
00143 
00144 
00145 //#####  evaluate_all_i_x  ####################################################
00146 
00147 void IIDNoiseKernel::evaluate_all_i_x(const Vec& x, const Vec& k_xi_x,
00148                                       real, int istart) const
00149 {
00150     // Noise component is ZERO between a test and any train example
00151     k_xi_x.fill(0.0);
00152     real kronecker_sigma = softplus(m_isp_kronecker_sigma);
00153     int i_max = min(istart + k_xi_x.size(), data->length());
00154     int j = 0;
00155     for (int i=istart ; i<i_max ; ++i, ++j) {
00156         Vec* train_row = dataRow(i);
00157         k_xi_x[j] = kronecker_sigma * inherited::evaluate(*train_row, x);
00158     }
00159 }
00160 
00161 
00162 //#####  computeGramMatrix  ###################################################
00163 
00164 void IIDNoiseKernel::computeGramMatrix(Mat K) const
00165 {
00166     if (!data)
00167         PLERROR("Kernel::computeGramMatrix: setDataForKernelMatrix not yet called");
00168     if (!is_symmetric)
00169         PLERROR("Kernel::computeGramMatrix: not supported for non-symmetric kernels");
00170     if (K.length() != data.length() || K.width() != data.length())
00171         PLERROR("Kernel::computeGramMatrix: the argument matrix K should be\n"
00172                 "of size %d x %d (currently of size %d x %d)",
00173                 data.length(), data.length(), K.length(), K.width());
00174                 
00175     // Compute Kronecker gram matrix. Multiply by kronecker sigma if there were
00176     // any Kronecker terms.
00177     inherited::computeGramMatrix(K);
00178     if (m_kronecker_indexes.size() > 0)
00179         K *= softplus(m_isp_kronecker_sigma);
00180     
00181     // Add iid noise contribution
00182     real noise_sigma = softplus(m_isp_noise_sigma);
00183     int  l   = data->length();
00184     int  m   = K.mod() + 1;               // Mind the +1 to go along diagonal
00185     real *Ki = K[0];
00186     
00187     for (int i=0 ; i<l ; ++i, Ki += m) {
00188         *Ki += noise_sigma;
00189     }
00190 
00191     if (cache_gram_matrix) {
00192         gram_matrix.resize(l,l);
00193         gram_matrix << K;
00194         gram_matrix_is_cached = true;
00195     }
00196 }
00197 
00198 
00199 //#####  makeDeepCopyFromShallowCopy  #########################################
00200 
00201 void IIDNoiseKernel::makeDeepCopyFromShallowCopy(CopiesMap& copies)
00202 {
00203     inherited::makeDeepCopyFromShallowCopy(copies);
00204 }
00205 
00206 
00207 //#####  computeGramMatrixDerivative  #########################################
00208 
00209 void IIDNoiseKernel::computeGramMatrixDerivative(Mat& KD, const string& kernel_param,
00210                                                  real epsilon) const
00211 {
00212     static const string INS("isp_noise_sigma");
00213     static const string IKS("isp_kronecker_sigma");
00214 
00215     if (kernel_param == INS) {
00216         if (!data)
00217             PLERROR("Kernel::computeGramMatrixDerivative should be called only after "
00218                     "setDataForKernelMatrix");
00219 
00220         int W = nExamples();
00221         KD.resize(W,W);
00222         KD.fill(0.0);
00223         real deriv = sigmoid(m_isp_noise_sigma);
00224         for (int i=0 ; i<W ; ++i)
00225             KD(i,i) = deriv;
00226     }
00227     else if (kernel_param == IKS) {
00228         computeGramMatrixDerivKronecker(KD);
00229     }
00230     else
00231         inherited::computeGramMatrixDerivative(KD, kernel_param, epsilon);
00232 }
00233 
00234 
00235 //#####  computeGramMatrixDerivKronecker  #####################################
00236 
00237 void IIDNoiseKernel::computeGramMatrixDerivKronecker(Mat& KD) const
00238 {
00239     // From the cached version of the Gram matrix, this function is easily
00240     // implemented: we first copy the Gram to the KD matrix, subtract the IID
00241     // noise contribution from the main diagonal, and multiply the remaining
00242     // matrix (made up of 0/1 elements) by the derivative of the kronecker
00243     // sigma hyperparameter.
00244 
00245     int l = data->length();
00246     KD.resize(l,l);
00247     PLASSERT_MSG(gram_matrix.width() == l && gram_matrix.length() == l,
00248                  "To compute the derivative with respect to 'isp_kronecker_sigma',\n"
00249                  "the Gram matrix must be precomputed and cached in IIDNoiseKernel.");
00250     
00251     KD << gram_matrix;
00252     real noise_sigma = softplus(m_isp_noise_sigma);
00253     for (int i=0 ; i<l ; ++i)
00254         KD(i,i) -= noise_sigma;
00255 
00256     KD *= sigmoid(m_isp_kronecker_sigma) / softplus(m_isp_kronecker_sigma);
00257 }
00258 
00259 } // end of namespace PLearn
00260 
00261 
00262 /*
00263   Local Variables:
00264   mode:c++
00265   c-basic-offset:4
00266   c-file-style:"stroustrup"
00267   c-file-offsets:((innamespace . 0)(inline-open . 0))
00268   indent-tabs-mode:nil
00269   fill-column:79
00270   End:
00271 */
00272 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=79 :
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines