exponenta event banner

ANOVA со случайными эффектами

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

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

Настройка модели

Загрузите образцы данных.

load mileage

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

Создайте массив, указывающий фабрику для каждого значения в пробеге. Этот массив равен 1 для первого столбца, 2 для второго и 3 для третьего.

factory  = repmat(1:3,6,1);

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

carmod = [ones(3,3); 2*ones(3,3)];

Превратите эти матрицы в векторы и отобразите их.

mileage = mileage(:);
factory = factory(:);
carmod = carmod(:);
[mileage factory carmod]
ans = 18×3

   33.3000    1.0000    1.0000
   33.4000    1.0000    1.0000
   32.9000    1.0000    1.0000
   32.6000    1.0000    2.0000
   32.5000    1.0000    2.0000
   33.0000    1.0000    2.0000
   34.5000    2.0000    1.0000
   34.8000    2.0000    1.0000
   33.8000    2.0000    1.0000
   33.4000    2.0000    2.0000
      ⋮

Подгонка модели случайных эффектов

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

[pvals,tbl,stats] = anovan(mileage, {factory carmod}, ... 
'model',2, 'random',1,'varnames',{'Factory' 'Car Model'});

Figure N-Way ANOVA contains objects of type uicontrol.

В версии с фиксированными эффектами этой посадки, которую можно получить, пропустив входные данные 'random',1 в предыдущем коде эффект модели автомобиля значителен, с p-значением 0,0039. Но в этом примере, который учитывает случайное изменение эффекта переменной 'Car Model' от одной фабрики к другой эффект все еще значительный, но с более высоким p-значением 0,0136.

F-Статистика для моделей со случайными эффектами

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

В примере, описанном в разделе Настройка модели, эффект переменной 'Factory' может варьироваться в зависимости от моделей автомобилей. В этом случае среднее квадрат взаимодействия занимает место среднего квадрата ошибки в F-статистике.

Найдите F-статистику.

F = 26.6756 / 0.02
F = 1.3338e+03

Степени свободы для статистики - это степени свободы для числителя (2) и знаменателя (2) средних квадратов.

Найдите значение p.

pval = 1 - fcdf(F,2,2)
pval = 7.4919e-04

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

Найдите коэффициенты этих линейных комбинаций.

stats.ems
ans = 4×4

    6.0000    0.0000    3.0000    1.0000
    0.0000    9.0000    3.0000    1.0000
    0.0000    0.0000    3.0000    1.0000
         0         0         0    1.0000

Это возвращает ems области stats структура.

Отображение текстовых представлений линейных комбинаций.

stats.txtems
ans = 4x1 cell
    {'6*V(Factory)+3*V(Factory*Car Model)+V(Error)'  }
    {'9*Q(Car Model)+3*V(Factory*Car Model)+V(Error)'}
    {'3*V(Factory*Car Model)+V(Error)'               }
    {'V(Error)'                                      }

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

В некоторых случаях не существует единого термина, ожидаемое значение которого соответствует значению, требуемому для знаменателя F-статистики. В этом случае знаменателем является линейная комбинация средних квадратов. Структура stats содержит поля, дающие определения знаменателей для каждой F-статистики. txtdenom поле, stats.txtdenom, содержит текстовое представление и denom содержит матрицу, которая определяет линейную комбинацию отклонений членов в модели. Для сбалансированных моделей, таких как эта, denom матрица, stats.denom, содержит нули и единицы, потому что знаменатель - это только один средний квадрат члена.

Отображение txtdenom поле.

stats.txtdenom
ans = 3x1 cell
    {'MS(Factory*Car Model)'}
    {'MS(Factory*Car Model)'}
    {'MS(Error)'            }

Отображение denom поле.

stats.denom
ans = 3×3

    0.0000    1.0000   -0.0000
    0.0000    1.0000   -0.0000
    0.0000   -0.0000    1.0000

Компоненты расхождения

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

Отображение имен случайных терминов.

stats.rtnames
ans = 3x1 cell
    {'Factory'          }
    {'Factory*Car Model'}
    {'Error'            }

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

Просмотрите оценку компонента отклонения для каждого члена.

stats.varest
ans = 3×1

    4.4426
   -0.0313
    0.1139

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

Создайте гистограмму компонентов.

bar(max(0,stats.varest))

Figure contains an axes. The axes contains an object of type bar.

gca.xtick = 1:3;
gca.xticklabel = stats.rtnames;

Можно также вычислить доверительные границы для оценки дисперсии. anovan функция делает это путем вычисления доверительных границ для ожидаемого среднего квадрата дисперсии и нахождения нижнего и верхнего пределов для каждой составляющей дисперсии, содержащей все эти границы. Эта процедура приводит к набору ограничений, которые являются консервативными для сбалансированных данных. (То есть 95% доверительные границы будут иметь вероятность по меньшей мере 95% содержать истинные дисперсии, если количество наблюдений для каждой комбинации переменных группировки одинаково.) Для несбалансированных данных это аппроксимации, которые не гарантированно являются консервативными.

Просмотрите оценки дисперсии и доверительные пределы для оценок дисперсии каждого компонента.

[{'Term' 'Estimate' 'Lower' 'Upper'};
 stats.rtnames, num2cell([stats.varest stats.varci])]
ans=4×4 cell array
    {'Term'             }    {'Estimate'}    {'Lower' }    {'Upper'   }
    {'Factory'          }    {[  4.4426]}    {[1.0736]}    {[175.6038]}
    {'Factory*Car Model'}    {[ -0.0313]}    {[   NaN]}    {[     NaN]}
    {'Error'            }    {[  0.1139]}    {[0.0586]}    {[  0.3103]}