1. Introduction

The least squares method is a core approach for solving overdetermined linear systems widely used in data fitting, signal processing, and machine learning. However, directly solving the normal equations may lead to numerical instability, especially when the matrix condition number is large [1] .

Matrix decomposition transforms the original matrix into a more manageable form, significantly improving computational efficiency and stability [2]. This paper focuses on QR decomposition, SVD decomposition, and Cholesky decomposition, analyzing their numerical properties in least squares problems and validating their applicability through experiments. Furthermore, this paper compares the performance of direct solving (without decomposition) with decomposition methods, highlighting the advantages of decomposition methods [3].

2. Theoretical background

2.1. Least squares problem (example)

The least squares solution  $$ {\mathrm{x}}^{\ast }$$  satisfies:

$$ {\mathrm{x}}^{\ast }=\mathrm{a}\mathrm{r}\mathrm{g}\mathrm{m}\mathrm{i}{\mathrm{n}}_{\mathrm{x}}\Vert \mathrm{A}\mathrm{x}-\mathrm{b}\Vert $$(1)

The analytical solution is

$$ \mathrm{x}=({{\mathrm{A}}^{\mathrm{T}}\mathrm{A})}^{-1}{\mathrm{A}}^{\mathrm{T}}\mathrm{b}$$(2)

but direct computation may amplify errors due to the squared condition number of  $$ {\mathrm{A}}^{\mathrm{T}}\mathrm{A}$$ 

2.2. Matrix decomposition methods

(1) QR Decomposition

Decompose  $$ \mathrm{A}$$  into an orthogonal matrix  $$ \mathrm{Q}$$  and an upper triangular matrix  $$ \mathrm{R}$$ :

$$ \mathrm{A}=\mathrm{Q}\mathrm{R}\stackrel{}{\Rightarrow }\mathrm{ }{\mathrm{x}}^{\ast }={\mathrm{R}}^{-1}{\mathrm{Q}}^{\mathrm{T}}\mathrm{b}$$(3)

Features: Avoids computing  $$ {\mathrm{A}}^{\mathrm{T}}\mathrm{A}$$ , preserving the condition number  $$ \mathrm{\kappa }\left(\mathrm{A}\right)$$ .

(2) SVD Decomposition

Decompose  $$ \mathrm{A}$$  into  $$ \mathrm{U}\mathrm{\Sigma }{\mathrm{V}}^{\mathrm{T}}$$ , with the least squares solution:

$$ {\mathrm{x}}^{\ast }=\mathrm{V}{\mathrm{\Sigma }}^{-1}{\mathrm{U}}^{\mathrm{T}}\mathrm{b}$$(4)

Features: Optimal stability, capable of handling rank-deficient matrices.

(3) Cholesky Decomposition

When is  $$ \mathrm{A}$$  full-rank, decompose  $$ {\mathrm{A}}^{\mathrm{T}}\mathrm{A}=\mathrm{L}{\mathrm{L}}^{\mathrm{T}}$$ , with the solution:

$$ {\mathrm{x}}^{\ast }={\mathrm{L}}^{-\mathrm{T}}{\mathrm{L}}^{-1}{\mathrm{A}}^{\mathrm{T}}\mathrm{b}$$(5)

Feature: High computational efficiency, but sensitive to ill-conditioned matrices, and can only be used in positive definite matrices

2.3. Condition number and numerical stability

(1) Condition Number

The matrix condition number

$$ \mathrm{\kappa }\left(\mathrm{A}\right)=\Vert \mathrm{A}\Vert \Vert {\mathrm{A}}^{-1}\Vert $$(6)

measures the sensitivity of the problem to the tiny perturbation of the initial values of the multiplication by matrices. For the least square problems specifically, the relatively error will satisfy the inequalities of

$$ \frac{\Vert \mathrm{\delta }\mathrm{x}\Vert }{\Vert \mathrm{x}\Vert }\le \mathrm{\kappa }\left(\mathrm{A}\right)(\frac{\Vert \mathrm{\delta }\mathrm{A}\Vert }{\Vert \mathrm{A}\Vert }+\frac{\Vert \mathrm{\delta }\mathrm{b}\Vert }{\Vert \mathrm{b}\Vert })+\mathrm{o}$$(7)

As a result, the error of the solutions can be magnified exponentially if the condition number  $$ \mathrm{\kappa }\left(\mathrm{A}\right)$$  is significantly large [4] .

For the Cholesky decomposition, direct calculation of  $$ {\mathrm{A}}^{\mathrm{T}}\mathrm{A}$$  is required, so that the condition number is magnified squarely by

$$ \mathrm{\kappa }\left({\mathrm{A}}^{\mathrm{T}}\mathrm{A}\right)=\mathrm{ }\mathrm{\kappa }{\left(\mathrm{A}\right)}^2$$(8)

