El problema de hipotesis
multiple
Edgar Acuna
Universidad de Puerto Rico en
Mayaguez
EL problema de pruebas de
hipotesis multiples
• Prueba simultanea de m hipotesis nulas,
una por cada gene j (j=1,…m)
Hj: No hay relacion entre el nivel de expresion
del gene j y las distintas condiciones (gene no
expresados)
• Debido a que los experimentos con
microarreglos monitorean
simultaneamente niveles de expresion de
miles de genes hay que hacer ajuste a los
p-values de las hipotesis individuales.
El p-value en hipotesis simples
Resultado de la
prueba (p>): No
se rechaza Ho
Resultado de la
prueba (p<). Se
rechaza Ho
La Ho realmente Especificidad
es verdadera
Error Tipo I()
La Ho no es
verdadera
Sensitividad
Error tipo II
(probabilidad )
(probabilidad )
El p-value (nivel de significacion observado) controla la
tasa de error tipo I. La probabilidad de rechazar una
hipotesis por error no es mayor que  (=.05). De 100
pruebas solo 5 serian significativas
Posibilidades en Hipotesis multiples
# no
Decision
rechazadas
Verdad
#
rechazadas
total
# H ciertas
U
V (F +)
m0
# H no
ciertas
T(F-)
S
m1
total
m-R
R
m
m es fijo, mo y m1 son fijas pero
desconocidas, R es aleatorio y observado, V
es aleatorio y no observado.
El problema de multiplicidad
Cuando miles de hipotesis son probadas
simultaneamente se incrementa la posibilidad de
los falsos positivos
Por ejemplo, en un microarray con 10,000 genes
al 5% de significacion se puede esperar que 500
de ellos sean identificados como diferencialmente
expresados ( es decir sus p-values serian
menores que .05) simplemente por el azar.
Se ha perdido control del error tipo I y los p-values
individuales no se pueden usar directamente para
concluir que el gen es expresado. Necesitan ser
ajustados (corregidos)
Tasas de error tipo I (Falsos Positivos)
a controlar
• Tasa de error por familia
PFER = E(V): numero esperado de falsos positivos
• Tasa de error por comparacion
PCER = E(V)/m: proporcion esperada de falsos positivos
• Family-wise Error Rate
FWER = p(V ≥ 1) (probabilidad de al menos un falso
positivo)
• False Discovery Rate:
FDR = E(V/R)P(R>0) Proporcion esperada de falsos
positivos entre las pruebas que fueron significativas
Comparacion de tasas de error Tipo I
• En general, dado un procedimiento de
prueba de hipotesis multiple
y
PCER  FWER  PFER,
FDR  FWER,
Tipos de control del error
en hipotesis multiples
• Control debil. Se controla el error tipo I
asumiendo
que
todas
las
hipotesis
nulas
m
C
(Ho  j1Hj ) son ciertas. Es decir, mo=m. No
bria genes diferncialemente expresados.
• Control fuerte. Se controla el error tipo
I asumiendo cualquier combinacion de
hipotesis nulas y falsas. Algunos genes
serian diferencialmente expresados
p-values ajustados (p*)
• Objetivo: dado una tasa de error tipo I ,
hay que usar un procedimiento para
selecionar el conjunto de genes significantes
de tal manera que error tipo I sea 
• Si el interes es controlar el FWER, entonces
el p-value ajustado para la hipotesis Hj is:
pj* = inf {’: Hj is rechazado al FWER }
• Hipotesis Hj es rechazada al FWER  si pj* 

