0.1 Obiettivo dell’analisi

L’obiettivo di questo lavoro è costruire un modello statistico in grado di prevedere il peso dei neonati alla nascita a partire da variabili cliniche e antropometriche raccolte presso tre ospedali.

L’analisi si articola in più fasi: controllo preliminare dei dati, analisi descrittiva, studio delle associazioni, costruzione del modello di regressione, selezione del modello finale, diagnostica e previsione in uno scenario clinico ipotetico.

1 Preparazione e controllo qualità dei dati

for(v in c("Sesso", "Tipo.parto", "Ospedale")){
  if(v %in% names(dati)) dati[[v]] <- factor(dati[[v]])
}

if("Fumatrici" %in% names(dati)){
  dati$Fumatrici <- factor(dati$Fumatrici, levels = c(0,1), labels = c("No","Si"))
}

if("N.gravidanze" %in% names(dati)){
  dati$N.gravidanze <- suppressWarnings(as.numeric(dati$N.gravidanze))
}

if("Anni.madre" %in% names(dati)){
  dati$Anni.madre <- suppressWarnings(as.numeric(dati$Anni.madre))
}

n <- nrow(dati)
p <- ncol(dati)

kable(data.frame(Osservazioni = n, Variabili = p),
      caption = "Dimensioni del dataset")
Dimensioni del dataset
Osservazioni Variabili
2500 10 
summary(dati)
##    Anni.madre     N.gravidanze     Fumatrici   Gestazione         Peso     
##  Min.   : 0.00   Min.   : 0.0000   No:2396   Min.   :25.00   Min.   : 830  
##  1st Qu.:25.00   1st Qu.: 0.0000   Si: 104   1st Qu.:38.00   1st Qu.:2990  
##  Median :28.00   Median : 1.0000             Median :39.00   Median :3300  
##  Mean   :28.16   Mean   : 0.9812             Mean   :38.98   Mean   :3284  
##  3rd Qu.:32.00   3rd Qu.: 1.0000             3rd Qu.:40.00   3rd Qu.:3620  
##  Max.   :46.00   Max.   :12.0000             Max.   :43.00   Max.   :4930  
##    Lunghezza         Cranio    Tipo.parto Ospedale   Sesso   
##  Min.   :310.0   Min.   :235   Ces: 728   osp1:816   F:1256  
##  1st Qu.:480.0   1st Qu.:330   Nat:1772   osp2:849   M:1244  
##  Median :500.0   Median :340              osp3:835           
##  Mean   :494.7   Mean   :340                                 
##  3rd Qu.:510.0   3rd Qu.:350                                 
##  Max.   :565.0   Max.   :390
struttura <- data.frame(
  Variabile = names(dati),
  Classe = sapply(dati, class),
  Valori_mancanti = sapply(dati, function(x) sum(is.na(x))),
  row.names = NULL
)

kable(struttura, caption = "Struttura delle variabili e valori mancanti")
Struttura delle variabili e valori mancanti
Variabile Classe Valori_mancanti
Anni.madre numeric 0
N.gravidanze numeric 0
Fumatrici factor 0
Gestazione integer 0
Peso integer 0
Lunghezza integer 0
Cranio integer 0
Tipo.parto factor 0
Ospedale factor 0
Sesso factor 0 
check_plaus <- list()

if("Anni.madre" %in% names(dati)){
  idx_eta_err <- which(dati$Anni.madre < 12 | dati$Anni.madre > 60)
  check_plaus$eta_madre_implausibile <- dati[idx_eta_err, , drop = FALSE]
}

if("Gestazione" %in% names(dati)){
  idx_gest_err <- which(dati$Gestazione < 20 | dati$Gestazione > 45)
  check_plaus$gestazione_implausibile <- dati[idx_gest_err, , drop = FALSE]
}

if("Peso" %in% names(dati)){
  idx_peso_err <- which(dati$Peso < 300 | dati$Peso > 7000)
  check_plaus$peso_implausibile <- dati[idx_peso_err, , drop = FALSE]
}

