Calculer la variance expliquée par chaque prédicteur en régression multiple en utilisant R

13

J'ai effectué une régression multiple dans laquelle le modèle dans son ensemble est significatif et explique environ 13% de la variance. Cependant, je dois trouver la quantité de variance expliquée par chaque prédicteur significatif. Comment puis-je faire cela en utilisant R?

Voici quelques exemples de données et de code:

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 )

Voici la sortie avec mes données réelles:

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

Cette question a été répondue ici , mais la réponse acceptée ne concerne que les prédicteurs non corrélés, et bien qu'il existe une réponse supplémentaire qui traite des prédicteurs corrélés, elle ne fournit qu'un indice général, pas une solution spécifique. Je voudrais savoir quoi faire si mes prédicteurs sont corrélés.

baixiwei
la source
2
Avez-vous regardé la réponse de Jeromy Anglim ici ?
Stat
Oui, c'est la réponse supplémentaire à laquelle je faisais référence. J'espérais quelque chose de plus spécifique et pas à pas. J'ai téléchargé ppcor mais je ne savais pas quoi faire de la sortie spcor. Aussi, je me demande s'il existe un moyen de le faire dans le noyau R? Cela semble être une tâche assez courante pour ne pas nécessiter de package spécial.
baixiwei
La réponse la plus courte à votre question sur les prédicteurs corrélés est que leur importance séparée ne peut pas être quantifiée, sans à tout le moins d'autres hypothèses et approximations. Considérez-le de cette façon: si c'est simple, pourquoi n'est-il pas facilement et facilement disponible, parce que de nombreux chercheurs pensent qu'ils le veulent?
Nick Cox
Je suggérerais de regarder dans le relaimpopaquet et son document d'accompagnement: jstatsoft.org/index.php/jss/article/view/v017i01/v17i01.pdf J'utilise fréquemment la méthode "LMG".
Phil

Réponses:

15

Le pourcentage expliqué dépend de la commande saisie.

Si vous spécifiez un ordre particulier, vous pouvez le calculer trivialement dans R (par exemple via les fonctions updateet anova, voir ci-dessous), mais un ordre d'entrée différent donnerait des réponses potentiellement très différentes.

[Une possibilité pourrait être de faire la moyenne de toutes les commandes ou quelque chose, mais cela deviendrait lourd et pourrait ne pas répondre à une question particulièrement utile.]

-

Comme le souligne Stat, avec un seul modèle, si vous recherchez une variable à la fois, vous pouvez simplement utiliser 'anova' pour produire la table des sommes incrémentielles des carrés. Cela découlerait de votre code:

 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        

-

Nous avons donc expliqué la variance incrémentielle; comment obtient-on la proportion?

Assez trivialement, mettez-les à l'échelle 1 divisé par leur somme. (Remplacez le 1 par 100 pour la variation en pourcentage expliquée.)

Ici, je l'ai affiché en tant que colonne ajoutée au tableau 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

-

Si vous décidez que vous voulez plusieurs ordres d'entrée particuliers, vous pouvez faire quelque chose d'encore plus général comme celui-ci (qui vous permet également d'entrer ou de supprimer des groupes de variables à la fois si vous le souhaitez):

 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

(Une telle approche peut également être automatisée, par exemple via des boucles et l'utilisation de get. Vous pouvez ajouter et supprimer des variables dans plusieurs commandes si nécessaire)

... puis évoluez en pourcentages comme précédemment.

(NB. Le fait que j'explique comment faire ces choses ne doit pas nécessairement être considéré comme le plaidoyer de tout ce que j'explique.)

Glen_b -Reinstate Monica
la source
2
R2anova(fit)m0m5
Cette réponse révisée est vraiment utile. Je pense que j'y arrive. Une question: si je calcule la proportion de variance expliquée pour iv5 (la dernière variable) de la manière que vous avez décrite, est-ce mathématiquement la même que la différence des valeurs R ^ 2 renvoyées par le résumé appliqué au modèle correspond avec et sans iv5? En fait, j'obtiens les mêmes valeurs et je voulais juste vérifier si ce sont conceptuellement la même chose.
baixiwei
Et encore une question: y a-t-il une raison pour laquelle je n'ai pas pu faire ce que je viens de décrire dans le commentaire précédent une fois pour chacun des deux IV différents? Serait-ce équivalent à votre deuxième méthode proposée impliquant différents ordres de saisie de variables?
baixiwei
R2summary.lm
2

J'ai prouvé que le pourcentage de variation expliqué par un prédicteur donné dans une régression linéaire multiple est le produit du coefficient de pente et de la corrélation du prédicteur avec les valeurs ajustées de la variable dépendante (en supposant que toutes les variables ont été normalisées pour avoir une moyenne nulle et variance un, qui est sans perte de généralité). Trouvez-le ici:

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

user128460
la source
3
Bienvenue à l'utilisateur 128460, mais il s'agit d'un site de questions et réponses, pas d'un site de questions et de liens vers des réponses.
Robert Long
N'est-ce pas le score de Pratt?
Brett
2

Vous pouvez utiliser la bibliothèque hier.part pour avoir des mesures d'ajustement pour les régressions d'une variable dépendante unique à toutes les combinaisons de N variables indépendantes

library(hier.part)
env <- D[,2:5]
all.regs(D$dv, env, fam = "gaussian", gof = "Rsqu",
     print.vars = TRUE)
MFR
la source