Test statistici in BioRender Graphing: metodi, controlli delle assunzioni e pacchetti R

Di seguito troverai i pacchetti specifici, le funzioni e i metodi utilizzati per ogni test. Per trasparenza, tutti i calcoli sono eseguiti con R (versione 4.5.1)

Indice

1. Elenco dei pacchetti R
2. Controlli delle assunzioni statistiche
3. Test t per un campione (e test non parametrici)
4. Test t per due campioni (e test non parametrici)
5. ANOVA a una via (e test non parametrici) e test di confronti multipli
6. ANOVA a due vie e test di confronti multipli
7. Correlazione
8. Regressione logistica
9. Curve dose-risposta
10. Analisi di sopravvivenza

Tutti i pacchetti R

• R versione 4.5.1
• afex versione 1.3-1
• dplyr versione 1.2.1
• effsize versione 0.8.1
• emmeans versione 1.10.1
• httr2 versione 1.0.0
• jsonlite versione 1.8.4
• lme4 versione 1.1-37
• logger versione 0.2.2
• multcomp versione 1.4-20
• nlraa versione 1.9.3
• OptimModel versione 2.0-3
• PMCMRplus versione 1.9.6
• readr versione 2.1.5
• rjson versione 0.2.21
• rstatix versione 0.7.1
• survival versione 3.8-3
• tictoc versione 1.1
• tidyr versione 1.3.0
• uuid versione 1.2-0

Controlli delle assunzioni statistiche

La nostra piattaforma esegue automaticamente test per verificare le assunzioni essenziali per ogni tipo di analisi. Suggeriamo opzioni di analisi pre-selezionate basate su questi risultati. Passando con il mouse sulle etichette "suggerito" nell'app vengono mostrati dettagli su queste raccomandazioni.

Test di Shapiro-Wilk per la normalità

I test parametrici t-test e ANOVA assumono che il dataset segua una distribuzione normale (curva a campana standard).

Eseguiamo in background il test di Shapiro-Wilk per verificare se il dataset viola questa assunzione.

Se il test di Shapiro-Wilk (p > 0.05) non può rifiutare l'ipotesi di normalità, raccomandiamo test parametrici (es. t-test).

Se il test di Shapiro-Wilk (p ≤ 0.05) indica che i dati non sono distribuiti normalmente, trasformiamo internamente i dati con il logaritmo e testiamo nuovamente. Se il test successivo suggerisce che i dati sono normalmente distribuiti nello spazio logaritmico, raccomandiamo test parametrici sui dati log-trasformati (es. t-test lognormale).

Se entrambi i test di Shapiro-Wilk rifiutano l'ipotesi di normalità, suggeriamo alternative non parametriche (es. test di Mann-Whitney U).

Questo controllo utilizza la funzione shapiro_test dal pacchetto rstatix.

Test di Levene per l'uguaglianza delle varianze

I test parametrici t-test e ANOVA assumono che la varianza delle osservazioni tra i diversi gruppi sia approssimativamente uguale.

Eseguiamo il test di Levene per verificare se il dataset viola questa assunzione.

Se le varianze differiscono, raccomandiamo la correzione di Welch per accomodare varianze disuguali nei t-test o ANOVA.

Questo controllo utilizza la funzione levene_test dal pacchetto rstatix.

Utilizziamo l'impostazione predefinita che applica la mediana come centro di ogni gruppo per un test più robusto.

Test t per un campione (e test non parametrici)

Un test t per un campione determina se la media di un singolo campione differisce significativamente da un valore noto o ipotizzato. Assume che i dati provengano da una popolazione distribuita normalmente.

Test t per un campione

