Comment comparer la moyenne de deux échantillons dont les données correspondent aux distributions exponentielles

10

J'ai deux échantillons de données, un échantillon de référence et un échantillon de traitement.

L'hypothèse est que l'échantillon de traitement a une moyenne plus élevée que l'échantillon de référence.

Les deux échantillons ont une forme exponentielle. Étant donné que les données sont assez volumineuses, je n'ai que la moyenne et le nombre d'éléments pour chaque échantillon au moment où je vais exécuter le test.

Comment puis-je tester cette hypothèse? Je suppose que c'est super facile, et je suis tombé sur plusieurs références à l'utilisation du F-Test, mais je ne sais pas comment les paramètres sont mappés.

Jonathan Dobbie
la source
2
Pourquoi n'avez-vous pas les données? Si les échantillons sont de très gros tests non paramétriques, cela devrait très bien fonctionner, mais il semble que vous essayez d'exécuter un test à partir des statistiques récapitulatives. Est-ce correct?
Mimshot
Les valeurs de référence et de traitement proviennent-elles du même ensemble de patients ou les deux groupes sont-ils indépendants?
Michael M
1
@Mimshot, les données sont en streaming, mais vous avez raison d'essayer d'exécuter un test à partir des statistiques récapitulatives. Cela fonctionne assez bien avec un test Z pour les données normales
Jonathan Dobbie
1
Dans ces circonstances, un test z approximatif est peut-être le meilleur que vous puissiez faire. Cependant, je me soucierais davantage de l'ampleur du véritable effet du traitement, et non de la signification statistique. Rappelez-vous qu'avec des échantillons suffisamment grands, tout petit effet vrai entraînera une petite valeur p.
Michael M
1
@janvier - bien que, si ses tailles d'échantillon sont suffisamment grandes, par le CLT, elles seront très proches de celles normalement distribuées. Dans l'hypothèse nulle, les variances seraient les mêmes (comme le sont les moyennes), donc, avec une taille d'échantillon suffisamment grande, un test t devrait fonctionner correctement; ce ne sera pas aussi bon que vous pouvez le faire avec toutes les données, mais ce serait quand même OK. , par exemple, serait plutôt bien. n1=n2=100
jbowman

Réponses:

14

