Esta página aún no se ha traducido para esta versión. Puede ver la versión más reciente de esta página en inglés.

Selección de un tamaño de muestra

Este ejemplo muestra cómo determinar el número de muestras u observaciones necesarias para llevar a cabo una prueba estadística. Ilustra los cálculos de tamaño de la muestra para un problema simple, a continuación, muestra cómo utilizar la función para calcular la potencia y el tamaño de la muestra para dos problemas más realistas.sampsizepwr Por último, ilustra el uso de las funciones estadísticas y machine learning Toolbox™ para calcular el tamaño de muestra necesario para una prueba que la función no admite.sampsizepwr

Probando una media normal con desviación estándar conocida, unilateral

Sólo para presentar algunos conceptos, vamos a considerar un ejemplo simple y poco realista donde queremos probar una media y conocemos la desviación estándar. Nuestros datos son continuos y se pueden modelar con la distribución normal. Necesitamos determinar un tamaño de muestra N para que podamos distinguir entre una media de 100 y una media de 110. Sabemos que la desviación estándar es 20.

Cuando realizamos una prueba estadística, generalmente probamos una contra una.null hypothesisalternative hypothesis Encontramos una estadística de prueba T, y miramos su distribución bajo la hipótesis nula. Si observamos un valor inusual, digamos uno que tiene menos de 5% de probabilidades de suceder si la hipótesis nula es verdadera, entonces rechazamos la hipótesis nula a favor de la alternativa. (La probabilidad del 5% se conoce como la de la prueba.)significance level Si el valor no es inusual, no rechazamos la hipótesis nula.

En este caso, la estadística de prueba T es la media de la muestra. Bajo la hipótesis nula tiene una media de 100 y una desviación estándar de 20/sqrt (N). Primero echemos un vistazo a un tamaño de muestra fijo, digamos N = 16. Rechazaremos la hipótesis nula si T está en la región sombreada, que es la cola superior de su distribución. Eso hace que esto sea una prueba unilateral, ya que no rechazaremos si T está en la cola inferior. El límite para esta región sombreada es 1,6 desviaciones estándar por encima de la media.

rng(0,'twister'); mu0 = 100; sig = 20; N = 16; alpha = 0.05; conf = 1-alpha; cutoff = norminv(conf, mu0, sig/sqrt(N)); x = [linspace(90,cutoff), linspace(cutoff,127)]; y = normpdf(x,mu0,sig/sqrt(N)); h1 = plot(x,y); xhi = [cutoff, x(x>=cutoff)]; yhi = [0, y(x>=cutoff)]; patch(xhi,yhi,'b'); title('Distribution of sample mean, N=16'); xlabel('Sample mean'); ylabel('Density'); text(96,.01,sprintf('Reject if mean>%.4g\nProb = 0.05',cutoff),'Color','b'); 

Así es como T se comporta bajo la hipótesis nula, pero ¿qué pasa con una alternativa? Nuestra distribución alternativa tiene una media de 110, representada por la curva roja.

mu1 = 110; y2 = normpdf(x,mu1,sig/sqrt(N)); h2 = line(x,y2,'Color','r'); yhi = [0, y2(x>=cutoff)]; patch(xhi,yhi,'r','FaceAlpha',0.25); P = 1 - normcdf(cutoff,mu1,sig/sqrt(N)); text(115,.06,sprintf('Reject if T>%.4g\nProb = %.2g',cutoff,P),'Color',[1 0 0]); legend([h1 h2],'Null hypothesis','Alternative hypothesis'); 

Existe una mayor probabilidad de rechazar la hipótesis nula si la alternativa es verdadera. Esto es justo lo que queremos. Es más fácil visualizar si miramos la función de distribución acumulativa (CDF) en lugar de la densidad (pdf). Podemos leer las probabilidades directamente desde este gráfico, en lugar de tener que calcular las áreas.

