Octave и мат.статистика

# Функция для построения эмпирической кривой
# плотности распределения:
function y = fEmp(x, Sx2, Xsr, n)
  tmp = 0;
  for i = 1:9
    tmp += (fb(sqrt(Sx2)*0.8*(i-4)+Xsr, n, x) - fb(sqrt(Sx2)*0.8*(i-5)+Xsr, n, x))*(xb(x-(sqrt(Sx2)*0.8*(i-5)+Xsr), n)-xb(x-(sqrt(Sx2)*0.8*(i-4)+Xsr), n))*(n/(0.8*sqrt(Sx2)));
  endfor
  y = tmp;
endfunction

-- Function File:  normpdf (X, M, S)
     For each element of X, compute the probability density function
     (PDF) at X of the normal distribution with mean M and standard
     deviation S.
     Default values are M = 0, S = 1

# TViMS.m - скрипт для выполнения курсовой работы
#
# Выполнил: Вытовтов П.К.
#           гр. № 4-78-1
#
# Среда разработки: GNU Octave, version 3.2.4
# -------------- [Первое задание] --------------
# Для выборок X, Y найти оценки: выборочное среднее
# выборочная дисперсия, оценка центрального момента
# порядка 3 и 4, коэффициент эксцесса, коэффициент
# ассимитрии, оценка корреляционного момента, оценка
# коэффициента корреляции, размах выборки X и Y
n = 180; # Объем выборки
m = 0; # Коэффициент сдвига
s = 1; # Коэффициент масштаба
a = 0.05; # Уровень значимости
# Формируем выборки с нормальным распределением:
X = normrnd(m, s, n, 1);
Y = normrnd(m, s, n, 1);
# Находим выборочное среднее:
Xsr = mean(X)
Ysr = mean(Y)
# Находим выборочную дисперсию:
Sx2 = 0;
Sy2 = 0;
for i = 1:n
  Sx2 += X(i) - Xsr;
  Sy2 += Y(i) - Ysr;
endfor
Sx2 /= n - 1
Sy2 /= n - 1
# Находим оценку центрального момента порядка 3 и 4:
M3X = 0;
M3Y = 0;
M4X = 0;
M4Y = 0;
for i = 1:n
  M3X += (X(i) - Xsr)^3;
  M3Y += (Y(i) - Ysr)^3;
  M4X += (X(i) - Xsr)^4;
  M4Y += (Y(i) - Ysr)^4;
endfor
M3X /= n - 1
M3Y /= n - 1
M4X /= n - 1
M4Y /= n - 1
# Находим коэффициент эксцесса:
Ex = (M4X / Sx2^2) - 3
Ey = (M4Y / Sy2^2) - 3
# Находим коэффициент ассимитрии:
Sx = M3X / Sx2^(3/2)
Sy = M3Y / Sy2^(3/2)
# Находим оценку корреляционного момента:
Rxy = 0;
for i = 1:n
  Rxy += (X(i) - Xsr) * (Y(i) - Ysr);
endfor
Rxy /= n - 1
# Находим оценку коэффициента корреляции:
ro = Rxy / (Sx2^(1/2) * Sy2^(1/2))
# Сортируем для нахождения размаха:
for i = 1:n
  for j = 1:n
    if(X(i)<X(j))
      tmp = X(i);
      X(i) = X(j);
      X(j) = tmp;
    endif
    if(Y(i)<Y(j))
      tmp = Y(i);
      Y(i) = Y(j);
      Y(j) = tmp;
    endif
  endfor
endfor
# Находим размах выборки:
wX = X(n) - X(1)
wY = Y(n) - Y(1)
# ----------------------------------------------
# -------------- [Второе задание] --------------
# Для оценок по выбоке X построить доверительные
# доверительные интервалы для математического
# ожидания (считая дисперсию известной и
# неизвестной) и дисперсии
# Находим доверительный интервал для математического ожидания
# в случае с известной дисперсией для выборки X:
dims1 = Xsr - ((Sx2^(1/2) / n^(1/2)) * normcdf(1-a/2, m, s))
dims2 = Xsr + ((Sx2^(1/2) / n^(1/2)) * normcdf(1-a/2, m, s))
# Находим доверительный интервал для математического ожидания
# в случае с неизвестной дисперсией для выборки X:
dimns1 = Xsr - ((Sx2^(1/2) / n^(1/2)) * tcdf(1-a/2, n-1))
dimns2 = Xsr + ((Sx2^(1/2) / n^(1/2)) * tcdf(1-a/2, n-1))
# Находим доверительный интервал для дисперсии
# для выборки X:
dis1 = (n-1)*Sx2/chi2cdf(1-a/2, n-1)
dis2 = (n-1)*Sx2/chi2cdf(a/2, n-1)
# ----------------------------------------------
# -------------- [Третье задание] --------------
# По выборкам X, Y построить доверительный интервал
# для коэффициента корреляции
# Определяем доверительный интервал по коэффициенту
# корреляции для выборок X, Y:
z = (1/2)*log((1+Rxy)/(1-Rxy));
z1 = z - tcdf(1-a/2, n) * (1/(n-3)^(1/2));
z2 = z + tcdf(1-a/2, n) * (1/(n-3)^(1/2));
ro1 = tanh(z1)
ro2 = tanh(z2)
# ----------------------------------------------
# -------------- [Четвертое задание] --------------
# Для X построить эмпирическую кривую плотности
# распределения, разбив интервал (X(1), X(n)) на 10
# интервалов. На этом же графике изобразить
# теоретическую кривую
# Строим теоретическую кривую плотности распределения:
Ypt = normpdf(X);
# Строим эмпирическую кривую плотности распределения:
j = 1;
Xe = [0 0 0 0 0 0 0 0 0 0];
for i = 1:18:n
  Xe(j) = X(i);
  j++;