if("Lunghezza" %in% names(dati)){
  idx_lung_err <- which(dati$Lunghezza < 200 | dati$Lunghezza > 700)
  check_plaus$lunghezza_implausibile <- dati[idx_lung_err, , drop = FALSE]
}

if("Cranio" %in% names(dati)){
  idx_cranio_err <- which(dati$Cranio < 150 | dati$Cranio > 450)
  check_plaus$cranio_implausibile <- dati[idx_cranio_err, , drop = FALSE]
}

# Mostro in tabella solo i casi davvero problematici

if(nrow(check_plaus$eta_madre_implausibile) > 0){
  kable(check_plaus$eta_madre_implausibile,
        caption = "Osservazioni con età materna implausibile")
}
Osservazioni con età materna implausibile
Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio Tipo.parto Ospedale Sesso
1152 1 1 No 41 3250 490 350 Nat osp2 F
1380 0 0 No 39 3060 490 330 Nat osp3 M 
# Trattamento valori errati dell' età della madre

dati$Anni.madre[dati$Anni.madre < 12 | dati$Anni.madre > 60] <- NA

kable(data.frame(
  Minimo_post_correzione = min(dati$Anni.madre, na.rm = TRUE),
  Massimo_post_correzione = max(dati$Anni.madre, na.rm = TRUE),
  NA_totali = sum(is.na(dati$Anni.madre))
), caption = "Controllo di Anni.madre dopo il trattamento dei valori non utilizzabili")
Controllo di Anni.madre dopo il trattamento dei valori non utilizzabili
Minimo_post_correzione Massimo_post_correzione NA_totali
13 46 2

Il summary() iniziale ha evidenziato valori implausibili nella variabile Anni.madre. Il controllo di plausibilità ha confermato la presenza di due osservazioni errate (0 e 1 anno), verosimilmente dovute a errori di inserimento. Tali valori sono stati pertanto trattati come mancanti (NA) prima delle analisi successive.

2 Analisi descrittiva

vars_num <- names(dati)[sapply(dati, is.numeric)]

desc_num <- do.call(rbind, lapply(vars_num, function(v){
  x <- dati[[v]]
  data.frame(
    Variabile = v,
    N = sum(!is.na(x)),
    Media = mean(x, na.rm = TRUE),
    Mediana = median(x, na.rm = TRUE),
    SD = sd(x, na.rm = TRUE),
    Min = min(x, na.rm = TRUE),
    Q1 = quantile(x, 0.25, na.rm = TRUE),
    Q3 = quantile(x, 0.75, na.rm = TRUE),
    Max = max(x, na.rm = TRUE),
    Skewness = skewness(x, na.rm = TRUE),
    Kurtosis_excess = kurtosis(x, na.rm = TRUE) - 3
  )
}))

kable(as.data.frame(lapply(desc_num, rd)),
      caption = "Statistiche descrittive delle variabili numeriche")
Statistiche descrittive delle variabili numeriche
Variabile N Media Mediana SD Min Q1 Q3 Max Skewness Kurtosis_excess
Anni.madre 2498 28.186 28 5.217 13 25 32 46 0.151 -0.106
N.gravidanze 2500 0.981 1 1.281 0 0 1 12 2.514 10.989
Gestazione 2500 38.980 39 1.869 25 38 40 43 -2.065 8.258
Peso 2500 3284.081 3300 525.039 830 2990 3620 4930 -0.647 2.032
Lunghezza 2500 494.692 500 26.319 310 480 510 565 -1.515 6.487
Cranio 2500 340.029 340 16.425 235 330 350 390 -0.785 2.946 
# 4.2 Focus sulla variabile risposta: Peso