ynull = normcdf(x,mu0,sig/sqrt(N)); yalt = normcdf(x,mu1,sig/sqrt(N)); h12 = plot(x,ynull,'b-',x,yalt,'r-'); zval = norminv(conf); cutoff = mu0 + zval * sig / sqrt(N); line([90,cutoff,cutoff],[conf, conf, 0],'LineStyle',':'); msg = sprintf(' Cutoff = 100 + %.2g \\times 20 / \\surd{n}',zval); text(cutoff,.15,msg,'Color','b'); text(min(x),conf,sprintf('   %g%% test',100*alpha),'Color','b',...                          'verticalalignment','top') palt = normcdf(cutoff,mu1,sig/sqrt(N)); line([90,cutoff],[palt,palt],'Color','r','LineStyle',':'); text(91,palt+.02,sprintf(' Power is 1-%.2g = %.2g',palt,1-palt),'Color',[1 0 0]); legend(h12,'Null hypothesis','Alternative hypothesis'); 

Este gráfico muestra la probabilidad de obtener una estadística significativa (rechazando la hipótesis nula) para dos valores MU diferentes cuando N = 16.

El se define como la probabilidad de rechazar la hipótesis nula cuando la alternativa es verdadera.power function Depende del valor de la alternativa y del tamaño de la muestra. Vamos a graficar la potencia (que es una menos la CDF) como una función de N, la fijación de la alternativa en 110. Seleccionaremos N para lograr una potencia de 80%. El gráfico muestra que necesitamos alrededor de N = 25.

DesiredPower = 0.80; Nvec = 1:30; cutoff = mu0 + norminv(conf)*sig./sqrt(Nvec); power = 1 - normcdf(cutoff, mu1, sig./sqrt(Nvec)); plot(Nvec,power,'bo-',[0 30],[DesiredPower DesiredPower],'k:'); xlabel('N = sample size') ylabel('Power') title('Power function for the alternative: \mu = 110') 

En este ejemplo excesivamente simple hay una fórmula para calcular el valor necesario directamente para obtener una potencia de 80%:

mudiff = (mu1 - mu0) / sig; N80 = ceil(((norminv(1-DesiredPower)-norminv(conf)) / mudiff)^2) 
 N80 =      25  

Para verificar que esto funcione, hagamos una simulación Monte Carlo y generemos 400 muestras de tamaño 25, tanto bajo la hipótesis nula con la media 100, como bajo la hipótesis alternativa con la media 110. Si Probamos cada muestra para ver si su media es 100, debemos esperar que el 5% del primer grupo y el 80% del segundo grupo sean significativos.

nsamples = 400; samplenum = 1:nsamples; N = 25; h0 = zeros(1,nsamples); h1 = h0; for j = 1:nsamples     Z0 = normrnd(mu0,sig,N,1);     h0(j) = ztest(Z0,mu0,sig,alpha,'right');     Z1 = normrnd(mu1,sig,N,1);     h1(j) = ztest(Z1,mu0,sig,alpha,'right'); end p0 = cumsum(h0) ./ samplenum; p1 = cumsum(h1) ./ samplenum; plot(samplenum,p0,'b-',samplenum,p1,'r-') xlabel('Sample number') ylabel('Proportion significant') title('Verification of power computation') legend('Null hypothesis','Alternative hypothesis','Location','East') 

Probando una media normal con desviación estándar desconocida, doble cara

Ahora Supongamos que no conocemos la desviación estándar, y queremos realizar una prueba bilateral, es decir, una que rechace la hipótesis nula si la media de la muestra es demasiado alta o demasiado baja.

La estadística de prueba es una estadística t, que es la diferencia entre la media de la muestra y la media que se está probando, dividida por el error estándar de la media. Bajo la hipótesis nula, esto tiene la distribución de Student t con N-1 grados de libertad. Bajo la hipótesis alternativa, la distribución es una distribución no central t con un parámetro de no centralidad igual a la diferencia normalizada entre la media verdadera y la media que se está probando.

