MATLAB和Python线性系统解中的病态和正则化

线性系统解决方案中的病态

病态线性系统是一种线性系统,它对右侧矩阵上的系数矩阵或矢量的较小摄动做出响应,而系统解却发生了较大变化。 看到这一点,将向要考虑的两个示例提供两种小扰动。 在第一个示例中,将引入系数矩阵一个分量的轻微变化,在第二个示例中,将在右侧的向量上进行摄动。

**示例1,**在本例中,我们考虑线性系统 F x = d \boldsymbol{Fx}=\boldsymbol{d} Fx=d,其中

A = ( 5.0 1.0 1.0 1.0 5.0 1.0 1.0 1.0 5.0 ) 和    b = ( 5.0 − 3.0 5.0 ) \boldsymbol{A}=\left( \begin{matrix} 5.0& 1.0& 1.0\\ 1.0& 5.0& 1.0\\ 1.0& 1.0& 5.0\\\end{matrix} \right) 和\,\,\boldsymbol{b}=\left( \begin{array}{c} 5.0\\ -3.0\\ 5.0\\\end{array} \right) A=5.01.01.01.05.01.01.01.05.0b=5.03.05.0

为了在MATLAB中求解此线性系统,我们使用以下命令:

>> A = [5., 1., 1.; 1., 5., 1.; 1., 1., 5.]
A =
5 1 1
1 5 1
1 1 5
>> b = [5.; -3.; 5.]
b =
5
-3
5
>> x = A\b
x =
1.0000
-1.0000
1.0000
In Python:
In [1]: import numpy as np
In [2]: A = np.array([[5., 1., 1.], [1., 5., 1.], [1., 1., 5.]])
In [3]: b = np.array([[5.], [-3.], [5.]])
In [4]: x = np.linalg.solve(A, b)
In [5]: print(’x = \n’, x)
x =
[[ 1.]
[-1.]
[ 1.]]

现在,我们考虑的第一个扰动将在分量 A 11 \boldsymbol{A}_{11} A11上,其中我们将 B = A B=A B=A设置为 B 11 = A 11 + 1 0 − 4 = A 11 + 0.0001 \boldsymbol{B}_{11}=\boldsymbol{A}_{11}+10^{-4}=\boldsymbol{A}_{11}+0.0001 B11=A11+104=A11+0.0001

>> B = A
B =
5 1 1
1 5 1
1 1 5
>> B(1, 1) = B(1, 1) + 0.0001
B =
5.0001 1.0000 1.0000
1.0000 5.0000 1.0000
1.0000 1.0000 5.0000
>> y = B\b
y =
1.0000
-1.0000
1.0000
>> disp(’x - y = ’), disp(num2str(x-y))
x - y =
2.1428e-05
-3.5714e-06
-3.5714e-06
>> disp([’|| x-y||_2 = ’ num2str(norm(x-y, 2))])
|| x-y||_2 = 2.2015e-05

Python代码

正在考虑的第二个扰动将在右侧,其中 b 3 \boldsymbol{b}_3 b3将被 b 3 + 1 0 − 4 \boldsymbol{b}_3+10^{-4} b3+104替换,从而在右侧给出向量 c c c。 原始系数矩阵 A A A将保持不变,并且将求解线性系统 A z = c \boldsymbol{Az}=\boldsymbol{c} Az=c

MATLAB代码

Python代码

从示例1中可以注意到,系数矩阵或右侧向量的某些分量的较小变化会导致解的较小变化。 因此,给定的线性系统对小扰动不敏感。 对小扰动不敏感的系统称为适定系统。

**示例2,**在本例中,考虑线性系统 F x = d Fx=d Fx=d,其中

F = ( 1001 − 999 999 1 1 − 1 1000 − 1000 1000 ) 和    d = ( 1001 1 1000 ) \boldsymbol{F}=\left( \begin{matrix} 1001& -999& 999\\ 1& 1& -1\\ 1000& -1000& 1000\\\end{matrix} \right) 和\,\,\boldsymbol{d}=\left( \begin{array}{c} 1001\\ 1\\ 1000\\\end{array} \right) F=1001110009991100099911000d=100111000

该示例的目的是显示系数矩阵F或向量 d d d的一个条目的微小变化如何导致线性系统 F x = d Fx=d Fx=d的解的急剧变化。

MATLAB用于通过以下命令求解线性系统 A x = b Ax=b Ax=b

>> F = [1001, -999, 999; 1, 1, -1; 1000, -1000, 1000] ;
>> d = [1001; 1; 1000] ;
>> x = F\d
Warning: Matrix is close to singular or badly scaled. Results
may be inaccurate.
RCOND = 1.067985e-16.
x =
1
0.4147
0.4147

警告表明矩阵接近奇异,因此结果可能不准确。 线性系统的MATLAB解决方案与所注意到的精确解决方案 [ 1.0 , 1.0 , 1.0 ] T \left[ 1.0,1.0,1.0 \right] ^{\boldsymbol{T}} [1.0,1.0,1.0]T相去甚远。

在Python中,可以使用Python命令解决上述线性系统:

