PLearn 0.1
|
00001 // -*- C++ -*- 00002 00003 // NonLocalManifoldParzen.cc 00004 // 00005 // Copyright (C) 2004 Yoshua Bengio & 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 /* ******************************************************* 00036 * $Id: NonLocalManifoldParzen.cc 8172 2007-10-10 22:41:33Z larocheh $ 00037 ******************************************************* */ 00038 00039 // Authors: Yoshua Bengio & Martin Monperrus 00040 00044 #include "NonLocalManifoldParzen.h" 00045 #include <plearn/display/DisplayUtils.h> 00046 #include <plearn/math/plapack.h> 00047 #include <plearn/var/AffineTransformVariable.h> 00048 #include <plearn/var/AffineTransformWeightPenalty.h> 00049 #include <plearn/var/ColumnSumVariable.h> 00050 #include <plearn/var/NllGeneralGaussianVariable.h> 00051 #include <plearn/var/NoBpropVariable.h> 00052 #include <plearn/var/ReshapeVariable.h> 00053 #include <plearn/var/SourceVariable.h> 00054 #include <plearn/var/SquareVariable.h> 00055 #include <plearn/var/SumOfVariable.h> 00056 #include <plearn/var/TanhVariable.h> 00057 #include <plearn/var/ThresholdBpropVariable.h> 00058 #include <plearn/var/Var_operators.h> 00059 #include <plearn/vmat/AppendNeighborsVMatrix.h> 00060 #include <plearn/vmat/ConcatColumnsVMatrix.h> 00061 00062 00063 namespace PLearn { 00064 using namespace std; 00065 00066 00067 NonLocalManifoldParzen::NonLocalManifoldParzen() 00068 : 00069 reference_set(0), 00070 ncomponents(1), 00071 nneighbors(5), 00072 nneighbors_density(-1), 00073 store_prediction(false), 00074 learn_mu(false), 00075 sigma_init(0.1), 00076 sigma_min(0.00001), 00077 mu_nneighbors(2), 00078 sigma_threshold_factor(-1), 00079 svd_threshold(1e-8), 00080 nhidden(10), 00081 weight_decay(0), 00082 penalty_type("L2_square"), 00083 batch_size(1) 00084 { 00085 } 00086 00087 PLEARN_IMPLEMENT_OBJECT(NonLocalManifoldParzen, 00088 "Non-Local version of Manifold Parzen Windows", 00089 "Manifold Parzen Windows density model, where the\n" 00090 "parameters of the kernel for each training point\n" 00091 "are predicted by a neural network.\n" 00092 ); 00093 00094 00095 void NonLocalManifoldParzen::declareOptions(OptionList& ol) 00096 { 00097 00098 declareOption(ol, "parameters", &NonLocalManifoldParzen::parameters, 00099 OptionBase::learntoption, 00100 "Parameters of the tangent_predictor function.\n" 00101 ); 00102 00103 declareOption(ol, "reference_set", &NonLocalManifoldParzen::reference_set, 00104 OptionBase::learntoption, 00105 "Reference points for density computation.\n" 00106 ); 00107 00108 declareOption(ol, "ncomponents", &NonLocalManifoldParzen::ncomponents, 00109 OptionBase::buildoption, 00110 "Number of \"principal components\" to predict\n" 00111 "for kernel parameters prediction.\n" 00112 ); 00113 00114 declareOption(ol, "nneighbors", &NonLocalManifoldParzen::nneighbors, 00115 OptionBase::buildoption, 00116 "Number of nearest neighbors to consider in training procedure.\n" 00117 ); 00118 00119 declareOption(ol, "nneighbors_density", 00120 &NonLocalManifoldParzen::nneighbors_density, 00121 OptionBase::buildoption, 00122 "Number of nearest neighbors to consider for\n" 00123 "p(x) density estimation.\n" 00124 ); 00125 00126 declareOption(ol, "store_prediction", 00127 &NonLocalManifoldParzen::store_prediction, 00128 OptionBase::buildoption, 00129 "Indication that the predicted parameters should be stored.\n" 00130 "This may make testing faster. Note that the predictions are\n" 00131 "stored after the last training stage\n" 00132 ); 00133 00134 00135 declareOption(ol, "paramsvalues", 00136 &NonLocalManifoldParzen::paramsvalues, 00137 OptionBase::learntoption, 00138 "The learned parameter vector.\n" 00139 ); 00140 00141 // ** Gaussian kernel options 00142 00143 declareOption(ol, "learn_mu", &NonLocalManifoldParzen::learn_mu, 00144 OptionBase::buildoption, 00145 "Indication that the deviation from the training point\n" 00146 "in a Gaussian kernel (called mu) should be learned.\n" 00147 ); 00148 00149 declareOption(ol, "sigma_init", &NonLocalManifoldParzen::sigma_init, 00150 OptionBase::buildoption, 00151 "Initial minimum value for sigma noise.\n" 00152 ); 00153 00154 declareOption(ol, "sigma_min", &NonLocalManifoldParzen::sigma_min, 00155 OptionBase::buildoption, 00156 "The minimum value for sigma noise.\n" 00157 ); 00158 00159 declareOption(ol, "mu_nneighbors", &NonLocalManifoldParzen::mu_nneighbors, 00160 OptionBase::buildoption, 00161 "Number of nearest neighbors to learn the mus \n" 00162 "(if < 0, mu_nneighbors = nneighbors).\n" 00163 ); 00164 00165 declareOption(ol, "sigma_threshold_factor", 00166 &NonLocalManifoldParzen::sigma_threshold_factor, 00167 OptionBase::buildoption, 00168 "Threshold factor of the gradient on the sigma noise\n" 00169 "parameter of the Gaussian kernel. If < 0, then\n" 00170 "no threshold is used." 00171 ); 00172 00173 declareOption(ol, "svd_threshold", 00174 &NonLocalManifoldParzen::svd_threshold, OptionBase::buildoption, 00175 "Threshold to accept singular values of F in solving for\n" 00176 "linear combination weights on tangent subspace.\n" 00177 ); 00178 00179 // ** Neural network predictor ** 00180 00181 declareOption(ol, "nhidden", 00182 &NonLocalManifoldParzen::nhidden, OptionBase::buildoption, 00183 "Number of hidden units of the neural network.\n" 00184 ); 00185 00186 declareOption(ol, "weight_decay", &NonLocalManifoldParzen::weight_decay, 00187 OptionBase::buildoption, 00188 "Global weight decay for all layers.\n"); 00189 00190 declareOption(ol, "penalty_type", &NonLocalManifoldParzen::penalty_type, 00191 OptionBase::buildoption, 00192 "Penalty to use on the weights (for weight and bias decay).\n" 00193 "Can be any of:\n" 00194 " - \"L1\": L1 norm,\n" 00195 " - \"L2_square\" (default): square of the L2 norm.\n"); 00196 00197 declareOption(ol, "optimizer", &NonLocalManifoldParzen::optimizer, 00198 OptionBase::buildoption, 00199 "Optimizer that optimizes the cost function.\n" 00200 ); 00201 00202 declareOption(ol, "batch_size", 00203 &NonLocalManifoldParzen::batch_size, OptionBase::buildoption, 00204 "How many samples to use to estimate the average gradient\n" 00205 "before updating the weights. If <= 0, is equivalent to\n" 00206 "specifying training_set->length() \n"); 00207 00208 00209 // ** Stored outputs of neural network 00210 00211 declareOption(ol, "mus", 00212 &NonLocalManifoldParzen::mus, OptionBase::learntoption, 00213 "The stored mu vectors for the reference set.\n" 00214 ); 00215 00216 declareOption(ol, "sns", &NonLocalManifoldParzen::sns, 00217 OptionBase::learntoption, 00218 "The stored sigma noise values for the reference set.\n" 00219 ); 00220 00221 declareOption(ol, "sms", &NonLocalManifoldParzen::sms, 00222 OptionBase::learntoption, 00223 "The stored sigma manifold values for the reference set.\n" 00224 ); 00225 00226 declareOption(ol, "Fs", &NonLocalManifoldParzen::Fs, OptionBase::learntoption, 00227 "The storaged \"principal components\" (F) values for\n" 00228 "the reference set.\n" 00229 ); 00230 00231 00232 // Now call the parent class' declareOptions 00233 inherited::declareOptions(ol); 00234 } 00235 00236 void NonLocalManifoldParzen::build_() 00237 { 00238 00239 if (inputsize_>0) 00240 { 00241 if (nhidden <= 0) 00242 PLERROR("NonLocalManifoldParzen::Number of hidden units " 00243 "should be positive, now %d\n",nhidden); 00244 00245 Var log_n_examples(1,1,"log(n_examples)"); 00246 if(train_set) 00247 { 00248 L = train_set->length(); 00249 reference_set = train_set; 00250 } 00251 00252 log_L= pl_log((real) L); 00253 parameters.resize(0); 00254 00255 // Neural network prediction of principal components 00256 00257 x = Var(inputsize_); 00258 x->setName("x"); 00259 00260 W = Var(nhidden+1,inputsize_,"W"); 00261 parameters.append(W); 00262 00263 Var a; // outputs of hidden layer 00264 a = affine_transform(x,W); 00265 a->setName("a"); 00266 00267 V = Var(ncomponents*(inputsize_+1),nhidden,"V"); 00268 parameters.append(V); 00269 00270 // TODO: instead, make NllGeneralGaussianVariable use vector... (DONE) 00271 //components = reshape(affine_transform(V,a),ncomponents,n); 00272 components = affine_transform(V,a); 00273 components->setName("components"); 00274 00275 // Gaussian kernel parameters prediction 00276 00277 muV = Var(inputsize_+1,nhidden,"muV"); 00278 snV = Var(2,nhidden,"snV"); 00279 00280 parameters.append(muV); 00281 parameters.append(snV); 00282 00283 if(learn_mu) 00284 mu = affine_transform(muV,a); 00285 else 00286 { 00287 mu = new SourceVariable(inputsize_,1); 00288 mu->value.clear(); 00289 } 00290 mu->setName("mu"); 00291 00292 min_sig = new SourceVariable(1,1); 00293 min_sig->value[0] = sigma_min; 00294 min_sig->setName("min_sig"); 00295 init_sig = Var(1,1); 00296 init_sig->setName("init_sig"); 00297 parameters.append(init_sig); 00298 00299 sn = square(affine_transform(snV,a)) + min_sig + square(init_sig); 00300 sn->setName("sn"); 00301 00302 if(sigma_threshold_factor > 0) 00303 sn = threshold_bprop(sn,sigma_threshold_factor); 00304 00305 predictor = Func(x, parameters , components & mu & sn ); 00306 00307 Var target_index = Var(1,1); 00308 target_index->setName("target_index"); 00309 Var neighbor_indexes = Var(nneighbors,1); 00310 neighbor_indexes->setName("neighbor_indexes"); 00311 00312 tangent_targets = Var(nneighbors,inputsize_); 00313 if(mu_nneighbors < 0 ) mu_nneighbors = nneighbors; 00314 00315 Var nll; 00316 nll = nll_general_gaussian(components, mu, sn, tangent_targets, 00317 log_L, learn_mu, mu_nneighbors); 00318 00319 Var knn = new SourceVariable(1,1); 00320 knn->setName("knn"); 00321 knn->value[0] = nneighbors; 00322 sum_nll = new ColumnSumVariable(nll) / knn; 00323 00324 // Weight decay penalty 00325 if(weight_decay > 0 ) 00326 { 00327 sum_nll += affine_transform_weight_penalty( 00328 W,weight_decay,0,penalty_type) + 00329 affine_transform_weight_penalty( 00330 V,weight_decay,0,penalty_type) + 00331 affine_transform_weight_penalty( 00332 muV,weight_decay,0,penalty_type) + 00333 affine_transform_weight_penalty( 00334 snV,weight_decay,0,penalty_type); 00335 } 00336 00337 cost_of_one_example = Func(x & tangent_targets & target_index & 00338 neighbor_indexes, parameters, sum_nll); 00339 00340 if(nneighbors_density >= L || nneighbors_density < 0) 00341 nneighbors_density = L; 00342 00343 // Output storage variables 00344 t_row.resize(inputsize_); 00345 Ut_svd.resize(inputsize_,inputsize_); 00346 V_svd.resize(ncomponents,ncomponents); 00347 F.resize(components->length(),components->width()); 00348 z.resize(inputsize_); 00349 x_minus_neighbor.resize(inputsize_); 00350 neighbor_row.resize(inputsize_); 00351 00352 // log_density and Kernel methods variables 00353 U_temp.resize(ncomponents,inputsize_); 00354 mu_temp.resize(inputsize_); 00355 sm_temp.resize(ncomponents); 00356 sn_temp.resize(1); 00357 diff.resize(inputsize_); 00358 00359 mus.resize(L, inputsize_); 00360 sns.resize(L); 00361 sms.resize(L,ncomponents); 00362 Fs.resize(L); 00363 for(int i=0; i<L; i++) 00364 { 00365 Fs[i].resize(ncomponents,inputsize_); 00366 } 00367 00368 if(paramsvalues.length() == parameters.nelems()) 00369 parameters << paramsvalues; 00370 else 00371 { 00372 paramsvalues.resize(parameters.nelems()); 00373 initializeParams(); 00374 if(optimizer) 00375 optimizer->reset(); 00376 } 00377 parameters.makeSharedValue(paramsvalues); 00378 } 00379 00380 } 00381 00382 void NonLocalManifoldParzen::knn(const VMat& vm, const Vec& x, const int& k, TVec<int>& neighbors, bool sortk) const 00383 { 00384 int n = vm->length(); 00385 distances.resize(n,2); 00386 distances.column(1) << Vec(0, n-1, 1); 00387 dk.setDataForKernelMatrix(vm); 00388 t_dist.resize(n); 00389 dk.evaluate_all_i_x(x, t_dist); 00390 distances.column(0) << t_dist; 00391 partialSortRows(distances, k, sortk); 00392 neighbors.resize(k); 00393 00394 for (int i=0; i < k && i<n; i++) 00395 { 00396 neighbors[i] = int(distances(i,1)); 00397 } 00398 00399 } 00400 00401 // ### Nothing to add here, simply calls build_ 00402 void NonLocalManifoldParzen::build() 00403 { 00404 inherited::build(); 00405 build_(); 00406 } 00407 00408 #ifdef __INTEL_COMPILER 00409 #pragma warning(disable:1419) // Get rid of compiler warning. 00410 #endif 00411 extern void varDeepCopyField(Var& field, CopiesMap& copies); 00412 #ifdef __INTEL_COMPILER 00413 #pragma warning(default:1419) 00414 #endif 00415 00416 void NonLocalManifoldParzen::makeDeepCopyFromShallowCopy(CopiesMap& copies) 00417 { 00418 inherited::makeDeepCopyFromShallowCopy(copies); 00419 00420 // Protected 00421 00422 deepCopyField(cost_of_one_example, copies); 00423 varDeepCopyField(x, copies); 00424 varDeepCopyField(W, copies); 00425 varDeepCopyField(V, copies); 00426 varDeepCopyField(muV, copies); 00427 varDeepCopyField(snV, copies); 00428 varDeepCopyField(tangent_targets, copies); 00429 varDeepCopyField(components, copies); 00430 varDeepCopyField(mu, copies); 00431 varDeepCopyField(sn, copies); 00432 varDeepCopyField(sum_nll, copies); 00433 varDeepCopyField(min_sig, copies); 00434 varDeepCopyField(init_sig, copies); 00435 deepCopyField(predictor, copies); 00436 deepCopyField(U_temp,copies); 00437 deepCopyField(F, copies); 00438 deepCopyField(distances,copies); 00439 deepCopyField(mu_temp,copies); 00440 deepCopyField(sm_temp,copies); 00441 deepCopyField(sn_temp,copies); 00442 deepCopyField(diff,copies); 00443 deepCopyField(z,copies); 00444 deepCopyField(x_minus_neighbor,copies); 00445 deepCopyField(t_row,copies); 00446 deepCopyField(neighbor_row,copies); 00447 deepCopyField(log_gauss,copies); 00448 deepCopyField(t_dist,copies); 00449 deepCopyField(t_nn,copies); 00450 deepCopyField(Ut_svd, copies); 00451 deepCopyField(V_svd, copies); 00452 deepCopyField(S_svd, copies); 00453 deepCopyField(mus, copies); 00454 deepCopyField(sns, copies); 00455 deepCopyField(sms, copies); 00456 deepCopyField(Fs, copies); 00457 deepCopyField(train_set_with_targets, copies); 00458 deepCopyField(targets_vmat, copies); 00459 varDeepCopyField(totalcost, copies); 00460 deepCopyField(paramsvalues, copies); 00461 00462 // Public 00463 00464 deepCopyField(parameters, copies); 00465 deepCopyField(reference_set,copies); 00466 deepCopyField(optimizer, copies); 00467 00468 } 00469 00470 00471 void NonLocalManifoldParzen::forget() 00472 { 00473 inherited::forget(); 00474 if (train_set) initializeParams(); 00475 if(optimizer) optimizer->reset(); 00476 stage = 0; 00477 } 00478 00479 void NonLocalManifoldParzen::train() 00480 { 00481 // Check whether gradient descent is going to be done 00482 // If not, then we don't need to store the parameters, 00483 // except for sn... 00484 bool flag = (nstages == stage); 00485 00486 // Update sigma_min, in case it was changed, 00487 // e.g. using an HyperLearner 00488 min_sig->value[0] = sigma_min; 00489 00490 // Set train_stats if not already done. 00491 if (!train_stats) 00492 train_stats = new VecStatsCollector(); 00493 00494 if (!cost_of_one_example) 00495 PLERROR("NonLocalManifoldParzen::train: build has not been run after setTrainingSet!"); 00496 00497 if(stage == 0) 00498 { 00499 targets_vmat = append_neighbors( 00500 train_set, nneighbors, true); 00501 nsamples = batch_size>0 ? batch_size : train_set->length(); 00502 00503 totalcost = meanOf(train_set_with_targets, cost_of_one_example, nsamples); 00504 00505 if(optimizer) 00506 { 00507 optimizer->setToOptimize(parameters, totalcost); 00508 optimizer->build(); 00509 } 00510 else PLERROR("NonLocalManifoldParzen::train can't train without setting an optimizer first!"); 00511 } 00512 00513 int optstage_per_lstage = train_set->length()/nsamples; 00514 00515 PP<ProgressBar> pb; 00516 if(report_progress>0) 00517 pb = new ProgressBar("Training NonLocalManifoldParzen from stage " + tostring(stage) + " to " + tostring(nstages), nstages-stage); 00518 00519 t_row.resize(train_set.width()); 00520 00521 int initial_stage = stage; 00522 bool early_stop=false; 00523 while(stage<nstages && !early_stop) 00524 { 00525 optimizer->nstages = optstage_per_lstage; 00526 train_stats->forget(); 00527 optimizer->early_stop = false; 00528 optimizer->optimizeN(*train_stats); 00529 train_stats->finalize(); 00530 if(verbosity>2) 00531 cout << "Epoch " << stage << " train objective: " << train_stats->getMean() << endl; 00532 ++stage; 00533 if(pb) 00534 pb->update(stage-initial_stage); 00535 } 00536 if(verbosity>1) 00537 cout << "EPOCH " << stage << " train objective: " << train_stats->getMean() << endl; 00538 00539 if(store_prediction && !flag) 00540 { 00541 for(int t=0; t<L;t++) 00542 { 00543 reference_set->getRow(t,neighbor_row); 00544 predictor->fprop(neighbor_row, F.toVec() & mus(t) & sns.subVec(t,1)); 00545 // N.B. this is the SVD of F' 00546 lapackSVD(F, Ut_svd, S_svd, V_svd,'A',1.5); 00547 for (int k=0;k<ncomponents;k++) 00548 { 00549 sms(t,k) = mypow(S_svd[k],2); 00550 Fs[t](k) << Ut_svd(k); 00551 } 00552 sns[t] += sigma_min - min_sig->value[0]; 00553 } 00554 } 00555 } 00556 00558 // initializeParams // 00560 void NonLocalManifoldParzen::initializeParams() 00561 { 00562 real delta = 1.0 / sqrt(real(inputsize_)); 00563 random_gen->fill_random_uniform(W->value, -delta, delta); 00564 delta = 1.0 / real(nhidden); 00565 random_gen->fill_random_uniform(V->matValue, -delta, delta); 00566 random_gen->fill_random_uniform(snV->matValue, -delta, delta); 00567 random_gen->fill_random_uniform(muV->matValue, -delta, delta); 00568 W->matValue(0).clear(); 00569 V->matValue(0).clear(); 00570 muV->matValue(0).clear(); 00571 snV->matValue(0).clear(); 00572 init_sig->value[0] = sqrt(sigma_init); 00573 } 00574 00576 // log_density // 00578 real NonLocalManifoldParzen::log_density(const Vec& x) const { 00579 // Compute log-density. 00580 real ret = 0; 00581 t_row << x; 00582 real mahal = 0; 00583 real norm_term = 0; 00584 00585 // Update sigma_min, in case it was changed, 00586 // e.g. using an HyperLearner 00587 00588 if(store_prediction && min_sig->value[0] != sigma_min) 00589 { 00590 for(int i=0; i<L; i++) 00591 { 00592 sns[i] += sigma_min - min_sig->value[0]; 00593 } 00594 } 00595 00596 min_sig->value[0] = sigma_min; 00597 00598 if(nneighbors_density != L) 00599 { 00600 // Fetching nearest neighbors for density estimation. 00601 knn(reference_set,x,nneighbors_density,t_nn,0); 00602 log_gauss.resize(t_nn.length()); 00603 for(int neighbor=0; neighbor<t_nn.length(); neighbor++) 00604 { 00605 reference_set->getRow(t_nn[neighbor],neighbor_row); 00606 if(!store_prediction) 00607 { 00608 predictor->fprop(neighbor_row, F.toVec() & mu_temp & sn_temp); 00609 // N.B. this is the SVD of F' 00610 lapackSVD(F, Ut_svd, S_svd, V_svd,'A',1.5); 00611 for (int k=0;k<ncomponents;k++) 00612 { 00613 sm_temp[k] = mypow(S_svd[k],2); 00614 U_temp(k) << Ut_svd(k); 00615 } 00616 } 00617 else 00618 { 00619 if(learn_mu) 00620 mu_temp << mus(t_nn[neighbor]); 00621 sn_temp[0] = sns[t_nn[neighbor]]; 00622 sm_temp << sms(t_nn[neighbor]); 00623 U_temp << Fs[t_nn[neighbor]]; 00624 } 00625 if(learn_mu) 00626 { 00627 substract(t_row,neighbor_row,x_minus_neighbor); 00628 substract(x_minus_neighbor,mu_temp,z); 00629 } 00630 else 00631 substract(t_row,neighbor_row,z); 00632 00633 mahal = -0.5*pownorm(z)/sn_temp[0]; 00634 norm_term = - inputsize_/2.0 * Log2Pi 00635 - log_L - 0.5*(inputsize_-ncomponents)*pl_log(sn_temp[0]); 00636 00637 00638 for(int k=0; k<ncomponents; k++) 00639 { 00640 mahal -= square(dot(z,U_temp(k)))*(0.5/(sm_temp[k]+sn_temp[0]) 00641 - 0.5/sn_temp[0]); 00642 norm_term -= 0.5*pl_log(sm_temp[k]+sn_temp[0]); 00643 } 00644 00645 log_gauss[neighbor] = mahal + norm_term; 00646 } 00647 } 00648 else 00649 { 00650 // Fetching nearest neighbors for density estimation. 00651 log_gauss.resize(L); 00652 for(int t=0; t<L;t++) 00653 { 00654 reference_set->getRow(t,neighbor_row); 00655 if(!store_prediction) 00656 { 00657 predictor->fprop(neighbor_row, F.toVec() & mu_temp & sn_temp); 00658 00659 // N.B. this is the SVD of F' 00660 lapackSVD(F, Ut_svd, S_svd, V_svd,'A',1.5); 00661 for (int k=0;k<ncomponents;k++) 00662 { 00663 sm_temp[k] = mypow(S_svd[k],2); 00664 U_temp(k) << Ut_svd(k); 00665 } 00666 } 00667 else 00668 { 00669 if(learn_mu) 00670 mu_temp << mus(t); 00671 sn_temp[0] = sns[t]; 00672 sm_temp << sms(t); 00673 U_temp << Fs[t]; 00674 } 00675 00676 if(learn_mu) 00677 { 00678 substract(t_row,neighbor_row,x_minus_neighbor); 00679 substract(x_minus_neighbor,mu_temp,z); 00680 } 00681 else 00682 substract(t_row,neighbor_row,z); 00683 00684 mahal = -0.5*pownorm(z)/sn_temp[0]; 00685 norm_term = - inputsize_/2.0 * Log2Pi - log_L 00686 - 0.5*(inputsize_-ncomponents)*pl_log(sn_temp[0]); 00687 00688 for(int k=0; k<ncomponents; k++) 00689 { 00690 mahal -= square(dot(z,U_temp(k)))*(0.5/(sm_temp[k]+sn_temp[0]) 00691 - 0.5/sn_temp[0]); 00692 norm_term -= 0.5*pl_log(sm_temp[k]+sn_temp[0]); 00693 } 00694 00695 log_gauss[t] = mahal + norm_term; 00696 } 00697 } 00698 ret = logadd(log_gauss); 00699 00700 return ret; 00701 } 00702 00703 00705 // computeOutput // 00707 void NonLocalManifoldParzen::computeOutput(const Vec& input, Vec& output) const 00708 { 00709 switch(outputs_def[0]) 00710 { 00711 /* 00712 case 'r': 00713 { 00714 string fsave = ""; 00715 VMat temp; 00716 real step_size = rw_size_step; 00717 real dp; 00718 t_row << input; 00719 Vec last_F(inputsize()); 00720 for(int s=0; s<rw_n_step;s++) 00721 { 00722 if(s == 0) 00723 { 00724 predictor->fprop(t_row, F.toVec() & mu_temp & sn_temp); 00725 last_F << F(rw_ith_component); 00726 } 00727 predictor->fprop(t_row, F.toVec() & mu_temp & sn_temp); 00728 00729 // N.B. this is the SVD of F' 00730 lapackSVD(F, Ut_svd, S_svd, V_svd,'A',1.5); 00731 F(rw_ith_component) << Ut_svd(rw_ith_component); 00732 00733 if(s % rw_save_every == 0) 00734 { 00735 fsave = rw_file_name + tostring(s) + ".amat"; 00736 temp = new MemoryVMatrix(t_row.toMat(1,t_row.length())); 00737 temp->saveAMAT(fsave,false,true); 00738 //PLearn::save(fsave,t_row); 00739 } 00740 dp = dot(last_F,F(rw_ith_component)); 00741 if(dp>0) dp = 1; 00742 else dp = -1; 00743 t_row += step_size*F(rw_ith_component)*abs(S_svd[rw_ith_component])*dp; 00744 last_F << dp*F(rw_ith_component); 00745 } 00746 output << t_row; 00747 00748 t_row << input; 00749 for(int s=0; s<rw_n_step;s++) 00750 { 00751 if(s == 0) 00752 { 00753 predictor->fprop(t_row, F.toVec() & mu_temp & sn_temp); 00754 last_F << (-1.0)*F(rw_ith_component); 00755 } 00756 00757 00758 predictor->fprop(t_row, F.toVec() & mu_temp & sn_temp); 00759 00760 // N.B. this is the SVD of F' 00761 lapackSVD(F, Ut_svd, S_svd, V_svd,'A',1.5); 00762 F(rw_ith_component) << Ut_svd(rw_ith_component); 00763 00764 if(s % rw_save_every == 0) 00765 { 00766 fsave = rw_file_name + tostring(-s) + ".amat"; 00767 temp = new MemoryVMatrix(t_row.toMat(1,t_row.length())); 00768 temp->saveAMAT(fsave,false,true); 00769 //PLearn::save(fsave,t_row); 00770 } 00771 dp = dot(last_F,F(rw_ith_component)); 00772 if(dp>0) dp = 1; 00773 else dp = -1; 00774 t_row += step_size*F(rw_ith_component)*abs(S_svd[rw_ith_component])*dp; 00775 last_F << dp*F(rw_ith_component); 00776 } 00777 break; 00778 } 00779 case 't': 00780 { 00781 predictor->fprop(input, F.toVec() & mu_temp & sn_temp); 00782 output << F.toVec(); 00783 break; 00784 } 00785 */ 00786 default: 00787 00788 inherited::computeOutput(input,output); 00789 } 00790 } 00791 00793 // outputsize // 00795 int NonLocalManifoldParzen::outputsize() const 00796 { 00797 switch(outputs_def[0]) 00798 { 00799 /* 00800 case 'm': 00801 return ncomponents; 00802 break; 00803 case 'r': 00804 return n; 00805 case 't': 00806 return ncomponents*n; 00807 */ 00808 default: 00809 return inherited::outputsize(); 00810 } 00811 } 00812 00813 } // end of namespace PLearn 00814 00815 00816 /* 00817 Local Variables: 00818 mode:c++ 00819 c-basic-offset:4 00820 c-file-style:"stroustrup" 00821 c-file-offsets:((innamespace . 0)(inline-open . 0)) 00822 indent-tabs-mode:nil 00823 fill-column:79 00824 End: 00825 */ 00826 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=79 :