Vous pouvez tester l'égalité des paramètres moyens par rapport à l'alternative selon laquelle les paramètres moyens sont inégaux avec un test de rapport de vraisemblance (test LR). (Cependant, si les paramètres moyens diffèrent et que la distribution est exponentielle, il s'agit d'un décalage d'échelle, pas d'un décalage d'emplacement.)

Pour un test unilatéral (mais seulement asymptotiquement dans le cas bilatéral), je crois que le test LR se révèle être équivalent au suivant (pour montrer que c'est en fait le même que le test LR pour le unilatéral) dans le cas où il faudrait montrer que la statistique LR était monotone en ):X¯/y¯

Disons que nous paramétrons la ème observation dans la première exponentielle comme ayant pdf 1 / μ x exp ( - x i / μ x ) et la j ème observation dans le deuxième échantillon comme ayant pdf 1 / μ y exp ( - y j / μ y ) (sur les domaines évidents pour les observations et les paramètres). (Pour être clair, nous travaillons ici sous forme de moyenne et non sous forme de taux; cela n'affectera pas le résultat des calculs.)je1/μxexp(xi/μx)j1/μyexp(yj/μy)

Puisque la distribution de est un cas particulier du gamma, Γ ( 1 , μ x ) , la distribution de la somme des X , S x est distribuée Γ ( n x , μ x ) ; de même que pour la somme des Y s, S y est Γ ( n y , μ y ) .XiΓ(1,μx)XSxΓ(nx,μx)YSyΓ(ny,μy)

En raison de la relation entre les distributions gamma et les distributions khi-deux, il s'avère que est distribué χ 2 2 n x . Le rapport de deux chi-carrés sur leurs degrés de liberté est F. D'où le rapport, μ y2/μXSXχ2nX2.μyμXSX/nXSy/nyF2nX,2ny

Dans l'hypothèse nulle d'égalité des moyennes, alors, , et selon l'alternative bilatérale, les valeurs pourraient avoir tendance à être plus petites ou plus grandes qu'une valeur de la distribution nulle , vous avez donc besoin d'un test bilatéral.X¯/y¯F2nX,2ny


Simulation pour vérifier que nous n'avons pas commis d'erreur simple dans l'algèbre:

Ici, j'ai simulé 1000 échantillons de taille 30 pour et 20 pour Y à partir d'une distribution exponentielle avec la même moyenne, et calculé la statistique du rapport des moyennes ci-dessus.XOui

Vous trouverez ci-dessous un histogramme de la distribution résultante ainsi qu'une courbe montrant la distribution nous avons calculée sous le nul:F

exemple de distribution simulée de la statistique de rapport sous la valeur nulle


Exemple, avec discussion du calcul des valeurs p bilatérales :

Pour illustrer le calcul, voici deux petits échantillons de distributions exponentielles. L'échantillon X a 14 observations d'une population avec une moyenne de 10, l'échantillon Y a 17 observations d'une population avec une moyenne de 15:

x: 12.173  3.148 33.873  0.160  3.054 11.579 13.491  7.048 48.836 
   16.478  3.323  3.520  7.113  5.358

y:  7.635  1.508 29.987 13.636  8.709 13.132 12.141  5.280 23.447 
   18.687 13.055 47.747  0.334  7.745 26.287 34.390  9.596

Les moyennes d'échantillon sont respectivement 12,082 et 16,077. Le rapport des moyennes est de 0,7515

La zone à gauche est simple, car elle se trouve dans la queue inférieure (calculée en R):

 > pf(r,28,34) 
 [1] 0.2210767

Nous avons besoin de la probabilité pour l'autre queue. Si la distribution était symétrique à l'inverse, il serait simple de le faire.

Une convention courante avec le rapport des variances F-test (qui est pareillement à deux queues) est simplement de doubler la valeur p à une queue (effectivement ce qui se passe comme ici ; c'est aussi ce qui semble être fait dans R, par exemple ); dans ce cas, il donne une valeur de p de 0,44.

α/2α

Glen_b -Reinstate Monica
la source
Je suppose que c'est juste moi qui est épais, mais d'où vient 0.7515?
Jonathan Dobbie
r = moyenne (x) / moyenne (y) = 0,7515 - c'est-à-dire "le rapport des moyennes"
Glen_b -Reinstate Monica
D'accord, génial. J'ai obtenu 0,67, mais cela est probablement dû à une erreur de saisie de données.
Jonathan Dobbie
1
J'ai fait la distinction entre les moyennes de la population et les moyennes de l'échantillon résultant plus clairement
Glen_b -Reinstate Monica
αα2
3

nXJournalnXXje+nyJournalnyyj-(nX+ny)JournalnX+nyXje+yj
nXJournal(nXny+1r)+nyJournal(nynX+r)+nXJournalnynX+ny+nyJournalnXnX+ny
r=X¯y¯r=1

Pour effectuer le test de rapport de vraisemblance approprié pour une alternative bilatérale, vous pouvez toujours utiliser la distribution F; il suffit de trouver l'autre valeur du rapport des moyennes d'échantillon rELRrobsPr(R>rELR)rELR=1,3272Pr(R>rELR)=0,21420,43520,4315

entrez la description de l'image ici

Mais doubler la valeur p unilatérale est peut-être le moyen le plus courant d'obtenir une valeur p bilatérale: cela revient à trouver la valeur du rapport des moyennes de l'échantillon rETPPr(R>rETP)Pr(R<robs)Pr(R>rETP)μX>μyμX<μyμX>μyμX<μy

entrez la description de l'image ici

Le code R suit:

x <- c(12.173, 3.148, 33.873, 0.160, 3.054, 11.579, 13.491, 7.048, 48.836,
       16.478, 3.323, 3.520, 7.113, 5.358)

y <- c(7.635, 1.508, 29.987, 13.636, 8.709, 13.132, 12.141, 5.280, 23.447, 
       18.687, 13.055, 47.747, 0.334,7.745, 26.287, 34.390, 9.596)

# observed ratio of sample means
r.obs <- mean(x)/mean(y)

# sample sizes
n.x <- length(x)
n.y <- length(y)

# define log likelihood ratio function
calc.llr <- function(r,n.x,n.y){
  n.x * log(n.x/n.y + 1/r) + n.y*log(n.y/n.x + r) + n.x*log(n.y/(n.x+n.y)) + n.y*log(n.x/(n.x+n.y))
}

# observed log likelihood ratio
calc.llr(r.obs,n.x, n.y) -> llr.obs

# p-value in lower tail
pf(r.obs,2*n.x,2*n.y) -> p.lo

# find the other ratio of sample means giving an LLR equal to that observed
uniroot(function(x) calc.llr(x,n.x,n.y)-llr.obs, lower=1.2, upper=1.4, tol=1e-6)$root -> r.hi

#p.value in upper tail
p.hi <- 1-pf(r.hi,2*n.x,2*n.y)

# overall p.value
p.value <- p.lo + p.hi

#approximate p.value
1-pchisq(2*llr.obs, 1)
Scortchi - Réintégrer Monica
la source