peso_stats <- data.frame(
  Media = mean(dati$Peso, na.rm = TRUE),
  Mediana = median(dati$Peso, na.rm = TRUE),
  SD = sd(dati$Peso, na.rm = TRUE),
  Skewness = skewness(dati$Peso, na.rm = TRUE),
  Kurtosis_excess = kurtosis(dati$Peso, na.rm = TRUE) - 3,
  Shapiro_W = shapiro.test(dati$Peso)$statistic,
  Shapiro_p = shapiro.test(dati$Peso)$p.value
)

kable(as.data.frame(lapply(peso_stats, rd)),
      caption = "Statistiche descrittive del peso")
Statistiche descrittive del peso
Media Mediana SD Skewness Kurtosis_excess Shapiro_W Shapiro_p
3284.081 3300 525.039 -0.647 2.032 0.971 0 
ggplot(dati, aes(x = Peso)) +
  geom_histogram(bins = 30) +
  theme_classic() +
  labs(title = "Distribuzione del peso neonatale",
       x = "Peso (g)",
       y = "Frequenza")

ggplot(dati, aes(sample = Peso)) +
  stat_qq() +
  stat_qq_line() +
  theme_classic() +
  labs(title = "QQ-plot del peso")

L’analisi descrittiva evidenzia valori medi plausibili per le principali variabili del dataset. L’età media delle madri è di circa 28 anni, la durata media della gestazione di circa 39 settimane e il peso medio alla nascita di circa 3284 grammi.

La distribuzione del peso neonatale risulta relativamente concentrata attorno alla media (deviazione standard circa 525 grammi). L’istogramma e il QQ-plot mostrano una distribuzione complessivamente vicina alla normalità; il test di Shapiro-Wilk risulta significativo.

3 Analisi esplorativa delle associazioni

# prendo solo variabili numeriche

dati_num <- dati[, sapply(dati, is.numeric)]

# elimino variabili con tutti NA o con deviazione standard zero

vars_ok <- sapply(dati_num, function(x){
  sum(!is.na(x)) > 1 & sd(x, na.rm = TRUE) > 0
})

dati_num <- dati_num[, vars_ok]

# matrice di correlazione

cor_mat <- cor(dati_num, use = "pairwise.complete.obs")

kable(as.data.frame(round(cor_mat, 2)),
      caption = "Matrice di correlazione")
Matrice di correlazione
Anni.madre N.gravidanze Gestazione Peso Lunghezza Cranio
Anni.madre 1.00 0.38 -0.13 -0.02 -0.06 0.02
N.gravidanze 0.38 1.00 -0.10 0.00 -0.06 0.04
Gestazione -0.13 -0.10 1.00 0.59 0.62 0.46
Peso -0.02 0.00 0.59 1.00 0.80 0.70
Lunghezza -0.06 -0.06 0.62 0.80 1.00 0.60
Cranio 0.02 0.04 0.46 0.70 0.60 1.00 
corrplot(cor_mat,
         method = "color",
         type = "upper",
         tl.col = "black",
         tl.srt = 45)

La matrice di correlazione mostra che il peso neonatale è positivamente associato soprattutto a lunghezza, diametro cranico e settimane di gestazione. Queste variabili appaiono quindi candidate naturali per la successiva modellazione del peso alla nascita.

4 Verifica delle ipotesi del progetto

In questa sezione vengono esaminate alcune relazioni tra le variabili del dataset mediante test statistici appropriati.

4.1 Associazione tra ospedale e tipo di parto

if(all(c("Ospedale","Tipo.parto") %in% names(dati))){
  tab_ot <- table(dati$Ospedale, dati$Tipo.parto)

  kable(as.data.frame.matrix(tab_ot),
        caption = "Tabella di contingenza Ospedale × Tipo di parto")

  chi_ot <- chisq.test(tab_ot)

  kable(data.frame(
    X2 = unname(chi_ot$statistic),
    df = unname(chi_ot$parameter),
    p_value = chi_ot$p.value
  ) |> as.data.frame() |> lapply(rd) |> as.data.frame(),
  caption = "Test chi-quadrato di associazione")
}
Test chi-quadrato di associazione
X2 df p_value
1.097 2 0.578