endfor
Ype = normpdf(Xe);
# Прорисовываем график:
#plot(X, Yt, Xe, Ye)
#grid("on")
# -------------------------------------------------
# -------------- [Пятое задание] --------------
# Для X построить эмпирическую интегральную функцию
# распределения и теоретическую.
function tmp = xb(u, n)
  if(u<0)
    tmp = 0;
  else
    tmp = 1 / n;
  endif
endfunction
function temp = fb(u, X, n)
  tmp = 0;
  for i = 1:n
    tmp += xb(u - X(i), n);
  endfor
  temp = tmp;
endfunction
# Строим эмпирическую кривую
for i = 1:n
  Yie(i) = fb(X(i), X, n);
endfor
# Строим теоретическую кривую
Yit = normcdf(X);
# Прописовываем график
#plot(X, Ye, X, Yt)
#grid("on")
# ---------------------------------------------
# -------------- [Шестое задание] --------------
# Проверить гипотезы: о величине среднего, дисперсии,
# о нормальном законе распределения по Пирсону
# и по Колмогорову
# Проверка гипотезы о величине среднего:
Sx = sqrt(Sx2/n);
t = (Xsr - a)/Sx;
tKrit = tcdf(1-(a/2), n-1);
disp("Выборка X (о величине среднего):");
if(tKrit>t)
  disp("Гипотеза отклоняется");
else
  disp("Гипотеза не отклоняется");
endif
Sy = sqrt(Sy2/n);
t = (Ysr - a)/Sy;
disp("Выборка Y (о величине среднего):");
if(tKrit>t)
  disp("Гипотеза отклоняется");
else
  disp("Гипотеза не отклоняется");
endif
# Проверка гипотезы о величине дисперсии:
HIx2 = n*Sx2/s;
HIy2 = n*Sy2/s;
HIv2 = chi2cdf(a, n-1);
disp("Выборка X (о величине дисперсии):");
if(HIx2>HIv2)
  disp("Гипотеза отклоняется");
else
  disp("Гипотеза не отклоняется");
endif
disp("Выборка Y (о величине дисперсии):");
if(HIy2>HIv2)
  disp("Гипотеза отклоняется");
else
  disp("Гипотеза не отклоняется");
endif
# Проверка гипотезы о нормальном законе
# распределения по Пирсону (хи-квадрат)
function tmp = Pi(Xsr, Sx2)
  tmp = ((((i+1)*18) - Xsr)/(sqrt(Sx2))) - (((i*18) - Xsr)/(sqrt(Sx2)));
endfunction
HI2Xe = 0;
HI2Ye = 0;
for i = 1:9
  HI2Xe += (18 - 18 * Pi(Xsr, Sx2))^2 / (18 * Pi(Xsr, Sx2));
  HI2Ye += (18 - 18 * Pi(Ysr, Sy2))^2 / (18 * Pi(Ysr, Sy2));
endfor
HI2T = chi2cdf(a, n-3);
disp("Выборка X (о нормальном законе по Пирсону):");
if(HI2Xe<HI2T)
  disp("Гипотеза не отклоняется");
else
  disp("Гипотеза отклоняется");
endif;
disp("Выборка Y (о нормальном законе по Пирсону):");
if(HI2Ye<HI2T)
  disp("Гипотеза не отклоняется");
else
  disp("Гипотеза отклоняется");
endif
# Проверка гипотезы о нормальном законе
# распределения по Колмогорову:
Dnx1 = [0 0 0 0 0 0 0 0 0];
Dny1 = [0 0 0 0 0 0 0 0 0];
for i = 1:9
  Dnx1(i) += normcdf(Xe(i)) - fb((Sx2 * 0.8 * (j - 4) + Xsr), X, n) - fb((Sx2 * 0.8 * (j - 5) + Xsr), X, n);
  Dny1(i) += normcdf(Xe(i)) - fb((Sy2 * 0.8 * (j - 4) + Ysr), Y, n) - fb((Sy2 * 0.8 * (j - 5) + Ysr), Y, n);
endfor
Dnx = max(Dnx1);
Dny = max(Dny1);
l = 1.36;
Dnx *= sqrt(n);
Dny *= sqrt(n);
disp("Выборка X (по Колмогорову):");
if(Dnx<l)
  disp("Гипотеза не отклоняется");
else
  disp("Гипотеза отклоняется");
endif
disp("Выборка Y (по Колмогорову):");
if(Dny<l)
  disp("Гипотеза не отклоняется");
else
  disp("Гипотеза отклоняется");
endif
# ----------------------------------------------
# -------------- [График к 4-ому и 5-ому заданию] --------------
plot(X, Ypt, Xe, Ype, X, Yie, X, Yit)
legend('boxon')
legend("4 - Teoreticheskaya", "4 - Empiricheskaya", "5 - Empiricheskaya", "5 - Teoreticheskaya", 2)
grid("on")
# --------------------------------------------------------------