leading to poor stability for Cholesky decomposition.

For the QR decomposition, the original matrix is decomposed into the multiplication of the orthogonal matrix  $$ \mathrm{Q}$$  and an upper triangular matrix  $$ \mathrm{R}$$ . Given that the orthogonal matrix  $$ \mathrm{Q}$$  has the norm equal to 0 for itself and its inverse  $$ {\mathrm{Q}}^{-1}$$ , the condition number remains the same that [5]_.

$$ \mathrm{\kappa }\left(\mathrm{A}\right)=\mathrm{\kappa }\left(\mathrm{R}\right)$$(9)

For the SVD, the orthogonal matrices act as the same functions, while the algorithm for SVD also contains the truncation of the small singular values, it can reduce the condition number implicitly and avoid the ill-conditioning problems.

(2) Floating-point Arithmetic

Due to the methods of the floating-point arithmetic, the computation of matrices will propagate errors inevitably [6].

For the Cholesky decomposition, the prerequisite is that the matrix should by positive definite, so the floating-point mechanism may damage the positive-definiteness and the occurrence of having square root of the negative values. The square root and division of calculation will also lead to the round-off error by calculating in the process of Cholesky decomposition.

$$ {\mathrm{L}}_{\mathrm{i}\mathrm{i}}=\sqrt{{\mathrm{A}}_{\mathrm{i}\mathrm{i}}-\sum _{\mathrm{k}=1}^{\mathrm{i}-1}{\mathrm{L}}_{\mathrm{i}\mathrm{k}}^2}$$(10)

For the QR decomposition, the Gram-Schmidt orthogonalization may loss orthogonality due to the round-off error by calculating

$$ {\mathrm{q}}_{\mathrm{i}}={\mathrm{a}}_{\mathrm{i}}-\sum _{\mathrm{j}=1}^{\mathrm{i}-1}\left({\mathrm{a}}_{\mathrm{i}}^{\mathrm{T}}{\mathrm{q}}_{\mathrm{j}}\right){\mathrm{q}}_{\mathrm{j}}$$(11)

Although improved methods such as Householder transformation and Givens rotation can maintain orthogonality, round - off errors still exist.

For SVD, the iteration algorithm will also produce round - off errors when calculating the singular values and singular vectors. Consequently, error propagation may lead to the loss of the matrix's orthogonality.

3. Experimental design and implementation

3.1. Experimental setup

In the numerical experiment, the matrix of a specific dimension will be generated randomly and be calculated by various decomposition approaches and the direct solving method. The final results will be represented by the errors between the accurate solution and the numerical solutions, and the computation time will also be considered as an index to exhibit the efficiency of the algorithm. Two categories of the matrices will be generated with well-conditioning being random and full-rank, while the ill-conditioned ones are constructed by adding small perturbation to approximate singularity and high condition number  $$ \mathrm{\kappa }\left(\mathrm{A}\right)$$  .

Given that the Cholesky decomposition can only be utilized in the scenario that the matrix is positive definite, two groups of matrices will be tested while the first group is the  $$ \mathrm{A}\in {\mathbb{ℝ}}^{100\times 50}$$  and the other group contains the matrices with  $$ \mathrm{A}\in {\mathbb{ℝ}}^{50\times 50}$$  and is positive definite.

Tools: Python + NumPy/SciPy.

3.2. Experimental results

3.3. Discussion

As shown in Table 1 and Table 2, from the above results, we can conclude that Cholesky decomposition has the smallest error and the highest efficiency when dealing with well - conditioned problems. However, it can only be applied to positive definite matrices under restricted conditions. Meanwhile, direct solving has a similar error but poorer stability.

For the ill-conditioned Matrix, SVD decomposition has the smallest error and the highest numerical stability, although with the longest calculating time, followed by QR decomposition which has moderate numerical stability. Cholesky decomposition and direct solving fail due to the instability in the circumstances of ill-conditioning problems.

4. Practical application case study

4.1. Structural engineering analysis: ill-conditioned stiffness matrix in bridge problems

(1) Background

In bridge structural mechanics analysis, the stiffness matrix  $$ \mathrm{K}$$  ,which is usually positive definite, may become ill-conditioned ( $$ \mathrm{\kappa }\left(\mathrm{K}\right)\gg {10}^6)$$  due to material parameter errors or redundant degrees of freedom such as the over-constraint problems. Careful selection of numerical methods is required when solving the equation  $$ \mathrm{K}\mathrm{u}=\mathrm{F}$$  to reach the most precision results without considering the limit of the efficiency of the algorithm, where  $$ \mathrm{K}$$  represents the characteristics of the stiffness of a specific engineering structure,  $$ \mathrm{u}$$  represents the displacement vector of nodes, and  $$ \mathrm{F}$$  is the load vector representing the force exerted on the nodes [7].