Il test chi-quadrato non evidenzia un’associazione statisticamente significativa tra ospedale e tipo di parto. Queste variabili vengono pertanto considerate solo a scopo descrittivo e non incluse nel modello predittivo finale, in quanto poco rilevanti in ottica previsionale.

4.2 Confronto con valori di riferimento

mu_peso <- 3300
mu_lung <- 500

t_peso <- t.test(dati$Peso, mu = mu_peso)
t_lung <- t.test(dati$Lunghezza, mu = mu_lung)

tab_test <- data.frame(
  Variabile = c("Peso", "Lunghezza"),
  Mu_pop = c(mu_peso, mu_lung),
  Media_camp = c(mean(dati$Peso, na.rm = TRUE),
                 mean(dati$Lunghezza, na.rm = TRUE)),
  t = c(unname(t_peso$statistic),
        unname(t_lung$statistic)),
  p_value = c(t_peso$p.value,
              t_lung$p.value)
)

kable(as.data.frame(lapply(tab_test, rd)),
      caption = "Test t a un campione")
Test t a un campione
Variabile Mu_pop Media_camp t p_value
Peso 3300 3284.081 -1.516 0.13
Lunghezza 500 494.692 -10.084 0.00

I test t a un campione confrontano le medie osservate nel dataset con valori di riferimento tratti da dati statistici sui nati in Italia pubblicati da ISTAT (2025). I risultati mostrano che il peso medio del campione non risulta significativamente diverso dal valore di riferimento, mentre la lunghezza media risulta significativamente inferiore.

4.3 Differenze per sesso nelle misure antropometriche

if("Sesso" %in% names(dati)){

  vars_sesso <- c("Peso", "Lunghezza", "Cranio")

  tab_sesso <- do.call(rbind, lapply(vars_sesso, function(v){

    test <- t.test(dati[[v]] ~ dati$Sesso, na.action = na.omit)

    medie <- tapply(dati[[v]], dati$Sesso, mean, na.rm = TRUE)

    data.frame(
      Variabile = v,
      Media_F = if("F" %in% names(medie)) medie["F"] else NA,
      Media_M = if("M" %in% names(medie)) medie["M"] else NA,
      t = unname(test$statistic),
      p_value = test$p.value
    )
  }))

  kable(as.data.frame(lapply(tab_sesso, rd)),
        caption = "Differenze per sesso nelle principali misure antropometriche")
}
Differenze per sesso nelle principali misure antropometriche
Variabile Media_F Media_M t p_value
Peso 3161.132 3408.215 -12.106 0
Lunghezza 489.764 499.667 -9.582 0
Cranio 337.633 342.449 -7.410 0

I test t mostrano che peso, lunghezza e diametro cranico differiscono significativamente tra maschi e femmine, con valori mediamente più elevati nei maschi. Questo suggerisce che il sesso possa rappresentare un potenziale predittore del peso neonatale e giustifica la sua considerazione nella fase di modellazione.

5 Costruzione del modello di regressione

5.1 Modello completo con sole variabili clinicamente significative

# Dataset per la modellazione: casi completi sulle variabili candidate

vars_modello <- c("Peso", "Anni.madre", "N.gravidanze", "Fumatrici",
                  "Gestazione", "Lunghezza", "Cranio", "Sesso")

dati_mod <- dati[complete.cases(dati[, vars_modello]), vars_modello]

kable(data.frame(
  Osservazioni_iniziali = nrow(dati),
  Osservazioni_usate_modello = nrow(dati_mod),
  Osservazioni_escluse = nrow(dati) - nrow(dati_mod)
), caption = "Osservazioni utilizzate per la modellazione")
Osservazioni utilizzate per la modellazione
Osservazioni_iniziali Osservazioni_usate_modello Osservazioni_escluse
2500 2498 2 
mod_full <- lm(Peso ~ Anni.madre + N.gravidanze + Fumatrici +
                 Gestazione + Lunghezza + Cranio + Sesso,
               data = dati_mod)

