# feature matrix

*Entry of the [Dictionary of Applied Machine Learning](https://dictionaryofml.org) (course edition), by Alexander Jung, Aalto University. CC BY 4.0. Updated 2026-11-11. [HTML](https://dictionaryofml.org/terms/featuremtx.html) · [typeset PDF](https://dictionaryofml.org/terms/featuremtx.pdf) · [LaTeX](https://dictionaryofml.org/terms/featuremtx.tex) · [Python demo](https://dictionaryofml.org/terms/featuremtx.py) · [notebook](https://dictionaryofml.org/terms/featuremtx.ipynb)*

Math is LaTeX with the book's macros expanded, so no notation package is needed to read it.

## Abstract

The feature matrix of a dataset collects the feature vectors of its data points as its rows, one row per data point and one column per feature. Its two most basic properties are the number of rows, which is the number of data points in the dataset, and the number of columns, which is the number of features of each data point. More properties of a dataset follow from the singular values of its feature matrix. The number of nonzero singular values is the rank, and a rank below the number of features leaves linear regression with infinitely many solutions. The ratio of the largest to the smallest nonzero singular value is the condition number of the feature matrix. The condition number determines the number of gradient descent steps required by linear regression and its robustness against perturbations of the features and labels in the training set.

## Definition

A weather station at Krems an der Donau records eight measurements on each of the $366$ days of the leap year $2024$: the lowest, highest, and mean air temperature, the humidity, the wind speed, the precipitation, the sunshine duration, and the air pressure. Writing the eight measurements of one day as a row, and placing the rows of the $366$ days under one another, gives a two-dimensional numeric array with $366$ rows and eight columns. That array is a matrix, the feature matrix of the dataset whose data points are the $366$ days, and each of its rows is the feature vector of one day. The feature matrix is the form in which machine learning software expects a dataset, and properties of the machine learning methods applied to the dataset, such as whether linear regression has a unique solution or how fast gradient descent converges, are determined by it.

![Figure 1 of the entry featuremtx](https://dictionaryofml.org/terms/fig/featuremtx-fig1.svg)

*Figure 1: The feature matrix of the $366$ days of $2024$ at Krems an der Donau: one row per day and one column per measurement, with the first three days and the last shown, and five of the eight measurements (the precipitation, the sunshine duration, and the air pressure are omitted). Temperatures are in degrees Celsius, the humidity in percent, and the wind speed in meters per second. The shaded row is the transposed feature vector of 2 January; the heavily framed column collects the lowest temperature of every day*

In general, consider a dataset $\mathcal{D}$ of $m$ data points with feature vectors ${\mathbf x}^{(1)}, \,\ldots, \,{\mathbf x}^{(m)}$, each in $\mathbb{R}^{d}$. Stacking the feature vectors row-wise yields the feature matrix:

$$
{\mathbf X} := \big({\mathbf x}^{(1)}, \,\ldots, \,{\mathbf x}^{(m)}\big)^{\top} =
\begin{pmatrix}
x^{(1)}_{1} & x^{(1)}_{2} & \cdots & x^{(1)}_{d} \\
x^{(2)}_{1} & x^{(2)}_{2} & \cdots & x^{(2)}_{d} \\
\vdots & \vdots & \ddots & \vdots \\
x^{(m)}_{1} & x^{(m)}_{2} & \cdots & x^{(m)}_{d}
\end{pmatrix}
\in \mathbb{R}^{m \times d} \text{.}
$$

The feature matrix is of size $m \times d$, i.e., it has $m$ rows and $d$ columns. Each row is the transposed feature vector of one data point, and each column holds the values of one feature for all data points (see Fig. 1). Many properties of a dataset, and of the learning tasks built from it, are computed from ${\mathbf X}$. The most basic ones are the numbers of rows and columns, $m$ and $d$.

The feature matrix also defines a linear map. A linear model with model parameters ${\mathbf w}$ from the parameter space $\mathbb{R}^{d}$ predicts ${\mathbf w}^{\top} {\mathbf x}^{(r)}$ for the $r$-th data point, so the $m$ predictions are the entries of ${\mathbf X} {\mathbf w}$. The feature matrix therefore maps each ${\mathbf w} \in \mathbb{R}^{d}$ to the vector ${\mathbf X} {\mathbf w} \in \mathbb{R}^{m}$, which is compared with the label vector ${\mathbf y} = \big(y^{(1)}, \,\ldots,
\,y^{(m)}\big)^{\top} \in \mathbb{R}^{m}$. Linear regression on a training set with feature matrix ${\mathbf X}$ and label vector ${\mathbf y}$ is the search for a vector ${\mathbf w}$ whose image ${\mathbf X} {\mathbf w}$ is close to ${\mathbf y}$ in Euclidean norm.

Linear regression has a unique solution only when ${\mathbf X}$ has full column rank, which requires $m \geq d$ (Golub and Loan, 2013, Sect. 5.3.1) (see normal equations). A rank below $d$ leaves infinitely many parameter vectors ${\mathbf w}$ with the same training error (Hastie et al., 2009, Sect. 3.2), and the one a method returns can have a large loss outside the training set. Using more features than data points (i.e., $d > m$) forces this case, since then ${\mathbf X}$ has fewer rows than columns, which rules out full column rank.

More properties of a dataset follow from the singular values of its feature matrix ${\mathbf X}$. The singular value decomposition can be written as (Golub and Loan, 2013, Sect. 2.4)

$$
{\mathbf X} = \sum_{j=1}^{\min\{m, d\}}
\lambda_{j} {\mathbf v}^{(j)}
\big({\mathbf u}^{(j)}\big)^{\top} \text{,}
$$

where the vectors ${\mathbf v}^{(1)}, \,\ldots, \,{\mathbf v}^{(\min\{m, d\})} \in \mathbb{R}^{m}$ are orthonormal among themselves, the vectors ${\mathbf u}^{(1)}, \,\ldots, \,{\mathbf u}^{(\min\{m, d\})} \in \mathbb{R}^{d}$ are orthonormal among themselves. The singular values are ordered as $\lambda_{1} \geq \ldots \geq \lambda_{\min\{m, d\}} \geq 0$. The number of nonzero singular values is the rank of ${\mathbf X}$ (Golub and Loan, 2013, Cor. 2.4.6). A rank below $d$ means that some feature is a linear combination of the others and adds no information beyond them. A small but nonzero singular value is the approximate form of such a dependence. When $m \geq d$, the $d$ right singular vectors form an orthonormal basis of $\mathbb{R}^{d}$, and the smallest singular value is the smallest Euclidean norm of a linear combination of the feature columns with unit-norm coefficients,

$$
\lambda_{d}
= \min_{{\mathbf u} \in \mathbb{R}^{d}, \, \Vert {{\mathbf u}} \Vert = 1}
\Vert {{\mathbf X} {\mathbf u}} \Vert \text{,}
$$

attained at the right singular vector ${\mathbf u}^{(d)}$. This holds since a unit-norm ${\mathbf u} = \sum_{j=1}^{d} c_{j}
{\mathbf u}^{(j)}$ gives $\Vert {{\mathbf X} {\mathbf u}} \Vert^{2} =
\sum_{j=1}^{d} c_{j}^{2} \lambda_{j}^{2}$, a weighted average of the squared singular values. When $d > m$, a unit-norm ${\mathbf u}$ orthogonal to all $m$ right singular vectors gives ${\mathbf X} {\mathbf u} = \mathbf{0}$, so the minimum is $0$. A small singular value therefore says that some combination of the feature columns with unit-norm coefficients is nearly the zero vector: the features are nearly collinear, and the right singular vector names the combination. For the weather year of Fig. 1, with every measurement scaled to zero sample mean and unit sample variance, the smallest singular value is $0.05$ against a largest of $36.4$. Its right singular vector places all of its weight on the lowest, the highest, and the mean temperature, with coefficients $-0.37$, $-0.46$, and $0.81$, because the mean temperature of a day is the midpoint of its lowest and its highest up to a residual with standard deviation $0.03$ degrees. Dropping the mean temperature raises the smallest singular value to $4.19$ (pythondemos/featuremtx.py). The size of a singular value depends on the units of the columns, which is why the scaling to unit sample variance comes first. The squares of the singular values are the eigenvalues of ${\mathbf X}^{\top} {\mathbf X}$, which for features with zero sample mean is $m$ times their sample covariance matrix.

The same eigenvalues and eigenvectors of ${\mathbf X}^{\top} {\mathbf X}$ determine principal component analysis. Its dimensionality reduction keeps the eigenvectors of the largest eigenvalues as the new features, and the sum over the data points of the squared errors of reconstructing the feature vectors from them is the sum of the eigenvalues dropped (Bishop, 2006, Sect. 12.1.2).

Another important property of a feature matrix is the condition number $\kappa\!\left({\mathbf X}\right) = \lambda_{1} / \lambda_{\operatorname{rank}({\mathbf X})}$, the ratio of the largest singular value to the smallest nonzero one. For full column rank, it is the ratio of the largest to the smallest singular value (Golub and Loan, 2013, Sect. 5.3.1). The condition number enters linear regression twice: through the number of gradient descent steps, and through the robustness of the learned model parameters against perturbations of the training set. The number that matters in both is its square $\kappa\!\left({\mathbf X}\right)^{2}$. The objective function $f({\mathbf w}) = (1/m) \Vert {{\mathbf X} {\mathbf w} -
{\mathbf y}} \Vert^{2}$ has the constant Hessian $(2/m) {\mathbf X}^{\top} {\mathbf X}$, whose eigenvalues are the squared singular values of ${\mathbf X}$ scaled by $2/m$. The condition number of ${\mathbf X}^{\top} {\mathbf X}$ is therefore $\kappa\!\left({\mathbf X}\right)^{2}$.

The square $\kappa\!\left({\mathbf X}\right)^{2}$ determines the convergence speed of gradient descent on the empirical risk minimization problem of linear regression. One step with an exact line search multiplies the gap between $f\big({\mathbf w}^{(t)}\big)$ and the minimum $f\big(\widehat{{\mathbf w}}\big)$ by at most $1 - 1/\kappa\!\left({\mathbf X}\right)^{2}$ (Boyd and Vandenberghe, 2004, Sect. 9.3.1). The number of steps that reaches a given gap therefore grows in proportion to $\kappa\!\left({\mathbf X}\right)^{2}$.

Fig. 2 shows this for two feature matrices built from the weather measurements of Fig. 1. The features are the relative humidity and the wind speed, the label is the highest temperature, and each has its sample mean removed. In the recorded units the condition number is $14.5$. After scaling each feature to unit sample variance, it is $1.45$. Gradient descent starts at the point at distance $\Vert {\widehat{{\mathbf w}}} \Vert$ from $\widehat{{\mathbf w}}$ for which it is slowest, where the gap shrinks by exactly $\big((\kappa\!\left({\mathbf X}\right)^{2} - 1) /
(\kappa\!\left({\mathbf X}\right)^{2} + 1)\big)^{2}$ per step (Boyd and Vandenberghe, 2004, Sect. 9.3.2); from no start does the gap shrink by less, since this factor bounds the shrinkage of every exact-line-search step on a quadratic objective function (Bertsekas, 2016, Prop. 1.3.1). A gap of $10^{-6}$ of the initial one then takes $731$ steps in the recorded units and $7$ after the scaling.

![Figure 2 of the entry featuremtx](https://dictionaryofml.org/terms/fig/featuremtx-fig2.svg)

*Figure 2: The convergence of gradient descent with an exact line search for the features relative humidity and wind speed of the $366$ days, in their recorded units (solid curve) and after scaling each to unit sample variance (dashed curve). Shown is the gap between $f\big({\mathbf w}^{(t)}\big)$ and the minimum $f^{\star} = \min_{{\mathbf w} \in \mathbb{R}^{2}} f({\mathbf w})$, relative to its initial value, over the first $50$ steps. Both runs start at the point at distance $\Vert {\widehat{{\mathbf w}}} \Vert$ from $\widehat{{\mathbf w}}$ for which gradient descent is slowest; the dashed curve ends where the gap falls below $10^{-12}$. Data generated by pythondemos/featuremtx.py*

The condition number $\kappa\!\left({\mathbf X}\right)$ also governs the robustness of the model parameters $\widehat{{\mathbf w}}$ delivered by linear regression against perturbations of the training set. Let ${\mathbf X}$ have full column rank. Linear regression delivers the model parameters $\widehat{{\mathbf w}} = {{\mathbf X}}^{+} {\mathbf y}$ with the pseudoinverse ${{\mathbf X}}^{+} = \big({\mathbf X}^{\top} {\mathbf X}\big)^{-1}
{\mathbf X}^{\top}$ (Golub and Loan, 2013, Sect. 5.5.2). They are the unique solution of the normal equations ${\mathbf X}^{\top} {\mathbf X} \widehat{{\mathbf w}} =
{\mathbf X}^{\top} {\mathbf y}$, which state that the gradient of $f$ vanishes at $\widehat{{\mathbf w}}$.

Consider perturbations of the training set that change ${\mathbf X}$ to ${\mathbf X} + \delta {\mathbf X}$ and ${\mathbf y}$ to ${\mathbf y} + \delta {\mathbf y}$, with ${\mathbf X} + \delta {\mathbf X}$ still of full column rank, and let $\widehat{{\mathbf w}} + \delta \widehat{{\mathbf w}}$ be the model parameters learned by linear regression on the perturbed training set. They solve the normal equations of the perturbed training set,

$$
({\mathbf X} + \delta {\mathbf X})^{\top}
({\mathbf X} + \delta {\mathbf X})
\big(\widehat{{\mathbf w}} + \delta \widehat{{\mathbf w}}\big)
= ({\mathbf X} + \delta {\mathbf X})^{\top}
({\mathbf y} + \delta {\mathbf y}) \text{.}
$$

Moving the terms with $\widehat{{\mathbf w}}$ to the right side and using ${\mathbf X}^{\top} {\mathbf r} = \mathbf{0}$ for the least squares residual ${\mathbf r} = {\mathbf y} - {\mathbf X} \widehat{{\mathbf w}}$, which is the normal equations of the unperturbed training set, gives the exact equation

$$
({\mathbf X} + \delta {\mathbf X})^{\top}
({\mathbf X} + \delta {\mathbf X}) \, \delta \widehat{{\mathbf w}}
= ({\mathbf X} + \delta {\mathbf X})^{\top}
\big(\delta {\mathbf y} - \delta {\mathbf X} \, \widehat{{\mathbf w}}\big)
+ \delta {\mathbf X}^{\top} {\mathbf r} \text{.}
$$

Solving for $\delta \widehat{{\mathbf w}}$ applies the pseudoinverse of ${\mathbf X} + \delta {\mathbf X}$ to the noise in the labels and to the noise in the features acting on $\widehat{{\mathbf w}}$, and the inverse of $({\mathbf X} + \delta {\mathbf X})^{\top}
({\mathbf X} + \delta {\mathbf X})$ to the noise in the features acting on the residual.

The largest singular value of ${{\mathbf X}}^{+}$ is $1 / \lambda_{d}$, and that of $\big({\mathbf X}^{\top} {\mathbf X}\big)^{-1}$ is $1 / \lambda_{d}^{2}$ (Golub and Loan, 2013, Sect. 5.3.6). The same holds for ${\mathbf X} + \delta {\mathbf X}$ with its own smallest singular value. Let the largest singular value of $\delta {\mathbf X}$ be at most $\varepsilon \lambda_{1}$ and the Euclidean norm of $\delta {\mathbf y}$ at most $\varepsilon \Vert {{\mathbf y}} \Vert$, with $\varepsilon \, \kappa\!\left({\mathbf X}\right) < 1$. The smallest singular value of ${\mathbf X} + \delta {\mathbf X}$ is then at least $\lambda_{d} - \varepsilon \lambda_{1}$ (Golub and Loan, 2013, Cor. 8.6.2), and taking Euclidean norms in the exact equation gives the perturbation bound

\begin{equation}
\tag{1}
\Vert {\delta \widehat{{\mathbf w}}} \Vert \leq
\frac{\varepsilon \, \kappa\!\left({\mathbf X}\right)}
{1 - \varepsilon \, \kappa\!\left({\mathbf X}\right)}
\Big( \frac{\Vert {{\mathbf y}} \Vert}{\lambda_{1}}
+ \Vert {\widehat{{\mathbf w}}} \Vert \Big)
+ \frac{\varepsilon \, \kappa\!\left({\mathbf X}\right)^{2}}
{\big(1 - \varepsilon \, \kappa\!\left({\mathbf X}\right)\big)^{2}}
\frac{\Vert {{\mathbf r}} \Vert}{\lambda_{1}} \text{.}
\end{equation}

The bound Eq. (1) holds for every perturbation with $\varepsilon \, \kappa\!\left({\mathbf X}\right) < 1$.

The bound Eq. (1) carries over to the prediction for any given feature vector ${\mathbf x}$, such as that of a test data point. The learned linear model predicts ${\mathbf x}^{\top} \widehat{{\mathbf w}}$ before and ${\mathbf x}^{\top} \big(\widehat{{\mathbf w}} + \delta \widehat{{\mathbf w}}\big)$ after the perturbation. The prediction therefore changes by ${\mathbf x}^{\top} \delta \widehat{{\mathbf w}}$, and the Cauchy-Schwarz inequality bounds this change by

\begin{equation}
\tag{2}
\big| {\mathbf x}^{\top} \delta \widehat{{\mathbf w}} \big|
\leq \Vert {{\mathbf x}} \Vert \, \Vert {\delta \widehat{{\mathbf w}}} \Vert
\leq \Vert {{\mathbf x}} \Vert \Bigg(
\frac{\varepsilon \, \kappa\!\left({\mathbf X}\right)}
{1 - \varepsilon \, \kappa\!\left({\mathbf X}\right)}
\Big( \frac{\Vert {{\mathbf y}} \Vert}{\lambda_{1}}
+ \Vert {\widehat{{\mathbf w}}} \Vert \Big)
+ \frac{\varepsilon \, \kappa\!\left({\mathbf X}\right)^{2}}
{\big(1 - \varepsilon \, \kappa\!\left({\mathbf X}\right)\big)^{2}}
\frac{\Vert {{\mathbf r}} \Vert}{\lambda_{1}} \Bigg) \text{.}
\end{equation}

The change of the prediction is therefore at most $\Vert {{\mathbf x}} \Vert$ times the right-hand side of the bound Eq. (1), and grows with the norm of the feature vector.

Fig. 3 shows the robustness of $\widehat{{\mathbf w}}$ against perturbations on a related but different and smaller dataset, the first days of the twelve months of $2024$, twelve of the rows of Fig. 1. Each day has the feature vector ${\mathbf x} = (1, x)^{\top}$ with $x$ its lowest temperature, and its highest temperature as the label. A hypothesis of the linear model then has intercept $w_{1}$ and slope $w_{2}$, and the feature matrix has a column of ones and the column of the lowest temperatures. Its condition number is $14.9$, and the least squares residual is $0.20$ times the norm of ${\mathbf y}$.

Two constructions of a perturbation of relative size $\varepsilon = 0.01$ are compared. The first is chosen adversarially, to make the change of the model parameters learned by linear regression as large as possible: $\delta {\mathbf X} = -\varepsilon \lambda_{1} {\mathbf v}^{(2)} \big({\mathbf u}^{(2)}\big)^{\top}$ and $\delta {\mathbf y} = \pm \varepsilon \Vert {{\mathbf y}} \Vert {\mathbf v}^{(2)}$, with the sign for which $\widehat{{\mathbf w}}$ moves more. It changes $\widehat{{\mathbf w}}$ by norm $2.23$, which is $64$ percent of the right-hand side of the bound Eq. (1). The second construction uses $100$ realizations of random $\delta {\mathbf X}$ and $\delta {\mathbf y}$. They are generated with independent and identically distributed Gaussian RVs as their entries and then normalized so that the largest singular value of $\delta {\mathbf X}$ is $\varepsilon \lambda_{1}$ and the Euclidean norm of $\delta {\mathbf y}$ is $\varepsilon \Vert {{\mathbf y}} \Vert$. The largest change among the $100$ realizations is $1.17$, which is $34$ percent of the right-hand side of the bound Eq. (1).

For a day with lowest temperature $x$ the feature vector $(1, x)^{\top}$ has norm $\sqrt{1 + x^{2}}$. By the prediction bound Eq. (2) the prediction therefore changes by at most $\sqrt{1 + x^{2}}$ times the right-hand side of the bound Eq. (1): $3.5$ degrees at $x = 0$ and $69$ degrees at $x = 20$. These are the two dotted curves in Fig. 3(a), which leave the plotted range. The hypotheses learned from the random perturbations stay within $1.6$ degrees of the learned hypothesis over the plotted range, and the hypothesis after the adversarial construction within $3.1$ degrees. Panel (b) repeats the comparison for $\varepsilon$ from $5 \cdot 10^{-4}$ to $5 \cdot 10^{-2}$. The adversarial construction reaches more than half of the right-hand side of the bound Eq. (1) up to $\varepsilon = 0.02$, the random construction less than half at every $\varepsilon$.

![Figure 3 of the entry featuremtx](https://dictionaryofml.org/terms/fig/featuremtx-fig3.svg)

*Figure 3: The perturbation bound for the features $(1, x)$, with $x$ the lowest temperature of a day, and the highest temperature as the label, over the first days of the twelve months of $2024$ (filled circles). (a) The least squares hypothesis (solid line), the hypothesis learned after the adversarial construction of a perturbation of relative size $\varepsilon = 0.01$ (dashed line), the range of the hypotheses learned after $100$ random constructions of the same relative size (gray band), and the bound Eq. (2) above and below the learned hypothesis (dotted curves, leaving the plotted range). (b) The change of the learned model parameters against $\varepsilon$: the right-hand side of the bound Eq. (1) (solid curve), the adversarial construction (dashed curve), and the largest change among $100$ random constructions (dotted curve with circles). Data generated by pythondemos/featuremtx.py*

**Synonyms:** design matrix, data matrix.

**See also:** [dataset](https://dictionaryofml.org/terms/dataset.html), [data point](https://dictionaryofml.org/terms/datapoint.html), [feature vector](https://dictionaryofml.org/terms/featurevec.html), [feature](https://dictionaryofml.org/terms/feature.html), array, [matrix](https://dictionaryofml.org/terms/matrix.html), design matrix, linear model, parameter space, [singular value decomposition](https://dictionaryofml.org/terms/svd.html), singular value, rank, condition number, sample covariance matrix, [principal component analysis](https://dictionaryofml.org/terms/pca.html), [dimensionality reduction](https://dictionaryofml.org/terms/dimred.html), column space.

## References

- Golub and Loan (2013). Matrix Computations. *The Johns Hopkins Univ. Press*.
- Hastie et al. (2009). The Elements of Statistical Learning: Data Mining, Inference, and Prediction. *Springer Science+Business Media*. <https://doi.org/10.1007/978-0-387-84858-7>
- Bishop (2006). Pattern Recognition and Machine Learning. *Springer Science+Business Media*. <https://doi.org/10.1007/978-0-387-45528-0>
- Boyd and Vandenberghe (2004). Convex Optimization. *Cambridge Univ. Press*. <https://doi.org/10.1017/CBO9780511804441>
- Bertsekas (2016). Nonlinear Programming. *Athena Scientific*.

## Cite this entry

```bibtex
@misc{dictml_featuremtx,
  author = {Jung, Alexander},
  editor = {Olioumtsevits, Konstantina and Schnoor, Ekkehard},
  title = {feature matrix},
  howpublished = {Dictionary of Applied Machine Learning (course edition)},
  year = {2026},
  doi = {10.5281/zenodo.21569296},
  note = {ISBN 978-952-64-3013-3, CC BY 4.0, retrieved 2026-11-11},
  url = {https://dictionaryofml.org/terms/featuremtx.html}
}
```