Advertisement
sashachca

lab2

May 21st, 2018
166
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Octave 11.34 KB | None | 0 0
  1. function MC_LR_2_KHMELEV
  2.   V = 17;
  3.   A = -1.7000;
  4.   B = 1;
  5.   SIGMA =  1.1700;
  6.   LAMBDA =  0.83000;
  7.  
  8.   P = X = W = XW = [];
  9.  
  10.   printf(' Variant is %i\n [a, b] = [%i, %i],\n sigma = %i, lambda = %i\n', V, A, B, SIGMA, LAMBDA);
  11.   for num = 1:3
  12.     XW = task(num);
  13.     X = XW(1, :);
  14.     W = XW(2, :);
  15.     W(end) = [];
  16.     len = length(W);
  17.    
  18.     % Расчеты для нормального распределения:
  19.     if num == 1
  20.       P = normcdf(X(2:end), A, SIGMA) - normcdf(X(1:end-1), A, SIGMA);
  21.     elseif num == 2
  22.     % Расчеты для показательного распределения:
  23.       P = expcdf(X(2:end), 1/LAMBDA) - expcdf(X(1:end-1), 1/LAMBDA);
  24.     elseif num == 3
  25.     % Расчеты для равномерного распределения:
  26.       P = unifcdf(X(2:end), A, B) - unifcdf(X(1:end-1),A, B);
  27.     endif
  28.     sum(P)
  29.    
  30.     WP = abs(W - P);
  31.    
  32.     printf('pi = { ');
  33.     for i = 1:len-1
  34.       printf('%i, ', myround(P(i), 5));
  35.     endfor
  36.     printf('%i }\n', myround(P(len), 5));
  37.    
  38.     printf('|wi - pi| = { ');
  39.     for i = 1:len-1
  40.       printf('%i, ', myround(WP(i), 5));
  41.     endfor
  42.     printf('%i }\n', myround(WP(len), 5));
  43.    
  44.     printf('max{|wi - pi|} = %i\n', myround(max(WP), 5));
  45.   endfor
  46.  
  47. endfunction
  48.  
  49.  
  50. function AW = task(number)
  51.  
  52.   printf('\n----------\n  TASK NUMBER %i\n', number);
  53.   N = 200;
  54.   a = n = x = [];
  55.  
  56.   xx = sort(oneOfArrays(number));
  57.   m = 1 + fix(log2(N)); % число интервалов по формуле Стерджеса (fix отбрасывает дробную часть)
  58.  
  59.   a(1) = xx(1); %а0
  60.   if number == 2
  61.     a(1) = 0; % для показательного распределения
  62.   endif
  63.   a(m+1) = xx(N); %аm
  64.   d = a(m+1) - a(1);
  65.   h = d/m; % поправка Шеппарда
  66.  
  67.   for k = 1:m-1
  68.     a(k+1) = d/m + a(k); %от а0 до аm или от а1 до а(m+1)
  69.   endfor
  70.  
  71.   [n, x] = hist(xx, m);
  72.   w = n/N;
  73.   printf('\n\n summa wi = %i\n\n\n', sum(w));
  74.  
  75.   printf('Interval series:\n\n[');
  76.   for i = 1:m
  77.     printf('%i, %i]\n n = %i, w = %i\n\n', a(i), a(i+1), n(i), w(i));
  78.     if i == m
  79.       printf('\n');
  80.       break
  81.     endif
  82.     printf('(');
  83.   endfor
  84.  
  85.   figure('Name','Emperical distribution function');
  86.   draw_edf(xx, N);
  87.  
  88.   figure('Name','Histogram of relative frequencies');
  89.   bar(x, w/h, 'hist');
  90.  
  91.  
  92.   X = 0; % выборочное среднее
  93.   X = sampleMoment(x, w, 1, X)
  94.  
  95.   D = 0; % выборочная дисперсия
  96.   D = sampleMoment(x, w, 2, X);
  97.   S = 0; % выборочная дисперсия c поправкой Шеппарда
  98.   S = D - h^2/12
  99.  
  100.   sigma = 0; %выборочное среднее квадратическое отклонение
  101.   sigma = sqrt(S)
  102.  
  103.   M0 = 0; %выборочная мода
  104.   M0 = Mode(a, w, m, h)
  105.  
  106.   Me = 0; %выборочная медиана
  107.   Me = Median(a, w, m, h)
  108.  
  109.   gamma1 = gamma2 = 0; % выборочные коэффициенты
  110.   if sigma ~= 0
  111.     num = sampleMoment(x, w, 3, X);
  112.     denum = sigma^3;
  113.     gamma1 = num/denum % выборочный коэффициент асимметрии
  114.    
  115.     num = sampleMoment(x, w, 4, X);
  116.     denum = sigma^4;
  117.     gamma2 = num/denum - 3 % выборочный коэффициент эксцесса
  118.   endif
  119.   w = [w, 0]; % для того, шоб можно было сделать матрицу с одинаковым количеством элементов в строке
  120.   AW = [a; w];
  121.  
  122. endfunction
  123.  
  124.  
  125. function Mu0 = sampleMoment(x, w, k, X) % выборочный (если X != 0, то еще и центральный) момент k-ого порядка
  126.   Mu0 = 0;
  127.  
  128.   Mu0 = sum(((x.-X).^k).*w);
  129.  
  130. endfunction
  131.  
  132.  
  133. function M0 = Mode(a, w, m, h) % выборочная мода
  134.   M0 = 0;
  135.   w_prev = 0;
  136.   w_next = 0;
  137.   [w_max, i] = max(w);
  138.  
  139.   if (i ~= m)
  140.     w_next = w(i+1);
  141.   endif
  142.   if (i ~= 1)
  143.     w_prev = w(i-1);
  144.   endif
  145.   M0 = a(i)+h*(w(i) - w_prev)/(2*w(i) - w_prev - w_next);
  146.  
  147. endfunction
  148.  
  149. function Me = Median(a, w, m, h) % выборочная медиана
  150.   Me = 0;
  151.  
  152.   sumw = 0; i = 1;
  153.   while (i <= m && sumw + w(i) < 0.5)
  154.     sumw += w(i);
  155.     i = i+1;
  156.   endwhile
  157.   Me = a(i+1);
  158.   if (sumw + w(i) > 0.5)
  159.     Me = a(i) + h/w(i)*(0.5 - sumw);
  160.   endif
  161.  
  162. endfunction
  163.  
  164.  
  165. function F = edf(y, xx, N) % emperical distribution function (эмпирическая функция распределения)
  166.   F = 0;
  167.  
  168.   for i = 1:N
  169.     if xx(i) <= y
  170.       F += 1/N;
  171.     endif
  172.   endfor
  173.  
  174. endfunction
  175.  
  176.  
  177. function draw_edf(xx, N)
  178.  
  179.   interval = xx(1)-1:0.01:xx(1);
  180.   line(interval, 0, "linewidth", 1.75);
  181.   for i = 2:N-1
  182.     F = edf(xx(i), xx, N);
  183.     interval = xx(i-1):0.01:xx(i);
  184.     line(interval, F, "linewidth", 1.75);
  185.   endfor
  186.   interval = xx(200):0.01:xx(200)+1;
  187.   line(interval, 1, "linewidth", 1.75);
  188.  
  189. endfunction
  190.  
  191.  
  192. function arr = oneOfArrays(a)
  193.   arr = [];
  194.  
  195.   if a == 1
  196.     arr = [ 0.945283, -2.06617, -2.71212, -1.95723, -2.82998, -1.81279, -2.07176, -0.710579, -1.47213, -1.48979, -0.517603, -0.855957, -1.65863, -2.84071, -3.32287, -1.60892, 0.147083, -3.61568, -1.86726, -1.08543, -3.14257, 0.31785, -0.473479, -0.345501, -2.80987, -1.88648, -1.74145, -1.92354, -4.11625, -3.24455, -3.23523, -0.70695, -0.79363, -2.14492, -2.58918, -3.96434, -1.68702, -0.134764, -1.12666, -0.269525, -2.75434, -2.6853, -1.99038, -2.17781, -3.15076, -1.58966, -2.17831, -2.23366, -0.0994655, -0.172237, -3.12701, -1.50827, 0.85855, -1.52703, -2.27682, -0.990274, -2.58586, -2.46557, -2.98211, 1.16515, -1.68262, -1.37891, -1.58473, -2.86909, -1.57624, -1.96669, -2.09494, -1.84911, -1.82433, -2.63872, -2.92327, 1.58335, -1.90007, 1.75789, -0.409697, -2.47554, -2.27445, -1.14097, -1.52075, -2.75879, -1.52013, -1.48938, -2.58281, -1.51727, -2.50423, -2.38588, -2.12786, -0.836045, -0.587623, -2.7854, -3.10065, -1.42406, -0.809771, -2.17878, -1.17799, -2.58562, -0.657522, -2.46965, -0.171438, -0.81828, -2.47422, -0.69163, -2.01337, -3.65011, -1.7803, -0.571835, -1.46896, -1.56939, -0.848872, -3.48812, 1.88424, -2.29347, 0.0223165, 1.4784, -3.64157, -3.49172, -1.7499, 0.26698, -2.71308, -1.71597, -1.6796, -2.62463, -2.3972, -1.72866, 0.694823, -2.52757, 0.869767, -3.05757, -0.678853, -1.39321, -4.2647, -4.11425, -1.20893, 0.606871, -2.52303, -1.43166, -1.83489, 0.581549, -0.556458, 1.17557, 0.503682, -1.19899, -1.84798, -1.63378, -2.82692, -2.75873, -1.97005, -0.577092, -1.72913, -0.661957, -2.35756, -3.52809, 0.576527, 0.567181, -2.08436, -2.39206, -1.84209, -0.575595, -1.15436, -1.81609, -0.992812, -0.536021, -1.97943, 0.480746, -2.01, -2.76832, -0.889094, -0.0145241, -3.1347, 0.638841, -0.784655, -2.10318, -1.50043, -1.52127, -1.2691, -1.81336, -2.74042, -1.64962, -0.961865, -0.767931, -1.55012, -3.14284, -3.12792, -2.37715, -3.36592, -1.07188, -1.47116, -2.55522, -2.86124, -0.737448, -2.13386, -1.29757, -1.39216, 0.286887, -2.22566, -2.11048, 0.437069, 0.069697, -0.759136, -1.01348 ];
  197.  
  198.   elseif a == 2    
  199.     arr = [ 0.70071, 0.552521, 2.42469, 0.385447, 2.2528, 0.961222, 0.132513, 1.09588, 2.0704, 0.38289, 0.00361903, 0.396812, 0.999958, 0.66701, 0.0678301, 0.834677, 2.49263, 0.405059, 0.776688, 0.825137, 1.38137, 2.38138, 0.477949, 1.35055, 0.778213, 1.25251, 1.91617, 0.273711, 1.17779, 0.777776, 0.25222, 0.40794, 0.94598, 0.00894616, 0.215276, 1.20029, 1.044, 3.13446, 1.76611, 0.568937, 3.27171, 0.017784, 2.92152, 0.306144, 0.440532, 0.603799, 1.13209, 1.06524, 0.110008, 0.0320471, 0.607945, 3.36284, 0.145036, 1.5464, 0.2367, 1.01822, 1.47319, 0.931407, 1.55182, 0.75239, 1.12243, 1.26739, 0.619438, 2.12079, 1.06426, 0.129365, 0.27428, 2.05554, 0.525062, 2.26446, 1.10673, 7.7782, 1.59408, 0.484774, 0.12352, 0.240475, 2.34159, 0.206745, 1.37961, 3.5023, 0.252737, 0.783518, 0.597136, 1.22923, 1.47077, 0.45055, 0.491048, 1.2428, 0.38804, 0.07631, 0.367377, 2.97136, 0.686305, 0.0963505, 1.40654, 0.618636, 0.241352, 0.72694, 0.34762, 0.115897, 0.277737, 3.70547, 1.82634, 0.517682, 0.694946, 0.120539, 2.31001, 3.62352, 0.870685, 4.28428, 0.0144237,0.018017, 0.485097, 0.72259, 1.66115, 0.251356, 0.366446, 4.57314, 1.7078, 1.83989, 0.478202, 0.112346, 3.48549, 1.08826, 0.150269, 0.325178, 1.0989, 0.719387, 2.29126, 1.82819, 1.18646, 1.67413, 1.21658, 0.0394879, 0.446564, 1.7921, 1.24164, 3.79219, 0.636219, 0.228781, 1.24732, 0.5403, 0.123129, 1.03522, 0.441509, 0.0556123, 0.330339, 0.274696, 1.89326, 0.437623, 3.16463, 1.07979, 0.0154662, 0.136733, 4.67929, 1.47281, 0.379262, 0.953914, 0.254726, 0.1286, 0.591864, 0.361418, 0.354341, 0.00280567, 1.65146, 1.42142, 3.61803, 0.613403, 1.82167, 0.0300214, 0.53848, 0.142309, 1.00385, 0.200341, 0.456653, 1.74318, 0.01088, 0.348768, 0.801765, 0.0162816, 1.68256, 0.0814346, 2.84202, 4.53684, 0.0215211, 0.195153, 0.133513, 0.0401515, 0.914902, 0.402067, 0.396135, 0.211113, 0.35011, 1.33248, 4.22435, 2.24393, 2.55864, 0.5586, 1.39324, 1.21052 ];
  200.  
  201.   elseif a == 3
  202.     arr = [ 0.476694, 0.0843134, -0.700786, 0.75856, -1.59373, -0.0624069, -0.633571, -0.368289, -0.518817, 0.854694, -1.26846, 0.462822, -1.66012, -0.131748, -1.01056, -0.130174, -1.56153, 0.283641, -1.46425, -0.856705, 0.00259962, -1.26732, -1.38335, 0.824514, -1.03893, 0.801303, -1.25802, -1.53605, -0.439681, -1.33648, -0.0665551, 0.96014, -0.486979, 0.971777, -1.56244, -0.58461, -1.60873, 0.77675, -0.400804, 0.61565, 0.216919, 0.463497, -0.137146,0.881553, -1.2856, 0.583729, -0.155113, -0.00713341, -1.12784, -1.07571, 0.289538, -0.521984, -0.31744, -0.727761, 0.147277, -0.248653, -1.67218, 0.315561, -0.643423, -1.22289, -0.188112, -0.951154, 0.657705, 0.712886, -0.390987, -1.2142, -1.44421, 0.454537, 0.875687, -1.20751, 0.731493, -1.1443, 0.39352, -1.63649, -0.524572, -1.55873, -1.546, -0.367047, 0.334602, -0.215564, -0.0881068, 0.555007, -1.59121, 0.0851864, 0.642726, 0.806017, 0.638842, -1.61453, -0.69139, -1.03321, -1.41483, -0.311866, 0.214791, -0.0981157, 0.361791, 0.977283, 0.318551, 0.280426, -1.1114, -0.747861, 0.851864, -1.01675,-1.20305, -0.926762, -1.26025, -0.00148996, -0.94144, -0.178128, 0.360537, -1.28936, -1.2505, -0.453266, -1.06159, 0.554422, -0.673038, -0.686217, -0.377933, 0.239516, 0.343865, 0.743257, -0.114357, -1.47893, 0.747183, 0.269477, -0.71139, 0.983945, -1.5337, -1.69058, 0.107548, -1.00136, -0.407417, -0.334595, -0.118671, -1.09543, -0.18074, -0.575652, 0.514438, 0.62945, -0.733271, 0.00731516, -1.60134, -0.157404, -1.23824, -0.400291, -1.06007, 0.96623, -0.261657, -0.243722, -0.626704, -0.474552, -0.21428, -0.176189, -1.56408, 0.667865, 0.248116, 0.935999, -1.52654, -1.17919, -0.0167115, 0.109524, -0.764633, -1.45357, -0.453395, 0.892928, 0.907744, -0.936402, 0.310306, -0.227736, -0.153756, 0.162371, -1.65934, 0.781652, -0.171568, -0.55599, 0.121784, 0.0579806, -1.47074, -0.00381734, -1.08758, 0.316077, -1.16503, 0.187724, 0.402, -1.02803, -1.50774, 0.679789, -1.03947, 0.520825, -1.01154, -0.111421, -1.20447, 0.0137564, -0.732928, -1.69531, -1.09774, 0.930254, -0.0619198, -0.696017, -1.13844, -0.520925 ];
  203.  
  204.   endif
  205.    
  206. endfunction
  207.  
  208. function R = myround(x, k)
  209. R = round(x*10^k)/10^k;
  210.  
  211. endfunction
  212.  
  213. %  n = x = zeros(1, m);
  214. %    
  215. %  jj = 1;
  216. %  for i = 1:m
  217. %    %printf('\na(%i) = %i\n', i, a(i));
  218. %    for j = jj:N
  219. %      if xx(j) >= a(i) && xx(j) <= a(i+1)
  220. %        %xx(j)
  221. %        n(i) += 1;
  222. %        jj = j;
  223. %      endif
  224. %    endfor
  225. %    %printf('a(%i) = %i\n', i+1, a(i+1));
  226. %    x(i) = (a(i) + a(i+1))/2;
  227. %  endfor
  228.  
  229. %  figure('Name','Histogram of frequencies');
  230. %  hist(xx, m);
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement