Стохастический образец с использованием случайного образца?

Question

Я хотел бы сделать стохастическое моделирование шага для P и I с randsample, как показано ниже.

P=zeros(1,5); I=zeros(1,5)

% легкий путь

for i=1:5
X=rand; dt=0.01; 

a=randi(50,1);  
b=randi(50,1);
c=randi(50,1);
d=randi(50,1);



if X<=a*dt, 
   P(i+1)=P(i+1)+1;
elseif X>a*dt && X<=(a+b)*dt
   P(i+1)=P(i)-1;
elseif X>(a+b)*dt && X<=(a+b+c)*dt
   I(i+1)=I(i)-1;
elseif X>(a+b+c)*dt && X<=(a+b+c+d)*dt
   I(i+1)=I(i)+1;       
else  %do nothing
  P(i+1)=P(i);
  I(i+1)=P(i);
end

% Использование randsample

    Pvec=[a b c d (some value for doing nothing)]*dt;
    Pvec=Pvec./sum(Pvec);
    s=randsample(1:5,1,'true',Pvec);

Это не правильно.Как бы вы сделали это эффективно?

Это то, что я пытаюсь сделать, но я не думаю, что это совершенно правильно ...

ОБНОВЛЕНИЕ с конкурирующими группами I и P code На основе этого набора уравнений.

theta_P=0.15;delta_P=0.01;alpha_I=0.4;gamma_I=0.01;delta_I=0.005;lambda_I=0.05;
m=100;   % # runs
time=10; % # Total time of simulation
dt=0.01; % # Time step
D=6000; T=10/dt;

P=zeros(m,time/dt); I=zeros(m,time/dt);

for i=1:m
 for j=1:time/dt         
  arrivalI=alpha_I+P(i,j)*lambda_I;
  lossI=I(i,j)*gamma_I+P(i,j)*I(i,j)*delta_I;

     if j<=T
        alpha_P=D/T;
     else
        alpha_P=0;
     end

  arrivalP=alpha_P+P(i,j)*theta_P;
  lossP=P(i,j)*I(i,j)*delta_P;


  X=rand;

Pvec=[arrivalI lossI arrivalP lossP]*dt;%
Pvec=Pvec./sum(Pvec);

s=randsample(1:4,1,'true',Pvec);


if s==1
    I(i,j+1)=I(i,j)+1;%;
    P(i,j+1)=P(i,j);
elseif s==2
    I(i,j+1)=I(i,j)-1;%
    P(i,j+1)=P(i,j);
elseif s==3

    P(i,j+1)=P(i,j)+1;%
    I(i,j+1)=I(i,j);
elseif s==4

    P(i,j+1)=P(i,j)-1;%;
    I(i,j+1)=I(i,j);
else

    P(i,j+1)=P(i,j);  %check
    I(i,j+1)=I(i,j);
end


end

   subplot(2,2,1:2)
%     
   if P(i,j)>5
   loglog(abs(P(i,:)),'-r')
%    
   else
    loglog(abs(P(i,:)),'-b')
%          
   end
  hold on 
  axis([1 1e3 1 1e4])
end

Answer 1

Я не думаю, что вы можете скопировать свой первый блок кода «простым способом» с помощью вызова randsample.

Первый блок кода генерирует P и I рекурсивно.

В то время как randsample создает выборки с или без замены совокупности: в данном случае 1:5.

Вы можете попробовать randseq (требуется Биоинформатика Инструментарий).С точки зрения эффективности, как правило, нет способа векторизовать рекурсивную операцию.