Utilizziamo un'implementazione personalizzata base di R che calcola direttamente la statistica t come (media_campione − μ) / (DS / √n) usando pt e qt dal pacchetto stats per valori p e valori critici rispettivamente. I gradi di libertà sono n − 1. Gli intervalli di confidenza sono calcolati sulla differenza (media_campione − μ): per un test a due code, diff ± qt(1 − α/2, df) × SE; per test a una coda, l'IC è [diff − qt(conf_level, df) × SE, Inf] oppure [−Inf, diff + qt(conf_level, df) × SE]. La dimensione dell'effetto Cohen's d è riportata come (media_campione − μ) / DS.

Test t per rapporto a un campione (lognormale)

Quando i dati si prevede seguano una distribuzione lognormale, il test t per rapporto verifica se la media geometrica differisce da un valore positivo ipotizzato. Il test trasforma in logaritmo sia i dati che il valore ipotizzato, esegue il test t personalizzato identico sulla scala logaritmica e riporta i risultati trasformati indietro. La media geometrica, la DS geometrica e il rapporto rispetto al valore teorico sono restituiti sulla scala originale, insieme agli intervalli di confidenza sia per il rapporto sia per la media geometrica.

Test di Wilcoxon per un campione

L'alternativa non parametrica per confronti su un campione. Utilizziamo wilcox.test dal pacchetto stats con exact = NULL, lasciando che base R selezioni il metodo esatto o asintotico. Il flag p_value_exact è impostato a TRUE quando il numero di differenze non nulle n_nz < 50 (strettamente) e non ci sono differenze pari, rispecchiando il comportamento previsto dalla logica interna di R. In caso contrario viene applicata un'approssimazione normale con correzione di continuità (correct = TRUE).

La pseudo-mediana di Hodges–Lehmann è riportata come stima di posizione: è la mediana di tutte le medie di Walsh a coppie (x_i + x_j) / 2 per i ≤ j (inclusi i casi i = j) calcolate sul vettore originale dei valori — questo è lo stimatore H-L della mediana della popolazione. Viene inoltre restituita la statistica V (somma dei ranghi firmati positivi) e le somme separate dei ranghi positivi e negativi.

Test t per due campioni (e test non parametrici)

I test per due campioni confrontano un risultato continuo tra due gruppi. BioRender offre test t parametrici (assumendo normalità) e alternative non parametriche che non fanno assunzioni sulla distribuzione, insieme a varianti lognormali per dati in cui i rapporti, piuttosto che le differenze, sono l'effetto rilevante.

Tutti i test parametrici per due campioni usano t_test dal pacchetto rstatix, e tutti i test non parametrici usano wilcox_test dallo stesso pacchetto. Le dimensioni dell'effetto per i test parametrici sono calcolate con cohen.d dal pacchetto effsize. Nota che t_test di rstatix calcola la differenza come livello1 − livello2 (ordine alfabetico); BioRender nega la statistica t e gli IC per riportare livello2 − livello1 (in accordo con la convenzione di GraphPad Prism), e i valori p a una coda sono ricalcolati dalla statistica negata usando pt per assicurare l'uso della coda corretta.

Test t per due campioni indipendenti (Student)

t_test(..., var.equal = TRUE). Assume varianze uguali; usa una deviazione standard pooled tra entrambi i gruppi. L'errore standard pooled è pooled_SD × √(1/n₁ + 1/n₂). Cohen's d è calcolato con cohen.d(y, x, paired = FALSE).

Test t per due campioni indipendenti con correzione di Welch (aka test t di Welch)

t_test(..., var.equal = FALSE). Non assume varianze uguali; usa l'approssimazione di Welch per i gradi di libertà (Satterthwaite). L'errore standard è √(s₁²/n₁ + s₂²/n₂). Cohen's d è calcolato con cohen.d(y, x, paired = paired) da effsize, dove paired è l'argomento della funzione (default FALSE).

Test di Mann-Whitney U (test di rango di Wilcoxon)

