Analisi del Peso dei Neonati

Questo documento presenta l’analisi del peso dei neonati utilizzando R. L’obiettivo è costruire un modello di regressione lineare multipla per stimare il peso del neonato in base a diverse variabili materne e neonatali.

1. Raccolta Dati e Struttura del Dataset

Il dataset “neonati.csv” è stato caricato in R. Le prime righe del dataset e la sua struttura sono state esaminate per comprendere le variabili e i loro tipi.

setwd("C:/Users/Admin/Desktop/Master AI") # Impostare la directory

dati <- read.csv("neonati.csv", stringsAsFactors = TRUE, sep = ",") # Caricare i dati

str(dati)
## 'data.frame':    2500 obs. of  10 variables:
##  $ Anni.madre  : int  26 21 34 28 20 32 26 25 22 23 ...
##  $ N.gravidanze: int  0 2 3 1 0 0 1 0 1 0 ...
##  $ Fumatrici   : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ Gestazione  : int  42 39 38 41 38 40 39 40 40 41 ...
##  $ Peso        : int  3380 3150 3640 3690 3700 3200 3100 3580 3670 3700 ...
##  $ Lunghezza   : int  490 490 500 515 480 495 480 510 500 510 ...
##  $ Cranio      : int  325 345 375 365 335 340 345 349 335 362 ...
##  $ Tipo.parto  : Factor w/ 2 levels "Ces","Nat": 2 2 2 2 2 2 2 2 1 1 ...
##  $ Ospedale    : Factor w/ 3 levels "osp1","osp2",..: 3 1 2 2 3 2 3 1 2 2 ...
##  $ Sesso       : Factor w/ 2 levels "F","M": 2 1 2 2 1 1 1 2 1 1 ...

Il dataset è composto da 2500 osservazioni su 10 variabili. Le prime sette variabili sono di tipo numerico intero e rappresentano caratteristiche materne e neonatali, mentre le ultime tre variabili sono categoriche.

Analizziamo le statistiche di queste variabili:

kable(summary(dati)) # Riepilogo statistico del dataset
Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio Tipo.parto Ospedale Sesso
Min. : 0.00 Min. : 0.0000 Min. :0.0000 Min. :25.00 Min. : 830 Min. :310.0 Min. :235 Ces: 728 osp1:816 F:1256
1st Qu.:25.00 1st Qu.: 0.0000 1st Qu.:0.0000 1st Qu.:38.00 1st Qu.:2990 1st Qu.:480.0 1st Qu.:330 Nat:1772 osp2:849 M:1244
Median :28.00 Median : 1.0000 Median :0.0000 Median :39.00 Median :3300 Median :500.0 Median :340 NA osp3:835 NA
Mean :28.16 Mean : 0.9812 Mean :0.0416 Mean :38.98 Mean :3284 Mean :494.7 Mean :340 NA NA NA
3rd Qu.:32.00 3rd Qu.: 1.0000 3rd Qu.:0.0000 3rd Qu.:40.00 3rd Qu.:3620 3rd Qu.:510.0 3rd Qu.:350 NA NA NA
Max. :46.00 Max. :12.0000 Max. :1.0000 Max. :43.00 Max. :4930 Max. :565.0 Max. :390 NA NA NA

Durante l’esplorazione iniziale dei dati, è stata rilevata un’anomalia nella variabile Anni.madre:

  • Min: 0.00

Per verificare la presenza di osservazioni con età materna sospetta, si esegue il seguente filtro:

kable(dati[dati$Anni.madre < 12 | dati$Anni.madre > 50, ])
Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio Tipo.parto Ospedale Sesso
1152 1 1 0 41 3250 490 350 Nat osp2 F
1380 0 0 0 39 3060 490 330 Nat osp3 M

Dall’output è emersa una osservazione con Anni.madre = 0 ed una osservazione con Anni.madre = 1. Poiché l’età pari a 0 o 1 rappresenta un chiaro errore di inserimento e non vi sono informazioni aggiuntive per correggere il valore, si è deciso di escludere i record dal dataset:

dati <- dati[dati$Anni.madre != 0 & dati$Anni.madre != 1, ]

summary(dati$Anni.madre)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   13.00   25.00   28.00   28.19   32.00   46.00
attach(dati)  # Rende le variabili accessibili

# Verifica asimmetria e curtosi della variabile Peso
skewness <- moments::skewness(Peso)
kurtosis <- moments::kurtosis(Peso) - 3

tab <- data.frame(
  Statistica = c("Skewness", "Kurtosis (excess)"),
  Valore = round(c(skewness, kurtosis), 3)
)

# Stampa tabella
kable(tab, caption = "Skewness e Kurtosi della variabile Peso")
Skewness e Kurtosi della variabile Peso
Statistica Valore
Skewness -0.647
Kurtosis (excess) 2.029

L’asimmetria (skewness) risulta pari a -0.647, indicando una lieve asimmetria negativa: la distribuzione del peso tende leggermente verso sinistra. Tuttavia, il valore è vicino a zero, quindi possiamo considerare la distribuzione quasi simmetrica.

La curtosi, pari a 2.03, è positiva e superiore a zero, suggerendo una distribuzione più appuntita rispetto a quella normale (leptocurtica).

Dunque, la distribuzione del peso si discosta leggermente dalla normalità, ma non in modo marcato.

2. Analisi e Modellizzazione