Para esta prueba bilateral tenemos que asignar el 5% de probabilidades de un error bajo la hipótesis nula por igual a ambas colas, y rechazar si la estadística de prueba es demasiado extrema en cualquier dirección. También tenemos que considerar ambas colas bajo cualquier alternativa.

N = 16; df = N-1; alpha = 0.05; conf = 1-alpha; cutoff1 = tinv(alpha/2,df); cutoff2 = tinv(1-alpha/2,df); x = [linspace(-5,cutoff1), linspace(cutoff1,cutoff2),linspace(cutoff2,5)]; y = tpdf(x,df); h1 = plot(x,y); xlo = [x(x<=cutoff1),cutoff1]; ylo = [y(x<=cutoff1),0]; xhi = [cutoff2,x(x>=cutoff2)]; yhi = [0, y(x>=cutoff2)]; patch(xlo,ylo,'b'); patch(xhi,yhi,'b'); title('Distribution of t statistic, N=16'); xlabel('t'); ylabel('Density'); text(2.5,.05,sprintf('Reject if t>%.4g\nProb = 0.025',cutoff2),'Color','b'); text(-4.5,.05,sprintf('Reject if t<%.4g\nProb = 0.025',cutoff1),'Color','b'); 

En lugar de examinar la función de potencia justo debajo de la hipótesis nula y un único valor alternativo de MU, podemos mirarla como una función de MU. La potencia aumenta a medida que Mu se aleja del valor de la hipótesis nula en cualquier dirección. Podemos utilizar la función para calcular la potencia.sampsizepwr Para un cálculo de potencia necesitamos especificar un valor para la desviación estándar, que sospechamos será aproximadamente 20. Aquí hay una imagen de la función de potencia para un tamaño de muestra N = 16.

N = 16; x = linspace(90,127); power = sampsizepwr('t',[100 20],x,[],N); plot(x,power); xlabel('True mean') ylabel('Power') title('Power function for N=16') 

Queremos una potencia de 80% cuando la media es 110. De acuerdo con este gráfico, nuestra potencia es inferior al 50% con un tamaño de muestra de N = 16. ¿Qué tamaño de muestra dará la potencia que queramos?

N = sampsizepwr('t',[100 20],110,0.8) 
 N =      34  

Necesitamos un tamaño de muestra de alrededor de 34. En comparación con el ejemplo anterior, tenemos que tomar nueve observaciones adicionales para permitir una prueba bilateral y para hacer para no conocer la verdadera desviación estándar.

Podemos hacer una trama de la función de potencia para varios valores de N.

Nvec = 2:40; power = sampsizepwr('t',[100 20],110,[],Nvec); plot(Nvec,power,'bo-',[0 40],[DesiredPower DesiredPower],'k:'); xlabel('N = sample size') ylabel('Power') title('Power function for the alternative: \mu = 110') 

Y podemos hacer una simulación similar a la anterior para verificar que tenemos el poder que necesitamos.

nsamples = 400; samplenum = 1:nsamples; N = 34; h0 = zeros(1,nsamples); h1 = h0; for j = 1:nsamples     Z0 = normrnd(mu0,sig,N,1);     h0(j) = ttest(Z0,mu0,alpha);     Z1 = normrnd(mu1,sig,N,1);     h1(j) = ttest(Z1,mu0,alpha); end p0 = cumsum(h0) ./ samplenum; p1 = cumsum(h1) ./ samplenum; plot(samplenum,p0,'b-',samplenum,p1,'r-') xlabel('Sample number') ylabel('Proportion significant') title('Verification of power computation') legend('Null hypothesis','Alternative hypothesis','Location','East') 

Supongamos que pudiéramos permitirnos tomar un tamaño de muestra de 50. Presumiblemente, nuestro poder para detectar el valor alternativo MU = 110 sería mayor que 80%. Si mantenemos la potencia al 80%, ¿qué alternativa podríamos detectar?

