Загрузите выборочные данные.

load flu

flu массив набора данных имеет Date переменная и 10 переменных, содержащих оцененные уровни гриппа (в 9 различных областях, оцененных от поисковых запросов Google®, плюс общенациональная оценка из Центра по контролю и профилактике заболеваний, CDC).

Чтобы подбирать линейно смешанную модель эффектов, ваши данные должны быть в правильно отформатированном массиве набора данных. Чтобы подбирать линейную модель смешанных эффектов с уровнями гриппа как ответы и область как переменный предиктор, объедините эти девять столбцов, соответствующих областям в массив. Новый массив набора данных, flu2, должен иметь переменную отклика, FluRate, номинальная переменная, Region, это показывает, какая область каждая оценка от, и сгруппированная переменная Date.

flu2 = stack(flu,2:10,'NewDataVarName','FluRate','IndVarName','Region');
flu2.Date = nominal(flu2.Date);

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

Область является категориальной переменной. Можно задать контрасты для категориальных переменных с помощью DummyVarCoding аргумент пары "имя-значение", подбирая модель. Когда вы не задаете контрасты, fitlme использует 'reference' контраст по умолчанию. Поскольку модель имеет прерывание, fitlme берет первую область, NE, как ссылка и создает восемь фиктивных переменных, представляющих другие восемь областей. Например, $$I[MidAtl]$$ фиктивная переменная, представляющая область MidAtl. Для получения дополнительной информации смотрите Фиктивные Переменные.

Соответствующая модель

где $$y_{im}$$ наблюдение $$i$$ для уровня $$m$$ из сгруппированной переменной Date, $${\beta }_j$$, $$j$$ = 0, 1..., 8, коэффициенты фиксированных эффектов, с $${\beta }_0$$ будучи коэффициентом для области NE. $$b_{0m}$$ случайный эффект для уровня $$m$$ из сгруппированной переменной Date, и $${\epsilon }_{im}$$ ошибка наблюдения для наблюдения $$i$$. Случайный эффект имеет предшествующее распределение, $$b_{0m}\sim N(0,{\sigma }_b^2)$$ и остаточный член имеет распределение, $${\epsilon }_{im}\sim N(0,{\sigma }^2)$$.

lme = fitlme(flu2,'FluRate ~ 1 + Region + (1|Date)')

lme = 
Linear mixed-effects model fit by ML

Model information:
    Number of observations             468
    Fixed effects coefficients           9
    Random effects coefficients         52
    Covariance parameters                2

Formula:
    FluRate ~ 1 + Region + (1 | Date)

Model fit statistics:
    AIC       BIC       LogLikelihood    Deviance
    318.71    364.35    -148.36          296.71  

Fixed effects coefficients (95% CIs):
    Name                        Estimate    SE          tStat      DF 
    {'(Intercept)'     }          1.2233    0.096678     12.654    459
    {'Region_MidAtl'   }        0.010192    0.052221    0.19518    459
    {'Region_ENCentral'}        0.051923    0.052221     0.9943    459
    {'Region_WNCentral'}         0.23687    0.052221     4.5359    459
    {'Region_SAtl'     }        0.075481    0.052221     1.4454    459
    {'Region_ESCentral'}         0.33917    0.052221      6.495    459
    {'Region_WSCentral'}           0.069    0.052221     1.3213    459
    {'Region_Mtn'      }        0.046673    0.052221    0.89377    459
    {'Region_Pac'      }        -0.16013    0.052221    -3.0665    459


    pValue        Lower        Upper    
     1.085e-31       1.0334       1.4133
       0.84534    -0.092429      0.11281
        0.3206    -0.050698      0.15454
    7.3324e-06      0.13424      0.33949
       0.14902     -0.02714       0.1781
    2.1623e-10      0.23655      0.44179
       0.18705    -0.033621      0.17162
       0.37191    -0.055948      0.14929
     0.0022936     -0.26276    -0.057514

Random effects covariance parameters (95% CIs):
Group: Date (52 Levels)
    Name1                  Name2                  Type           Estimate
    {'(Intercept)'}        {'(Intercept)'}        {'std'}        0.6443  


    Lower     Upper  
    0.5297    0.78368

Group: Error
    Name               Estimate    Lower      Upper
    {'Res Std'}        0.26627     0.24878    0.285

$$p$$- значения 7.3324e-06 и 2.1623e-10 соответственно показывают что фиксированные эффекты уровней гриппа в областях WNCentral и ESCentral существенно отличаются относительно уровней гриппа в области NE.

Пределы достоверности для стандартного отклонения термина случайных эффектов, $${\sigma }_b$$, не включайте 0 (0.5297, 0.78368), который указывает, что термин случайных эффектов является значительным. Можно также протестировать значение условий случайных эффектов с помощью compare метод.

Условное выражение соответствовало, ответ из модели при данном наблюдении включает вклады от фиксированных и случайных эффектов. Например, предполагаемый лучше всего линейный несмещенный предиктор (BLUP) уровня гриппа для области WNCentral на неделе 10/9/2005

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

beta = fixedEffects(lme);
[~,~,STATS] = randomEffects(lme); % Compute the random-effects statistics (STATS)
STATS.Level = nominal(STATS.Level);
y_hat = beta(1) + beta(4) + STATS.Estimate(STATS.Level=='10/9/2005')

y_hat = 1.2884

В предыдущем вычислении, beta(1) соответствует оценке для $${\beta }_0$$ и beta(4) соответствует оценке для $${\beta }_3$$. Можно просто отобразить подходящее значение с помощью fitted метод.

F = fitted(lme);
F(flu2.Date == '10/9/2005' & flu2.Region == 'WNCentral')

ans = 1.2884

Предполагаемый крайний ответ для области WNCentral на неделе 10/9/2005

Вычислите подходящий крайний ответ.

F = fitted(lme,'Conditional',false);
F(flu2.Date == '10/9/2005' & flu2.Region == 'WNCentral')

ans = 1.4602

Загрузите выборочные данные.

load(fullfile(matlabroot,'examples','stats','weight.mat'));

weight содержит данные из продольного исследования, где 20 предметов случайным образом присвоены 4 программам подготовки, и их потеря веса зарегистрирована более чем шесть 2-недельных периодов времени. Это - симулированные данные.

Храните данные в таблице. Задайте Subject и Program как категориальные переменные.

tbl = table(InitialWeight,Program,Subject,Week,y);
tbl.Subject = nominal(tbl.Subject);
tbl.Program = nominal(tbl.Program);

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

lme = fitlme(tbl,'y ~ InitialWeight + Program*Week + (Week|Subject)');

Вычислите подходящие значения и необработанные остаточные значения.

F = fitted(lme);
R = residuals(lme);

Постройте остаточные значения по сравнению с подходящими значениями.

plot(F,R,'bx')
xlabel('Fitted Values')
ylabel('Residuals')

Теперь постройте остаточные значения по сравнению с подходящими значениями, сгруппированными программой.

figure()
gscatter(F,R,Program)