OpenCV  4.1.2
Open Source Computer Vision
Introduction to Principal Component Analysis (PCA)

Goal

In this tutorial you will learn how to:

  • Use the OpenCV class cv::PCA to calculate the orientation of an object.

What is PCA?

Principal Component Analysis (PCA) is a statistical procedure that extracts the most important features of a dataset.

Consider that you have a set of 2D points as it is shown in the figure above. Each dimension corresponds to a feature you are interested in. Here some could argue that the points are set in a random order. However, if you have a better look you will see that there is a linear pattern (indicated by the blue line) which is hard to dismiss. A key point of PCA is the Dimensionality Reduction. Dimensionality Reduction is the process of reducing the number of the dimensions of the given dataset. For example, in the above case it is possible to approximate the set of points to a single line and therefore, reduce the dimensionality of the given points from 2D to 1D.

Moreover, you could also see that the points vary the most along the blue line, more than they vary along the Feature 1 or Feature 2 axes. This means that if you know the position of a point along the blue line you have more information about the point than if you only knew where it was on Feature 1 axis or Feature 2 axis.

Hence, PCA allows us to find the direction along which our data varies the most. In fact, the result of running PCA on the set of points in the diagram consist of 2 vectors called eigenvectors which are the principal components of the data set.

The size of each eigenvector is encoded in the corresponding eigenvalue and indicates how much the data vary along the principal component. The beginning of the eigenvectors is the center of all points in the data set. Applying PCA to N-dimensional data set yields N N-dimensional eigenvectors, N eigenvalues and 1 N-dimensional center point. Enough theory, let’s see how we can put these ideas into code.

How are the eigenvectors and eigenvalues computed?

The goal is to transform a given data set X of dimension p to an alternative data set Y of smaller dimension L. Equivalently, we are seeking to find the matrix Y, where Y is the Karhunen–Loève transform (KLT) of matrix X:

\[ \mathbf{Y} = \mathbb{K} \mathbb{L} \mathbb{T} \{\mathbf{X}\} \]

Organize the data set

Suppose you have data comprising a set of observations of p variables, and you want to reduce the data so that each observation can be described with only L variables, L < p. Suppose further, that the data are arranged as a set of n data vectors \( x_1...x_n \) with each \( x_i \) representing a single grouped observation of the p variables.

  • Write \( x_1...x_n \) as row vectors, each of which has p columns.
  • Place the row vectors into a single matrix X of dimensions \( n\times p \).

Calculate the empirical mean

  • Find the empirical mean along each dimension \( j = 1, ..., p \).
  • Place the calculated mean values into an empirical mean vector u of dimensions \( p\times 1 \).

    \[ \mathbf{u[j]} = \frac{1}{n}\sum_{i=1}^{n}\mathbf{X[i,j]} \]

Calculate the deviations from the mean

Mean subtraction is an integral part of the solution towards finding a principal component basis that minimizes the mean square error of approximating the data. Hence, we proceed by centering the data as follows:

  • Subtract the empirical mean vector u from each row of the data matrix X.
  • Store mean-subtracted data in the \( n\times p \) matrix B.

    \[ \mathbf{B} = \mathbf{X} - \mathbf{h}\mathbf{u^{T}} \]

    where h is an \( n\times 1 \) column vector of all 1s:

    \[ h[i] = 1, i = 1, ..., n \]

Find the covariance matrix

  • Find the \( p\times p \) empirical covariance matrix C from the outer product of matrix B with itself:

    \[ \mathbf{C} = \frac{1}{n-1} \mathbf{B^{*}} \cdot \mathbf{B} \]

    where * is the conjugate transpose operator. Note that if B consists entirely of real numbers, which is the case in many applications, the "conjugate transpose" is the same as the regular transpose.

