给定ℝ上无解线性方程组
A
x
=
b
boldsymbol{Ax}=boldsymbol{b}
Ax=b,构造
A
T
A
boldsymbol{A}^text{T}boldsymbol{A}
ATA及
A
T
b
boldsymbol{A}^text{T}boldsymbol{b}
ATb,然后调用博文《线性方程组的通解》定义的mySolve函数,解方程组
A
T
A
x
=
A
T
b
boldsymbol{A}^text{T}boldsymbol{Ax}=boldsymbol{A}^text{T}boldsymbol{b}
ATAx=ATb。取任一特解
x
0
boldsymbol{x}_0
x0即为解线性方程组
A
x
=
b
boldsymbol{Ax}=boldsymbol{b}
Ax=b的一个最小二乘解(即
A
boldsymbol{A}
A的列向量生成空间中距
b
boldsymbol{b}
b最近的向量)。
例1 用Python计算无解方程组
{
4
x
1
+
2
x
2
−
x
3
=
2
3
x
1
−
x
2
+
2
x
3
=
10
11
x
1
+
3
x
2
=
8
begin{cases}4x_1+2x_2-x_3=2\3x_1-x_2+2x_3=10\11x_1+3x_2quadquad=8end{cases}
⎩⎪⎨⎪⎧4x1+2x2−x3=23x1−x2+2x3=1011x1+3x2=8的最小二乘解。
import numpy as np #导入numpy
from fractions import Fraction as F #导入Fraction
np.set_printoptions(formatter= #设置输出数据格式
{'all':lambda x:str(F(x).limit_denominator())})
A=np.array([[4,2,-1], #设置系数矩阵A
[3,-1,2],
[11,3,0]],dtype='float')
b=np.array([2,10,8]) #常数项量b
B=np.matmul(A.T,A) #A的转置与A的积
c=np.matmul(A.T,b.reshape(3,1)) #A的转置与b的积
X=mySolve(B,c) #解最小二乘方程组
print(X[:,0])
程序的第5~7行设置原方程组的系数矩阵A,第8行设置原方程组的常数向量b。第9行调用numpy的matmul函数计算 A T A boldsymbol{A}^text{T}boldsymbol{A} ATA,存于B。第10行计算 A T A b boldsymbol{A}^text{T}boldsymbol{Ab} ATAb,存于c。第11行调用函数mySolve(见博文《线性方程组的通解》)解方程组 A T A x = A T b boldsymbol{A}^text{T}boldsymbol{Ax}=boldsymbol{A}^text{T}boldsymbol{b} ATAx=ATb,解集记为X。注意,X中第1列(X[:,0])存储的是方程组的特解。运行程序,输出
[9/5 -18/5 0]
即原方程组
A
x
=
b
boldsymbol{Ax}=boldsymbol{b}
Ax=b的最小二乘解为
x
0
=
(
9
5
−
18
5
0
)
boldsymbol{x}_0=begin{pmatrix}frac{9}{5}\-frac{18}{5}\0end{pmatrix}
x0=⎝⎛59−5180⎠⎞。它是
α
1
=
(
4
3
1
)
,
α
2
=
(
2
−
1
3
)
,
α
3
=
(
1
2
0
)
boldsymbol{alpha}_1=begin{pmatrix}4\3\1end{pmatrix},boldsymbol{alpha}_2=begin{pmatrix}2\-1\3end{pmatrix},boldsymbol{alpha}_3=begin{pmatrix}1\2\0end{pmatrix}
α1=⎝⎛431⎠⎞,α2=⎝⎛2−13⎠⎞,α3=⎝⎛120⎠⎞的生成子空间
W
W
W中距
b
=
(
2
10
8
)
boldsymbol{b}=begin{pmatrix}2\10\8end{pmatrix}
b=⎝⎛2108⎠⎞最近的向量。
写博不易,敬请支持:
如果阅读本文于您有所获,敬请点赞、评论、收藏,谢谢大家的支持!