Notacion
• Para hipotesis Hj, j = 1, …, m
prueba estadistica calculada: tj
p-value no ajustado: pj
• Ordernamiento de tj (absoluto)
observado: {rj}
tal que |tr1|  |tr2|  …  |trm|
• Ordernamiento de pj observado: {rj}
tal que |pr1|  |pr2|  …  |prm|
Control de la FWER
• Metodo de Sidak
Rechazar Hj con pj1-(1-)1/m, p-value
ajustado pj* = 1-(1-pj)m
• Metodo de Bonferroni
Rechazar Hj con pj/m, p-value ajustado
pj* = min (mpj, 1)
Los metodos de Sidak y Bonferroni son faciles de calcular pero son
muy conservativos. Para un numero grande de genes, el nivel
individual de significancia por gene se vuelve bien pequeno bien
rapidamente.
Estos metodos son llamados de un solo paso porque usan el mismo
ajuste de multiplicidad para todas las hipotesis y no usan el
ordenamiento de los p-values observados
Los p-values ajustados pueden ser obtenidos usando la funcion
mt.rawp2adjp de la libreria multtest
Control de la FWER
• Metodo step-down de Holm (1979)
Los p-values ajustados se definen como
p*r1=mpr1
p*ri=max(p*ri-1,(m-i+1)pri) para 2im
con p*ri  1 tomados como 1.
O sea si j*=min{j:Prj >/(m-j+1)}, rechazar Hj para j=1,…j*-1.
Los p-values ajustados estan dados por:
prj* = maxk = 1…i{min ((m-k+1)prk, 1)}
Los p-values ajustados pueden ser obtenidos usando la funcion
mt.rawp2adjp de la libreria multtest
Control de la FWER
Por ejemplo supongamos que los 4 primeros p-values
ordenados son: 0.0006950, 0.0008659,
0.00149758, 0.00417016 de un total de m=1000 pvalues. Entonces los p-values ajustados
correspondientes seran: p*1=.6950,
p*2=max(.6950, 999*0.0008659)=0.86510,
p*3=max(.86510,998*0.00194758)=1.494585. Luego
se toma p*3 =1 pues el p-value no debe exceder a 1.
P*4=max(1,997*0.004170)==1 y todos los que siguen
tambien se toman como 1.
Control de la FWER
• Metodo step up de Hochberg (1988)
Los p-values ajustados se definen como
p*rm=prm
p*ri=min(p*ri+1,(m-i+1)pri) para i=m-1,…1
con p-values  1 tomados como 1.
O sea si j*=min{j:Prj/(m-j+1)}, rechazar Hj para j=1,…j*.
Los p-values ajustados vienen dados por
prj* = mink = j..m {min ((m-k+1)prk, 1) }
Los p-values ajustados pueden ser obtenidos usando la funcion
mt.rawp2adjp de la libreria multtest
Control de la FWER
P-value
Sidak
ordenado
Bonferroni Holm
pr1
1-(1-)1/m
/m
/m
/m
pr2
1-(1-)1/m
/m
/(m-1)
/(m-1)
….
………
…….
prj
1-(1-)1/m
/m
/(m-j+1)
/(m-j+1)
….
………..
……
prm
1-(1-)1/m
/m



Hochberg

Ejemplo; Aplicacion a datos de Golub
• 38 pacientes de Leucemia, 27 con
leucemia aguda linfoblastica (ALL) y
11 con leucemia aguda meloide (AML).
Inicialmente hay 6817 genes que se
reducen a 3051 aplicando cierto
criterios de exclusion para valores de
expresion
Ejemplo (cont)
Histogram of pvals
600
200
400
Frequency
300
200
0
100
0
Frequency
400
800
500
1000
histograma de valores t
-5
0
5
pruebash$t
10
0.0
0.4
0.8
pvals
1045 genes expresados con p-values <.05
Ejemplo (cont)
histograma p-values:Holm
2000
1000
1500
Frequency
1500
1000
0
500
500
0
Frequency
2000
2500
2500
histograma p-values:Sidak
0.0
0.2
0.4
0.6
a1$adjp[, 2]
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
a3$adjp[, 2]
En ambos casos se detectan 98 genes expresados
Top 10 genes
t
wilcox sidak
bonferroni
Holm
Hochberg
829
896
829
829
829
829
378
2124
378
378
378
378
2124
829
2124
2124
2124
2124
808
2670
808
808
808
808
2489
2939
2489
2489
2489
2489
394
394
394
394
394
394
2670
766
2670
2670
2670
2670
1009
808
1009
1009
1009
1009
1995
1834
1995
1995
1995
1995
Pruebas basadas en permutaciones
Estimar la distribucion conjunta de los estadisticos de prueba T1,..,Tm ,
donde m es el número de genes, mediante permutaciones de las columnas de
la matriz de expression genética
Labels originales de los grupos estadistico t
1.45
Labels permutados de los grupos
1.34
1.89
-1.33
-0.78
2.17
Algoritmo de Permutación para pvalues-no ajustados
• Calcular las pruebas tj para cada una de las hipótesis Hj.
• Hacer B permutaciones (B aprox 1000) de las columnas
de la matriz de expresión genética
• Para la b-ésima permutacion b=1,2,….B
A) Calcular las pruebas estadisticas t1,b,……tm,b
B) El p-value estimado por permutacion para la prueba de
hipotesis Hj está dado por la proporcion de |tib| ‘s que
son mayores que |ti|
1
p

