Discussion :

[CupOfAppleTea]

Bonjour,
Je dispose de 2 variables sur lesquelles je souhaite appliquer une régression linéaire. Vous trouverez en pièce jointe un extrait de mes données pour 2 mois, Août et Juin (sur 4 au total).

Voici le code pour l'analyse de ces 2 mois de données:

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
a=read.table("StatsDescript.txt", header=TRUE)
summary(a)
 
#JUIN
aJuin=subset(a, mois=="Juin")
modJuin=lm(aJuin$V2~aJuin$V1)
modJuin=aov(aJuin$V2~aJuin$V1)
plot(modJuin)
summary(modJuin)
shapiro.test(residuals(modJuin))
hist(residuals(modJuin))
bartlett.test(residuals(modJuin)~aJuin$V1)
 
#AOÛT
aAout=subset(a, mois=="Aout")
modAout=lm(aAout$V2~aAout$V1)
plot(modAout)
summary(modAout)
shapiro.test(residuals(modAout))
hist(residuals(modAout))
bartlett.test(residuals(modAout)~aAout$V1)

Mon problème réside dans le fait que lorsque je souhaite vérifier l'hétéroscédasticité des résidus du modèle linéaire via la commande suivante:

Code :

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

bartlett.test(residuals(modJuin)~aJuin$V1)

J'obtiens ce message d'erreur:

Code :

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

1
2
3
> bartlett.test(residuals(modJuin)~aJuin$V1)
Error in model.frame.default(formula = residuals(modJuin) ~ aJuin$V1) : 
  variable lengths differ (found for 'aJuin$V1')

Or, quand je verifie la longueur de "residuals(modJuin)" et "aJuin$V1":

Code :

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

1
2
3
4
> length(residuals(modJuin))
[1] 353
> length(aJuin$V1)
[1] 355

La longueur diffère en effet.

Mais via cette commande:

Code :

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

1
2
residuals(modJuin)
aJuin$V1

les 2 vecteurs contiennent bien 355 valeurs.

Pourriez-vous m'éclairer? Sur ce premier point?

En attendant, j'ai donc opté pour une évaluation visuelle de la normalité et de l'hétéroscédasticité via la fonction "plot" comme vous pouvez le voir dans le premier script que j'ai copié dans mon message. D'après moi, la normalité n'est pas respectée pour les valeurs extrêmes (voir un peu plus largement) mais l'hétéroscédasticité reste dans les clous. Le mois d'Aout est le cas le plus extrême que j'ai rencontré pour cette dernière hypothèse. Que feriez-vous à ma place? J'aurais tendance à valider ces 2 hypothèses et travailler en paramétrique. L'équivalent non-paramétrique serait le coefficient de Pearson (ou Spearman).

Merci pour votre aide!

[Yoan73]

Salut,

souvent les problème liés à des différences de longueur entre deux vecteurs proviennent des valeurs manquantes.
Généralement, par défaut, le logiciel enlève les NA avec la fonction na.omit qui soustrait des lignes du vecteur.
Il est possible de vérifier la commande par défaut avec cette ligne :

Code :

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

1
2
3
4
 
> options("na.action")
$na.action
[1] "na.omit"

Cependant, dans votre cas, cela pose un problème pour une possible comparaison. Pour conserver ces lignes avec la mention NA, il faut alors plutôt utiliser na.exclude. Pour vérifier rapidement la différence, regardez sur cet exemple tiré de : http://www.ats.ucla.edu/stat/r/faq/missing.htm

Code :

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

1
2
3
4
5
6
7
8
9
 
resid(model.omit)
##      4      5      6      7      8      9     10     11 
##  0.727  1.575 -0.799 -0.743 -1.553 -0.425  2.190 -0.971 
resid(model.exclude)
##      1      2      3      4      5      6      7      8      9     10 
##     NA     NA     NA  0.727  1.575 -0.799 -0.743 -1.553 -0.425  2.190 
##     11 
## -0.971

Le nombre de lignes n'est plus le même selon l'action employée.
Pour vous, cela devrais ce résumer à ces quelques changements :

Code :

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

1
2
3
 
modJuin=lm(aJuin$V2~aJuin$V1, na.action="na.exclude")
modJuin=aov(aJuin$V2~aJuin$V1, na.action="na.exclude")

Pour ensuite obtenir :

Code :

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

1
2
3
4
5
6
7
 
> bartlett.test(residuals(modJuin)~aJuin$V1)
 
	Bartlett test of homogeneity of variances
 
data:  residuals(modJuin) by aJuin$V1
Bartlett's K-squared = 83.1399, df = 61, p-value = 0.03131

Voila qui pourra peut être vous aider sur le 1er point.

Bon courage,

Yoan

[CupOfAppleTea]

Merci en effet ça résout ce premier problème!