UN Distribuzione t di Student non è altro che una distribuzione gaussiana con code più pesanti. In altre parole, possiamo dire che la distribuzione gaussiana è un caso speciale della distribuzione t di Student. La distribuzione gaussiana è definita dalla media (μ) e dalla deviazione standard (μ). La distribuzione t di Student, d'altro canto, aggiunge un ulteriore parametro, i gradi di libertà (df), che controlla lo “spessore” della distribuzione. Questo parametro assegna maggiore probabilità agli eventi più lontani dalla media. Questa caratteristica è particolarmente utile per campioni di piccole dimensioni, come in biomedicina, dove il presupposto della normalità è discutibile. Si noti che all'aumentare dei gradi di libertà, la distribuzione t di Student si avvicina alla distribuzione gaussiana. Possiamo visualizzarlo utilizzando i grafici di densità:
# 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 dentro Figura 1 che la collina attorno alla media diventa più piccola man mano che i gradi di libertà diminuiscono a causa della massa di probabilità che va verso le code, che sono più spesse. Questa proprietà è ciò che conferisce alla distribuzione t di Student una sensibilità ridotta ai valori anomali. Per maggiori dettagli su questo argomento, puoi controllare Questo blog.
Carichiamo le librerie richieste:
library(ggplot2)
library(brms)
library(ggdist)
library(easystats)
library(dplyr)
library(tibble)
library(ghibli)
Quindi, tralasciamo le simulazioni di dati e facciamo sul serio. Lavoreremo con dati reali che ho acquisito dai topi che eseguono il test rotarod.
Per prima cosa carichiamo il set di dati nel nostro ambiente e impostiamo i livelli di fattore corrispondenti. Il set di dati contiene gli ID degli animali, una variabile di brancolamento (Genotipo), un indicatore per due giorni diversi in cui è stato eseguito il test (giorno) e diverse prove per lo stesso giorno. Per questo articolo, modelliamo solo una delle prove (Trial3). Salveremo le altre prove per un futuro articolo sulla variazione della modellazione.
Come implica la gestione dei dati, la nostra strategia di modellazione sarà basata sul genotipo e sul giorno come predittori categorici della distribuzione di Trial3.
Nella scienza biomedica, i predittori categorici o i fattori di raggruppamento sono più comuni dei predittori continui. Agli scienziati in questo campo piace dividere i loro campioni in gruppi o condizioni e applicare trattamenti diversi.
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)
Diamo una prima visione dei dati utilizzando Trame di nuvole di pioggia come mostrato da Guilherme A. Franchi, PhD In Questo ottimo post sul 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 sembra diverso dall'originale di Guilherme A. Franchi, PhD perché stiamo tracciando due fattori invece di uno. Tuttavia, la natura della trama è la stessa. Presta attenzione ai punti rossi, sono quelli che possono essere considerati osservazioni estreme che inclinano le misure di tendenza centrale (soprattutto la media) verso una direzione. Osserviamo anche che le varianze sono diverse, quindi anche la modellazione sigma può fornire stime migliori. Il nostro compito ora è modellare l'output utilizzando il file brms pacchetto.
Fonte: towardsdatascience.com
