Архив рубрики: Наука

Презентация с доклада по защите кандидатской диссертации

28 октября защитил диссертацию!)

Как и обещал, выкладываю презентацию с защиты кандидатской. Работа, кстати, называется «Методы и алгоритмы выбора композиции веб-сервисов в системах с сервисно-ориентированной архитектурой» по специальности 05.13.15 «Вычислительные машины, комплексы и компьютерные сети».

Кому интересны детали, смотрите дополнительные слайды в конце. Основные слайды я старался сократить, как мог, потому что суть доклада – это общая характеристика работ, на детали в 15 минутном выступлении времени нет. Доп. слайды оказались полезны во время ответов на вопросы.

PS: Возможно, позже напишу резюме по всему процессу подготовки диссертации, об аспирантуре и пр. Пока надо морально от этого отдалиться.)

Диссертация и автореферат

28 октября защищаю тут кандидатскую, так что все интересующиеся могут

скачать автореферат моей диссертации, скачать саму диссертацию.

Думаю, такие вещи из диссертации, как:

  • пример решения задачи машинного обучения (определение класса чувствительности),
  • решение задачи многокритериального выбора,
  • архитектура СОА систем,
  • способ определения времени обработки запроса с наименьшим разбросом значений
  • и некоторые другие

могут быть кому-то интересны =)

Позднее выложу презентацию с доклада.

Mendeley: как присвоить выбранным источникам русский язык

Если вы работаете в Latex и управляете библиографией с помощью Mendeley, вы, безусловно, счастливый человек.

Однако радость бытия может омрачить необходимость установить свойство language = russian, чтобы в списке литературе русские ссылки на источники были оформлены с использованием русских слов, типа «С., Т., под ред.» и т.д.

Оказалось, что в Mendeley это сделать для большого числа источников очень просто:

  1. Выделяем русские работы в All Documents
  2. Идем в Preferences -> Document Details -> Show fields:
  3. Ставим галочку у поля Language
  4. Закрываем настройки
  5. Теперь в правом меню для всех выбранных источников можно установить значение поля Language в russian!

Mindmap: Верификация программного обеспечения

Когда-то надо было изучить эту тему. С тех пор осталась весьма подробная «карта знаний», или mindmap.

Карта знаний по верификации в формате .xmind.

Кстати, Xmind – очень удобное open-source приложение.

Карта памяти по теме верификации ПО

 

В качестве литературы в основном использовал:

Кулямин В.В. Методы верификации программного обеспечения / Всероссийский конкурсный отбор обзорно-аналитических статей по приоритетному направлению «Информационно-телекоммуникационные системы», 2008. — 117 с.

Font ftmr6a at not found или самая быстрая установка pscyr в mactex 2012 (2013)

Если вы не можете установить шрифт PsCyr в MacTex 2012, то скорей всего вы просто запускаете все команды, начиная с sudo!

Ошибки, которые вы можете встретить:

mktexpk --mfmode / --bdpi 600 --mag 1+240/600 --dpi 840 ftmr6a
что-то вроде ... bitmap error ... 

Font ftmr6a at 840 not found

Чтобы установить этот шрифт достаточно:

  1. Скачать его, например, отсюда в папке helpers.
  2. (Опционально) Настроить права доступа: sudo chmod -R 777 /usr/local/texlive/
  3. Скопировать все файлы в /usr/local/texlive/texmf-local/
  4. Создать папку /usr/local/texlive/texmf-local/fonts/map
  5. Туда скопировать файл pscyr.map
  6. Ввести последовательно команды:
    mktexlsr
    updmap --enable Map pscyr.map
    
  7. profit!!!

Самая распространенная ошибка в мануалах – это совет все команды вводить от рута, т.е. начиная с sudo, что в корне неверно и ведет к бессмысленно потраченному времени в поисках ошибки.

Обсуждение на форуме, которое помогло решить проблему – http://tug.org/pipermail/macostex-archives/2012-August/049748.html

Наука: наш доклад на АКТО-2012

Вчера были в Казани – пожалуй, единственном городе в России (из тех, что я видел), где архитекторам не ограничивают фантазию!)

Впрочем, мы были там по делу – выступали на конференции «Авиакосмические технологии, современные материалы и оборудование. Казань-2012» в секции 7 – «Современные требования к подготовке кадров для работы в авиакосмической отрасли. Новые образовательные технологии в подготовке и переподготовке кадров.»

Доклад наш назывался «Применение инструментов символьных вычислений для проверки решений задач из курсов для конструкторов и технологов в интеллектуальной обучающей веб-системе «Волга»«. Презентацию с доклада вы можете посмотреть ниже:

Основная заслуга по праву принадлежит Наталии Смирновой – главному теоретическому и практическому исследователю и разработчику системы.) Я в основном консультирую и помогаю по технической стороне вопроса, именно поэтому проект разрабатывается на Django 😉

Разберемся с измерением случайных величин (ч. 4/4)

Доверительный интервал

Доверительный интервал – интервал, который с заданной вероятностью накроет неизвестное значение оцениваемого параметра распределения (вообще говоря, это относится к любому измеримому значению).

Оценка какого-либо параметра распределения с помощью доверительного интервала называют интервальной оценкой.

Главное преимущество интервальной оценки – возможность характеризовать уверенность в вычисленных параметрах распределения. Например, можно будет сказать: «С вероятностью 95% среднее значение \mu нормального распределения для моих данных будет лежать в интервале от 100 до 102, или, как еще можно это записать, \mu=100\pm1.