同样,尽管Python在找到线性系统解上比在MATLAB中具有更高的准确性,但在求解线性系统上的误差仍然不小。

同样,小的扰动将被引入系数矩阵 F F F和向量 d d d。 首先,将扰动 1 0 − 5 10^{-5} 105加到 F 11 \boldsymbol{F}_{11} F11,得到矩阵 G G G。然后,计算线性系统 G y = d Gy=d Gy=d ∥ y − x ∥ 2 \left\| \boldsymbol{y}-\boldsymbol{x} \right\| _2 yx2的解:

matlab

在Python中,将使用以下命令:

其次,向向量 d d d的第三分量添加一个较小的 1 0 − 5 10^{-5} 105以获得向量 g g g并求解线性系统 F z = g Fz=g Fz=g。 在MATLAB中,使用以下命令:

matlab

Python代码

Python不会发布一条消息来抱怨矩阵的奇异性,但是确实矩阵接近奇异。

示例2中的结果表明,系数矩阵的较小变化或右侧矢量的变化都将导致线性系统解的较大变化。 对小变化的敏感性表明是病态线性系统。

不适定系统的更多示例

在本节中,将考虑两个最小二乘近似问题。 在第一个问题中,将构建一个范德蒙矩阵以找到最适合表格数据的回归系数。 在第二个示例中,将构建希尔伯特矩阵,以找到特定度数的多项式系数来近似给定函数。

示例3在此示例中,考虑了表2.1中给出的数据。 它由九个人的九个身高和体重数据对组成。 目的是找到回归系数,使得 W = ∑ j = 0 9 α j H j \boldsymbol{W}=\sum\nolimits_{\boldsymbol{j}=0}^9{\boldsymbol{\alpha }_{\boldsymbol{j}}\boldsymbol{H}^{\boldsymbol{j}}} W=j=09αjHj是最小二乘解。 在MATLAB中,我们按如下方式计算回归系数:

python

现在,我们在高度矢量的第三个分量上呈现一个小的变化,并观察该变化将如何影响最终的回归系数。

python

条件编号和病态矩阵

将条件编号链接到与矩阵相关的特征值

不适定系统的进一步分析

线性系统解的正则化

截断的SVD方法

Tikhonov正则化

L曲线法

差异原则

详情请参考 - 亚图跨际

  • 0
    点赞
  • 5
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
假设我们有一个线性方程组 Ax=b,我们可以通过 L1 和 L2 正则化组合的方法求解。 首先,我们可以将问题转化为一个最小化问题: min ||Ax-b||^2 + λ1||x||1 + λ2||x||2^2 其中,λ1 和 λ2 是两个正则化参数,||x||1 和 ||x||2^2 分别表示 L1 和 L2 正则化项。这个问题可以通过坐标下降算法求解。 下面是 MATLAB 代码示例: ``` % 生成数据 n = 100; % 变量数 m = 50; % 方程数 A = rand(m,n); % 系数矩阵 b = rand(m,1); % 右侧向量 % 求解线性方程组 x0 = rand(n,1); % 初始解 lambda1 = 0.01; % L1 正则化参数 lambda2 = 0.1; % L2 正则化参数 max_iter = 1000; % 最大迭代次数 tol = 1e-6; % 收敛精度 x = l1l2_solve(A,b,x0,lambda1,lambda2,max_iter,tol); % 输出结果 disp(x); % 定义 L1 和 L2 正则化组合求解函数 function x = l1l2_solve(A,b,x0,lambda1,lambda2,max_iter,tol) n = length(x0); x = x0; for iter=1:max_iter for i=1:n % 按照坐标轴顺序更新变量 x(i) = l1l2_shrinkage(A,b,x,lambda1,lambda2,i); end % 判断是否收敛 if norm(A*x-b) < tol break; end end end % 定义 L1 和 L2 正则化项收缩函数 function y = l1l2_shrinkage(A,b,x,lambda1,lambda2,i) % 计算梯度和 Hessian 矩阵 [G,H] = l1l2_grad_hess(A,b,x,i); % 计算收缩系数 if lambda1 == 0 alpha = -1/H; elseif lambda2 == 0 alpha = -G/(H+eps); else alpha = max((abs(G)-lambda1)/((1+2*lambda2)*H),0); end % 应用收缩操作 y = sign(G)*max(abs(G)-alpha*lambda1,0)/(H+alpha*lambda2); end % 定义 L1 和 L2 正则化项的梯度和 Hessian 矩阵计算函数 function [G,H] = l1l2_grad_hess(A,b,x,i) G = 2*sum(A(:,i).*(A*x-b)); % 梯度 H = 2*sum(A(:,i).^2); % Hessian 矩阵 end ``` 上述代码中,我们首先生成 100 个变量和 50 个方程的随机线性方程组,然后使用 L1 和 L2 正则化组合的坐标下降算法求解。 其中,l1l2_solve 函数用于求解线性方程组,l1l2_shrinkage 函数用于进行 L1 和 L2 正则化项的收紧操作,l1l2_grad_hess 函数用于计算梯度和 Hessian 矩阵。在收紧操作中,我们使用了 LARS 算法中的步长计算方法,详见《The Elements of Statistical Learning》一书。

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值