wilcox_test(..., paired = FALSE). Confronto non parametrico di due campioni indipendenti. Le statistiche U₁ e U₂ sono calcolate dalle somme dei ranghi: U₁ = R₁ − n₁(n₁+1)/2 e U₂ = R₂ − n₂(n₂+1)/2, dove R₁ e R₂ sono le somme dei ranghi per ciascun gruppo. La statistica U riportata è min(U₁, U₂).

Test t accoppiato

t_test(..., paired = TRUE). I soggetti sono abbinati tra le due condizioni; il test è eseguito sulle differenze all'interno delle coppie. Cohen's d per disegni accoppiati usa la formula within-subjects: |mean(x − y)| / √(SD₁² + SD₂² − 2r·SD₁·SD₂), dove r è la correlazione di Pearson tra le due condizioni. Vengono inoltre restituiti il coefficiente di correlazione di Pearson e il suo valore p.

Test di Wilcoxon per ranghi firmati accoppiati

wilcox_test(..., paired = TRUE). Equivalente non parametrico del test t accoppiato. Sono restituiti la somma dei ranghi firmati positivi, la somma dei ranghi firmati negativi e la somma complessiva dei ranghi firmati W. Il flag p_value_exact è impostato a TRUE quando il numero di coppie è inferiore a 50 e non ci sono differenze pari, riflettendo il comportamento previsto dalla logica interna di rstatix. Vengono inoltre restituiti il coefficiente di correlazione di rango di Spearman e il suo valore p tra le due condizioni.

Varianti lognormali (rapporto)

Quando i dati seguono una distribuzione lognormale, l'interesse è nel rapporto delle medie geometriche piuttosto che nella differenza delle medie aritmetiche. Per i test non accoppiati e di Welch, tutti i valori sono trasformati in logaritmo naturale e viene eseguito il corrispondente t_test sulla scala logaritmica. Le medie geometriche sono trasformate indietro come exp(mean(log(values))). Seguendo la convenzione livello2 − livello1, geometric_mean1 corrisponde al secondo gruppo in ordine alfabetico e geometric_mean2 al primo; il rapporto riportato è geometric_mean1 / geometric_mean2. Gli IC sul rapporto si ottengono invertendo e scambiando i limiti degli IC sulla scala logaritmica: conf_low_ratio = exp(−conf.high_log), conf_high_ratio = exp(−conf.low_log). Per il test lognormale accoppiato, Cohen's d è calcolato anche sulla scala logaritmica; sono inoltre restituiti media, DS e SEM delle differenze logaritmiche e la correlazione di Pearson sui dati logaritmici.

ANOVA a una via (e test non parametrici) e test di confronti multipli

L'ANOVA a una via verifica se le medie di tre o più gruppi differiscono significativamente. BioRender supporta misure indipendenti, misure ripetute e varianti a effetti misti, oltre ad alternative non parametriche e ANOVA lognormale.

ANOVA classica a una via

Implementata tramite aov_car dal pacchetto afex. I contrasti sono impostati a somma zero (contr.sum) per tutta la durata dell'analisi e ripristinati al termine. Tipo di somma dei quadrati: la SS di tipo II è usata solo per disegni bilanciati indipendenti a una via (tutti i gruppi hanno n uguale e nessun fattore a misure ripetute). La SS di tipo III è usata in tutti gli altri casi: disegni sbilanciati, misure ripetute e tutti i disegni multifattoriali indipendentemente dal bilanciamento. La dimensione dell'effetto è riportata come eta-quadrato generalizzato (ges) dalla tabella ANOVA di afex.

ANOVA a una via con correzione di Welch (aka ANOVA di Welch a una via)

welch_anova_test da rstatix. Un test F che non assume varianze uguali tra i gruppi, usando un aggiustamento di Welch–Satterthwaite ai gradi di libertà del denominatore. Riporta la statistica W e i gradi di libertà numeratore/denominatore.

ANOVA lognormale

I valori sono trasformati in log₁₀ e viene eseguita la classic_anova sui dati trasformati. Le statistiche F, i valori p e i gradi di libertà sono invarianti alla base del logaritmo; SS e MS sono riportati sulla scala log₁₀ per corrispondere a GraphPad Prism.

