PRINCIPIOS DE DISEÑO EXPERIMENTAL: COMPARACIONES MÚLTIPLES

ISADORE NABI

PRINCIPIOS DE DISEÑO EXPERIMENTAL: ANÁLISIS DE VARIANZA DE UNA VÍA

ISADORE NABI

ASPECTOS CONCEPTUALES GENERALES DEL DISEÑO EXPERIMENTAL POR BLOQUES

ISADORE NABI

Como se señala en (Dey, 2010, págs. 1-2), en determinadas situaciones experimentales, puede haber variaciones sistemáticas presentes entre las unidades experimentales[1]. Por ejemplo, en un experimento de campo, las unidades experimentales suelen ser parcelas de tierra. En un experimento de este tipo, puede haber un gradiente de fecundidad tal que las parcelas del mismo nivel de fecundidad sean más homogéneas que las que tienen distintos niveles de fecundidad. En experimentos con lechones como unidades experimentales, es muy plausible que los lechones pertenecientes a la misma camada estén genéticamente más cercanos entre sí (naciendo del mismo par de padres) que los que pertenecen a diferentes camadas. De manera similar, en experimentos con ganado, pueden estar involucradas diferentes razas (o diferentes edades) y se espera que los animales que pertenecen a la misma raza sean más parecidos que los que pertenecen a diferentes razas. En el contexto de los ensayos clínicos con pacientes que forman las unidades experimentales, el ensayo puede realizarse en diferentes centros (principalmente para obtener un número suficiente de observaciones) y los pacientes del mismo centro pueden ser más parecidos que los de diferentes centros debido a las diferencias en el tratamiento. prácticas y/o procedimientos de gestión seguidos en los diferentes centros. Los ejemplos anteriores, que son meramente ilustrativos y de ninguna manera exhaustivos, demuestran que en muchas situaciones existe una variación sistemática entre las unidades experimentales. En tales situaciones, el uso de un diseño completamente aleatorio no es apropiado. Más bien, se debe aprovechar la información a priori sobre esta variación sistemática mientras se diseña el experimento en el sentido de que esta información se debe utilizar durante el diseño para eliminar el efecto de dicha variabilidad. El impacto de este esfuerzo se verá reflejado en un error reducido, aumentando así la sensibilidad del experimento. Las consideraciones anteriores llevaron a la noción de control o bloqueo local. Los grupos de unidades experimentales relativamente homogéneas se denominan bloques. Cuando el bloqueo se realiza de acuerdo con un atributo, obtenemos un diseño de bloque. En un diseño de bloques, los tratamientos se aplican aleatoriamente a las unidades experimentales dentro de un bloque, y la asignación aleatoria de tratamientos a las unidades experimentales dentro de un bloque se realiza de forma independiente en cada bloque. El más simple entre los diseños de bloques es el diseño de bloques completos al azar.

Adicionalmente, (Batabyal, Sarkar, & Mandal, 2015, pág. 19) señalan que el experimento de gradiente de fertilidad (específicamente el del experimento por ellos analizado) se realizó antes del experimento del cultivo de prueba según la metodología inductiva propuesta por Ramamoorthy et al (1967), durante el verano de 2008-09, dividiendo el campo experimental en tres franjas rectangulares a lo largo del ancho. Los gradientes de fertilidad se crearon aplicando dosis graduadas de fertilizante N, P y K en las tiras como se muestra en la Tabla 1. El maíz forrajero se cultivó exhaustivamente para ayudar a que los fertilizantes se transformaran en el suelo por la planta y los microbios.

Figura 1

Fuente: (Batabyal, Sarkar, & Mandal, 2015, pág. 19).

La referencia anterior permite comprender conceptualmente el concepto gradiente al que se refirió Aloke Day en penúltima referencia realizada, así como también generalizar conceptualmente lo expuesto por este autor. Así, expresando de forma abstracta lo anterior, puede afirmarse que, en ciertas condiciones experimentales, pueden presentarse variaciones sistemáticas entre las unidades experimentales. En tales experimentos, existe variabilidad diferenciada en la distribución de los datos muestrales en las subregiones del espacio de muestra (en la teoría del diseño muestral estas subregiones son conocidas como bloques) a causa de un conjunto de factores subyacentes (por ello se considera la variabilidad de carácter sistemático) y esa variabilidad diferenciada por regiones se expresa matemáticamente como un gradiente, es decir, como una matriz en cuyo interior se contienen las derivadas parciales de primer orden de la función objetivo (la que explica la propagación diferenciada de la variabilidad) evaluadas en las subregiones pertinentes. Esta es la forma usual en que en el contexto de la teoría del diseño de experimentos y los ensayos clínicos se maneja el problema de volatilidad diferenciada de la varianza. El concepto bloque formaliza la noción de control local e implica cierta homogeneidad mínima entre los elementos de cada grupo, tiene como objetivo diseñar el experimento de tal forma que se elimine el efecto de esta variabilidad sistemática. Por supuesto, en otros contextos aplicados distintos de los ensayos clínicos puede desearse analizar el comportamiento del fenómeno estudiado considerando los efectos que la variabilidad diferenciada de la varianza y por ello existen modelos como el de heterocedasticidad condicional autorregresiva. Cuando los bloques han sido organizados alrededor de un atributo se está en presencia de un diseño por bloques y, en tal escenario, las variables explicativas (los tratamientos experimentales, para el caso aplicado) son consideradas como aleatorias (se aplican aleatoriamente sobre las unidades experimentales -pacientes humanos o de otra especie- al interior de un bloque). Esto se establece bajo el supuesto de que la aplicación de tales tratamientos a los elementos de cada bloque es linealmente independiente o, lo que es lo mismo, que al realizar la operación producto vectorial (exterior) entre los vectores que contienen las variables consideradas como estocásticas se genera un sistema de ecuaciones homogéneo de solución no nula.