Зачем вообще нужна интервальная оценка? Предположим, что мы хотим узнать время обработки запроса каким-либо веб-сервисом. Для этого мы производим несколько измерений: мы последовательно посылаем запросы, получаем ответы на запросы, записываем время между отправкой запроса и получением ответа (время отлика). Как вы уже могли заметить из предыдущих записей, в ходе измерения мы получаем множество различных времен отклика. Какие-то значения встречаются чаще, какие-то – реже. Всё это говорит о том, что мы имеем дело со случайной величиной. У случайной величины есть некий закон распределения, который задается рядом параметров. Параметры распределения мы не знаем, но можем их оценить на основе имеющейся выборки (результатов измерений). Например, оценкой дисперсии нормального распределения будет выборочная дисперсия.

Однако, наша оценка – это тоже случайная величина. В этом легко убедиться, если запустить программу сбора данных из первой записи еще раз. Вы получите примерно те же данные, но тем не менее они будут другими, а значит и выборочная дисперсия будет другой. Другими словами эта оценка (или, как еще говорят, статистика) зависит от конкретных выборок.

Согласно центральной предельной теореме при большом количестве измерений наша оценка параметра распределения имеет распределение, близкое к нормальному. Поэтому, чем больше у нас будет выборка, тем ближе выборочная дисперсия будет к дисперсии (параметру) нормального закона распределения! Но когда остановится? Понятно, что если мы сделаем 10 000 измерений, то мы достаточно (с уверенностью в 99,9999%) точно сможем оценить истинные значения параметров распределения, но часто такой возможности нет. На помощь приходит интервальная оценка.

Мы можем сказать так: «Мы не хотим делать 10 000 измерений, но нам достаточно быть уверенным в наших оценках параметров распределения на 95%». Проведя вычисления доверительного интервала, вполне может оказаться, что для 95% уверенности достаточно всего лишь 50 измерений. Существенная экономия сил и ресурсов!)

Обозначим за (1-\alpha)*100\% доверительную вероятность, или надежность, или доверительный уровень (англ. level of confidence) – вероятность того, что истинное значение измеряемой окажется в доверительном интервале. Саму величину \alpha называют уровнем значимости. Математически 1-\alpha можно записать так:
P(a\le \theta \le b)=1-\alpha,
где \theta – это истинное значение оцениваемого интервала, a – нижнее значение интервала, b – верхнее.

Интересно, что \alpha уже встречался ранее, также назывался уровнем значимости, но обозначал вероятность отвержения верной нулевой гипотезы.

Итак, чем больше \alpha, тем меньше доверительная вероятность и тем уже доверительный интервал. Обычно достаточно \alpha=0.05\%, тогда доверительная вероятность будет равна 95%.

У разных распределений разные параметры и, соответственно, построение интервальных оценок для этих параметров разное.

Доверительные интервалы нормального распределения

Посмотрим, как выглядит доверительный интервал для математического ожидания нормального распределения при известном стандартном отклонение \sigma (напоминаю, стандартное отклонение – это квадратный корень от дисперсии). Возникает вопрос, а откуда нам может быть известно стандартное отклонение? Эта характеристика распределения случайной величины может быть задана, например, при производстве какого-либо изделия технический процесс должен обеспечит заданное стандартное отклонение. В этой ситуации доверительный интервал, построенные на базе случайной выборки, поможет определить – соответствует ли изготовленная продукция стандартам тех. процесса.

\bar{X}\pm Z_{1-\alpha}\frac{\sigma}{\sqrt{n}}, где n – объем выборки, \overline{X} – выборочное среднее и, самое интересное, Z_{1-\alpha} – квантиль стандартного нормального распределения (нормальное распределение с мат. ожиданием = 0 и стандартном отклонением = 1) для заданной вероятности 1-\alpha.

Не будем приводить все выкладки, почему именно так определяется доверительный интервал. Подробнее вы можете почитать в [1] (стр. 50-55).

Поясню еще раз про квантиль. Квантиль – это такое значение x, которое случайная величина не превысит с заданной вероятностью 1-\alpha. Квантиль – это значение обратной функции распределения, т.е. 1-\alpha=F(Z).

Статистика малых выборок

В случаях, когда стандартное отклонение неизвестно, а также объем выборки не превышает 10-30 значений, используется t-статистика для оценки значения параметра распределения. t-статистика вычисляется как квантиль t-распределения (распределение Стьюдента, см. первую часть этой серии записей) и задается только одним параметром df – степенью свободы. Для нормального распределения df=n-1, где n – объем выборки.

Тогда интервальная оценка мат. ожидания будет выглядеть следующим образом:

\bar{X}\pm t_{\alpha,df} \frac{S_x^2}{\sqrt{n}}, где t_{\alpha,df} – коэффициента Стьюдента (значение квантили статистики t порядка P=1-\alpha/2 для числа степеней свободы df = n – 1), S_x^2 – выборочная дисперсия.

Доверительные интервалы экспоненциального распределения

После предыдущей записи, я попробовал еще одно распределение для имеющихся данных – экспоненциальное. Оказалось, что оно лучше подходит, чем распределение Эрланга. Результаты измерений вы найдете на иллюстрации ниже и в прилагаемом файле ipython notebook. Видимо, это произошло из-за того, что распределение этих конкретных данных не имеет левого хвоста.

Поэтому, продолжая наш пример, найдем доверительный интервал для параметра \lambda экспоненциального распределения. Кстати, если бы мы использовали экспоненциальное распределение в расчетах в рамках теории массового обслуживания, то параметр \lambda называли «интенсивностью поступления заявок» и он обозначал бы среднее время поступления новых заявок за единицу времени.

Оценкой \lambda является обратное значение среднего арифметического по выборке, т.е.
\hat\lambda=1/\bar{x}.

