粒子群-相关向量机.rar

%% PSO-RVM clear clc close all format compact %% load data % data=xlsread('数据(15).xlsx','A2:E219'); % input=data(:,1:4); % output=data(:,end); % save data input output load data %% divide trainset and test set rand('seed',0) [m n]=sort(rand(1,size(input,1))); m=100; x = input(n(1:m),:); t = output(n(1:m),1); xts = input(n(m+1:end),:); tts = output(n(m+1:end),:); data.x=x;data.t=t;data.xts=xts;data.tts=tts; %% d=1;%对径向基核函数的权重alpha寻优 [alpha,trace]=psoforrvm(data,d); disp('优化后高斯核函数的权重为：') alpha figure plot(trace) xlabel('迭代数') ylabel('适应度值/均方差') %% N=size(x,1); % number of train set Nts=size(xts,1); % number of test set %% RVM var1=1000;%核参数 K1=RBFfun1(x,x,var1); %生成K矩阵：N*(N+1) var2=1; K2=PLOYfun(x,x,var2); K=alpha*K1+(1-alpha)*K2; % 随机初始化 alpha（系数w的先验方差倒数） 和 % beta（p(t|w,x,1/beta）=N(t|y(x),beta)=N(t|w*K,1/beta) m=size(K,2); alp=rand(1,m); beta=rand(); ee = zeros(100,1) ; best = inf; for ii=1:100 sig=pinv(diag(alp)+beta*(K'*K)); % 系数w的后验方差矩阵Sigma mu=sig*K'*t*beta; % 系数w的后验均值mu/u gamma=1-alp.*diag(sig)'; alp_old=alp; beta_old=beta; idx=abs(alp)<1e3; % 部分alp会趋向于无穷大，对应的mu会趋向于会向于0，对于这部分alp不再更新 alp(idx)=gamma(idx)./(mu(idx)'.^2); beta=(N-sum(gamma))/((t-K*mu)'*(t-K*mu)); % 判断收敛则退出循环 tmp_err=max(abs(alp(idx)-alp_old(idx))./abs(alp(idx)))+abs(beta-beta_old)/abs(beta); if tmp_err<0.1, break; end; if best>tmp_err best = tmp_err; end ee(ii) = best; end; % 计算并呈现训练集的预测结果 tpred=K(:,idx)*mu(idx); pso_mse_train=mse(tpred-t) % 训练结果 figure subplot(121) hold on; plot(t,'-k*'); plot(tpred,'-ro'); xlabel('训练集样本编号') ylabel('标签类别') title('训练集分类') legend('期望输出','实际输出') %% 测试 xtsK1=RBFfun1(xts,x,var1); xtsK2=PLOYfun(xts,x,var2); xtsK=alpha*xtsK1+(1-alpha)*xtsK2; ttspred=xtsK(:,idx)*mu(idx); % 测试集精度 % 绝对误差 er2=abs(ttspred-tts); % 相对误差 er3=er2./tts*100; result=[tts ttspred]; pso_mse_test=mse(ttspred-tts) subplot(1,2,2) plot(tts,'-*'); hold on; plot(ttspred,'-o'); hold off title('测试图'); xlabel('测试集样本编号') ylabel('标签类别') legend('期望数据','预测数据');