Глава 5. Задачи о выборках: двухфакторный дисперсионный анализ

Дальше

Назад

Начало

Конец

Список

5.3. Теория и практика двухфакторного дисперсионного анализа

Формулировка задачи

Пусть имеется выборка из m значений X₁, X₂, …, X_m некоторого гидробиологического показателя, измеренного в количественной шкале. Изучается влияние, которое оказывают на отклик два качественных признака – фактор A, который имеет k уровней (градаций) A₁, ..., A_k , и фактор B, разбитый на n уровней B₁, ..., B_n

Необходимо:

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

если нулевая гипотеза отвергается, то оценить степень влияния каждого фактора;

выделить все пары групп {i – j}, i = 1,2,…,k× n, j = 1,2,…,k× n, i ¹ j, которые имеют между собой статистически значимые различия.

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

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

Рекомендуемая литература та же, что и для раздела 5.2.

Математический лист

Предполагается, что измеряемая величина X есть результат действия факторов A и B и составляющей e , которая является независимой, нормально N(0,s ²) распределенной случайной величиной: .

Принимается аддитивная и независимая модель действия факторов:

, причем , . (5.43)

Величины a_j и b_i называются вкладами факторов. Последние два условия всегда можно выполнить масштабированием величин a_j и b_i за счет изменения величины c.

Для каждого наблюдения из рассматриваемой совокупности справедливо уравнение:

x_ij = c + a_j+ b_i + e _ij, i =1, ..., n; j =1, ..., k. (5.44)

Обычно наблюдения представляют структурной таблицей статистического комплекса. Приведем простейший двухфакторный комплекс, в которой каждому сочетанию (A_j, B_j) уровней (градаций) факторов, т.е. одной клетке таблицы, соответствует одно наблюдение (в таблице сочетание символов “( )^” обозначает статистическую оценку групповых средних):

Фактор B	Фактор A A₁ A₂ ... A_k	Средние по строкам (оценки вкладов B)
B₁ B₂ ... B_n	x₁₁ x₁₂ ... x_1k x₂₁ x₂₂ ... x_2k ... x_n1 x_n2 ... x_nk	x_1·=(c+b₁)^ x_2·=(c+b₂)^ ... x_n·=(c+b_n)^
Средние по столбцам (оценки вкладов A)	x· ₁= x· ₂= x· _k= (c+a₁)^ (c+a₂)^ (c+a_k)^	x· · =c^

Основное тождество дисперсионного анализа

Оценки c, b_i, a_j могут быть получены с помощью метода наименьших квадратов (МНК) минимизацией суммы (5.45)

Основываясь на МНК-оценках

, ,

, (5.46)

введем следующие обозначения:

для сумм квадратов отклонений под влиянием k уровней фактора А и n уровней фактора В:

, ;(5.47)

для остаточной суммы квадратов:

;(5.48)

для полной суммы квадратов наблюдений относительно общего среднего :

.(5.49)

Тогда справедливо следующее соотношение:

, (5.50)

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

Проверка нулевых гипотез

По имеющимся наблюдениям требуется проверить предположение об отсутствии влияния фактора A (или B) на результат измерения, т.е. проверить гипотезы

H_A: a₁ = a₂ = . . . = a_k = 0

H_B: b₁ = b₂ = . . . = b_n = 0 .

Основой процедуры проверки гипотезы является сравнение двух статистически независимых оценок дисперсии s ².

Одна из них, s ^2* оценивает дисперсию вне зависимости от того, верна или нет гипотеза H_A (или H_В) и основана на сумме квадратов случайных отклонений:

. (5.51)

Другая, s ^2** оценивает дисперсию, если H_A (или H_В) верна. Для гипотезы H_A эта дисперсия основана на сумме квадратов разностей “между столбцами”, т.е. по уровням фактора A:

. (5.52)

Если гипотеза H_A верна, то отношение

(5.53)

имеет F-распределение Фишера с (k -1) и r степенями свободы. Если

F_A ³ F_1-a , (5.54)

где F_1-a – квантиль этого распределения порядка 1- a , a – выбранный уровень значимости, то гипотеза H_A отклоняется.

Вместо (5.54) можно использовать эквивалентную процедуру: гипотеза H_A отклоняется, если

P{ F ³ F_A } £ a, (5.55)

где P{ F ³ F_A } – вероятность при справедливости H_A получить значение, большее, чем F_A; F – случайная величина, имеющая распределение Фишера.

Для проверки гипотезы H_В используют сумму квадратов разностей "между строками", то есть по уровням фактора B: .(5.56)

Аналогичным образом, если отношение велико, то гипотеза H_B отклоняется.

