精华内容
下载资源
问答
  • 本篇解析放大器共模抑制比参数定义与其影响的评估方法,以及结合一个实际 讨论影响电路共模抑制的因素。 在讨论共模抑制比之前,先认识两个专有名词,差模增益Ad、共模增益Ac。 如图2.42(a),差模增益定义为...
  • \qquad 现在我给大家介绍一种不需要知道内部原理也可以调整PID参数的方法,而且调整的思路是所需即所得——你需要的参数是什么样子的,这个算法就会往这个方向调整PID参数。如果这个系统是多输入多输出的,我需要让...


    阅读前必看:

    1. 本代码基于MATLAB2017a版本,如果版本不同可能会报错
    2. 请从set_para.m文件开始运行,其他M文件(+下载的资源包里面的slx文件)放在一个文件夹
    3. 每次迭代都会打印出来,如果运行时间过长怀疑程序死机可以观察迭代次数是否有变化
    4. 最后会输出3幅图(figure)
    5. slx文件下载链接:https://download.csdn.net/download/weixin_44044411/11609837

    混合粒子群算法链接:https://blog.csdn.net/weixin_44044411/article/details/103638611?spm=1001.2014.3001.5502
    MATALB的自带的粒子群优化器来了:particleswarm

    0.背景

    \qquad 调整PID参数时,我们常常需要按照各种方法查各种表,还要根据建立的数学模型和之前的结果做人工调整。这个过程常常是漫长的,而且最终得到的很可能也不是最优的参数。若系统的参数有改动,又要人工再调一次。而3个参数6个方向的调整,如果输入的变量很多,那么调节的可能性组合就更多了,这个工作常常不是非专业人士所能胜任的。
    \qquad 现在我给大家介绍一种不需要知道内部原理也可以调整PID参数的方法,而且调整的思路是所需即所得——你需要的参数是什么样子的,这个算法就会往这个方向调整PID参数。如果这个系统是多输入多输出的,我需要让A输出的超调小一些,而B输出的调节时间可以尽量大一些,只需要调整我们的评价函数就可以让算法自动调参以满足客户需求了。那这样的算法使用起来是不是很容易呢?

    1.粒子群算法

    1.1.算法简介

    \qquad 粒子群算法(PSO)于1995年提出,和遗传算法一样,也是一种群体迭代算法,不同的是,粒子群算法需要整定的参数更少,不存在交叉和变异过程,所以收敛速度更快。
    \qquad 由于是群体迭代算法,因此该算法一般做成并行的。如果做成串行的,也是可以运行的,就是会随着种群规模的增大越来越慢。对于Python运用熟练的朋友可以利用Python做成并行的,对于不熟悉并行算法的朋友,直接使用MATLAB就好了,因为MATLAB的矩阵运算已经集成了并行算法。
    \qquad 为了让大家理解,我们做一个很形象的比喻。将我们的粒子群比作一群可以互相通讯的鸟,每只鸟不知道食物的具体位置,但是可以根据嗅觉知道离食物的距离和自己在地图中的位置。种群中的任意2只鸟均可以知道对方离食物的距离和对方的具体位置。假设每只鸟的飞行都是有惯性的,而每只鸟都会根据自己距离食物的距离和位置及其他鸟距离食物的距离和位置选择下一次的飞行方向。最终,所有鸟都会在信息交换的基础上,汇聚在食物地点。
    \qquad 这只是对粒子群算法的一个形象的比喻,不必过多计较该模型的现实性,闲言少叙,让我们进入 算 法 步 骤 \color{red}算法步骤

    1.2.算法步骤

    【step1】确定参数维度 N \color{blue}{N} N、参数范围(即每只鸟的初始位置),确定惯性系数 c 1 , c 2 , w c_1,c_2,w c1,c2,w,确定种群规模m,迭代次数n。
    每个粒子的初始位置是随机的,设输入向量 x = ( x 1 , x 2 , . . . , x N ) T x=(x_1,x_2,...,x_N)^T x=(x1,x2,...,xN)T必须满足参数范围限制为:
    { x m i n ( 1 ) < x 1 < x m a x ( 1 ) x m i n ( 2 ) < x 2 < x m a x ( 2 ) . . . x m i n ( N ) < x 1 < x m a x ( N ) \begin{cases}x_{min}^{(1)}<x_1<x_{max}^{(1)} \\x_{min}^{(2)}<x_2<x_{max}^{(2)} \\... \\x_{min}^{(N)}<x_1<x_{max}^{(N)} \end{cases} xmin(1)<x1<xmax(1)xmin(2)<x2<xmax(2)...xmin(N)<x1<xmax(N)
    X m i n = ( x m i n ( 1 ) , x m i n ( 2 ) , . . . x m i n ( N ) ) , X m a x = ( x m a x ( 1 ) , x m a x ( 2 ) , . . . x m a x ( N ) ) \color{blue}X_{min}=(x_{min}^{(1)},x_{min}^{(2)},...x_{min}^{(N)}),X_{max}=(x_{max}^{(1)},x_{max}^{(2)},...x_{max}^{(N)}) Xmin=(xmin(1),xmin(2),...xmin(N)),Xmax=(xmax(1),xmax(2),...xmax(N))
    则输入向量 x x x应满足 X m i n < x < X m a x X_{min}<x<X_{max} Xmin<x<Xmax
    每个粒子的初速度定为0,即 v 0 = 0 v_0=0 v0=0,第 j j j个粒子( 1 ≤ j ≤ m 1≤j≤m 1jm)的下一次迭代的速度 v ( j ) v^{(j)} v(j)由三部分组成:
    v ( j ) = w ⋅ v 0 + c 1 ⋅ r a n d ⋅ ( P ( j ) − X ( j ) ) + c 2 ⋅ r a n d ⋅ ( P G − X ( j ) ) v^{(j)}=w\cdot v_0+c_1\cdot rand\cdot (P^{(j)}-X^{(j)})+c_2\cdot rand\cdot(P_G-X^{(j)}) v(j)=wv0+c1rand(P(j)X(j))+c2rand(PGX(j))
    注:rand是(0,1)的随机数, v 0 v_0 v0代表上一次粒子的速度。
    第一部分为自身惯性因子,因为下一次的迭代次数保留了上一次的速度信息;
    第二个部分为自身最优因子, P ( j ) \color{blue}{P^{(j)}} P(j)为第 j j j个因子在之前所有迭代次数中自适应度最高的位置,可以理解为历史上自身离食物最近的位置。
    第三部分为社会因子, P G \color{blue}{P_G} PG为种群在之前所有迭代次数中自适应度最高的位置,可以理解为历史上所有粒子离食物最近的位置中的最近的一个。
    一般情况下,取 w = 1 , c 1 = c 2 = 2 w=1,c_1=c_2=2 w=1,c1=c2=2,当种群规模较大时,可以让 w w w的值随迭代次数减小以增加收敛速度。
    【step2】按照step1的速度公式计算下一次的速度,并计算下一次的粒子位置。对于第 j j j个粒子,第 k + 1 k+1 k+1次迭代(第 k + 1 k+1 k+1代)的位置 X k + 1 ( j ) \color{blue}{X_{k+1}^{(j)}} Xk+1(j)与第 k k k次迭代的位置 X K ( j ) \color{blue}{X_K^{(j)}} XK(j)与速度 v k ( k + 1 ) \color{blue}{v_k^{(k+1)}} vk(k+1)关系为:
    X k + 1 ( j ) = X k ( j ) + v k ( j + 1 ) ⋅ d t X^{(j)}_{k+1}=X^{(j)}_{k}+v_k^{(j+1)}\cdot dt Xk+1(j)=Xk(j)+vk(j+1)dt
    其中 d t \color{blue}{dt} dt是仿真间隔,一般情况下可以取1,如果希望仿真得慢一点,搜索平滑一点,可以适当减小 d t \color{blue}{dt} dt
    【step3】计算每个粒子的自适应度 F k + 1 ( j ) \color{blue}{F^{(j)}_{k+1}} Fk+1(j),为了计算出step1公式中的 P ( j ) 和 P G \color{blue}{P^{(j)}和P_G} P(j)PG。为了方便起见,我们记前k次计算得到了的 P G P_G PG P G ( k ) P_G^{(k)} PG(k),则第k+1次迭代中我们将适应度最高的粒子位置记为 P G ( k + 1 ) P_G^{(k+1)} PG(k+1),则最终的 P G P_G PG为:
    P G = { P G ( k ) , F ( P G ( k ) ) > F ( P G ( k + 1 ) ) P G ( k + 1 ) , F ( P G ( k ) ) < F ( P G ( k + 1 ) ) P_G=\begin{cases}P_G^{(k)},\qquad F(P_G^{(k)})>F(P_G^{(k+1)}) \\[2ex]P_G^{(k+1)},\quad F(P_G^{(k)})<F(P_G^{(k+1)}) \end{cases} PG=PG(k)F(PG(k))>F(PG(k+1))PG(k+1)F(PG(k))<F(PG(k+1))
    同样,我们记前k次计算得到的第 j j j个粒子的位置为 P k ( j ) P^{(j)}_{k} Pk(j),第k+1次计算得到的第 j j j个粒子的位置为 P k + 1 ( j ) P^{(j)}_{k+1} Pk+1(j),则最终的第 j j j个粒子的历史最优解 P ( j ) P^{(j)} P(j)为:
    P ( j ) = { P k ( j ) , F ( P k ( j ) ) > F ( P k + 1 ( j ) ) P k + 1 ( j ) , F ( P k ( j ) ) < F ( P k + 1 ( j ) ) P^{(j)}=\begin{cases}P_{k}^{(j)},\quad F(P_{k}^{(j)})>F(P_{k+1}^{(j)}) \\[2ex]P_{k+1}^{(j)},\quad F(P_{k}^{(j)})<F(P_{k+1}^{(j)}) \end{cases} P(j)=Pk(j)F(Pk(j))>F(Pk+1(j))Pk+1(j)F(Pk(j))<F(Pk+1(j))
    【step4】更新每个粒子的信息。
    由于我们在step2的位置迭代中已经更新过粒子的位置信息,在step1中的速度迭代中已经更新过粒子的速度信息,而在step3中又更新了每个粒子的历史最优位置 P ( j ) P^{(j)} P(j)及种群最优位置 P G P_G PG,因此实际上如果仅需要知道最优解的情况下我们不需要做这一步。
    但如果需要作出参数随迭代次数的改变的图的话,每次迭代产生最优解 P G ( k ) P_G^{(k)} PG(k)及最优适应度 F ( P G ( k ) ) F(P_G^{(k)}) F(PG(k))需要用数组保存下来。

    1.3.算法举例

    \qquad 根据以上的算法步骤,本人编写了一段粒子群算法的MATLAB代码,这个代码很简单,是无约束问题下的全局粒子群算法。全局粒子群算法收敛速度快,占用内存小,但是容易陷入局部最优解,每次迭代的最优解可能存在细微差异,当然如果这个无约束问题本身是单解的话,这个差异会随着种群规模的扩大和迭代次数增加而逐渐消除。

    function [Pg,fmin]=PSO(f,dimension,n,m,xmax,xmin,vmax,vmin)
    %全局粒子群算法,f为目标函数,dimension为维度,n为代数,m为种群规模
    %位置限幅为[xmin,xmax],速度限幅为[vmin,vmax]
        Savef=zeros(n+1,1);
        SaveData=zeros(m,dimension,10);%记录11代种群的位置
        w=1;c1=2;c2=2;%速度惯性系数
        dt=0.3;%位移仿真间隔
        v=zeros(m,dimension);%初始速度为0
        X=(xmax-xmin)*rand(m,dimension)+xmin;%初始位置满足(-xmax,xmax)内均匀分布
        P=X;%P为每个粒子每代的最优位置
        last_f=f(X);
        [fmin,min_i]=min(last_f);%Pg为所有代中的最优位置 
        Pg=X(min_i,:);
        Savef(1)=fmin;
        N=0;
        for i=1:n
            v=w*v+c1*rand*(P-X)+c2*rand*(ones(m,1)*Pg-X);
            v=(v>vmax).*vmax+(v>=vmin & v<=vmax).*v+(v<vmin).*vmin;
            X=X+v*dt;
            X=(X>xmax).*xmax+(X>=xmin & X<=xmax).*X+(X<xmin).*xmin;
            new_f=f(X);%新的目标函数值
            update_j=find(new_f<last_f);
            P(update_j,:)=X(update_j,:);%修正每个粒子的历史最优值
            [new_fmin,min_i]=min(new_f);
            new_Pg=X(min_i,:);
            Pg=(new_fmin<fmin)*new_Pg+(new_fmin>=fmin)*Pg;
            last_f=new_f;%保存当前的函数值
            fmin=min(new_fmin,fmin);%更新函数最小值
             Savef(i)=fmin;
             if mod(i,floor(n/10))==0%10代记录一次种群参数
                 N=N+1;
                SaveData(:,:,N)=X;
             end
    %         w=w-i/n*0.7*w;
        end
        for j=1:10
            figure(j)
            plot(SaveData(:,1,j),SaveData(:,2,j),'o');
            xlim([-xmax,xmax])
            ylim([-xmax,xmax])
            axis tight
        end
        figure
        plot(Savef','b-')
        disp(Pg)
        disp(fmin)
    end
    

    这是一个粒子群算法的函数,可以求函数在一定范围内的极小值,下面我们来定义一个测试函数测试一下:
    f = ( 1 − x 1 ) 2 + ( x 1 2 − 2 x 2 ) 2 f=(1-x_1)^2+(x_1^2-2x_2)^2 f=(1x1)2+(x122x2)2
    我们很容易知道这个函数的极小值为 x ∗ = ( 1 , 0.5 ) T x^*=(1,0.5)^T x=(1,0.5)T
    我们将X和Y方向上的限制都设定为[-40,40],种群规模200,仿真代数100次

    f=@(x)(1-x(:,1)).^2+(x(:,1).^2-2*x(:,2)).^2;
    [Pg,fmin]=PSO(f,2,100,200,40,-40,40,-40)
    

    运行结果为:
    Alt

    序号粒子群位置随迭代次数的变化
    1Alt
    2在这里插入图片描述
    3在这里插入图片描述
    4在这里插入图片描述
    5在这里插入图片描述
    6在这里插入图片描述

    \qquad 读者们可以发现粒子群算法在迭代的过程中部分粒子甚至有反向运动的趋势,这是由于粒子的惯性速度过大导致的,在迭代点附近,惯性因子过大,会导致很长时间无法收敛到一个点,同时也很难找到相比上一次迭代更优的解。
    \qquad 我们只需要随着迭代次数增加减少我们的惯性因子 w w w即可,在这里作者给出的表达式为
    w k + 1 = w k − i / n ∗ 0.7 w k w_{k+1}=w_k-i/n*0.7w_k wk+1=wki/n0.7wk
    i i i为当前迭代代数, n n n为总代数。
    调用同样的测试函数,种群规模、迭代次数及其他参数不变,我们得到结果如下:
    在这里插入图片描述
    和上一次对比我们可以发现fmin小了很多,即在相同的计算规模下我们的最优解更加接近于极小值,可谓一个突破性的改变,我们再看看粒子随迭代次数的分布:

    序号粒子群位置随迭代次数的变化
    1在这里插入图片描述
    2在这里插入图片描述
    3在这里插入图片描述

    \qquad 可以发现粒子群的分布相比之前稳定得多,几乎全部收敛到了中心点,这便是减小惯性因子的作用。

    2.PID自整定

    2.1.基于M文件编写的PID参数自整定

    \qquad 利用M文件可以制作一个简单的PID引擎,用于测试我们的PSO算法是否有效,书写简单PID引擎的代码可以参考如下链接:
    【MATLAB】仿真PID控制
    \qquad 我们仅需要在此基础上做一些改变并做成并行算法即可,当然,我们还需要设置一个自适应度函数,在自适应度函数的作用下,PSO算法会选择到我们满意的参数。在调整PID参数时,我们通常会关注响应曲线的超调、上升时间、调节时间、峰值时间等等。利用这些参数组合为评价函数,让评价函数越小,仿真得到的参数越好,以此作为解的好坏的对比依据(相当于鸟离食物的距离)为了简便起见,我们只选取了调节时间 t s t_s ts和超调量 σ \sigma σ作为我们的比较依据,我们将评价函数设定为:
    f a c c e s s = l n ( t s 5 ⋅ 1 0 − 2 + 1 ) + l n ( σ 1 0 − 4 + 1 ) f_{access}=ln(\frac{t_s}{5\cdot 10^{-2}}+1)+ln(\frac{\sigma}{10^{-4}}+1) faccess=ln(5102ts+1)+ln(104σ+1)
    如果对该函数求导就可以发现:
    ∂ f ∂ t s ∣ t s = 5 ⋅ 1 0 − 2 = 1 2 ⋅ 5 ⋅ 1 0 − 2 \frac{\partial f}{\partial t_s}_{|t_s=5\cdot 10^{-2}}=\frac{1}{2}\cdot 5\cdot 10^{-2} tsfts=5102=215102
    ∂ f ∂ σ ∣ σ = 1 0 − 4 = 1 2 ⋅ 1 0 − 4 \frac{\partial f}{\partial \sigma}_{|\sigma=10^{-4}}=\frac{1}{2}\cdot 10^{-4} σfσ=104=21104
    \qquad 因此 f a c c e s s f_{access} faccess t s = 5 ⋅ 1 0 − 2 t_s=5\cdot 10^{-2} ts=5102的调整率为 2.5 × 1 0 − 2 2.5×10^{-2} 2.5×102,在 σ = 1 0 − 4 \sigma=10^{-4} σ=104时的调整率为 0.5 × 1 0 − 4 0.5×10^{-4} 0.5×104如果将期望参数定为 t s ∗ = 5 ⋅ 1 0 − 2 , σ ∗ = 1 0 − 4 t_s^*=5\cdot 10^{-2},\sigma ^*=10^{-4} ts=5102,σ=104则评价函数在小于任意一个参数的情况下就会下降,但是下降的速度会越来越大,而大于任意一个参数就会上升,但上升速率会越来越慢。由于两个参数均为正数,因此下降的速率不会大于 ∂ f ∂ t s ∣ t s = 0 = 1 \frac{\partial f}{\partial t_s}_{|t_s=0}=1 tsfts=0=1 ∂ f ∂ σ ∣ σ = 0 = 1 \frac{\partial f}{\partial \sigma}_{|\sigma=0}=1 σfσ=0=1,而上升的速率可以无限小。若将评价函数小视为自适应度高,则调整过后的系统参数总是在期望值附近调整,并尽量满足$所有参数都小于期望值的情况。因为在期望值附近上升速率是越来越慢而下降速率是越来越快的,而参数最小只能取到0。因此,有多个参数大于期望值而有少个参数小于期望值的和会比有少个参数大于期望值而多个参数小于期望值的和要大得多。
    \qquad 写好评价函数之后,我们可以设自适应度函数 F = − f a c c e s s F=-f_{access} F=faccess,因为我们的评价函数是越小越好的,而自适应度是越高越好。
    \qquad 根据1.2的算法步骤我们可以编写出粒子群算法,而PID并行仿真引擎的编写需要参考上文提供的连接。(simulink的PID仿真请参考下一小节2.2.)粒子群算法每迭代出一代的参数(PSO.m),就把这一代里面最好的参数 P G P_G PG返还给PID引擎进行仿真(PID_sim.m),再编写一个参数自识别的程序(parameter_cal.m),将需要识别的指标(调节时间、超调)计算出来,传递给 f a c c e s s f_{access} faccess进行计算(PID_access.m),将计算得到的参数

    【PID_sim.m】(PID并行仿真引擎)

    function [f_infty,tp,ts,sigma]=PID_sim(kp,ki,kd,debug)
        %kp,ki,kd为PID参数,T0为采样时间,total_t为仿真时间
        Tt=5e-4;%仿真采样周期
        T0=1e-2;%控制采样周期
        Tf=1e-3;%微分环节的滤波器系数
        alp=Tf/(Tt+Tf);
        total_t=1;
        N=floor(total_t/Tt);%样本总数
        M=numel(kp);%仿真序列数
        kp=reshape(kp,M,1);
        ki=reshape(ki,M,1);
        kd=reshape(kd,M,1);
        sys=tf(0.998,[0.021,1]);   %建立被控对象传递函数,即式4.1
        dsys=c2d(sys,Tt,'z');      %离散化
        [num,den]=tfdata(dsys,'v');   %
        e_1=zeros(M,1);      %前一时刻的偏差      
        Ee=zeros(M,1);       %累积偏差
        u_1=zeros(M,1);    %前一时刻的控制量
        y_1=zeros(M,1);       %前一时刻的输出
        ud_1=zeros(M,1);     %前一时刻的微分输出
        %预分配内存
        r=1500;
        y=zeros(M,N);
        u=zeros(M,N);
        for k=1:1:N
            y(:,k)=-1*den(2).*y_1+num(2)*u_1+num(1).*u(:,k);%系统响应输出序列
            e0=r-y(:,k);   %误差信号
            ud=alp.*ud_1+(1-alp)/T0.*kd.*(e0-e_1);
            u(:,k)=kp.*e0+T0*ki.*Ee+ud; %系统PID控制器输出序列
            Ee=Ee+e0;    %误差的累加和
            u_1=u(:,k);    	%前一个的控制器输出值
            y_1=y(:,k);    	%前一个的系统响应输出值
            e_1=e0;		%前一个误差信号的值
            ud_1=ud;
        end
        if debug %非调试模式下不显示也不打印图像
            plot(linspace(0,total_t,N),y)
        end
        [f_infty,tp,ts,sigma]=parameter_cal(y,Tt,2e-2,debug);%取得阶跃响应指标
    end
    

    【parameter_cal.m】(计算仿真曲线参数)

    function [f_infty,tp,ts,sigma]=parameter_cal(y,T0,delt_err,debug)
    %y是系统输出序列
    %T0是采样时间,可以结合时间序列点序号算出实际时间
    %delt_err是调节时间的误差精度
        M=size(y,1);N=size(y,2);%M为运算维度,N为时间序列长度
        f_infty=y(:,N);%稳态值序列
        err=y-f_infty*ones(1,N);%通过稳态值计算误差序列
        ferr=fliplr(abs(err));%倒序并取绝对值
        [~,ts_i]=max(ferr>delt_err*f_infty,[],2);
        ts_i=N*ones(M,1)-ts_i;
        ts=ts_i*T0;%调节时间
        [fp,tp]=max(y,[],2);%峰值和函数峰值
        tp=tp*T0;
        tp(abs(fp-f_infty)<1e-5)=NaN; %过阻尼无超调,没有峰值时间
        sigma=(fp-f_infty)./f_infty;
        if debug && M==1 %非调试模式下不显示
            disp(['系统稳态值为',num2str(f_infty)])
            disp(['系统超调量为',num2str(sigma*100),'%'])
            if isnan(tp)
                disp('系统不存在峰值时间')
            else
                disp(['系统峰值时间为',num2str(tp),'s'])
            end
            disp(['系统的调节时间为',num2str(ts),'s'])
        end
    end
    

    【PID_access.m】(根据参数得出评价度的函数,程序中为评价度越低越好)

    function y=PID_access(para)
    %PID调节性能的指标参数
    kp=para(:,1);ki=para(:,2);kd=para(:,3);
    [~,~,ts,sigma]=PID_sim(kp,ki,kd,false);
    y=log(ts/5e-2+1)+log(sigma/1e-2+1);
    end
    

    【PSO.m】(PSO算法优化函数)

    function [Pg,fmin]=PSO(dimension,n,m)
    %全局粒子群算法,f为目标函数,dimension为维度,n为代数,m为种群规模
        w=1;c1=2;c2=2;%速度惯性系数
        sigma_data=zeros(1,n);
        ts_data=zeros(1,n);
        dt=0.3;%位移仿真间隔
        vmax=1;%速度限幅
        xmax=[15,50,2];%位置限幅
        xmin=[0.2,0,0];
        v=zeros(m,dimension);%初始速度为0
        X=(xmax-xmin).*rand(m,dimension)+xmin;%初始位置满足(xmin,xmax)内均匀分布
        P=X;%P为每个粒子每代的最优位置
        last_f=PID_access(X);
        [fmin,min_i]=min(last_f);%Pg为所有代中的最优位置 
        Pg=X(min_i,:);
        N=0;
        for i=1:n
            v=w*v+c1*rand*(P-X)+c2*rand*(ones(m,1)*Pg-X);
            v=(v>vmax).*vmax+(v>=-vmax & v<=vmax).*v+(v<-vmax).*(-vmax);
            X=X+v*dt;
            X=(X>xmax).*xmax+(X>=xmin & X<=xmax).*X+(X<xmin).*(xmin);
            new_f=PID_access(X);%新的目标函数值
            update_j=find(new_f<last_f);
            P(update_j,:)=X(update_j,:);%修正每个粒子的历史最优值
            [new_fmin,min_i]=min(new_f);
            new_Pg=X(min_i,:);
            Pg=(new_fmin<fmin)*new_Pg+(new_fmin>=fmin)*Pg;
            last_f=new_f;%保存当前的函数值
            fmin=min(new_fmin,fmin);%更新函数最小值
            [~,~,ts,sigma]=PID_sim(Pg(1),Pg(2),Pg(3),true)
            sigma_data(1,i)=sigma;
            ts_data(1,i)=ts;
            hold on
        end
        legend('1','2','3','4','5','6','7','8','9','10','11','12','13','14','15')
        title('响应曲线随迭代次数的变化')
        figure
        plot(ts_data)
        title('调节时间随迭代次数的变化')
        figure
        plot(sigma_data)
        title('超调量随迭代次数的变化')
    end
    

    【main.m】(主函数,从这里开始)

    dimension=3;n=20;m=200;%PID为3位参数,n是迭代次数,m为种群规模
    [Pg,fmin]=PSO(dimension,n,m)
    

    仿真结果如下:
    在这里插入图片描述
    注:响应曲线可能要在左上角局部放大才能获得该图像。
    在这里插入图片描述
    在这里插入图片描述
    \qquad 这是一个简单的一阶惯性系统的PID调参仿真,用于模拟直流电机的转速模型,可以看出直流电机的转速响应是越来越快的(调节时间越来越小),超调量在第一次就达到了0,因此系统不存在超调。

    *2.2.复杂系统的PID自整定(基于simulink仿真)

    \qquad 我们以双闭环直流调速系统为例来说明PSO对PID自整定参数的作用。直流调速系统由外环(转速环)和内环(电流环)构成,转速环是一个PI调节器,调节转速为期望转速,用于抵抗负载扰动、电机参数变化带来的扰动;电流环也是一个PI调节器,与转速环不同的是,电流调节器的目的是在于让电流快速跟随给定值,以获得更好的动态特性。我们以电机的启动过程为例仿真这个系统,直流电机启动时,电流首先达到容许最大值,电机以最大加速度进行加速,到达期望转速时让电流迅速减小并达到负载转矩所需电流,电机维持期望转速不变。

    2.2.1.PSO优化PID的过程详解

    \qquad 下面我们设定我们需要达成的调节指标——由于启动时的电机电流会迅速上升到过载电流(允许的最大值),是一个恒定值,故会有一个电流平台,记上升到电流平台产生的超调为电流超调 σ i \sigma_i σi,从超调回到电流平台所需的时间为电流的调节时间 t s i t_{si} tsi,最终转速稳定前的超调记为转速超调 σ n \sigma_n σn,转速到达稳定值区域的时间为调节时间 t s n t_{sn} tsn
    【搭建simulink模型】
    双闭环调速系统的结构框图如下:
    在这里插入图片描述
    搭建simulink仿真模型如下:
    在这里插入图片描述
    注:上图可能不够清楚,根据以下参数搭建simulink模型,并且命名为double_circle.slx(和PSO的所有程序同名)
    电机的参数如下:

    • 额定电压/电流/功率:800V/2164A 1500KW
    • 额定转速:44r/min
    • 过载系数: λ i = 1.5 \lambda_i=1.5 λi=1.5
    • 励磁:110V/191A
    • 电枢电阻 R a R_a Ra:0.1 Ω \Omega Ω
    • 主回路电阻R:0.1014 Ω \Omega Ω
    • 电枢(主回路)电感L:3.25mH
    • 变流器放大倍数 K s K_s Ks=80
    • 变流器通态压降:3V
    • 系统飞轮惯量: 8.6 × 1 0 5 N ⋅ m 2 8.6×10^5N\cdot m^2 8.6×105Nm2

    \qquad 根据工程设计方法当然可以得到一套电流和转速的PI参数,计算的相关参数为: [ K n p , K n I , K i p , K i I ] E N = [ 7.3660 , 1.5658 , 0.8257 , 30.2773 ] [K_{np},K_{nI},K_{ip},K_{iI}]_{EN}=[7.3660 ,1.5658 ,0.8257, 30.2773] [Knp,KnI,Kip,KiI]EN=[7.3660,1.5658,0.8257,30.2773]
    通过仿真,我们得到的转速和电流曲线如下:

    额定负载启动空载启动
    在这里插入图片描述在这里插入图片描述

    调用我们的参数计算程序计算一下参数可以得到: σ n = 4.07 × 1 0 − 4 % , t s n = 0.7875 , σ i = 20.38 % , t s i = 0.0172 \sigma_n=4.07×10^{-4}\%,t_{sn}=0.7875,\sigma_i=20.38\%,t_{si}=0.0172 σn=4.07×104%,tsn=0.7875,σi=20.38%,tsi=0.0172
    \qquad 这个参数当然是比较理想的参数,如果我们不知道工程设计方法或者说工程设计方法不适用于我们的模型时,那么就可以使用PSO整定PID参数的方法。
    【修改simulink模型以适应PSO的M文件调用】
    \qquad 这里需要调整的参数为转速调节器ASR的 K n p , K n I K_{np},K_{nI} KnpKnI和电流调节器ACR的 K i p , K i I K_{ip},K_{iI} Kip,KiI,如果我们希望能实时调整simulink的参数进行动态仿真,我们需要直接将图中的几个比例环节依次设为符号( K n p , K n I , K i p , K i I K_{np},K_{nI},K_{ip},K_{iI} Knp,KnI,Kip,KiI),并在工作区实时修改这4个参数再进行仿真。如果我们要将实时调整PI参数的代码写在M文件中仿真,必须要将变量声明为全局(global),此时变量会在工作区里面,simulink就可以调用了。另存为PSO_double_circle.slx
    【编写辅助的M文件】
    cal_n.m:单组PID参数仿真结果的参数计算

    function [n_sigma,ts_n,I_sigma,ts_transient_i,n_err]=cal_n(parameters)
        global Kpn KIn Kpi KIi
        %将simulink所需参数导入工作区,以便仿真使用
        obj_slx='PSO_double_circle.slx';
        Kpn=parameters(1);
        KIn=parameters(2);
        Kpi=parameters(3);
        KIi=parameters(4);
        sim(obj_slx);
        I_data=I.signals.values;n_data=n.signals.values;time=I.time;
        n_infty=mean(n_data(numel(n_data)-10:numel(n_data)));
        n_err=abs(n_infty-10/0.227);%期望与输出的转速误差
        [i_max,tr_index]=max(I_data);
        diffI=[0;diff(I_data)];%电流微分
        pos=find(abs(diffI)<1);%找到电流波动较小的部分
        pos=pos(pos>tr_index);%排除启动时电流不变的情况
        min_pos=min(pos);
        I_transient=mean(I_data(min_pos+10:min_pos+20));%恒流升速阶段的暂稳态电流值
        I_sigma=(i_max-I_transient)/I_transient;
        ts_transient_i=time(min_pos-tr_index);%电流暂稳态的调节时间
        n_max=max(n_data);
        n_sigma=(n_max-n_infty)/n_infty;
        ts_n_index=numel(n_data)-find(abs(flip(n_data)-n_infty)/n_infty>2e-2,1);
        ts_n=time(ts_n_index);
    end
    

    PID_n_access.m:PID参数的仿真,以代为单位的

    function [y,ts_n,n_sigma,ts_transient_i,I_sigma]=PID_n_access(parameters_list)
        M=size(parameters_list,1);%计算需要仿真的组数(simulink无法并行仿真)
        ts_n0=0.8;n_sigma0=1e-2;I_sigma0=0.02;ts_i0=1e-4;%期望值
        %将simulink所需参数导入工作区,以便仿真使用
        y=zeros(M,1);
        n_sigma=zeros(M,1);ts_n=zeros(M,1);
        for sim_i=1:M
            parameters=parameters_list(sim_i,:);
            [n_sigma,ts_n,I_sigma,ts_transient_i,n_err]=cal_n(parameters);
            y(sim_i)=log(ts_n./ts_n0+1)+log(n_sigma./n_sigma0+1)+log(I_sigma/I_sigma0+1)+log(ts_transient_i/ts_i0+1)+exp(n_err*10);
        end
    end
    

    para_cal.m:参数计算,以代为单位的

    obj_slx='double_circle.slx';
    sim(obj_slx);
    I_data=I.signals.values;n_data=n.signals.values;time=I.time;
    I_infty=mean(I_data(numel(I_data)-10:numel(I_data)));
    n_infty=mean(n_data(numel(n_data)-10:numel(n_data)));
    [i_max,tr_index]=max(I_data);
    diffI=[0;diff(I_data)];%电流微分
    pos=find(abs(diffI)<1);%找到电流波动较小的部分
    pos=pos(pos>tr_index);%排除启动时电流不变的情况
    min_pos=min(pos);
    I_transient=mean(I_data(min_pos+10:min_pos+20));%恒流升速阶段的暂稳态电流值
    I_sigma=(i_max-I_transient)/I_transient;
    if I_sigma<1e-4
        tr_i=NaN;    
    else
        tr_i=time(tr_index);%电流的上升时间
    end
    ts_transient_i=time(min_pos-tr_index);%电流暂稳态的调节时间
    ts_i_index=numel(I_data)-find(abs(flip(I_data)-I_infty)/I_infty>2e-2,1);
    ts_i=time(ts_i_index);%转速稳态时电流的调节时间
    [n_max,tr_index]=max(n_data);
    n_sigma=(n_max-n_infty)/n_infty;
    %参数展示
    %转速参数展示
    if n_sigma<1e-4
        tr_n=NaN;
    else
        tr_n=time(tr_index);
    end
    ts_n_index=numel(n_data)-find(abs(flip(n_data)-n_infty)/n_infty>2e-2,1);
    ts_n=time(ts_n_index);
    if isnan(tr_n)
        disp('转速无超调')
    else
        disp(['转速上升时间为',num2str(tr_n),'s'])
        disp(['转速超调量为',num2str(100*n_sigma),'%'])
    end
    disp(['转速调节时间为',num2str(ts_n),'s'])
    %电流参数展示
    if isnan(tr_i)
        disp('电流无超调')
    else
        disp(['电流上升时间为',num2str(tr_i),'s'])
        disp(['电流超调量为',num2str(100*I_sigma),'%'])
        disp(['电流暂稳态调节时间为',num2str(ts_transient_i),'s'])
    end
    disp(['电流的最终调节时间为',num2str(ts_i),'s'])
    

    draw.m:仅根据一组PID参数画出转速波形图

    function draw(parameters)
        global Kpn Kin KpI KiI
        Kpn=parameters(1);
        Kin=parameters(2);
        KpI=parameters(3);
        KiI=parameters(4);
        sim('PSO_double_circle')
        n_data=n.signals.values;time=n.time;%获取转速序列及对应时刻
        plot(time,n_data)
    

    draw_n_I.m:根据一组PID参数绘制电流和转速曲线

    function draw_n_I(parameters)
        global Kpn KIn Kpi KIi
        Kpn=parameters(1);
        KIn=parameters(2);
        Kpi=parameters(3);
        KIi=parameters(4);
        sim('PSO_double_circle')
        n_data=n.signals.values;time=n.time;%获取转速序列及对应时刻
        I_data=I.signals.values;
        figure
        subplot(211)
        plot(time,I_data,'r-','LineWidth',1.5)
        title('电流波形')
        subplot(212)
        plot(time,n_data,'b-','LineWidth',1.5)
        title('转速波形')
    end
    

    PSO.m:PSO优化程序

    function [Pg,fmin]=PSO(n,m,xmax,xmin,cal_f,f)
    %全局粒子群算法,f为目标函数,dimension为维度,n为代数,m为种群规模
        dimension=numel(xmax);
        w=1;c1=2;c2=2;%速度惯性系数
        sigma_data=zeros(1,n);
        ts_data=zeros(1,n);
        dt=0.3;%位移仿真间隔
        vmax=0.1*(xmax-xmin);%速度限幅
        v=zeros(m,dimension);%初始速度为0
        X=(xmax-xmin).*rand(m,dimension)+xmin;%初始位置满足(xmin,xmax)内均匀分布
        P=X;%P为每个粒子每代的最优位置
        last_f=f(X);
        [fmin,min_i]=min(last_f);%Pg为所有代中的最优位置 
        Pg=X(min_i,:);
        last_Pg=Pg;
        legend_str=cell(0);
        legend_i=1;
        figure(1)
        legend_str{legend_i}=num2str(1);
        draw(Pg)
        for i=1:n
            v=w*v+c1*rand*(P-X)+c2*rand*(ones(m,1)*Pg-X);
            v=(v>vmax).*vmax+(v>=-vmax & v<=vmax).*v+(v<-vmax).*(-vmax);
            X=X+v*dt;
            X=(X>xmax).*xmax+(X>=xmin & X<=xmax).*X+(X<xmin).*(xmin);
            new_f=f(X);%新的目标函数值
            update_j=find(new_f<last_f);
            P(update_j,:)=X(update_j,:);%修正每个粒子的历史最优值
            [new_fmin,min_i]=min(new_f);
            new_Pg=X(min_i,:);
            Pg=(new_fmin<fmin)*new_Pg+(new_fmin>=fmin)*Pg;
            last_f=new_f;%保存当前的函数值
            fmin=min(new_fmin,fmin);%更新函数最小值
            [sigma,ts,~,~]=cal_f(Pg);
            sigma_data(1,i)=sigma;
            ts_data(1,i)=ts;
            if last_Pg~=Pg
                legend_i=legend_i+1;
                figure(1)
                legend_str{legend_i}=num2str(i);
                draw(Pg)
                hold on
            end
            last_Pg=Pg;
            disp(['迭代次数:',num2str(i)])
        end
        figure(1)
        legend(legend_str)
        title('响应曲线随迭代次数的变化')
        figure(2)
        plot(ts_data)
        title('调节时间随迭代次数的变化')
        figure(3)
        plot(sigma_data)
        title('超调量随迭代次数的变化')
    end
    

    set_para.m(设置参数程序,从这里开始)

    clear
    best_parameters=[8.4725,97,72,1.794,56.06];
    global Kpn KIn Kpi KIi
    Kpn=0;KIn=0;Kpi=0;KIi=0;
    n=10;m=100;xmax=[20,120,10,60];xmin=[0,0,0,0];cal_f=@cal_n;f=@PID_n_access;
    [Pg,fmin]=PSO(n,m,xmax,xmin,cal_f,f);
    
    

    我们将期望的参数设定为
    σ n ∗ = 0 , 1 % , t s n ∗ = 1 , σ i ∗ = 20 % , t s i ∗ = 0.005 \sigma_n^*=0,1\%,t_{sn}^*=1,\sigma_i^*=20\%,t_{si}^*=0.005 σn=0,1%,tsn=1,σi=20%,tsi=0.005利用PSO算法构造的评价函数为
    f a c c e s s = l n ( t s n 1 + 1 ) + l n ( σ n 0.001 + 1 ) + l n ( σ i 0.2 + 1 ) + l n ( t s i 0.005 + 1 ) ; f_{access}=ln(\frac{ts_n}{1}+1)+ln(\frac{\sigma_n}{0.001}+1)+ln(\frac{\sigma_i}{0.2}+1)+ln(\frac{ts_i}{0.005}+1); faccess=ln(1tsn+1)+ln(0.001σn+1)+ln(0.2σi+1)+ln(0.005tsi+1);
    得到的最佳参数为:
    [ K n p , K n I , K i p , K i I ] P S O = [ 18.2613 , 73.0965 , 1.6843 , 16.5703 ] P S O [K_{np},K_{nI},K_{ip},K_{iI}]_{PSO}=[18.2613,73.0965 ,1.6843,16.5703]_{PSO} [Knp,KnI,Kip,KiI]PSO=[18.2613,73.0965,1.6843,16.5703]PSO
    下的仿真曲线如下:
    超调量和调节时间随迭代次数变化
    在这里插入图片描述在这里插入图片描述
    PSO与工程设计方法电机启动曲线对比(额定负载启动)

    工程设计方法PSO方法
    在这里插入图片描述在这里插入图片描述

    通过计算程序得到的参数为: σ n = 6.7556 × 1 0 − 4 % , t s n = 0.7862 , σ i = 18.83 % , t s i = 0.01341 \sigma_n=6.7556×10^{-4}\%,t_{sn}=0.7862,\sigma_i=18.83\%,t_{si}=0.01341 σn=6.7556×104%,tsn=0.7862,σi=18.83%,tsi=0.01341
    工程设计方法的参数为:
    σ n = 4.07 × 1 0 − 4 % , t s n = 0.7875 , σ i = 20.38 % , t s i = 0.0172 \sigma_n=4.07×10^{-4}\%,t_{sn}=0.7875,\sigma_i=20.38\%,t_{si}=0.0172 σn=4.07×104%,tsn=0.7875,σi=20.38%,tsi=0.0172
    \qquad 可以发现除了最后一项指标之外,其余指标均达到了期望值要求,并且和工程设计方法得出的指标相比甚至还略优于它。

    2.2.2.在PSO优化过程中修改参数价值权重

    \qquad 转速的超调和调节时间,电流的超调和调节时间,我们最后评估的指标除了转速无静差的基本要求外,这4项指标就是我们的评价标准,如果我们更希望牺牲某一个或多个参数去优化另一个或另一些参数时,该怎么调节我们的算法呢?
    \qquad 答案也非常简单,调节的方法有2个,一个是期望值,一个是评价函数的各项系数。这里我们用调节期望值的方法做一个示范:
    \qquad 之前的期望值为:
    σ n ∗ = 0.1 % , t s n ∗ = 1 , σ i ∗ = 20 % , t s i ∗ = 0.005 \sigma_n^*=0.1\%,t_{sn}^*=1,\sigma_i^*=20\%,t_{si}^*=0.005 σn=0.1%,tsn=1,σi=20%,tsi=0.005
    修改期望值为:
    σ n ∗ = 10 % , t s n ∗ = 0.8 , σ i ∗ = 2 % , t s i ∗ = 0.0001 \sigma_n^*=10\%,t_{sn}^*=0.8,\sigma_i^*=2\%,t_{si}^*=0.0001 σn=10%,tsn=0.8,σi=2%,tsi=0.0001
    重新仿真,对比二者的结果如下:

    PSO1PSO2
    在这里插入图片描述在这里插入图片描述

    \qquad 可以看出转速超调量略微增大了一些,但是转速和电流的调节时间明显缩短,电流的超调量基本不变。这一点也表明多项优化指标是相互耦合的,可以牺牲一者满足另一者,只要参数在这方面是可调的。
    \qquad 需要注意的是,simulink中无法使用串行的算法,因此代码中使用时是串行的,可能仿真比较慢,请大家耐心等待(每迭代一次都会打印一次)
    \qquad 希望本文对您有帮助,谢谢阅读!

    展开全文
  • 现实中的金融时间序列存在严重的波动集聚性,...在此基础上,利用Holder不等式和Jensen不等式及相关理论评估了带有指数为α>0定期变化分布扰动的GARCH(1,1)模型参数,最后进一步证明所构造评估的可信性和渐近正态分布性。
  • 一方面它们要足够大,才足够评估参数、模型。 另一方面,如果它们太大,则会浪费数据(验证集和训练集的数据无法用于训练)。 在 k-fold 交叉验证中:先将所有数据拆分成 k 份,然后其中 1 份...

    机器学习入门系列(2)–如何构建一个完整的机器学习项目,第十一篇!

    该系列的前 10 篇文章:

    上一篇文章介绍了性能评估标准,但如何进行模型评估呢,如何对数据集进行划分出训练集、验证集和测试集呢?如何应对可能的过拟合和欠拟合问题,还有超参数的调优,如何更好更快找到最优的参数呢?

    本文会一一介绍上述的问题和解决方法。


    2. 模型评估的方法

    2.1 泛化能力

    1. 泛化能力:指模型对未知的、新鲜的数据的预测能力,通常是根据测试误差来衡量模型的泛化能力,测试误差越小,模型能力越强;
    2. 统计理论表明:如果训练集和测试集中的样本都是独立同分布产生的,则有 模型的训练误差的期望等于模型的测试误差的期望
    3. 机器学习的“没有免费的午餐定理”表明:在所有可能的数据生成分布上,没有一个机器学习算法总是比其他的要好。
      • 该结论仅在考虑所有可能的数据分布时才成立。
      • 现实中特定任务的数据分布往往满足某类假设,从而可以设计在这类分布上效果更好的学习算法。
      • 这意味着机器学习并不需要寻找一个通用的学习算法,而是寻找一个在关心的数据分布上效果最好的算法。
    4. 正则化是对学习算法做的一个修改,这种修改趋向于降低泛化误差(而不是降低训练误差)。
      • 正则化是机器学习领域的中心问题之一。
      • 没有免费的午餐定理说明了没有最优的学习算法,因此也没有最优的正则化形式。

    2.2 泛化能力的评估

    常用的对模型泛化能力的评估方法有以下几种,主要区别就是如何划分测试集。

    • 留出法(Holdout)
    • k-fold 交叉验证(Cross Validation)
    • 留一法(Leave One Out, LOO)
    • 自助法(bootstrapping)
    2.2.1 留出法(Holdout)

    留出法是最简单也是最直接的验证方法,它就是将数据集随机划分为两个互斥的集合,即训练集和测试集,比如按照 7:3 的比例划分,70% 的数据作为训练集,30% 的数据作为测试集。也可以划分为三个互斥的集合,此时就增加一个验证集,用于调试参数和选择模型

    直接采用 sklearn 库的 train_test_split 函数即可实现,一个简单的示例代码如下,这里简单调用 knn 算法,采用 Iris 数据集。

    from sklearn.model_selection import train_test_split
    from sklearn.datasets import load_iris
    from sklearn.neighbors import KNeighborsClassifier
    
    # 加载 Iris 数据集
    dataset = load_iris()
    # 划分训练集和测试集
    (trainX, testX, trainY, testY) = train_test_split(dataset.data, dataset.target, random_state=3, test_size=0.3)
    # 建立模型
    knn = KNeighborsClassifier()
    # 训练模型
    knn.fit(trainX, trainY)
    # 将准确率打印
    print('hold_out, score:', knn.score(testX, testY))
    

    留出法的使用需要注意:

    1. 数据集的划分要尽可能保持数据分布的一致性,避免因为数据划分过程引入额外的偏差而对最终结果产生影响。比如训练、验证和测试集的类别比例差别很大,则误差估计将由于三个集合数据分布的差异而产生偏差。

      因此,分类任务中必须保持每个集合中的类别比例相似。从采样的角度看数据集的划分过程,这种保留类别比例的采样方式称为“分层采样”。

    2. 即便确定了训练、验证、测试集的比例,还是有多种划分方式,比如排序后划分、随机划分等等,这些不同的划分方式导致单次留出法得到的估计结果往往不够稳定可靠。因此,使用留出法的时候,往往采用若干次随机划分、重复进行实验后,取平均值作为最终评估结果

    分层采样的简单代码实现如下所示,主要是调用了 sklearn.model_selection 中的 StratifiedKFold

    from sklearn.datasets import load_iris
    from sklearn.neighbors import KNeighborsClassifier
    from sklearn.model_selection import StratifiedKFold
    from sklearn.base import clone
    
    def StratifiedKFold_method(n_splits=3):
        '''
        分层采样
        :return:
        '''
        # 加载 Iris 数据集
        dataset = load_iris()
        data = dataset.data
        label = dataset.target
        # 建立模型
        knn = KNeighborsClassifier()
        print('use StratifiedKFold')
        skfolds = StratifiedKFold(n_splits=n_splits, random_state=42)
        scores = 0.
        for train_index, test_index in skfolds.split(data, label):
            clone_clf = clone(knn)
            X_train_folds = data[train_index]
            y_train_folds = (label[train_index])
            X_test_fold = data[test_index]
            y_test_fold = (label[test_index])
            clone_clf.fit(X_train_folds, y_train_folds)
            y_pred = clone_clf.predict(X_test_fold)
            n_correct = sum(y_pred == y_test_fold)
            print(n_correct / len(y_pred))
            scores += n_correct / len(y_pred)
        print('mean scores:', scores / n_splits)
    

    留出法也存在以下的缺点:

    1. 在验证集或者测试集上的评估结果和划分方式有关系,这也就是为什么需要多次实验,取平均值;
    2. 我们希望评估的是在原始数据集上训练得到的模型的能力,但留出法在划分两个或者三个集合后,训练模型仅使用了原始数据集的一部分,这会降低评估结果的保真性。但这个问题没有完美的解决方法,常见做法是将大约 2/3 ~ 4/5 的样本作为训练集,剩余的作为验证集和测试集。
    2.2.2 k-fold 交叉验证(Cross Validation)

    k-fold 交叉验证 的工作流程:

    1. 将原始数据集划分为 k 个大小相等且互斥的子集;
    2. 选择 k-1 个子集作为训练集,剩余作为验证集进行模型的训练和评估,重复 k 次(每次采用不同子集作为验证集);
    3. k 次实验评估指标的平均值作为最终的评估结果。

    通常,k 取 10。

    但和留出法类似,同样存在多种划分 k 个子集的方法,所以依然需要随时使用不同方式划分 p 次,每次得到 k 个子集。

    同样,采用 sklearn.cross_validationcross_val_score 库可以快速实现 k-fold 交叉验证法,示例如下:

    from sklearn.datasets import load_iris
    from sklearn.neighbors import KNeighborsClassifier
    from sklearn.cross_validation import cross_val_score
    # 加载 Iris 数据集
    dataset = load_iris()
    data = dataset.data
    label = dataset.target
    # 建立模型
    knn = KNeighborsClassifier()
    # 使用K折交叉验证模块
    scores = cross_val_score(knn, data, label, cv=10, scoring='accuracy')
    # 将每次的预测准确率打印出
    print(scores)
    # 将预测准确平均率打印出
    print(scores.mean())
    
    2.2.3 留一法

    留一法是 k-fold 交叉验证的一个特例情况,即让 k=N, 其中 N 是原始数据集的样本数量,这样每个子集就只有一个样本,这就是留一法

    留一法的优点就是训练数据更接近原始数据集了,仅仅相差一个样本而已,通过这种方法训练的模型,几乎可以认为就是在原始数据集上训练得到的模型 。

    但缺点也比较明显,计算速度会大大降低特别是原始数据集非常大的时候,训练 N 个模型的计算量和计算时间都很大,因此一般实际应用中很少采用这种方法。

    2.2.4 自助法

    在留出法和 k-fold 交叉验证法中,由于保留了一部分样本用于测试,因此实际训练模型使用的训练集比初始数据集小,这必然会引入一些因为训练样本规模不同而导致的估计偏差

    留一法受训练样本规模变化的影响较小,但是计算复杂度太高

    自助法是一个以自助采样法(bootstrap sampling)为基础的比较好的解决方案。同时,它也是随机森林算法中用到的方法。

    它的做法就是对样本数量为 N 的数据集进行 N 次有放回的随机采样,得到一个大小是 N 的训练集。

    在这个过程中将会有一部分数据是没有被采样得到的,一个样本始终没有被采样出来的概率是 ( 1 − 1 N ) N (1-\frac{1}{N})^N (1N1)N,根据极限可以计算得到:
    l i m N → + ∞ ( 1 − 1 N ) N = 1 e ≈ 0.368 lim_{N\rightarrow +\infty}(1-\frac{1}{N})^N=\frac{1}{e}\approx 0.368 limN+(1N1)N=e10.368
    也就是采用自助法,会有 36.8% 的样本不会出现在训练集中,使用这部分样本作为测试集。这种方法也被称为包外估计。

    自助法的优点有:

    • 数据集比较小、难以有效划分训练/测试集时很有用:

    • 能从初始数据集中产生多个不同的训练集,这对集成学习等方法而言有很大好处。

    但也存在如下缺点:

    • 产生的数据集改变了初始数据集的分布,这会引入估计偏差。因此在初始数据量足够时,留出法和折交叉验证法更常用

    2.3 训练集、验证集、测试集

    简单介绍下训练集、验证集和测试集各自的作用:

    1. 训练集:主要就是训练模型,理论上越大越好;
    2. 验证集:用于模型调试超参数。通常要求验证集比较大,避免模型会对验证集过拟合;
    3. 测试集:用于评估模型的泛化能力。理论上,测试集越大,评估结果就约精准。另外,测试集必须不包含训练样本,否则会影响对模型泛化能力的评估。

    验证集和测试集的对比:

    • 测试集通常用于对模型的预测能力进行评估,它是提供模型预测能力的无偏估计;如果不需要对模型预测能力的无偏估计,可以不需要测试集;
    • 验证集主要是用于超参数的选择。

    2.4 划分数据集的比例选择方法

    那么一般如何选择划分训练、验证和测试集的比例呢?通常可以按照如下做法:

    1. 对于小批量数据,数据的拆分的常见比例为:
      • 如果未设置验证集,则将数据三七分:70% 的数据用作训练集、30% 的数据用作测试集。
      • 如果设置验证集,则将数据划分为:60% 的数据用作训练集、20%的数据用过验证集、20% 的数据用作测试集。
    2. 对于大批量数据,验证集和测试集占总数据的比例会更小
      • 对于百万级别的数据,其中 1 万条作为验证集、1 万条作为测试集即可。
      • 验证集的目的就是验证不同的超参数;测试集的目的就是比较不同的模型。
        • 一方面它们要足够大,才足够评估超参数、模型。
        • 另一方面,如果它们太大,则会浪费数据(验证集和训练集的数据无法用于训练)。
    3. k-fold 交叉验证中:先将所有数据拆分成 k 份,然后其中 1 份作为测试集,其他 k-1 份作为训练集。
      • 这里并没有验证集来做超参数的选择。所有测试集的测试误差的均值作为模型的预测能力的一个估计。
      • 使用 k-fold 交叉的原因是:样本集太小。如果选择一部分数据来训练,则有两个问题:
        • 训练数据的分布可能与真实的分布有偏离k-fold 交叉让所有的数据参与训练,会使得这种偏离得到一定程度的修正。
        • 训练数据太少,容易陷入过拟合k-fold 交叉让所有数据参与训练,会一定程度上缓解过拟合。

    2.5 分布不匹配

    深度学习时代,经常会发生:训练集和验证集、测试集的数据分布不同

    如:训练集的数据可能是从网上下载的高清图片,测试集的数据可能是用户上传的、低像素的手机照片。

    • 必须保证验证集、测试集的分布一致,它们都要很好的代表你的真实应用场景中的数据分布。
    • 训练数据可以与真实应用场景中的数据分布不一致,因为最终关心的是在模型真实应用场景中的表现。

    如果发生了数据不匹配问题,则可以想办法让训练集的分布更接近验证集

    • 一种做法是:收集更多的、分布接近验证集的数据作为训练集合
    • 另一种做法是:人工合成训练数据,使得它更接近验证集。该策略有一个潜在问题:你可能只是模拟了全部数据空间中的一小部分。导致你的模型对这一小部分过拟合。

    当训练集和验证集、测试集的数据分布不同时,有以下经验原则:

    • 确保验证集和测试集的数据来自同一分布

      因为需要使用验证集来优化超参数,而优化的最终目标是希望模型在测试集上表现更好。

    • 确保验证集和测试集能够反映未来得到的数据,或者最关注的数据

    • 确保数据被随机分配到验证集和测试集上

    当训练集和验证集、测试集的数据分布不同时,分析偏差和方差的方式有所不同

    • 如果训练集和验证集的分布一致,那么当训练误差和验证误差相差较大时,我们认为存在很大的方差问题

    • 如果训练集和验证集的分布不一致,那么当训练误差和验证误差相差较大时,有两种原因:

      • 第一个原因:模型只见过训练集数据,没有见过验证集的数据导致的,是数据不匹配的问题
      • 第二个原因:模型本来就存在较大的方差

      为了弄清楚原因,需要将训练集再随机划分为:训练-训练集训练-验证集。这时候,训练-训练集训练-验证集 是同一分布的。

      • 模型在训练-训练集训练-验证集 上的误差的差距代表了模型的方差
      • 模型在训练-验证集 和 验证集上的误差的差距代表了数据不匹配问题的程度

    3. 过拟合、欠拟合

    机器学习的两个主要挑战是过拟合和欠拟合

    过拟合(overfitting)指算法模型在训练集上的性能非常好,但是泛化能力很差,泛化误差很大,即在测试集上的效果却很糟糕的情况

    • 过拟合的原因:将训练样本本身的一些特点当作了所有潜在样本都具有的一般性质,这会造成泛化能力下降;另一个原因是模型可能学到训练集中的噪声,并基于噪声进行了预测
    • 过拟合无法避免,只能缓解。因为机器学习的问题通常是 NP 难甚至更难的,而有效的学习算法必然是在多项式时间内运行完成。如果可以避免过拟合,这就意味着构造性的证明了 P=NP

    欠拟合(underfitting)模型的性能非常差,在训练数据和测试数据上的性能都不好,训练误差和泛化误差都很大。其原因就是模型的学习能力比较差。

    一般可以通过挑战模型的容量来缓解过拟合和欠拟合问题。模型的容量是指其拟合各种函数的能力

    • 容量低的模型容易发生欠拟合,模型拟合能力太弱。
    • 容量高的模型容易发生过拟合,模型拟合能力太强。

    一般解决过拟合的方法有:

    • 简化模型,这包括了采用简单点的模型、减少特征数量,比如神经网络中减少网络层数或者权重参数,决策树模型中降低树的深度、采用剪枝等;
    • 增加训练数据,采用数据增强的方法,比如人工合成训练数据等;
    • 早停,当验证集上的误差没有进一步改善,训练提前终止;
    • 正则化,常用 L1 或者 L2 正则化。
    • 集成学习方法,训练多个模型,并以每个模型的平均输出作为结果,降低单一模型的过拟合风险,常用方法有 baggingboostingdropout(深度学习中的方法)等;
    • 噪声注入:包括输入噪声注入、输出噪声注入、权重噪声注入。将噪声分别注入到输入/输出/权重参数中,虽然噪声可能是模型过拟合的一个原因,但第一可以通过交叉验证来避免;第二就是没有噪声的完美数据也是很有可能发生过拟合;第三可以选择在特征、权值参数加入噪声,而非直接在数据加入噪声。

    解决欠拟合的方法有:

    • 选择一个更强大的模型,带有更多参数
    • 更好的特征训练学习算法(特征工程)
    • 减小对模型的限制(比如,减小正则化超参数)

    4. 超参数调优

    超参数调优是一件非常头疼的事情,很多时候都需要一些先验知识来选择合理的参数值,但如果没有这部分先验知识,要找到最优的参数值是很困难,非常耗费时间和精力。但超参数调优确实又可以让模型性能变得更加的好。

    在选择超参数调优算法前,需要明确以下几个要素:

    • 目标函数。算法需要最大化/最小化的目标;
    • 搜索范围。一般通过上下限来确定;
    • 算法的其他参数,比如搜索步长。

    4.1 搜索策略

    常用的几种超参数搜索策略如下:

    • 手动搜索:需要较好的先验知识经验
    • 网格搜索:超参数的数据相对较少的时候,这个方法比较实用
    • 随机搜索:通常推荐这种方式
    • 贝叶斯优化算法:基于模型的搜索方法,利用了历史搜索结果
    4.1.1 手动搜索
    1. 手动选择超参数需要了解超参数做了些什么,以及机器学习模型如何才能取得良好的泛化

    2. 手动搜索超参数的任务是:在给定运行时间和内存预算范围的条件下,最小化泛化误差

    3. 手动调整超参数时不要忘记最终目标:提升测试集性能

      • 加入正则化只是实现这个目标的一种方法。

      • 如果训练误差很低,也可以通过收集更多的训练数据来减少泛化误差。

        如果训练误差太大,则收集更多的训练数据就没有意义。

      • 实践中的一种暴力方法是:不断提高模型容量和训练集的大小

        这种方法增加了计算代价,只有在拥有充足的计算资源时才可行

    4.1.2 网格搜索

    网格搜索可能是最简单也是应用最广泛的超参数搜索算法了。它的几种做法如下:

    • 采用较大的搜索范围和较小的搜索步长,很大概率会搜索到全局最优值,但十分耗费计算资源和时间,特别是超参数比较多的时候;
    • 先采用较大搜索范围和较大步长,寻找全局最优的可能位置,然后逐渐缩小搜索范围和步长,来确定更精确的最优值。可以降低所需要的计算时间和计算量,但由于目标函数一般都是非凸的,可能会错过全局最优值。

    网格搜索也可以借助 sklearn 实现,简单的示例代码如下:

    from sklearn.model_selection import	GridSearchCV
    from sklearn.ensemble import RandomForestClassifier
    param_grid = [
        {'n_estimators': [3, 10, 30], 'max_features': [2, 4, 6, 8]},
        {'bootstrap': [False], 'n_estimators': [3, 10], 'max_features': [2, 3, 4]},
    ]
    forest_reg = RandomForestRegressor()
    grid_search = GridSearchCV(forest_reg, param_grid, cv=5,
                               scoring='neg_mean_squared_error')
    grid_search.fit(data, labels)
    
    4.1.3 随机搜索

    随机搜索是一种可以替代网格搜索的方法,它编程简单、使用方便、能更快收敛到超参数的良好取值。

    • 首先为每个超参数定义一个边缘分布,如伯努利分布(对应着二元超参数)或者对数尺度上的均匀分布(对应着正实值超参数)。
    • 然后假设超参数之间相互独立,从各分布中抽样出一组超参数。
    • 使用这组超参数训练模型。
    • 经过多次抽样 -> 训练过程,挑选验证集误差最小的超参数作为最好的超参数。

    随机搜索的优点如下:

    • 不需要离散化超参数的值,也不需要限定超参数的取值范围。这允许我们在一个更大的集合上进行搜索。
    • 当某些超参数对于性能没有显著影响时,随机搜索相比于网格搜索指数级地高效,它能更快的减小验证集误差

    随机搜索比网格搜索更快的找到良好超参数的原因是:没有浪费的实验

    • 在网格搜索中,两次实验之间只会改变一个超参数 (假设为 m)的值,而其他超参数的值保持不变。如果这个超参数 m 的值对于验证集误差没有明显区别,那么网格搜索相当于进行了两个重复的实验。
    • 在随机搜索中,两次实验之间,所有的超参数值都不会相等,因为每个超参数的值都是从它们的分布函数中随机采样而来。因此不大可能会出现两个重复的实验。
    • 如果 m 超参数与泛化误差无关,那么不同的 m 值:
      • 在网格搜索中,不同 m 值、相同的其他超参数值,会导致大量的重复实验
      • 在随机搜索中,其他超参数值每次也都不同,因此不大可能出现两个重复的实验(除非所有的超参数都与泛化误差无关)。

    随机搜索可以采用 sklearn.model_selection 中的 RandomizedSearchCV 方法。

    4.1.4 贝叶斯优化方法

    贝叶斯优化方法是基于模型的参数搜索算法的一种比较常见的算法。它相比于前面的网格搜索和随机搜索,最大的不同就是利用历史的搜索结果进行优化搜索。主要是由四部分组成的:

    1. 目标函数。大部分情况是模型验证集上的损失;
    2. 搜索空间。各类待搜索的超参数;
    3. 优化策略。建立的概率模型和选择超参数的方式;
    4. 历史的搜索结果。

    贝叶斯优化算法的步骤如下:

    1. 根据先验分布,假设一个搜索函数;
    2. 然后,每一次采用新的采样点来测试目标函数时,利用这个信息更新目标函数的先验分布;
    3. 最后,算法测试由后验分布给出的全局最优最可能出现的位置的点。

    需要特别注意的是,贝叶斯优化算法容易陷入局部最优值:它在找到一个局部最优值后,会不断在该区域进行采样。

    因此,贝叶斯优化算法会在探索和利用之间找到一个平衡点,探索是在还未取样的区域获取采样点,利用则是根据后验分布在最可能出现全局最优的区域进行采样。

    4.2 调整原则

    1. 通常先对超参数进行粗调,然后在粗调中表现良好的超参数区域进行精调

    2. 超参数随机搜索,并不意味着是在有效范围内随机均匀取值。需要选择合适的缩放来进行随机选取。

    3. 通常情况下,建议至少每隔几个月重新评估或者修改超参数。因为随着时间的变化,真实场景的数据会逐渐发生改变:

      • 可能是由于用户的行为、偏好发生了改变。
      • 可能是采样的方式发生了改变。
      • 也可能仅仅是由于数据中心更新了服务器。

      由于这些变化,原来设定的超参数可能不再适用。

    4. 有两种超参数调整策略:

      • 如果数据足够大且没有足够的计算资源,此时只能一次完成一个试验。

        可以每天观察模型的表现,实时的、动态的调整超参数

      • 如果数据不大,有足够的计算资源可以同一时间完成大量的试验,则可以设置多组超参数设定,然后选择其中表现最好的那个


    小结

    关于模型评估方面的内容就介绍这么多,文章有些长,而且内容也比较多。

    关于如何构建一个机器学习项目的内容,基本到本文就介绍完毕了,从开始的评估问题,获取数据,到数据预处理、特征工程,然后就是各种常见机器学习算法的评估,最后就是模型评估部分的内容了。

    当然了,本系列的文章还是偏向于理论,代码比较少,主要也是整理和总结书本以及网上文章的知识点。

    所以下一篇文章会是介绍一篇手把手教你运用机器学习算法来做分类的文章,来自国外一个大神的博客文章,主要是面向机器学习的初学者。


    参考:

    欢迎关注我的微信公众号–机器学习与计算机视觉,或者扫描下方的二维码,大家一起交流,学习和进步!

    往期精彩推荐

    机器学习系列
    Github项目 & 资源教程推荐
    展开全文
  • yolov--10--目标检测模型的参数评估指标详解、概念解析 yolov--11--YOLO v3的原版训练记录、mAP、AP、recall、precision、time等评价指标计算 yolov--12--YOLOv3的原理深度剖析和关键点讲解 1、目标检测中...

    Yolov-1-TX2上用YOLOv3训练自己数据集的流程(VOC2007-TX2-GPU)

    Yolov--2--一文全面了解深度学习性能优化加速引擎---TensorRT

    Yolov--3--TensorRT中yolov3性能优化加速(基于caffe)

    yolov-5-目标检测:YOLOv2算法原理详解

    yolov--8--Tensorflow实现YOLO v3

    yolov--9--YOLO v3的剪枝优化

    yolov--10--目标检测模型的参数评估指标详解、概念解析

    yolov--11--YOLO v3的原版训练记录、mAP、AP、recall、precision、time等评价指标计算

    yolov--12--YOLOv3的原理深度剖析和关键点讲解


    1、目标检测中的mAP是什么含义?

    • 关于Ground Truth:ground truth包括图像中物体的类别以及该图像中每个物体的真实边界框
    • mAP含义及计算:如何量化呢? 采用IoU(Intersection over Union),IoU是预测框与ground truth的交集和并集的比值。这个量也被称为Jaccard指数,并于20世纪初由Paul Jaccard首次提出。

    • 鉴别正确的检测结果并计算precision和recall

    为了计算precision和recall,与所有机器学习问题一样,我们必须鉴别出True Positives(真正例)、False Positives(假正例)、True Negatives(真负例)和 False Negatives(假负例)。


    TP TN FP FN是什么意思?

    T或者N代表的是该样本是否被分类分对,P或者N代表的是该样本被分为什么

    TP(True Positives)意思我们倒着来翻译就是“被分为正样本,并且分对了”,TN(True Negatives)意思是“被分为负样本,而且分对了”,FP(False Positives)意思是“被分为正样本,但是分错了”,FN(False Negatives)意思是“被分为负样本,但是分错了”。

    按下图来解释,左半矩形是正样本,右半矩形是负样本。一个2分类器,在图上画了个圆,分类器认为圆内是正样本,圆外是负样本。那么左半圆分类器认为是正样本,同时它确实是正样本,那么就是“被分为正样本,并且分对了”即TP,左半矩形扣除左半圆的部分就是分类器认为它是负样本,但是它本身却是正样本,就是“被分为负样本,但是分错了”即FN。右半圆分类器认为它是正样本,但是本身却是负样本,那么就是“被分为正样本,但是分错了”即FP。右半矩形扣除右半圆的部分就是分类器认为它是负样本,同时它本身确实是负样本,那么就是“被分为负样本,而且分对了”即TN

    有了上面TP TN FP FN的概念,这个Precision和Recall的概念一张图就能说明。

     

    ,翻译成中文就是“分类器认为是正类并且确实是正类的部分占所有分类器认为是正类的比例”,衡量的是一个分类器分出来的正类的确是正类的概率。两种极端情况就是,如果精度是100%,就代表所有分类器分出来的正类确实都是正类。如果精度是0%,就代表分类器分出来的正类没一个是正类。光是精度还不能衡量分类器的好坏程度,比如50个正样本和50个负样本,我的分类器把49个正样本和50个负样本都分为负样本,剩下一个正样本分为正样本,这样我的精度也是100%,但是傻子也知道这个分类器很垃圾。

    ,翻译成中文就是“分类器认为是正类并且确实是正类的部分占所有确实是正类的比例”,衡量的是一个分类能把所有的正类都找出来的能力。两种极端情况,如果召回率是100%,就代表所有的正类都被分类器分为正类。如果召回率是0%,就代表没一个正类被分为正类。

    原文:https://blog.csdn.net/hsqyc/article/details/81702437 


    mAP定义及相关概念

    • mAP: mean Average Precision, 即各类别AP的平均值
    • AP: PR曲线下面积,后文会详细讲解
    • PR曲线: Precision-Recall曲线
    • Precision: TP / (TP + FP)
    • Recall: TP / (TP + FN)
    • TP: IoU>0.5的检测框数量(同一Ground Truth只计算一次)
    • FP: IoU<=0.5的检测框,或者是检测到同一个GT的多余检测框的数量
    • FN: 没有检测到的GT的数量,IoU=0

    mAP如何计算?

    由前面定义,我们可以知道,要计算mAP必须先绘出各类别PR曲线,计算出AP。而如何采样PR曲线,VOC采用过两种不同方法。参见:The PASCAL Visual Object Classes Challenge 2012 (VOC2012) Development Kit

    在VOC2010以前,只需要选取当Recall >= 0, 0.1, 0.2, ..., 1共11个点时的Precision最大值,然后AP就是这11个Precision的平均值。

    在VOC2010及以后,需要针对每一个不同的Recall值(包括0和1),选取其大于等于这些Recall值时的Precision最大值,然后计算PR曲线下面积作为AP值。

    mAP计算示例

    假设,对于Aeroplane类别,我们网络有以下输出(BB表示BoundingBox序号,IoU>0.5时GT=1):

    BB  | confidence | GT
    ----------------------
    BB1 |  0.9       | 1
    ----------------------
    BB2 |  0.9       | 1
    ----------------------
    BB1 |  0.8       | 1
    ----------------------
    BB3 |  0.7       | 0
    ----------------------
    BB4 |  0.7       | 0
    ----------------------
    BB5 |  0.7       | 1
    ----------------------
    BB6 |  0.7       | 0
    ----------------------
    BB7 |  0.7       | 0
    ----------------------
    BB8 |  0.7       | 1
    ----------------------
    BB9 |  0.7       | 1
    ----------------------
    

    因此,我们有 TP=5 (BB1, BB2, BB5, BB8, BB9), FP=5 (重复检测到的BB1也算FP)。除了表里检测到的5个GT以外,我们还有2个GT没被检测到,因此: FN = 2. 这时我们就可以按照Confidence的顺序给出各处的PR值,如下:

    rank=1  precision=1.00 and recall=0.14
    ----------
    rank=2  precision=1.00 and recall=0.29
    ----------
    rank=3  precision=0.66 and recall=0.29
    ----------
    rank=4  precision=0.50 and recall=0.29
    ----------
    rank=5  precision=0.40 and recall=0.29
    ----------
    rank=6  precision=0.50 and recall=0.43
    ----------
    rank=7  precision=0.43 and recall=0.43
    ----------
    rank=8  precision=0.38 and recall=0.43
    ----------
    rank=9  precision=0.44 and recall=0.57
    ----------
    rank=10 precision=0.50 and recall=0.71
    ----------
    

    对于上述PR值,如果我们采用:

    1. VOC2010之前的方法,我们选取Recall >= 0, 0.1, ..., 1的11处Percision的最大值:1, 1, 1, 0.5, 0.5, 0.5, 0.5, 0.5, 0, 0, 0。此时Aeroplane类别的 AP = 5.5 / 11 = 0.5
    2. VOC2010及以后的方法,对于Recall >= 0, 0.14, 0.29, 0.43, 0.57, 0.71, 1,我们选取此时Percision的最大值:1, 1, 1, 0.5, 0.5, 0.5, 0。此时Aeroplane类别的 AP = (0.14-0)*1 + (0.29-0.14)*1 + (0.43-0.29)*0.5 + (0.57-0.43)*0.5 + (0.71-0.57)*0.5 + (1-0.71)*0 = 0.5

    mAP就是对每一个类别都计算出AP然后再计算AP平均值就好了

    更多信息

    建议参考GluonCV库里面的voc_detection.py实现了两种mAP计算方式,思路清晰:

    https://github.com/dmlc/gluon-cv/blob/master/gluoncv/utils/metrics/voc_detection.py​github.com

    https://www.zhihu.com/search?type=content&q=%E7%9B%AE%E6%A0%87%E6%A3%80%E6%B5%8B%E4%B8%AD%E7%9A%84mAP%E6%98%AF%E4%BB%80%E4%B9%88


    precision、Recall如何计算?

    为了获得True Positives and False Positives,我们需要使用IoU。计算IoU,我们从而确定一个检测结果(Positive)是正确的(True)还是错误的(False)。

    • 最常用的阈值是0.5,即如果IoU> 0.5,则认为它是True Positive,否则认为是False Positive。而COCO数据集的评估指标建议对不同的IoU阈值进行计算,但为简单起见,我们这里仅讨论一个阈值0.5,这是PASCAL VOC数据集所用的指标
    • 为了计算Recall,我们需要Negatives的数量。由于图片中我们没有预测到物体的每个部分都被视为Negative,因此计算True Negatives比较难办。但是我们可以只计算False Negatives,即我们模型所漏检的物体
    • 另外一个需要考虑的因素是模型所给出的各个检测结果的置信度。通过改变置信度阈值,我们可以改变一个预测框是Positive还是 Negative,即改变预测值的正负性(不是box的真实正负性,是预测正负性)。基本上,阈值以上的所有预测(Box + Class)都被认为是Positives,并且低于该值的都是Negatives

    对于每一个图片,ground truth数据会给出该图片中各个类别的实际物体数量。我们可以计算每个Positive预测框与ground truth的IoU值,并取最大的IoU值,认为该预测框检测到了那个IoU最大的ground truth。然后根据IoU阈值,我们可以计算出一张图片中各个类别的正确检测值(True Positives, TP)数量以及错误检测值数量(False Positives, FP)。据此,可以计算出各个类别的precision:
    \\precision = \frac{TP}{TP+FP}

    既然我们已经得到了正确的预测值数量(True Positives),也很容易计算出漏检的物体数(False Negatives, FN)。据此可以计算出Recall(其实分母可以用ground truth总数):
    \\recall = \frac{TP}{TP+FN}


     计算mAP的详细流程?

    在目标检测中,mAP的定义首先出现在PASCAL Visual Objects Classes(VOC)竞赛中,这个大赛包含许多图像处理任务,详情可以参考这个paper(里面包含各个比赛的介绍以及评估等)。

    前面我们已经讲述了如何计算Precision和Recall,但是,正如前面所述,至少有两个变量会影响Precision和Recall,即IoU和置信度阈值。IoU是一个简单的几何度量,可以很容易标准化,比如在PASCAL VOC竞赛中采用的IoU阈值为0.5,而COCO竞赛中在计算mAP较复杂,其计算了一系列IoU阈值(0.05至0.95)下的mAP。但是置信度却在不同模型会差异较大,可能在我的模型中置信度采用0.5却等价于在其它模型中采用0.8置信度,这会导致precision-recall曲线变化。为此,PASCAL VOC组织者想到了一种方法来解决这个问题,即要采用一种可以用于任何模型的评估指标。在paper中,他们推荐使用如下方式计算Average Precision(AP):

    For a given task and class, the precision/recall curve is computed from a method’s ranked output. Recall is defined as the proportion of all positive examples ranked above a given rank. Precision is the proportion of all examples above that rank which are from the positive class. The AP summarises the shape of the precision/recall curve, and is defined as the mean precision at a set of eleven equally spaced recall levels [0,0.1,...,1]:

    可以看到,为了得到precision-recall曲线,首先要对模型预测结果进行排序(ranked output,按照各个预测值置信度降序排列)。那么给定一个rank,Recall和Precision仅在高于该rank值的预测结果中计算,改变rank值会改变recall值。这里共选择11个不同的recall([0, 0.1, ..., 0.9, 1.0]),可以认为是选择了11个rank,由于按照置信度排序,所以实际上等于选择了11个不同的置信度阈值。那么,AP就定义为在这11个recall下precision的平均值,其可以表征整个precision-recall曲线(曲线下面积)。

    另外,在计算precision时采用一种插值方法(interpolate):

    The precision at each recall level r is interpolated by taking the maximum precision measured for a method for which the corresponding recall exceeds r:
    The intention in interpolating the precision/recall curve in this way is to reduce the impact of the “wiggles” in the precision/recall curve, caused by small variations in the ranking of examples.

    及对于某个recall值r,precision值取所有recall>=r中的最大值(这样保证了p-r曲线是单调递减的,避免曲线出现摇摆):

    不过这里VOC数据集在2007年提出的mAP计算方法,而在2010之后却使用了所有数据点,而不是仅使用11个recall值来计算AP(详细参考这篇paper):

    Up until 2009 interpolated average precision (Salton and Mcgill 1986) was used to evaluate both classification and detection. However, from 2010 onwards the method of computing AP changed to use all data points rather than TREC-style sampling (which only sampled the monotonically decreasing curve at a fixed set of uniformly-spaced recall values 0, 0.1, 0.2,..., 1). The intention in interpolating the precision–recall curve was to reduce the impact of the ‘wiggles’ in the precision–recall curve, caused by small variations in the ranking of examples. However, the downside of this interpolation was that the evaluation was too crude to discriminate between the methods at low AP.

    对于各个类别,分别按照上述方式计算AP,取所有类别的AP平均值就是mAP。这就是在目标检测问题中mAP的计算方法。可能有时会发生些许变化,如COCO数据集采用的计算方式更严格,其计算了不同IoU阈值和物体大小下的AP(详情参考COCO Detection Evaluation)。

    当比较mAP值,记住以下要点:

    1. mAP通常是在一个数据集上计算得到的。
    2. 虽然解释模型输出的绝对量化并不容易,但mAP作为一个相对较好的度量指标可以帮助我们。 当我们在流行的公共数据集上计算这个度量时,该度量可以很容易地用来比较目标检测问题的新旧方法。
    3. 根据训练数据中各个类的分布情况,mAP值可能在某些类(具有良好的训练数据)非常高,而其他类(具有较少/不良数据)却比较低。所以你的mAP可能是中等的,但是你的模型可能对某些类非常好,对某些类非常不好。因此,建议在分析模型结果时查看各个类的AP值。这些值也许暗示你需要添加更多的训练样本。

    代码实现

    Facebook开源的Detectron包含VOC数据集的mAP计算,这里贴出其核心实现,以对mAP的计算有更深入的理解。首先是precision和recall的计算:

    # 按照置信度降序排序
    sorted_ind = np.argsort(-confidence)
    BB = BB[sorted_ind, :]   # 预测框坐标
    image_ids = [image_ids[x] for x in sorted_ind] # 各个预测框的对应图片id
    
    # 便利预测框,并统计TPs和FPs
    nd = len(image_ids)
    tp = np.zeros(nd)
    fp = np.zeros(nd)
    for d in range(nd):
        R = class_recs[image_ids[d]]
        bb = BB[d, :].astype(float)
        ovmax = -np.inf
        BBGT = R['bbox'].astype(float)  # ground truth
    
        if BBGT.size > 0:
            # 计算IoU
            # intersection
            ixmin = np.maximum(BBGT[:, 0], bb[0])
            iymin = np.maximum(BBGT[:, 1], bb[1])
            ixmax = np.minimum(BBGT[:, 2], bb[2])
            iymax = np.minimum(BBGT[:, 3], bb[3])
            iw = np.maximum(ixmax - ixmin + 1., 0.)
            ih = np.maximum(iymax - iymin + 1., 0.)
            inters = iw * ih
    
            # union
            uni = ((bb[2] - bb[0] + 1.) * (bb[3] - bb[1] + 1.) +
                   (BBGT[:, 2] - BBGT[:, 0] + 1.) *
                   (BBGT[:, 3] - BBGT[:, 1] + 1.) - inters)
    
            overlaps = inters / uni
            ovmax = np.max(overlaps)
            jmax = np.argmax(overlaps)
        # 取最大的IoU
        if ovmax > ovthresh:  # 是否大于阈值
            if not R['difficult'][jmax]:  # 非difficult物体
                if not R['det'][jmax]:    # 未被检测
                    tp[d] = 1.
                    R['det'][jmax] = 1    # 标记已被检测
                else:
                    fp[d] = 1.
        else:
            fp[d] = 1.
    
    # 计算precision recall
    fp = np.cumsum(fp)
    tp = np.cumsum(tp)
    rec = tp / float(npos)
    # avoid divide by zero in case the first detection matches a difficult
    # ground truth
    prec = tp / np.maximum(tp + fp, np.finfo(np.float64).eps)

    这里最终得到一系列的precision和recall值,并且这些值是按照置信度降低排列统计的,可以认为是取不同的置信度阈值(或者rank值)得到的。然后据此可以计算AP:

    def voc_ap(rec, prec, use_07_metric=False):
        """Compute VOC AP given precision and recall. If use_07_metric is true, uses
        the VOC 07 11-point method (default:False).
        """
        if use_07_metric:  # 使用07年方法
            # 11 个点
            ap = 0.
            for t in np.arange(0., 1.1, 0.1):
                if np.sum(rec >= t) == 0:
                    p = 0
                else:
                    p = np.max(prec[rec >= t])  # 插值
                ap = ap + p / 11.
        else:  # 新方式,计算所有点
            # correct AP calculation
            # first append sentinel values at the end
            mrec = np.concatenate(([0.], rec, [1.]))
            mpre = np.concatenate(([0.], prec, [0.]))
    
            # compute the precision 曲线值(也用了插值)
            for i in range(mpre.size - 1, 0, -1):
                mpre[i - 1] = np.maximum(mpre[i - 1], mpre[i])
    
            # to calculate area under PR curve, look for points
            # where X axis (recall) changes value
            i = np.where(mrec[1:] != mrec[:-1])[0]
    
            # and sum (\Delta recall) * prec
            ap = np.sum((mrec[i + 1] - mrec[i]) * mpre[i + 1])
        return ap

    计算各个类别的AP值后,取平均值就可以得到最终的mAP值了。但是对于COCO数据集相对比较复杂,不过其提供了计算的API,感兴趣可以看一下cocodataset/cocoapi

    注:译文相比原文略有删改。

    https://zhuanlan.zhihu.com/p/37910324

    https://www.zhihu.com/search?type=content&q=%E7%9B%AE%E6%A0%87%E6%A3%80%E6%B5%8B%E4%B8%AD%E7%9A%84mAP%E6%98%AF%E4%BB%80%E4%B9%88


    举例计算mAP


    有3张图如下,要求算法找出face。蓝色框代表标签label,绿色框代表算法给出的结果pre,旁边的红色小字代表置信度。设定第一张图的检出框叫pre1,第一张的标签框叫label1。第二张、第三张同理。

    1.根据IOU计算TP,FP
    首先我们计算每张图的pre和label的IOU,根据IOU是否大于0.5来判断该pre是属于TP还是属于FP。显而易见,pre1是TP,pre2是FP,pre3是TP。

    2.根据置信度阈值排序
    根据每个pre的置信度阈值进行从高到低排序,这里pre1、pre2、pre3置信度刚好就是从高到低。

    3.在不同置信度阈值下获得Precision和Recall
    首先,设置阈值为0.9,无视所有小于0.9的pre。那么检测器检出的所有框pre即TP+FP=1,并且pre1是TP,那么Precision=1/1。因为所有的label=3,所以Recall=1/3。这样就得到一组P、R值。

    然后,设置阈值为0.8,无视所有小于0.8的pre。那么检测器检出的所有框pre即TP+FP=2,因为pre1是TP,pre2是FP,那么Precision=1/2=0.5。因为所有的label=3,所以Recall=1/3=0.33。这样就又得到一组P、R值。

    再然后,设置阈值为0.7,无视所有小于0.7的pre。那么检测器检出的所有框pre即TP+FP=3,因为pre1是TP,pre2是FP,pre3是TP,那么Precision=2/3=0.67。因为所有的label=3,所以Recall=2/3=0.67。这样就又得到一组P、R值。

    4.绘制PR曲线并计算AP值
    根据上面3组PR值绘制PR曲线如下。然后每个“峰值点”往左画一条线段直到与上一个峰值点的垂直线相交。这样画出来的红色线段与坐标轴围起来的面积就是AP值。

    在这里

    5.计算mAP
    AP衡量的是对一个类检测好坏,mAP就是对多个类的检测好坏。就是简单粗暴的把所有类的AP值取平均就好了。比如有两类,类A的AP值是0.5,类B的AP值是0.2,那么mAP=(0.5+0.2)/2=0.35

    Github
    在这里提供一个计算mAP的github地址,上面那张AP图就是它生成的。我觉得还蛮好用的。

    https://github.com/Cartucho/mAP


    若加微信请备注下姓名_公司/学校,相遇即缘分,感谢您的支持,愿真诚交流,共同进步,谢谢~ 

     


    原文:https://blog.csdn.net/hsqyc/article/details/81702437 
    https://pan.baidu.com/s/1b19ewu

    展开全文
  • 通过 OPC 算法 [1] 评估具有 1、2 或 3 个参数的 Mittag-Leffler (ML) 函数。 该例程评估 ML 函数 E 的近似值 Et,使得 |E-Et|/(1+|E|) 大约为 1.0e-15 E = ML(z,alpha) 用一个参数 alpha 为 z 的对应元素计算 ML ...
  • 模型评估参数调优

    千次阅读 2018-09-11 17:54:33
    基于流水线的工作流 scikit-learn 中的 Pipeline 类 ...在本例中,需要验证的是参数 C,即定义在 scikit-learn pipeline 中逻辑回归分类器的正则化参数,将其记为 ‘clf__C’,并通过 param_range 参数...

    基于流水线的工作流

    scikit-learn 中的 Pipeline 类。它使得我们可以拟合包含任意多个处理步骤的模型,并模型用于新数据的预测。

    案例1:威斯康星乳腺癌数据集

    使用Breast Cancer Wisconsin 数据集,此数据集共包含569个恶性或良性肿瘤细胞样本。数据集的前两列分别存储了样本唯一的ID以及对样本的诊断结果(M代表恶性,B代表良性)。数据集的3~32列包含了30个从细胞核照片中提取、用实数值标识的特征,它们可以用于构建判定模型,对肿瘤是良性还是恶性做出预测。

    1.数据加载
    import pandas as pd 
    df = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-\
    databases/breast-cancer-wisconsin/wdbc.data', header=None)
    
    from sklearn.preprocessing import LabelEncoder
    X = df.loc[:, 2:].values
    y = df.loc[:, 1].values
    le = LabelEncoder()
    y = le.fit_transform(y)
    
    from sklearn.model_selection import train_test_split
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2,
        random_state=1)
    

    使用 scikit-learn 中的 LabelEncoder 类,可以将类标从原始字符串表示(M 或 B)转换为整数。

    2.在pipeline中集成数据转换及评估操作

    通过 Pipeline 类将 StandardScaler、PCA、LogisticRegression 对象串联起来,无需在训练数据集和测试数据集上分别进行模型拟合。

    from sklearn.preprocessing import StandardScaler
    from sklearn.decomposition import PCA
    from sklearn.linear_model import LogisticRegression
    from sklearn.pipeline import Pipeline
    pipe_lr = Pipeline([
        ('sc1', StandardScaler()),
        ('pca', PCA(n_components=2)),
        ('clf', LogisticRegression(random_state=1)),
    ])
    pipe_lr.fit(X_train, y_train)
    print('Test Accuracy: %.3f' % (pip_lr.score(X_test, y_test)))

    执行结果如下:

    Test Accuracy: 0.947

    Pipeline 对象采用元组的序列作为输入,其中每个元组中的第一个值为一个字符串,可以是任意的标识符,我们通过它来访问 pipeline 中的元素,而元组的第二个值则为 scikit-learn 中的一个转换器或分类器。

    注意:pipeline 中没有设定中间步骤的数量。pipeline 的工作方式可用下图来描述:
    pipeline工作方式

    k-fold 交叉验证评估模型性能

    构建机器学习模型的一个关键步骤就是在新数据上对模型性能进行评估。如果一个模型过于简单,将会面临欠拟合(高偏差)的问题,而模型基于训练数据构造得过于复杂,则会到时过拟合(高方差)问题。为了在偏差和方差之间找到可接受的折中方案,需要对模型惊醒评估。目前,常用的交叉验证技术有:holdout 交叉验证和k-fold交叉验证。

    holdout法

    通过 holdout 将原始数据集划分为训练数据集和测试数据集,为进一步提高模型在预测未知数据上的性能,需要对不同参数设置进行调优和比较。

    但是,如果再模型选择过程中不断重复使用相同的测试数据,它们可看作训练数据的一部分,模型更易于过拟合。

    使用 holdout 进行模型选择更好的方法是将数据划分为三个部分:训练数据集、验证数据集和测试数据集。训练数据集用于不同模型的拟合,模型在验证数据集上的性能表现作为模型选择的标准。交叉验证中,在使用不同参数值对模型进行训练后,使用验证数据集反复进行模型性能的评估。一旦参数优化获得较为满意的结果,就可以在测试数据集上对模型的泛化误差进行评估。
    holdout法

    holdout 法的一个缺点在于:模型性能的评估对训练数据集划分为训练及验证子集的方法是敏感的,评价结果会随着样本的不同而发生变化。

    k-fold交叉验证法

    在 k-fold 交叉验证中,将训练数据集划分为 k 个,其中 k - 1 个用于模型的训练,剩余1个用于测试。重复此过程 k 次,就得到了 k 个模型及对模型性能的评价。

    与 holdout 方法相比,k-fold 得到的结果对数据划分方法的敏感性相对较低。通过情况下,将 k-fold 交叉验证用于模型的调优,一旦找到了满意的超参值,就可以在全部的训练数据上重新训练模型,并使用独立的测试数据集对模型性能做出最终评价。

    由于 k-fold 交叉验证使用了无重复抽样技术,该方法的优势在于(每次迭代过程中)每个样本点只有一次被划入训练数据集或测试数据集的机会,与 holdout 法相比,这将使得模型性能的评估具有较小的方差。
    k-fold交叉验证法

    k-fold 交叉验证法的一个特例就是留一交叉验证法(leave-one-out, LOO),LOO将训练数据集划分为等同于样本数的子集,这样每次只有一个样本用于测试,当数据集非常小时,建议使用此方法进行验证。

    分层 k-fold 交叉验证法对标准 k-fold 交叉验证法进行了改进,它可以获得偏差和方差都较低的评价,特别是类别比例相差较大时。在分层交叉验证中,类别比列在每一个分块中得以保持,这使得每个分块中的类别比例与训练数据的整体比例一致。

    from sklearn.cross_validation import StratifiedKFold
    import numpy as np 
    kfold = StratifiedKFold(y=y_train, n_folds=10, random_state=1)
    scores = []
    for k, (train, test) in enumerate(kfold):
        pipe_lr.fit(X_train[train], y_train[train])
        score = pipe_lr.score(X_train[test], y_train[test])
        scores.append(score)
        print('Fold: %s, Class dist: %s, Acc: %.3f' % (k + 1, np.bincount(y_train[train]), score))
    print('CV accuracy: %.3f +/- %.3f' % (np.mean(scores), np.std(scores)))

    执行结果如下:

    Fold: 1, Class dist: [256 153], Acc: 0.891
    Fold: 2, Class dist: [256 153], Acc: 0.978
    Fold: 3, Class dist: [256 153], Acc: 0.978
    Fold: 4, Class dist: [256 153], Acc: 0.913
    Fold: 5, Class dist: [256 153], Acc: 0.935
    Fold: 6, Class dist: [257 153], Acc: 0.978
    Fold: 7, Class dist: [257 153], Acc: 0.933
    Fold: 8, Class dist: [257 153], Acc: 0.956
    Fold: 9, Class dist: [257 153], Acc: 0.978
    Fold: 10, Class dist: [257 153], Acc: 0.956
    CV accuracy: 0.950 +/- 0.029

    scikit-learn 中同样也实现了 k-fold 交叉验证评分的计算,这使得我们可以更加高效地使用分层 k-fold 交叉验证对模型进行评估:

    from sklearn.cross_validation import cross_val_score
    scores = cross_val_score(estimator=pipe_lr, X=X_train, y=y_train, cv=10, n_jobs=1)
    print('CV accuracy scores: %s' % (scores))
    print('CV accuracy: %.3f +/- %.3f' % (np.mean(scores), np.std(scores)))
    

    执行结果如下:

    CV accuracy scores: [0.89130435 0.97826087 0.97826087 0.91304348 0.93478261 0.97777778
     0.93333333 0.95555556 0.97777778 0.95555556]
    CV accuracy: 0.950 +/- 0.029

    通过学习及验证曲线来调试算法

    偏差与方差
    左上方图像显示的是一个高偏差模型,此模型的训练准确率和交叉验证准确率都很低,这表明此模型未能很好地拟合数据,解决此问题的常用方法是增加模型中参数的数量。
    右上方图像的模型面临高方差问题,表明训练准确度与交叉验证准确度之间有很大差距,针对此类过拟合问题,可以收集更多的训练数据或降低模型的复杂度。

    from sklearn.preprocessing import StandardScaler
    from sklearn.pipeline import Pipeline
    from sklearn.linear_model import LogisticRegression
    from sklearn.learning_curve import learning_curve
    import matplotlib.pyplot as plt 
    
    pipe_lr = Pipeline([
        ('sc1', StandardScaler()),
        ('clf', LogisticRegression(penalty='l2', random_state=0))
    ])
    train_sizes, train_scores, test_scores = learning_curve(
        estimator=pipe_lr,
        X=X_train,
        y=y_train,
        train_sizes=np.linspace(0.1, 1.0, 10),
        cv=10,
        n_jobs=1
    )
    train_mean = np.mean(train_scores, axis=1)
    train_std = np.std(train_scores, axis=1)
    test_mean = np.mean(test_scores, axis=1)
    test_std = np.std(test_scores, axis=1)
    plt.plot(train_sizes, train_mean, color='blue', marker='o', markersize=5, label='training accuracy')
    plt.plot(train_sizes, train_mean + train_std, train_mean - train_std, alpha=0.15, color='blue')
    plt.plot(train_sizes, test_mean, color='green', linestyle='--', marker='s', markersize=5, label='validation accuracy')
    plt.fill_between(train_sizes, test_mean + test_std, test_mean - test_std, alpha=0.15, color='green')
    plt.grid()
    plt.xlabel('Number of training samples')
    plt.ylabel('Accuracy')
    plt.legend(loc='lower right')
    plt.ylim([0.8, 1.0])
    plt.show()

    学习曲线评估模型

    由上图学习曲线图像可见,模型咋测试数据集上表现良好。但是,在训练准确率曲线与交叉验证准确率曲线之间,存在着相对较小差距,说明模型对训练数据有轻微的过拟合

    通过验证曲线来判定过拟合与欠拟合

    验证曲线是一种通过定位过拟合或欠拟合等诸多问题所在,来帮助提高模型性能的有效工具。验证曲线与学习曲线相似,不过绘制的不是样本大小与训练准确率、测试准确率之间的函数关系,而是准确率与模型参数之间的关系。

    from sklearn.learning_curve import validation_curve
    import matplotlib.pyplot as plt
    param_range = [0.001, 0.01, 0.1, 1.0, 10.0, 100.0]
    train_scores, test_scores = validation_curve(estimator=pipe_lr, X=X_train, y=y_train, param_name='clf__C', 
    param_range=param_range, cv=10)
    train_mean = np.mean(train_scores, axis=1)
    train_std = np.std(train_scores, axis=1)
    test_mean = np.mean(test_scores, axis=1)
    test_std = np.std(test_scores, axis=1)
    plt.plot(param_range, train_mean, color='blue', marker='o', markersize=5, label='training accuracy')
    plt.fill_between(param_range, train_mean + train_std, train_mean - train_std, alpha=0.15, color='blue')
    plt.plot(param_range, test_mean, color='green', linestyle='--', marker='s', markersize=5, label='validation accuracy')
    plt.fill_between(param_range, test_mean + test_std, test_mean - test_std, alpha=0.15, color='green')
    plt.grid()
    plt.xscale('log')
    plt.legend(loc='lower right')
    plt.xlabel('Parameter C')
    plt.ylabel('Accuracy')
    plt.ylim([0.8, 1.0])
    plt.show()

    验证曲线

    与learning_curve函数相似,如果使用的是分类算法,则validation_curve函数默认使用分层 k-fold 交叉验证来评价模型的性能。在 validation_curve 函数内,可以指定想要验证的参数。在本例中,需要验证的是参数 C,即定义在 scikit-learn pipeline 中逻辑回归分类器的正则化参数,将其记为 ‘clf__C’,并通过 param_range 参数来设定其值的范围。

    如果加大正则化强度(较小的 C 值)会导致模型轻微的欠拟合;如果增加 C 值,意味着降低正则化的强度,因此模型会趋于过拟合。

    使用网格搜索调优机器学习模型

    在机器学习中,有两类参数:通过训练数据学习得到的参数,如逻辑回归的回归系数;以及学习算法中需要单独进行优化的参数,即超参数。

    网格搜索(grid search)是一种功能强大的超参数优化技巧:它通过寻找最优的超参值的组合以进一步提高模型的性能。

    网格搜索法通过对指定的不同超参列表惊醒暴力穷举搜索,并计算评估每个组合对模型性能的影响,以获得参数的最优组合。

    from sklearn.grid_search import GridSearchCV
    from sklearn.svm import SVC
    pipe_svc = Pipeline([
        ('sc1', StandardScaler()),
        ('clf', SVC(random_state=1))
    ])
    param_range = [0.0001, 0.001, 0.01, 0.1, 1.0, 10.0, 100.0, 1000.0]
    param_grid = [
        {'clf__C': param_range, 'clf__kernel': ['linear']},
        {'clf__C': param_range, 'clf__gamma': param_range, 'clf__kernel': ['rbf']}
    ]
    gs = GridSearchCV(estimator=pipe_svc, param_grid=param_grid, scoring='accuracy', cv=10, n_jobs=-1)
    gs = gs.fit(X_train, y_train)
    print(gs.best_score_)
    print(gs.best_params_)

    执行结果如下:

    0.978021978021978
    {'clf__C': 0.1, 'clf__kernel': 'linear'}

    最后,使用独立的测试数据集,通过 GridSearchCV 对象的 best_estimator_ 属性对最优模型进行性能评估:

    clf = gs.best_estimator_
    clf.fit(X_train, y_train)
    print('Test accuracy: %.3f' % (clf.score(X_test, y_test)))

    执行结果如下:

    Test accuracy: 0.965

    虽然网格搜索是寻找最优参数集合的一种功能强大的方法,但评估所有参数组合的计算成本是相当昂贵的。使用 scikit-learn 抽取不同参数组合的另一种方法是随机搜索,借助于 RandomizedSearchCV 类,可以以特定的代价从抽样分布中抽取随机的参数组合。

    通过嵌套交叉验证选择算法

    如果要在不同机器学习算法中做出选择,则采用嵌套交叉验证。
    在嵌套交叉验证的外部循环中,将数据划分为训练块和测试块;而在用于模型选择的内部循环中,则基于这些训练块使用 k-fold 交叉验证。在完成模型的选择后,测试块用于模型性能的评估。
    嵌套交叉验证
    可以通过如下方式使用嵌套交叉验证:

    gs = GridSearchCV(estimator=pipe_svc, param_grid=param_grid, scoring='accuracy', cv=10, n_jobs=-1)
    scores = cross_val_score(gs, X=X_train, y=y_train, scoring='accuracy', cv=5)
    print('CV accuracy: %.3f +/- %.3f' % (np.mean(scores), np.std(scores)))

    执行结果如下:

    CV accuracy: 0.971 +/- 0.018

    使用嵌套交叉验证方法比较 SVM 模型与简单决策树分类器:

    from sklearn.tree import DecisionTreeClassifier
    gs = GridSearchCV(estimator=DecisionTreeClassifier(random_state=0), param_grid=[
        {'max_depth': [1, 2, 3, 4, 5, 6, 7, None]}
    ], scoring='accuracy', cv=5)
    scores = cross_val_score(gs, X=X_train, y=y_train, cv=5, scoring='accuracy')
    print('CV accuracy: %.3f +/- %.3f' % (np.mean(scores), np.std(scores)))

    执行结果如下:

    CV accuracy: 0.908 +/- 0.045

    由此可见,嵌套交叉验证对 SVM 模型性能的评分远高于决策树,由此,可以预期:SVM 是用于对此数据集未知数据进行分类的一个好的选择。

    不同性能评价指标

    混淆矩阵
    混淆矩阵

    from sklearn.metrics import confusion_matrix
    pipe_svc.fit(X_train, y_train)
    y_pred = pipe_svc.predict(X_test)
    confmat = confusion_matrix(y_true=y_test, y_pred=y_pred)
    print(confmat)

    执行结果如下:

    [[71  1]
     [ 2 40]]

    使用 matplotlib matshow 函数绘制混淆矩阵:

    import matplotlib.pyplot as plt 
    fig, ax = plt.subplots(figsize=(2.5, 2.5))
    ax.matshow(confmat, cmap=plt.cm.Blues, alpha=0.3)
    for i in range(confmat.shape[0]):
        for j in range(confmat.shape[1]):
            ax.text(x=j, y=i, s=confmat[i, j], va='center', ha='center')
    plt.xlabel('predicted label')
    plt.ylabel('true label')
    plt.show()

    混淆矩阵

    误差(ERR)
    预测错误样本的数据与所有被预测样本数量的比值。

    准确率(ACC)
    正确预测样本的数量与所有被预测样本数量的比值。

    精准率(PRE)和召回率(REC)
    PRE与REC

    F1分数
    F1分数

    from sklearn.metrics import precision_score
    from sklearn.metrics import recall_score, f1_score
    print('Precision: %.3f' % (precision_score(y_true=y_test, y_pred=y_pred)))
    print('Recall: %.3f' % (recall_score(y_true=y_test, y_pred=y_pred)))
    print('F1: %.3f' % (f1_score(y_true=y_test, y_pred=y_pred)))

    执行结果如下:

    Precision: 0.976
    Recall: 0.952
    F1: 0.964

    scikit-learn中将正类标识为1, 如果要指定一个不同的正类类标,可通过 make_scorer 函数来构建:

    from sklearn.metrics import make_scorer, f1_score
    scorer = make_scorer(f1_score, pos_label=0)

    ROC曲线

    受试者工作特征曲线(receiver operator characteristic, ROC)是基于模型假正率和真正率等性能指标进行分类模型选择的工具,假正率和真正率可通过移动分类器的分类阈值来计算

    ROC对角线可以理解为随机猜测,如果分类器性能曲线在对角线以下,那么其性能就比随机猜测还差。对于完美分类器来说,其真正率为1,假正率为0,这时的 ROC 曲线即为横轴0与纵轴1组成的折线。

    计算 ROC 线下区域(area under the curve, AUC),用来刻画分类模型的性能。

    from sklearn.metrics import roc_curve, auc
    from scipy import interp
    X_train2 = X_train[:, [4, 14]]
    cv = StratifiedKFold(y_train, n_folds=3, random_state=1)
    fig = plt.figure(figsize=(7, 5))
    mean_tpr = 0.0
    mean_fpr = np.linspace(0, 1, 100)
    all_ptr = []
    for i, (train, test) in enumerate(cv):
        probas = pip_lr.fit(X_train2[train], y_train[train]).predict_proba(X_train2[test])
        fpr, tpr, thresholds= roc_curve(y_train[test], probas[:, 1], pos_label=1)
        mean_tpr += interp(mean_fpr, fpr, tpr)
        mean_tpr[0] = 0.0
        roc_auc = auc(fpr, tpr)
        plt.plot(fpr, tpr, lw=1, label='ROC fold %d (area = %.2f)' % (i + 1, roc_auc))
    plt.plot([0, 1], [0, 1], linestyle='--', color=(0.6, 0.6, 0.6), label='random guessing')
    mean_tpr /= len(cv)
    mean_tpr[-1] = 0.0
    mean_auc = auc(mean_fpr, mean_tpr)
    plt.plot(mean_fpr, mean_tpr, 'k--', label='mean ROC (area = %.2f)' % (mean_auc), lw=2)
    plt.plot([0, 0, 1], [0, 1, 1], lw=2, linestyle=':', color='black', label='perfect performance')
    plt.xlim([-0.05, 1.05])
    plt.ylim([-0.05, 1.05])
    plt.xlabel('false positive rate')
    plt.ylabel('true positive rate')
    plt.title('Receiver  Operator Characteristic')
    plt.legend(loc='lower right')
    plt.show()

    ROC曲线

    多类别分类评价标准

    scikit-learn 实现了 macro 及 micro 均值方法,旨在通过一对多的方式将评分标准扩展到了多分类别分类问题。

    micro 均值是通过系统的真正、真负、假正、假负来计算的。
    micro均值

    macro 均值仅计算不同系统的评价分值:
    macro均值

    当等同看待每个实例或每次预测时,微均值是有用的,而宏均值是等同看待各个类别,将其用于评估分类器针对最频繁类标(样本数量最多的类)的整体性能。

    参考文献:
    Python Machine Learning/Chapter 6: Learning Practices for Model Evaluation and Hyperparameter Tuning

    展开全文
  • 遥感图像分类

    万次阅读 多人点赞 2019-01-28 21:36:02
    训练样本的选取和评估需花费较多的人力,时间等。 非监督分类: 优点:无需对分类区域有广泛的了解,仅需一定的知识来解释分类出集群组;人为误差小;独特的,覆盖量小的类别均能被识别;简单,速度快等。 缺点...
  • 1.模型选择与评估 1.1 模型参数、超参数 了解什么是模型参数,什么是超参数是很必要的,不同类型的参数需要在不同的阶段进行调节(参数——训练集;超参数——验证集) 模型参数:模型内部的配置变量,可以通过数据...
  • 软件测试_笔记(完整版)

    万次阅读 多人点赞 2018-07-02 08:51:28
    广义的软件测试定义:人工或自动地运行或测定某系统的过程,目的在于检验它是否满足规定的需求或弄清预期结果和实际结果间的差别 为什么要做软件测试 发现软件缺陷 功能错 功能遗漏 超出需求部分(画蛇添足) ...
  • 可将参数组合评估为网格(每个参数组合)或列表(在行中尝试几种已定义参数组合)。 PyExperimentSuite还支持继续执行中断执行(例如,电源故障,进程被终止等)时中断的所有实验。 内置的Python界面可以通过多种...
  • 模型评估评估指标

    千次阅读 2019-04-27 16:56:20
    1 模型评价 2 评价指标 2.1 一级指标 序号 一级指标 真实值 模型预测值 说明 1 TP True ...假负,即真实值为正,预测值为负,表示预测错误,若正表示通过,该参数可表示误拒...
  • 线性回归

    万次阅读 多人点赞 2016-06-03 10:22:43
    线性回归模型的基本特性就是:模型是参数的线性函数。最简单的线性回归模型当然是模型是参数的线性函数的同时,也是输入变量的线性函数,或者叫做线性组合。如果我们想要获得更为强大的线性模型,
  • SAP MM 评估类型 评估类别

    千次阅读 2019-04-23 11:26:00
    (2)定义评估类型”(Valuation Type)(如本文中提到“自制品”、“外购品”两种评估类型,并选择相应帐户分类参考参数(与评估类相关),评估类型是评估类别的细分); (3)定义评估类别”(Valuation ...
  • WPF开发教程

    万次阅读 多人点赞 2019-07-02 23:13:20
    第二,如果评估动画系统,您将看到它几乎是完全声明性的。无需要求开发人员计算下一位置或下一颜色,您可以将动画表示为动画对象上的一组属性。于是,这些动画可以表示开发人员或设计人员的意图(在 5 秒内将此按钮...
  • 提出了一个新的三个层次的结构,最高层次为不可估计模型,最低层次为可评估的模型,中间层次为部分评估模型,后者链接最高层次的特征和最低层次的几何模型,并且也提供了最高层次和最低层次几何模型中的特征定义的...
  • SAP MM分割评估

    万次阅读 2017-08-07 16:34:55
    同一物料的使用,既有“自制品”,又有“外购品”,并且其来源不同,如同一外购品由不同的供应商提供,价格也不相同。也就是说:同一物料有不同的价值指派,即在不同的...SAP系统提供的“分割评估”(Split Valuation)
  • 软件测试入门知识了解

    万次阅读 多人点赞 2018-09-05 14:59:58
    质量评估 3.软件测试过程: 需求评审和设计评审是验证软件产品的需求定义和设计实现,验证所定义的产品特性是否符合客户的期望、系统的设计是否合理、是否具有可测试性以及满足非功能质量特性的要求。这个阶段...
  • 模型评估

    千次阅读 2019-05-30 22:24:53
    知道每种评估指标的精确定义、有针对性地选择合适的评估指标、根据评估指标的反馈进行模型调整,这些都是机器学习在模型评估阶段的关键问题。 首先,我们先来了解一下关于模型评估的基础概念。 【误差(error)】:...
  • 【转】评估类型 评估类别 评估

    千次阅读 2020-07-03 15:35:36
    评估类:valuation class  是连接移动类型与财务科目。与物料类型,移动类型一起确定唯一的记账科目。  因为SAP系统的以财务以核心及高度集成等特点,若能在后勤操作的同时自动记账而不是由财务人员手动控制,
  • 模型的选择、评估和优化-上

    千次阅读 2018-07-21 15:20:45
    对于一个模型而言,我们也有很多模型参数需要人工选择,本章将对模型的评估选择和优化进行详细介绍。 概念介绍 过拟合和欠拟合 在机器学习中,我们期望通过训练集来得到在新样本上表现的很好的学习器,找出...
  • 图像质量评估指标 SSIM / PSNR / MSE

    万次阅读 多人点赞 2017-12-09 00:25:21
    图像质量评估指标 SSIM / PSNR / MSE
  • 运维业务相关知识总结

    千次阅读 2019-08-20 17:55:26
    巡检即从各服务器中输出重要的运行日志参数以大概评估服务的整体运行情况,以发现服务隐患;故障处理指的是运维工程师能对服务出现的异常进行及时的处理,当然在服务上线时也需要给出一系列可能的故障处理预案,以在...
  • 建模的评估一般可以分为回归、分类和聚类的评估,本文主要介绍回归和分类的模型评估: 一、回归模型的评估 主要有以下方法: 指标 描述 metrics方法 Mean Absolute Error(MAE) 平均绝对误差 from ...
  • java方法可选参数_Java可选参数

    千次阅读 2020-07-05 23:15:52
    java方法可选参数 在Java类中设计方法时,某些参数对于其执行而言可能是可选的。 无论是在DTO,胖模型域对象还是简单的无状态服务类... 我们将停下来看看Java 8 Optional,并评估它是否符合我们的需求。 让我们开...
  • - 纵轴:真正类率(true postive rate TPR)灵敏度,Sensitivity(正类覆盖率),和recall定义相同。 因此一个好的ROC曲线是横轴还很小的时候,纵轴就很大,即已经预测了n个案例中,把对当对的越多越好,把错当对的越...
  • JedisPool参数说明资源设置与使用相关参数空闲资源检测相关参数关键参数设置建议maxTotal(最大连接数)maxIdle与minIdle使用监控获取合理值常见问题资源不足预热JedisPool Pre Redis Version : 5.0.3 maven依赖:...
  • Pytorch实战总结篇之模型训练、评估与使用

    千次阅读 多人点赞 2020-10-02 09:24:01
    4. 模型的训练设置 这里定义模型的评估指标, 损失函数, 优化方法等。 import datetime import numpy as np import pandas as pd from sklearn.metrics import accuracy_score def accuracy(y_pred,y_true): y_pred...
  • 机器学习之模型评估与选择

    千次阅读 2018-10-06 19:34:16
    模型评估与选择 目录 模型评估与选择 一、经验误差与过拟合 二、评估方法 三、性能度量 四、比较检验   一、经验误差与过拟合 误差:模型实际预测输出与样本的真实输出之间的差异 错误率:分类错误的样本...
  • 根据定义,吸力的绝对大小取决于有效应力参数的大小和矩阵吸力本身。 因此,通过建立地震状态下的兰金关系,检验了无支撑板桩侧面的主动面和被动面的组成方法,并通过类似的方法评估了土性参数对这些面的影响。 ...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 127,182
精华内容 50,872
关键字:

参数评估定义