PLearn 0.1
Grapher.cc
Go to the documentation of this file.
00001 
00002 // -*- C++ -*-
00003 
00004 // Grapher.cc
00005 //
00006 // Copyright (C) 2003  Pascal Vincent 
00007 // 
00008 // Redistribution and use in source and binary forms, with or without
00009 // modification, are permitted provided that the following conditions are met:
00010 // 
00011 //  1. Redistributions of source code must retain the above copyright
00012 //     notice, this list of conditions and the following disclaimer.
00013 // 
00014 //  2. Redistributions in binary form must reproduce the above copyright
00015 //     notice, this list of conditions and the following disclaimer in the
00016 //     documentation and/or other materials provided with the distribution.
00017 // 
00018 //  3. The name of the authors may not be used to endorse or promote
00019 //     products derived from this software without specific prior written
00020 //     permission.
00021 // 
00022 // THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR
00023 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
00024 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN
00025 // NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
00026 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
00027 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00028 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00029 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00030 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00031 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00032 // 
00033 // This file is part of the PLearn library. For more information on the PLearn
00034 // library, go to the PLearn Web site at www.plearn.org
00035 
00036 /* *******************************************************      
00037  * $Id: Grapher.cc 8621 2008-03-03 19:43:36Z tihocan $ 
00038  ******************************************************* */
00039 
00041 #include "Grapher.h"
00042 #include <plearn/io/load_and_save.h>
00043 #include <plearn/math/VecStatsCollector.h>
00044 #include <plearn/vmat/VMat_basic_stats.h>
00045 #include <plearn/display/GhostScript.h>
00046 #include <plearn/display/Gnuplot.h>
00047 #include <plearn/vmat/RegularGridVMatrix.h>
00048 #include <plearn/base/stringutils.h>
00049 #include <plearn/math/TMat_maths.h>
00050 
00051 namespace PLearn {
00052 using namespace std;
00053 
00054 
00059 void DX_write_2D_fields(ostream& out, const string& basename, TVec<Mat> fields, real x0, real y0, real deltax, real deltay, 
00060                         TVec<string> fieldnames=TVec<string>())
00061 {
00062     int nfields = fields.length();
00063     int nx = fields[0].length();
00064     int ny = fields[0].width();
00065 
00066     string posname = string("\"") + basename + "_gridpos\"";
00067 
00068     out << "object " << posname << " class gridpositions counts " << nx << " " << ny << "\n"
00069         << "origin  " << x0 << " " << y0 << "\n"
00070         << "delta   " << deltax << " 0 \n"
00071         << "delta    0 " << deltay << " \n\n\n";
00072 
00073     string conname = string("\"") + basename + "_gridcon\"";
00074 
00075     out << "object " << conname << " class gridconnections counts " << nx << " " << ny << "\n"
00076         //      << "attribute \"element type\" string \"cubes\" \n"
00077         << "attribute \"ref\" string \"positions\" \n\n\n";
00078 
00079     for(int k=0; k<nfields; k++)
00080     {
00081         Mat& m = fields[k];
00082         string fieldname = tostring(k);
00083         if(fieldnames)
00084             fieldname = fieldnames[k];
00085 
00086         string dataname = string("\"") + basename + "_" + fieldname + "_data\"";
00087 
00088         out << "object " << dataname << " class array type float rank 0 items " << nx*ny << " data follows \n";
00089         for(int i=0; i<nx; i++)
00090         {
00091             for(int j=0; j<ny; j++)
00092                 out << m(i,j) << " ";
00093             out << "\n";
00094         }
00095         out << "attribute \"dep\" string \"positions\" \n\n\n";
00096 
00097         out << "object \"" << fieldname << "\" class field \n"
00098             << "component \"positions\" " << posname << " \n"
00099             << "component \"connections\" " << conname << " \n"
00100             << "component \"data\" " << dataname << " \n\n\n";
00101     }
00102 }
00103 
00104 
00105 void DX_write_2D_fields(ostream& out, const string& basename, Vec X, Vec Y, TVec<Mat> fields)
00106 {
00107     int nfields = fields.length();
00108     int nx = fields[0].length();
00109     int ny = fields[0].width();
00110 
00111     /*
00112       out << "object \"" << basename << "_X\" class array type float rank 0 items " << nx << " data follows \n";
00113       for(int i=0; i<nx; i++)
00114       out << X[i] << "\n";
00115       out << "\n\n";
00116     
00117       out << "object \"" << basename << "_Y\" class array type float rank 0 items " << ny << " data follows \n";
00118       for(int i=0; i<ny; i++)
00119       out << Y[i] << "\n";
00120     */
00121 
00122     string posname = string("\"") + basename + "_gridpos\"";
00123     out << "object " << posname << " class array type float rank 1 shape 2 items " << nx*ny << " data follows\n";
00124     for(int i=0; i<nx; i++)
00125         for(int j=0; j<ny; j++)
00126             out << X[i] << " " << Y[j] << "\n";
00127     out << "\n\n";
00128 
00129     string conname = string("\"") + basename + "_gridcon\"";
00130     out << "object " << conname << " class gridconnections counts " << nx << " " << ny << "\n"
00131         //      << "attribute \"element type\" string \"cubes\" \n"
00132         << "attribute \"ref\" string \"positions\" \n\n\n";
00133 
00134     for(int k=0; k<nfields; k++)
00135     {
00136         Mat& m = fields[k];
00137         string fieldname = "output" + tostring(k);
00138         string dataname = string("\"") + basename + "_" + fieldname + "_data\"";
00139 
00140         out << "object " << dataname << " class array type float rank 0 items " << nx*ny << " data follows \n";
00141         for(int i=0; i<nx; i++)
00142         {
00143             for(int j=0; j<ny; j++)
00144                 out << m(i,j) << " ";
00145             out << "\n";
00146         }
00147         out << "attribute \"dep\" string \"positions\" \n\n\n";
00148 
00149         out << "object \"" << fieldname << "\" class field \n"
00150             << "component \"positions\" " << posname << " \n"
00151             << "component \"connections\" " << conname << " \n"
00152             << "component \"data\" " << dataname << " \n\n\n";
00153     }
00154 }
00155 
00156 
00157 TVec<Mat> computeOutputFields(PP<PLearner> learner, Vec X, Vec Y)
00158 {
00159     int noutputs = learner->outputsize();
00160 
00161     int nx = X.length();
00162     int ny = Y.length();
00163     int nfields = noutputs;
00164     TVec<Mat> fields(nfields);
00165 
00166     for(int k=0; k<nfields; k++)
00167         fields[k].resize(nx,ny);
00168 
00169     Vec input(2);
00170     Vec output(noutputs);
00171 
00172     ProgressBar pb("Computing " + tostring(nx) + " x " + tostring(ny) + " output field",nx*ny);
00173   
00174     for(int i=0; i<nx; i++)
00175         for(int j=0; j<ny; j++)
00176         {
00177             input[0] = X[i];
00178             input[1] = Y[j];
00179             learner->computeOutput(input,output);
00180             // cerr << "in: " << input << " out: " << output << endl;
00181             for(int k=0; k<noutputs; k++)
00182                 fields[k](i,j) = output[k];
00183             pb.update(i*nx+j);
00184         }
00185 
00186     return fields;
00187 }
00188 
00189 
00190 TVec<Mat> computeOutputFields(PP<PLearner> learner, int nx, int ny, real x0, real y0, real deltax, real deltay)
00191 {
00192     int noutputs = learner->outputsize();
00193     int nfields = noutputs;
00194 
00195     TVec<Mat> fields(nfields);
00196     for(int k=0; k<nfields; k++)
00197         fields[k].resize(nx,ny);
00198 
00199     Vec input(2);
00200     Vec output(noutputs);
00201 
00202     ProgressBar pb("Computing " + tostring(nx) + " x " + tostring(ny) + " output field",nx*ny);
00203 
00204     real x = x0;
00205     real y = y0;
00206     for(int i=0; i<nx; i++, x+=deltax)
00207         for(int j=0; j<ny; j++, y+=deltay)
00208         {
00209             input[0] = x;
00210             input[1] = y;
00211             learner->computeOutput(input,output);
00212             // cerr << "in: " << input << " out: " << output << endl;
00213             for(int k=0; k<noutputs; k++)
00214                 fields[k](i,j) = output[k];
00215             pb.update(i*nx+j);
00216         }
00217 
00218     return fields;
00219 }
00220 
00221 // Finds appropriate x0, y0, deltax, deltay from the dataset range, computes the fields and returns them
00222 // extraspace of .10 means we'll look 10% beyond the data range on every side
00223 TVec<Mat> computeOutputFieldsAutoRange(PP<PLearner> learner, VMat dataset, int nx, int ny, 
00224                                        real& x0, real& y0, real& deltax, real& deltay, real extraspace=.10)
00225 {
00226     Vec minv(2);
00227     Vec maxv(2);
00228     computeRange(dataset.subMatColumns(0,2), minv, maxv);
00229     real extrax = (maxv[0]-minv[0])*extraspace;
00230     x0 = minv[0]-extrax;
00231     deltax = (maxv[0]+extrax-x0)/nx;
00232     real extray = (maxv[1]-minv[1])*extraspace;
00233     y0 = minv[1]-extray;
00234     deltay = (maxv[1]+extray-y0)/ny;
00235     return computeOutputFields(learner, nx, ny, x0, y0, deltax, deltay);
00236 }
00237 
00238 
00239 void computeXYPositions(VMat dataset, int nx, int ny, Vec& X, Vec& Y, real extraspace=.10)
00240 {
00241     Vec minv(2);
00242     Vec maxv(2);
00243     computeRange(dataset.subMatColumns(0,2), minv, maxv);
00244     real extrax = (maxv[0]-minv[0])*extraspace;
00245     real x0 = minv[0]-extrax;
00246     real deltax = (maxv[0]+extrax-x0)/nx;
00247     real extray = (maxv[1]-minv[1])*extraspace;
00248     real y0 = minv[1]-extray;
00249     real deltay = (maxv[1]+extray-y0)/ny;
00250 
00251     set<real> xpos;
00252     set<real> ypos;
00253     int l = dataset.length();
00254     Vec datapoint(2);
00255     for(int i=0; i<l; i++)
00256     {
00257         dataset->getRow(i,datapoint);
00258         xpos.insert(datapoint[0]);
00259         ypos.insert(datapoint[1]);
00260     }
00261     real x = x0;
00262     for(int i=0; i<nx; i++, x+=deltax)
00263         xpos.insert(x);
00264     real y = y0;
00265     for(int j=0; j<ny; j++, y+=deltay)
00266         ypos.insert(y);
00267     set<real>::iterator it;
00268     X.resize((int)xpos.size());
00269     real* xptr = X.data();
00270     it = xpos.begin();
00271     while(it!=xpos.end())
00272         *xptr++ = *it++;
00273     Y.resize((int)ypos.size());
00274     real* yptr = Y.data();
00275     it = ypos.begin();
00276     while(it!=ypos.end())
00277         *yptr++ = *it++;
00278 }
00279 
00280 
00281 
00284 void DX_create_dataset_outputs_file(const string& filename, PP<PLearner> learner, VMat dataset)
00285 {
00286     ofstream out(filename.c_str());
00287 
00288     int l = dataset.length();
00289     int inputsize = learner->inputsize();
00290     int targetsize = learner->targetsize();
00291     int outputsize = learner->outputsize();
00292 
00293     // First write data points (input -> target, output)
00294     Vec input(inputsize);
00295     Vec target(targetsize);
00296     real weight;
00297     Vec output(outputsize);
00298 
00299     // write 2D positions
00300     out << "object \"dset_pos\" class array type float rank 1 shape " << inputsize << " items " << l << " data follows \n";
00301     for(int i=0; i<l; i++)
00302     {
00303         dataset->getExample(i,input,target,weight);
00304         for(int j=0; j<inputsize; j++)
00305             out << input[j] << " ";
00306         out << "\n";
00307     }
00308     out << "\n\n\n";
00309 
00310     // Now write data for those positions (target and output)
00311     if(targetsize+outputsize>0)
00312     {
00313         ProgressBar pb("Computing outputs for dataset points",l);
00314         out << "object \"dset_value\" class array type float rank 1 shape " << targetsize+outputsize << " items " << l << " data follows \n";
00315         for(int i=0; i<l; i++)
00316         {
00317             dataset->getExample(i,input,target,weight);
00318             for(int j=0; j<targetsize; j++)
00319                 out << target[j] << " ";
00320             learner->computeOutput(input, output);
00321             for(int j=0; j<outputsize; j++)
00322                 out << output[j] << " ";
00323             out << "\n";
00324             pb.update(i);
00325         }
00326         out << "attribute \"dep\" string \"positions\" \n\n\n";
00327     }
00328 
00329     // Field is created with two components: "positions" and "data"
00330     out << "object \"dset\" class field \n"
00331         << "component \"positions\" \"dset_pos\" \n";
00332     if(targetsize+outputsize>0)
00333         out << "component \"data\" \"dset_value\" \n";
00334     out << "\n\n\n";
00335 
00336 
00337   
00338     out << "end" << endl;
00339 }
00340 
00341 
00348 
00349 void DX_create_grid_outputs_file(const string& filename, PP<PLearner> learner, VMat dataset, 
00350                                  int nx, int ny, bool include_datapoint_grid=false, 
00351                                  real xmin=MISSING_VALUE, real xmax=MISSING_VALUE, 
00352                                  real ymin=MISSING_VALUE, real ymax=MISSING_VALUE,
00353                                  real extraspace=.10)
00354 {
00355     ofstream out(filename.c_str());
00356 
00357     double logsum = -FLT_MAX;
00358 
00359     int l = dataset.length();
00360     int inputsize = learner->inputsize();
00361     int targetsize = learner->targetsize();
00362     int outputsize = learner->outputsize();
00363 
00364     Vec input(inputsize);
00365     Vec target(targetsize);
00366     real weight;
00367     Vec output(outputsize);
00368 
00369     // Create the grid field
00370 
00371     set<real> xpos;
00372     set<real> ypos;
00373 
00374     // First the regular grid coordinates
00375     Vec minv(2);
00376     Vec maxv(2);
00377     computeRange(dataset.subMatColumns(0,2), minv, maxv);
00378     real extrax = (maxv[0]-minv[0])*extraspace;
00379     real extray = (maxv[1]-minv[1])*extraspace;
00380     if(is_missing(xmin))
00381         xmin = minv[0]-extrax;
00382     if(is_missing(xmax))
00383         xmax = maxv[0]+extrax;
00384     if(is_missing(ymin))
00385         ymin = minv[1]-extray;
00386     if(is_missing(ymax))
00387         ymax = maxv[1]+extray;
00388     real deltax = (xmax-xmin)/nx;
00389     real deltay = (ymax-ymin)/ny;
00390 
00391     real x = xmin;
00392     for(int i=0; i<nx; i++, x+=deltax)
00393         xpos.insert(x);
00394     real y = ymin;
00395     for(int j=0; j<ny; j++, y+=deltay)
00396         ypos.insert(y);
00397 
00398     // also include irregular grid coordinates based on coordinates of dataset points?
00399     if(include_datapoint_grid) 
00400     {
00401         for(int i=0; i<l; i++)
00402         {
00403             dataset->getExample(i,input,target,weight);
00404             x = input[0];
00405             y = input[1];
00406             if(x>xmin && x<xmax)
00407                 xpos.insert(x);
00408             if(y>ymin && y<ymax)
00409                 ypos.insert(y);
00410         }
00411     }
00412 
00413     nx = (int)xpos.size();
00414     ny = (int)ypos.size();
00415     set<real>::iterator itx;
00416     set<real>::iterator ity;
00417 
00418     out << "object \"outputs_gridpos\" class array type float rank 1 shape 2 items " << nx*ny << " data follows\n";
00419     for(itx=xpos.begin(); itx!=xpos.end(); ++itx)
00420         for(ity=ypos.begin(); ity!=ypos.end(); ++ity)
00421             out << *itx << " " << *ity << "\n";
00422     out << "\n\n";
00423 
00424     out << "object \"outputs_gridcon\" class gridconnections counts " << nx << " " << ny << "\n"
00425         //      << "attribute \"element type\" string \"cubes\" \n"
00426         << "attribute \"ref\" string \"positions\" \n\n\n";
00427 
00428     out << "object \"outputs_values\" class array type float rank 1 shape " << outputsize << " items " << nx*ny << " data follows \n";
00429   
00430     ProgressBar pb("Computing outputs for grid positions: " + tostring(nx)+"x"+tostring(ny), nx*ny);
00431     int n = 0;
00432     for(itx=xpos.begin(); itx!=xpos.end(); ++itx)
00433     {
00434         input[0] = *itx;
00435         for(ity=ypos.begin(); ity!=ypos.end(); ++ity)
00436         {
00437             input[1] = *ity;
00438             learner->computeOutput(input, output);
00439             for(int j=0; j<outputsize; j++)
00440                 out << output[j] << " ";
00441             out << "\n";
00442             if(is_equal(logsum, -FLT_MAX))
00443                 logsum = output[0];
00444             else 
00445                 logsum = logadd(logsum, double(output[0]));
00446             pb.update(n++);
00447         }
00448     }
00449     pb.close();
00450     out << "attribute \"dep\" string \"positions\" \n\n\n";
00451 
00452     out << "object \"outputs\" class field \n"
00453         << "component \"positions\" \"outputs_gridpos\" \n"
00454         << "component \"connections\" \"outputs_gridcon\" \n"
00455         << "component \"data\" \"outputs_values\" \n\n\n";
00456   
00457     out << "end" << endl;
00458 
00459     double surfelem = deltax*deltay;
00460     double surfintegral = exp(logsum)*surfelem;
00461     cerr << "Estimated integral over sampled domain: " << surfintegral << endl;
00462 }
00463 
00464 void Grapher::computeAutoGridrange()
00465 {
00466     int d = trainset->inputsize();
00467     gridrange.resize(d);
00468     for(int j=0; j<d; j++)
00469         gridrange[j] = pair<real,real>(FLT_MAX,-FLT_MAX);
00470     Vec input;
00471     Vec target;
00472     real weight;
00473     int l = trainset.length();
00474     for(int i=0; i<l; i++)
00475     {
00476         trainset->getExample(i, input, target, weight);
00477         for(int j=0; j<d; j++)
00478         {
00479             real x_j = input[j];
00480             if(x_j<gridrange[j].first)
00481                 gridrange[j].first = x_j;
00482             if(x_j>gridrange[j].second)
00483                 gridrange[j].second = x_j;
00484         }      
00485     }
00486 
00487     // Now add extra 10%
00488     real extra = .10;
00489     for(int j=0; j<d; j++)
00490     {
00491         real extent = extra*(gridrange[j].second-gridrange[j].first);
00492         gridrange[j].first -= extent;
00493         gridrange[j].second += extent;
00494     }
00495 }
00496 
00497 
00498 Grapher::Grapher() 
00499     :basename("dxplot"), task(""), class1_threshold(0.5),
00500      radius(-0.01), bw(false)
00501 {
00502 }
00503 
00504 PLEARN_IMPLEMENT_OBJECT(Grapher, "ONE LINE DESCR", "NO HELP");
00505 
00506 void Grapher::declareOptions(OptionList& ol)
00507 {
00508     // ### Declare all of this object's options here
00509     // ### For the "flags" of each option, you should typically specify  
00510     // ### one of OptionBase::buildoption, OptionBase::learntoption or 
00511     // ### OptionBase::tuningoption. Another possible flag to be combined with
00512     // ### is OptionBase::nosave
00513 
00514     declareOption(ol, "basename", &Grapher::basename, OptionBase::buildoption,
00515                   "Base name of the .dx data file to generate. Running this class will generate\n"
00516                   "files basename_dset.dx containing targets and outputs for the given dataset positions\n"
00517                   "and basename_outputs.dx containing outputs computed at grid positions\n");
00518     declareOption(ol, "task", &Grapher::task, OptionBase::buildoption,
00519                   "Desired plotting task. Can be \"1D regression\",\n"
00520                   "\"2D clustering\", \"2D density\", \"2D classification\",\n"
00521                   "\"2D regression\"");
00522     declareOption(ol, "class1_threshold", &Grapher::class1_threshold, OptionBase::buildoption,
00523                   "In the case of 1 output 2D classification, the output threshold to\n"
00524                   "have class=1 (below the threshold, the class = 0).  The default\n"
00525                   "is 0.5, which appropriate for sigmoidal-output learners, but use 0\n"
00526                   "if the learner outputs -1/1, such as for SVM's");
00527     declareOption(ol, "learner", &Grapher::learner, OptionBase::buildoption,
00528                   "The learner to train/test");
00529     declareOption(ol, "trainset", &Grapher::trainset, OptionBase::buildoption,
00530                   "The training set to train the learner on\n");
00531     declareOption(ol, "gridrange", &Grapher::gridrange, OptionBase::buildoption,
00532                   "A vector of low:high pairs with as many dimensions as the input space\n"
00533                   "ex for 2D: [ -10:10 -3:4 ] \n"
00534                   "If empty, it will be automatically inferred from the range of the\n"
00535                   "trainset inputs (with an extra 10%)");
00536     declareOption(ol, "griddim", &Grapher::griddim, OptionBase::buildoption,
00537                   "A vector of integers giving the number of sample coordinates\n"
00538                   "for each dimension of the grid. Ex for 2D: [ 100 100 ]\n");
00539     declareOption(ol, "radius", &Grapher::radius, OptionBase::buildoption,
00540                   "The radius of the discs around data points.\n"
00541                   "(If negative, it's considered to be expressed as a percentage of the x range)\n");
00542     declareOption(ol, "bw", &Grapher::bw, OptionBase::buildoption,
00543                   "Set this to true if you want to generate black and white eps");
00544     declareOption(ol, "save_learner_as", &Grapher::save_learner_as, OptionBase::buildoption,
00545                   "(Optionally) save trained learner in this file (.psave)");
00546 
00547     // Now call the parent class' declareOptions
00548     inherited::declareOptions(ol);
00549 }
00550 
00551 void Grapher::build_()
00552 {
00553     // ### This method should do the real building of the object,
00554     // ### according to set 'options', in *any* situation. 
00555     // ### Typical situations include:
00556     // ###  - Initial building of an object from a few user-specified options
00557     // ###  - Building of a "reloaded" object: i.e. from the complete set of all serialised options.
00558     // ###  - Updating or "re-building" of an object after a few "tuning" options have been modified.
00559     // ### You should assume that the parent class' build_() has already been called.
00560 }
00561 
00562 real color(int colornum, real lightness)
00563 {
00564     real col = 0;
00565     switch(colornum)
00566     {
00567     case 0:
00568         col = rgb2real(lightness,1,1);
00569         break;
00570     case 1:
00571         col = rgb2real(1,lightness,1);
00572         break;
00573     case 2:
00574         col = rgb2real(1,1,lightness);
00575         break;
00576     case 3:
00577         col = rgb2real(1,lightness,lightness);
00578         break;
00579     case 4:
00580         col = rgb2real(lightness,1,lightness);
00581     case 5:
00582         col = rgb2real(lightness,lightness,1);
00583     default:
00584         PLERROR("No color setting for colornum %d", colornum);
00585     }
00586     return col;
00587 }
00588 
00589 void Grapher::plot_1D_regression(string basename, VMat trainset, 
00590                                  TVec<int> griddim, TVec< pair<real,real> > gridrange, 
00591                                  VMat gridoutputs, VMat trainoutputs, bool bw)
00592 {
00593     if(griddim.size()!=1)
00594         PLERROR("In Grapher::plot_1D_regression, not a 1D grid!");  
00595   
00596     int nx = griddim[0];
00597     real x = gridrange[0].first;
00598     real w = gridrange[0].second - x;
00599     real dx = w/(nx-1);
00600     int l = trainset.length();
00601     Mat curve(nx+l,2);
00602     Mat points(l,2);
00603     for(int i=0; i<nx; i++)
00604     {
00605         curve(i,0) = x;
00606         curve(i,1) = gridoutputs(i,0);
00607         x += dx;
00608     }
00609 
00610     Vec input(1);
00611     Vec target(1);
00612     Vec output(1);
00613     real weight;
00614     for(int i=0; i<l; i++)
00615     {
00616         trainset->getExample(i, input, target, weight);
00617         trainoutputs->getRow(i, output);
00618         points(i,0) = input[0];
00619         points(i,1) = target[0];
00620         curve(nx+i,0) = input[0];
00621         curve(nx+i,1) = output[0];
00622     }
00623 
00624     sortRows(curve);
00625 
00626     saveAscii(basename+"_points.amat", points);
00627     saveAscii(basename+"_curve.amat", curve);
00628     Gnuplot gp;
00629     gp << "plot '" << basename+"_points.amat" << "' with points, '" 
00630        << basename + "_curve.amat" << "' with lines" << endl;
00631 
00632     pgetline(cin);
00633 }
00634 
00635 void Grapher::plot_2D_classification(string epsfname, VMat trainset, 
00636                                      TVec<int> griddim, TVec< pair<real,real> > gridrange,
00637                                      VMat gridoutputs, real radius, bool bw)
00638 {
00639     cerr << "Plotting 2D classification result" << endl;
00640     if(griddim.size()!=2 || gridrange.size()!=2)
00641         PLERROR("In Grapher::plot_2D_classification griddim and gridrange must be of size 2");
00642     int nx = griddim[0];
00643     int ny = griddim[1];
00644     int nclasses = gridoutputs.width();
00645     real x = gridrange[0].first;
00646     real y = gridrange[1].first;
00647     real w = gridrange[0].second - x;
00648     real h = gridrange[1].second - y;
00649     if(nclasses<2)
00650         PLERROR("In Grapher::plot_2D_classification number of classes (width of gridoutputs) must be at least 2: it is currently %d",nclasses);
00651 
00652     real gswidth = 600;
00653     real gsheight = gswidth/w*h;
00654     GhostScript gs(epsfname, 0, 0, gswidth, gsheight);
00655     gs.mapping(x, y, w, h, 0, 0, gswidth, gsheight);
00656 
00657 
00658     Mat image(ny,nx);
00659     Vec output(nclasses);
00660     for(int i=0; i<ny; i++)
00661         for(int j=0; j<nx; j++)
00662         {
00663             gridoutputs->getRow((ny-i-1)+ny*j,output);
00664             int winner = argmax(output);
00665             // cout << i << " " << j << " " << output << endl;
00666             // real winnerval = output[winner];
00667             if(bw) // grayscale
00668                 image(i,j) = (winner==0 ?0.7 :0.3);
00669             else // color
00670                 image(i,j) = color(winner,0.7);
00671         }
00672 
00673     // cout << "IMAGE:" << endl << image << endl;
00674   
00675     if(bw)
00676         gs.displayGray(image, x, y, w, h);
00677     else
00678         gs.displayRGB(image, x, y, w, h);
00679 
00680     int l = trainset.length();
00681     Vec input;
00682     Vec target;
00683     real weight;
00684     gs.setgray(0); // black
00685     for(int i=0; i<l; i++)
00686     {
00687         trainset->getExample(i,input,target,weight);
00688         if(target.length()==1)
00689         {
00690             if(bw)
00691                 gs.setgray(fast_exact_is_equal(target[0], 0) ?0 :1);
00692             else
00693                 gs.setcolor(color(int(target[0]),0.2));
00694         }
00695         else
00696         {
00697             if(bw)
00698                 gs.setgray(fast_exact_is_equal(target[0], 0) ?0 :1);
00699             else
00700                 gs.setcolor(color(argmax(target),0.2));
00701         }
00702         gs.fillCircle(input[0],input[1],radius);
00703     }
00704 
00705 }
00706 
00707 /*
00708   void Grapher::plot_2D_density(Mat gridoutputs) const
00709   {
00710   
00711   }
00712 
00713   void Grapher::plot_1D_regression(Mat gridoutputs, Mat trainoutputs) const
00714   {
00715   
00716   }
00717 */
00718 
00720 void Grapher::run()
00721 {
00722     int l = trainset->length();
00723     PP<VecStatsCollector> statscol = new VecStatsCollector();
00724     learner->setTrainStatsCollector(statscol);
00725     learner->setTrainingSet(trainset);
00726 
00727     cerr << "*** Training learner on trainset of length " << l << " ... ***" << endl; 
00728     learner->train();
00729     cerr << "Final traincosts: " << statscol->getMean() << endl;
00730 
00731     if(save_learner_as!="")
00732     {
00733         cerr << "Saving trained learner in file " << save_learner_as << endl;
00734         PLearn::save(save_learner_as, *learner);
00735     }
00736 
00737 
00738     cerr << "*** Computing outputs on trainset inputs... ***" << endl; 
00739     Mat trainoutputs(l, learner->outputsize());
00740     learner->use(trainset, trainoutputs);
00741 
00742     // Try to determine the task if not specified
00743     if (task == "") {
00744         if(trainset->inputsize()==1 &&
00745            trainset->targetsize()==1 && learner->outputsize()==1)
00746             task = "1D regression";
00747         else if(trainset->inputsize()==2)
00748         {
00749             switch(trainset->targetsize())
00750             {
00751             case 0:  // density estimation or clustering
00752                 if(learner->outputsize()>1)
00753                     task = "2D clustering";
00754                 else
00755                     task = "2D density";
00756                 break;
00757             case 1: // classif or regression
00758                 if(learner->outputsize()>1)
00759                     task = "2D classification";
00760                 else
00761                     task = "2D regression";
00762                 break;
00763             default:
00764                 PLERROR("Tasks with targetsize > 1 (multi-regression) not supported");
00765             }
00766         }
00767         else
00768             PLERROR("Task wih inputsize=%d, targetsize=%d, outputsize=%d not supported", 
00769                     trainset->inputsize(), trainset->targetsize(), learner->outputsize());
00770     }
00771     
00772     // Now compute outputs on grid
00773     if(gridrange.isEmpty())
00774         computeAutoGridrange();
00775     cerr << "*** Computing outputs on " << griddim << " grid... ***" << endl; 
00776     VMat gridinputs = new RegularGridVMatrix(griddim, gridrange);
00777     Mat gridoutputs(gridinputs->length(),learner->outputsize());
00778     learner->use(gridinputs, gridoutputs);
00779 
00780     if (task == "2D classification" && learner->outputsize() == 1) {
00781         // Transform a one-class output into a two-class one...
00782         Mat newgridoutputs(gridinputs->length(), 2);
00783         for (int i=0; i<gridinputs->length(); ++i) {
00784             newgridoutputs(i,0) = gridoutputs(i,0) <= class1_threshold;
00785             newgridoutputs(i,1) = gridoutputs(i,0) > class1_threshold;
00786         }
00787         gridoutputs = newgridoutputs;
00788     }
00789 
00790     cerr << ">>> TASK PERFORMED: " << task << endl;
00791 
00792     string epsfname = basename+".eps";
00793     cerr << "Creating file: " << epsfname << endl;
00794 
00795     if(radius<0)
00796         radius = fabs(radius * (gridrange[0].second-gridrange[0].first));
00797 
00798     if(task=="2D classification" || task=="2D clustering")
00799         plot_2D_classification(epsfname, trainset, griddim, gridrange, gridoutputs, radius, bw);
00800     else if(task=="1D regression")
00801         plot_1D_regression(basename, trainset, griddim, gridrange, gridoutputs, trainoutputs, bw);
00802     //  else if(task=="2D regression")
00803     //  plot_2D-regression(basename, 
00804     /*
00805       else if(task=="2D density")
00806       plot_2D_density(gridoutputs);
00807     */
00808 
00809     // Old DX stuff
00810     /*
00811       string dset_fname = basename+"_dset.dx";
00812       cerr << "Computing and writing trainset output field to file " << dset_fname << endl;
00813       DX_create_dataset_outputs_file(dset_fname, learner, trainset);
00814 
00815       string outputs_fname = basename+"_outputs.dx";
00816       cerr << "Computing and writing grid output field to file " << outputs_fname << endl; 
00817       DX_create_grid_outputs_file(outputs_fname, learner, trainset, nx, ny, 
00818       include_datapoint_grid, 
00819       xmin, xmax, ymin, ymax);
00820       cerr << "You can now view those files with OpenDX." << endl;
00821     */
00822 }
00823 
00824 
00825 // ### Nothing to add here, simply calls build_
00826 void Grapher::build()
00827 {
00828     inherited::build();
00829     build_();
00830 }
00831 
00832 
00833 void Grapher::makeDeepCopyFromShallowCopy(CopiesMap& copies)
00834 {
00835     inherited::makeDeepCopyFromShallowCopy(copies);
00836 }
00837 
00838 } // end of namespace PLearn
00839 
00840 
00841 /*
00842   Local Variables:
00843   mode:c++
00844   c-basic-offset:4
00845   c-file-style:"stroustrup"
00846   c-file-offsets:((innamespace . 0)(inline-open . 0))
00847   indent-tabs-mode:nil
00848   fill-column:79
00849   End:
00850 */
00851 // 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