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

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

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

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

> Y: =[10.63, 26.04, 6.09, 23.42, 5.25, 24.87, 3.24, 6.24, 4.96, 13.74,

13.25, 21.71, 20.96, 34.72, 8.71, 9.06, 19.12, 20.02, 8.58, 34.52,

14.29, 32.13, 13.40, 26.62, 20.13, 6.48, 30.30, 9.16, 12.39, 21.48,

5.28, 13.82, 21.77, 32.26, 21.70, 7.87, 29.74, 21.11, 17.79, 17.67,

27.76, 27.34, 5.87, 5.02, 12.32, 25.43, 31.07, 24.85, 15.14, 25.85,

7.14, 12.78, 24.99, 27.51, 22.59, 29.00, 34.62, 17.65, 9.02, 21.51,

11.24, 22.13, 10.48, 13.20, 12.34, 25.25, 31.73, 28.72, 14.11, 9.62,

17.54, 12.87, 27.15, 18.08, 19.94, 29.86, 30.53, 10.30, 33.13, 23.41];

Определим объём выборки (подсчитаем количество значений в выборке) и рассчитаем количество интервалов разбиения 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);

 

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

> m: =mean(Y);

> S: =variance(Y);

> s: =standarddeviation(Y1);

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

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

> s1: =sqrt(S1);

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

> a: = m-sqrt(3.0)*s1; b: = m+sqrt(3.0)*s1;

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

> p[1]: =(Z[2]-a)/(b-a);

> p[k]: =(b-Z[k])/(b-a);

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

> for j from 2 to 6 do p[j]: =(Z[j+1]-Z[j])/(b-a) od;

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

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

Так как все np_i > 5, то пересчёт не делаем, число интервалов остаётся прежним: k = 7.

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

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

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

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

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

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

> f: =piecewise(x< a, 0, x> =a and x< =b, 1/(b-a), x> b, 0);

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

> plots[display](Hist, f1);

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

и построим её график.

> F: =piecewise(x< a, 0, x> =a and x< =b, (x-a)/(b-a), x> b, 1);

> F1: =plot(F, x=0..ymax+10):

> plots[display](F1);