Test di Kruskal-Wallis

kruskal_test da rstatix. Analogo non parametrico dell'ANOVA a una via. Riporta la statistica H e il valore p.

Test post-hoc di confronti multipli

I test post-hoc sono eseguiti dopo un'ANOVA significativa per identificare quali coppie di gruppi differiscono.

• Test di confronti multipli di Tukey usa le medie marginali stimate da emmeans e contrasti a coppie aggiustati con Tukey tramite pairs(..., adjust = "tukey"). Sono restituiti intervalli di confidenza e statistiche t per ogni confronto a coppie.
• Test di confronti multipli di Bonferroni usa lo stesso framework emmeans con pairs(..., adjust = "bonferroni"). Il valore critico t per l'IC è qt(0.05 / (2 × k), df), dove k è il numero di confronti.
• Test di Dunnett (un controllo vs tutti gli altri) usa glht da multcomp con mcp(categories = "Dunnet"). Il gruppo di controllo è reimpostato come riferimento prima dell'adattamento; gli intervalli di confidenza sono estratti tramite confint(glht_result).
• Test di Games–Howell (post-hoc per varianze disuguali) usa games_howell_test da rstatix. La funzione rstatix applica internamente un fattore √0.5 all'errore standard pooled; BioRender lo rimuove dividendo per √0.5 (equivalente a moltiplicare per √2) per corrispondere alla formulazione convenzionale.
• Test di Dunnett T3 (per varianze disuguali) usa dunnettT3Test da PMCMRplus. I gradi di libertà per coppia sono calcolati con la formula di Welch–Satterthwaite applicata a ogni singola coppia. I valori critici sono derivati dalla distribuzione t multivariata tramite qmvt.
• Test di Dunn (non parametrico, segue Kruskal–Wallis) usa dunn_test da rstatix con aggiustamento p-value Bonferroni per tutti i confronti a coppie. Quando è specificato un gruppo di controllo, sono restituiti solo i confronti che coinvolgono il controllo, e la correzione di Bonferroni usa la dimensione della famiglia di Dunnett (k − 1 confronti) invece di k(k−1)/2.

ANOVA a due vie e test di confronti multipli

L'ANOVA a due vie testa gli effetti di due fattori categorici e della loro interazione su un risultato continuo. BioRender usa la stessa funzione classic_anova dell'ANOVA a una via, estesa per includere una seconda colonna fattore nel dataframe.

• Somma dei quadrati: Il tipo III è sempre usato per l'ANOVA a due vie (bilanciata o no), che valuta ogni effetto dopo aver considerato tutti gli altri effetti. I contrasti sono impostati a contr.sum (somma zero) come richiesto per la SS di tipo III.
• Dimensione dell'effetto è riportata come eta-quadrato generalizzato (ges) per ogni effetto principale e per il termine di interazione.
• Misure ripetute: L'ANOVA a due vie con misure ripetute è supportata tramite la stessa infrastruttura aov_car / afex. Viene calcolato il test di sfericità di Mauchly per i fattori within-subject. Può essere applicata la correzione di Greenhouse–Geisser agli effetti che violano la sfericità.
• Fallback a modelli a effetti misti: Quando i dati a misure ripetute presentano osservazioni mancanti, si usa un modello a effetti misti (tramite lme4::lmer), sempre o solo in presenza di dati mancanti, a seconda dell'opzione missing_value_handling.

I test post-hoc per l'ANOVA a due vie usano le stesse funzioni dell'ANOVA a una via (Tukey, Bonferroni, Dunnett, Dunnett T3, Games–Howell) ma operano sulle medie delle celle di interazione tramite emmeans. I prefissi di livello fattore che emmeans aggiunge ai livelli fattore numerici sono rimossi dalle etichette dei contrasti prima di restituire i risultati, e i nomi delle colonne dei fattori within-subject preceduti da afex (es. Xfactor1) sono corretti ai nomi originali.

