runtime/ComputeAverageRMS.cpp

   1 #include "TrainingDatabase.hpp"
   2 #include <string>
   3 #include <fcntl.h>
   4 #include "Utilities.hxx"
   5 #include <mysql++/mysql++.h>
   6
   7 class CompAvgException
   8   : public std::exception
   9 {
  10 public:
  11   CompAvgException(const std::string& aWhat)
  12     : mWhat(aWhat)
  13   {
  14   }
  15
  16   ~CompAvgException() throw()
  17   {
  18   }
  19
  20   const char* what() const throw()
  21   {
  22     return mWhat.c_str();
  23   }
  24
  25 private:
  26   std::string mWhat;
  27 };
  28
  29 class AverageComputer
  30 {
  31 public:
  32   AverageComputer(mysqlpp::Connection& aConn, const char* aFilename)
  33     : mConn(aConn), mAnswers(NULL), mSum(NULL), mTemp(NULL)
  34   {
  35     mFile = open(aFilename, O_RDONLY | O_LARGEFILE);
  36     if (mFile == -1)
  37       throw CompAvgException("Can't open training file.");
  38
  39     if (read(mFile, &mTdh, sizeof(TrainingDatabaseHeader)) != sizeof(TrainingDatabaseHeader))
  40       throw CompAvgException("Can't read the training file.");
  41
  42     mJump = mTdh.ninput_signals * sizeof(double);
  43
  44     mAnswers = new double[mTdh.noutput_signals];
  45     mSum = new double[mTdh.noutput_signals];
  46     mTemp = new double[mTdh.noutput_signals];
  47   }
  48
  49   ~AverageComputer()
  50   {
  51     if (mFile != -1)
  52       close(mFile);
  53     if (mAnswers)
  54       delete [] mAnswers;
  55     if (mTemp)
  56       delete [] mTemp;
  57     if (mSum)
  58       delete [] mSum;
  59   }
  60
  61   void
  62   computeRMSForAllVariants()
  63   {
  64     uint32_t v = mTdh.first_variant_id;
  65     uint32_t lv = v + mTdh.nvariants;
  66
  67     for (; v < lv; v++)
  68     {
  69       double RMS = computeRMSFromAverageExcluding(v);
  70       mysqlpp::Query q(mConn.query());
  71       q << "UPDATE crossval_rss SET average_rss='" << RMS
  72         << "' WHERE leftout_variant_id='" << v << "'";
  73       q.execute();
  74     }
  75   }
  76
  77   double
  78   computeRMSFromAverageExcluding(uint32_t aExclude)
  79   {
  80     if (aExclude < mTdh.first_variant_id)
  81       throw CompAvgException("Excluded variant out of range.");
  82
  83     aExclude -= mTdh.first_variant_id;
  84     if (aExclude >= mTdh.nvariants)
  85       throw CompAvgException("Excluded variant out of range.");
  86
  87     lseek64(mFile, sizeof(mTdh), SEEK_SET);
  88     memset(mSum, 0, sizeof(double) * mTdh.noutput_signals);
  89
  90     uint32_t i(0);
  91     for (; i < mTdh.nvariants; i++)
  92     {
  93       // Skip the inputs...
  94       if (mJump != 0)
  95         lseek64(mFile, mJump, SEEK_CUR);
  96
  97       double * into = (aExclude == i) ? mAnswers : mTemp;
  98
  99       if (read(mFile, into, sizeof(double) * mTdh.noutput_signals) !=
 100           sizeof(double) * mTdh.noutput_signals)
 101         throw CompAvgException("Couldn't load signals.");
 102
 103       if (into == mTemp)
 104       {
 105         uint32_t j;
 106         for (j = 0; j < mTdh.noutput_signals; j++)
 107           mSum[j] += mTemp[j];
 108       }
 109     }
 110
 111     double factor = 1.0 / (mTdh.nvariants - 1);
 112
 113     // RSS isn't robust against a few outliers, so use residual median of
 114     // squares instead.
 115     std::vector<double> residualsSq(mTdh.noutput_signals);
 116     for (i = 0; i < mTdh.noutput_signals; i++)
 117     {
 118       double t = mAnswers[i] - factor * mSum[i];
 119       residualsSq[i] = t * t;
 120     }
 121     std::sort(residualsSq.begin(), residualsSq.end());
 122
 123     double rms;
 124
 125     if ((mTdh.noutput_signals) & 1 == 1)
 126     {
 127       rms = residualsSq[mTdh.noutput_signals / 2];
 128     }
 129     else
 130     {
 131       uint32_t ofs = mTdh.noutput_signals / 2;
 132       rms = (residualsSq[ofs - 1] + residualsSq[ofs]) / 2;
 133     }
 134
 135     return rms;
 136   }
 137
 138 private:
 139   int mFile;
 140   TrainingDatabaseHeader mTdh;
 141   mysqlpp::Connection& mConn;
 142   double* mAnswers, * mSum, * mTemp;
 143   off64_t mJump;
 144 };
 145
 146 int
 147 main(int argc, char** argv)
 148 {
 149   if (argc < 6)
 150   {
 151     printf("Usage: ComputeAverageRMS db_host db_user db_password db_database trainingfile\n");
 152     // printf("Set variant=-1 if RMS is sum for all variants.\n");
 153     return 1;
 154   }
 155
 156   try
 157   {
 158     mysqlpp::Connection conn(mysqlpp::use_exceptions);
 159     conn.connect(argv[4], argv[1], argv[2], argv[3]);
 160
 161     AverageComputer ac(conn, argv[5]);
 162     ac.computeRMSForAllVariants();
 163   }
 164   catch (std::exception& e)
 165   {
 166     printf("Error: %s\n", e.what());
 167   }
 168 }