Примеры выполнения работы

Проверка гипотезы о нормальном распределении

Загрузим пакет stats и подпакеты transform, describe.

> restart: with(stats): with(transform): with(describe):

Вводим реализацию выборки (см. данные своего варианта):

> Y: =[15.41, 13.32, 14.28, 12.26, 12.70, 13.97, 10.89, 13.46, 12.79,

13.96, 15.83, 13.27, 14.19, 14.78, 13.35, 16.56, 14.22, 13.26, 13.46,

14.98, 14.30, 14.23, 14.99, 11.90, 15.34, 13.80, 12.13, 13.06, 13.37,

13.69, 12.15, 14.50, 13.34, 13.37, 14.06, 15.82, 11.85, 12.30, 11.86,

12.86, 13.87, 16.39, 12.49, 13.93, 15.33, 14.44, 13.96, 14.74, 16.09,

12.65, 13.40, 13.44, 14.54, 13.23, 12.86, 15.91, 14.54, 12.16, 14.42,

14.76, 13.60, 12.86, 13.60, 13.58, 13.91, 13.49, 13.82, 15.51, 13.92,

15.59, 12.44, 15.70, 14.71, 15.61, 12.88, 11.79, 13.23, 11.79, 16.06,

12.29];

Определим объём выборки (подсчитаем количество значений в выборке) и рассчитаем количество интервалов разбиения k:

> n: =count(Y); k: =round(1+1.4*ln(n));

Проведём сортировку выборки (варианты расположим в порядке возрастания):

> Y1: =statsort(Y);

Находим минимальное и максимальное значения выборки и длину интервала разбиения:

> ymin: =Y1[1]; ymax: =Y1[n]; h: =(ymax-ymin)/k;

Вычислим границы интервалов разбиения:

> Y2: =[seq(ymin+(i-1)*(h+0.0001)..ymin+i*(h+0.0001), i=1..k)];

 

Находим вектор точек разбиения:

> Z: =[seq(ymin+(i-1)*(h+0.0001), i=1..k+1)];

Составляем интервальный ряд частот Y3 (каждому интервалу поставим в соответствие частоту n_i, т.е. число элементов выборки, попадающих в данный интервал) и вектор частот Y3f:

> Y3: =statsort(transform[tallyinto](Y1, Y2));

> Y3f: =transform[frequency](Y3);

Получим интервальный ряд относительных частот (каждому интервалу поставим в соответствие относительную частоту, т.е. частоту, делённую на объём выборки):

> Y4: =transform[scaleweight[1/n]](Y3);

 

Строим гистограмму относительных частот:

> Hist: =statplots[histogram](Y4, color=green):

 

> plots[display](Hist);

 

По виду гистограммы выдвигаем гипотезу о нормальном распределении генеральной совокупности.

Находим накопленные частоты Y5 (накопленная частота показывает, сколько наблюдалось значений, меньших заданного x) и относительные накопленные частоты Y6:

> Y5: =transform[cumulativefrequency](Y3);

> Y6: =transform[cumulativefrequency](Y4);

.

Строим график эмпирической функции распределения:

> p: =[seq(plot(Y6[i], Y2[i], color=blue), i=1..k)]: plots[display](p);

 

Находим точечные оценки математического ожидания a (выборочное среднее значение), дисперсии S и среднего квадратического отклонения s:

> a: =mean(Y);

> S: =variance(Y);

> s: =standarddeviation(Y1);

.

Находим исправленные оценки дисперсии (несмещённая оценка дисперсии) и среднего квадратического отклонения:

> S1: =S*n/(n-1);

> s1: =sqrt(S1);

.

Вычислим вероятности попадания значения случайной величины в первый и последний (k- ый) интервалы:

> p[1]: =evalf(1/(sqrt(2*Pi)*s1)*int(exp(-(t-a)^2/(2*S1)), t=-infinity..Z[2]));

.

> p[k]: =evalf(1/(sqrt(2*Pi)*s1)*int(exp(-(t-a)^2/(2*S1)), t=Z[k]..infinity));

.

Вычислим вероятности попадания значения случайной величины во 2, 3, …, k -1 интервалы по формулам , где :

> for i from 2 to k-1 do p[i]: =evalf(1/(sqrt(2*Pi)*s1)*int(exp(-(t-a)^2/(2*S1)), t=Z[i]..Z[i+1])) od;

Находим теоретические частоты np_i:

> for i from 1 to k do n*p[i] od;

Так как на первом и последнем интервалах np_i < 5, то объединим 1-й со 2-м и 6-й с 7-м интервалы и пересчитаем соответствующие вероятности и частоты:

> p[2]: =p[1]+p[2]; Y3f[2]: =Y3f[1]+Y3f[2]; p[6]: =p[6]+p[7]; Y3f[6]: = Y3f[6] +Y3f[7];

.

Сравним эмпирические n_i и теоретические np_i частоты, для этого находим наблюдаемое значение по формуле , где i = 2, 3, …, 6, так как два первых и два последних интервала объединили.

> chi2: =sum((Y3f[j]-n*p[j])^2/(n*p[j]), j=2..6);

.

По таблице критических точек распределения , по заданномууровню значимости a и числу степеней свободы ν = s-l-1 (s – число интервалов после пересчёта, l – число параметров в гипотетической функции распределения) находят критическую точку . В нашем случае a = 0, 01(см. задание), s = 5, l = 2, т.е. ν = 5-2-1=2, тогда .

Так как , то гипотеза о нормальном распределении генеральной совокупности принимается.

Запишем гипотетическую функцию плотности распределения и построим на одном рисунке гистограмму относительных частот и график плотности гипотетического распределения.

> f: =evalf(1/(sqrt(2*Pi)*s1)*exp(-(x-a)^2/(2*S1)));

> f1: =plot(f, x=ymin-2..ymax+2):

> plots[display](Hist, f1);

 

Запишем гипотетическую функцию распределения и построим её график.

> F: =evalf(1/(sqrt(2*Pi)*s1))*Int(exp(-(t-a)^2/(2*S1)), t=-infinity..x);

> F1: =plot(F, x=ymin-2..ymax+2):

> plots[display](F1);