I
(|
t
|

|t
|
)
j
j
,
b
j
b
B
*
donde I (.) representa la función indicadora
La libreria multtest tiene la funcion mt.sample.teststat que
calcula el test estadistico por permutaciones pero lo
hace vector por vector
Westfall & Young (1993) pvalues ajustados
• En cada paso hace pequenos ajustes
• Toma en cuenta la distribucion conjunta
(dependencia ) de la pruebas estadisticas
• Menos conservativo que los anteriores
metodos de ajustar p-values.
• Puede ser estimado por remuestreo
resampling pero tarda bastante
(especialmente la version minP)
Metodo maxT de Westfall & Young
• Ordenar los t-values observados
• Para la b-esima (b=1…B) permutacion de las columnas
del conjunto de datos calcular
a) los tj,b values para cada Hj (j=1,…m)
b) Calcular los maximos consecutivos de la prueba
estadistica.
um,b=|trm,b|
uj,b=max(uj+1,b,|trj,b|) para j=m-1,……1
Calcular los pvalues1
ajustados usando
*
p

I(|
u
|
|tr |)
r
j,b
j
j
B
b
Metodo maxT de Westfall & Young
Tambien es conocido como step-down maxT
Un formula equivalentemente para los p-values
ajustados es
prj* = maxk = 1…j { p(maxl{rk…rm} |Tl| ≥ |trk| H0C )}
Tl es el valor de la prueba estadistica
correspondiente a la l-esima hipotesis.
Ejemplo del metodo MaxT
gene |t|
Gene
|tb| |ub| I(ub>|t|)

