PLearn 0.1
ManifoldKNNDistribution.cc
Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 // ManifoldKNNDistribution.cc
00004 //
00005 // Copyright (C) 2007 Hugo Larochelle
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: Hugo Larochelle
00036 
00040 #include "ManifoldKNNDistribution.h"
00041 #include <plearn/math/pl_erf.h>
00042 #include <plearn/math/plapack.h>
00043 
00044 namespace PLearn {
00045 using namespace std;
00046 
00047 PLEARN_IMPLEMENT_OBJECT(
00048     ManifoldKNNDistribution,
00049     "K nearest neighbors density estimator locally taking into account the manifold",
00050     ""
00051     );
00052 
00054 // ManifoldKNNDistribution //
00056 ManifoldKNNDistribution::ManifoldKNNDistribution()
00057     : manifold_dimensionality(5),
00058       min_sigma_square(1e-5),
00059       center_around_manifold_neighbors(false),
00060       use_gaussian_distribution(false)
00061 {}
00062 
00064 // declareOptions //
00066 void ManifoldKNNDistribution::declareOptions(OptionList& ol)
00067 {
00068     declareOption(ol, "knn_manifold", &ManifoldKNNDistribution::knn_manifold,
00069                   OptionBase::buildoption,
00070                   "Nearest neighbors search algorithms for local manifold structure"
00071                   "estimation.");
00072 
00073     declareOption(ol, "knn_density", &ManifoldKNNDistribution::knn_density,
00074                   OptionBase::buildoption,
00075                   "Nearest neighbors search algorithms for density "
00076                   "estimation from ellipsoid\n"
00077                   "volume.");
00078 
00079     declareOption(ol, "manifold_dimensionality", 
00080                   &ManifoldKNNDistribution::manifold_dimensionality,
00081                   OptionBase::buildoption,
00082                   "Dimensionality of the manifold.");
00083 
00084     declareOption(ol, "min_sigma_square", 
00085                   &ManifoldKNNDistribution::min_sigma_square,
00086                   OptionBase::buildoption,
00087                   "Minimum variance in all directions on the manifold. This value"
00088                   "is added\n"
00089                   "to the estimated covariance matrix.");
00090 
00091     declareOption(ol, "center_around_manifold_neighbors", 
00092                   &ManifoldKNNDistribution::center_around_manifold_neighbors,
00093                   OptionBase::buildoption,
00094                   "Indication that the estimation of the manifold tangent vectors\n"
00095                   "should be made around the knn_manifold neighbors' mean vector,\n"
00096                   "not around the test point."
00097                   );
00098 
00099     declareOption(ol, "use_gaussian_distribution", 
00100                   &ManifoldKNNDistribution::use_gaussian_distribution,
00101                   OptionBase::buildoption,
00102                   "Indication that a Gaussian distribution should be used as the\n"
00103                   "knn_manifold nearest neighbors distribution, instead of the\n"
00104                   "uniform in the ellipsoid."
00105                   );
00106 
00107     declareOption(ol, "density_learner", 
00108                   &ManifoldKNNDistribution::density_learner,
00109                   OptionBase::buildoption,
00110                   "Generic density learner for knn_manifold nearest neighbors."
00111                   );
00112 
00113     // Now call the parent class' declareOptions().
00114     inherited::declareOptions(ol);
00115 }
00116 
00118 // build //
00120 void ManifoldKNNDistribution::build()
00121 {
00122     // ### Nothing to add here, simply calls build_().
00123     inherited::build();
00124     build_();
00125 }
00126 
00128 // build_ //
00130 void ManifoldKNNDistribution::build_()
00131 {
00132     if(min_sigma_square < 0)
00133         PLERROR("In ManifoldKNNDistribution::build_(): min_sigma_square should be"
00134                 " >= 0.");
00135 
00136     if(!knn_manifold)
00137         PLERROR("In ManifoldKNNDistribution::build_(): knn_manifold must be"
00138                 " provided.");
00139 
00140     if(!knn_density)
00141         PLERROR("In ManifoldKNNDistribution::build_(): knn_density must be"
00142                 " provided.");
00143 
00144     if(inputsize_ > 0)
00145     {
00146         if(manifold_dimensionality > inputsize_)
00147             manifold_dimensionality = inputsize_;
00148         if(manifold_dimensionality < 1)
00149             PLERROR("In ManifoldKNNDistribution::build_(): manifold_dimensionality"
00150                     " should be > 0.");
00151         eig_vectors.resize(manifold_dimensionality,inputsize_);
00152         eig_values.resize(manifold_dimensionality);
00153         Ut.resize(inputsize_,inputsize_);
00154         V.resize(knn_manifold->num_neighbors,knn_manifold->num_neighbors);
00155         eig_vectors_projection.resize(manifold_dimensionality);
00156         neighbors_mean.resize(inputsize_);
00157     }
00158 
00159     if(train_set)
00160     {
00161         knn_manifold->setTrainingSet(train_set,true);
00162         knn_density->setTrainingSet(train_set,true);
00163         
00164         knn_manifold->train();
00165         knn_density->train();
00166     }
00167 
00168     if(use_gaussian_distribution && !center_around_manifold_neighbors)
00169     {
00170         PLWARNING("In ManifoldKNNDistribution::build_(): when using "
00171                   "use_gaussian_distribution=true, center_around_manifold_neighbors"
00172                   "must be true too. Setting center_around_manifold_neighbors=true...");
00173         center_around_manifold_neighbors = true;
00174     }
00175 }
00176 
00178 // cdf //
00180 real ManifoldKNNDistribution::cdf(const Vec& y) const
00181 {
00182     PLERROR("cdf not implemented for ManifoldKNNDistribution"); return 0;
00183 }
00184 
00186 // expectation //
00188 void ManifoldKNNDistribution::expectation(Vec& mu) const
00189 {
00190     PLERROR("expectation not implemented for ManifoldKNNDistribution");
00191 }
00192 
00193 // ### Remove this method if your distribution does not implement it.
00195 // forget //
00197 void ManifoldKNNDistribution::forget()
00198 {}
00199 
00201 // generate //
00203 void ManifoldKNNDistribution::generate(Vec& y) const
00204 {
00205     PLERROR("generate not implemented for ManifoldKNNDistribution");
00206 }
00207 
00209 // log_density //
00211 real ManifoldKNNDistribution::log_density(const Vec& y) const
00212 {
00213 
00214     real ret = 0;
00215     if(density_learner)
00216     {
00217         knn_manifold->computeOutput(y,nearest_neighbors_manifold_vec);
00218         nearest_neighbors_manifold = 
00219             nearest_neighbors_manifold_vec.toMat(
00220                 knn_manifold->num_neighbors,inputsize_);
00221         density_learner_train_set = VMat( nearest_neighbors_manifold );
00222         density_learner_train_set->defineSizes(inputsize_,0);
00223         density_learner->setTrainingSet(density_learner_train_set,true);
00224         density_learner->train();
00225         density_learner->log_density(y);
00226         ret = density_learner->log_density(y) + 
00227             pl_log((real)knn_manifold->num_neighbors)-pl_log((real)n_examples);
00228     }
00229     else if(use_gaussian_distribution)
00230     {
00231         computeLocalPrincipalComponents(y,eig_values,eig_vectors);
00232         if(!center_around_manifold_neighbors)
00233             columnMean(nearest_neighbors_manifold,neighbors_mean);
00234 
00235         // Compute log-normalization constant
00236         ret = - inputsize_ *0.5 *Log2Pi;
00237         for(int i=0; i<manifold_dimensionality; i++)
00238             ret -= 0.5 * pl_log(eig_values[i]+min_sigma_square);
00239         ret -= (inputsize_-manifold_dimensionality)*0.5*pl_log(min_sigma_square);
00240         substract(y,neighbors_mean,test_minus_mean);
00241         product(eig_vectors_projection,eig_vectors,test_minus_mean);
00242         for(int j=0; j<eig_values.length(); j++)
00243             ret -= mypow(eig_vectors_projection[j],2) * 
00244                 (1/(eig_values[j]+min_sigma_square) 
00245                  - 1/min_sigma_square) ;
00246         ret -= pownorm(test_minus_mean,2)
00247             /min_sigma_square;
00248         ret += pl_log((real)knn_manifold->num_neighbors)-pl_log((real)n_examples);
00249     }
00250     else
00251     {
00252         computeLocalPrincipalComponents(y,eig_values,eig_vectors);
00253 
00254         // Find volume of ellipsoid defined by eig_values, eig_vectors and
00255         // min_sigma_square that covers all the nearest_neighbors found by knn_density
00256         knn_density->computeOutput(y,nearest_neighbors_density_vec);
00257         nearest_neighbors_density = 
00258             nearest_neighbors_density_vec.toMat(knn_density->num_neighbors,inputsize_);
00259         nearest_neighbors_density -= y;
00260         real max = -1;
00261         real scaled_projection=0;
00262         for(int i=0; i<nearest_neighbors_density.length(); i++)
00263         {
00264             scaled_projection = 0;
00265             product(eig_vectors_projection,eig_vectors,nearest_neighbors_density(i));
00266             for(int j=0; j<eig_values.length(); j++)
00267                 scaled_projection += mypow(eig_vectors_projection[j],2) * 
00268                     (1/(eig_values[j]+min_sigma_square) 
00269                      - 1/min_sigma_square) ;
00270             scaled_projection += pownorm(nearest_neighbors_density(i),2)
00271                 /min_sigma_square;
00272             if(max < scaled_projection)
00273                 max = scaled_projection;
00274         }
00275         
00276         // Compute log-volume of the ellipsoid: pi
00277         real log_vol = 0.5 * inputsize_ * pl_log(scaled_projection);
00278         for(int i=0; i<manifold_dimensionality; i++)
00279             log_vol += 0.5 * pl_log(eig_values[i]+min_sigma_square);
00280         log_vol += (inputsize_-manifold_dimensionality)*0.5*pl_log(min_sigma_square);
00281         log_vol += 0.5*inputsize_*pl_log(Pi) - pl_gammln(0.5*inputsize_+1);
00282         
00283         ret = pl_log((real)knn_density->num_neighbors)-pl_log((real)n_examples)-log_vol;
00284     }
00285     return ret;
00286 }
00287 
00289 // makeDeepCopyFromShallowCopy //
00291 void ManifoldKNNDistribution::makeDeepCopyFromShallowCopy(CopiesMap& copies)
00292 {
00293     inherited::makeDeepCopyFromShallowCopy(copies);
00294 
00295     deepCopyField(knn_manifold, copies);
00296     deepCopyField(knn_density, copies);
00297     deepCopyField(density_learner, copies);
00298     deepCopyField(nearest_neighbors_manifold, copies);
00299     deepCopyField(nearest_neighbors_manifold_vec, copies);
00300     deepCopyField(nearest_neighbors_density, copies);
00301     deepCopyField(nearest_neighbors_density_vec, copies);
00302     deepCopyField(eig_vectors, copies);
00303     deepCopyField(eig_values, copies);
00304     deepCopyField(Ut, copies);
00305     deepCopyField(V, copies);
00306     deepCopyField(S, copies);
00307     deepCopyField(eig_vectors_projection, copies);
00308     deepCopyField(neighbors_mean,copies);
00309     deepCopyField(test_minus_mean,copies);
00310     deepCopyField(density_learner_train_set,copies);
00311 }
00312 
00314 // resetGenerator //
00316 void ManifoldKNNDistribution::resetGenerator(long g_seed) const
00317 {
00318     PLERROR("resetGenerator not implemented for ManifoldKNNDistribution");
00319 }
00320 
00322 // survival_fn //
00324 real ManifoldKNNDistribution::survival_fn(const Vec& y) const
00325 {
00326     PLERROR("survival_fn not implemented for ManifoldKNNDistribution"); return 0;
00327 }
00328 
00329 // ### Remove this method, if your distribution does not implement it.
00331 // train //
00333 void ManifoldKNNDistribution::train()
00334 {}
00335 
00337 // variance //
00339 void ManifoldKNNDistribution::variance(Mat& covar) const
00340 {
00341     PLERROR("variance not implemented for ManifoldKNNDistribution");
00342 }
00343 
00344 void ManifoldKNNDistribution::computeLocalPrincipalComponents(const Vec& x, 
00345                                          Vec& eig_values, Mat& eig_vectors) const
00346 {
00347     knn_manifold->computeOutput(x,nearest_neighbors_manifold_vec);
00348     nearest_neighbors_manifold = 
00349         nearest_neighbors_manifold_vec.toMat(knn_manifold->num_neighbors,inputsize_);
00350 
00351     if(center_around_manifold_neighbors)
00352     {
00353         columnMean(nearest_neighbors_manifold,neighbors_mean);
00354         nearest_neighbors_manifold -= neighbors_mean;
00355     }
00356     else
00357         nearest_neighbors_manifold -= x;
00358     
00359     // Compute principal components
00360     // N.B. this is the SVD of F'
00361     lapackSVD(nearest_neighbors_manifold, Ut, S, V,'A',1.5);
00362     for (int k=0;k<manifold_dimensionality;k++)
00363     {
00364         eig_values[k] = mypow(S[k],2);
00365         eig_vectors(k) << Ut(k);
00366     }  
00367 }
00368 
00369 } // end of namespace PLearn
00370 
00371 
00372 /*
00373   Local Variables:
00374   mode:c++
00375   c-basic-offset:4
00376   c-file-style:"stroustrup"
00377   c-file-offsets:((innamespace . 0)(inline-open . 0))
00378   indent-tabs-mode:nil
00379   fill-column:79
00380   End:
00381 */
00382 // 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