Доверительный интервал выглядит следующим образом:
$latex \frac{2n}{\widehat{\lambda}\chi^2_{1-\alpha/2,2n}} < \frac{1}{\lambda} < \frac{2n}{\widehat{\lambda}\chi^2_{\alpha/2,2n}}$, где $latex \chi^2_{\alpha/2,2n}$ – квантиль распределения хи-квадрат при заданных доверительной вероятностью и степенью свободы $latex k=2n$.

Код ниже проводит сначала графическую оценку параметра распределения, потом точечную с помощью критерия Пирсона и далее интервальную.

#=========================================
# 1) Графический метод
#=========================================

hist(latencies, bins_number, normed=1)  # Построим нормализованную гистограмму

print(ss.expon.fit(latencies))
exp_loc, exp_scale = ss.expon.fit(latencies)        # Оценим положение распределения

# Наберем больше промежуточных точек для плавности графика
x = np.linspace(np.min(latencies), np.max(latencies) + 10, number_of_samples * 2)

# Получим соответствующие значение по функции плотности вероятности для распределения Эрланга
# В качестве начала построения (параметр loc) выберем минимальное значение  из выборки,
# т.е. сместим распределение так, чтобы оно начиналось с нуля.
# В scipy параметр scale = 1 / lambda, где lambda - параметр экспоненциального распределения.
distr_data = ss.expon.pdf(x, loc=np.min(latencies), scale=4.8)

# Построим график
plot(x, distr_data, 'r--')


#=========================================
# 2) Точечная оценка параметра распределения
#=========================================


# Оценим с помощью хи-критерия соответствие распределения
counts, bins = np.histogram(latencies, bins=bins_number)
p_observed = np.array(counts, dtype=float)

tmp = ss.expon.cdf(bins, loc=exp_loc-0.6, scale=4.8)
#print(tmp)

# Now get probabilities for intervals
p_expected = []
for idx, p in enumerate(tmp):
    if idx == len(tmp) - 1: break
    p_expected.append(tmp[idx+1] - tmp[idx])

p_expected = number_of_samples*np.array(p_expected)

# Calculate using scipy built in chi-square test
chi_test_value, probability_fit = ss.chisquare(p_observed, np.array(p_expected, dtype=float))
print('Вероятность того, что экспоненциальное распределение подходит к полученным данным: {0:f}'.format(probability_fit))

# Calculate by ourselved
chi_star = np.sum((p_observed-p_expected)**2/p_expected)
print("chi_star = {}".format(chi_star))
conf_interval = 0.95
df = len(bins) - 3
chi = ss.chi2.ppf(conf_interval, df)  # chi-square quntile for alpha = conf-interval, degrees of freedom = df
print("chi = {}".format(chi))

Получим следующее:

Следующий код проводит интервальную оценку параметра:

#==================================================
# 3) Интервальная оценка параметра распределения
#==================================================


n = size(latencies)
alpha = 0.05

# Сдвинем данные так, чтобы они начинались с 0
data = np.array(latencies) - np.min(latencies)

# Посчитаем оценку лямбда
lambda_estimate = 1/np.mean(data)
print("Оценка 1/lambda на основе выборки: {}".format(1/lambda_estimate))

# Посчитаем левую границу интервала
a = 2 * n / (lambda_estimate * ss.chi2.ppf(1 - alpha / 2, 2 * n))

# Правую
b = 2 * n / (lambda_estimate * ss.chi2.ppf(alpha / 2, 2 * n))
             
print("95% доверительный интервал: {} < 1/lambda < {}".format(a, b))


ci_width = b - a
print("Ширина доверительного интервала: {}".format(ci_width))

Получим следующее:

    Оценка lambda на основе выборки: 0.249376558603
    Оценка 1/lambda на основе выборки: 4.01
    95% доверительный интервал: 3.32700158323 < 1/lambda < 4.92847012339

Я везде привожу \frac{1}{\lambda}, т.к. это значение использовалось в графическом методе вместо просто \lambda.

Определения необходимого объема выборки для нормального распределения при известном стандартном отклонении

А теперь разберемся с самым интересным! Вот ради чего всё и затевалось.)

До этого мы заранее фиксировали объем выборки, не учитывая доверительную вероятность. В действительности же сначала задается, насколько точной должна быть оценка, и только потом определяют объем выборки.

Опять начнем с простой ситуации, когда известно стандартное отклонение. Такой вариант маловероятен, но тем не менее возможен. Вспомним запись мат. ожидания с доверительным интервалом:
\bar{X}\pm Z_{1-\alpha} \frac{\sigma}{\sqrt{n}}
Величина, добавляемая и вычитаемая из \bar{X}, называется ошибкой выборочного исследования и обозначается следующим образом:
e=Z_{1-\alpha} \frac{\sigma}{\sqrt{n}}

Из формулы e можно определить объем выборки n:
n=\frac{Z^2 \sigma^2}{e^2}

Таким образом, для определения объема выборки необходимо знать три параметра.
1. Требуемый доверительный уровень, который влияет на величину Z, являющуюся критическим значением стандартизованного нормального распределения.
2. Приемлемую ошибку выборочного исследования e. Задается в измеряемых единицах. Например, если бы у нас было нормальное распределение e можно было бы взять равным 5 (мс).
3. Стандартное отклонение ?.

Определения необходимого объема выборки для нормального распределения при неизвестном стандартном отклонении

Здесь, как и при определении доверительного интервала, воспользуемся t-распределением и выборочной дисперсией S_x^2.