tab_full <- as.data.frame(summary(mod_full)$coefficients)
tab_full$Termine <- rownames(tab_full)
rownames(tab_full) <- NULL
tab_full <- tab_full[, c("Termine","Estimate","Std. Error","t value","Pr(>|t|)")]

kable(as.data.frame(lapply(tab_full, rd)),
      caption = "Coefficienti del modello completo")
Coefficienti del modello completo
Termine Estimate Std..Error t.value Pr…t..
(Intercept) -6712.240 141.334 -47.492 0.000
Anni.madre 0.880 1.149 0.766 0.444
N.gravidanze 11.379 4.677 2.433 0.015
FumatriciSi -30.396 27.608 -1.101 0.271
Gestazione 32.947 3.829 8.605 0.000
Lunghezza 10.232 0.301 33.979 0.000
Cranio 10.520 0.427 24.633 0.000
SessoM 78.079 11.213 6.963 0.000 
kable(data.frame(
  R2 = summary(mod_full)$r.squared,
  Adj_R2 = summary(mod_full)$adj.r.squared,
  AIC = AIC(mod_full),
  BIC = BIC(mod_full)
) |> as.data.frame() |> lapply(rd) |> as.data.frame(),
caption = "Indicatori di adattamento del modello completo")
Indicatori di adattamento del modello completo
R2 Adj_R2 AIC BIC
0.727 0.726 35155.07 35207.48

Il modello completo è stato costruito includendo solo variabili con plausibile rilevanza clinica e predittiva: età materna, numero di gravidanze, fumo materno, settimane di gestazione, lunghezza, diametro cranico e sesso del neonato. Sono state invece escluse a priori variabili come ospedale e tipo di parto, ritenute poco informative in ottica previsiva.

5.2 Verifica della multicollinearità

v <- car::vif(mod_full)

if (is.matrix(v)) {
  vif_tab <- data.frame(
    Variabile = rownames(v),
    GVIF = round(v[, "GVIF"], 2),
    Df = v[, "Df"],
    GVIF_adj = round(v[, "GVIF^(1/(2*Df))"], 2),
    row.names = NULL
  )
  kable(vif_tab, caption = "Multicollinearità (GVIF) - modello completo")
} else {
  vif_tab <- data.frame(
    Variabile = names(v),
    VIF = round(as.numeric(v), 2)
  )
  kable(vif_tab, caption = "Multicollinearità (VIF) - modello completo")
}
Multicollinearità (VIF) - modello completo
Variabile VIF
Anni.madre 1.19
N.gravidanze 1.19
Fumatrici 1.01
Gestazione 1.69
Lunghezza 2.08
Cranio 1.63
Sesso 1.04

I valori di VIF/GVIF risultano contenuti e non segnalano criticità rilevanti di multicollinearità. Le stime del modello possono quindi essere interpretate con sufficiente affidabilità.

5.3 Scrematura delle variabili tramite criterio AIC

mod_step <- MASS::stepAIC(mod_full, direction = "both", trace = FALSE)

kable(data.frame(
  Modello = c("Completo", "Selezionato_AIC"),
  AIC = c(AIC(mod_full), AIC(mod_step)),
  BIC = c(BIC(mod_full), BIC(mod_step)),
  R2 = c(summary(mod_full)$r.squared, summary(mod_step)$r.squared),
  Adj_R2 = c(summary(mod_full)$adj.r.squared, summary(mod_step)$adj.r.squared)
) |> as.data.frame() |> lapply(rd) |> as.data.frame(),
caption = "Confronto tra modello completo e modello selezionato")
Confronto tra modello completo e modello selezionato
Modello AIC BIC R2 Adj_R2
Completo 35155.07 35207.48 0.727 0.726
Selezionato_AIC 35152.89 35193.65 0.727 0.726