Correlazione

L'analisi di correlazione quantifica l'associazione lineare (Pearson) o monotona (Spearman) tra due variabili continue e testa se tale associazione differisce da zero.

Usiamo cor.test dal pacchetto stats. I valori mancanti sono gestiti tramite cancellazione pairwise. Sono richieste almeno quattro osservazioni complete.

Correlazione di Pearson

cor.test(x, y, method = "pearson") restituisce il coefficiente di correlazione di Pearson r, una statistica t con n − 2 gradi di libertà, un valore p a due code e intervalli di confidenza esatti al 95% (tramite trasformazione z di Fisher, applicata internamente da R).

Correlazione di Spearman

cor.test(x, y, method = "spearman") restituisce il coefficiente di correlazione di rango di Spearman ρ (rho). La distribuzione esatta per permutazione è usata quando non ci sono valori pari e n ≤ 50 (exact = TRUE); altrimenti si usa l'approssimazione normale asintotica. Poiché cor.test di R non restituisce intervalli di confidenza per la correlazione di Spearman, gli IC al 95% sono approssimati tramite la trasformazione z di Fisher: z = 0.5 × ln((1 + ρ) / (1 − ρ)), SE = 1 / √(n − 3), e i limiti degli IC sono trasformati indietro alla scala di correlazione.

Visualizzazione

Indipendentemente dalla selezione di correlazione di Pearson o Spearman, la linea di previsione e la banda di confidenza al 95% visualizzate sul grafico a dispersione derivano da una regressione lineare semplice (lm) adattata ai dati, su 100 punti equispaziati nell'intervallo osservato di x, tramite predict(model, interval = "confidence", level = 0.95).

Regressione lineare semplice

La regressione lineare semplice modella la relazione tra un singolo predittore continuo e un risultato continuo come una linea retta.

Usiamo lm dal pacchetto stats: lm(dependent ~ independent). È disponibile un'opzione per vincolare l'intercetta a zero tramite lm(dependent ~ 0 + independent).

I coefficienti di pendenza e intercetta, i loro errori standard, le statistiche t e i valori p sono presi direttamente da summary(model)$coefficients. L'intervallo di confidenza al 95% su ogni coefficiente è stima ± qt(0.975, df.residual(model)) × SE. L'intercetta sull'asse x è calcolata come −β₀ / β₁.

L'errore standard della stima (Sy.x) è √(SS_residuals / (n − 2)) per il modello non vincolato, o √(SS_residuals / (n − 1)) quando l'intercetta è forzata a zero, corrispondente alla convenzione di GraphPad Prism.

La significatività complessiva del modello è valutata con un test F: pf(F, df1, df2, lower.tail = FALSE) usando la statistica F da summary(model)$fstatistic.

Una curva adattata con banda di confidenza al 95% è restituita su 250 punti equispaziati che coprono l'intervallo osservato del predittore, tramite predict(model, interval = "confidence").

Regressione logistica

La regressione logistica modella il log-odds di un risultato binario come funzione lineare di uno o più predittori. BioRender offre tre varianti: semplice (un predittore), multipla (più predittori) e regressione logistica penalizzata di Firth per campioni piccoli o squilibrati.

Tutti i modelli di regressione logistica usano glm dal pacchetto stats con family = binomial(link = "logit") e maxit = 100 (il default di 25 può restituire silenziosamente un iterato non convergente su dati con quasi-separazione). Prima dell'adattamento, il servizio verifica la quasi-separazione completa (probabilità adattate numericamente 0 o 1) e valuta se è raccomandata la penalizzazione di Firth basata sulla rilevazione di separazione, basso numero di eventi per variabile (EPV) o tassi di eventi rari.

Regressione logistica semplice

