diff options
Diffstat (limited to 'ml/dlib/dlib/test/cca.cpp')
-rw-r--r-- | ml/dlib/dlib/test/cca.cpp | 460 |
1 files changed, 460 insertions, 0 deletions
diff --git a/ml/dlib/dlib/test/cca.cpp b/ml/dlib/dlib/test/cca.cpp new file mode 100644 index 000000000..a2014121d --- /dev/null +++ b/ml/dlib/dlib/test/cca.cpp @@ -0,0 +1,460 @@ +// Copyright (C) 2013 Davis E. King (davis@dlib.net) +// License: Boost Software License See LICENSE.txt for the full license. + +#include <dlib/statistics.h> +#include <dlib/sparse_vector.h> +#include <dlib/timing.h> +#include <map> + +#include "tester.h" + +namespace +{ + using namespace test; + using namespace dlib; + using namespace std; + + logger dlog("test.cca"); + + dlib::rand rnd; +// ---------------------------------------------------------------------------------------- + + /* + std::vector<std::map<unsigned long, double> > make_really_big_test_matrix ( + ) + { + std::vector<std::map<unsigned long,double> > temp(30000); + for (unsigned long i = 0; i < temp.size(); ++i) + { + for (int k = 0; k < 30; ++k) + temp[i][rnd.get_random_32bit_number()%10000] = 1; + } + return temp; + } + */ + + template <typename T> + std::vector<std::map<unsigned long, T> > mat_to_sparse ( + const matrix<T>& A + ) + { + std::vector<std::map<unsigned long,T> > temp(A.nr()); + for (long r = 0; r < A.nr(); ++r) + { + for (long c = 0; c < A.nc(); ++c) + { + temp[r][c] = A(r,c); + } + } + return temp; + } + +// ---------------------------------------------------------------------------------------- + + template <typename EXP> + matrix<typename EXP::type> rm_zeros ( + const matrix_exp<EXP>& m + ) + { + // Do this to avoid trying to correlate super small numbers that are really just + // zero. Doing this avoids some potential false alarms in the unit tests below. + return round_zeros(m, max(abs(m))*1e-14); + } + +// ---------------------------------------------------------------------------------------- + + /* + void check_correlation ( + matrix<double> L, + matrix<double> R, + const matrix<double>& Ltrans, + const matrix<double>& Rtrans, + const matrix<double,0,1>& correlations + ) + { + // apply the transforms + L = L*Ltrans; + R = R*Rtrans; + + // compute the real correlation values. Store them in A. + matrix<double> A = compute_correlations(L, R); + + for (long i = 0; i < correlations.size(); ++i) + { + // compare what the measured correlation values are (in A) to the + // predicted values. + cout << "error: "<< A(i) - correlations(i); + } + } + */ + +// ---------------------------------------------------------------------------------------- + + void test_cca3() + { + print_spinner(); + const unsigned long rank = rnd.get_random_32bit_number()%10 + 1; + const unsigned long m = rank + rnd.get_random_32bit_number()%15; + const unsigned long n = rank + rnd.get_random_32bit_number()%15; + const unsigned long n2 = rank + rnd.get_random_32bit_number()%15; + const unsigned long rank2 = rank + rnd.get_random_32bit_number()%5; + + dlog << LINFO << "m: " << m; + dlog << LINFO << "n: " << n; + dlog << LINFO << "n2: " << n2; + dlog << LINFO << "rank: " << rank; + dlog << LINFO << "rank2: " << rank2; + + + matrix<double> L = randm(m,rank, rnd)*randm(rank,n, rnd); + //matrix<double> R = randm(m,rank, rnd)*randm(rank,n2, rnd); + matrix<double> R = L*randm(n,n2, rnd); + //matrix<double> L = randm(m,n, rnd); + //matrix<double> R = randm(m,n2, rnd); + + matrix<double> Ltrans, Rtrans; + matrix<double,0,1> correlations; + + { + correlations = cca(L, R, Ltrans, Rtrans, min(m,n), max(n,n2)); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + dlog << LINFO << "correlations: "<< trans(correlations); + + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST_MSG(corr_error < 1e-13, Ltrans << "\n\n" << Rtrans); + + const double trans_error = max(abs(L*Ltrans - R*Rtrans)); + dlog << LINFO << "trans_error: "<< trans_error; + DLIB_TEST_MSG(trans_error < 1e-9, trans_error); + } + { + correlations = cca(mat_to_sparse(L), mat_to_sparse(R), Ltrans, Rtrans, min(m,n), max(n,n2)+6, 4); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + dlog << LINFO << "correlations: "<< trans(correlations); + dlog << LINFO << "computed cors: " << trans(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans))); + + const double trans_error = max(abs(L*Ltrans - R*Rtrans)); + dlog << LINFO << "trans_error: "<< trans_error; + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST_MSG(corr_error < 1e-13, Ltrans << "\n\n" << Rtrans); + + DLIB_TEST(trans_error < 2e-9); + } + + dlog << LINFO << "*****************************************************"; + } + + void test_cca2() + { + print_spinner(); + const unsigned long rank = rnd.get_random_32bit_number()%10 + 1; + const unsigned long m = rank + rnd.get_random_32bit_number()%15; + const unsigned long n = rank + rnd.get_random_32bit_number()%15; + const unsigned long n2 = rank + rnd.get_random_32bit_number()%15; + + dlog << LINFO << "m: " << m; + dlog << LINFO << "n: " << n; + dlog << LINFO << "n2: " << n2; + dlog << LINFO << "rank: " << rank; + + + matrix<double> L = randm(m,n, rnd); + matrix<double> R = randm(m,n2, rnd); + + matrix<double> Ltrans, Rtrans; + matrix<double,0,1> correlations; + + { + correlations = cca(L, R, Ltrans, Rtrans, min(n,n2), max(n,n2)-min(n,n2)); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + dlog << LINFO << "correlations: "<< trans(correlations); + + if (Ltrans.nc() > 1) + { + // The CCA projection directions are supposed to be uncorrelated for + // non-matching pairs of projections. + const double corr_rot1_error = max(abs(compute_correlations(rm_zeros(L*rotate<0,1>(Ltrans)), rm_zeros(R*Rtrans)))); + dlog << LINFO << "corr_rot1_error: "<< corr_rot1_error; + DLIB_TEST(std::abs(corr_rot1_error) < 1e-10); + } + // Matching projection directions should be correlated with the amount of + // correlation indicated by the return value of cca(). + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST(corr_error < 1e-13); + } + { + correlations = cca(mat_to_sparse(L), mat_to_sparse(R), Ltrans, Rtrans, min(n,n2), max(n,n2)-min(n,n2)); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + dlog << LINFO << "correlations: "<< trans(correlations); + + if (Ltrans.nc() > 1) + { + // The CCA projection directions are supposed to be uncorrelated for + // non-matching pairs of projections. + const double corr_rot1_error = max(abs(compute_correlations(rm_zeros(L*rotate<0,1>(Ltrans)), rm_zeros(R*Rtrans)))); + dlog << LINFO << "corr_rot1_error: "<< corr_rot1_error; + DLIB_TEST(std::abs(corr_rot1_error) < 1e-10); + } + // Matching projection directions should be correlated with the amount of + // correlation indicated by the return value of cca(). + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST(corr_error < 1e-13); + } + + dlog << LINFO << "*****************************************************"; + } + + void test_cca1() + { + print_spinner(); + const unsigned long rank = rnd.get_random_32bit_number()%10 + 1; + const unsigned long m = rank + rnd.get_random_32bit_number()%15; + const unsigned long n = rank + rnd.get_random_32bit_number()%15; + + dlog << LINFO << "m: " << m; + dlog << LINFO << "n: " << n; + dlog << LINFO << "rank: " << rank; + + matrix<double> T = randm(n,n, rnd); + + matrix<double> L = randm(m,rank, rnd)*randm(rank,n, rnd); + //matrix<double> L = randm(m,n, rnd); + matrix<double> R = L*T; + + matrix<double> Ltrans, Rtrans; + matrix<double,0,1> correlations; + + { + correlations = cca(L, R, Ltrans, Rtrans, rank); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + if (Ltrans.nc() > 1) + { + // The CCA projection directions are supposed to be uncorrelated for + // non-matching pairs of projections. + const double corr_rot1_error = max(abs(compute_correlations(rm_zeros(L*rotate<0,1>(Ltrans)), rm_zeros(R*Rtrans)))); + dlog << LINFO << "corr_rot1_error: "<< corr_rot1_error; + DLIB_TEST(std::abs(corr_rot1_error) < 2e-9); + } + // Matching projection directions should be correlated with the amount of + // correlation indicated by the return value of cca(). + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST(corr_error < 1e-13); + + const double trans_error = max(abs(L*Ltrans - R*Rtrans)); + dlog << LINFO << "trans_error: "<< trans_error; + DLIB_TEST(trans_error < 2e-9); + + dlog << LINFO << "correlations: "<< trans(correlations); + } + { + correlations = cca(mat_to_sparse(L), mat_to_sparse(R), Ltrans, Rtrans, rank); + DLIB_TEST(Ltrans.nc() == Rtrans.nc()); + if (Ltrans.nc() > 1) + { + // The CCA projection directions are supposed to be uncorrelated for + // non-matching pairs of projections. + const double corr_rot1_error = max(abs(compute_correlations(rm_zeros(L*rotate<0,1>(Ltrans)), rm_zeros(R*Rtrans)))); + dlog << LINFO << "corr_rot1_error: "<< corr_rot1_error; + DLIB_TEST(std::abs(corr_rot1_error) < 2e-9); + } + // Matching projection directions should be correlated with the amount of + // correlation indicated by the return value of cca(). + const double corr_error = max(abs(compute_correlations(rm_zeros(L*Ltrans), rm_zeros(R*Rtrans)) - correlations)); + dlog << LINFO << "correlation error: "<< corr_error; + DLIB_TEST(corr_error < 1e-13); + + const double trans_error = max(abs(L*Ltrans - R*Rtrans)); + dlog << LINFO << "trans_error: "<< trans_error; + DLIB_TEST(trans_error < 2e-9); + + dlog << LINFO << "correlations: "<< trans(correlations); + } + + dlog << LINFO << "*****************************************************"; + } + +// ---------------------------------------------------------------------------------------- + + void test_svd_fast( + long rank, + long m, + long n + ) + { + print_spinner(); + matrix<double> A = randm(m,rank,rnd)*randm(rank,n,rnd); + matrix<double> u,v; + matrix<double,0,1> w; + + dlog << LINFO << "rank: "<< rank; + dlog << LINFO << "m: "<< m; + dlog << LINFO << "n: "<< n; + + svd_fast(A, u, w, v, rank, 2); + DLIB_TEST(u.nr() == m); + DLIB_TEST(u.nc() == rank); + DLIB_TEST(w.nr() == rank); + DLIB_TEST(w.nc() == 1); + DLIB_TEST(v.nr() == n); + DLIB_TEST(v.nc() == rank); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-13); + svd_fast(mat_to_sparse(A), u, w, v, rank, 2); + DLIB_TEST(u.nr() == m); + DLIB_TEST(u.nc() == rank); + DLIB_TEST(w.nr() == rank); + DLIB_TEST(w.nc() == 1); + DLIB_TEST(v.nr() == n); + DLIB_TEST(v.nc() == rank); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-13); + + svd_fast(A, u, w, v, rank, 0); + DLIB_TEST(u.nr() == m); + DLIB_TEST(u.nc() == rank); + DLIB_TEST(w.nr() == rank); + DLIB_TEST(w.nc() == 1); + DLIB_TEST(v.nr() == n); + DLIB_TEST(v.nc() == rank); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST_MSG(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-9,max(abs(tmp(A - u*diagm(w)*trans(v))))); + svd_fast(mat_to_sparse(A), u, w, v, rank, 0); + DLIB_TEST(u.nr() == m); + DLIB_TEST(u.nc() == rank); + DLIB_TEST(w.nr() == rank); + DLIB_TEST(w.nc() == 1); + DLIB_TEST(v.nr() == n); + DLIB_TEST(v.nc() == rank); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-10); + + svd_fast(A, u, w, v, rank+5, 0); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-11); + svd_fast(mat_to_sparse(A), u, w, v, rank+5, 0); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-11); + svd_fast(A, u, w, v, rank+5, 1); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-12); + svd_fast(mat_to_sparse(A), u, w, v, rank+5, 1); + DLIB_TEST(max(abs(trans(u)*u - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(trans(v)*v - identity_matrix<double>(u.nc()))) < 1e-13); + DLIB_TEST(max(abs(tmp(A - u*diagm(w)*trans(v)))) < 1e-12); + } + + void test_svd_fast() + { + for (int iter = 0; iter < 1000; ++iter) + { + const unsigned long rank = rnd.get_random_32bit_number()%10 + 1; + const unsigned long m = rank + rnd.get_random_32bit_number()%10; + const unsigned long n = rank + rnd.get_random_32bit_number()%10; + + test_svd_fast(rank, m, n); + + } + test_svd_fast(1, 1, 1); + test_svd_fast(1, 2, 2); + test_svd_fast(1, 1, 2); + test_svd_fast(1, 2, 1); + } + +// ---------------------------------------------------------------------------------------- + + /* + typedef std::vector<std::pair<unsigned int, float>> sv; + sv rand_sparse_vector() + { + static dlib::rand rnd; + sv v; + for (int i = 0; i < 50; ++i) + v.push_back(make_pair(rnd.get_integer(400000), rnd.get_random_gaussian()*100)); + + make_sparse_vector_inplace(v); + return v; + } + + sv rand_basis_combo(const std::vector<sv>& basis) + { + static dlib::rand rnd; + sv result; + + for (int i = 0; i < 5; ++i) + { + sv temp = basis[rnd.get_integer(basis.size())]; + scale_by(temp, rnd.get_random_gaussian()); + result = add(result,temp); + } + return result; + } + + void big_sparse_speed_test() + { + cout << "making A" << endl; + std::vector<sv> basis; + for (int i = 0; i < 100; ++i) + basis.emplace_back(rand_sparse_vector()); + + std::vector<sv> A; + for (int i = 0; i < 500000; ++i) + A.emplace_back(rand_basis_combo(basis)); + + cout << "done making A" << endl; + + matrix<float> u,v; + matrix<float,0,1> w; + { + timing::block aosijdf(0,"call it"); + svd_fast(A, u,w,v, 100, 5); + } + + timing::print(); + } + */ + +// ---------------------------------------------------------------------------------------- + + class test_cca : public tester + { + public: + test_cca ( + ) : + tester ("test_cca", + "Runs tests on the cca() and svd_fast() routines.") + {} + + void perform_test ( + ) + { + //big_sparse_speed_test(); + for (int i = 0; i < 200; ++i) + { + test_cca1(); + test_cca2(); + test_cca3(); + } + test_svd_fast(); + } + } a; + + + +} + + + + |