P*=/B
1
0.1
tr5
1
1.3
1.3
1
935
.935
4
0.2
tr4
4
0.8
1.3
1
876
.876
5
2.8
tr3
5
3.0
3.0
1
138
.138
2
3.4
tr2
2
2.1
3.0
0
145
.145
3
7.1
tr1
3
1.8
3.0
0
48
.048
Genes ordenados
segun su valor t
Valores u para la b-esima
permutacion
P-Values ajustados
Explicacion de los valores de la tabla
De la primera tabla hay que guardar los
indices de los ordenamientos. Esto es,
r1=3, r2=2,r3=5,r4=4,r5=1. Luego, de la
segunda tabla ur5=|tr5,b|=1.3, Tambien
Ur4=max(ur5,|tr4|)=max(1.3,0.8)=1.3
Ur3=max(ur4,|tr3|)=max(1.3,3.0)=3.0
Ur2=max(ur3,|tr2|)=max(3.0,2.1)=3.0 y
Ur1=max(ur2,|tr1|)=max(3.0,1.8)=2.0
Metodo minP de Westfall & Young
• Tambien llamado step-down minP
En este caso los p-values ajustados vienen dados
por:
p*r1=p(minl{r1…rm} Pl  pr1 H0C )}
prj* = max(prj-1,p(minl{rj…rm} Pl  prj H0C )} para
j=2,…m
Pl es la variable aleatoria para el p-value de la lesima hipotesis. Por ejemplo, Pl ~U(0,1)
maxT vs. minP
• Los p-values ajustados por maxT y minP
son los mismos cuando las pruebas
estadisticas son identicamente
distribuidas (id)
• maxT es mas rapida computacionalmente
que minP
• maxT es mas poderosa en los casos en que
el numero de genes m es grande y numero
de arreglos n es pequeno.
Top 10 genes
maxT minP
2124
2124
829
829
896
896
766
766
2600
2600
2939
2939
1995
1995
2386
2386
717
717
• El criterio de controlar el FWER es
demasiado conservativo esto significa que
muchos genes que son diferencialmente
expresados podrian no ser detectados.
• Para el ejemplo de Golub solo se detectan
98 genes como expresados con los ajusten
de bonferoni y Holm. El maxT de Westfall
detecta 91, pero el minP no detecta
ninguno.
• El criterio FDR trata de resolver este
problema
Control de la FDR
• Benjamini & Hochberg (1995): step-up
Asumiendo que el FDR es de nivel , se rechaza Hj para
j=1,…j* , donde j*=max{j: prj<=(j/m)}.
El p-value ajustado esta dado por
prj* = mink = j…m { min ((m/k) prk, 1) }
• Benjamini & Yuketieli (2001): conservative step-up.
Se rechaza Hj para j=1,…j* , donde j*=max{j: prj<=j/
mj=1m[1/j]}. Con p-value ajustado dado por
prj* = mink = j…m { min (mj=1m[1/j]/k] prk, 1) }
Los p-values ajustados pueden ser obtenidos usando la
funcion mt.rawp2adjp de la libreria multtest
Top 10 genes
BH
BY
829
829
378
378
2124
2124
808
808
2489
2489
394
394
2670
2670
1009
1009
1995
1995
Bejamini y Hocberg detecta 681
genes expresados y Bejnjamini y
Yuketilei detecta 269.
Tambien se puede ajustar los pvalues obtenidos con maxT o minP
res1=mt.maxT(golub,golub.cl)
Rawp=res1$rawp[order(res1$index)]
# Permutation adjusted p-values for
simple multiple testing procedures
Res2=mt.rawp2adjp(rawp,”BH”)
Otras propuestas para
multiple testing
• ‘Significance Analysis of Microarrays (SAM)’
(2 versions)
– Tusher et al. (2001)
– Efron et al. (2001), basada en empirical Bayes
• SAM tambien estima el ‘FDR’, pero este es
definido como E(V|H0C)/R y no como E(V/R)
• La libreria siggenes de Bioconductor encuentra
los genes diferenciados por SAM
Analisis de significancia de
Microarrays (SAM)
• Tusher et al. (2001)
• Does not assume normal distribution
– Instead, p-values computed via values
computed via permutation
• Test statistic: similar to t Test
statistic: but with modified with
modified variance estimate
– Improved for small experiments.
Analisis de significancia de
Microarrays (SAM)
•
First, compute test statistic per gene for the observed data.
t
xT  xC
sp  so
donde s*p es un estimado de la desviacion estandar combinada y so es el
coeficiente de variacion minimo de la expresion de todos los genes.
Compute average of test statistics distribution, over all statistics
distribution, over all permutations
– This gives estimate of distribution, if treatment has no effect.
• Un gen es considerado significamente diferenciado si su distancia
con respecto a la media de su distribucin excede un threshold 
• SAM es facil de implementar
• El estimado usa todos los genes, luego si uno deellos se afecta por
el tratamiento tambie se afectara el estimado.
• No es confiable en experimentos pequenos.
*
Ejemplo de SAM
Row
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
d.value
829 8.165222
2124 7.964784
2600 6.102371
2664 5.975750
766 5.970848
2489 -5.726212
717 -5.704438
1995 -5.696514
2939 -5.576921
2663 5.547021
378 5.408458
1778 5.336856
1911 5.170084
1413 5.168875
808 5.139462
stdev
rawp q.value
0.2958251
0
0.1778697
0
0.1911219
0
0.3918749 0
0.1731333 0
0.2154975 0
0.2068956 0
0.1933259 0
0.1650727 0
0.4178283 0
0.3024200 0
0.2215924 0
0.1897508 0
0.2809704 0
0.1819399 0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
R.fold
7.2771792
3.3953035
2.6686992
4.7229540
2.4972299
0.3449532
0.3450401
0.3735687
0.4132925
5.4551523
4.4230278
2.7817426
2.3574207
3.2926668
2.4278709
SAM usando empirical Bayes
(EBAM)
• Aqui se estima la probabilidad
posterior p1(Z) =1-pofo(Z)/f(Z) de que
un gen con score Z sea expresado. La
razon fo(Z)/f(Z) es estimada usando
una regresion logistica basada en las
densidades relativas de los scores Zi
Ejemplo de SAM (empirical Bayes)
# El valor optimo del factor a0 is determinado, donde
# los posibles valores de a0 son 0 y los quantiles 0, 0.05 y 0.1
# de la desviacion estandard de los genes. Hacer rand=123
find.out=find.a0(golub,golub.cl,alpha=c(0,0.05,0.1),rand=123)
# Una vez que se establece el valor optimo de a0,se efectua
#un analisis por Empirical Bayes.
>ebam.out=ebam(find.out,gene.names=golub.gnames[,3])
>cat("\n el numero de genes diferenciados
es",length(ebam.out$row.sig.genes),"\n")
el numero de genes diferenciados es 714
>
cat("\n los top 10 genes son\n")
los top 10 genes son
>
ebam.out$row.sig.genes[1:10]
[1] 2489 394 1995 2939 717 1042 2702 523 1811 849
References
• Alizadeh et al. (2000) Distinct types of diffuse large
B-cell lymphoma identified by gene expression profiling.
Nature 403: 503-511
• Benjamini and Hochberg (1995) Controlling the false
discovery rate: a practical and powerful approach to
multiple testing. JRSSB 57: 289-200
• Benjamini and Yuketieli (2001) The control of false
discovery rate in multiple hypothesis testing under
dependency. Annals of Statistics
• Efron et al. (2000) Microarrays and their use in a
comparative experiment. Tech report, Stats, Stanford
• Golub et al. (1999) Molecular classification of cancer.
Science 286: 531-537
References
• Hochberg (1988) A sharper Bonferroni procedure for
multiple tests of significance. Biometrika 75: 800-802
• Holm (1979) A simple sequentially rejective multiple
testing procedure. Scand. J Statistics 6: 65-70
• Tusher et al. (2001) Significance analysis of
microarrays applied to transcriptional responses to
ionizing radiation. PNAS 98: 5116 -5121
• Westfall and Young (1993) Resampling-based multiple
testing: Examples and methods for p-value adjustment.
New York: Wiley
• Yuketieli and Benjamini (1999) Resampling based false
discovery rate controlling multiple test procedures for
correlated test statistics. J Stat Plan Inf 82: 171-196
Descargar

Document