Глава 3. Функция вероятности

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

Для представления ФВ мы будем использовать объекты класса Pmf, содержащего набор значений и их вероятности. Объекты Pmf окажутся полезными для вычисления среднего значения и дисперсии распределения, а также коэффициента асимметрии, показывающего, скошено распределение влево или вправо. И наконец, мы рассмотрим «парадокс проверки» (inspection paradox), проявляющийся в том, что выборка может дать необъективную оценку распределения.

Функция вероятности

Объект класса Pmf похож на FreqTab, только содержит вероятности, а не частоты. Поэтому один из способов создать Pmf — это задействовать объект FreqTab. Например, рассмотрим FreqTab, представляющий распределение значений в короткой последовательности:

from empiricaldist import FreqTab

 

ftab = FreqTab.from_seq([1, 2, 2, 3, 5])

ftab

	частота
1	1
2	2
3	1
5	1

Сумма частот равна размеру исходного ряда:

n = ftab.sum()

n

 

5

Если мы разделим частоту на n, то полученное значение будет представлять собой пропорцию, а не количество:

pmf = ftab / n

pmf

	вероятность
1	0.2
2	0.4
3	0.2
5	0.2

Результат указывает на то, что 20 % значений в последовательности равны 1, 40 % — 2 и т.д.

Эти пропорции также можно считать вероятностями в следующем смысле. Если из исходной последовательности выбрать случайное значение, то вероятность получить значение 1 равна 0.2, значение 2 — 0.4 и т.д.

Поскольку мы делили все частоты на n, то сумма вероятностей равна 1, это означает, что полученное распределение нормированное (normalized):

pmf.sum()

1.0

Нормированный объект FreqTab представляет собой функцию вероятности, названную так потому, что вероятности, связанные с дискретными значениями, также называют «вероятностными мерами».

В библиотеке empiricaldist реализован класс Pmf, представляющий функцию вероятности. Поэтому вместо создания объекта FreqTab и последующей его нормировки можно сразу создать объект Pmf:

from empiricaldist import Pmf

 

pmf = Pmf.from_seq([1, 2, 2, 3, 5])

pmf

	вероятность
1	0.2
2	0.4
3	0.2
5	0.2

Объект Pmf нормированный, поэтому суммарная вероятность равна 1:

pmf.sum()

1.0

Объекты Pmf и FreqTab во многом похожи. Чтобы найти вероятность, связанную с каким-либо значением, можно использовать оператор квадратных скобок:

pmf[2]

0.4

Или можно использовать круглые скобки для вызова Pmf как функции:

pmf(2)

0.4

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

pmf[2] = 0.2

pmf(2)

0.2

Вы можете изменить объект Pmf, увеличив вероятность, связанную с некоторым значением:

pmf[2] += 0.3

pmf[2]

0.5

Или можете умножить вероятность на коэффициент:

pmf[2] *= 0.5

pmf[2]

0.25

После изменения результирующий объект Pmf может оказаться ненормированным, то есть вероятности в сумме могут больше не равняться 1:

pmf.sum()

0.8500000000000001

Метод normalize перенормирует Pmf путем деления на сумму вероятностей и возвращает эту сумму:

pmf.normalize()

0.8500000000000001

Объекты Pmf имеют метод копирования copy, поэтому можно создавать и изменять копию, не меняя оригинальный объект:

pmf.copy()

	вероятность
1	0.235294
2	0.294118
3	0.235294
5	0.235294

Как и FreqTab, объект Pmf имеет атрибут qs для доступа к величинам и атрибут ps для доступа к вероятностям.

В нем также есть метод bar, отображающий Pmf в виде столбчатой диаграммы, и метод plot, изображающий распределение в виде линейного графика.

Сводные характеристики функции вероятности

В разделе «Сводные статистические показатели» на с. 28 мы вычисляли среднее значение выборки, складывая все элементы и деля сумму на количество элементов. Вот простой пример:

seq = [1, 2, 2, 3, 5]

 

n = len(seq)

mean = np.sum(seq) / n

mean

