A Distribución t de Student no es más que una distribución gaussiana con colas más pesadas. En otras palabras, podemos decir que la distribución gaussiana es un caso especial de la distribución t de Student. La distribución gaussiana está definida por la media (μ) y la desviación estándar (σ). La distribución t de Student, por otro lado, agrega un parámetro adicional, los grados de libertad (gl), que controla el “espesor” de la distribución. Este parámetro asigna mayor probabilidad a eventos más alejados de la media. Esta característica es particularmente útil para muestras de tamaño pequeño, como en biomedicina, donde la suposición de normalidad es cuestionable. Tenga en cuenta que a medida que aumentan los grados de libertad, la distribución t de Student se acerca a la distribución gaussiana. Podemos visualizar esto usando gráficos de densidad:
# Load necessary libraries
library(ggplot2)# Set seed for reproducibility
set.seed(123)
# Define the distributions
x <- seq(-4, 4, length.out = 200)
y_gaussian <- dnorm(x)
y_t3 <- dt(x, df = 3)
y_t10 <- dt(x, df = 10)
y_t30 <- dt(x, df = 30)
# Create a data frame for plotting
df <- data.frame(x, y_gaussian, y_t3, y_t10, y_t30)
# Plot the distributions
ggplot(df, aes(x)) +
geom_line(aes(y = y_gaussian, color = "Gaussian")) +
geom_line(aes(y = y_t3, color = "t, df=3")) +
geom_line(aes(y = y_t10, color = "t, df=10")) +
geom_line(aes(y = y_t30, color = "t, df=30")) +
labs(title = "Comparison of Gaussian and Student t-Distributions",
x = "Value",
y = "Density") +
scale_color_manual(values = c("Gaussian" = "blue", "t, df=3" = "red", "t, df=10" = "green", "t, df=30" = "purple")) +
theme_classic()
Nota en Figura 1 que la colina alrededor de la media se hace más pequeña a medida que los grados de libertad disminuyen como resultado de que la masa de probabilidad va hacia las colas, que son más gruesas. Esta propiedad es lo que le da a la distribución t de Student una sensibilidad reducida a los valores atípicos. Para más detalles sobre este asunto, puedes consultar este Blog.
Cargamos las bibliotecas requeridas:
library(ggplot2)
library(brms)
library(ggdist)
library(easystats)
library(dplyr)
library(tibble)
library(ghibli)
Entonces, saltemos las simulaciones de datos y pongámonos serios. Trabajaremos con datos reales que adquirí de ratones que realizaron la prueba del rotarod.
Primero, cargamos el conjunto de datos en nuestro entorno y configuramos los niveles de factor correspondientes. El conjunto de datos contiene identificaciones de los animales, una variable de tanteo (genotipo), un indicador para dos días diferentes en los que se realizó la prueba (día) y diferentes ensayos para el mismo día. Para este artículo, modelamos solo uno de los ensayos (Prueba3). Guardaremos las otras pruebas para un artículo futuro sobre la variación del modelado.
Como implica el manejo de datos, nuestra estrategia de modelado se basará en Genotipo y Día como predictores categóricos de la distribución de Trial3.
En la ciencia biomédica, los predictores categóricos o factores de agrupación son más comunes que los predictores continuos. A los científicos de este campo les gusta dividir sus muestras en grupos o condiciones y aplicar diferentes tratamientos.
data <- read.csv("Data/Rotarod.csv")
data$Day <- factor(data$Day, levels = c("1", "2"))
data$Genotype <- factor(data$Genotype, levels = c("WT", "KO"))
head(data)
Tengamos una vista inicial de los datos usando Parcelas de nubes de lluvia como se muestra Guilherme A. Franchi, PhD en este gran publicación de blog.
edv <- ggplot(data, aes(x = Day, y = Trial3, fill=Genotype)) +
scale_fill_ghibli_d("SpiritedMedium", direction = -1) +
geom_boxplot(width = 0.1,
outlier.color = "red") +
xlab('Day') +
ylab('Time (s)') +
ggtitle("Rorarod performance") +
theme_classic(base_size=18, base_family="serif")+
theme(text = element_text(size=18),
axis.text.x = element_text(angle=0, hjust=.1, vjust = 0.5, color = "black"),
axis.text.y = element_text(color = "black"),
plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5),
legend.position="bottom")+
scale_y_continuous(breaks = seq(0, 100, by=20),
limits=c(0,100)) +
# Line below adds dot plots from {ggdist} package
stat_dots(side = "left",
justification = 1.12,
binwidth = 1.9) +
# Line below adds half-violin from {ggdist} package
stat_halfeye(adjust = .5,
width = .6,
justification = -.2,
.width = 0,
point_colour = NA)
edv
Figura 2 se ve diferente del original por Guilherme A. Franchi, PhD porque estamos trazando dos factores en lugar de uno. Sin embargo, la naturaleza de la trama es la misma. Presta atención a los puntos rojos, estos son los que pueden considerarse observaciones extremas que inclinan las medidas de tendencia central (especialmente la media) hacia una dirección. También observamos que las varianzas son diferentes, por lo que modelar también sigma puede dar mejores estimaciones. Nuestra tarea ahora es modelar la salida usando el brms paquete.