n\ge t^2_{\alpha,df}  (\frac{S_x}{\delta})^2=t^2_{\alpha,df} \epsilon^2,
где \delta – требуемая точность. Обычно требуемую точность берут в размерах выборочной дисперсии. Например, если нам нужна точность измерений в 2 раза выше выборочной дисперсии, то \epsilon=\frac{S_x}{0.5S_x}=2.

Стоит отметить, что по сути выборочная дисперсия зависит от объема выборки, поэтому вычисляемый объем выборки необходимо решать методом последовательных приближений. В качестве первого значения выборки можно взять, например, df=10. Посчитать n. Получим, например, 8. Снова n, но при df=8-1=7. Получим n=7. И т.д., пока значение n станет приблизительно неизменным.

Определения необходимого объема выборки для экспоненциального распределения

Вернемся к нашему примеру. Мы установили, что время отлика сайта имеет экспоненциальное распределение. Посчитали доверительный интервал для нашей выборки в 100 значений.

Поставим следующую задачу: сколько необходимо измерений, чтобы быть уверенным в том, что оценка параметра \lambda попадет в доверительный интервал с вероятностью 95%?

Если честно ни русской, ни в иностранной литературе ответа на этот вопрос я не нашел, поэтому решил вывести формулу сам по типу того, как это делали с нормальным распределением.

Из неравенства доверительного интервала выведем формулу расчета необходимого объема выборки. Значение ширины доверительного интервала рассчитываем как разницу между верхней и нижней интервальными оценками и примем за \epsilon (эпсилон). Это значение мы зафиксируем.

\frac{1}{n} = \frac{2}{\epsilon \widehat{\lambda}} (\frac{1}{\chi^2_{\alpha/2,2n}} - \frac{1}{\chi^2_{1-\alpha/2,2n}})

Код ниже просто рассчитывает приведенную формулу:

#==================================================
# 4) Определение необходимого объема выборки 
# при фиксированной ширине доверительного интервала
#==================================================

n = 100
epsilon = ci_width
alpha = 0.1

n_ci = (2 / (epsilon * lambda_estimate)) * (1/ss.chi2.ppf(alpha / 2, 2 * n) - 1/ss.chi2.ppf(1 - alpha / 2, 2 * n))
n_ci = 1 / n_ci
print('Необходимый объем выборки = {} для попадания истинного значения параметра распределения' \
     ' в интервал шириной = {} с доверительной вероятностью = {}%'.format(round(n_ci), epsilon, (1-alpha)*100))

Результ: «Необходимый объем выборки = 120.0 для попадания истинного значения параметра в интервал шириной = 1.60146854016 с доверительной вероятностью = 90.0%».

Нахождение объема выборки при фиксированной доверительной вероятности и ширине доверительного интервала

Дальше идут весьма сомнительные рассуждения, поэтому не рекомендую их читать (тем более всерьез). =) Задача в заголовке раздела вполне может появится, но весьма редко. В ходе решения предыдущей задачи (когда мы могли изменять доверительную вероятность) я неспециально решил задачу в заголовке.

Как в случае с нормальным распределением, формула вышла рекуррентной. Оптимальное значение можно найти с помощью метода последовательных приближений или с помощью линейной регрессии. Для использования метода линейной регрессии необходимо вычислить частную производную по n, что в нашем случае не так-то просто, поэтому ограничимся простым методом приближения.

Код (который вы также найдете в файле ipython notebook в конце поста):

#==================================================
# 5*) Определение необходимого объема выборки 
# при фиксированной ширине доверительного интервала и доверительной вероятности
#==================================================
alpha = 0.05
epsilon =  0.5 * ci_width
optimum_n = 0
N = [0]

for i in range(0, 500, 10):
    n = i
    
    for k in range(5):
        n_ci = (2 / (epsilon * lambda_estimate)) * (1/ss.chi2.ppf(alpha / 2, 2 * n) - 1/ss.chi2.ppf(1 - alpha / 2, 2 * n))
        n_ci = 1 / n_ci
        n = n_ci
        N.append(n_ci)
        #print("{}".format(n_ci))
    
    #print('---')
    
    # Замечено, что при увеличении первоначального значения n
    # в определенной итерации метода приближения вычисленный объем выборки
    # перестает падать или расти, а колеблется около одного и того же числа.
    # Это число можно считать оптимальным объемом выборки при заданной ширине доверительного интервала.
    # Что требует доказательства.
    
    # Ниже мы грубо оцениваем размах колебания значения n_ci в этой итерации.
    # Если оно мало, значит это значение n можно считать оптимальным.
    if np.std(N[-5:]) < 5:
        optimum_n = n_ci

print('Эверстически оптимальный объем выборки: {}'.format(math.ceil(optimum_n)))
print('При заданной ширине доверительного интервала: {}'.format(epsilon))
print('И доверительной вероятности: {}%'.format((1-alpha)*100))
plot(N)

Если использовать в качестве начального значения n изначальный объем выборки (в нашем случае это 100), то каждое следующее приближение будет уменьшать вычисленной значение n_ci до нуля, однако, я заметил, что при увеличении первоначального значения n (больше первоначального объем) в определенной итерации метода приближения вычисленный объем выборки перестает падать или расти, а колеблется около одного и того же числа. Это число можно считать оптимальным объемом выборки при заданной ширине доверительного интервала. Это, конечно, требует доказательства и является полнейшей эвристикой.

Полученный с помощью этой эвристики результаты хорошо согласуются с предполагаемым поведением этой функции: чем шире значение доверительного интервала, тем меньше нужен объем выборки и наоборот. Однако сам алгоритм нестабилен и не выдает результатов при \epsilon на порядок меньше исходной ширины интервала. При больших значениях он соответственнно выдаст 0 (что можно понимать как выборку с 1 значением), т.к. если интервал достаточно широк, то попасть в него будет очень просто с любой выборкой.