Adatta un singolo predittore continuo o binario: glm(dependent ~ independent, family = binomial(link = "logit")).

Gli intervalli di confidenza per i coefficienti sono calcolati sulla scala del log-odds usando la verosimiglianza profilo tramite confint, con fallback automatico agli IC di Wald se il profiling fallisce. Gli odds ratio e i loro IC al 95% si ottengono esponenziando le stime dei coefficienti e i limiti di confidenza.

Il valore di X alla probabilità prevista del 50% (dove il log-odds è zero, cioè X = −β₀/β₁) è riportato con un IC al 95% derivato dal metodo delta.

L'adattamento del modello è valutato tramite un test di verosimiglianza rispetto al modello nullo (solo intercetta). Le statistiche di bontà di adattamento riportate includono R² di Tjur, R² di Cox–Snell, l'AUC della curva ROC (calcolata tramite un algoritmo vettoriale basato sul rango) e AICc. Sono restituite probabilità previste e una banda di confidenza sull'intervallo osservato del predittore (250 punti griglia) per la visualizzazione; la banda è calcolata con IC di Wald sulla scala del logit e trasformata indietro tramite plogis.

Regressione logistica penalizzata di Firth

Quando si rileva separazione perfetta o quasi completa, o quando EPV è basso (meno di ~10 eventi per predittore), le stime di massima verosimiglianza standard sono infinite o fortemente distorte. In questi casi, BioRender adatta il modello di verosimiglianza penalizzata di Firth (Firth 1993; Heinze & Schemper 2002), che aggiunge una penalizzazione prior di Jeffreys alla log-verosimiglianza per produrre stime finite e meno distorte.

Il modello di Firth è adattato usando un'implementazione IRLS personalizzata. La convergenza è valutata congiuntamente su tre criteri: max|Δβ| < 1e-8 (variazione parametro L∞), max|U*(β)| < 1e-5 (valore assoluto massimo del vettore score modificato, cioè norma score L∞), e |ΔL*| < 1e-5 (variazione della log-verosimiglianza penalizzata).

Gli intervalli di confidenza sono calcolati tramite verosimiglianza profilo penalizzata usando una funzione firth_profile_ci() personalizzata. Per ogni parametro, gli altri parametri sono riottimizzati a ogni valore candidato β tramite IRLS vincolata. Un'espansione dell'intervallo di ricerca con step di raddoppio (partendo da max(SE×3, 1.5), raddoppiando fino a un limite massimo di 50 sulla scala log-odds) individua un intervallo di cambio di segno, quindi uniroot() (metodo di Brent) trova la radice. Se non si trova un intervallo, il limite è riportato come NA e ci_fallback_used = TRUE.

Le statistiche di bontà di adattamento includono R² di Tjur, R² di Cox–Snell, R² di McFadden, devianza penalizzata (-2 × L*), devianza penalizzata nulla e AICc (calcolato dalla log-verosimiglianza non penalizzata).

Curve dose-risposta

Le curve dose-risposta modellano la relazione sigmoide tra uno stimolo (concentrazione, dose) e una risposta (inibizione, stimolazione). BioRender adatta il modello logistico a quattro parametri (4PL) di default; le varianti a tre e cinque parametri sono casi particolari.

Modello a 4 parametri (4PL)

L'equazione logistica a quattro parametri è:

• Inibizione (log x): y = bottom + (top − bottom) / (1 + 10^((log_IC50 − x) × hill_slope))
• Inibizione (x lineare): y = bottom + (top − bottom) / (1 + (IC50 / x)^hill_slope)
• Stimolazione (x lineare): y = bottom + (top − bottom) × x^hill_slope / (x^hill_slope + EC50^hill_slope)

I quattro parametri liberi sono: bottom (risposta minima), top (risposta massima), hill_slope (coefficiente di Hill) e EC50/IC50 (concentrazione al punto medio tra bottom e top).

Modello a 3 parametri (3PL)