mu1 = sampsizepwr('t',[100 20],[],.8,50) 
 mu1 =    108.0837  

Probando una proporción

Ahora vamos a recurrir al problema de determinar el tamaño de la muestra necesario para distinguir entre dos proporciones. Supongamos que estamos probando una población en la que aproximadamente el 30% favorece a algún candidato, y queremos muestrear suficientes personas para que podamos distinguir este valor del 33%.

La idea aquí es la misma que antes. Aquí podemos usar el recuento de muestras como nuestro estadístico de prueba. Este recuento tiene una distribución binomial. Para cualquier tamaño de muestra N podemos calcular el límite para rechazar la hipótesis nula P = 0.30. Para N = 100, por ejemplo, rechazamos la hipótesis nula si el recuento de muestras es mayor que un valor de corte calculado de la siguiente manera:

N = 100; alpha = 0.05; p0 = 0.30; p1 = 0.33; cutoff = binoinv(1-alpha, N, p0) 
 cutoff =      38  

Una complicación con la distribución binomial viene porque es discreta. La probabilidad de exceder el valor límite no es exactamente del 5%:

1 - binocdf(cutoff, N, p0) 
 ans =      0.0340  

Una vez más, calculemos la potencia contra la alternativa P = 0.33 para una gama de tamaños de muestra. Usaremos una prueba unilateral (de cola derecha) porque solo estamos interesados en valores alternativos superiores al 30%.

Nvec = 50:50:2000; power = sampsizepwr('p',p0,p1,[],Nvec,'tail','right'); plot(Nvec,power,'bo-',[0 2000],[DesiredPower DesiredPower],'k:'); xlabel('N = sample size') ylabel('Power') title('Power function for the alternative: p = 0.33') 

Podemos utilizar la función para solicitar el tamaño de muestra requerido para una potencia de 80%.sampsizepwr

approxN = sampsizepwr('p',p0,p1,0.80,[],'tail','right') 
Warning: Values N>200 are approximate.  Plotting the power as a function of N may reveal lower N values that have the required power.   approxN =          1500  

Un mensaje de advertencia nos informa que la respuesta es simplemente aproximada. Si miramos la función de potencia para diferentes tamaños de muestra, podemos ver que la función es generalmente creciente, pero irregular porque la distribución binomial es discreta. Echemos un vistazo a la probabilidad de rechazar la hipótesis nula para p = 0.30 y p = 0.33 en el rango de tamaños de muestras de 1470 a 1480.

subplot(3,1,1); Nvec = 1470:1480; power = sampsizepwr('p',p0,p1,[],Nvec,'tail','right'); plot(Nvec,power,'ro-',[min(Nvec),max(Nvec)],[DesiredPower DesiredPower],'k:'); ylabel(sprintf('Prob[T>cutoff]\nif p=0.33')) h_gca = gca; h_gca.XTickLabel = ''; ylim([.78, .82]);  subplot(3,1,2); alf = sampsizepwr('p',p0,p0,[],Nvec,'tail','right'); plot(Nvec,alf,'bo-',[min(Nvec),max(Nvec)],[alpha alpha],'k:'); ylabel(sprintf('Prob[T>cutoff]\nif p=0.30')) h_gca = gca; h_gca.XTickLabel = ''; ylim([.04, .06]);  subplot(3,1,3); cutoff = binoinv(1-alpha, Nvec, p0); plot(Nvec,cutoff,'go-'); xlabel('N = sample size') ylabel('Cutoff') 

Esta gráfica revela que la curva de la función de potencia (parcela superior) no sólo es irregular, sino que también disminuye en algunos tamaños de muestra. Estos son los tamaños de muestra para los que es necesario aumentar el valor de corte (gráfico inferior) con el fin de mantener el nivel de significancia (trama intermedia) no mayor que 5%. Podemos encontrar un tamaño de muestra más pequeño dentro de este rango que da la potencia deseada de 80%:

