mirror of
https://gitlab.com/libeigen/eigen.git
synced 2025-06-04 18:54:00 +08:00
kill the LeastSquares module.
I didn't even put it in Eigen2Support because it requires several other modules. But if you want we can always create a new module, Eigen2Support_LeastSquares...
This commit is contained in:
Benoit Jacob
parent
664f2d4508
commit
38facbd55b
@ -4,4 +4,4 @@
|
||||
#include "QR"
|
||||
#include "SVD"
|
||||
#include "Geometry"
|
||||
#include "LeastSquares"
|
||||
#include "Eigenvalues"
|
||||
@ -1,28 +0,0 @@
|
||||
#ifndef EIGEN_REGRESSION_MODULE_H
|
||||
#define EIGEN_REGRESSION_MODULE_H
|
||||
|
||||
#include "Core"
|
||||
|
||||
#include "src/Core/util/DisableMSVCWarnings.h"
|
||||
|
||||
#include "Eigenvalues"
|
||||
#include "Geometry"
|
||||
|
||||
namespace Eigen {
|
||||
|
||||
/** \defgroup LeastSquares_Module LeastSquares module
|
||||
* This module provides linear regression and related features.
|
||||
*
|
||||
* \code
|
||||
* #include <Eigen/LeastSquares>
|
||||
* \endcode
|
||||
*/
|
||||
|
||||
#include "src/LeastSquares/LeastSquares.h"
|
||||
|
||||
} // namespace Eigen
|
||||
|
||||
#include "src/Core/util/EnableMSVCWarnings.h"
|
||||
|
||||
#endif // EIGEN_REGRESSION_MODULE_H
|
||||
/* vim: set filetype=cpp et sw=2 ts=2 ai: */
|
||||
@ -5,7 +5,6 @@ ADD_SUBDIRECTORY(SVD)
|
||||
ADD_SUBDIRECTORY(Cholesky)
|
||||
ADD_SUBDIRECTORY(Array)
|
||||
ADD_SUBDIRECTORY(Geometry)
|
||||
ADD_SUBDIRECTORY(LeastSquares)
|
||||
ADD_SUBDIRECTORY(Sparse)
|
||||
ADD_SUBDIRECTORY(Jacobi)
|
||||
ADD_SUBDIRECTORY(Householder)
|
||||
|
||||
@ -1,6 +0,0 @@
|
||||
FILE(GLOB Eigen_LeastSquares_SRCS "*.h")
|
||||
|
||||
INSTALL(FILES
|
||||
${Eigen_LeastSquares_SRCS}
|
||||
DESTINATION ${INCLUDE_INSTALL_DIR}/Eigen/src/LeastSquares COMPONENT Devel
|
||||
)
|
||||
@ -1,181 +0,0 @@
|
||||
// This file is part of Eigen, a lightweight C++ template library
|
||||
// for linear algebra.
|
||||
//
|
||||
// Copyright (C) 2006-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
|
||||
//
|
||||
// Eigen is free software; you can redistribute it and/or
|
||||
// modify it under the terms of the GNU Lesser General Public
|
||||
// License as published by the Free Software Foundation; either
|
||||
// version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Alternatively, you can redistribute it and/or
|
||||
// modify it under the terms of the GNU General Public License as
|
||||
// published by the Free Software Foundation; either version 2 of
|
||||
// the License, or (at your option) any later version.
|
||||
//
|
||||
// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
|
||||
// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
|
||||
// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
|
||||
// GNU General Public License for more details.
|
||||
//
|
||||
// You should have received a copy of the GNU Lesser General Public
|
||||
// License and a copy of the GNU General Public License along with
|
||||
// Eigen. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
#ifndef EIGEN_LEASTSQUARES_H
|
||||
#define EIGEN_LEASTSQUARES_H
|
||||
|
||||
/** \ingroup LeastSquares_Module
|
||||
*
|
||||
* \leastsquares_module
|
||||
*
|
||||
* For a set of points, this function tries to express
|
||||
* one of the coords as a linear (affine) function of the other coords.
|
||||
*
|
||||
* This is best explained by an example. This function works in full
|
||||
* generality, for points in a space of arbitrary dimension, and also over
|
||||
* the complex numbers, but for this example we will work in dimension 3
|
||||
* over the real numbers (doubles).
|
||||
*
|
||||
* So let us work with the following set of 5 points given by their
|
||||
* \f$(x,y,z)\f$ coordinates:
|
||||
* @code
|
||||
Vector3d points[5];
|
||||
points[0] = Vector3d( 3.02, 6.89, -4.32 );
|
||||
points[1] = Vector3d( 2.01, 5.39, -3.79 );
|
||||
points[2] = Vector3d( 2.41, 6.01, -4.01 );
|
||||
points[3] = Vector3d( 2.09, 5.55, -3.86 );
|
||||
points[4] = Vector3d( 2.58, 6.32, -4.10 );
|
||||
* @endcode
|
||||
* Suppose that we want to express the second coordinate (\f$y\f$) as a linear
|
||||
* expression in \f$x\f$ and \f$z\f$, that is,
|
||||
* \f[ y=ax+bz+c \f]
|
||||
* for some constants \f$a,b,c\f$. Thus, we want to find the best possible
|
||||
* constants \f$a,b,c\f$ so that the plane of equation \f$y=ax+bz+c\f$ fits
|
||||
* best the five above points. To do that, call this function as follows:
|
||||
* @code
|
||||
Vector3d coeffs; // will store the coefficients a, b, c
|
||||
linearRegression(
|
||||
5,
|
||||
&points,
|
||||
&coeffs,
|
||||
1 // the coord to express as a function of
|
||||
// the other ones. 0 means x, 1 means y, 2 means z.
|
||||
);
|
||||
* @endcode
|
||||
* Now the vector \a coeffs is approximately
|
||||
* \f$( 0.495 , -1.927 , -2.906 )\f$.
|
||||
* Thus, we get \f$a=0.495, b = -1.927, c = -2.906\f$. Let us check for
|
||||
* instance how near points[0] is from the plane of equation \f$y=ax+bz+c\f$.
|
||||
* Looking at the coords of points[0], we see that:
|
||||
* \f[ax+bz+c = 0.495 * 3.02 + (-1.927) * (-4.32) + (-2.906) = 6.91.\f]
|
||||
* On the other hand, we have \f$y=6.89\f$. We see that the values
|
||||
* \f$6.91\f$ and \f$6.89\f$
|
||||
* are near, so points[0] is very near the plane of equation \f$y=ax+bz+c\f$.
|
||||
*
|
||||
* Let's now describe precisely the parameters:
|
||||
* @param numPoints the number of points
|
||||
* @param points the array of pointers to the points on which to perform the linear regression
|
||||
* @param result pointer to the vector in which to store the result.
|
||||
This vector must be of the same type and size as the
|
||||
data points. The meaning of its coords is as follows.
|
||||
For brevity, let \f$n=Size\f$,
|
||||
\f$r_i=result[i]\f$,
|
||||
and \f$f=funcOfOthers\f$. Denote by
|
||||
\f$x_0,\ldots,x_{n-1}\f$
|
||||
the n coordinates in the n-dimensional space.
|
||||
Then the resulting equation is:
|
||||
\f[ x_f = r_0 x_0 + \cdots + r_{f-1}x_{f-1}
|
||||
+ r_{f+1}x_{f+1} + \cdots + r_{n-1}x_{n-1} + r_n. \f]
|
||||
* @param funcOfOthers Determines which coord to express as a function of the
|
||||
others. Coords are numbered starting from 0, so that a
|
||||
value of 0 means \f$x\f$, 1 means \f$y\f$,
|
||||
2 means \f$z\f$, ...
|
||||
*
|
||||
* \sa fitHyperplane()
|
||||
*/
|
||||
template<typename VectorType>
|
||||
void linearRegression(int numPoints,
|
||||
VectorType **points,
|
||||
VectorType *result,
|
||||
int funcOfOthers )
|
||||
{
|
||||
typedef typename VectorType::Scalar Scalar;
|
||||
typedef Hyperplane<Scalar, VectorType::SizeAtCompileTime> HyperplaneType;
|
||||
const int size = points[0]->size();
|
||||
result->resize(size);
|
||||
HyperplaneType h(size);
|
||||
fitHyperplane(numPoints, points, &h);
|
||||
for(int i = 0; i < funcOfOthers; i++)
|
||||
result->coeffRef(i) = - h.coeffs()[i] / h.coeffs()[funcOfOthers];
|
||||
for(int i = funcOfOthers; i < size; i++)
|
||||
result->coeffRef(i) = - h.coeffs()[i+1] / h.coeffs()[funcOfOthers];
|
||||
}
|
||||
|
||||
/** \ingroup LeastSquares_Module
|
||||
*
|
||||
* \leastsquares_module
|
||||
*
|
||||
* This function is quite similar to linearRegression(), so we refer to the
|
||||
* documentation of this function and only list here the differences.
|
||||
*
|
||||
* The main difference from linearRegression() is that this function doesn't
|
||||
* take a \a funcOfOthers argument. Instead, it finds a general equation
|
||||
* of the form
|
||||
* \f[ r_0 x_0 + \cdots + r_{n-1}x_{n-1} + r_n = 0, \f]
|
||||
* where \f$n=Size\f$, \f$r_i=retCoefficients[i]\f$, and we denote by
|
||||
* \f$x_0,\ldots,x_{n-1}\f$ the n coordinates in the n-dimensional space.
|
||||
*
|
||||
* Thus, the vector \a retCoefficients has size \f$n+1\f$, which is another
|
||||
* difference from linearRegression().
|
||||
*
|
||||
* In practice, this function performs an hyper-plane fit in a total least square sense
|
||||
* via the following steps:
|
||||
* 1 - center the data to the mean
|
||||
* 2 - compute the covariance matrix
|
||||
* 3 - pick the eigenvector corresponding to the smallest eigenvalue of the covariance matrix
|
||||
* The ratio of the smallest eigenvalue and the second one gives us a hint about the relevance
|
||||
* of the solution. This value is optionally returned in \a soundness.
|
||||
*
|
||||
* \sa linearRegression()
|
||||
*/
|
||||
template<typename VectorType, typename HyperplaneType>
|
||||
void fitHyperplane(int numPoints,
|
||||
VectorType **points,
|
||||
HyperplaneType *result,
|
||||
typename NumTraits<typename VectorType::Scalar>::Real* soundness = 0)
|
||||
{
|
||||
typedef typename VectorType::Scalar Scalar;
|
||||
typedef Matrix<Scalar,VectorType::SizeAtCompileTime,VectorType::SizeAtCompileTime> CovMatrixType;
|
||||
EIGEN_STATIC_ASSERT_VECTOR_ONLY(VectorType)
|
||||
ei_assert(numPoints >= 1);
|
||||
int size = points[0]->size();
|
||||
ei_assert(size+1 == result->coeffs().size());
|
||||
|
||||
// compute the mean of the data
|
||||
VectorType mean = VectorType::Zero(size);
|
||||
for(int i = 0; i < numPoints; ++i)
|
||||
mean += *(points[i]);
|
||||
mean /= numPoints;
|
||||
|
||||
// compute the covariance matrix
|
||||
CovMatrixType covMat = CovMatrixType::Zero(size, size);
|
||||
for(int i = 0; i < numPoints; ++i)
|
||||
{
|
||||
VectorType diff = (*(points[i]) - mean).conjugate();
|
||||
covMat += diff * diff.adjoint();
|
||||
}
|
||||
|
||||
// now we just have to pick the eigen vector with smallest eigen value
|
||||
SelfAdjointEigenSolver<CovMatrixType> eig(covMat);
|
||||
result->normal() = eig.eigenvectors().col(0);
|
||||
if (soundness)
|
||||
*soundness = eig.eigenvalues().coeff(0)/eig.eigenvalues().coeff(1);
|
||||
|
||||
// let's compute the constant coefficient such that the
|
||||
// plane pass trough the mean point:
|
||||
result->offset() = - (result->normal().cwiseProduct(mean)).sum();
|
||||
}
|
||||
|
||||
|
||||
#endif // EIGEN_LEASTSQUARES_H
|
||||
@ -202,7 +202,6 @@ ALIASES = "only_for_vectors=This is only for vectors (either row-
|
||||
"geometry_module=This is defined in the %Geometry module. \code #include <Eigen/Geometry> \endcode" \
|
||||
"householder_module=This is defined in the %Householder module. \code #include <Eigen/Householder> \endcode" \
|
||||
"jacobi_module=This is defined in the %Jacobi module. \code #include <Eigen/Jacobi> \endcode" \
|
||||
"leastsquares_module=This is defined in the %LeastSquares module. \code #include <Eigen/LeastSquares> \endcode" \
|
||||
"lu_module=This is defined in the %LU module. \code #include <Eigen/LU> \endcode" \
|
||||
"qr_module=This is defined in the %QR module. \code #include <Eigen/QR> \endcode" \
|
||||
"svd_module=This is defined in the %SVD module. \code #include <Eigen/SVD> \endcode" \
|
||||
|
||||
@ -156,7 +156,6 @@ ei_add_test(geo_eulerangles)
|
||||
ei_add_test(geo_hyperplane)
|
||||
ei_add_test(geo_parametrizedline)
|
||||
ei_add_test(geo_alignedbox)
|
||||
ei_add_test(regression)
|
||||
ei_add_test(stdvector)
|
||||
ei_add_test(stdvector_overload)
|
||||
ei_add_test(stdlist)
|
||||
|
||||
@ -1,153 +0,0 @@
|
||||
// This file is part of Eigen, a lightweight C++ template library
|
||||
// for linear algebra.
|
||||
//
|
||||
// Copyright (C) 2008 Benoit Jacob <jacob.benoit.1@gmail.com>
|
||||
//
|
||||
// Eigen is free software; you can redistribute it and/or
|
||||
// modify it under the terms of the GNU Lesser General Public
|
||||
// License as published by the Free Software Foundation; either
|
||||
// version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Alternatively, you can redistribute it and/or
|
||||
// modify it under the terms of the GNU General Public License as
|
||||
// published by the Free Software Foundation; either version 2 of
|
||||
// the License, or (at your option) any later version.
|
||||
//
|
||||
// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
|
||||
// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
|
||||
// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
|
||||
// GNU General Public License for more details.
|
||||
//
|
||||
// You should have received a copy of the GNU Lesser General Public
|
||||
// License and a copy of the GNU General Public License along with
|
||||
// Eigen. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
#include "main.h"
|
||||
#include <Eigen/LeastSquares>
|
||||
|
||||
template<typename VectorType,
|
||||
typename HyperplaneType>
|
||||
void makeNoisyCohyperplanarPoints(int numPoints,
|
||||
VectorType **points,
|
||||
HyperplaneType *hyperplane,
|
||||
typename VectorType::Scalar noiseAmplitude)
|
||||
{
|
||||
typedef typename VectorType::Scalar Scalar;
|
||||
const int size = points[0]->size();
|
||||
// pick a random hyperplane, store the coefficients of its equation
|
||||
hyperplane->coeffs().resize(size + 1);
|
||||
for(int j = 0; j < size + 1; j++)
|
||||
{
|
||||
do {
|
||||
hyperplane->coeffs().coeffRef(j) = ei_random<Scalar>();
|
||||
} while(ei_abs(hyperplane->coeffs().coeff(j)) < 0.5);
|
||||
}
|
||||
|
||||
// now pick numPoints random points on this hyperplane
|
||||
for(int i = 0; i < numPoints; i++)
|
||||
{
|
||||
VectorType& cur_point = *(points[i]);
|
||||
do
|
||||
{
|
||||
cur_point = VectorType::Random(size)/*.normalized()*/;
|
||||
// project cur_point onto the hyperplane
|
||||
Scalar x = - (hyperplane->coeffs().head(size).cwiseProduct(cur_point)).sum();
|
||||
cur_point *= hyperplane->coeffs().coeff(size) / x;
|
||||
} while( cur_point.norm() < 0.5
|
||||
|| cur_point.norm() > 2.0 );
|
||||
}
|
||||
|
||||
// add some noise to these points
|
||||
for(int i = 0; i < numPoints; i++ )
|
||||
*(points[i]) += noiseAmplitude * VectorType::Random(size);
|
||||
}
|
||||
|
||||
template<typename VectorType>
|
||||
void check_linearRegression(int numPoints,
|
||||
VectorType **points,
|
||||
const VectorType& original,
|
||||
typename VectorType::Scalar tolerance)
|
||||
{
|
||||
int size = points[0]->size();
|
||||
assert(size==2);
|
||||
VectorType result(size);
|
||||
linearRegression(numPoints, points, &result, 1);
|
||||
typename VectorType::Scalar error = (result - original).norm() / original.norm();
|
||||
VERIFY(ei_abs(error) < ei_abs(tolerance));
|
||||
}
|
||||
|
||||
template<typename VectorType,
|
||||
typename HyperplaneType>
|
||||
void check_fitHyperplane(int numPoints,
|
||||
VectorType **points,
|
||||
const HyperplaneType& original,
|
||||
typename VectorType::Scalar tolerance)
|
||||
{
|
||||
int size = points[0]->size();
|
||||
HyperplaneType result(size);
|
||||
fitHyperplane(numPoints, points, &result);
|
||||
result.coeffs() *= original.coeffs().coeff(size)/result.coeffs().coeff(size);
|
||||
typename VectorType::Scalar error = (result.coeffs() - original.coeffs()).norm() / original.coeffs().norm();
|
||||
VERIFY(ei_abs(error) < ei_abs(tolerance));
|
||||
}
|
||||
|
||||
void test_regression()
|
||||
{
|
||||
for(int i = 0; i < g_repeat; i++)
|
||||
{
|
||||
#ifdef EIGEN_TEST_PART_1
|
||||
{
|
||||
Vector2f points2f [1000];
|
||||
Vector2f *points2f_ptrs [1000];
|
||||
for(int i = 0; i < 1000; i++) points2f_ptrs[i] = &(points2f[i]);
|
||||
Vector2f coeffs2f;
|
||||
Hyperplane<float,2> coeffs3f;
|
||||
makeNoisyCohyperplanarPoints(1000, points2f_ptrs, &coeffs3f, 0.01f);
|
||||
coeffs2f[0] = -coeffs3f.coeffs()[0]/coeffs3f.coeffs()[1];
|
||||
coeffs2f[1] = -coeffs3f.coeffs()[2]/coeffs3f.coeffs()[1];
|
||||
CALL_SUBTEST(check_linearRegression(10, points2f_ptrs, coeffs2f, 0.05f));
|
||||
CALL_SUBTEST(check_linearRegression(100, points2f_ptrs, coeffs2f, 0.01f));
|
||||
CALL_SUBTEST(check_linearRegression(1000, points2f_ptrs, coeffs2f, 0.002f));
|
||||
}
|
||||
#endif
|
||||
|
||||
#ifdef EIGEN_TEST_PART_2
|
||||
{
|
||||
Vector2f points2f [1000];
|
||||
Vector2f *points2f_ptrs [1000];
|
||||
for(int i = 0; i < 1000; i++) points2f_ptrs[i] = &(points2f[i]);
|
||||
Hyperplane<float,2> coeffs3f;
|
||||
makeNoisyCohyperplanarPoints(1000, points2f_ptrs, &coeffs3f, 0.01f);
|
||||
CALL_SUBTEST(check_fitHyperplane(10, points2f_ptrs, coeffs3f, 0.05f));
|
||||
CALL_SUBTEST(check_fitHyperplane(100, points2f_ptrs, coeffs3f, 0.01f));
|
||||
CALL_SUBTEST(check_fitHyperplane(1000, points2f_ptrs, coeffs3f, 0.002f));
|
||||
}
|
||||
#endif
|
||||
|
||||
#ifdef EIGEN_TEST_PART_3
|
||||
{
|
||||
Vector4d points4d [1000];
|
||||
Vector4d *points4d_ptrs [1000];
|
||||
for(int i = 0; i < 1000; i++) points4d_ptrs[i] = &(points4d[i]);
|
||||
Hyperplane<double,4> coeffs5d;
|
||||
makeNoisyCohyperplanarPoints(1000, points4d_ptrs, &coeffs5d, 0.01);
|
||||
CALL_SUBTEST(check_fitHyperplane(10, points4d_ptrs, coeffs5d, 0.05));
|
||||
CALL_SUBTEST(check_fitHyperplane(100, points4d_ptrs, coeffs5d, 0.01));
|
||||
CALL_SUBTEST(check_fitHyperplane(1000, points4d_ptrs, coeffs5d, 0.002));
|
||||
}
|
||||
#endif
|
||||
|
||||
#ifdef EIGEN_TEST_PART_4
|
||||
{
|
||||
VectorXcd *points11cd_ptrs[1000];
|
||||
for(int i = 0; i < 1000; i++) points11cd_ptrs[i] = new VectorXcd(11);
|
||||
Hyperplane<std::complex<double>,Dynamic> *coeffs12cd = new Hyperplane<std::complex<double>,Dynamic>(11);
|
||||
makeNoisyCohyperplanarPoints(1000, points11cd_ptrs, coeffs12cd, 0.01);
|
||||
CALL_SUBTEST(check_fitHyperplane(100, points11cd_ptrs, *coeffs12cd, 0.025));
|
||||
CALL_SUBTEST(check_fitHyperplane(1000, points11cd_ptrs, *coeffs12cd, 0.006));
|
||||
delete coeffs12cd;
|
||||
for(int i = 0; i < 1000; i++) delete points11cd_ptrs[i];
|
||||
}
|
||||
#endif
|
||||
}
|
||||
}
|
||||
@ -202,7 +202,6 @@ ALIASES = "only_for_vectors=This is only for vectors (either row-
|
||||
"qr_module=This is defined in the %QR module. \code #include <Eigen/QR> \endcode" \
|
||||
"svd_module=This is defined in the %SVD module. \code #include <Eigen/SVD> \endcode" \
|
||||
"geometry_module=This is defined in the %Geometry module. \code #include <Eigen/Geometry> \endcode" \
|
||||
"leastsquares_module=This is defined in the %LeastSquares module. \code #include <Eigen/LeastSquares> \endcode" \
|
||||
"label=\bug" \
|
||||
"redstar=<a href='#warningarraymodule' style='color:red;text-decoration: none;'>*</a>" \
|
||||
"nonstableyet=\warning This is not considered to be part of the stable public API yet. Changes may happen in future releases. See \ref Experimental \"Experimental parts of Eigen\""
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user