基于粒子群算法的三态马尔可夫过程数据包服务模型

根据所描述的数据包到达和服务过程，我们可以将MMOO改写为三态的马尔可夫过程，其中g为1态的概率，m为2态的概率，b为3态的概率。依据概率转移矩阵，我们可以将数据包的服务过程表示为：/n/n$$/begin{bmatrix}g&m&b/end{bmatrix}/begin{bmatrix}1-p&p&(1-q)//a*(1-p)&ap&q(1-a)//0&0&1-b/end{bmatrix}^{n}/begin{bmatrix}1//0//0/end{bmatrix}$$其中，$n$表示数据包的离去次数，即数据包的服务次数。根据这个公式，我们可以将'Threshold_MMOO_01'函数改写为'Threshold_Markov_3state'函数，代码如下：matlabfunction [threshold03]=Threshold_Markov_3state(p,q,a,Dmax,x,y)% g: 1态的概率，m: 2态的概率，b: 3态的概率theta0 = [0.1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5];f = @(theta) largedeviationMarkov_3state(p,q,a,x,y,theta);for j = 1:length(theta0) theta_tmp(j) = fsolve(f,[theta0(j)],optimset('Display','off')); end max=theta_tmp(1);for i=2:length(theta0) if (theta_tmp(i)>max) max=theta_tmp(i); endendtheta1=max;% 计算1态和2态的概率T=[1-p p (1-q);a*(1-p) ap q(1-a);0 0 1-b];[T_eigenvector,T_eigenvalue] = eig(T);[v,w] = size(T_eigenvalue);max_eigenvalue=T_eigenvalue(1,1);max_index=1;for i=2:v if (T_eigenvalue(i,i)>max_eigenvalue) max_eigenvalue=T_eigenvalue(i,i); max_index=i; endendstate_prob=T_eigenvector(:,max_index);% 计算阈值g=state_prob(1);m=state_prob(2);b=state_prob(3);threshold03=(gmexp(-theta1Dmax))/(1-g-b+gexp(-theta1Dmax));end其中，'largedeviationMarkov_3state'函数的代码如下：matlabfunction [fx]=largedeviationMarkov_3state(p,q,a,x,y,theta) % 计算三态马尔可夫过程中的大偏差概率T=[1-p p (1-q);a(1-p) ap q(1-a);0 0 1];T_theta=T.exp(theta);[T_theta_eigenvector,T_theta_eigenvalue] = eig(T_theta);[v,w] = size(T_theta_eigenvalue);max_eigenvalue=T_theta_eigenvalue(1,1);max_index=1;for i=2:v if (T_theta_eigenvalue(i,i)>max_eigenvalue) max_eigenvalue=T_theta_eigenvalue(i,i); max_index=i; endendstate_prob=T_theta_eigenvector(:,max_index);fx=(state_prob(1)state_prob(2)exp(-thetaDmax))/(1-state_prob(1)-state_prob(3)+state_prob(1)exp(-thetaDmax))-1; end 而主函数的代码则需要根据新的'Threshold_Markov_3state'函数进行修改，具体修改如下：matlabnum_particle=10; c1=2;c2=2;for j = 1:1:M arr_rate(1,j)=a(1-q)/(2-p-q); v=zeros(num_particle,3); position=zeros(num_particle,3); pbestf=zeros(num_particle,1); fitness=zeros(num_particle,1); record=zeros(itermax,1); position(:,1)=rand(num_particle,1); position(:,2)=rand(num_particle,1)(1-position(:,1)); position(:,3)=1-position(:,1)-position(:,2); echo=0; pbestp=position; for i=1:1:num_particle pbestf(i,1)=Threshold_Markov_3state(p,q,a,Dmax,position(i,1),position(i,2)); end gbestf(1,j)=pbestf(1,1); id=1; for i=2:1:num_particle if abs(threshold-pbestf(i,1))<abs(threshold-gbestf(1,j)) gbestf(1,j)=pbestf(i,1); id=i; end end gbestp=position(id,:); for echo=1:1:itermax echo w=0.9-0.5*echo/itermax; for i=1:1:num_particle v(i,:)=w.*v(i,:)+c1.rand().(gbestp-position(i,:))+c2.rand().(pbestp(i,:)-position(i,:)); position(i,:)=position(i,:)+v(i,:); if position(i,1)<0 position(i,1)=0; elseif position(i,1)>1 position(i,1)=1; end if position(i,2)<0 position(i,2)=0; elseif position(i,2)>1-position(i,1) position(i,2)=1-position(i,1); end position(i,3)=1-position(i,1)-position(i,2); fitness(i,1)=Threshold_Markov_3state(p,q,a,Dmax,position(i,1),position(i,2)); if abs(threshold-fitness(i,1))<abs(threshold-pbestf(i,1)) pbestf(i,1)=fitness(i,1); pbestp(i,:)=position(i,:); end if abs(threshold-fitness(i,1))<abs(threshold-gbestf(1,j)) gbestf(1,j)=fitness(i,1); gbestp=position(i,:); end end record(echo)=gbestf(1,j); end x(1,j)=gbestp(1,1); y(1,j)=gbestp(1,2); g(1,j)=1-x(1,j)-y(1,j); arr_ser(1,j)=(g(1,j)+m(1,j))/(2-x(1,j)-y(1,j)); rho=arr_rate(1,j)/arr_ser(1,j);end其中，我们新增了g、m、b三个参数，分别表示三态马尔可夫过程中1态、2态、3态的概率。在粒子群算法中，初始的'position'也需要修改为三个参数的随机分布，并且更新时需要分别更新三个参数。最后，我们计算出了参数g、m、b和数据包到达服从IBP过程得到的参数p、q、a所对应的数据包服务概率'arr_ser'，并计算出了负