2.6

Теперь предположим, что у нас есть функция вероятности значений в последовательности:

pmf = Pmf.from_seq(seq)

Имея объект Pmf, мы по-прежнему можем вычислить среднее значение, но процедура будет другой: нужно перемножить вероятности и соответствующие величины, а потом сложить полученные результаты:

mean = np.sum(pmf.ps * pmf.qs)

mean

2.6

Обратите внимание, что здесь нам не нужно делить на n, потому что мы уже сделали это, когда нормировали Pmf. У объектов Pmf есть метод mean, который делает то же самое:

pmf.mean()

2.6

Имея под рукой Pmf, мы можем вычислить дисперсию, посчитав отклонение (deviation) каждой величины от среднего значения:

deviations = pmf.qs - mean

Затем умножаем квадраты отклонений на вероятности и складываем полученные результаты:

var = np.sum(pmf.ps * deviations**2)

var

1.84

Метод var дает тот же результат:

pmf.var()

1.84

Зная дисперсию, можно вычислить стандартное отклонение обычным способом:

np.sqrt(var)

1.3564659966250536

Или воспользоваться методом std с тем же результатом:

pmf.std()

1.3564659966250536

Класс Pmf также предоставляет метод mode, находящий значение, связанное с наибольшей вероятностью:

pmf.mode()

2

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

Парадокс размера группы

В качестве примера того, для чего могут пригодиться объекты Pmf, рассмотрим феномен, который я называю «парадоксом размера группы» (class size paradox).

Во многих американских колледжах и университетах соотношение студентов и преподавателей равняется примерно 10:1. Но студентов часто удивляет то, что во многих группах обучается больше 10 человек, иногда заметно больше. Такое несоответствие объясняется следующими двумя причинами:

• Студенты обычно выбирают четыре или пять курсов в семестр, но преподаватели зачастую ведут только один или два курса.

• Количество студентов в группе малочисленного курса невелико, а в группе многочисленного — велико.

Первый факт очевиден — по крайней мере, когда на него обращают внимание; второй требует пояснений. Рассмотрим пример. Предположим, что колледж предлагает 65 курсов в течение одного семестра, и нам задано количество курсов в каждом из следующих диапазонов размера группы:

ranges = pd.interval_range(start=5, end=50, freq=5, closed="left")

ranges.name = "class size"

 

data = pd.DataFrame(index=ranges)

data["count"] = [8, 8, 14, 4, 6, 12, 8, 3, 2]

data

class size (размер группы)	count (количество)
[5, 10)	8
[10, 15)	8
[15, 20)	14
[20, 25)	4
[25, 30)	6
[30, 35)	12
[35, 40)	8
[40, 45)	3
[45, 50)	2

Функция Pandas interval_range создает объект Index, где каждая метка представляет диапазон значений. Запись [5, 10) означает, что 5 включено в интервал, а 10 — нет. Поскольку мы не знаем размеры групп в каждом интервале, то предположим, что каждая группа находится в середине своего диапазона:

sizes = ranges.left + 2

sizes

Index([7, 12, 17, 22, 27, 32, 37, 42, 47], dtype='int64')

Теперь создадим объект Pmf, представляющий распределение размера группы. Поскольку мы знаем размеры и их частоту, то можем создать Pmf напрямую, передав в качестве аргументов количество (counts), размеры (sizes) и название (name). Результатом нормирования нового объекта Pmf будет сумма значений:

counts = data["count"]

actual_pmf = Pmf(counts, sizes, name="фактическое")

actual_pmf.normalize()

65

Если вы спросите сотрудников колледжа о среднем размере группы, то получите в ответ среднее значение этого распределения, которое составляет 23.7:

actual_pmf.mean()

23.692307692307693

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

Следующая функция принимает на входе объект Pmf фактического размера групп и создает новый объект Pmf, соответствующий размеру групп в представлении студентов. Численности групп в обоих распределениях одинаковы, но вероятности во втором распределении умножаются на размер, потому что в группе размера x присутствует x студентов-наблюдателей. Таким образом, вероятность наблюдения группы пропорциональна ее размеру:

def bias(pmf, name):

    # умножить каждую вероятность на размер группы

    ps = pmf.ps * pmf.qs

 

    # создать новый объект Pmf и нормировать его

    new_pmf = Pmf(ps, pmf.qs, name=name)

    new_pmf.normalize()

    return new_pmf

Теперь можно посчитать наблюдаемое студентами смещенное (biased) распределение Pmf:

observed_pmf = bias(actual_pmf, name="наблюдаемое")

Вот как выглядят эти два распределения:

from thinkstats import two_bar_plots

 

two_bar_plots(actual_pmf, observed_pmf, width=2)

decorate(xlabel="Размер группы", ylabel="Вероятность")

В наблюдаемом распределении малочисленных групп меньше, а многочисленных — больше. Смещенное (biased) среднее значение равно 29.1, что почти на 25 % выше фактического среднего:

observed_pmf.mean()

29.123376623376622

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

Результат будет смещенным по причинам, которые мы только что рассмотрели, но его можно использовать для оценки фактического распределения. Вот функция, которая устраняет смещение Pmf путем деления вероятности на размер группы:

def unbias(pmf, name):

    # разделить каждую вероятность на размер группы

    ps = pmf.ps / pmf.qs

 

    new_pmf = Pmf(ps, pmf.qs, name=name)

    new_pmf.normalize()

    return new_pmf

Вот что мы получаем:

debiased_pmf = unbias(observed_pmf, "debiased")

debiased_pmf.mean()

23.692307692307693

Среднее значение Pmf после устранения смещения совпадает со средним значением фактического распределения, с которого мы начинали.

Если вас заинтересовал этот пример, то вам также может понравиться глава 2 из моей книги Probably Overthinking It (University of Chicago Press, 2023), в которой приводится этот и несколько других примеров так называемого парадокса проверки.

Данные NSFG

В предыдущей главе мы строили частотные диаграммы продолжительности беременности при рождении первых и последующих детей. Но размеры этих групп различаются, поэтому мы не можем сравнивать такие таблицы частот непосредственно. Но, поскольку функции вероятности (ФВ) нормированные, их можно подвергать сравнению. Поэтому снова загрузим данные исследования NSFG и, чтобы построить распределения продолжительности беременности, создадим объекты Pmf.

В модуле nsfg реализована функция read_nsfg_groups, которая считывает данные, выбирает строки, представляющие живорождение (live), а затем разбивает их на соответствующие первенцам (firsts) и последующим (others) детям строки. Функция возвращает три DataFrame-объекта:

from nsfg import get_nsfg_groups

 

live, firsts, others = get_nsfg_groups()

Теперь можно воспользоваться переменными firsts и others, чтобы создать объекты Pmf для ФВ продолжительности беременности в каждой группе:

first_pmf = Pmf.from_seq(firsts["prglngth"], name="первые")

other_pmf = Pmf.from_seq(others["prglngth"], name="последующие")

Вот столбчатые диаграммы ФВ для первых и последующих детей:

two_bar_plots(first_pmf, other_pmf)

decorate(xlabel="Недели", ylabel="Вероятность", xlim=[20, 50])

Построив диаграммы ФВ вместо диаграмм частот, мы можем сравнить два распределения, не будучи введенными в заблуждение разницей в размерах выборок. Если посмотреть на рисунок, может показаться, что первые дети с меньшей вероятностью, чем последующие, появляются на свет вовремя (39-я неделя) и с большей вероятностью — позже срока (41-я и 42-я недели).

Дополнительная визуализация

Таблицы частот и функции вероятности полезны на этапе разведки данных для выявления закономерностей и взаимосвязей. Когда у вас сформируется понимание происходящего, создайте такую визуализацию, которая представит выявленные вами закономерности наиболее отчетливо.

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

Когда объект Pmf вызывается как функция, можно задать на входе последовательность величин и получить на выходе последовательность вероятностей:

weeks = range(35, 46)

first_pmf(weeks)