min(Nvec(power>=0.80)) 
 ans =          1478  

Probando una correlación

En los ejemplos que hemos considerado hasta ahora, hemos sido capaces de averiguar el límite para un estadístico de prueba para lograr un cierto nivel de significancia, y para calcular la probabilidad de exceder ese corte bajo una hipótesis alternativa. Para nuestro ejemplo final, consideremos un problema donde eso no es tan fácil.

Imagine que podemos tomar muestras de dos variables X e y, y queremos saber qué tamaño de la muestra tendríamos que probar si no están correlacionadas frente a la alternativa de que su correlación es tan alta como 0,4. Aunque es posible trabajar la distribución de la correlación de muestra transformándola en una distribución t, vamos a usar un método que podemos usar incluso en problemas en los que no podemos trabajar la distribución de la estadística de prueba.

Para un tamaño de muestra determinado, podemos utilizar la simulación de Montecarlo para determinar un valor de corte aproximado para una prueba de la correlación. Hagamos una gran ejecución de simulación para que podamos obtener este valor con precisión. Empezaremos con un tamaño de muestra de 25.

nsamples = 10000; N = 25; alpha = 0.05; conf = 1-alpha; r = zeros(1,nsamples); for j = 1:nsamples     xy = normrnd(0,1,N,2);     r(j) = corr(xy(:,1),xy(:,2)); end cutoff = quantile(r,conf) 
 cutoff =      0.3372  

Entonces podemos generar muestras bajo la hipótesis alternativa, y estimar la potencia de la prueba.

nsamples = 1000; mu = [0; 0]; sig = [1 0.4; 0.4 1]; r = zeros(1,nsamples); for j = 1:nsamples     xy = mvnrnd(mu,sig,N);     r(j) = corr(xy(:,1),xy(:,2)); end [power,powerci] = binofit(sum(r>cutoff),nsamples) 
 power =      0.6470   powerci =      0.6165    0.6767  

Estimamos que el poder es del 65%, y tenemos un 95% de confianza en que el valor real está entre el 62% y el 68%. Para obtener una potencia de 80%, necesitamos un tamaño de muestra más grande. Podríamos intentar aumentar N a 50, estimar el valor de corte para este tamaño de muestra y repetir la simulación de potencia.

nsamples = 10000; N = 50; alpha = 0.05; conf = 1-alpha; r = zeros(1,nsamples); for j = 1:nsamples     xy = normrnd(0,1,N,2);     r(j) = corr(xy(:,1),xy(:,2)); end cutoff = quantile(r,conf)  nsamples = 1000; mu = [0; 0]; sig = [1 0.4; 0.4 1]; r = zeros(1,nsamples); for j = 1:nsamples     xy = mvnrnd(mu,sig,N);     r(j) = corr(xy(:,1),xy(:,2)); end [power,powerci] = binofit(sum(r>cutoff),nsamples) 
 cutoff =      0.2315   power =      0.8990   powerci =      0.8786    0.9170  

Este tamaño de muestra da una potencia mejor que nuestro objetivo de 80%. Podríamos seguir experimentando de esta manera, tratando de encontrar un tamaño de muestra menor que 50 que satisfaga nuestros requisitos.

Conclusión

Las funciones de probabilidad en el cuadro de herramientas de estadísticas y aprendizaje automático se pueden utilizar para determinar el tamaño de muestra necesario para lograr un nivel de potencia deseado en una prueba de hipótesis. En algunos problemas, el tamaño de la muestra se puede calcular directamente; en otros, es necesario buscar en un rango de tamaños de muestra hasta que se encuentra el valor correcto. Los generadores de números aleatorios pueden ayudar a verificar que se cumple la potencia deseada, y también pueden utilizarse para estudiar la potencia de una prueba específica en condiciones alternativas.