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

13

Я провел множественную регрессию, в которой модель в целом значима и объясняет около 13% дисперсии. Тем не менее, мне нужно найти величину дисперсии, объясняемой каждым значимым предиктором. Как я могу сделать это с помощью R?

Вот некоторые примеры данных и кода:

D = data.frame(
    dv = c( 0.75, 1.00, 1.00, 0.75, 0.50, 0.75, 1.00, 1.00, 0.75, 0.50 ),
    iv1 = c( 0.75, 1.00, 1.00, 0.75, 0.75, 1.00, 0.50, 0.50, 0.75, 0.25 ),
    iv2 = c( 0.882, 0.867, 0.900, 0.333, 0.875, 0.500, 0.882, 0.875, 0.778, 0.867 ),
    iv3 = c( 1.000, 0.067, 1.000, 0.933, 0.875, 0.500, 0.588, 0.875, 1.000, 0.467 ),
    iv4 = c( 0.889, 1.000, 0.905, 0.938, 0.833, 0.882, 0.444, 0.588, 0.895, 0.812 ),
    iv5 = c( 18, 16, 21, 16, 18, 17, 18, 17, 19, 16 ) )
fit = lm( dv ~ iv1 + iv2 + iv3 + iv4 + iv5, data=D )
summary( fit )

Вот вывод с моими фактическими данными:

Call: lm(formula = posttestScore ~ pretestScore + probCategorySame + 
    probDataRelated + practiceAccuracy + practiceNumTrials, data = D)

Residuals:
    Min      1Q  Median      3Q     Max 
-0.6881 -0.1185  0.0516  0.1359  0.3690 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)
 (Intercept)        0.77364    0.10603    7.30  8.5e-13 ***
 iv1                0.29267    0.03091    9.47  < 2e-16 ***
 iv2                0.06354    0.02456    2.59   0.0099 **
 iv3                0.00553    0.02637    0.21   0.8340
 iv4               -0.02642    0.06505   -0.41   0.6847
 iv5               -0.00941    0.00501   -1.88   0.0607 .  
--- Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Residual standard error: 0.18 on 665 degrees of freedom
 Multiple R-squared:  0.13,      Adjusted R-squared:  0.123
 F-statistic: 19.8 on 5 and 665 DF,  p-value: <2e-16

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

baixiwei
источник
2
Вы смотрели на ответ Джерома Англима здесь ?
Стат
Да, это был дополнительный ответ, о котором я говорил. Я надеялся на что-то более конкретное и пошаговое. Я загрузил ppcor, но не знал, что делать с выходом spcor. Кроме того, мне интересно, есть ли способ сделать это в ядре R? Кажется, это достаточно распространенная задача, для которой не требуется специальный пакет.
baixiwei
Кратчайший ответ на ваш вопрос о коррелированных предикторах заключается в том, что их отдельная важность не может быть определена количественно, по крайней мере, без дополнительных предположений и приближений. Подумайте об этом так: если это просто, почему это не легко и легко доступно, потому что многие исследователи думают, что они этого хотят?
Ник Кокс
Я бы посоветовал взглянуть на relaimpoпакет и сопровождающий его документ: jstatsoft.org/index.php/jss/article/view/v017i01/v17i01.pdf Я часто использую метод «LMG».
Фил

Ответы:

15

Объясненный процент зависит от введенного заказа.

Если указать конкретный заказ, вы можете вычислить это тривиально в R (например , через updateи anovaфункции, смотри ниже), а другой порядок въезда даст потенциально очень разные ответы.

[Одной из возможностей может быть усреднение по всем заказам или что-то в этом роде, но это может стать громоздким и может не дать ответа на особенно полезный вопрос.]

-

Как указывает Стат, в одной модели, если вы хотите использовать одну переменную за раз, вы можете просто использовать 'anova' для создания таблицы добавочных сумм квадратов. Это будет следовать из вашего кода:

 anova(fit)
Analysis of Variance Table

Response: dv
          Df   Sum Sq  Mean Sq F value Pr(>F)
iv1        1 0.033989 0.033989  0.7762 0.4281
iv2        1 0.022435 0.022435  0.5123 0.5137
iv3        1 0.003048 0.003048  0.0696 0.8050
iv4        1 0.115143 0.115143  2.6294 0.1802
iv5        1 0.000220 0.000220  0.0050 0.9469
Residuals  4 0.175166 0.043791        