2.1 Analisi Preliminare

Variabili Numeriche

Si definisce una funzione analisi_descrittiva() per calcolare statistiche descrittive (media, mediana, deviazione standard, quartili, asimmetria, curtosi) per le variabili numeriche.

analisi_descrittiva <- function(x, nome_var = "Variabile") {
  # Calcolo delle statistiche
  media <- mean(x, na.rm = TRUE)
  mediana <- median(x, na.rm = TRUE)
  deviazione <- sd(x, na.rm = TRUE)
  quartili <- quantile(x, probs = c(0.25, 0.75), na.rm = TRUE)
  asimmetria <- moments::skewness(x, na.rm = TRUE)
  curtosi <- moments::kurtosis(x, na.rm = TRUE) - 3
  
  # Crea un data frame con i risultati
  risultati <- data.frame(
    Variabile = nome_var,
    Media = media,
    Mediana = mediana,
    DeviazioneStandard = deviazione,
    Quartile_25 = unname(quartili[1]),  # Usa unname per rimuovere i nomi
    Quartile_75 = unname(quartili[2]),  # Usa unname per rimuovere i nomi
    Asimmetria = asimmetria,
    Curtosi = curtosi
  )
  
  return(risultati)
}

# Calcola le statistiche per ogni variabile e combina i risultati
risultati_peso <- analisi_descrittiva(dati$Peso, "Peso")
risultati_anni_madre <- analisi_descrittiva(dati$Anni.madre, "Anni.madre")
risultati_n_gravidanze <- analisi_descrittiva(dati$N.gravidanze, "Numero Gravidanze")
risultati_gestazione <- analisi_descrittiva(dati$Gestazione, "Gestazione")
risultati_lunghezza <- analisi_descrittiva(dati$Lunghezza, "Lunghezza")
risultati_cranio <- analisi_descrittiva(dati$Cranio, "Diametro Cranio")

# Combina tutti i risultati in un unico data frame
analisi_completa <- rbind(risultati_peso, risultati_anni_madre, risultati_n_gravidanze,
                          risultati_gestazione, risultati_lunghezza, risultati_cranio)

# Usa kable per creare la tabella
kable(analisi_completa)
Variabile Media Mediana DeviazioneStandard Quartile_25 Quartile_75 Asimmetria Curtosi
Peso 3284.1841473 3300 525.229374 2990 3620 -0.6474036 2.0287531
Anni.madre 28.1861489 28 5.217206 25 32 0.1510624 -0.1056061
Numero Gravidanze 0.9815853 1 1.280949 0 1 2.5134123 10.9816269
Gestazione 38.9795837 39 1.868950 38 40 -2.0651309 8.2555158
Lunghezza 494.6957566 500 26.328847 480 510 -1.5145746 6.4809304
Diametro Cranio 340.0292234 340 16.429469 330 350 -0.7850906 2.9448704
  • Le variabili mostrano diverse caratteristiche di distribuzione, con alcune che presentano skewness e curtosi significative, suggerendo la presenza di outliers o distribuzioni non normali.
  • La variabilità all’interno dei dati è ben rappresentata dalle deviazioni standard e dai quartili.
  • Le analisi di skewness e curtosi possono essere utili per ulteriori indagini statistiche e per comprendere la forma delle distribuzioni delle variabili.

Questi risultati forniscono un quadro chiaro e dettagliato delle caratteristiche statistiche delle variabili analizzate, utile per interpretazioni più approfondite e per guidare ulteriori analisi dei dati.

Variabili Categoriche

analisi_categorica <- function(x, nome_var = "Variabile") {
  # Frequenze assolute e relative
  freq_ass <- table(x)
  freq_rel <- prop.table(freq_ass)

  # Crea il dataframe per la tabella
  freq_df <- data.frame(
    Categoria = names(freq_ass),
    Fr_Assoluta = as.integer(freq_ass),
    Fr_Relativa = round(as.numeric(freq_rel), 3)
  )

  # Formatta i nomi delle colonne
  colnames(freq_df) <- c(nome_var, "Fr. Assoluta", "Fr. Relativa")
  
    # Stampa la tabella
  kable(freq_df, caption = paste("Distribuzione della variabile", nome_var), align = "lcc", escape = FALSE)

}

analisi_categorica(dati$Sesso, "Sesso")
Distribuzione della variabile Sesso
Sesso Fr. Assoluta Fr. Relativa
F 1255 0.502
M 1243 0.498 
analisi_categorica(dati$Fumatrici, "Fumatrici")
Distribuzione della variabile Fumatrici
Fumatrici Fr. Assoluta Fr. Relativa
0 2394 0.958
1 104 0.042 
analisi_categorica(dati$Tipo.parto, "Tipo_Parto")
Distribuzione della variabile Tipo_Parto
Tipo_Parto Fr. Assoluta Fr. Relativa
Ces 728 0.291
Nat 1770 0.709 
analisi_categorica(dati$Ospedale, "Ospedale")
Distribuzione della variabile Ospedale
Ospedale Fr. Assoluta Fr. Relativa
osp1 816 0.327
osp2 848 0.339
osp3 834 0.334

Quindi possiamo osservare: - Sesso: quasi perfettamente bilanciato tra Femmine (50.2%) e Maschi (49.8%). - Fumatrici: fortemente sbilanciato, con la maggior parte non fumatrice (95.8%). - Tipo.parto: prevalentemente naturale (70.9%). - Ospedale: molto uniforme tra i tre ospedali (circa 33% ciascuno).