Equivalente al 4PL con la pendenza di Hill vincolata: hill_slope = −1 per inibizione, hill_slope = +1 per stimolazione. Tutti gli altri parametri rimangono liberi.

Logistica a 5 parametri (5PL)

Estende il 4PL aggiungendo un fattore di asimmetria S:

• Log x: y = bottom + (top − bottom) / (1 + 10^((log_EC50 − x) × hill_slope))^S
• x lineare (inibizione): y = bottom + (top − bottom) / (1 + (EC50/x)^hill_slope)^S

Quando S = 1 il 5PL si riduce al 4PL. S è inizializzato a 1 e stimato liberamente. Nota: il parametro log_EC50 nella formula 5PL è il punto di inflessione della curva. Quando S ≠ 1, il vero EC50 differisce ed è calcolato applicando la correzione: log_EC50_actual = log_inflection − log10(2^(1/S) − 1) / hill_slope.

Procedura di fitting

Tutti i modelli sono adattati usando nls dal pacchetto stats. I valori iniziali derivano dai dati: bottom e top sono inizializzati ai valori minimi e massimi osservati delle medie di gruppo; EC50/IC50 è inizializzato tramite interpolazione lineare al punto medio y; hill_slope è inizializzato a ±1 in base alla direzione dei dati.

L'algoritmo PORT (nl2sol) è provato per primo (algorithm = "port"), poiché supporta vincoli di limite ed evita il fallimento "minFactor" che l'algoritmo Gauss-Newton predefinito incontra sull'asse del profilo EC50. Se PORT fallisce e non sono attivi vincoli forniti dall'utente, si tenta un fallback con nls() senza algoritmo esplicito (default Gauss-Newton). Qualsiasi vincolo sui parametri specificato dall'utente è passato come vettori lower / upper a nls; la pendenza di Hill è inoltre limitata a ±50 per prevenire divergenze, con il segno fissato per corrispondere alla direzione dei dati.

Gli intervalli di confidenza per i parametri di adattamento sono calcolati tramite verosimiglianza profilo usando confint(fit). Se il profiling fallisce (es. gradiente singolare), si usa un'approssimazione di Wald: stima ± qt(0.975, df_residual) × SE. EC50/IC50 e il suo IC sono riportati sia sulla scala log₁₀ sia sulla scala originale di concentrazione.

