tests multiples - procédure de Storey (PFDR)

Bonjour,

Je galère méchamment à appliquer la procédure corrective de Storey sur cet exemple:


Ca fait maintenant un moment que je cherche alors je poste l'algorithme que je suis et mes calculs dans l'espoir de trouver enfin la méthode à suivre...

La page SAS sur cette procédure présente l'algorithme en question: https://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_multtest_sect014.htm (cf en bas de la page à Positive False Discovery)


Les paramètres sont:
- lambda = 0.5 => N(lambda) = 3 (en effet p1, p2, p3 sont inférieurs à 0.5)
- option finite activitée donc f = 1
- m = 5 hypothèses nulles

NB: N(.) = nombre de p-valeurs inférieures ou égales à ".".

Calcul de pi0_lambda = pi0 (pour faire plus simple) = (m - N(lambda) + f) / ((1 - lambda) * m):


Calcul de FDR_lambda (p) = FDR(p) (pour faire plus simple) = (pi0 * p) / (max(N(p),1) / m):

- Comme on a choisi p1 = "<0.0001" on en déduit que FDR(p1) = pFDR(p1) = "<0.0001"
- FDR(p2) = (1.2 * 0.0004) / (max(2,1) * 5) = 0.0012
- FDR(p3) = (1.2 * 0.0007) / (max(3,1) * 5) = 0.0014
- FDR(p4) = 1 car p4 = 0.6969 > lambda = 0.5
- FDR(p5) = 1 car p5 = 0.8129 > lambda = 0.5

Calcul de pFDR_lambda (p) = pFDR(p) (pour faire simple) = FDR(p) / (1 - (1 - p) ^m):

- pFDR(p1) = "<0.0001"
- pFDR(p2) = 0.0012 / (1 - (1 - 0.0004) ^5) = 0.6004802
- pFDR(p3) = 0.0014 / (1 - (1 - 0.0007) ^5) = 0.4005604
- pFDR(p4) = 1 / (1 - (1 - 0.6969) ^5) = 1.002565
- pFDR(p5) = 1 / (1 - (1 - 0.8129) ^5) = 1.000229

Enfin on peut calculer les p* (j'imagine que dans le lien sur SAS il y a une erreur, enfin l'oublie d'un p dans la formule finale): pi * = q_lambda (pi) = inf_{p >= pi} pFDR (p) quelque soit i = 1, ..., 5:
- p1* = "<0.0001"
- p2* = 0.4005604
- p3* = 0.4005604
- p4* = 1.000229 réajusté à 1
- p5* = 1.000229 réajusté à 1

Soit,


quand la solution est, on le rappel,


En gros, à 3 résultats prés sur 5, j'ai tout faux... J'ai consulté un max de doc, la formule change un coup sur deux, j'ai même consulté le papier de Storey qui parle de bootrstrap quand c'est évoqué nul part... bref je suis perdu.


Merci d'avance pour le coup de main car là je bloque totalement et comble du comble c'est la seule procédure pour laquelle je n'arrive pas à retrouver la solution sur cet exemple... même la procédure par combinatoire de Fisher j'ai fini par y mettre la main dessus et ce fut bien chiant...

Salut,

j'ai regardé un peu et tes calculs m'ont l'air correct (tu as juste une formule que tu as mal recopiée "(max(2,1) * 5" est en fait "(max(2,1) / 5)" mais tu avais bien calculé le bon)
As-tu un jeu de données à partager et le code SAS qui va avec pour que je rentre mieux dans le sujet quand j'aurai le temps?

Le bootstrap peut servir apparemment, pour trouver une valeur de lambda optimale.

Niaboc

Salut Niaboc et merci pour ton égard à mon topic

Bon j'ai finalement fini par trouver la solution. En fait pi0 = min(formule indiquée, 1). A partir de là ça marche!

Après ce qui me gène sur la procédure de Storey est principalement qu'il semble y avoir plusieurs méthodos. Si je compare les résultats de SAS et ceux de la fonction qvalue (je ne me souviens plus du package) de R j'ai des différences car N(lambda) n'est plus le nombre de p-valeurs inférieurs à lambda mais cette fois-ci supérieur. Dans la biblio que j'ai pu monter, un coup c'est le premier cas un coup l'autre tout en mentionnant à chaque fois procédure de Storey et non une énième version du pFDR.

De toute manière ça me casse les nouilles, c'est vraiment THE critique que je fais dans notre domaine d'application, on s'accorde pas toujours sur les formules d'usages et la plupart du temps les mecs qui publient sur le net balance des papelards simplement théorique et sans jamais l'algorithme d'application concret. Les trois exemples les plus concrets que j'ai pu trouver sont:
- l'analyse discriminante de Fisher, j'ai fini par mettre la main sur l'algorithme pour calculer les coefficients et écrire la formule décisionnelle mais bon sang que ce fut laborieux,
- les réseaux de neurones et les svm, mais là j'ai toujours pas mis la main sur les dit algorithmes

Il m'arrive parfois de me dire que les mecs qui publient sur le net sont encore plus paumés que moi la plupart du temps mais veulent se la jouer grand théoricien, je pari que 80% d'entre eux ne seraient pas capable de coder dans un logiciel de leur choix à partir de la théorie qu'ils écrivent.

C'était le coup de gueule du lapinou du fofo.