Результаты расчетов

Формирование речных сообществ, их видовой состав и продуктивность гидробионтов находятся в постоянной динамике под воздействием большого количества экологических факторов. Будем оценивать временнỳю изменчивость сообществ, когда набор видов и количество особей претерпевают ежедневные, сезонные и многолетние колебания под воздействием температурного режима, мгновенной гидродинамики водотоков, сбросов токсикантов и т.д. Для прогнозирования динамики экологических систем методами статистического анализа существует обширный и специализированный математический аппарат, который остается за рамками настоящего изложения. Подробно теоретические и прикладные вопросы анализа хронологических трендов экосистем с использованием функциональных предикторов временных рядов рассматривались нами ранее [Розенберг с соавт., 1994]. Здесь мы ограничимся изложением частной методики дисперсионного анализа для оценки влияния многолетней и сезонной составляющих на общие показатели обилия зообентоса.

На основе данных из базы по малым рекам Самарской области сформируем выборку из суммарных значений численности N_s (экз./м²) и биомассы В_s (г/м²) зообентоса и индекса разнообразия Шеннона H для каждой пробы наблюдений. Данные по численности и биомассе предварительно прологарифмируем. Каждому значению варьируемой переменной поставим в соответствие три фактора: водоток, из которого взяты пробы, год исследования (с 1988-99 гг.) и порядковый номер месяца (с 5 по 9) отбора пробы.

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

Таблица 5.8

Результаты проверки гипотезы о влиянии места отбора проб на логарифм численности зообентоса по F-критерию

Источник вариации	Сумма квадратов	Степеней свободы	Средние квадраты	F-отношение	p-значениe
Между группами	85.58	32	2.67	0.99	0.49 (фактор не значим)
Внутри групп	1357.12	500	2.71
Итого	1442.71	532

Оговоримся, что проблемы пространственной изменчивости, отражающей распределение показателей обилия в зависимости от географических координат поверхности, неизмеримо сложнее описываемого примера, поэтому в данном контексте речь идет не об оценке влияния места отбора пробы вообще, а о конкретном разбиении влияющего фактора на градации. Во-первых, проблематично само понятие “река”, как средство обобщения гидробиологических данных, т.к. вариации данных между станциями одной реки, как правило, превышают межрегиональную вариацию (см. фиг. “а” рис. 5.10). Во-вторых, сама природа дисперсионного анализа предполагает тенденцию к недооценке влияния фактора при числе градаций больше 10. Однако выводы табл. 5.8 дают нам формальные основания осуществить дисперсионный анализ в градациях остальных двух факторов: “год” – “месяц”, основные результаты которого приведены в табл. 5.9.

Рис 5.10 Диаграммы изменения групповых средних по результатам дисперсионного анализа показателей обилия и биоразнообразия зообентоса

Таблица 5.9

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

Выборки для анализа	Источник вариации	Сумма квадратов	Степеней свободы	Средние квадраты	F-отноше-ния	p-значения
1. Численность зообентоса в пробе (логарифм *N_s*)	Фактор “Месяц”	24.23	4	6.06	2.58	0.0365
	Фактор “Год”	208.12	13	16.01	6.82	0.0
	Остаточная	1208.84	515	2.35
	Всего	1442.71	532
2. Численность Chironomidae в пробе (логарифм *N_s*)	Фактор “Месяц”	25.65	4	6.41	2.23	0.0651
	Фактор “Год”	80.38	13	6.18	2.15	0.0108
	Остаточная	1341.48	446	2.87
	Всего	1455.08	463
3. Индекс Шеннона H (по численности зообентоса)	Фактор “Месяц”	11.93	4	2.98	3.95	0.0036
	Фактор “Год”	46.6	13	3.58	4.74	0.0
	Остаточная	388.6	514	0.756
	Всего	442.16	531
4. Биомасса зообентоса в пробе (логарифм *B_s*)	Фактор “Месяц”	40.08	4	10.02	3.09	0.0157
	Фактор “Год”	203.16	13	15.63	4.82	0.0
	Остаточная	1451.89	448	3.24
	Всего	1694.25	465
5. Индекс плотности населения (*N_s***B_s*)^1/2	Фактор “Месяц”	145189	4	36297	1.28	0.27
	Фактор “Год”	598507	13	46039	1.62	0.0762
	Остаточная	12727600	448	28409
	Всего	13472500	465
6. Логарифм индекса плотности населения ln((*N_s***B_s*)^1/2)	Фактор “Месяц”	22.49	4	5.62	2.09	0.081
	Фактор “Год”	119.1	13	9.16	3.40	0.0001
	Остаточная	1206.8	448	2.69
	Всего	1356.6	465