La procedura di selezione tramite AIC consente di ottenere un modello più parsimonioso, mantenendo una buona capacità esplicativa. Il modello selezionato rappresenta il punto di partenza per le verifiche successive su eventuali effetti di secondo ordine.

5.4 Verifica di effetti di secondo ordine

# Modello base per valutare effetti di secondo ordine

mod_base_secondo <- lm(Peso ~ Anni.madre + N.gravidanze + Fumatrici +
                         Gestazione + Lunghezza + Cranio + Sesso,
                       data = dati_mod)

# Interazione tra gestazione e fumo

mod_int <- update(mod_base_secondo, . ~ . + Gestazione:Fumatrici)

# Termine quadratico della gestazione

mod_quad <- update(mod_base_secondo, . ~ . + I(Gestazione^2))

cmp_secondo <- data.frame(
  Modello = c("Base", "Interazione", "Quadratico"),
  AIC = c(AIC(mod_base_secondo), AIC(mod_int), AIC(mod_quad)),
  BIC = c(BIC(mod_base_secondo), BIC(mod_int), BIC(mod_quad)),
  R2 = c(summary(mod_base_secondo)$r.squared,
         summary(mod_int)$r.squared,
         summary(mod_quad)$r.squared)
)

kable(as.data.frame(lapply(cmp_secondo, rd)),
      caption = "Confronto tra modelli con effetti di secondo ordine")
Confronto tra modelli con effetti di secondo ordine
Modello AIC BIC R2
Base 35155.07 35207.48 0.727
Interazione 35155.90 35214.13 0.727
Quadratico 35151.79 35210.02 0.728 
anova_int <- anova(mod_base_secondo, mod_int)

kable(as.data.frame(lapply(as.data.frame(anova_int), rd)),
      caption = "Test formale dell'interazione Gestazione × Fumatrici")
Test formale dell’interazione Gestazione × Fumatrici
Res.Df RSS Df Sum.of.Sq F Pr..F.
2490 187905214 NA NA NA NA
2489 187817498 1 87715.43 1.162 0.281

Sono stati esplorati due possibili effetti di secondo ordine: un’interazione tra gestazione e fumo materno e un termine quadratico per la gestazione. In assenza di un miglioramento sostanziale e/o di evidenza statistica convincente, si privilegia il modello più parsimonioso.

5.5 Modello finale

mod_final <- mod_step

tab_final <- as.data.frame(summary(mod_final)$coefficients)
tab_final$Termine <- rownames(tab_final)
rownames(tab_final) <- NULL
tab_final <- tab_final[, c("Termine","Estimate","Std. Error","t value","Pr(>|t|)")]

kable(as.data.frame(lapply(tab_final, rd)),
      caption = "Coefficienti del modello finale")
Coefficienti del modello finale
Termine Estimate Std..Error t.value Pr…t..
(Intercept) -6681.725 135.804 -49.201 0.000
N.gravidanze 12.455 4.342 2.869 0.004
Gestazione 32.383 3.801 8.520 0.000
Lunghezza 10.245 0.301 34.059 0.000
Cranio 10.541 0.426 24.717 0.000
SessoM 77.981 11.211 6.956 0.000 
kable(data.frame(
  R2 = summary(mod_final)$r.squared,
  Adj_R2 = summary(mod_final)$adj.r.squared,
  RMSE = sqrt(mean(residuals(mod_final)^2, na.rm = TRUE))
) |> as.data.frame() |> lapply(rd) |> as.data.frame(),
caption = "Prestazioni del modello finale")
Prestazioni del modello finale
R2 Adj_R2 RMSE
0.727 0.726 274.367

Il modello finale evidenzia che il peso alla nascita è significativamente associato a diverse variabili cliniche e antropometriche. In particolare, la durata della gestazione mostra un effetto positivo sul peso neonatale: a parità delle altre variabili incluse nel modello, un aumento di una settimana di gestazione è associato a un incremento medio del peso alla nascita pari al coefficiente stimato per la variabile Gestazione.