Ниже график, на котором показано вычисление объема выборки методом последовательных приближений при различных начальных значениях n. Хорошо виден интервал, в котором вычисленное значение колеблется, но заметно не растет и не падает.

В общем, если что-то подобное пригодится в будущем, я спрошу совета у продвинутых математиков.) Если тот, кто читает эту запись, предложит другой способ определения объема выборки для экспоненциального распределения, я с удовольствием его послушаю!

Скачать файл с ipython notebook раздела

Заключение

Занимаясь написанием этой череды постов, я узнал много нового про распределения и методы работы с ними, надеюсь и вы почерпнули для себя что-нибудь интересное и полезное. =)

Список литературы

  1. Спирин H.A. Методы планирования и обработки результатов инженерного эксперимента : Конспект лекций. / H.A. Спирин, В.В. Лавров. Екатеринбург : ГОУ ВПО УГТУ-УПИ, 2004. — 257 с.
  2. Статистика для менеджеров с использованием Microsoft Office Excel 4-e издание. Дэвид М. Левин и др. – отлично рассказано в 7 главе про практическое использование доверительных интервалов. И именно 7-я глава есть в открытом доступе.
  3. Confidence Intervals, Guy Lebanon – отличная работа, рассказывающая о доверительных интервалах.
  4. http://en.wikipedia.org/wiki/Exponential_distribution#Maximum_likelihood – про оценку параметра экспоненциального распределения.

Разберемся с измерением случайных величин (ч. 3/4)

Подбор подходящего распределения на основе экспериментальных данных

На английском такую операцию называют «distribution fitting«.

Перед тем как углубляться в рассмотрение этих способов введем понятие гипотезы, который вы практически наверняка встретите в обсуждении темы подбора распределения.

Статистические гипотезы

Статистической гипотезой называется любое предположение о виде или параметрах неизвестного закона распределения.

Проверяемую гипотезу обычно называют нулевой (или основной) и обозначают H_0. Наряду с нулевой гипотезой H_0 рассматривают альтернативную, или конкурирующую, гипотезу H_1, являющуюся логическим отрицанием H_0. Нулевая и альтернативная гипотезы представляют собой две возможности выбора, осуществляемого в задачах проверки статистических гипотез. Правило, по которому гипотеза H_0 отвергается или принимается, называется статистическим критерием или статистическим тестом.

Алгоритм выбора подходящего распределений

Существует множество способов (статистических критериев) определения подходящего распределения. Объединим их в общий алгоритм выбора распределений, начиная с самого простого и очевидного и заканчивая наиболее точным и надежным. В качестве экспериментальных данных будем использовать пример из первой части с распределением времен открытия веб-страниц.

Организационный момент: все вычисления будут проводится в веб-версии ipython – ipython notebook, которую вы можете запустить из консоли соответствующей командой: «ipython notebook». Это весьма удобная система для подобных задач, и, кстати, по общей концепции работы она похожа на рабочую тетрадь из дорогущей мат. системы Wolfram Mathematica, слава open-source!)

В конце я выложу файл с рабочей тетрадью (*.ipnb), в котором вы найдете все исходные коды из данной записи и сможете одним кликом запустить их у себя в своем ipython notebook.

1.1. Графический метод (ищем возможные подходящие распределения)

Суть графического метода в построении гистограммы полученных экспериментальных данных и последующим сравнении внешнего вида полученного графика с графиками плотности вероятностей (p.d.f.) различных законов распределения. И вправду, по сути гистограмма очень похожа на функцию плотности вероятности, только гистограмма весьма груба, дискретна и, значит, неточна, но для начала нам только это нужно для формирования гипотезы о действительном распределении нашей случайной величины.

Итак, запустим ipython с поддержкой построения графиков в самой рабочей тетради (что мега-круто!)) с помощью следующей команды в консоли:

ipython notebook --pylab inline

1. Получим необходимые данные с помощью следующего кода:

import httplib
import time

latencies = []
number_of_samples = 100
for i in range(number_of_samples):
    conn = httplib.HTTPConnection("livejournal.com")
    conn.request("GET", "/")
    start_time = time.time()
    r1 = conn.getresponse()
    latency = round((time.time() - start_time) * 1000)
    latencies.append(latency)
    print("{1}\\{2}. Задержка: {0}мс.".format(latency, i+1, number_of_samples))

2. Определим число участков разбиения и построим гистограмму времен открытия веб-страницы

Число участков я определял эмпирически, чтобы с одной стороны участков разбиения было по-больше, и с другой – чтобы каждый участок включал по-больше значение. По-хорошему стоило бы воспользоваться либо формулой Стэджесса

k=\lfloor1+3.322*lg(n)\rfloor

либо просто

k=\sqrt{n},

где  k – число участков разбиения в гистограмме, n – число измерений.

В этот раз рассмотрим популярный сайт livejournal.com, он же ЖЖ.

import math
bins_number = round(10+math.sqrt(len(latencies)))
print(u"Число участков разбиения в гистограмме:{0}".format(bins_number))
<p>hist(latencies, bins_number)
title("Histogram of latencies of web-page load")
xlabel("Latency")
ylabel("Frequency")
show()

Получим следующее

3. Проанализируем гистограмму

Для начала отметим, что распределение несимметрично относительно значения с максимальной частотой появления (назовем такое значение «центром«). Правая часть распределения (или как говорят «хвост«) весьма растянута, т.е. вероятность большего значения случайного значения не так уж и мала. Такой хвост так и называют «длинный хвост» (long-tail). Если вероятность удаленного значения относительно велика, то такую ситуацию называют «тяжелый хвост» (heavy-tail).

