Discussion :

[Peio64]

Bonjour,

Je suis débutant en R et je souhaiterais écrire un script me permettant des faire des ANOVA en routine sur l'ensemble de mes variables, ainsi que de récupérer certains résultats sous forme de tableau à la fin.


Mon modèle est par exemple le suivant:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

> modele <- lm(taille ~ id_genotype + ht05cov:id_genotype)

j'obtiens les résultats de l'anova type I comme ceci:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
> anova(modele)
Analysis of Variance Table
 
Response: taille
                    Df  Sum Sq Mean Sq F value   Pr(>F)   
id_genotype          7  544.84   77.83  2.3197 0.042340 * 
id_genotype:ht05cov  7  730.36  104.34  3.1096 0.009592 **
Residuals           43 1442.81   33.55                    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Je veux aussi les résultats de type III, j'ai trouvé ceci:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
> drop1(modele,~.,test="F")
Single term deletions
 
Model:
taille ~ id_genotype + ht05cov:id_genotype
                    Df Sum of Sq     RSS     AIC F value    Pr(F)   
<none>                           1442.81  216.41                    
id_genotype          6    779.38 2222.19  229.46  3.8713 0.003551 **
id_genotype:ht05cov  7    730.36 2173.17  226.16  3.1096 0.009592 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

J'aimerai donc écrire un script
- qui réaliserait ces opérations à la suite pour toutes les variables de mes données
- qui récupère ensuite certains de ces résultats et qui les liste sous cette forme (le R² étant le rapport entre le SS du facteur sur le SS total):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
		Type	Df	R²modèle	R²id_génotype	F	p-val
variable1	I	xx	xx		xx		xx	xx
variable1	III	xx	xx		xx		xx	xx
variable2	I	xx	xx		xx		xx	xx
variable2	III	xx	xx		xx		xx	xx
variable3	I	xx	xx		xx		xx	xx
variable3	III	xx	xx		xx		xx	xx

Auriez-vous quelques pistes pour m'aiguiller?

Merci bcp

Pierre

[pitipoisson]

Bonjour et bienvenue,

Je peux te donner quelques pistes, en espèrant que la réponse n'est pas trop tardive.

Pour ce qui est de faire des modèles en routine, en supposant que tes noms de variable (réponse ?) soient dans un vecteur variables, tu peux faire une boucle :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
for (i in variables)
{
    formula <- as.formula(paste(variables[i], " ~ id_genotype + ht05cov:id_genotype", sep=""))
    resI <- lm(formula)
    resIII <- drop1(resI,~.,test="F")
 
}

À adapter évidemment, si ce ne sont pas les variables dépendantes qui changent.

Pour pouvoir récupérer les résultats (toujours dans la boucle), un exemple de la méthode sera sans doute plus parlant :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
> x <- 1:10
 
> y <- 2 * x - 0.5 +rnorm(10)
 
> reslm <- lm(y ~ x)
 
> reslm
 
Call:
lm(formula = y ~ x)
 
Coefficients:
(Intercept)            x  
     -1.134        2.072  
 
> summary(reslm)
 
Call:
lm(formula = y ~ x)
 
Residuals:
    Min      1Q  Median      3Q     Max 
-1.8994 -1.0569  0.5238  0.7240  1.2219 
 
Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -1.1343     0.8549  -1.327    0.221    
x             2.0720     0.1378  15.038 3.78e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
 
Residual standard error: 1.252 on 8 degrees of freedom
Multiple R-squared: 0.9658,	Adjusted R-squared: 0.9616 
F-statistic: 226.1 on 1 and 8 DF,  p-value: 3.779e-07 
 
> names(summary(reslm))
 [1] "call"          "terms"         "residuals"     "coefficients" 
 [5] "aliased"       "sigma"         "df"            "r.squared"    
 [9] "adj.r.squared" "fstatistic"    "cov.unscaled" 
 
> summary(reslm)$adj.r.squared
[1] 0.961561
 
## Pareil pour anova()
> names(anova(reslm))
[1] "Df"      "Sum Sq"  "Mean Sq" "F value" "Pr(>F)" 
 
> anova(reslm)$Df
[1] 1 8