test/distance.cc

Code
Comments
Other
Rev Date Author Line
889 25 Sep 07 markus 1 // $Id$
889 25 Sep 07 markus 2
999 23 Dec 07 jari 3 /*
2119 12 Dec 09 peter 4   Copyright (C) 2007 Jari Häkkinen, Markus Ringnér
2119 12 Dec 09 peter 5   Copyright (C) 2008 Jari Häkkinen, Peter Johansson, Markus Ringnér
3550 03 Jan 17 peter 6   Copyright (C) 2009, 2010, 2012, 2016 Peter Johansson
999 23 Dec 07 jari 7
1437 25 Aug 08 peter 8   This file is part of the yat library, http://dev.thep.lu.se/yat
999 23 Dec 07 jari 9
999 23 Dec 07 jari 10   The yat library is free software; you can redistribute it and/or
999 23 Dec 07 jari 11   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 12   published by the Free Software Foundation; either version 3 of the
999 23 Dec 07 jari 13   License, or (at your option) any later version.
999 23 Dec 07 jari 14
999 23 Dec 07 jari 15   The yat library is distributed in the hope that it will be useful,
999 23 Dec 07 jari 16   but WITHOUT ANY WARRANTY; without even the implied warranty of
999 23 Dec 07 jari 17   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
999 23 Dec 07 jari 18   General Public License for more details.
999 23 Dec 07 jari 19
999 23 Dec 07 jari 20   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 21   along with yat. If not, see <http://www.gnu.org/licenses/>.
999 23 Dec 07 jari 22 */
999 23 Dec 07 jari 23
2881 18 Nov 12 peter 24 #include <config.h>
2881 18 Nov 12 peter 25
1233 15 Mar 08 peter 26 #include "Suite.h"
1233 15 Mar 08 peter 27
890 25 Sep 07 markus 28 #include "yat/classifier/DataLookupWeighted1D.h"
890 25 Sep 07 markus 29 #include "yat/classifier/MatrixLookupWeighted.h"
1052 07 Feb 08 peter 30 #include "yat/statistics/EuclideanDistance.h"
1052 07 Feb 08 peter 31 #include "yat/statistics/PearsonDistance.h"
2334 15 Oct 10 peter 32 #include "yat/utility/concept_check.h"
1582 15 Oct 08 peter 33 #include "yat/utility/DataIterator.h"
1121 22 Feb 08 peter 34 #include "yat/utility/Matrix.h"
1483 09 Sep 08 peter 35 #include "yat/utility/MatrixWeighted.h"
1120 21 Feb 08 peter 36 #include "yat/utility/Vector.h"
1582 15 Oct 08 peter 37 #include "yat/utility/WeightIterator.h"
889 25 Sep 07 markus 38
2202 21 Feb 10 peter 39 #include <boost/concept_archetype.hpp>
2334 15 Oct 10 peter 40 #include <boost/concept_check.hpp>
2202 21 Feb 10 peter 41
889 25 Sep 07 markus 42 #include <cassert>
889 25 Sep 07 markus 43 #include <fstream>
889 25 Sep 07 markus 44 #include <iostream>
1247 17 Mar 08 peter 45 #include <limits>
900 27 Sep 07 markus 46 #include <list>
900 27 Sep 07 markus 47 #include <vector>
889 25 Sep 07 markus 48
889 25 Sep 07 markus 49
889 25 Sep 07 markus 50 using namespace theplu::yat;
889 25 Sep 07 markus 51
3537 22 Dec 16 peter 52 void check_equality(double, double, test::Suite&, const std::string&,
1273 10 Apr 08 jari 53                     unsigned long int N=1);
1234 15 Mar 08 peter 54 utility::Matrix data(void);
1582 15 Oct 08 peter 55 utility::MatrixWeighted data_weighted(void);
1234 15 Mar 08 peter 56
1234 15 Mar 08 peter 57 template<class Distance>
1247 17 Mar 08 peter 58 void test_distance(Distance, theplu::yat::test::Suite&, unsigned long int N=1);
1234 15 Mar 08 peter 59
1234 15 Mar 08 peter 60 template<class Distance>
1247 17 Mar 08 peter 61 void test_duplicate(Distance, theplu::yat::test::Suite&, unsigned long int N=1);
1234 15 Mar 08 peter 62
1234 15 Mar 08 peter 63 template<class Distance>
1247 17 Mar 08 peter 64 void test_rescaling(Distance, theplu::yat::test::Suite&, unsigned long int N=1);
1234 15 Mar 08 peter 65
1234 15 Mar 08 peter 66 template<class Distance>
3537 22 Dec 16 peter 67 void test_unity_weights(Distance, theplu::yat::test::Suite&,
1247 17 Mar 08 peter 68                         unsigned long int N=1);
1234 15 Mar 08 peter 69
1234 15 Mar 08 peter 70 template<class Distance>
3537 22 Dec 16 peter 71 void test_self_distance(Distance, theplu::yat::test::Suite&,
1247 17 Mar 08 peter 72                         unsigned long int N=1);
1234 15 Mar 08 peter 73
1234 15 Mar 08 peter 74 template<class Distance>
1247 17 Mar 08 peter 75 void test_symmetry(Distance, theplu::yat::test::Suite&, unsigned long int N=1);
1234 15 Mar 08 peter 76
1234 15 Mar 08 peter 77 template<class Distance>
3537 22 Dec 16 peter 78 void test_zero_weight(Distance, theplu::yat::test::Suite&,
1247 17 Mar 08 peter 79                       unsigned long int N=1);
1234 15 Mar 08 peter 80
1234 15 Mar 08 peter 81 utility::Matrix weight(void);
1234 15 Mar 08 peter 82
1233 15 Mar 08 peter 83 int main(int argc, char* argv[])
3537 22 Dec 16 peter 84 {
1233 15 Mar 08 peter 85   theplu::yat::test::Suite suite(argc, argv);
1233 15 Mar 08 peter 86   suite.err() << "testing distance" << std::endl;
3537 22 Dec 16 peter 87
1120 21 Feb 08 peter 88   utility::Vector a(3,1);
889 25 Sep 07 markus 89   a(1) = 2;
1120 21 Feb 08 peter 90   utility::Vector b(3,0);
889 25 Sep 07 markus 91   b(2) = 1;
3537 22 Dec 16 peter 92
1050 07 Feb 08 peter 93   statistics::EuclideanDistance eucl_dist;
3534 21 Dec 16 peter 94   BOOST_CONCEPT_ASSERT((utility::DistanceConcept<statistics::EuclideanDistance>));
1234 15 Mar 08 peter 95   suite.err() << "testing EuclideanDistance" << std::endl;
1247 17 Mar 08 peter 96   test_distance(eucl_dist, suite, 100);
1050 07 Feb 08 peter 97   double dist=eucl_dist(a.begin(),a.end(),b.begin());
1665 20 Dec 08 peter 98   if (!suite.equal_fix(dist, 2.23606797749978988178, 1e-16)) {
1233 15 Mar 08 peter 99     suite.err() << "Error in unweighted Euclidean distance " << std::endl;
1233 15 Mar 08 peter 100     suite.add(false);
890 25 Sep 07 markus 101   }
3537 22 Dec 16 peter 102
1050 07 Feb 08 peter 103   statistics::PearsonDistance pear_dist;
1234 15 Mar 08 peter 104   suite.err() << "testing PearsonDistance" << std::endl;
1247 17 Mar 08 peter 105   test_distance(pear_dist, suite, 1000);
1050 07 Feb 08 peter 106   dist=pear_dist(a.begin(),a.end(),b.begin());
1582 15 Oct 08 peter 107   if (!suite.equal(dist, 1.5)) {
1233 15 Mar 08 peter 108     suite.err() << "Error in unweighted Pearson distance " << std::endl;
1233 15 Mar 08 peter 109     suite.add(false);
890 25 Sep 07 markus 110   }
3537 22 Dec 16 peter 111
3537 22 Dec 16 peter 112
890 25 Sep 07 markus 113   // Testing weighted versions
1121 22 Feb 08 peter 114   utility::Matrix m(2,3,1);
890 25 Sep 07 markus 115   m(0,1)=2;
890 25 Sep 07 markus 116   m(1,0)=0;
890 25 Sep 07 markus 117   m(1,1)=0;
1582 15 Oct 08 peter 118   utility::MatrixWeighted m_w(m);
1582 15 Oct 08 peter 119   m_w(0,0).weight()=0;
1582 15 Oct 08 peter 120   classifier::MatrixLookupWeighted mw(m_w);
890 25 Sep 07 markus 121   classifier::DataLookupWeighted1D aw(mw,0,true);
890 25 Sep 07 markus 122   classifier::DataLookupWeighted1D bw(mw,1,true);
3537 22 Dec 16 peter 123
1050 07 Feb 08 peter 124   dist=eucl_dist(aw.begin(),aw.end(),bw.begin());
3537 22 Dec 16 peter 125
1582 15 Oct 08 peter 126   if (!suite.equal_sqrt(dist, sqrt(6))) {
1233 15 Mar 08 peter 127     suite.err() << "Error in weighted Euclidean distance " << std::endl;
1233 15 Mar 08 peter 128     suite.add(false);
890 25 Sep 07 markus 129   }
3537 22 Dec 16 peter 130
1050 07 Feb 08 peter 131   dist=pear_dist(aw.begin(),aw.end(),bw.begin());
3537 22 Dec 16 peter 132
1582 15 Oct 08 peter 133   if (!suite.equal(dist, 2)) {
1233 15 Mar 08 peter 134     suite.err() << "Error in weighted Pearson distance " << std::endl;
1233 15 Mar 08 peter 135     suite.add(false);
890 25 Sep 07 markus 136   }
3537 22 Dec 16 peter 137
3537 22 Dec 16 peter 138
898 26 Sep 07 markus 139   // Test with std::vectors
898 26 Sep 07 markus 140   std::vector<double> sa(3,1);
898 26 Sep 07 markus 141   sa[1] = 2;
898 26 Sep 07 markus 142   std::vector<double> sb(3,0);
898 26 Sep 07 markus 143   sb[2] = 1;
3537 22 Dec 16 peter 144
1704 08 Jan 09 peter 145   double tolerance=1e-4;
3537 22 Dec 16 peter 146   dist=eucl_dist(sa.begin(),sa.end(),sb.begin());
1704 08 Jan 09 peter 147   if(!suite.equal_fix(dist, 2.23607,tolerance)) {
1233 15 Mar 08 peter 148     suite.err() << "Error in distance for std::vector " << std::endl;
1233 15 Mar 08 peter 149     suite.add(false);
898 26 Sep 07 markus 150   }
3537 22 Dec 16 peter 151
900 27 Sep 07 markus 152   // Test for a std::list and a std::vector
900 27 Sep 07 markus 153   std::list<double> la;
900 27 Sep 07 markus 154   std::copy(sa.begin(),sa.end(),std::back_inserter<std::list<double> >(la));
1050 07 Feb 08 peter 155   dist=eucl_dist(la.begin(),la.end(),sb.begin());
1704 08 Jan 09 peter 156   if(!suite.equal_fix(dist, 2.23607, tolerance) ) {
1233 15 Mar 08 peter 157     suite.err() << "Error in distance for std::list " << std::endl;
1233 15 Mar 08 peter 158     suite.add(false);
898 26 Sep 07 markus 159   }
1234 15 Mar 08 peter 160
3537 22 Dec 16 peter 161
1233 15 Mar 08 peter 162   return suite.return_value();
889 25 Sep 07 markus 163 }
889 25 Sep 07 markus 164
889 25 Sep 07 markus 165
3537 22 Dec 16 peter 166 void check_equality(double dist1, double dist2, test::Suite& suite,
1273 10 Apr 08 jari 167                     const std::string& msg, unsigned long int N)
1234 15 Mar 08 peter 168 {
1234 15 Mar 08 peter 169   if (!suite.equal(dist1, dist2, N)) {
1234 15 Mar 08 peter 170     suite.err() << "Error: " << msg << " failed.\n";
1234 15 Mar 08 peter 171     suite.add(false);
1234 15 Mar 08 peter 172   }
1234 15 Mar 08 peter 173 }
1234 15 Mar 08 peter 174
1234 15 Mar 08 peter 175
1234 15 Mar 08 peter 176 utility::Matrix data(void)
1234 15 Mar 08 peter 177 {
1234 15 Mar 08 peter 178   utility::Matrix res(2,10);
1234 15 Mar 08 peter 179   for (size_t i = 0; i<res.columns(); ++i){
1234 15 Mar 08 peter 180     res(0,i) = i*i+1;
1234 15 Mar 08 peter 181     res(1,i) = 2*i+3;
1234 15 Mar 08 peter 182   }
1234 15 Mar 08 peter 183   return res;
1234 15 Mar 08 peter 184 }
1234 15 Mar 08 peter 185
1582 15 Oct 08 peter 186 utility::MatrixWeighted data_weighted(void)
1582 15 Oct 08 peter 187 {
1582 15 Oct 08 peter 188   utility::Matrix x = data();
1582 15 Oct 08 peter 189   utility::Matrix w = weight();
1582 15 Oct 08 peter 190   utility::MatrixWeighted res(x.rows(), x.columns());
1582 15 Oct 08 peter 191   std::copy(x.begin(), x.end(), utility::data_iterator(res.begin()));
1582 15 Oct 08 peter 192   std::copy(w.begin(), w.end(), utility::weight_iterator(res.begin()));
1582 15 Oct 08 peter 193   return res;
1582 15 Oct 08 peter 194 }
1582 15 Oct 08 peter 195
1234 15 Mar 08 peter 196 template<class Distance>
1247 17 Mar 08 peter 197 void test_distance(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 198                    unsigned int long N)
1234 15 Mar 08 peter 199 {
2334 15 Oct 10 peter 200   BOOST_CONCEPT_ASSERT((utility::DistanceConcept<Distance>));
1247 17 Mar 08 peter 201   test_duplicate(dist, suite, N);
1247 17 Mar 08 peter 202   test_rescaling(dist, suite, N);
1247 17 Mar 08 peter 203   test_unity_weights(dist, suite, N);
1247 17 Mar 08 peter 204   test_self_distance(dist, suite, N);
1247 17 Mar 08 peter 205   test_symmetry(dist, suite, N);
1247 17 Mar 08 peter 206   test_zero_weight(dist, suite, N);
2202 21 Feb 10 peter 207
2202 21 Feb 10 peter 208   // do not run compiler test
2202 21 Feb 10 peter 209   if (false) {
2202 21 Feb 10 peter 210     boost::forward_iterator_archetype<double> iter;
2202 21 Feb 10 peter 211     boost::forward_iterator_archetype<utility::DataWeight> witer;
2202 21 Feb 10 peter 212     dist(iter, iter, iter);
2202 21 Feb 10 peter 213     dist(iter, iter, witer);
2202 21 Feb 10 peter 214     dist(witer, witer, iter);
2202 21 Feb 10 peter 215     dist(witer, witer, witer);
2202 21 Feb 10 peter 216   }
1234 15 Mar 08 peter 217 }
1234 15 Mar 08 peter 218
1234 15 Mar 08 peter 219 template<class Distance>
3537 22 Dec 16 peter 220 void test_duplicate(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 221                     unsigned long int N)
1234 15 Mar 08 peter 222 {
1582 15 Oct 08 peter 223   utility::MatrixWeighted x(data_weighted());
1582 15 Oct 08 peter 224   utility::MatrixWeighted mw(x.rows(), 2*x.columns());
1234 15 Mar 08 peter 225   for (size_t i=0; i<x.rows(); ++i){
1582 15 Oct 08 peter 226     std::copy(x.begin_row(i), x.end_row(i), mw.begin_row(i));
1582 15 Oct 08 peter 227     std::copy(x.begin_row(i), x.end_row(i), mw.begin_row(i)+x.columns());
1234 15 Mar 08 peter 228   }
1582 15 Oct 08 peter 229   double dist1 = dist(mw.begin_row(0), mw.end_row(0), mw.begin_row(1));
1582 15 Oct 08 peter 230   for (size_t i=0; i<x.columns(); ++i)
1582 15 Oct 08 peter 231     mw(0,i).weight()=0.0;
1582 15 Oct 08 peter 232   double dist2 = dist(mw.begin_row(0), mw.end_row(0), mw.begin_row(1));
1247 17 Mar 08 peter 233   check_equality(dist1, dist2, suite, "duplicate property", N);
1234 15 Mar 08 peter 234 }
1234 15 Mar 08 peter 235
1234 15 Mar 08 peter 236 template<class Distance>
3537 22 Dec 16 peter 237 void test_rescaling(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 238                     unsigned long int N)
1234 15 Mar 08 peter 239 {
1582 15 Oct 08 peter 240   utility::MatrixWeighted wx=data_weighted();
1582 15 Oct 08 peter 241   double dist1 = dist(wx.begin_row(0), wx.end_row(0), wx.begin_row(1));
3537 22 Dec 16 peter 242   // rescale weights
1582 15 Oct 08 peter 243   for (size_t i=0; i<wx.rows(); ++i)
1582 15 Oct 08 peter 244     for (size_t j=0; j<wx.columns(); ++j)
1582 15 Oct 08 peter 245       wx(i,j).weight() *= 2.13;
1582 15 Oct 08 peter 246   double dist2 = dist(wx.begin_row(0), wx.end_row(0), wx.begin_row(1));
1247 17 Mar 08 peter 247   check_equality(dist1, dist2, suite, "rescaling", N);
1234 15 Mar 08 peter 248 }
1234 15 Mar 08 peter 249
1234 15 Mar 08 peter 250 template<class Distance>
3537 22 Dec 16 peter 251 void test_unity_weights(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 252                         unsigned long int N)
1234 15 Mar 08 peter 253 {
1234 15 Mar 08 peter 254   utility::Matrix x=data();
1483 09 Sep 08 peter 255   utility::MatrixWeighted mw(x);
1483 09 Sep 08 peter 256   double dist1 = dist(mw.begin_row(0), mw.end_row(0), mw.begin_row(1));
1234 15 Mar 08 peter 257   double dist2 = dist(x.begin_row(0), x.end_row(0), x.begin_row(1));
1247 17 Mar 08 peter 258   check_equality(dist1, dist2, suite, "unity weights", N);
1234 15 Mar 08 peter 259 }
1234 15 Mar 08 peter 260
1234 15 Mar 08 peter 261 template<class Distance>
3537 22 Dec 16 peter 262 void test_self_distance(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 263                         unsigned long int N)
1234 15 Mar 08 peter 264 {
1234 15 Mar 08 peter 265   utility::Matrix x = data();
1234 15 Mar 08 peter 266   double self = dist(x.begin(), x.end(), x.begin());
1704 08 Jan 09 peter 267   if (!suite.equal_fix(self,0, N*std::numeric_limits<double>().epsilon()) ) {
1234 15 Mar 08 peter 268     suite.err() << "error: self distance is " << self << "\n"
1234 15 Mar 08 peter 269                 << "supposed to be zero.\n";
1234 15 Mar 08 peter 270     suite.add(false);
1234 15 Mar 08 peter 271   }
1234 15 Mar 08 peter 272 }
1234 15 Mar 08 peter 273
1234 15 Mar 08 peter 274
1234 15 Mar 08 peter 275 template<class Distance>
3537 22 Dec 16 peter 276 void test_symmetry(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 277                    unsigned long int N)
1234 15 Mar 08 peter 278 {
1234 15 Mar 08 peter 279   utility::Matrix x = data();
1234 15 Mar 08 peter 280   double distab = dist(x.begin_row(0), x.end_row(0), x.begin_row(1));
1234 15 Mar 08 peter 281   double distba = dist(x.begin_row(1), x.end_row(1), x.begin_row(0));
1247 17 Mar 08 peter 282   check_equality(distab, distba, suite, "symmetry test", N);
1234 15 Mar 08 peter 283 }
1234 15 Mar 08 peter 284
1234 15 Mar 08 peter 285
1234 15 Mar 08 peter 286 template<class Distance>
3537 22 Dec 16 peter 287 void test_zero_weight(Distance dist, theplu::yat::test::Suite& suite,
1247 17 Mar 08 peter 288                       unsigned long int N)
1234 15 Mar 08 peter 289 {
1582 15 Oct 08 peter 290   utility::MatrixWeighted wx=data_weighted();
1582 15 Oct 08 peter 291   wx(0,0).weight() = 0.0;
1582 15 Oct 08 peter 292   double dist1 = dist(wx.begin_row(0), wx.end_row(0), wx.begin_row(1));
1582 15 Oct 08 peter 293   wx(0,0).weight() = 100*std::numeric_limits<double>().epsilon();
1582 15 Oct 08 peter 294   double dist2 = dist(wx.begin_row(0), wx.end_row(0), wx.begin_row(1));
1247 17 Mar 08 peter 295   check_equality(dist1, dist2, suite, "zero weight", N);
1234 15 Mar 08 peter 296 }
1234 15 Mar 08 peter 297
1234 15 Mar 08 peter 298 utility::Matrix weight(void)
1234 15 Mar 08 peter 299 {
1234 15 Mar 08 peter 300   utility::Matrix res(2,10);
1234 15 Mar 08 peter 301   for (size_t i = 0; i<res.columns(); ++i){
1234 15 Mar 08 peter 302     res(0,i) = 1.0/(1+i);
1234 15 Mar 08 peter 303     res(1,i) = 1.0-0.1*i;
1234 15 Mar 08 peter 304   }
1234 15 Mar 08 peter 305   return res;
1234 15 Mar 08 peter 306 }