Руководствуясь изображениями распределений, например, отсюда, можно предположить, что неплохо подходит частный случай гамма-распределенияраспределение Эрланга. Функция плотности этого распределения выглядит так:

4. Подберем параметры распределения и построим график распределения с наложением на гистограмму

И опять на помощь приходит numpy/scipy. В этой библиотеке есть порядка 80 распределений с функциями вероятности, плотности распределения и подбора параметров распределения на основе имеющихся данных.

Следующий код сделает ровно то, о чем гласит подзаголовок этого пункта.

import scipy.stats as ss
import scipy as sp
import numpy as np

hist(latencies, bins_number, normed=1)  # Построим нормализованную гистограмму

k = 1.5
print(ss.erlang.fit(latencies, k))
shape, loc, lambd = ss.erlang.fit(latencies, k)        # Оценим положение распределения

# Наберем больше промежуточных точек для плавности графика
x = np.linspace(np.min(latencies), np.max(latencies) + 10, number_of_samples * 2)

# Получим соответствующие значение по функции плотности вероятности для распределения Эрланга
distr_data = ss.erlang.pdf(x, k, loc=loc, scale=6)

# Построим график
plot(x, distr_data, 'r--')

Вы получите примерно следующее:

Левая хвост не особо вписался в наши данные, зато правый хвост подошел.

Любопытно, что если поменять параметр \lambda на 0.5, то мы получим распределение хи-квадрат с k-степенями свободы. Впрочем, это нам ничем здесь не поможет.)

В качестве альтернативного распределения выберем распределение Стьюдента, или t-распределение. По форме оно очень похоже на нормальное распределение и при значении единственного параметра, называемого степенью свободы (degree of freedom, df) близкого к 100-120 практически идентичного нормальному распределению. Главное его достоинство – более тяжелые хвосты, что полезно в ситуации малых выборок (меньше 30). На рисунке ниже голубая линия – нормальное распределение, красная – t-распределение.

Следующий код строит t-распределение для наших данных. Значения параметров распределения я подобрал вручную, поскольку функция fit() для t-распределения хорошо определяет только центр распределения, с другими параметрами у неё беда.

hist(latencies, bins_number, normed=1)  # Построим нормализованную гистограмму

df = 100

print(ss.t.fit(latencies, df))
tmp, loc, scale = ss.t.fit(latencies, df)        # Оценим положение распределения

# Наберем больше промежуточных точек для плавности графика
x = np.linspace(np.min(latencies)-10, np.max(latencies) + 10, number_of_samples * 2)

# Получим соответствующие значение по функции плотности вероятности для t-распределения
distr_data = ss.t.pdf(x, df, loc=loc, scale=8)

# Построим график
plot(x, distr_data, 'r--')

На глаз получилось не так уж плохо.)

Еще пачка важный обозначений

Перед рассказом о методе оценки коэффициента вариации введем пару обозначений, которые понадобятся нам в будущем.

Статистика – функция результатов наблюдений, используемая для оценки параметров распределения и (или) для проверки статистических гипотез.

По выборке невозможно найти параметры распределения случайной величины (поскольку для этого требуется бесконечное количество результатов наблюдений – изучение всей генеральной совокупности), поэтому, имея в своем распоряжении всегда ограниченный объем экспериментальных данных, исследователю остается довольствоваться только лишь получением некоторых оценок.

Оценка – статистика, являющаяся основой для оценивания неизвестного параметра распределения.

Оценка может быть точечной и интервальной.

Точечной оценкой неизвестного параметра называют число (точку на числовой оси), которое приблизительно равно оцениваемому параметру и может заменить его с достаточной степенью точности в статистических расчетах.

Интервальной называют оценку, которая определяется двумя числами – концами отрезка и вероятностью попадания истинного значения оцениваемого параметра в этот отрезок..

В таблице ниже приведены оценки различных параметров распределения.

Данная таблица разделена на две группы: обозначения, относящиеся к выборочным данным (результатам измерений), и обозначения, относящиеся к случайной величине.

Выборка (стат. совокупность)Случайная величина (генеральная совокупность)
Обозначение: {x_1,x_2,...,x_n}Обозначение: X
Характеристика: n – объем выборки, число элементов в выборкеХарактеристика: Закон распределения
Среднее арифметическоеArithmetic mean, mean, average

\overline{x}=\frac{1}{n}\sum_{i=1}^{n}x_i

Является оценкой мат. ожидания.

Математическое ожидание

Expected value, expectation, mathematical expectation, mean, first moment

\mu, E[X]

Cтандартное отклонение, стандартный разброс, выборочное стандартное отклонение
Sample standart deviation

S=\sqrt{\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\overline{x})^2}

Является несмещённой оценкой стандартного отклонения случайной величины,т.е. мат. ожидание этой оценки равно оцениваемому параметру.

Стандартное отклонение

Standart deviation

\sigma=\sqrt{E[(X-\mu^2)]}

\sigma=\sqrt{E[X^2]-(E[X])^2}

Среднеквадратическое отклонение, среднеквадратичное отклонение, квадратичное отклонение
Standard deviation of the sample
S_x
S_N=\sqrt{\frac{1}{n}\sum_{i=1}^{n}(x_i-\overline{x})^2}

Также является оценкой стандартного отклонения распределения. В отличие от предыдущей оценки является смещенной оценкой.

Несмещенная (исправленная) выборочная дисперсия

Sample variance

S^2=\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\overline{x})^2

Является несмещенной оценкой дисперсии.

Дисперсия

Variance

\sigma^2, \sigma_X^2, Var(X), D[X]

Выборочная дисперсия

Sample variance