Anche le misure antropometriche del neonato risultano fortemente associate al peso: aumenti della lunghezza e del diametro cranico sono associati a incrementi del peso alla nascita. Le variabili lunghezza e diametro cranico rappresentano dimensioni corporee correlate; tuttavia, l’analisi della multicollinearità non evidenzia criticità rilevanti, consentendo di mantenerle entrambe nel modello.

Il sesso del neonato e il numero di gravidanze risultano inoltre variabili statisticamente significative. In particolare, i neonati maschi presentano mediamente un peso alla nascita superiore rispetto alle femmine, a parità delle altre variabili incluse nel modello.

5.6 Diagnostica dei residui

RMSE_final <- sqrt(mean(residuals(mod_final)^2, na.rm = TRUE))
kable(data.frame(RMSE = round(RMSE_final, 2)),
      caption = "RMSE del modello finale")
RMSE del modello finale
RMSE
274.37 
op <- par(mfrow = c(2,2))
plot(mod_final)

par(op)

tab_diag_final <- data.frame(
  Test = c("Shapiro-Wilk", "Breusch-Pagan", "Durbin-Watson"),
  Stat = c(shapiro.test(residuals(mod_final))$statistic,
           bptest(mod_final)$statistic,
           dwtest(mod_final)$statistic),
  p_value = c(shapiro.test(residuals(mod_final))$p.value,
              bptest(mod_final)$p.value,
              dwtest(mod_final)$p.value)
)

kable(as.data.frame(lapply(tab_diag_final, rd)),
      caption = "Test diagnostici sui residui del modello finale")
Test diagnostici sui residui del modello finale
Test Stat p_value
Shapiro-Wilk 0.974 0.000
Breusch-Pagan 90.297 0.000
Durbin-Watson 1.953 0.121

I grafici diagnostici e i test formali consentono di valutare l’adeguatezza del modello. Eventuali deviazioni dalla normalità o segnali di eteroschedasticità vanno interpretati con cautela, soprattutto considerando la numerosità campionaria elevata.

5.7 Studio degli outlier e dei punti influenti

cook_f <- cooks.distance(mod_final)
lev_f <- hatvalues(mod_final)
rstud_f <- rstudent(mod_final)

soglia_cook <- 4 / nrow(model.frame(mod_final))

tab_infl <- data.frame(
  Osservazione = seq_along(cook_f),
  Cook = cook_f,
  Leverage = lev_f,
  Rstudent = rstud_f
)

influenti <- tab_infl %>%
  dplyr::filter(Cook > soglia_cook | abs(Rstudent) > 3)

kable(head(as.data.frame(lapply(influenti, rd)), 15),
      caption = "Osservazioni potenzialmente influenti")
Osservazioni potenzialmente influenti
Osservazione Cook Leverage Rstudent
13 0.002 0.006 1.320
34 0.002 0.007 -1.372
119 0.005 0.003 -3.060
130 0.003 0.002 3.199
134 0.004 0.008 1.714
146 0.002 0.002 2.547
155 0.030 0.007 5.025
161 0.003 0.020 -0.850
204 0.002 0.014 0.901
220 0.004 0.007 1.858
295 0.004 0.004 -2.525
310 0.068 0.029 -3.705
312 0.002 0.013 0.870
329 0.003 0.003 2.699
348 0.002 0.005 1.553 
if(nrow(influenti) > 0){
  righe_infl <- dati_mod[influenti$Osservazione, , drop = FALSE]
  righe_infl$Osservazione <- influenti$Osservazione
  kable(head(righe_infl, 15),
        caption = "Righe del dataset corrispondenti alle osservazioni influenti")
}
Righe del dataset corrispondenti alle osservazioni influenti
Peso Anni.madre N.gravidanze Fumatrici Gestazione Lunghezza Cranio Sesso Osservazione
13 3060 36 5 No 38 455 325 F 13
34 3150 27 0 No 39 480 382 F 34
119 3410 31 0 No 40 550 372 M 119
130 4240 30 2 No 39 485 352 M 130
134 3950 38 6 No 37 500 350 M 134
146 3820 24 1 No 40 500 320 F 146
155 3610 30 0 No 36 410 330 M 155
161 3760 35 9 No 42 540 348 F 161
204 3850 30 8 No 40 518 340 F 204
220 3520 23 1 No 40 445 363 F 220
295 1850 18 0 No 40 460 305 F 295
310 1560 40 3 No 28 420 379 F 310
312 1280 26 1 No 32 360 276 M 312
329 4560 25 1 No 40 540 340 M 329
348 3560 32 0 No 38 460 360 M 348