По существу рассматриваемой проблемы выполненные расчеты влияния факторов позволяют сделать следующие выводы:

Во всех случаях (выборки 1-4, 6 в табл. 5.9) отчетливо проявляется влияние многолетнего фактора. Динамика изменения групповых средних по годам для численности зообентоса представлена на рис. 5.10, фиг. “б”. Однако вряд ли следует понимать это влияние слишком буквально, как результат изменения гидрологических или радиационно-климатических условий. Можно, например, предложить следующую версию влияния многолетнего фактора. В 1987-91 гг. в коллектив исследователей пришла группа ученых, владеющих современными методами биологического контроля, что существенно расширило диапазон исследований и создало эффект роста численности определяемого материала. Начиная с 1992 г., наблюдалось резкое снижение уровня финансирования академических исследований, в результате чего сузилась широта географического охвата, и внимание исследователей сконцентрировалось на водоемах, подверженных экстремальным антропогенным воздействиям. Рост обилия гидробионтов с 1997 г. можно объяснить как расширением региона исследований на более "чистые" зоны, так и эффектом постепенного самоочищения водотоков после экономической стагнации 1991-97 гг. В любом случае, многолетний фактор является сложным конгломератом трудно идентифицируемых "человеческих", экономических и гидрометеорологических составляющих.

В очередной раз подтверждается необходимость обращать внимание на закон распределения зависимого показателя: влияние факторов на индекс плотности населения (N_s*B_s)^1/2 не было обнаружено, но проявилось после логарифмирования отклика.

С чисто формальных позиций, влияние сезонного фактора, выраженного календарным месяцем отбора пробы, является значимым для индекса Шеннона, численности и биомассы зообентоса (выборки 1, 3, 4), но незначимо для численности хирономид и индекса плотности населения (выборки 2, 6). Более детальные содержательные сведения можно получить из анализа сдвигов групповых средних, представленных на рис. 5.10. Месяц май характеризуют высокие значения численности N_sи биоразнообразия по Шеннону H, обусловленные развитием личинок амфибиотических насекомых весенней генерации. Начиная с июня, происходит постепенный захват экологических ниш видами-доминантами, в результате чего общая численность и, особенно резко, биоразнообразие Н начинают падать. В то же время, для июня характерна максимальная продуктивность и энергообмен в трофических цепях, выраженные через биомассу. Результатом майско-июньской конкуренции между N_s и B_sявилось сглаживание контрастов индекса плотности населения (ИПН). Увеличение индекса Шеннона в сентябре обусловлено развитием осенних генераций амфибиотических насекомых в бентосе рек на фоне небольших показателей их обилия в конце лета.

Представляет интерес сравнить результаты двухфакторного и однофакторного анализа, например, для выборки 2 (численность хирономид). Если оценить по Фишеру локальное влияние сезонного фактора, то его вклад в общий статистический разброс численности следует трактовать как значимый (статистика Фишера F = 2.80, а соответствующее ему значение p = 0.0256). При переходе к двухфакторной схеме оказалось, что часть дисперсии, приписываемой в однофакторной схеме сезонному фактору, на самом деле объясняется многолетней составляющей (например, в разные годы экспедиции проводились во время разных месяцев).

Хотя у нас нет формальных оснований отвергнуть нулевую гипотезу о равенстве групповых дисперсий по тестам Кохрена и Бартлетта для данного примера, однако значения уровней значимости для обоих критериев находились в слишком опасной близости от порогового значения (a = 0.05):

по тесту Кохрена: G-критерий = 0.259, p-значение = 0.07;

по тесту Бартлетта: М-критерий = 1.019, p-значение = 0.06.

Поскольку в рамках схемы двухфакторного дисперсионного анализа оценка различий групповых средних представляет собой нетривиальную задачу, выполним на основе однофакторного анализа с использованием метода Шеффе выделение в течение сезона гомогенных групп значений фактора, т.е. комбинаций месяцев, между которыми нет значимых отличий в средних значениях (см. табл. 5.10). Например, для численности хирономид (выборка 2 табл. 5.9) можно составить три последовательности месяцев, отмеченных затененными клетками таблицы 5.10, внутри которых данные можно считать статистически однородными: а) с мая по июль, б) с июня по август, в) с июля по сентябрь. Таким образом, не отвергнув нулевую гипотезу об отсутствии влияния фактора сезонности для этого показателя, мы можем ставить вопрос об объединении в гомогенные группы некоторые подмножества уровней фактора.

Таблица 5.10

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

Дальше

Назад

Начало

Конец

Список