S_n^2=\frac{1}{n}\sum_{i=1}^{n}(x_i-\overline{x})^2

Является смещенной оценкой дисперсии.

Среднее арифметического значение по полученным данным, среднее значение, выборочное среднее арифметическое, выборочное среднее – это простое среднее арифметическое, взятое из полученных данных. С терминологией здесь опять кто в лес, кто по дрова. Встречаются обозначения m*, \overline{M}, \overline{x}.

VarianceДисперсия – мера рассеивания значений относительно среднего значения. Малое значение дисперсии указывает на то, что разброс значений невелик относительно общего среднего значения.  Тут в литературе могут иметь в виду сразу два понятия:

Mean, expected value, expectationМатематическое ожидание – параметр закона распределения, по смыслу означает самое вероятное значение случайной величины. На графиках плотности распределения обычно на пике кривой. Есть у всех распределений (на сколько мне известно. Обозначается как M_x, m_x, M(x), \mu, E[X]. Считается для всех распределений по-разному. Для нормального распределения считается как среднее арифметическое по равномерно взятым значениям нормального распределения.

Фу-ух, с обозначениями пока всё. Теперь перейдем непосредственно к

1.2. Оценка коэффициента вариации

В случае, если вы уверены, что данные имеют, например, нормальное распределение, то в качестве простой альтернативы можно использовать метод оценки коэффициента вариации. Также его можно использовать дополнительно к графическому методу.

Коэффициент вариации V = \frac{\S_n}{\overline{x}} – это мера относительного разброса случайной величины; показывает, какую долю среднего значения этой величины составляет её средний разброс.

Полученное значение сравнивается со следующей таблицей и делается вывод, какое распределение лучше подходит:

ПРЕДЕЛЫ ИЗМЕНЕНИЯ КОЭФФИЦИЕНТА ВАРИАЦИИ VXЗАКОН РАСПРЕДЕЛЕНИЯ СЛУЧАЙНОЙ ВЕЛИЧИНЫ Х
Vx <= 0.3Нормальный
0.3 <= Vx < 0.4Гамма-распределение
0.4 <= Vx < 1Вейбулла
Vx = 1 Экспоненциальный

Попробуем метод на практике!

Тут я должен извиниться, я посеял данные с прошлого раза, поэтому собрал их снова, и теперь они выглядят немного по-другому:

Общая тенденция осталась примерно та же, но левый хвост практически исчез. Теперь первоначальная подгонка распределения Эрланга по данным выглядит следующим образом:

Получилось даже лучше, чем раньше.)

Ну а теперь, наконец, посчитаем коэффициент вариации!

Данные получили следующие:

[189.0, 189.0,

186.0, 189.0, 189.0, 191.0, 190.0, 186.0, 188.0, 186.0, 187.0, 188.0,

186.0, 189.0, 208.0, 186.0, 194.0, 187.0, 186.0, 193.0, 190.0, 194.0,

186.0, 187.0, 201.0, 188.0, 188.0, 187.0, 186.0, 188.0, 189.0, 188.0,

185.0, 188.0, 191.0, 193.0, 194.0, 191.0, 188.0, 192.0, 187.0, 190.0,

189.0, 188.0, 193.0, 190.0, 188.0, 193.0, 187.0, 185.0, 186.0, 189.0,

190.0, 190.0, 189.0, 186.0, 186.0, 186.0, 187.0, 191.0, 186.0, 189.0,

188.0, 188.0, 194.0, 186.0, 195.0, 187.0, 185.0, 190.0, 187.0, 187.0,

188.0, 187.0, 188.0, 190.0, 188.0, 189.0, 189.0, 193.0, 187.0, 193.0,

186.0, 189.0, 188.0, 190.0, 186.0, 190.0, 188.0, 186.0, 185.0, 189.0,

189.0, 189.0, 194.0, 187.0, 192.0, 187.0, 190.0, 191.0]

Посчитаем среднее \overline{x}, стандартное отклонение S и коэф. вариации V_x

n = number_of_samples
mean = np.sum(latencies)/n
print("mean = {0}".format(mean))
S = np.sqrt(np.sum((latencies-np.mean(latencies))**2)/n-1)
print("S = {0}".format(S))
V_x = S/mean
print("V_x = {0}".format(V_x))

Получим \overline{x} = 189.01, S = 3.16384260038 и V_x = 0.016739022276

Если судить по таблице, то наши данные имеют нормальный закон распределения, что не так. В нашем примере анализ коэффициента вариации не помог.

Критерии согласия

Критерий согласия – статистический критерий для проверки гипотезы о согласии (равенстве) распределения случайной величины исследуемой совокупности с теоретическим распределением или гипотезы о согласии распределений в двух и более совокупностях.

Кстати, любой статистический критерий является случайной величиной, поэтому статистические гипотезы носят вероятностный характер. Уровень значимости \alpha – это вероятность отвержения верной нулевой гипотезы.

\beta – вероятность принятия неверной нулевой гипотезы. Мощность критерия 1-\beta – вероятность отвержения нулевой гипотезы, когда она неверна.