Le osservazioni influenti non vengono rimosse automaticamente. In accordo con un corretto approccio statistico, tali casi vengono esaminati per valutarne la plausibilità clinica e l’eventuale impatto sulle stime, ma restano inclusi nel modello salvo evidenza di errore manifesto nei dati.

5.8 Forecasting

pred_vars <- all.vars(delete.response(terms(mod_final)))

scenario <- list(
  Lunghezza = 500,
  Cranio = 340,
  Gestazione = 39,
  N.gravidanze = 3,
  Anni.madre = 30,
  Fumatrici = "No",
  Sesso = "F"
)

newd <- as.data.frame(setNames(replicate(length(pred_vars), NA, simplify = FALSE), pred_vars))

for (v in pred_vars) {
  if (v %in% names(scenario)) {
    if (is.factor(dati_mod[[v]])) {
      newd[[v]] <- factor(scenario[[v]], levels = levels(dati_mod[[v]]))
    } else {
      newd[[v]] <- scenario[[v]]
    }
  } else if (is.numeric(dati_mod[[v]])) {
    newd[[v]] <- median(dati_mod[[v]], na.rm = TRUE)
  } else if (is.factor(dati_mod[[v]])) {
    mode_lv <- names(sort(table(dati_mod[[v]]), decreasing = TRUE))[1]
    newd[[v]] <- factor(mode_lv, levels = levels(dati_mod[[v]]))
  }
}

predF <- predict(mod_final, newdata = newd, interval = "prediction")

kable(data.frame(
  Peso_previsto = round(predF[1], 1),
  PI_low = round(predF[2], 1),
  PI_high = round(predF[3], 1)
), caption = "Previsione del peso neonatale in uno scenario clinico ipotetico")
Previsione del peso neonatale in uno scenario clinico ipotetico
Peso_previsto PI_low PI_high
3325.2 2786 3864.4

Il modello è stato utilizzato per stimare il peso alla nascita in uno scenario clinico ipotetico caratterizzato da una gestazione di 39 settimane, lunghezza del neonato pari a 500 mm, diametro cranico di 340 mm, madre di 30 anni, tre gravidanze pregresse, madre non fumatrice e neonato di sesso femminile. In questo scenario il peso previsto alla nascita risulta pari a circa 3325 grammi.

6 Conclusioni

L’analisi ha permesso di costruire un modello di regressione lineare multipla per la previsione del peso neonatale, seguendo un percorso metodologico strutturato che ha incluso il controllo preliminare dei dati, l’analisi descrittiva, lo studio delle associazioni tra variabili e la successiva costruzione e selezione del modello.

Il modello finale evidenzia che il peso alla nascita è principalmente associato alla durata della gestazione e alle principali misure antropometriche del neonato, in particolare lunghezza e diametro cranico. Anche il sesso del neonato e il numero di gravidanze risultano variabili statisticamente significative, con i neonati maschi che tendono ad avere un peso mediamente superiore rispetto alle femmine.

Dal punto di vista quantitativo, il modello è in grado di spiegare circa il 72% della variabilità del peso neonatale, con un errore medio di previsione pari a circa 274 grammi. Questi risultati indicano una buona capacità esplicativa del modello e suggeriscono che le variabili considerate rappresentano efficacemente i principali fattori associati al peso alla nascita.