(2) Experimental Design

Generate Ill-conditioned Stiffness Matrix: Construct a matrix  $$ \mathrm{K}\in {\mathbb{ℝ}}^{50\times 50}$$  with a condition number  $$ \mathrm{\kappa }\left(\mathrm{K}\right)>{10}^{12}$$  by adding perturbation, Generate a load vector  $$ \mathrm{F}$$  with added Gaussian noise.

Solution Methods: Use Cholesky decomposition, QR decomposition, and SVD decomposition to solve for  $$ \mathrm{u}$$ 

Evaluation Metrics: Relative error  $$ \Vert \mathrm{K}\mathrm{u}-\mathrm{F}\Vert /\Vert \mathrm{F}\Vert $$ .

Computation time.

(3) Results

As shown in Table 3, in engineering problems, efficiency can be negligible since instantaneity is not required, while precision is a key factor to consider. Ill - conditioned problems prevent the application of Cholesky Decomposition due to its numerical instability. As a result, SVD can be the most suitable matrix decomposition for this problem.

4.2. Finance modeling: well-conditioned covariance matrix in investment combination

(1) Background

In the optimization of investment combination, the covariance matrix  $$ \mathrm{\Sigma }$$  is positive definite and well-conditioned ( $$ \mathrm{\kappa }\left(\mathrm{K}\right)<{10}^3$$ ) due to the feature of the data from the financial markets. The aim is to calculate the optimal combination of the investment with least risk or highest profits by minimize or maximize  $$ {\mathrm{\omega }}^{\mathrm{T}}\mathrm{\Sigma }\mathrm{\omega }$$ , where  $$ \mathrm{\Sigma }$$  represents for the covariance of the assets and  $$ \mathrm{\omega }$$  represents for the weight of the assets

(2) Experimental Design

Step1: Generate Well-conditioning Covariance Matrix:

Construct a matrix  $$ \mathrm{\Sigma }\in {\mathbb{ℝ}}^{100\times 100}$$  by using the data from the history yield rate in 1000 days.

Step2:Solution Methods:

Use Cholesky decomposition, QR decomposition to get the answers.

Step3: Evaluation Metrics:

The difference of the objective function and the optimal solution theoretically.

Computation time.

(3) Results

As shown in Table 4, the financial problems have a great demand in the efficiency of the algorithm in order to capture the opportunities to maximize the profits and avoid loss and need a flexible analysis of the trends of the market. Considering that the covariance matrices are usually well-conditioning and positive definite, it is an optimal choice to utilize the Cholesky Decomposition as an approach.

4.3. Signal processing: moderate condition number designing matrix in filter

(1) Background

In signal denoising, the design matrix  $$ \mathrm{A}$$  may form moderate ill-condition problems due to the redundant configuration of the sensors. The target is to compute the filter coefficient  $$ \mathrm{x}$$  requires both the numerical stability and the instantaneity. The equation containing  $$ \mathrm{A}\mathrm{x}=\mathrm{b}$$  where  $$ \mathrm{A}$$  is the design matrix and  $$ \mathrm{b}$$  is the observed signal containing the noise [8,9]_.

(2) Experimental Design

Step1: Generate Moderate-conditioning Matrix:

Construct a matrix  $$ \mathrm{A}\in {\mathbb{ℝ}}^{100\times 20}$$  by using the data to simulate the signal received by the sensors

Step2:Solution Methods:

Use QR decomposition and SVD decomposition to get the answers.

Step3: Evaluation Metrics:

The SNR improvement by the calculation.

Computation time.

(3) Results

As shown in Table 5, for signal processing problems, there are small differences in error, improvement, and computation times between the two approaches, and no numerical instability issues occur for them. As a result, both of the algorithm can be applied in solving the equations of signal processing.

5. Conclusion

In terms of numerical stability and computational efficiency, SVD decomposition performs best for ill - conditioned problems and has the highest precision in solving least - square problems. However, it has higher computational costs and is unable to obtain the answer instantaneously for high - dimensional matrices because of its high computational complexity. QR decomposition offers a balance between stability and efficiency. It can provide moderate numerical stability when facing ill - conditioned problems and also has relatively higher computational efficiency. Cholesky decomposition is the fastest for well-conditioned matrices with low computational complexity but is sensitive to condition numbers and can only be restricted in solving the positive definite matrices.

As for the limitations of direct solving, Direct solving produces large errors for ill-conditioned matrices and is not recommended due to the lack of orthogonality that can reduce the condition number and easily to propagate significant errors in the multiplication and solving inverse.

This paper puts forward application recommendations, hich can meet the requirements of high-precision : Prioritize SVD decomposition;For real-time computation: Use Cholesky decomposition for well-conditioned problems and QR decomposition for ill-conditioned problems.