Часто используются следующие критерии согласия:

  1. Критерий Пирсона, или критерий хи-квадрат (\chi^2) (англ. chi-square test, который по сути вычисляет меру различия между гистограммой полученных данных и гистограммой теоретического распределения, а далее сравнивается полученное значение с предполагаемым значением (квантили можно взять из таблицы или посчитать самим) при заданных уровне значимости и числе степеней свободы df. Если полученное число меньше предполагаемого, то распределение подходит.

    df вычисляется как число интервалов минус число оцениваемых параметров распределения и минус 1. Например, если сравниваем с нормальным распределением, то число оцениваемых параметров равно двум (дисперсия и среднее), значит df=n-2-1=n-3

    \chi^2* = n*\sum_{i=1}^{m}\frac{(p_i*-p_i)^2}{p_i}

    Есть соответствующие функции в scipy, matlab.

  2. Критерий Колмогорова-Смирнова (англ. KS-test) – похож на критерий Пирсона, только анализируется разность между функцией распределения (cdf) предполагаемого распределения и имеющегося распределения данных измерений.

    Есть соответствующие функции в scipy, matlab.

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

Примем за нулевую гипотезу распределение Эрланга. Проверим эту гипотезу используя критерий согласия Пирсона. В коде ниже происходит следующее:

  1. В переменную p_observed записываем частоты появления интервалов. В переменную bins – границы интервалов.
  2. В tmp записываем значения функции распределения Эрланга для наших интервалов. В качестве параметров распределения берем значения из предыдущего шага с графическим методом. Тут я немного сместил центр распределения (параметр loc) для более точного подбора распределения.
  3. В цикле считаем вероятности значений в рамках наших интервалов. Далее получаем частоты появления интервалов согласно предполагаемому распределению Эрланга
  4. С помощью функции scipy получаем вероятность того, что нулевая гипотеза верна. Получилось около 10%.
  5. Далее делаем всё тоже самое, что делает функция scipy, только сами.

В прилагаемом файле ipython notebook вы также найдете некоторые комментарии.

counts, bins = np.histogram(latencies, bins=bins_number)
p_observed = np.array(counts, dtype=float)
#print(p_observed)

tmp = ss.erlang.cdf(bins, k, loc=erlang_loc-0.55, scale=lambd)
#print(tmp)

# Now get probabilities for intervals
p_expected = []
for idx, p in enumerate(tmp):
    if idx == len(tmp) - 1: break
    p_expected.append(tmp[idx+1] - tmp[idx])

p_expected = n*np.array(p_expected)

# Calculate using scipy built in chi-square test
chi_test_value, probability_fit = ss.chisquare(p_observed, np.array(p_expected, dtype=float))
print('Вероятность того, что распределение Эрланга подходит к полученным данным: {0:f}'.format(probability_fit))

# Calculate by ourselved
chi_star = np.sum((p_observed-p_expected)**2/p_expected)
print("chi_star = {}".format(chi_star))
conf_interval = 0.95
df = len(bins) - 3
chi = ss.chi2.ppf(conf_interval, df)  # chi-square quntile for alpha = conf-interval, degrees of freedom = df
print("chi = {}".format(chi))

    Вероятность того, что распределение Эрланга подходит к полученным данным: 0.098665
    chi_star = 27.2634383493
    chi = 28.8692994304

\chi>\chi* – это указывает на то, что распределение Эрланга подходит для наших данных.

Однако, стоит отметить, что разница между величинами небольшая.

В коде также есть вычисление квантиля распределения (функция ppf, percent point function). Квантиль – это значение обратной функции распределения. Если в обычной функции распределения F(x) в качестве x мы принимаем какое-либо значение и получаем вероятность того, что случайная величина X с этим распределением примет значения от минус бесконечности до x, то в обратной функции мы решаем обратную задачу: при заданном p мы получаем x, то есть получаем такое значение x, которое случайная величина не превысит с заданной вероятностью p.

Если вы предполагаете, что у вас нормальное распределение, вы можете воспользоваться еще одним простым способом проверки – методом моментов. В нём просто сравниваются выборочное среднее значение (среднее значение по полученным данным) с математическим ожиданием предполагаемого распределения и тоже самое делается для дисперсии. Чем меньше разница – тем лучше.

Итог

Сегодня мы так и не добрались до доверительных интервалов, однако изучили много нового: основные определения статистики c распространенными синонимами, два основных подхода выбора подходящего распределения для экспериментальных данных: графический и статистический, а также попробовали всё это на практике, используя превосходных open-source пакет NumPy/SciPy для Python.

Скачать файл ipython notebook

В нём уже содержится весь код из записи и результаты его выполнения. Вы можете сами менять параметры и распределения. Я уверен, что немного поразбиравшись в теме, вы найдете более подходящие распределения и получите впечатляющие результаты критериев согласия!)

Используемая литература

  1. http://www.statsoft.com/textbook/distribution-fitting/ – много занимательных gif-ок, в которых показывают различные распределения и измерения кривых распределения и плотности в зависимости от различных значений параметров распределений. Полезная информация для графического подбора распределения.
  2. http://iglin.exponenta.ru/All/ContData/ContData.html – полезная информация о проверки стат. гипотез.
  3. http://sernam.ru/book_tp.php?id=19
  4. http://vfkomd.ru/docs/lections/ms/18_ocenka%20raspred.pdf – глава из неизвестного учебника про подбор распределения.
  5. http://www.hr-portal.ru/statistica/gl3/gl3.php – просто и понятно про статистику
  6. http://www.mathworks.com/... – немало информации по теме от mathworks.
  7. http://www.slideshare.net/enthought/numpyscipy-statistics – отличная презентация о возможностях nympy/scipy.
  8. http://www.grandars.ru/student/statistika/generalnaya-sovokupnost.html – хорошее практическое объяснение всех терминов.
  9. http://statanaliz.info/teoriya-i-praktika/10-variatsiya/15-dispersiya-standartnoe-otklonenie-koeffitsient-variatsii.html – статистика простыми словами
  10. http://esotericpl.narod.ru/Ps_enc/Ps_enc9_4_0.html – хорошо рассказано про статистику малых выборок
  11. http://www.enviroliteracy.org/pdf/materials/1210.pdf – хорошая работа по критерию Пирсона.