-

Таким образом, у нас есть объяснение возрастающей дисперсии; как мы получаем пропорцию?

Довольно тривиально, масштабируйте их на 1, деленное на их сумму. (Замените 1 на 100 для объяснения процентной дисперсии.)

Здесь я отобразил его как добавленный столбец к таблице anova:

 af <- anova(fit)
 afss <- af$"Sum Sq"
 print(cbind(af,PctExp=afss/sum(afss)*100))
          Df       Sum Sq      Mean Sq    F value    Pr(>F)      PctExp
iv1        1 0.0339887640 0.0339887640 0.77615140 0.4280748  9.71107544
iv2        1 0.0224346357 0.0224346357 0.51230677 0.5137026  6.40989591
iv3        1 0.0030477233 0.0030477233 0.06959637 0.8049589  0.87077807
iv4        1 0.1151432643 0.1151432643 2.62935731 0.1802223 32.89807550
iv5        1 0.0002199726 0.0002199726 0.00502319 0.9468997  0.06284931
Residuals  4 0.1751656402 0.0437914100         NA        NA 50.04732577

-

Если вы решите, что вам нужно несколько конкретных порядков ввода, вы можете сделать что-то еще более общее, например, это (что также позволяет вам вводить или удалять группы переменных одновременно, если вы хотите):

 m5 = fit
 m4 = update(m5, ~ . - iv5)
 m3 = update(m4, ~ . - iv4)
 m2 = update(m3, ~ . - iv3)
 m1 = update(m2, ~ . - iv2)
 m0 = update(m1, ~ . - iv1)

 anova(m0,m1,m2,m3,m4,m5)
Analysis of Variance Table

Model 1: dv ~ 1
Model 2: dv ~ iv1
Model 3: dv ~ iv1 + iv2
Model 4: dv ~ iv1 + iv2 + iv3
Model 5: dv ~ iv1 + iv2 + iv3 + iv4
Model 6: dv ~ iv1 + iv2 + iv3 + iv4 + iv5
  Res.Df     RSS Df Sum of Sq      F Pr(>F)
1      9 0.35000                           
2      8 0.31601  1  0.033989 0.7762 0.4281
3      7 0.29358  1  0.022435 0.5123 0.5137
4      6 0.29053  1  0.003048 0.0696 0.8050
5      5 0.17539  1  0.115143 2.6294 0.1802
6      4 0.17517  1  0.000220 0.0050 0.9469

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

... а затем масштабировать до процентов, как и раньше.

(NB. Тот факт, что я объясняю, как делать эти вещи, не обязательно должен восприниматься как пропаганда всего, что я объясняю.)

Glen_b - Восстановить Монику
источник
2
р2anova(fit)м0м5
Этот пересмотренный ответ действительно полезен. Я думаю, что я получаю там. Один вопрос: если я вычислю пропорцию дисперсии, объясненную для iv5 (последней переменной) описанным вами способом, будет ли это математически таким же, как разница в значениях R ^ 2, полученных в результате суммирования, применяемых к модели, соответствует и без iv5? Я на самом деле получаю одинаковые значения и просто хотел проверить, являются ли они концептуально одинаковыми.
baixiwei
И еще один вопрос: есть ли причина, по которой я не смог сделать то, что только что описал в предыдущем комментарии, один раз для каждой из двух разных iv? Будет ли это эквивалентно вашему второму предложенному методу, включающему различные порядки ввода переменных?
baixiwei
baixiwei - да, эти два различия в р2будет идентичным Да, вы можете сделать «подходящий последний» вклад для каждого. (Фактически, t-критерий в таблице регрессии summary.lmдает вам тест значимости каждой переменной, которая уже установлена ​​последней).
Glen_b -Reinstate Monica
2

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

https://www.researchgate.net/publication/306347340_A_Natural_Decomposition_of_R2_in_Multiple_Linear_Regression

user128460
источник
3
user128460 приветствуется, но это сайт вопросов и ответов, а не сайт вопросов и ссылок.
Роберт Лонг
Разве это не оценка Пратта?
Бретт
2

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

library(hier.part)
env <- D[,2:5]
all.regs(D$dv, env, fam = "gaussian", gof = "Rsqu",
     print.vars = TRUE)
MFR
источник