Find the eigenvectors and eigenvalues of the covariance matrix

  • Compute the matrix V of eigenvectors which diagonalizes the covariance matrix C:

    \[ \mathbf{V^{-1}} \mathbf{C} \mathbf{V} = \mathbf{D} \]

    where D is the diagonal matrix of eigenvalues of C.

  • Matrix D will take the form of an \( p \times p \) diagonal matrix:

    \[ D[k,l] = \left\{\begin{matrix} \lambda_k, k = l \\ 0, k \neq l \end{matrix}\right. \]

    here, \( \lambda_j \) is the j-th eigenvalue of the covariance matrix C

  • Matrix V, also of dimension p x p, contains p column vectors, each of length p, which represent the p eigenvectors of the covariance matrix C.
  • The eigenvalues and eigenvectors are ordered and paired. The j th eigenvalue corresponds to the j th eigenvector.
Note
sources [1], [2] and special thanks to Svetlin Penkov for the original tutorial.

Source Code

Note
Another example using PCA for dimensionality reduction while maintaining an amount of variance can be found at opencv_source_code/samples/cpp/pca.cpp

Explanation

  • Read image and convert it to binary

Here we apply the necessary pre-processing procedures in order to be able to detect the objects of interest.

  • Extract objects of interest

Then find and filter contours by size and obtain the orientation of the remaining ones.

  • Extract orientation

Orientation is extracted by the call of getOrientation() function, which performs all the PCA procedure.

First the data need to be arranged in a matrix with size n x 2, where n is the number of data points we have. Then we can perform that PCA analysis. The calculated mean (i.e. center of mass) is stored in the cntr variable and the eigenvectors and eigenvalues are stored in the corresponding std::vector’s.

  • Visualize result

The final result is visualized through the drawAxis() function, where the principal components are drawn in lines, and each eigenvector is multiplied by its eigenvalue and translated to the mean position.

Results

The code opens an image, finds the orientation of the detected objects of interest and then visualizes the result by drawing the contours of the detected objects of interest, the center point, and the x-axis, y-axis regarding the extracted orientation.

cv::THRESH_OTSU
flag, use Otsu algorithm to choose the optimal threshold value
Definition: imgproc.hpp:323
cv::String
std::string String
Definition: cvstd.hpp:150
cv::Point_< int >
cv::Mat::clone
Mat clone() const CV_NODISCARD
Creates a full copy of the array and the underlying data.
cv::aruco::drawAxis
void drawAxis(InputOutputArray image, InputArray cameraMatrix, InputArray distCoeffs, InputArray rvec, InputArray tvec, float length)
Draw coordinate system axis from pose estimation.
cv::cvtColor
void cvtColor(InputArray src, OutputArray dst, int code, int dstCn=0)
Converts an image from one color space to another.
cv::THRESH_BINARY
Definition: imgproc.hpp:317
cv::sin
softdouble sin(const softdouble &a)
Sine.
cv::RETR_LIST
Definition: imgproc.hpp:413
cv::samples::findFile
cv::String findFile(const cv::String &relative_path, bool required=true, bool silentMode=false)
Try to find requested data file.
cv::threshold
double threshold(InputArray src, OutputArray dst, double thresh, double maxval, int type)
Applies a fixed-level threshold to each array element.
cv::waitKey
int waitKey(int delay=0)
Waits for a pressed key.
cv::sqrt
softfloat sqrt(const softfloat &a)
Square root.
cv::Point_::y
_Tp y
y coordinate of the point
Definition: types.hpp:187
cv::Point_::x
_Tp x
x coordinate of the point
Definition: types.hpp:186
highgui.hpp
core.hpp
cv::Scalar_< double >
cv::quality::quality_utils::scale
void scale(cv::Mat &mat, const cv::Mat &range, const T min, const T max)
Definition: quality_utils.hpp:90
cv::findContours
void findContours(InputArray image, OutputArrayOfArrays contours, OutputArray hierarchy, int mode, int method, Point offset=Point())
Finds contours in a binary image.
cv::cudev::atan2
__device__ __forceinline__ float1 atan2(const uchar1 &a, const uchar1 &b)
Definition: vec_math.hpp:812
cv::line
void line(InputOutputArray img, Point pt1, Point pt2, const Scalar &color, int thickness=1, int lineType=LINE_8, int shift=0)
Draws a line segment connecting two points.
cv::imread
Mat imread(const String &filename, int flags=IMREAD_COLOR)
Loads an image from a file.
cv::Mat::empty
bool empty() const
Returns true if the array has no elements.
cv::contourArea
double contourArea(InputArray contour, bool oriented=false)
Calculates a contour area.
cv::cos
softdouble cos(const softdouble &a)
Cosine.
cv::dnn::print
static void print(const MatShape &shape, const String &name="")
Definition: shape_utils.hpp:188
cv::findContours
void findContours(InputArray image, OutputArrayOfArrays contours, int mode, int method, Point offset=Point())
CV_64F
#define CV_64F
Definition: interface.h:79
cv::imshow
void imshow(const String &winname, InputArray mat)
Displays an image in the specified window.
cv::Scalar
Scalar_< double > Scalar
Definition: types.hpp:669
cv::drawContours
void drawContours(InputOutputArray image, InputArrayOfArrays contours, int contourIdx, const Scalar &color, int thickness=1, int lineType=LINE_8, InputArray hierarchy=noArray(), int maxLevel=INT_MAX, Point offset=Point())
Draws contours outlines or filled contours.
cv::PCA
Principal Component Analysis.
Definition: core.hpp:2392
cv::Point
Point2i Point
Definition: types.hpp:194
cv::Mat
n-dimensional dense array class
Definition: mat.hpp:791
cv::Point2d
Point_< double > Point2d
Definition: types.hpp:193
cv::imshow
void imshow(const String &winname, const ogl::Texture2D &tex)
Displays OpenGL 2D texture in the specified window.
cv::CommandLineParser
Designed for command line parsing.
Definition: utility.hpp:831
cv::CHAIN_APPROX_NONE
Definition: imgproc.hpp:429
cv::COLOR_BGR2GRAY
convert between RGB/BGR and grayscale, color conversions
Definition: imgproc.hpp:540
cv
"black box" representation of the file storage associated with a file on disk.
Definition: affine.hpp:51
imgproc.hpp
CV_PI
#define CV_PI
Definition: cvdef.h:320
cv::mean
Scalar mean(InputArray src, InputArray mask=noArray())
Calculates an average (mean) of array elements.
cv::Vec::channels
Definition: matx.hpp:337
cv::datasets::circle
Definition: gr_skig.hpp:62
cv::LINE_AA
antialiased line
Definition: imgproc.hpp:805
cv::circle
void circle(InputOutputArray img, Point center, int radius, const Scalar &color, int thickness=1, int lineType=LINE_8, int shift=0)
Draws a circle.