Test d'hypothèse: numérique vs classé

Sep 11 2020

Je pensais que le test d'hypothèse le plus puissant pour juger si un seul échantillon provient$N(0,1)$ ou de $N(1,1)$utilise la valeur moyenne comme statistique de test . Ainsi, j'ai calculé la taille de l'échantillon à l'aide d'une simulation de Monte Carlo et comparé deux méthodes: (a) en utilisant la valeur moyenne numérique comme statistique de test, et (b) en transformant d'abord les données en rangs puis en calculant la valeur moyenne. Ce que j'ai trouvé, c'est que les données classées donnent une taille d'échantillon plus petite. Comment est-ce possible?

Ce qui m'intéresse enfin, c'est de déterminer la taille de l'échantillon pour une future expérience. J'aime prendre$n$des échantillons provenant d'un processus industriel avec une distribution connue, analyser les échantillons et montrer que le processus est bien dans ses limites de spécification. Cela me permet de fabriquer le produit sans mesurer chaque pièce.

Le test d'hypothèse comporte de nombreuses difficultés subtiles. C'est pourquoi j'ai configuré le modèle de jouet suivant pour concentrer l'attention sur la question «numérique vs classée» décrite ci-dessus.

set.seed(2020)
transformToRank = FALSE
alpha           = 0.05
betaTarget      = 0.20

T = matrix(nrow = 1e4, ncol = 2) # alloc space for test statistic
for ( nSample in 3:100 ) {
    for ( i in 1:1e4 ){
        # combine data:
        random = c(rnorm(nSample,0,1), rnorm(nSample,1,1))
        if ( transformToRank ){
            # rank data (combined):
            random = rank(random)
        } 
        # Calc test statistic:
        T[i,] = c(mean(random[1:nSample]), mean(random[(nSample+1):(2*nSample)]) );
    }
    cutoff1 = quantile(T[, 1], alpha/2)
    cutoff2 = quantile(T[, 1], 1-alpha/2)
    beta    = sum(cutoff1 < T[, 2] & T[, 2] < cutoff2) / 1e4
    if ( beta <= betaTarget ){
        break 
    } else {
        nSample = nSample + 1
    }
}
print(beta)
print(nSample)

La casse numérique peut être vérifiée en utilisant le code suivant:

library(BSDA)
set.seed(2020)
pValue = replicate(1e4,z.test(rnorm(8,1,1), 
                    alternative="two.sided", mu=0, sigma.x=1)$p.value); 
power = mean(pValue <= 0.05)
beta  = 1 - power
print(beta)

Si j'utilise la valeur numérique, la taille de l'échantillon $ n = 8 $ satisfait la condition $ \ beta $ -risk. En revanche, si je transforme les nombres aléatoires en rangs, je n'ai besoin que de $ n = 4 $ . Ainsi, si je suis prêt à utiliser un ensemble de données généré aléatoirement pour analyser (classer) les données expérimentales, la puissance du test augmente considérablement. Cette logique s'applique également si j'échantillonne à partir d'une version à l'échelle de l'emplacement de la distribution $ t $ . Qu'est-ce que je rate?

Ma question clé est la suivante: pourquoi la version classée est-elle si puissante en termes de puissance? Je pensais que j'obtiendrais approximativement la même taille d'échantillon, car on sait que la transformation de rang fournit un pont vers les tests d'hypothèse non paramétriques, voir par exemple Ref1 ou Ref2 .

Réponses

2 BruceET Sep 12 2020 at 05:24

Voici des simulations comparant deux échantillons de taille 15 de $\mathsf{Norm}(0,1)$ et $\mathsf{Norm}(1,1),$respectivement. Ma simulation montre que le test t groupé a une meilleure puissance que le test de Wilcoxon à deux échantillons, qui est bien connu, et qu'aucun des tests n'a de puissance$0.8.$

set.seed(2020)
pv = replicate(10^4, t.test(rnorm(15,0,1),
                            rnorm(15,1,1), var.eq=T)$p.val) mean(pv <= 0.05) [1] 0.7525 set.seed(911) pv = replicate(10^4, wilcox.test(rnorm(15,0,1), rnorm(15,1,1))$p.val)
mean(pv <= 0.05)
[1] 0.7118

Il semble que j'ai mal compris ce que vous faites, que votre code de simulation est erroné, ou les deux. Il peut être utile d'avoir une explication plus claire de ce que vous faites avec les rangs et de voir la boucle interne de votre programme où vous calculez la puissance.

[Cela n'a aucun sens de prendre les moyennes des rangs pour les deux échantillons séparément: par exemple, si$n=15,$ alors les deux ensembles de rangs iraient de 1 à 15 et les deux ensembles de rangs seraient toujours totalisés à 120. Vous voudrez peut-être regarder ce que fait le test de la somme des rangs de Wilcoxon avec les rangs des deux échantillons.]

Voici des simulations avec des tailles d'échantillon $n=25$ et différence $0.5$en termes de population. Dans les deux cas, la puissance est proche de 80%.

set.seed(1066)
pv = replicate(10^4, t.test(rnorm(25,0,1),
                            rnorm(25,.5,1), var.eq=T)$p.val) mean(pv <= 0.05) [1] 0.3978 set.seed(1776) pv = replicate(10^4, wilcox.test(rnorm(25,0,1), rnorm(25,.5,1))$p.val)
mean(pv <= 0.05)
[1] 0.3867

Remarque: Pour les tests t groupés, voici une calculatrice en ligne de «puissance et taille d'échantillon», qui fonctionne pour des paramètres raisonnables.