array([0.03602991, 0.03897575, 0.04713347, 0.06163608, 0.4790392,

       0.12145932, 0.08157716, 0.04645366, 0.01971448, 0.00521187,

       0.00135962])

other_pmf(weeks)

array([0.03210137, 0.03146779, 0.05216473, 0.07074974, 0.54466737,

       0.12249208, 0.04794087, 0.02597677, 0.01288279, 0.00485744,

       0.00084477])

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

diffs = first_pmf(weeks) — other_pmf(weeks)

diffs

array([0.00392854, 0.00750796, -0.00503126, -0.00911366, -0.06562817,

       -0.00103276, 0.03363629, 0.02047689, 0.00683169, 0.00035443,

       0.00051485])

Вот как они выглядят после умножения на 100, чтобы выразить разницу в процентах:

plt.bar(weeks, diffs * 100)

decorate(xlabel="Недели", ylabel="Разница (%)")

Этот график делает картину более ясной: первые дети реже рождаются на 39-й неделе и несколько чаще — на 41-й и 42-й.

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

Глоссарий

В этой главе не так много новых терминов, как в предыдущих главах:

Нормированный (normalized)

Множество вероятностей нормированное, если эти вероятности в сумме дают 1.

Функция вероятности (probability mass function, PMF)

Функция, описывающая распределение вероятностей путем установления соответствия каждой величины с ее вероятностью.

Упражнения

Для упражнений в этой главе используйте файл респондентов NSFG, в котором каждая строка соответствует одной респондентке. Инструкции по загрузке данных приведены в Jupyter-блокноте для этой главы.

В модуле nsfg.py реализована функция, которая считывает файл респондентов и возвращает объект DataFrame:

from nsfg import read_fem_resp

 

resp = read_fem_resp()

resp.shape

(7643, 3092)

Этот DataFrame содержит 7643 строки и 3092 столбца.

Упражнение 3.1

Выберите столбец numbabes, в котором записано «количество младенцев, родившихся живыми» у каждой респондентки. Создайте объект FreqTab и постройте частотную диаграмму значений в этом столбце. Убедитесь, что они соответствуют частотам, указанным в кодбуке. Есть ли какие-то специальные значения, которые следует заменить на NaN?

Теперь создайте объект Pmf и отобразите его в виде столбчатой диаграммы. Является ли распределение симметричным, скошенным влево или вправо?

Упражнение 3.2

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

Для имеющейся выборки асимметрию можно вычислить, сложив возведенные в куб отклонения и разделив полученную сумму на стандартное отклонение в кубе. Например, вот как рассчитать асимметрию переменной numbabes:

numbabes = resp["numbabes"].replace(97, np.nan)

deviations = numbabes - numbabes.mean()

skewness = np.mean(deviations**3) / numbabes.std(ddof=0) ** 3

skewness

1.7018914266755958

Положительное значение указывает на скошенность распределения вправо, а отрицательное — на скошенность влево.

Если у вас вместо последовательности значений есть объект Pmf, то асимметрию можно вычислить следующим образом:

1. Вычислите отклонение каждой величины в Pmf от среднего значения.

2. Возведите отклонения в куб, умножьте на вероятности, записанные в Pmf, и сложите полученные результаты.

3. Разделите полученную сумму на стандартное отклонение в третьей степени.

Напишите функцию pmf_skewness, которая принимает на входе объект Pmf и возвращает его асимметрию. Используйте свою функцию и объект Pmf переменной numbabes для вычисления асимметрии. Убедитесь, что полученный результат соответствует значению, вычисленному выше.

Упражнение 3.3

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

В датафрейме resp выберите столбец numkdhh, содержащий количество детей в возрасте до 18 лет в домохозяйстве каждой респондентки. На основе значений этого столбца создайте объект Pmf.

Используйте функцию bias для вычисления распределения, которое бы мы увидели, если бы у детей респонденток спросили, сколько детей младше 18 лет (включая самих опрошенных) проживает в их семьях.

Постройте диаграммы фактического и смещенного распределений, вычислите их средние значения.

---