PLearn 0.1
|
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 :