Può essere calcolata anche una concentrazione personalizzata (EC_x o IC_x a una percentuale di risposta specificata dall'utente). Per inibizione, EC_x = EC50 × ((1/frazione − 1))^(1/hill_slope); per stimolazione, EC_x = EC50 × (frazione / (1 − frazione))^(1/hill_slope). L'IC al 95% sulla concentrazione personalizzata è approssimato propagando l'errore relativo dall'IC di EC50.

Le statistiche di bontà di adattamento riportate sono R² (1 − SS_res / SS_tot), somma dei quadrati residui, Sy.x (errore standard residuo) e gradi di libertà. Le bande di confidenza su griglie di concentrazione lineari e logaritmiche (250 punti ciascuna) sono calcolate con il metodo delta numerico: la formula del modello è derivata numericamente rispetto a ogni parametro per costruire la matrice Jacobiana J, quindi SE(ŷ) = √(J' × vcov(fit) × J) e la banda è ŷ ± t_{0.975, df} × SE(ŷ). Se il metodo delta fallisce, la banda ricade sulle stime puntuali.

Analisi di sopravvivenza

L'analisi di sopravvivenza modella il tempo fino a un evento di interesse (es. morte, recidiva) tenendo conto delle osservazioni censurate — soggetti che non hanno ancora sperimentato l'evento alla fine del follow-up.

Tutte le analisi di sopravvivenza usano funzioni dal pacchetto survival.

Curve di Kaplan–Meier

Le curve di sopravvivenza Kaplan–Meier sono stimate con survfit usando conf.type = "log-log" (trasformazione log-log complementare della funzione di sopravvivenza), metodo usato da GraphPad Prism e più adatto per probabilità di sopravvivenza vicine a 0 o 1 rispetto al metodo plain-log predefinito. I gruppi sono fattorizzati nell'ordine della prima apparizione nei dati per garantire che l'ordinamento delle curve corrisponda all'input.

Il tempo mediano di sopravvivenza e il suo IC al 95% sono estratti con quantile(km_fit, probs = 0.5). La curva è estesa al tempo massimo osservato in ogni gruppo (anche se quel tempo ha solo osservazioni censurate) affinché la funzione a gradini tracciata non sia troncata prematuramente.

Test log-rank

Il test log-rank confronta le distribuzioni di sopravvivenza tra gruppi sotto l'ipotesi nulla di rischi proporzionali. Usiamo survdiff con rho = 0. La statistica chi-quadrato è valutata su (numero di gruppi − 1) gradi di libertà.

Test di Gehan–Breslow–Wilcoxon

Alternativa al test log-rank che dà maggiore peso ai tempi iniziali. Questo test è implementato come una statistica chi-quadrato pesata personalizzata (non tramite survdiff). Ad ogni tempo di evento t, si applica un peso w_t = N_t (numero a rischio). Il vettore score è Z_j = Σ_t N_t (d_jt − e_jt) e la matrice di varianza-covarianza usa la varianza ipergeometrica con i termini N_t² applicati. La statistica del test χ² = Z' V⁻¹ Z è valutata su K−1 gradi di libertà (decomposizione di Cholesky; pseudoinversa QR come fallback).

Rapporti di rischio

Sono forniti due metodi di stima del rapporto di rischio per confronti a coppie tra gruppi.

Metodo logrank: Usa il rapporto osservato/sperato HR = (O₁/E₁) / (O₂/E₂), dove O ed E sono ottenuti dall'output log-rank di survdiff. L'errore standard di ln(HR) è √(1/E₁ + 1/E₂), dando IC al 95% di exp(ln(HR) ± 1.96 × SE). Questo corrisponde al rapporto di rischio "logrank" di GraphPad Prism.

Metodo Mantel–Haenszel: Usa la formula L = (O₁ − E₁) / V, HR = exp(L), con IC al 95% di exp(L ± 1.96 / √V), dove V è la varianza di (O₁ − E₁) sotto l'ipotesi nulla, presa da survdiff$var[1,1] di un test log-rank a coppie. Questo corrisponde al rapporto di rischio "Mantel–Haenszel" di GraphPad Prism (Machin, Cheung & Parmar; Vaeth).

Regressione proporzionale dei rischi di Cox

La regressione di Cox PH è adattata con coxph. In modalità auto-rilevamento, il metodo di gestione dei legami è selezionato automaticamente: "exact" quando non ci sono tempi di evento legati, e "efron" quando ci sono; il metodo dei legami può anche essere specificato esplicitamente. Gli intervalli di confidenza per i rapporti di rischio sono calcolati tramite una funzione di verosimiglianza profilo personalizzata (compute_profile_ci_cox()), che traccia la superficie di verosimiglianza profilo per trovare il confine dell'IC invece di affidarsi all'approssimazione di Wald.

L'assunzione di rischi proporzionali è testata usando cox.zph, che adatta i residui scalati di Schoenfeld nel tempo e riporta un valore p globale e per covariata; un valore p piccolo suggerisce che il rapporto di rischio cambia nel tempo, violando l'assunzione PH.

Correzioni per confronti multipli

Quando si confrontano tre o più gruppi e viene selezionato un metodo di correzione, i valori p dei test a coppie sono corretti per confronti multipli. I metodi disponibili sono Bonferroni e Benjamini–Hochberg (tasso di falsi positivi), entrambi applicati tramite p.adjust, e Holm–Šídák, implementato con una procedura stepwise personalizzata.