Test di Ipotesi

# Calcolo del test
tab_parto_ospedale <- table(dati$Tipo.parto, dati$Ospedale)
chi_test <- chisq.test(tab_parto_ospedale)

Risultato Test del Chi-quadrato

  • Statistica Chi-quadrato: 1.083
  • Gradi di libertà: 2
  • p-value: 0.58188

Il test del Chi-quadrato non ha evidenziato un’associazione significativa tra il tipo di parto e l’ospedale (p = 0.578). Pertanto, non ci sono differenze statisticamente significative nella distribuzione del tipo di parto tra i diversi ospedali.

Tramite ricerca sul web e consultando gli Standard di Crescita dell’Organizzazione Mondiale della Sanità (WHO Child Growth Standards) (https://www.cdc.gov/growthcharts/who-growth-charts.htm) è stato possibile stabilire dei valori di riferimento per la media del peso e della lunghezza neonatale da utilizzare come confronto con i nostri dati. A tale scopo, abbiamo utilizzato un valore di peso medio di 3300 grammi (3.3 kg) e una lunghezza media di 500 millimetri (50 cm), che rappresentano le mediane tipiche per i neonati secondo queste fonti.

Test t per confronto con la popolazione

# Test t per Peso
pop_media_peso <- 3300
test_peso <- t.test(dati$Peso, mu = pop_media_peso, alternative = "two.sided")

# Test t per Lunghezza
pop_media_lunghezza <- 500
test_lunghezza <- t.test(dati$Lunghezza, mu = pop_media_lunghezza, alternative = "two.sided")

Peso alla nascita:

  • Media campionaria: 3284.2
  • Valore t: -1.505n
  • Intervallo di confidenza al 95%: [3263.6, 3304.8]
  • Valore p: 0.1324

Lunghezza alla nascita:

  • Media campionaria: 494.7 mm
  • Valore t: -10.069
  • Intervallo di confidenza al 95%: [493.7, 495.7]
  • Valore p: < 1e-04

La media del peso dei neonati non differisce significativamente dalla media teorica della popolazione, con p = 0.13. Non si rifiuta l’ipotesi nulla.

La lunghezza media dei neonati è significativamente inferiore rispetto alla media ipotizzata di 500 mm (p < 0.001). Si rifiuta l’ipotesi nulla.

Test t per confronto tra sessi

# Test t tra sessi
test_peso_sesso <- t.test(Peso ~ Sesso, data = dati)
test_lung_sesso <- t.test(Lunghezza ~ Sesso, data = dati)
test_cranio_sesso <- t.test(Cranio ~ Sesso, data = dati)

# Estrazione e pulizia risultati
estrai_info <- function(test) {
  list(
    media_F = round(test$estimate["mean in group F"], 1),
    media_M = round(test$estimate["mean in group M"], 1),
    t = round(test$statistic, 3),
    p = format.pval(test$p.value, digits = 4, eps = 0.0001),
    ci = round(test$conf.int, 2)
  )
}

peso <- estrai_info(test_peso_sesso)
lung <- estrai_info(test_lung_sesso)
cranio <- estrai_info(test_cranio_sesso)

Peso alla nascita:

  • Media femmine: 3161.1 g
  • Media maschi: 3408.5 g
  • Valore t: -12.115
  • Intervallo di confidenza al 95%: [-287.48, -207.38]
  • Valore p: < 1e-04

Lunghezza alla nascita:

  • Media femmine: 489.8 mm
  • Media maschi: 499.7 mm
  • Valore t: -9.582
  • Intervallo di confidenza al 95%: [-11.94, -7.88]
  • Valore p: < 1e-04

Circonferenza cranica:

  • Media femmine: 337.6 mm
  • Media maschi: 342.5 mm
  • Valore t: -7.437
  • Intervallo di confidenza al 95%: [-6.11, -3.56]
  • Valore p: < 1e-04

I maschi hanno un peso medio significativamente maggiore rispetto alle femmine (p < 0.001). Differenza media: circa 247 g.

Anche la lunghezza media è maggiore nei maschi rispetto alle femmine, con una differenza significativa (p < 0.001).

La circonferenza cranica media è maggiore nei maschi, con differenza significativa (p < 0.001).

Gestione dei Dati Mancanti

È stato verificato il numero di dati mancanti per ciascuna variabile. In questo caso, non è stato riscontrato nessun dato mancante.

kable(colSums(is.na(dati))) # Numero di dati mancanti per colonna
x
Anni.madre 0
N.gravidanze 0
Fumatrici 0
Gestazione 0
Peso 0
Lunghezza 0
Cranio 0
Tipo.parto 0
Ospedale 0
Sesso 0

2.2 Creazione del Modello di Regressione

È stato costruito un modello di regressione lineare multipla per stimare il peso del neonato.

Correlazione: È stata visualizzata la matrice di correlazione tra le variabili per identificare eventuali forti correlazioni.

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

panel.cor <- function(x, y, digits = 2, prefix = "", cex.cor, ...) {
  par(usr = c(0, 1, 0, 1))
  r <- cor(x, y)
  txt <- format(c(r, 0.123456789), digits = digits)[1]
  txt <- paste0(prefix, txt)
  if (missing(cex.cor)) cex.cor <- 0.8 / strwidth(txt)
  text(0.5, 0.5, txt, cex = cex.cor, col = "red") # Aggiunto colore rosso
}

pairs(dati_numeric, lower.panel = panel.cor, upper.panel = panel.smooth)

La matrice di correlazione permette di valutare la relazione lineare tra le variabili numeriche del dataset. Dall’analisi si nota che il peso del neonato è moderatamente correlato sia con la lunghezza che con il diametro cranico, suggerendo che neonati più lunghi e con cranio più ampio tendono ad avere un peso maggiore. Anche la durata della gestazione mostra una correlazione positiva, seppur più debole, con il peso, indicando che un maggior numero di settimane di gestazione è associato a un lieve aumento del peso alla nascita. Le relazioni tra le altre variabili numeriche risultano contenute e non segnalano particolari problemi di multicollinearità. Le variabili categoriali non sono incluse nella matrice in quanto la correlazione di Pearson si applica esclusivamente a variabili numeriche continue.

Verifica di Normalità: È stato effettuato uno Shapiro-Wilk test per verificare la normalità della variabile Peso.

shapiro.test(Peso) #non Normale
## 
##  Shapiro-Wilk normality test
## 
## data:  Peso
## W = 0.97068, p-value < 2.2e-16

Il test ha evidenziato che la variabile non è normalmente distribuita.

Costruzione e Selezione del Modello: Sono stati creati e confrontati diversi modelli utilizzando ANOVA e BIC per selezionare il modello più parsimonioso e con buona capacità predittiva. P.S. Qui ho deciso di lasciare l’output della console per semplicità di codice

mod1 <- lm(Peso ~ .-Tipo.parto -Ospedale, data = dati)
summary(mod1)
## 
## Call:
## lm(formula = Peso ~ . - Tipo.parto - Ospedale, data = dati)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1160.6  -181.3   -15.7   163.6  2630.7 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -6712.2405   141.3339 -47.492  < 2e-16 ***
## Anni.madre       0.8803     1.1491   0.766    0.444    
## N.gravidanze    11.3789     4.6767   2.433    0.015 *  
## Fumatrici      -30.3958    27.6080  -1.101    0.271    
## Gestazione      32.9472     3.8288   8.605  < 2e-16 ***
## Lunghezza       10.2316     0.3011  33.979  < 2e-16 ***
## Cranio          10.5198     0.4271  24.633  < 2e-16 ***
## SessoM          78.0787    11.2132   6.963 4.24e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 274.7 on 2490 degrees of freedom
## Multiple R-squared:  0.7272, Adjusted R-squared:  0.7264 
## F-statistic: 948.3 on 7 and 2490 DF,  p-value: < 2.2e-16
car::vif(mod1) # Verifica multicollinearità
##   Anni.madre N.gravidanze    Fumatrici   Gestazione    Lunghezza       Cranio 
##     1.189264     1.187447     1.006692     1.694331     2.079749     1.628987 
##        Sesso 
##     1.040493

Il modello iniziale mod1 è stato costruito per analizzare l’effetto di alcune variabili esplicative sul peso neonatale (Peso), escludendo fin da subito le variabili Tipo.parto e Ospedale, che non hanno valore predittivo diretto nella pratica e potrebbero introdurre variabilità non generalizzabile.

Coefficiente di Determinazione:

  • R² = 0.7272 e Adjusted R² = 0.7264 indicano un buon livello di spiegazione della variabilità del peso neonatale da parte del modello.

  • L’F-statistic è molto elevato (948.3), con un p-value < 2.2e-16, confermando che il modello nel suo complesso è altamente significativo.

Significatività delle Variabili: Dal summary del modello, risultano statisticamente significative (p < 0.05) le seguenti variabili:

  • N.gravidanze (p = 0.015)

  • Gestazione, Lunghezza, Cranio e SessoM (p << 0.001)

Queste variabili mostrano una relazione lineare importante con il peso del neonato:

  • Ogni settimana in più di gestazione è associata a un aumento medio del peso.

  • Lunghezza e circonferenza cranica sono entrambe predittori molto forti del peso alla nascita.

  • I neonati maschi tendono ad avere un peso medio superiore rispetto alle femmine.

  • Anche il numero di gravidanze precedenti risulta avere un’influenza, seppur meno marcata.

Le variabili Anni.madre e Fumatrici non risultano significative (p > 0.05), suggerendo che in questo dataset non hanno un impatto rilevante sul peso neonatale, se controlliamo per le altre covariate.

Verifica della Multicollinearità: I valori di VIF sono tutti inferiori a 2.1, indicando assenza di multicollinearità preoccupante. Questo rafforza l’affidabilità delle stime dei coefficienti.

2.3 Selezione del Modello Ottimale

mod2 <- update(mod1, ~ . -Anni.madre -Fumatrici )
summary(mod2)
## 
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio + 
##     Sesso, data = dati)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1149.37  -180.98   -15.57   163.69  2639.09 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -6681.7251   135.8036 -49.201  < 2e-16 ***
## N.gravidanze    12.4554     4.3416   2.869  0.00415 ** 
## Gestazione      32.3827     3.8008   8.520  < 2e-16 ***
## Lunghezza       10.2455     0.3008  34.059  < 2e-16 ***
## Cranio          10.5410     0.4265  24.717  < 2e-16 ***
## SessoM          77.9807    11.2111   6.956 4.47e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 274.7 on 2492 degrees of freedom
## Multiple R-squared:  0.727,  Adjusted R-squared:  0.7265 
## F-statistic:  1327 on 5 and 2492 DF,  p-value: < 2.2e-16
car::vif(mod2)
## N.gravidanze   Gestazione    Lunghezza       Cranio        Sesso 
##     1.023462     1.669779     2.075747     1.624568     1.040184
#confronto tra mod1 e mod2
anova(mod1, mod2) 
## Analysis of Variance Table
## 
## Model 1: Peso ~ (Anni.madre + N.gravidanze + Fumatrici + Gestazione + 
##     Lunghezza + Cranio + Tipo.parto + Ospedale + Sesso) - Tipo.parto - 
##     Ospedale
## Model 2: Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio + Sesso
##   Res.Df       RSS Df Sum of Sq      F Pr(>F)
## 1   2490 187905214                           
## 2   2492 188042054 -2   -136840 0.9067  0.404
BIC(mod1, mod2)
##      df      BIC
## mod1  9 35207.48
## mod2  7 35193.65

Confronto tra mod1 e mod2

Il modello mod2 è stato ottenuto semplificando mod1 tramite la rimozione delle variabili non significative Anni.madre e Fumatrici, con l’obiettivo di mantenere una buona capacità predittiva riducendo la complessità del modello.

Variabili significative in mod2:

Tutte le variabili rimaste (N.gravidanze, Gestazione, Lunghezza, Cranio, Sesso) risultano altamente significative (p < 0.01), suggerendo un’influenza rilevante sul peso neonatale.
I coefficienti mantengono segno e dimensioni simili rispetto a mod1, a conferma della robustezza del modello anche dopo la rimozione delle variabili non informative.

Bontà di adattamento:

  • R²: passa da 0.7272 (mod1) a 0.7270 (mod2)

  • R² aggiustato: da 0.7264 a 0.7265
    → La capacità esplicativa rimane sostanzialmente invariata.

  • F-statistic: pari a 1327 con p < 2.2e-16 → modello globalmente molto significativo.

Multicollinearità:

L’indice VIF per tutte le variabili in mod2 è < 2.1, ben al di sotto della soglia critica di 5, indicando assenza di collinearità tra i predittori.
La rimozione delle variabili superflue potrebbe aver contribuito a migliorare la stabilità del modello.

Confronto tra modelli:

  • ANOVA (anova(mod1, mod2))
    p-value = 0.404 → la differenza tra i due modelli non è statisticamente significativa.
    La semplificazione non comporta perdita sostanziale di informazione.

  • BIC (BIC(mod1, mod2))
    mod1: 35207.48
    mod2: 35193.65
    → Il valore più basso per mod2 suggerisce una preferenza per il modello più semplice.

Conclusione:

mod2 rappresenta una versione più parsimoniosa di mod1, senza perdita rilevante di capacità predittiva.
I criteri di selezione del modello (R², BIC, ANOVA) indicano che la rimozione di variabili non significative è giustificata.
Pertanto, mod2 risulta preferibile per interpretazione e predizione.

Analisi dei residui

L’analisi dei residui del modello mod2 è stata condotta per verificarne le assunzioni.

# Analisi dei residui
par(mfrow=c(2,2))
plot(mod2)

I grafici plot(mod2) mostrano:

Residuals vs Fitted: Si osserva una tendenza alla dispersione crescente dei residui all’aumentare dei valori fitted, indicando una potenziale eteroschedasticità, ulteriormente suggerita dalla curvatura. Sono presenti anche alcuni punti estremi.

Q-Q residuals: I residui si discostano leggermente dalla linea ideale, segnalando una possibile violazione dell’assunzione di normalità. Questo è confermato dal test di Shapiro-Wilk (W = 0.974, p < 2.2e-16), che rigetta l’ipotesi di normalità dei residui.

Scale-Location: La dispersione dei residui standardizzati rispetto ai valori fitted non è del tutto omogenea, suggerendo una lieve eteroschedasticità.

Residuals vs Leverage: Alcune osservazioni mostrano valori di leverage elevati, in particolare la numero 1551, che potrebbe influenzare in modo rilevante il modello.

shapiro.test(residuals(mod2))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(mod2)
## W = 0.97414, p-value < 2.2e-16
lmtest::bptest(mod2)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod2
## BP = 90.297, df = 5, p-value < 2.2e-16
lmtest::dwtest(mod2)
## 
##  Durbin-Watson test
## 
## data:  mod2
## DW = 1.9532, p-value = 0.1209
## alternative hypothesis: true autocorrelation is greater than 0

Il test di Shapiro–Wilk restituisce W = 0.97 e p‑value < 2.2e‑16, indicando che i residui si discostano significativamente da una distribuzione normale. La deviazione è evidente anche nella QQ‑plot, soprattutto nelle code estreme.

Il test di Breusch-Pagan (bptest, BP = 90.3, df = 5, p < 2.2e-16) indica una significativa eteroschedasticità.

Il test di Durbin-Watson (dwtest, DW = 1.9532, p = 0.1209) non fornisce evidenza significativa di autocorrelazione dei residui.

#leverage
par(mfrow=c(2,2))
lev<-hatvalues(mod2)
plot(lev)
p<-sum(lev)
n<-length(lev)
soglia=2*p/n
abline(h=soglia,col=2)
lev[lev>soglia]

#outliers
plot(rstudent(mod2))
abline(h=c(-2,2))

#distanza di cook
cook<-cooks.distance(mod2)
plot(cook,ylim = c(0,1)) 

car::outlierTest(mod2)
##       rstudent unadjusted p-value Bonferroni p
## 1551 10.046230         2.6345e-23   6.5810e-20
## 155   5.025345         5.3818e-07   1.3444e-03
## 1306  4.824963         1.4848e-06   3.7092e-03

Punti di Leva (Leverage): Il grafico dei leverage plot(lev) mostra alcuni punti con valori di leverage superiori alla soglia (2*p/n ≈ 0.004). Le osservazioni con leverage elevato sono state identificate e includono, tra le altre, 15, 155, e 161.

Outlier: Il grafico degli outliers mostra che l’osservazione 1551 supera ampiamente la soglia di ±2, configurandosi come outlier significativo, come confermato anche dal test di Bonferroni (t = 10.05, p < 2e-16). Anche le osservazioni 155 e 1306 risultano potenzialmente influenti, con p-value corretti significativi.

Distanza di Cook: Il grafico della distanza di Cook plot(cook) mostra che alcune osservazioni, in particolare la 1551, esercitano un’influenza considerevole sul modello.

Conclusioni: L’analisi dei residui suggerisce alcune violazioni delle assunzioni del modello lineare, in particolare la normalità e l’omoschedasticità. Sono stati identificati diversi punti con elevato leverage e outlier influenti (in primis l’osservazione 1551), che potrebbero influenzare le stime dei coefficienti e la validità delle inferenze del modello.

kable(dati[1551, ])
Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio Tipo.parto Ospedale Sesso
1553 30 4 0 35 4520 520 360 Nat osp2 F

È stato verificato il dato grezzo relativo all’osservazione n. 1551 e, sebbene i valori risultino biologicamente plausibili, la combinazione risulta molto rara e distante rispetto alla distribuzione generale del dataset. Questa anomalia ha generato un residuo molto elevato, indicando una forte deviazione dal modello atteso.

Considerata l’influenza rilevante di questo dato sull’adattamento del modello, con effetti distorsivi sulle stime dei coefficienti e violazioni delle assunzioni di normalità e omoschedasticità, si è deciso di escludere temporaneamente l’osservazione dalle analisi successive. Tale scelta è stata fatta per ottenere un modello più stabile e rappresentativo della popolazione generale.

dati_no1551 <- dati[rownames(dati) != "1551", ]

mod2_no1551 <- lm(Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio + Sesso, data = dati_no1551)

summary(mod2_no1551)
## 
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio + 
##     Sesso, data = dati_no1551)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1165.68  -179.74   -12.42   162.92  1410.68 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -6683.8326   133.1602 -50.194  < 2e-16 ***
## N.gravidanze    13.1448     4.2576   3.087  0.00204 ** 
## Gestazione      29.6341     3.7369   7.930 3.27e-15 ***
## Lunghezza       10.8899     0.3019  36.077  < 2e-16 ***
## Cranio           9.9192     0.4227  23.465  < 2e-16 ***
## SessoM          78.1376    10.9929   7.108 1.53e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 269.3 on 2491 degrees of freedom
## Multiple R-squared:  0.7372, Adjusted R-squared:  0.7367 
## F-statistic:  1398 on 5 and 2491 DF,  p-value: < 2.2e-16
par(mfrow=c(2,2))
plot(mod2_no1551)

shapiro.test(residuals(mod2_no1551))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(mod2_no1551)
## W = 0.98891, p-value = 5.23e-13
lmtest::bptest(mod2_no1551)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod2_no1551
## BP = 11.343, df = 5, p-value = 0.04499
lmtest::dwtest(mod2_no1551)
## 
##  Durbin-Watson test
## 
## data:  mod2_no1551
## DW = 1.9536, p-value = 0.1227
## alternative hypothesis: true autocorrelation is greater than 0
#leverage
par(mfrow=c(2,2))
lev<-hatvalues(mod2_no1551)
plot(lev)
p<-sum(lev)
n<-length(lev)
soglia=2*p/n
abline(h=soglia,col=2)
lev[lev>soglia]

#outliers
plot(rstudent(mod2_no1551))
abline(h=c(-2,2))
car::outlierTest(mod2_no1551)

#distanza di cook
cook<-cooks.distance(mod2_no1551)
plot(cook,ylim = c(0,1))

La rimozione dell’osservazione 1551 ha eliminato l’outlier più estremo e influente dal dataset, come confermato dall’assenza del punto ‘1551’ nei nuovi grafici diagnostici di Cook’s distance e nell’output di outlierTest. Questo ha portato a un modello marginalmente più stabile (lieve riduzione dell’errore standard dei residui). Tuttavia, nonostante questa pulizia del dato, il modello lineare standard mod2_no1551 continua a mostrare significative violazioni delle assunzioni fondamentali: i residui non sono distribuiti normalmente (p-value di Shapiro-Wilk << 0.05). È ancora presente eteroschedasticità (p-value di Breusch-Pagan < 0.05). Dato il persistere dell’eteroschedasticità e la non normalità dei residui nel mod2_no1551, e considerando le relazioni non lineari osservate tra Peso e variabili come Gestazione e Lunghezza, si è deciso di affinare il modello provando a trasformare la variabile x Lunghezza e Gestazione come suggerito dai grafici:

mod2_no1551_quad <- update(mod2_no1551, ~ . + I(Lunghezza^2) + I(Gestazione^2))

summary(mod2_no1551_quad)
## 
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio + 
##     Sesso + I(Lunghezza^2) + I(Gestazione^2), data = dati_no1551)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1187.6  -180.9   -12.4   164.2  1319.9 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -2.811e+03  9.024e+02  -3.115 0.001860 ** 
## N.gravidanze     1.438e+01  4.219e+00   3.408 0.000664 ***
## Gestazione       2.000e+02  6.621e+01   3.020 0.002552 ** 
## Lunghezza       -1.925e+01  4.537e+00  -4.244 2.28e-05 ***
## Cranio           1.010e+01  4.204e-01  24.011  < 2e-16 ***
## SessoM           7.339e+01  1.093e+01   6.715 2.32e-11 ***
## I(Lunghezza^2)   3.090e-02  4.621e-03   6.687 2.80e-11 ***
## I(Gestazione^2) -2.130e+00  8.678e-01  -2.454 0.014198 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 266.7 on 2489 degrees of freedom
## Multiple R-squared:  0.7426, Adjusted R-squared:  0.7419 
## F-statistic:  1026 on 7 and 2489 DF,  p-value: < 2.2e-16
BIC(mod2_no1551,mod2_no1551_quad)
##                  df      BIC
## mod2_no1551       7 35081.40
## mod2_no1551_quad  9 35045.21
par(mfrow=c(2,2))
plot(mod2_no1551_quad)

shapiro.test(residuals(mod2_no1551_quad))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(mod2_no1551_quad)
## W = 0.99048, p-value = 8.393e-12
lmtest::bptest(mod2_no1551_quad)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod2_no1551_quad
## BP = 17.125, df = 7, p-value = 0.01661
lmtest::dwtest(mod2_no1551_quad)
## 
##  Durbin-Watson test
## 
## data:  mod2_no1551_quad
## DW = 1.9496, p-value = 0.104
## alternative hypothesis: true autocorrelation is greater than 0
#leverage
par(mfrow=c(2,2))
lev<-hatvalues(mod2_no1551_quad)
plot(lev)
p<-sum(lev)
n<-length(lev)
soglia=2*p/n
abline(h=soglia,col=2)
lev[lev>soglia]

#outliers
plot(rstudent(mod2_no1551_quad))
abline(h=c(-2,2))
car::outlierTest(mod2_no1551_quad)

#distanza di cook
cook<-cooks.distance(mod2_no1551_quad)
plot(cook,ylim = c(0,1))

Il modello mod2_no1551_quad è stato sviluppato includendo termini quadratici per Lunghezza e Gestazione. Questa estensione mira a catturare potenziali relazioni non lineari tra il peso alla nascita e le variabili esplicative.

Il modello mostra una buona capacità esplicativa, con un R-quadro Multiplo pari a 0.7426 e un R-quadro Aggiustato di 0.7419, indicando che circa il 74.2% della variabilità del peso alla nascita è spiegata dai predittori inclusi. Il test F complessivo (F-statistic: 1026, p-value < 2.2e-16) risulta altamente significativo, a conferma della validità statistica del modello nel suo insieme.

L’aggiunta dei termini quadratici ha migliorato il valore del BIC (da 35081.40 a 35045.21), suggerendo che mod2_no1551_quad rappresenta una versione più parsimoniosa e meglio adattata rispetto al modello lineare senza termini quadratici.

Tutti i predittori inclusi nel modello sono risultati statisticamente significativi (p-values molto bassi). Questo evidenzia l’importanza di ciascuna variabile nel predire il peso alla nascita. La significatività dei termini quadratici I(Lunghezza^2) e I(Gestazione^2) conferma l’esistenza di relazioni non lineari.

Diagnostica dei Residui e Verifica delle Assunzioni:

L’analisi dei residui è stata condotta per valutare l’aderenza del modello alle assunzioni della regressione lineare:

Normalità dei Residui: Il Q-Q plot dei residui mostra deviazioni dalla linea teorica nelle code, e il test di Shapiro-Wilk (p-value = 8.393e-12) indica che i residui non sono normalmente distribuiti.

Omoschedasticità: I grafici “Residui vs Fitted” e “Scale-Location” mostrano un chiaro fenomeno di eteroschedasticità (aumento della dispersione dei residui all’aumentare dei valori predetti). Questa violazione è statisticamente supportata dal test di Breusch-Pagan (p-value = 0.01661).

Autocorrelazione: Il test di Durbin-Watson (p-value = 0.104) suggerisce l’assenza di autocorrelazione significativa nei residui, essendo il valore DW prossimo a 2.

Outlier e Punti Influenzanti: Nonostante la rimozione precedente di un outlier significativo, i grafici della Distanza di Cook e del Leverage indicano la presenza di alcuni punti che, pur non essendo estremi, esercitano una certa influenza sul modello.

Conclusioni sul Modello: Il modello mod2_no1551_quad è altamente predittivo per il peso alla nascita e incorpora efficacemente le relazioni non lineari. Tuttavia, le violazioni delle assunzioni di normalità e, in particolare, di omoschedasticità dei residui sono significative. Sebbene il modello sia robusto dal punto di vista descrittivo e predittivo, queste violazioni potrebbero compromettere l’affidabilità delle inferenze statistiche e degli standard error dei coefficienti.

3. Previsioni e risultati

Per effettuare una previsione pratica, si è considerato lo scenario di una neonata (Sesso Femminile) alla terza gravidanza della madre, con una gestazione di 39 settimane. Poiché i valori di Lunghezza e Diametro Cranio non erano forniti nella richiesta, sono state utilizzate le medie di queste variabili, calcolate precedentemente tramite la funzione analisi_descrittiva, per ottenere una stima del peso neonatale.

# Estrazione delle medie di lunghezza e cranio
media_lunghezza_pratica <- analisi_completa[analisi_completa$Variabile == "Lunghezza", "Media"]
media_cranio_pratica <- analisi_completa[analisi_completa$Variabile == "Diametro Cranio", "Media"]

## Creazione di un nuovo dataframe con le caratteristiche specificate
nuovi_dati_previsione_pratica_medie <- data.frame(
  N.gravidanze = 3,
  Gestazione = 39,
  Lunghezza = media_lunghezza_pratica, # Usa la media calcolata
  Cranio = media_cranio_pratica,       # Usa la media calcolata
  Sesso = factor("F", levels = c("F", "M")) # Assicurati che i livelli siano corretti
)

# Codifica della variabile Sesso per corrispondere al modello (se il modello usa SessoM)
nuovi_dati_previsione_pratica_medie$SessoM <- ifelse(nuovi_dati_previsione_pratica_medie$Sesso == "M", 1, 0)


# Previsione con il modello mod2_no1551_quad, inclusi gli intervalli di confidenza
previsione_risultati_completa <- predict(mod2_no1551_quad, newdata = nuovi_dati_previsione_pratica_medie,
                                         interval = "prediction", level = 0.95)

Risultato della Previsione: Il peso previsto per una neonata con madre alla terza gravidanza, che partorisce alla 39esima settimana, con una lunghezza media (494.7 mm) e un diametro cranico medio (340 mm) del nostro dataset, è di circa **3262.7 grammi con un intervallo di confidenza al 95% = [2739.2g, 3786.1g].

4. Visualizzazioni

ggplot(data = dati_no1551) +
  geom_point(aes(x = Gestazione,
                 y = Peso,
                 col = Fumatrici), position = "jitter") +
  geom_smooth(aes(x = Gestazione,
                  y = Peso,
                  col = Fumatrici), se = FALSE, method = "lm") +
  labs(title = "Relazione tra Durata della Gestazione e Peso alla Nascita",
       x = "Durata della Gestazione (settimane)",
       y = "Peso alla Nascita (grammi)",
       color = "Fumo") +
  theme_minimal()

ggplot(data = dati_no1551) +
  geom_point(aes(x = Gestazione,
                 y = Peso,
                 col = Sesso), position = "jitter") +
  geom_smooth(aes(x = Gestazione,
                  y = Peso,
                  col = Sesso), se = FALSE, method = "lm") +
  labs(title = "Relazione tra Durata della Gestazione e Peso alla Nascita",
       x = "Durata della Gestazione (settimane)",
       y = "Peso alla Nascita (grammi)",
       color = "Sesso") +
  theme_minimal()

ggplot(data = dati_no1551) +
  geom_point(aes(x = N.gravidanze,
                 y = Peso,
                 col = Fumatrici), position = "jitter") +
  geom_smooth(aes(x = N.gravidanze,
                  y = Peso,
                  col = Fumatrici), se = FALSE, method = "lm") +
  labs(title = "Relazione tra N.gravidanze e Peso alla Nascita",
       x = "N.gravidanze",
       y = "Peso alla Nascita (grammi)",
       color = "Fumo") +
  theme_minimal()

ggplot(data = dati_no1551) +
  geom_point(aes(x = N.gravidanze,
                 y = Peso,
                 col = Sesso), position = "jitter") +
  geom_smooth(aes(x = N.gravidanze,
                  y = Peso,
                  col = Sesso), se = FALSE, method = "lm") +
  labs(title = "Relazione tra N.gravidanze e Peso alla Nascita",
       x = "N.gravidanze",
       y = "Peso alla Nascita (grammi)",
       color = "Sesso") +
  theme_minimal()

5. Conclusioni

In questo studio, è stato sviluppato un modello di regressione lineare mod2_no1551_quad con l’obiettivo di predire il peso alla nascita. Il modello si è dimostrato altamente performante, spiegando circa il 74% della varianza del peso e identificando predittori significativi, inclusi effetti non lineari per la lunghezza e la gestazione.

Nonostante l’ottima capacità predittiva e il miglioramento del modello grazie all’introduzione di termini quadratici, l’analisi diagnostica ha rivelato la presenza di eteroschedasticità e la non-normalità dei residui, violazioni che possono influenzare l’affidabilità delle inferenze statistiche.

In futuro, per affrontare queste limitazioni e migliorare l’affidabilità del modello, sarebbe opportuno esplorare l’applicazione di un Modello Lineare Generalizzato (GLM), che permetterebbe di gestire in modo più appropriato la distribuzione e la varianza dei dati. Nonostante queste sfide metodologiche, il modello attuale fornisce già uno strumento pratico e robusto per stimare il peso alla nascita in base a diverse caratteristiche chiave.