anfis matlab

clear;close all;
gamma=0.75;%设定惯性因子

eps1=0.005;%设定停止训练的条件参数
m1=8;%设定隶属函数个数
m2=8;

a=-1;b=1;
w0=a+（b-a）*rand（1m1*m2）;%初始化权值阵

for i=1:2
    switch i
        case 1beta=0.75;%设定学习率
        otherwisebeta=0.25;
    end
    c=[2/7*（0:m1-1）-1;2/7*（0:m2-1）-1];%初始化Cij
    sigma=0.1213*ones（2m1）;%初始化σij　
    w=w0;
    mu=zeros（2m1）;%初始化隶属度层
    alpha=zeros（1m1*m2）;%初始化规则层
    alpha_=zeros（1m1*m2）;%初始化归一化层

    delta2=zeros（2m1）;%初始化δ2
    dw=zeros（1m1*m2）;%初始化权值变化量
    dc=zeros（2m1）;%初始化c变化量
    dsigma=zeros（2m1）;%初始化σij变化量

%----------------训练过程----------------%
    err=1;%初始化误差
    er=[];
    counter=0;%统计训练步数
    while（err>=eps1）
        Par_E_w=zeros（1m1*m2）;%误差对权值阵的偏导数
        Par_E_c=zeros（2m1）;%误差对c的偏导数
        Par_E_sigma=zeros（2m1）;%误差对σ的偏导数
        E=0;
        for x1=-1:2/19:1
            for x2=-1:2/19:1
                yd=sin（pi*x1）*cos（pi*x2）;%期望输出，导师信号
                %正向传播
                mu（1:）=exp（-（x1-c（1:））.^2./sigma（1:）.^2）;%计算隶属度层
                mu（2:）=exp（-（x2-c（2:））.^2./sigma（2:）.^2）;
                s=zeros（2m1m1*m2）;%初始化Sij
                for m=1:m1
                    for n=1:m2
                        alpha（（m-1）*m2+n）=min（mu（1m）mu（2n））;%计算规则层（取小运算）
                        if mu（1m）<=mu（2n）
                            s（1m（m-1）*m2+n）=1;
                        end
                        if mu（1m）>=mu（2n）
                            s（2n（m-1）*m2+n）=1;
                        end
                    end
                end
                alpha_=alpha/sum（alpha）;%计算归一化层
            
                y=alpha_*w.‘;%计算网络输出
                E=E+1/2*（yd-y）^2;%计算误差
                %反向传播
                delta5=yd-y;
                delta4=delta5*w;

                for k=1:m1*m2
                    delta3（k）=delta4（k）*（sum（alpha）-alpha（k））./sum（alpha）^2;
                end
                for m=1:2
                    for n=1:m1
                        delta2（mn）=0;
                        for l=1:m1*m2
                            delta2（mn）=delta2（mn）+delta3（l）*s（mnl）*mu（mn）;
                        end
                    end
                end
            
                Par_E_w=Par_E_w-delta5*alpha_;%计算偏导数
                Par_E_c（1:）=Par_E_c（1:）-2*delta2（1:）.*（x1-c（1:））./sigma（1:）.^2;
                Par_E_c（2:）=Par_E_c（2:）-2*delta2（2:）.*（x2-c（2:））./sigma（2:）.^2;
                Par_E_sigma（1:）=Par_E_sigma（1:）-2*delta2（1:）.*（x1-c（1:））.^2./sigma（1:）.^3;
                Par_E_sigma（2:）=Par_E_sigma（2:）-2*delta2（2:）.*（x2-c（2:））.^2./sigma（2:）.^3;
            end
        end
  
        num=20*20;
        Par_E_w=Par_E_w/num;
        Par_E_c=Par_E_c/num;
        Par_E_si