PSEUDO BAYESIAN STABILIZATION ALGORITHM IN MATLAB


/ Published in: MatLab
Save to your folder(s)

PSEUDO BAYESIAN STABILIZATION ALGORITHM IN MATLAB, ESTIMATION OF BACKLOGGED NODES


Copy this code and paste it in your HTML
  1. %Student: Santiago Pagola
  2. %LiU-ID & Personal Number: 930302-T092, sanpa993
  3.  
  4.  
  5. %clear all variables from workspace, close all figures and clear everything
  6. %in command line
  7.  
  8.  
  9. %IMPORTANT NOTE: AN EXTERNAL FUNCTION WAS CREATED IN ORDER TO WRAP THE SUMS
  10. %OF THE VECTORS: a = sumk(vector,k). It will be described in the report,
  11. %but basically sums the first k elements of an array called 'vector'.
  12.  
  13. m = 100; %100 nodes on no-buffering assumption
  14. t = 1000; %Time: 1000 slots of time
  15. backlog_array = zeros(size(1:t)); %Array of the values of the system backlog
  16. backlog_estimate = zeros(size(1:t)); %Array of the estimated values of the system backlog
  17. slots = 1:t; %Slots array
  18. n = 0; %Initially, all nodes are unbacklogged
  19. %lambda = 1/exp(1); %Overal arrival rate will be set to 1/exp(1) (max. stable throughput)
  20. n_estimated = 0; %Estimated backlog initialization. This will be updated during runtime
  21. q_r = 0; %Value of q_r is now stabilized by this algorithm
  22. packets_leaving = 1:t; %Array which will count #packets leaving the system
  23. packets_arriving = 1:t; %Array which will count #packets entering the system
  24. state_probs = zeros(size(1,m)); %State probabilitites array
  25. W = zeros(size(1,7));%Array definition for the P.B. delay
  26. D = zeros(size(1,7));%Array definition for LIttle Theorem delay
  27. v = 1; %Auxiliary variable
  28. %main loop
  29. for lambda2 = 0.05:0.05:0.35
  30. qa = 1-exp(-lambda2/m); %qa is the probability that each unbacklogged node has to transmit in the next slot
  31. for i = 1:t
  32. a = rand(1); %realization of Pr_x
  33. b = rand(1);
  34. c = rand(1);
  35.  
  36. %Estimation of q_r starting from the estimated backlog
  37. if n_estimated >= 0 && n_estimated < 1
  38. q_r = 1;
  39. else
  40. q_r = 1/n_estimated;
  41. end
  42.  
  43. %Now we are going to create an array with the Poisson pmf's taking into
  44. %account the number of transmitting nodes
  45. Pr = zeros(size(1:11));
  46. for j = 0:10
  47. Pr(j+1) = ((lambda2)^j/factorial(j))*exp(-lambda2);
  48. end
  49.  
  50. if 0 <= a && a <= Pr(1)
  51. n;
  52. elseif sumk(Pr,1) < a && a <= sumk(Pr,2)
  53. n = n+1;
  54. elseif sumk(Pr,2) < a && a <= sumk(Pr,3)
  55. n = n+2;
  56. elseif sumk(Pr,3) < a && a <= sumk(Pr,4)
  57. n = n+3;
  58. elseif sumk(Pr,4) < a && a <= sumk(Pr,5)
  59. n = n+4;
  60. elseif sumk(Pr,5) < a && a <= sumk(Pr,6)
  61. n = n+5;
  62. elseif sumk(Pr,6) < a && a <= sumk(Pr,7)
  63. n = n+6;
  64. elseif sumk(Pr,7) < a && a <= sumk(Pr,8)
  65. n = n+7;
  66. elseif sumk(Pr,8) < a && a <= sumk(Pr,9)
  67. n = n+8;
  68. elseif sumk(Pr,9) < a && a <= sumk(Pr,10)
  69. n = n+9;
  70. elseif sumk(Pr,10) < a && a <= 1
  71. n = n+10;
  72. end
  73.  
  74. %Now we are going to create an array with the Qr's (%Probabilities that up to 10 backlogged nodes
  75. %retransmit (up to 10 new arrivals))
  76. Qr = zeros(size(1:11));
  77. for j = 0:10
  78. Qr(j+1) = binopdf(j,n,q_r);
  79. end
  80.  
  81. if 0 <= n && n < m %Check if the backlog is less than m = 100
  82.  
  83. if 0 <= b && b <= sumk(Qr,1) %CASE 1: IDLE SLOT, NO NEW ARRIVALS, FEEDBACK 0
  84. n_estimated=max(lambda2, n_estimated + lambda2 - 1);
  85. packets_arriving(i) = 0;
  86. packets_leaving(i) = 0;
  87.  
  88. elseif sumk(Qr,1) < b && b <= sumk(Qr,2) %CASE 2: NO RETRASMISSION + 1 UNBACKLOGGED NODE TRANSMISION -> SUCCESS! , FEEDBACK 1
  89. n = n-1;
  90. n_estimated = max(lambda2, n_estimated + lambda2 - 1);
  91. packets_arriving(i) = 1;
  92. packets_leaving(i) = 1;
  93.  
  94. elseif sumk(Qr,2) < b && c <= 1 % CASE 3: COLLISION: RETRANSMISSION OF MORE THAN ...
  95.  
  96. if sumk(Qr,2) < b && b <= sumk(Qr,3) %... 2 PACKETS
  97. x = 2;
  98. elseif sumk(Qr,3) < b && b <= sumk(Qr,4) %...3 PACKETS
  99. x = 3;
  100. elseif sumk(Qr,4) < b && b <= sumk(Qr,5) %...AND SO ON
  101. x = 4;
  102. elseif sumk(Qr,5) < b && b <= sumk(Qr,6)
  103. x = 5;
  104. elseif sumk(Qr,6) < b && b <= sumk(Qr,7)
  105. x = 6;
  106. elseif sumk(Qr,7) < b && b <= sumk(Qr,8)
  107. x = 7;
  108. elseif sumk(Qr,8) < b && b <= sumk(Qr,9)
  109. x = 8;
  110. elseif sumk(Qr,9) < b && b <= sumk(Qr,10)
  111. x = 9;
  112. elseif sumk(Qr,10) < b && b <= 1
  113. x = 10;
  114. end
  115.  
  116. %estimation of n for FEEDBACK e
  117. n_estimated = n_estimated + lambda2 + (exp(1)-2)^-1;
  118. packets_arriving(i) = x;
  119. packets_leaving(i) = 0;
  120.  
  121. end
  122. end
  123. backlog_array(i) = n;
  124. backlog_estimate(i) = n_estimated; %Vector that saves the estimated backlog after each slot
  125.  
  126. end
  127. N = mean(backlog_array);
  128. nelements = hist(backlog_array,m); %In each bin, it counts how many times element i was seen ("numerical plot of the hist function")
  129. D(v) = N/lambda2;
  130. W(v) = ((exp(1)-0.5)/(1-lambda2*exp(1)))-(((exp(1)-1)*((exp(1)^lambda2)-1))/(lambda2*(1-((exp(1)-1)*((exp(1)^lambda2)-1)))));
  131. v = v + 1;
  132. end
  133. figure(1) %Setting up the plotting environment for the backlog of the system
  134. xlabel('Slot number, n')
  135. ylabel('Backlogged packets')
  136. plot(slots, backlog_array)
  137. hold on
  138. plot(slots,backlog_estimate,'g')
  139. title('Real backlog (blue) vs. Pseudo Bayesian Estimated backlog (green)')
  140.  
  141. figure(2) %Setting up the plotting environment for the packets of the system
  142. plot(slots,packets_arriving)
  143. hold on
  144. plot(slots,packets_leaving,'k')
  145. grid on
  146. xlabel('Slot number')
  147. ylabel('Number of packets')
  148. title('Number of packets entering (blue) or leaving (black) the system vs. Slot number')
  149.  
  150. %Question 8
  151. max = max(backlog_array);
  152.  
  153. lambda2array = [0.05 0.1 0.15 0.2 0.25 0.3 0.35];
  154.  
  155. plot(lambda2array,D)
  156. xlabel('Variation of lambda')
  157. ylabel('Average delay')
  158. title('Average delay according to Little Theorem and P-B Stabilization')
  159. hold on
  160. plot(lambda2array,W,'r')
  161. legend('Little Theorem','P-B Stabilization')

Report this snippet


Comments

RSS Icon Subscribe to comments

You need to login to post a comment.