Como señala (Dey, 2010, pág. 3), el más simple entre los diseños de bloques es el diseño de bloques completos al azar. En tal diseño, se requiere que cada bloque tenga tantas unidades experimentales como el número de tratamientos a evaluar, es decir, el tamaño del bloque sea igual al número de tratamientos. Sin embargo, no siempre es posible adoptar un diseño de bloques completos al azar en cada situación experimental. En primer lugar, si se supone que la variación intrabloque depende directamente del tamaño del bloque, entonces es preferible la adopción de un diseño con bloques de tamaños pequeños a uno que tenga tamaños de bloque grandes. Esto restringe el uso de diseños de bloques completos al azar en situaciones donde el número de tratamientos es grande. Por ejemplo, en los experimentos agronómicos, el experimentador generalmente elige un bloque de tamaño 10-12 y, si se acepta, no se puede adoptar un diseño de bloque completo al azar en situaciones en las que, digamos, se van a comparar 20 tratamientos. Además, en muchas situaciones experimentales, el tamaño del bloque está determinado por la naturaleza del experimento. Por ejemplo, con algunos experimentos en psicología, es bastante común considerar a los dos miembros de un par de gemelos como unidades experimentales de un bloque. En ese caso, claramente no se puede realizar diseño de bloques completos al azar si el número de tratamientos es mayor que dos (puesto que por definición en cada bloque existirá únicamente una observación -. De manera similar, es razonable tomar a los compañeros de camada (por ejemplo, ratones) como unidades de un bloque y el tamaño de la camada puede no ser adecuado para acomodar todos los tratamientos bajo prueba.

Los pocos ejemplos considerados anteriormente muestran claramente que, en muchas situaciones, no se puede adoptar un diseño de bloques completos al azar y, por lo tanto, es necesario buscar diseños en los que no todos los tratamientos aparezcan en cada bloque. Estos diseños se denominan diseños de bloques incompletos. Como señala (Dey, 2010, págs. 3-4), el tipo de diseño más importante del conjunto de diseños balanceados es el diseño de bloques balanceado e incompleto (BIB, por su nombre en inglés) y sobre ellos debe decirse que estos todavía se encuentran útiles en el diseño de experimentos en diversos campos y en los últimos años se han encontrado aplicaciones más nuevas de estos diseños, por ejemplo, en criptografía visual (véase, por ejemplo, Bose y Mukerjee (2006), Adhikary, Bose, Kumar y Roy (2007) y las referencias allí citadas).

Como señala (Wikipedia, 2021), el diseño por bloques es una estructura de incidencia[2] consistente en un conjunto de elementos expresados en familias denominadas bloques, escogidos tal que las frecuencias de los elementos satisfacen ciertas condiciones que permiten que la colección de bloques exhiba simetría (balance de bloques). Si no se dan más especificaciones, usualmente por “diseño de bloques” se hace referencia a un diseño de bloques balanceado e incompleto. Se dice que un diseño está balanceado hasta τ si todos los subconjuntos τ del conjunto original se presentan (como evento estadístico) en la misma cantidad de bloques λ. Cuando τ no está especificado, generalmente se puede suponer que es 2, lo que significa que cada par de elementos se encuentra en el mismo número de bloques y el diseño está “balanceado por pares”. Cualquier diseño equilibrado hasta τ también está equilibrado en todos los valores más bajos de τ (aunque con diferentes valores para λ). Por ejemplo, un diseño balanceado por pares (τ=2) es también regular (τ=1). Cuando falla el requisito de equilibrio, un diseño puede estar parcialmente equilibrado si los subconjuntos τ se pueden dividir en n-ésimas clases, cada una con su correspondiente (y diferente) valor para λ.

Así, señala última fuente referida que, para el caso de τ=2, estos diseños por bloques se conocen como PBIB(n), cuyas clases forman un esquema de asociación[3]. La teoría de los esquemas de asociación generaliza la teoría de los caracteres de representación lineal de grupos (y, por consiguiente, los esquemas de asociación generalizan la noción de grupos). Por lo general, se asume que los diseños están incompletos, lo que significa que ningún bloque contiene todos los elementos del conjunto, descartando así un diseño trivial (esta es otra forma en que se expresa la creencia de la estadística matemática clásica de que, si en un sistema de ecuaciones una ecuación es linealmente dependiente de otra u otras, entonces la variable que es descrita mediante tal ecuación no aporta información relevante). Los diseños por bloques pueden tener (o no) bloques repetidos. Cuando no tienen bloques repetidos, se denominan simples, en cuyo caso la familia de bloques es un conjunto en lugar de un multiconjunto. En estadística, el diseño de bloques se extiende a diseños de bloques no binarios los cuales pueden contener múltiples copias de un elemento de X, lo que implica que un diseño es regular sólo si es también binario. La matriz de incidencia de un diseño no binario (véase más adelante) enlista el número de veces que cada elemento de repite en cada bloque.

Adicionalmente, como señala (Dey, 2010, pág. 4), existen generalizaciones de los diseños BIB. Los diseños BIB son los únicos diseños en la clase de diseños de experimentos binarios[4], equirreplicados[5] y propios[6] que son balanceados (según se definió antes) tanto en varianza como en eficiencia; sin embargo, es posible encontrar otros diseños con equilibrio de varianza y eficiencia si uno expande la clase de diseños a diseños no binarios, no equirreplicados o no apropiados. Los métodos de construcción de diseños balanceados en varianza y eficiencia con replicaciones posiblemente desiguales y tamaños de bloques desiguales son el estado del arte más avanzado en el ámbito de los diseños balanceados.

Además, como se señala en la última fuente referida, existen los diseños de experimentos por bloques parcialmente balanceados, dentro de los cuales los más estudiados y aplicados empíricamente son los diseños de bloques parcialmente balanceados (PBIB, por su nombre en inglés). Los diseños de PBIB se introducen formalmente mediante la noción de esquema de asociación antes definida. Existen diseños con dos o más clases asociadas, así como también otros diseños parcialmente balanceados que no son necesariamente diseños PBIB. Estos incluyen los diseños conocidos como de celosía, cíclico, bloque enlazado, diseños en C y diseños α.

REFERENCIAS

Batabyal, K., Sarkar, D., & Mandal, B. (2015). Fertilizer-prescription equations for targeted yield in radish under integrated nutrient management system. Journal of Horticultural Sciences, X(1), 18-23. Obtenido de blob:resource://pdf.js/782dc541-51e4-4535-9551-8b7db5f35d1b

Dey, A. (2010). Incomplete Block Design. Tob Tuck Link, Singapore: World Scientific Publishing Co. Pte. Ltd.

Gupta, S. C., & Jones, B. (Agosto de 1983). Equireplicate Balanced Block Designs with Unequal Block Sizes. Biometrika, LXX(2), 433-440.

Shah, K. R., & Ashish, D. (Septiembre de 1992). Binary Designs Are Not Always the Bes. The Canadian Journal of Statistics / La Revue Canadienne de Statistique, XX(3), 347-351.

Weisstein, E. W. (19 de Septiembre de 2021). Monoid. Obtenido de MathWorld – A Wolfram Web Resource: https://mathworld.wolfram.com/Monoid.html

Wikipedia. (6 de Julio de 2021). Block design. Obtenido de Design of experiments: https://en.wikipedia.org/wiki/Block_design

Wikipedia. (27 de Septiembre de 2021). Incidence structure. Obtenido de Incidence geometry: https://en.wikipedia.org/wiki/Incidence_structure


[1] La variabilidad diferenciada de la varianza se encuentra en la literatura bajo el nombre de homogeneidad de varianza u homocedasticidad. Sin embargo, debe señalarse que no necesariamente son equivalentes metodológicamente, lo cual se explica por el hecho de que filosófica e históricamente no lo son. La heterocedasticidad tiene su génesis conceptual en el contexto de las estructuras de datos conocidas como series temporales, mientras que la homogeneidad de varianza tiene su génesis histórica en las estructuras de datos de sección cruzada. La posibilidad de su divergencia metodológica puede verificarse si alguna prueba para varianza no se puede hacer para sección cruzada y sí se puede hacer para series temporales (y/o para datos de panel, lo cual sería necesario investigar) o, por supuesto, a la inversa, que se pueda realizar para datos de sección cruzada y/o de panel, mientras que no sea posible para datos de series temporales (o bien, que sea posible para una estructura de datos, mientras que para la otra y su estructura combinada -datos de panel- no sea posible). Si no son diferentes, resulta técnicamente adecuado (con el fin de evitar ambigüedades) hablar de variabilidad diferenciada de la varianza, en lugar de hablar de “homogeneidad de varianza” o “heterocedasticidad” (puesto que los diferenciaría en un contexto en que son equivalentes). Así, aunque la posibilidad de su equivalencia o divergencia conceptual sea una cuestión fundamentalmente filosófica, su diferenciación se encuentra en las minucias metodológicas de las distintas técnicas estadísticas para medir y clasificar la variabilidad de la varianza. Por supuesto, sus diferencias filosóficas están basada en un hecho histórico-técnico concreto innegable: ambas estructuras de datos son distintas.

[2] Como señala (Wikipedia, 2021), una estructura de incidencia es un sistema abstracto consistente en dos tipos de objeto y una relación única entre ellos que se conoce como estructura de incidencia Se consideran una generalización del concepto de plano. Por su definición, son una estructura métrica vinculada a una estructura algebraica.

[3] Un esquema de asociación es un concepto algebraico que generaliza la noción de grupo. Un grupo es un monoide en el que además se cumple que sus elementos son invertibles. Un monoide es un conjunto cerrado (el equivalente matemático de autocontenido) bajo una operación asociativa binaria y con elemento identidad I que pertenece a S tal que para todo elemento a que pertenece a S se cumple que I*a = a*I = a y se diferencia de un grupo en el sentido de que no exige que sus elementos sean invertibles bajo alguna operación. Véase (Weisstein, 2021).

[4] Como se señala en (Shah & Ashish, 1992, pág. 347), un diseño en el que cada tratamiento aparece como máximo una vez en cualquier bloque en particular.

[5] Como se señala en (Gupta & Jones, 1983, pág. 433), un diseño por bloques equirreplicado es aquel en el que las variables independientes (en el contexto de la bioestadística y la psicometría usualmente son los tipos de tratamiento) se repiten en cada bloque la misma cantidad de veces.  

[6] Como se señala en (Wikipedia, 2021), un diseño por bloques es propio cuando todos los bloques tienen el mismo tamaño. También, como se señala en la fuente referida, se estudian también diseños por bloques que no son necesariamente uniformes; para τ=2 se conocen en la bibliografía bajo el nombre general de diseños equilibrados por pares, en donde cada par de elementos de X (cada par de elementos el conjunto de variables independientes) está contenido en exactamente en λ subconjuntos o bloques, en donde λ pertenece a los números naturales.

Advertisement

FUNDAMENTOS GENERALES DEL PROCESO DE ESTIMACIÓN Y PRUEBA DE HIPÓTESIS EN R STUDIO. PARTE II, CÓDIGO EN R STUDIO CON COMENTARIOS

ISADORE NABI

##ESTABLECER EL DIRECTORIO DE TRABAJO

setwd(“(…)”)

##LEER EL ARCHIVO DE DATOS. EN ESTE CASO, SUPÓNGASE QUE LOS DATOS SON DE UNA MUESTRA ALEATORIA DE 21 TIENDAS UBICADAS EN DIFERENTES PARTES DEL PAÍS Y A LAS CUALES SE LES REALIZÓ VARIOS ESTUDIOS. PARA ELLO SE MIDIERON ALGUNAS VARIABLES QUE SE PRESENTAN A CONTINUACIÓN

###- menor16= es un indicador de limpieza del lugar, a mayor número más limpio. 

###- ipc= es un indice de producto reparado con defecto, indica el % de producto que se pudo reparar y posteriormente comercializar.

###- ventas= la cantidad de productos vendidos en el último mes.

read.table(“estudios.txt”)

## CREAR EL ARCHIVO Y AGREGAR NOMBRES A LAS COLUMNAS

estudios = read.table(“estudios.txt”, col.names=c(“menor16″,”ipc”,”ventas”))

names(estudios)

nrow(estudios)

ncol(estudios)

dim(estudios)

## REVISAR LA ESTRUCTURA DEL ARCHIVO Y CALCULAR LA MEDIA, LA DESVIACIÓN ESTÁNDAR Y LOS CUANTILES PARA LAS VARIABLES DE ESTUDIO Y, ADICIONALMENTE, CONSTRÚYASE UN HISTOGRAMA DE FRECUENCIAS PARA LA VARIABLE “VENTAS”

str(estudios)

attach(estudios)

ventas

###Nota: la función “attach” sirve para adjuntar la base de datos a la ruta de búsqueda R. Esto significa que R busca en la base de datos al evaluar una variable, por lo que se puede acceder a los objetos de la base de datos simplemente dando sus nombres.

###Nota: Al poner el comando “attach”, la base de datos se adjunta a la dirección de búsqueda de R. Entonces ahora pueden llamarse las columnas de la base de datos por su nombre sin necesidad de hacer referencia a la base de datos ventas es una columna -i.e., una variable- de la tabla estudios). Así, al escribrlo, se imprime (i.e., se genera visualmente para la lectura ocular)

## CALCULAR LOS ESTADÍSTICOS POR VARIABLE Y EN CONJUNTO

mean(ventas)

sd(ventas)

var(ventas)

apply(estudios,2,mean)

apply(estudios,2,sd)

###Nota: la función “apply” sirve para aplicar otra función a las filas o columnas de una tabla de datos

###Nota: Si en “apply” se pone un “1” significa que aplicará la función indicada sobre las filas y si se pone un “2” sobre las columnas

## APLICAR LA FUNCIÓN “quantile”.

quantile(ventas) ###El cuantil de función genérica produce cuantiles de muestra correspondientes a las probabilidades dadas. La observación más pequeña corresponde a una probabilidad de 0 y la más grande a una probabilidad de 1.

apply(estudios,2,quantile)

###Nótese que para aplicar la función “apply” debe haberse primero “llamado” (i.e., escrito en una línea de código) antes la función que se aplicará (en este caso es la función “quantile”).

(qv = quantile(ventas,probs=c(0.025,0.975)))

###Aquí se está creando un vector de valores correspondientes a determinada probabilidad (las ventas, en este caso), que para este ejemplo son probabilidades de 0.025 y 0.975 de probabilidad, que expresan determinada proporción de la unidad de estudio que cumple con una determinada característica (que en este ejemplo esta proporción es el porcentaje de tiendas que tienen determinado nivel de ventas -donde la característica es el nivel de ventas-).

## GENERAR UN HISTOGRAMA DE FRECUENCIAS PARA LA COLUMNA “ventas”

hist(ventas)

abline(v=qv,col=2)

###Aquí se indica con “v” el conjunto de valores x para los cuales se graficará una línea. Como se remite a “qv” (que es un vector numérico de dos valores, 141 y 243) en el eje de las x, entonces graficará dos líneas color rojo (una en 141 y otra en 243).

###Aquí “col” es la sintaxis conocida como parte de los “parámetros gráficos” que sirve para especificar el color de las líneas

hist(ventas, breaks=7, col=”red”, xlab=”Ventas”, ylab=”Frecuencia”,

     main=”Gráfico

   Histograma de las ventas”)

detach(estudios)

###”breaks” es la indicación de cuántas particiones tendrá la gráfica (número de rectángulos, para este caso).

## GENERAR UNA DISTRIBUCIÓN N(35,4) CON NÚMEROS PSEUDOALEATORIOS PARA UN TAMAÑO DE MUESTRA n=1000

y = rnorm(1000,35,2)

hist(y)

qy = quantile(y,probs=c(0.025,0.975))

hist(y,freq=F)

abline(v=qy,col=2)

lines(density(y),col=2) #”lines” es una función genérica que toma coordenadas dadas de varias formas y une los puntos correspondientes con segmentos de línea.

## GENERAR UNA FUNCIÓN CON LAS VARIABLES n (CANTIDAD DE DATOS), m (MEDIA MUESTRAL) y  s (DESVIACIÓN ESTÁNDAR MUESTRAL) QUE ESTIME Y GRAFIQUE, ADEMÁS DE LOS CÁLCULOS DEL INCISO ANTERIOR, LA MEDIA.

plot.m = function(n,m,s) {

  y = rnorm(n,m,s)

  qy = quantile(y,probs=c(0.025,0.975))

  hist(y,freq=F)

  abline(v=qy,col=2)

  lines(density(y),col=2) ###Aquí se agrega una densidad teórica (una curva que dibuja una distribución de probabilidad -de masa o densidad- de referencia), la cual aparece en color rojo.

  mean(y)

}

## OBTENER UNA MUESTRA DE TAMAÑO n=10 DE N(100, 15^2)

plot.m(10000,100,15)

###Nótese que formalmente la distribución normal se caracteriza siempre por su media y varianza, aunque en la sintaxis “rnorm” de R se introduzca su media y la raíz de su varianza (la desviación estándar muestral)

##Generar mil repeticiones e ingresarlas en un vector. Compárense sus medias y desviaciones estándar.

n=10000; m=100;s=15

I = 1000 ###”I” son las iteraciones

medias = numeric(I)

for(i in 1:I)           {#”for” es un bucle (sintaxis usada usualmente para crear funciones personalizadas)

  sam=rnorm(n,m,s) ###Aquí se crea una variable llamada “sam” (de “sample”, i.e., muestra) que contiene una la distribución normal creada con números pseudoaleatorios.

  medias[i]=mean(sam)   } ###”sam” se almacena en la i-ésima posición la i-ésima media generada con “rnorm” que le corresponde dentro del vector numérico de iteraciones (el que contiene las medias de cada iteración) medias[i] (que contiene los elementos generados con la función “mean(sam)”).

###Un bucle es una interrupción repetida del flujo regular de un programa; pueden concebirse como órbitas (en el contexto de los sistemas dinámicos) computacionales. Un programa está diseñado para ejecutar cada línea ordenadamente (una a una) de forma secuencial 1,2,3,…,n. En la línea m el programa entiende que tiene que ejecutar todo lo que esté entre la línea n y la línea m y repetirlo, en orden secuencial, una cantidad x de veces. Entonces el flujo del programa sería, para el caso de un flujo regular  1,2,3,(4,5,…,m),(4,5,…,m),…*x,m+1,m+2,…,n.

## UTILIZAR LA VARIABLE “medias[i]” GENERADA EN EL INCISO ANTERIOR PARA DETERMINAR LA DESVIACIÓN ESTÁNDAR DE ESE CONJUNTO DE MEDIAS (ALMACENADO EN “medias[i]”) Y DETERMINAR SU EQUIVALENCIA CON EL ERROR ESTÁNDAR DE LA MEDIA (e.e.)

###Lo anterior evidentemente implica que se está construyendo sintéticamente (a través de bucles computacionales) lo que, por ejemplo, en un laboratorio botánico se registra a nivel de datos (como en el que Karl Pearson y Student hacían sus experimentos y los registraban estadísticamente) y luego se analiza en términos de los métodos de la estadística descriptiva e inferencial (puesto que a esos dominios pertenece el e.e.).

sd(medias)     ### desviación de la distribución de las medias

(ee = s/sqrt(n)  )### equivalencia teórica

## COMPARAR LA DISTRIBUCIÓN DE MEDIAS

m

mean(medias)

## GRAFICAR LA DISTRIBUCIÓN DE MEDIAS GENERADA EN EL INCISO ANTERIOR

hist(medias)

qm = quantile(medias,probs=c(0.025,0.975))

hist(medias,freq=F)

abline(v=qm,col=2)

lines(density(medias),col=2)

## GENERAR UN INTERVALO DE CONFIANZA CON UN NIVEL DE 0.95 PARA LA MEDIA DE LAS VARIABLES SUJETAS A ESTUDIO

attach(estudios)

### Percentil 0.975 de la distribución t-student para 95% de área bajo la curva

n = length(ventas) ###Cardinalidad o módulo del conjunto de datos

t = qt(0.975,n-1) ###valor t de la distribución t de student correspondiente a un nivel de probabilidad y n-1 gl

###Se denominan pruebas t porque todos los resultados de la prueba se basan en valores t. Los valores T son un ejemplo de lo que los estadísticos llaman estadísticas de prueba. Una estadística de prueba es un valor estandarizado que se calcula a partir de datos de muestra durante una prueba de hipótesis. El procedimiento que calcula la estadística de prueba compara sus datos con lo que se espera bajo la hipótesis nula (fuente: https://blog.minitab.com/en/adventures-in-statistics-2/understanding-t-tests-t-values-and-t-distributions).

###”qt” es la sintaxis que especifica un valor t determinado de la variable aleatoria de manera que la probabilidad de que esta variable sea menor o igual a este determinado valor t es igual a la probabilidad dada (que en la sintaxis de R se designa como p)

###Para más información véase https://marxianstatistics.com/2021/09/05/analisis-teorico-de-la-funcion-cuantil-en-r-studio/

###”n-1″ son los grados de libertad de la distribución t de student.

#### Error Estándar

ee = sd(ventas)/sqrt(n)

### Intervalo

mean(ventas)-t*ee

mean(ventas)+t*ee

mean(ventas)+c(-1,1)*t*ee ###c(-1,1) es un vector que se introduce artificialmente para poder construir el intervalo de confianza al 95% (u a otro nivel de confianza deseado) en una sola línea de código.

## ELABORAR UNA FUNCIÓN QUE PERMITA CONSTRUIR UN INTERVALO DE CONFIANZA AL P% DE NIVEL DE CONFIANZA PARA LA VARIABLE X

ic = function(x,p) {

  n = length(x)

  t = qt(p+((1-p)/2),n-1)

  ee = sd(x)/sqrt(n)

  mean(x)+c(-1,1)*t*ee

}

###Intervalo de 95% confianza para ventas

ic(ventas,0.95)

ic(ventas,0.99)

###El nivel de confianza hace que el intervalo de confianza sea más grande pues esto implica que los estadísticos de prueba (las versiones muestrales de los parámetros poblacionales) son más estadísticamente más robustos, por lo que su vecindario de aplicación es más amplio.

ic(ipc,0.95)

ic(menor16,0.95)

## REALIZAR LA PRUEBA DE HIPÓTESIS (PARA UNA MUESTRA) DENTRO DEL INTERVALO DE CONFIANZA GENERADO AL P% DE NIVEL DE CONFIANZA

t.test(ventas,mu=180) ###Por defecto, salvo que se cambie tal configuración, R realiza esta prueba a un nivel de confianza de 0.95.

### Realizando manualmente el cálculo anterior:

(t2=(mean(ventas)-180)/ee) ###Aquí se calcula el valor t por separado (puesto que la sintaxis “t.test” lo estima por defecto, como puede verificarse en la consola tras correr el código). Se denota con “t2” porque anteriormente se había definido en la línea de código 106 t = qt(0.975,n-1) para la construcción manual de los intervalos de confianza.

2*(1-pt(t2,20)) ###Aquí se calcula manualmente el valor p. Se multiplica por dos para tener la probabilidad acumulada total (considerando ambas colas) al valor t (t2, siendo más precisos) definido, pues esta es la definición de valor p. Esto se justifica por el hecho de la simetría geométrica de la distribución normal, la cual hace que la probabilidad acumulada (dentro de un intervalo de igual longitud) a un lado de la media sea igual a la acumulada (bajo la condición especificada antes) a la derecha de la media.

2*(pt(-t2,20)) ###Si el signo resultante de t fuese negativo. Además, 20 es debido a n-1 = 21-1 = 20.

###La sintaxis “pt” calcula el valor de la función de densidad acumulada (cdf) de la distribución t de Student dada una determinada variable aleatoria x y grados de libertad df (degrees of freedom, equivalente a gl en español), véase https://www.statology.org/working-with-the-student-t-distribution-in-r-dt-qt-pt-rt/

## CREAR UNA VARIABLE QUE PERMITA SEPARAR ESPACIALMENTE (AL INTERIOR DE LA GRÁFICA QUE LOS REPRESENTA) AQUELLOS ipc MENORES A UN VALOR h (h=117) DE AQUELLOS QUE SON IGUALES O MAYORES QUE h (h=117)

(ipc1 = 1*(ipc<17)+2*(ipc>=17))

ipc2=factor(ipc1,levels=c(1,2),labels=c(“uno”,”dos”))

plot(ipc2,ipc)

abline(h=17,col=2)

## GENERAR GRÁFICO DE DIAMENTE CON LOS INTERVALOS DE CONFIANZA AL 0.95 DE NdC CENTRADOS EN LAS MEDIAS DE CADA GRUPO CREADO ALREDEDOR DE 17 Y UN BOX-PLOT

library(gplots)

plotmeans(ventas~ipc2) ###Intervalos del 95% alrededor de la media (GRÁFICO DE DIMANTES)

boxplot(ventas~ipc2)

## REALIZAR LA PRUEBA DE HIPÓTESIS DE QUE LA MEDIA ES LA MISMA PARA LOS DOS GRUPOS GENERADOS ALREDEDOR DE h=17

(med = tapply(ventas,ipc1,mean))

(dev = tapply(ventas,ipc1,sd))

(var = tapply(ventas,ipc1,var))

(n   = table(ipc1))

dif=med[1]-med[2]

###La sintaxis “tapply” aplica una función a cada celda de una matriz irregular (una matriz es irregular si la cantidad de elementos de cada fila varía), es decir, a cada grupo (no vacío) de valores dados por una combinación única de los niveles de ciertos factores.

### PRUEBA DE HIPÓTESIS EN ESCENARIO 1: ASUMIENDO VARIANZAS IGUALES (SUPUESTO QUE EN ESCENARIOS REALES DEBERÁ VERIFICARSE CON ANTELACIÓN)

varpond= ((n[1]-1)*var[1] + (n[2]-1)*var[2])/(n[1]+n[2]-2) ###Aquí se usa una varianza muestral ponderada como medida más precisa (dado que el tamaño de los grupos difiere) de una varianza muestral común entre los dos grupos construidos alrededor de h=17

e.e=sqrt((varpond/n[1])+(varpond/n[2]))

dif/e.e

t.test(ventas~ipc1,var.equal=T)

t.test(ventas~ipc1)  #Por defecto la sintaxis “t.test” considera las varianzas iguales, por lo que en un escenario de diferentes varianzas deberá ajustarse esto como se muestra a continuación.

### PRUEBA DE HIPÓTESIS EN ESCENARIO 2: ASUMIENDO VARIANZAS DESIGUALES (AL IGUAL QUE ANTES, ESTO DEBE VERIFICARSE)

e.e2=sqrt((var[1]/n[1])+(var[2]/n[2]))

dif/e.e2

a=((var[1]/n[1]) + (var[2]/n[2]))^2

b=(((var[1]/n[1])^2)/(n[1]-1)) +(((var[2]/n[2])^2)/(n[2]-1))

(glmod=a/b)

t.test(ventas~ipc1,var.equal=F)

###Para aceptar o rechazar la hipótesis nula el intervalo debe contener al cero (porque la Ho afirma que la verdadera diferencia en las medias -i.e., su significancia estadística- es nula).

###Conceptualmente hablando, una diferencia estadísticamente significativa expresa una variación significativa en el patrón geométrico que describe al conjunto de datos. Véase https://marxianstatistics.com/2021/08/27/modelos-lineales-generalizados/. Lo que define si una determinada variación es significativa o no está condicionado por el contexto en que se realiza la investigación y la naturaleza misma del fenómeno estudiado.

## REALIZAR PRUEBA F PARA COMPARAR LA VARIANZA DE LOS GRUPOS Y LA PROBABILIDAD ASOCIADA

(razon.2 = var[1]/var[2]) ###Ratio de varianzas (asumiendo que las varianzas poblacionales son equivalentes a la unidad, en otro caso su estimación sería matemáticamente diferente; véase https://sphweb.bumc.bu.edu/otlt/mph-modules/bs/bs704_power/bs704_power_print.html y https://stattrek.com/online-calculator/f-distribution.aspx).

pf(razon.2,n[1]-1,n[2]-1) ###Al igual que “pt” (para el caso de la t de Student que compara medias de dos grupos o muestras), “pf” en el contexto de la prueba F (que compara la varianza de dos grupos o muestras) calcula la probabilidad acumulada que existe hasta determinado valor.

###La forma general mínima (más sintética) de la sintaxis “pf” es “pf(x, df1, df2)”, en donde “x” es el vector numérico (en este caso, de un elemento), df1 son los gl del numerador y df2 son los grados de libertad del denominador de la distribución F (cuya forma matemática puede verificarse en la documentación de R; véase https://stat.ethz.ch/R-manual/R-devel/library/stats/html/Fdist.html).

(2*pf(razon.2,n[1]-1,n[2]-1)) ###Aquí se calcula el valor p manualmente.

###Realizando de forma automatizada el procedimiento anterior:

var.test(ventas~ipc1)

detach(estudios)

###Para aceptar o rechazar la hipótesis nula el intervalo debe contener al 1 porque la Ho afirma que la varianza de ambas muestras es igual (lo que implica que su cociente o razón debe ser 1), lo que equivale a afirmar que la diferencia real entre desviaciones (la significancia estadística de esta diferencia) es nula.

## EN EL ESCENARIO DEL ANÁLISIS DE MUESTRAS PAREADAS, ANALIZAR LOS DATOS SOBRE EL EFECTO DE DOS DROGAS EN LAS HORAS DE SUEÑO DE UN GRUPO DE PACIENTES (CONTENIDOS EN EL ARCHIVO “sleep” DE R)

attach(sleep) ###”sleep” es un archivo de datos nativo de R, por ello puede “llamarse” sin especificaciones de algún tipo.

plot(extra ~ group)

plotmeans(extra ~ group,connect=F)  ###Intervalos del 95% alrededor de la media. El primer insumo (entrada) de la aplicación “plotmeans” es cualquier expresión simbólica que especifique la variable dependiente o de respuesta (continuo) y la variable independiente o de agrupación (factor). En el contexto de una función lineal, como la función “lm()” que es empleada por “plotmeans” para graficar (véase la documentación de R sobre “plotmeans”), sirve para separar la variable dependiente de la o las variables independientes, las cuales en este caso de aplicación son los factores o variables de agrupación (puesto que se está en el contexto de casos clínicos y, en este contexto, las variables independientes son las variables que sirven de criterio para determinar la forma de agrupación interna del conjunto de datos; este conjunto de datos contiene las observaciones relativas al efecto de dos drogas diferentes sobre las horas de sueño del conjunto de pacientes-).

A = sleep[sleep$group == 1,] ###El símbolo “$” sirve para acceder a una variable (columna) de la matriz de datos, en este caso la número 1 (por ello el “1”).

B = sleep[sleep$group == 2,]

plot(1:10,A$extra,type=”l”,col=”red”,ylim=c(-2,7),main=”Gráfico 1

Horas de sueño entre pacientes con el tratamiento A y B”,ylab=”Horas”,xlab=”Numero de paciente”,cex.main=0.8)

lines(B$extra,col=”blue”)

legend(1,6,legend=c(“A”,”B”),col=c(“red”,”blue”),lwd=1,box.col=”black”,cex=1)

t.test(A$extra,B$extra)

t.test(A$extra,B$extra,paired=T)

t.test(A$extra-B$extra,mu=0)

###Una variable de agrupación (también llamada variable de codificación, variable de grupo o simplemente variable) clasifica las observaciones dentro de los archivos de datos en categorías o grupos. Le dice al sistema informático (sea cual fuere) cómo el usuario ha clasificado los datos en grupos. Las variables de agrupación pueden ser categóricas, binarias o numéricas.

###Cuando se desea realizar un comando dentro del texto (en un contexto de formato Rmd) se utiliza así,por ejemplo se podría decir que la media del sueño extra es `r mean(sleep$extra)` y la cantidad de datos son `r length(sleep$extra)`

## ESTIMACIÓN DE LA POTENCIA DE UNA PRUEBA DE HIPÓTESIS (PROBABILIDAD BETA DE COMETER ERROR TIPO II)

library(pwr) ###”pwr” es una base de datos nativa de R

delta=3 ###Nivel de Resolución de la prueba. Para un valor beta (probabilidad de cometer error tipo II) establecido el nivel de resolución es la distancia mínima que se desea que la prueba sea capaz de detectar, es decir, que si existe una distancia entre los promedios tal que la prueba muy probablemente rechace la hipótesis nula Ho. Para el cálculo manual de la probabilidad beta véase el complemento de este documento (FUNDAMENTOS GENERALES DEL PROCESO DE ESTIMACIÓN Y PRUEBA DE HIPÓTESIS EN R STUDIO. PARTE I, TEORÍA ESTADÍSTICA)

s=10.2 ###Desviación estándar muestral

(d=delta/s) #Tamano del efecto.

pwr.t.test(n=NULL,d=d,power =0.9,type=”one.sample”)

## ESTIMAR CON EL VALOR ÓPTIMO PARA EL NIVEL DE RESOLUCIÓN, PARTIENDO DE n=40 Y MANTENIENDO LA POTENCIA DE 0.9

(potencia=pwr.t.test(n=40,d=NULL,power =0.9,type=”one.sample”))

potencia$d*s  #Delta

## GRAFICAR LAS DIFERENTES COMBINACIONES DE TAMAÑO DE MUESTRA Y NIVEL DE RESOLUCIÓN PARA UNA POTENCIA DE LA PRUEBA FIJA

s=10.2

deltas=seq(2,6,length=30)

n=numeric(30)

for(i in 1:30) {

  (d[i]=deltas[i]/s)

  w=pwr.t.test(n=NULL,d=d[i],power =0.9,type=”one.sample”)

  n[i]=w$n

}

plot(deltas,n,type=”l”)

## SUPÓNGASE QUE SE QUIERE PROBAR SI DOS GRUPOS PRESENTAN DIFERENCIAS ESTADÍSTICAMENTE SIGNIFICATIVAS EN LOS NIVELES PROMEDIO DE AMILASA, PARA LO CUAL SE CONSIDERA IMPORTANTE DETECTAR DIFERENCIAS DE 15 UNIDADES/ML O MÁS ENTRE LOS PROMEDIOS

s2p=290.9  ###Varianza ponderada de los dos grupos

(sp=sqrt(s2p)) ###Desviación estándar ponderada de los dos grupos

delta=15

(d=delta/sp)

pwr.t.test(n=NULL,d=d,power =0.9,type=”two.sample”)