Modello Previsione Peso Neonatale
Luca Bassem Abdul Hay
2025-11-09
Modello Previsione Peso Neonatale
Introduzione ed analisi descrittiva del dataset
Questo progetto tratterà la formulazione di un affidabile modello di previsione del peso neonatale analizzando dati provenienti da tre ospedali. In particolare, il dataset offre informazioni dettagliate riguardo: età della madre al parto, numero di gravidanze, se sono fumatrici o no, settimane di gestazione, peso del neonato, lunghezza del cranio, diametro del cranio, tipo di parto, ospedale di nascita (1, 2, 3), sesso del neonato.
Analizziamo sommariamente le variabili:
## Anni.madre N.gravidanze Fumatrici Gestazione
## Min. : 0.00 Min. : 0.0000 Min. :0.0000 Min. :25.00
## 1st Qu.:25.00 1st Qu.: 0.0000 1st Qu.:0.0000 1st Qu.:38.00
## Median :28.00 Median : 1.0000 Median :0.0000 Median :39.00
## Mean :28.16 Mean : 0.9812 Mean :0.0416 Mean :38.98
## 3rd Qu.:32.00 3rd Qu.: 1.0000 3rd Qu.:0.0000 3rd Qu.:40.00
## Max. :46.00 Max. :12.0000 Max. :1.0000 Max. :43.00
## Peso Lunghezza Cranio Tipo.parto
## Min. : 830 Min. :310.0 Min. :235 Length:2500
## 1st Qu.:2990 1st Qu.:480.0 1st Qu.:330 Class :character
## Median :3300 Median :500.0 Median :340 Mode :character
## Mean :3284 Mean :494.7 Mean :340
## 3rd Qu.:3620 3rd Qu.:510.0 3rd Qu.:350
## Max. :4930 Max. :565.0 Max. :390
## Ospedale Sesso
## Length:2500 Length:2500
## Class :character Class :character
## Mode :character Mode :character
##
##
## Frequenze relative delle variabili qualitative:
## Ospedale
## osp1 osp2 osp3
## 0.3264 0.3396 0.3340## Fumatrici
## 0 1
## 0.9584 0.0416## Tipo.parto
## Ces Nat
## 0.2912 0.7088## Sesso
## F M
## 0.5024 0.4976Ad uno sguardo iniziale, i dati sembrano piuttosto completi e ben
distribuiti: i parti provengono in proporzione simile dai tre diversi
ospedali, i parti coinvolgono bambini maschi e bambine femmine in egual
misura, come è da aspettarsi. Noto che, tuttavia, un paio di elementi
che mi paiono degni di attenzione:
- l’età minima della madre è 0;
- la variabile N.gravidanze ha un valore massimo di 12 che potrebbe
rappresentare un outlier.
L’analisi del boxplot di Anni.madre conferma i miei sospetti. Le osservazioni 1152 e 1380 paiono francamente errori e posteriormente proverò ad eliminarle dall’analisi e a testarne gli effetti.
## Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio
## 1152 1 1 0 41 3250 490 350
## 1380 0 0 0 39 3060 490 330
## Tipo.parto Ospedale Sesso
## 1152 Nat osp2 F
## 1380 Nat osp3 MTest delle ipotesi
Testiamo, quindi, le seguenti ipotesi come ulteriore controllo della
affidabilità del campione:
1. In alcuni ospedali si praticano più
parti cesarei.
2. La media del peso e della lunghezza di questo
campione di neonati sono significativamente uguali a quelle della
popolazione
3. Le misure antropometriche sono significativamente
diverse tra i due sessi
Test 1
1. In alcuni ospedali si praticano più parti cesarei.
Per testare questa possibilità, ho bisogno di definire il tipo di test che dipende dal tipo di misura che devo confrontare. In questo caso, voglio saggiare l’ipotesi che il numero (o meglio la proporzione) di parti cesarei dipenda in qualche modo dall’ospedale in cui il parto è avvenuto. Per questo tipo di confronti, si può usare il test del Chi quadrato che verifica l’eventuale associazione di una mutabile (o variabile qualitative, o quantitativa divisa in classi) con un’altra mutabile.
Nel caso di specie, andremo a verificare se l’ospedale è associato ad un maggiore/minore numero di parti cesarei nel campione in esame. Per fare ciò definisco l’ipotesi nulla del test come “la percentuale di parti cesarei è indipendente dall’ospedale dove il parto è avvenuto”.
Innanzitutto costruiamo una tabella di contingenza tra i vari ospedali e il tipo di parto:
## Tipo.parto
## Ospedale Ces Nat
## osp1 242 574
## osp2 254 595
## osp3 232 603Procedo quindi ad eseguire un test Chi-quadrato su questa tabella di contingenza:
##
## Pearson's Chi-squared test
##
## data: table(Ospedale, Tipo.parto)
## X-squared = 1.0972, df = 2, p-value = 0.5778Il p-value ha un valore di 0.58 che non mi permette di rifiutare l’ipotesi nulla.
Posso concludere che nel campione non ci sono elementi che dimostrano una maggiore propensione di alcuni ospedali a praticare parti cesarei.
Test 2
2. La media del peso e della lunghezza di questo campione di neonati sono significativamente uguali a quelle della popolazione
Per poter testare l’enunciato di questo quesito, avrei bisogno di avere un valore di riferimento per la popolazione. Cercando online, ho trovato che uno studio sardo con dati del 2010 che riporta un peso medio alla nascita di 3147 g e una lunghezza media di 49.3cm. Assumendo che tali valori non siano significativamente diversi dai dati a livello nazionale (che posso considerare come la nostra popolazione di riferimento), li adopererò come parametri di confronto per il peso medio campionario e la lunghezza media campionaria.
A questo proposito, eseguirò un t test per saggiare l’ipotesi che peso e lunghezza del mio campione non siano diversi da quelli della popolazione di riferimento.
Test peso:
##
## One Sample t-test
##
## data: Peso
## t = 13.054, df = 2499, p-value < 2.2e-16
## alternative hypothesis: true mean is not equal to 3147
## 95 percent confidence interval:
## 3263.490 3304.672
## sample estimates:
## mean of x
## 3284.081Test lunghezza:
##
## One Sample t-test
##
## data: Lunghezza
## t = 3.2145, df = 2499, p-value = 0.001324
## alternative hypothesis: true mean is not equal to 493
## 95 percent confidence interval:
## 493.6598 495.7242
## sample estimates:
## mean of x
## 494.692I due test t eseguiti indicano che entrambi i valori campionari di peso e lunghezza sono significativamente diversi da quelli da me riportati per la popolazione generale ad un livello di significatività di <0.01. Ciò detto, pur non essendo un esperto di neonatologia, esiste una differenza di magnitudine tra lo scostamento del peso campionario rispetto al peso della popolazione e quello della lunghezza campionaria rispetto alla lunghezza della popolazione.
Infatti, il peso medio campionario è ca. 140g maggiore rispetto a quello della popolazione, con un intervallo di confidenza al 95% che va dai 3263g ai 3305g, indicando quindi una differenza di magnitudine oltre che statisticamente significativa. Al contrario, la lunghezza campionaria, pur significativamente inferiore, ha uno scostamento di circa 1 mm rispetto al valore della popolazione. Tale differenza potrebbe essere trascurabile.
Test 3
Le misure antropometriche sono significativamente diverse tra i due sessi
Il terzo quesito invita a saggiare le differenze antropometriche alla
nascita tra neonati maschi e neonate femmine. A questo proposito,
testerò se esiste una differenza significativa tra i due gruppi per le
seguenti variabili quantitative:
- Peso
- Lunghezza del
neonato
- Diametro del cranio
In questa fase vengono deliberatamente ignorati eventuali fattori confondenti (ad esempio settimane di gestazione, età materna, abitudini di fumo), che verranno invece considerati nella modellizzazione successiva.
Innanzitutto, procedo ad una rapida analisi delle distribuzioni delle tre variabili in esame per scegliere il test più adatto:
Prima di scegliere il test statistico, è stata eseguita una rapida analisi esplorativa delle distribuzioni delle variabili. Poiché il campione è numeroso, si applica il t-test di Welch, predefinito in R, che non assume omogeneità delle varianze. Anche qualora le varianze fossero simili, l’uso del t-test di Welch risulterebbe comunque appropriato e non comporterebbe conseguenze rilevanti nelle inferenze.
##
## Welch Two Sample t-test
##
## data: Lunghezza by Sesso
## t = -9.582, df = 2459.3, p-value < 2.2e-16
## alternative hypothesis: true difference in means between group F and group M is not equal to 0
## 95 percent confidence interval:
## -11.929470 -7.876273
## sample estimates:
## mean in group F mean in group M
## 489.7643 499.6672Test peso
##
## Welch Two Sample t-test
##
## data: Peso by Sesso
## t = -12.106, df = 2490.7, p-value < 2.2e-16
## alternative hypothesis: true difference in means between group F and group M is not equal to 0
## 95 percent confidence interval:
## -287.1051 -207.0615
## sample estimates:
## mean in group F mean in group M
## 3161.132 3408.215Test Diametro Cranio
##
## Welch Two Sample t-test
##
## data: Cranio by Sesso
## t = -7.4102, df = 2491.4, p-value = 1.718e-13
## alternative hypothesis: true difference in means between group F and group M is not equal to 0
## 95 percent confidence interval:
## -6.089912 -3.541270
## sample estimates:
## mean in group F mean in group M
## 337.6330 342.4486Per tutte e tre le misure antropometriche presenti nel dataset, esiste una differenza significativa tra il gruppo dei maschi e quello delle femmine.
Creazione del modello di regressione
Il modello di regressione multipla aiuterà a definire quali sono le variabili che influenzano maggiormente il peso dei neonati alla nascita.
Come primo esperimento, svilupperò un modello che contenga tutte le variabili contenute nel dataset. Procederò poi a ottimizzare il modello, analizzando l’eventuale trade-off tra semplicità ed accuratezza con gli indici preposti a tal fine ( \(R^2\), \(RMSE\), \(AIC\), \(BIC\)).
Ecco il modello completo:
##
## Call:
## lm(formula = Peso ~ Anni.madre + N.gravidanze + Fumatrici + Gestazione +
## Lunghezza + Cranio + Tipo.parto + Ospedale + Sesso)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1124.40 -181.66 -14.42 160.91 2611.89
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6738.4762 141.3087 -47.686 < 2e-16 ***
## Anni.madre 0.8921 1.1323 0.788 0.4308
## N.gravidanze 11.2665 4.6608 2.417 0.0157 *
## Fumatrici -30.1631 27.5386 -1.095 0.2735
## Gestazione 32.5696 3.8187 8.529 < 2e-16 ***
## Lunghezza 10.2945 0.3007 34.236 < 2e-16 ***
## Cranio 10.4707 0.4260 24.578 < 2e-16 ***
## Tipo.partoNat 29.5254 12.0844 2.443 0.0146 *
## Ospedaleosp2 -11.2095 13.4379 -0.834 0.4043
## Ospedaleosp3 28.0958 13.4957 2.082 0.0375 *
## SessoM 77.5409 11.1776 6.937 5.08e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 273.9 on 2489 degrees of freedom
## Multiple R-squared: 0.7289, Adjusted R-squared: 0.7278
## F-statistic: 669.2 on 10 and 2489 DF, p-value: < 2.2e-16Questo primo sviluppo del modello suggerisce che potremmo essere di fronte ad un caso di overfit. Infatti, diverse variabili non sembrano avere un effetto significativo sul peso dei neonati (es. età della madre, madre fumatrice, ospedale di nascita). Un altro possibile problema è quello della multicollinearità ovvero un’elevata correlazione tra regressori. A questo proposito, adopererò la funzione VIF (Variance Inflation Factor) del pacchetto car di R:
## Caricamento del pacchetto richiesto: carData## GVIF Df GVIF^(1/(2*Df))
## Anni.madre 1.187454 1 1.089704
## N.gravidanze 1.186428 1 1.089233
## Fumatrici 1.007392 1 1.003689
## Gestazione 1.695810 1 1.302233
## Lunghezza 2.085755 1 1.444214
## Cranio 1.630796 1 1.277026
## Tipo.parto 1.004242 1 1.002119
## Ospedale 1.004071 2 1.001016
## Sesso 1.040643 1 1.020119Come forse non ci si aspettava, il VIF del modello sembra sotto controllo (nessun regressore si avvicina la soglia limite di 5). Si può quindi concludere che tale modello non contiene variabili altamente correlate tra loro che potrebbero creare problemi di interpretabilità del modello stesso.
Selezione del modello ottimale
Il modello completo non è necessariamente il migliore. Si è infatti sottolineato come alcune variabili non abbiano un effetto significativo sulla variabile obiettivo e potrebbero pertanto essere rimosse. Esistono, però, dei metodi di valutazione della bontà di un modello che sono estrememamente utili a discriminare tra modelli concorrenti. Procederò, quindi, ad usare una procedura stepwise per provare varie versioni del modello che poi verranno confrontate usando il criterio AIC (Akaike Information Criterion) o il criterio BIC (Bayesian Information Criterion).
AIC:
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio +
## Tipo.parto + Ospedale + Sesso)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1113.18 -181.16 -16.58 161.01 2620.19
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6707.4293 135.9438 -49.340 < 2e-16 ***
## N.gravidanze 12.3619 4.3325 2.853 0.00436 **
## Gestazione 31.9909 3.7896 8.442 < 2e-16 ***
## Lunghezza 10.3086 0.3004 34.316 < 2e-16 ***
## Cranio 10.4922 0.4254 24.661 < 2e-16 ***
## Tipo.partoNat 29.2803 12.0817 2.424 0.01544 *
## Ospedaleosp2 -11.0227 13.4363 -0.820 0.41209
## Ospedaleosp3 28.6408 13.4886 2.123 0.03382 *
## SessoM 77.4412 11.1756 6.930 5.36e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 273.9 on 2491 degrees of freedom
## Multiple R-squared: 0.7287, Adjusted R-squared: 0.7278
## F-statistic: 836.3 on 8 and 2491 DF, p-value: < 2.2e-16BIC:
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio +
## Sesso)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1149.44 -180.81 -15.58 163.64 2639.72
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6681.1445 135.7229 -49.226 < 2e-16 ***
## N.gravidanze 12.4750 4.3396 2.875 0.00408 **
## Gestazione 32.3321 3.7980 8.513 < 2e-16 ***
## Lunghezza 10.2486 0.3006 34.090 < 2e-16 ***
## Cranio 10.5402 0.4262 24.728 < 2e-16 ***
## SessoM 77.9927 11.2021 6.962 4.26e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.6 on 2494 degrees of freedom
## Multiple R-squared: 0.727, Adjusted R-squared: 0.7265
## F-statistic: 1328 on 5 and 2494 DF, p-value: < 2.2e-16Il confronto tra i due modelli ottimizzati, rispettivamente, secondo il criterio AIC e quello BIC evidenzia delle differenze sostanziali. In particolare, la selezione del modello tramite stepwise BIC ha evidenziato che le variabili significativamente predittive del peso alla nascita sono il numero di gravidanze precedenti, la durata della gestazione, la lunghezza, il diametro craniale del neonato e il sesso. Il modello spiega circa il 72.6% della variabilità del peso (\(R^2\) aggiustato = 0.7265). Variabili come l’ospedale e il tipo di parto, pur incluse nel modello AIC, non apportano un contributo significativo e sono state eliminate dal modello BIC, più parsimonioso e interpretabile.
Interazioni e relazioni non lineari
Selezionato il modello ottimale secondo i criteri scelti, procedo a testare interazioni tra le variabili e relazioni non lineari. Un’interazione che potrebbe avere senso è quella tra gestazione (in settimane) e fumo materno (ovvero che una maggiore gestazione in donne fumatrici potrebbe avere un effetto minore rispetto alle non-fumatrici):
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione * Fumatrici + Lunghezza +
## Cranio + Sesso, data = neonati)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1149.9 -181.1 -17.1 163.6 2636.3
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6699.0708 136.6478 -49.024 < 2e-16 ***
## N.gravidanze 12.7641 4.3451 2.938 0.00334 **
## Gestazione 33.1472 3.8389 8.635 < 2e-16 ***
## Fumatrici 794.5870 757.2739 1.049 0.29415
## Lunghezza 10.2285 0.3009 33.988 < 2e-16 ***
## Cranio 10.5305 0.4263 24.704 < 2e-16 ***
## SessoM 78.7548 11.2151 7.022 2.81e-12 ***
## Gestazione:Fumatrici -21.0157 19.2765 -1.090 0.27572
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.6 on 2492 degrees of freedom
## Multiple R-squared: 0.7273, Adjusted R-squared: 0.7265
## F-statistic: 949.3 on 7 and 2492 DF, p-value: < 2.2e-16Tale effetto dell’interazione tra fumo e settimane di gestazione sul peso del neonato pare non esistere.
Un’altra interazione sensata potrebbe essere tra la Lunghezza e Sesso (ovvero l’aumento di peso per unità di lunghezza potrebbe essere diverso tra maschi e femmine)
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Fumatrici + Lunghezza *
## Sesso + Cranio, data = neonati)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1162.99 -179.80 -15.61 163.00 2570.05
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6493.7850 162.8062 -39.887 < 2e-16 ***
## N.gravidanze 12.9883 4.3441 2.990 0.00282 **
## Gestazione 32.8853 3.8051 8.642 < 2e-16 ***
## Fumatrici -31.7322 27.5830 -1.150 0.25008
## Lunghezza 9.8476 0.3532 27.882 < 2e-16 ***
## SessoM -365.8205 213.1202 -1.716 0.08620 .
## Cranio 10.5056 0.4262 24.649 < 2e-16 ***
## Lunghezza:SessoM 0.8962 0.4296 2.086 0.03706 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.4 on 2492 degrees of freedom
## Multiple R-squared: 0.7276, Adjusted R-squared: 0.7268
## F-statistic: 951 on 7 and 2492 DF, p-value: < 2.2e-16In questo caso, invece, il modello con l’interazione pare avere catturato un diverso effetto della lunghezza sul peso in funzione del sesso. L’interazione, infatti, è risultata significativa (\(β\) = 0.896, p = 0.037), indicando che per ogni unità di lunghezza in più, il peso dei maschi aumenta leggermente più delle femmine, confermando un effetto combinato di sesso e lunghezza.
Per quanto riguarda eventuali effetti non lineari, un candidato potrebbe essere la variabile Gestazione. È infatti plausibile che all’aumentare delle settimane di gestazione la crescita del neonato vari.
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) +
## Fumatrici + Lunghezza + Sesso + Cranio, data = neonati)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1145.00 -180.96 -13.12 165.00 2659.10
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4669.1027 898.3246 -5.198 2.18e-07 ***
## N.gravidanze 12.7990 4.3416 2.948 0.00323 **
## Gestazione -79.7741 49.7260 -1.604 0.10878
## I(Gestazione^2) 1.4999 0.6618 2.266 0.02352 *
## Fumatrici -29.2442 27.5772 -1.060 0.28904
## Lunghezza 10.3386 0.3042 33.989 < 2e-16 ***
## SessoM 75.9814 11.2351 6.763 1.68e-11 ***
## Cranio 10.6312 0.4280 24.841 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.4 on 2492 degrees of freedom
## Multiple R-squared: 0.7277, Adjusted R-squared: 0.7269
## F-statistic: 951.4 on 7 and 2492 DF, p-value: < 2.2e-16Includendo un termine quadratico per la gestazione nel modello lineare si cattura la non linearità della crescita del peso alla nascita. Il termine lineare è negativo (\(𝛽_1\) = - 79.77) mentre quello quadratico è positivo (\(𝛽_2\) = 1.50) indicando che l’effetto di ogni settimana aggiuntiva di gestazione aumenta con l’aumentare della durata della gravidanza. L’effetto marginale è quindi funzione della settimana di gestazione, confermando un aumento progressivo del peso per nascite più prossime al termine.
Se applico una procedura stepwise al modello con il fattore quadratico:
stepwise_BIC_quad <- MASS::stepAIC(mod_quad,
direction = "both",
k = log(nrow(neonati)))## Start: AIC=28126.86
## Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) + Fumatrici +
## Lunghezza + Sesso + Cranio
##
## Df Sum of Sq RSS AIC
## - Fumatrici 1 84652 187671672 28120
## - Gestazione 1 193737 187780757 28122
## - I(Gestazione^2) 1 386634 187973654 28124
## <none> 187587020 28127
## - N.gravidanze 1 654199 188241219 28128
## - Sesso 1 3442817 191029837 28165
## - Cranio 1 46451021 234038041 28672
## - Lunghezza 1 86964157 274551177 29071
##
## Step: AIC=28120.16
## Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) + Lunghezza +
## Sesso + Cranio
##
## Df Sum of Sq RSS AIC
## - Gestazione 1 200094 187871765 28115
## - I(Gestazione^2) 1 393875 188065546 28118
## <none> 187671672 28120
## - N.gravidanze 1 632219 188303891 28121
## + Fumatrici 1 84652 187587020 28127
## - Sesso 1 3426378 191098049 28158
## - Cranio 1 46500703 234172375 28666
## - Lunghezza 1 87399773 275071444 29068
##
## Step: AIC=28115
## Peso ~ N.gravidanze + I(Gestazione^2) + Lunghezza + Sesso + Cranio
##
## Df Sum of Sq RSS AIC
## <none> 187871765 28115
## - N.gravidanze 1 630973 188502738 28116
## + Gestazione 1 200094 187671672 28120
## + Fumatrici 1 91008 187780757 28122
## - Sesso 1 3592323 191464088 28155
## - I(Gestazione^2) 1 5658634 193530399 28181
## - Cranio 1 46389582 234261347 28659
## - Lunghezza 1 88923843 276795608 29076summary(stepwise_BIC_quad)##
## Call:
## lm(formula = Peso ~ N.gravidanze + I(Gestazione^2) + Lunghezza +
## Sesso + Cranio, data = neonati)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1146.87 -181.09 -15.12 165.40 2643.89
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.101e+03 1.235e+02 -49.380 < 2e-16 ***
## N.gravidanze 1.255e+01 4.338e+00 2.894 0.00383 **
## I(Gestazione^2) 4.379e-01 5.053e-02 8.667 < 2e-16 ***
## Lunghezza 1.026e+01 2.987e-01 34.358 < 2e-16 ***
## SessoM 7.733e+01 1.120e+01 6.906 6.32e-12 ***
## Cranio 1.056e+01 4.255e-01 24.816 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.5 on 2494 degrees of freedom
## Multiple R-squared: 0.7273, Adjusted R-squared: 0.7267
## F-statistic: 1330 on 5 and 2494 DF, p-value: < 2.2e-16Ottengo un modello senza la variabile Fumatrici e la variabile Gestazione base che seppur corretto da un punto di vista procedurale, rende il modello meno interpretabile. Scelgo quindi un modello senza la variabile Fumatrici ma includendo la variabile Gestazione.
mod_final <- lm(Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) +
Lunghezza + Sesso + Cranio, data = neonati)
summary(mod_final)##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) +
## Lunghezza + Sesso + Cranio, data = neonati)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1144.11 -181.17 -13.16 165.73 2662.56
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4650.2268 898.1706 -5.177 2.43e-07 ***
## N.gravidanze 12.5660 4.3361 2.898 0.00379 **
## Gestazione -81.0486 49.7127 -1.630 0.10316
## I(Gestazione^2) 1.5136 0.6617 2.287 0.02226 *
## Lunghezza 10.3534 0.3039 34.074 < 2e-16 ***
## SessoM 75.7900 11.2340 6.747 1.88e-11 ***
## Cranio 10.6363 0.4280 24.854 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274.4 on 2493 degrees of freedom
## Multiple R-squared: 0.7276, Adjusted R-squared: 0.7269
## F-statistic: 1110 on 6 and 2493 DF, p-value: < 2.2e-16Come ulteriore analisi, proietto su un grafico la relazione tra peso e gestazione e pare, in effetti, una relazione non lineare.
Decido, quindi, di includere nel modello finale il termine quadratico per la gestazione (\(R^2\) aggiustato = 0.7269).
Analisi della qualità del modello
Per validare definitivamente il modello prescelto, è necessario fare un’analisi accurata della capacità predittiva del modello. Ho già introdotto l’ \(R^2\), in questa fase si calcolerà anche il \(RMSE\)
## [1] 273.9866Questo valore indica qual è lo scostamento medio della previsione del modello rispetto al valore osservato. Chiaramente, più basso è, migliore è la capacità predittiva del modello.
Diagnostica dei residui
Ancora più importante, però, è la diagnostica dei residui per identificare eventuali osservazioni limite che potrebbero distorcere le previsioni.
Da questa analisi emerge un’osservazione in particolare, la numero 1551, che pare essere problematica.
Se approfondisco l’analisi sui valori di leva e gli outliers:
Valori di leva:
## [1] 118Ci sono ben 118 osservazioni con un valore di leva oltre la soglia.
Outliers:
Dal grafico, si evince che ci sono varie osservazioni che hanno dei valori outliers per la variabile esplicativa. Procedo quindi con un outlier.test per ottenere delle stime più precise.
## rstudent unadjusted p-value Bonferroni p
## 1551 10.158855 8.7351e-24 2.1838e-20
## 155 5.095784 3.7333e-07 9.3334e-04
## 1306 4.780859 1.8471e-06 4.6177e-03Le osservazioni effettivamente identificate come outliers sono: 1551, 155 e 1306.
La verifica finale consiste nel calcolo della distanza di Cook:
## 1551
## 0.7468715L’unica osservazione oltre la soglia limite di Cook (0.5) è la famigerata osservazione 1551.
Procedo con i test di normalità, omoschedasticità e …
Il test Breusch-Pagan di omoschedasticità:
## Caricamento del pacchetto richiesto: zoo##
## Caricamento pacchetto: 'zoo'## I seguenti oggetti sono mascherati da 'package:base':
##
## as.Date, as.Date.numeric##
## studentized Breusch-Pagan test
##
## data: mod_final
## BP = 95.737, df = 6, p-value < 2.2e-16Il test Durbin-Watson di autocorrelazione dei residui:
##
## Durbin-Watson test
##
## data: mod_final
## DW = 1.9524, p-value = 0.1171
## alternative hypothesis: true autocorrelation is greater than 0Il test Shapiro-Wilk di normalità dei residui:
##
## Shapiro-Wilk normality test
##
## data: residuals(mod_final)
## W = 0.97362, p-value < 2.2e-16
Proviamo ora a rimuovere quest’osservazione dal dataset e a ricalcolare
il modello:
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + I(Gestazione^2) +
## Lunghezza + Sesso + Cranio, data = neonati_corrected)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1159.71 -176.98 -11.72 165.54 1430.66
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4322.4420 880.8986 -4.907 9.85e-07 ***
## N.gravidanze 13.2781 4.2505 3.124 0.00181 **
## Gestazione -102.2464 48.7686 -2.097 0.03613 *
## I(Gestazione^2) 1.7596 0.6490 2.711 0.00675 **
## Lunghezza 11.0211 0.3050 36.137 < 2e-16 ***
## SessoM 75.5755 11.0105 6.864 8.43e-12 ***
## Cranio 10.0242 0.4237 23.656 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 268.9 on 2492 degrees of freedom
## Multiple R-squared: 0.738, Adjusted R-squared: 0.7373
## F-statistic: 1170 on 6 and 2492 DF, p-value: < 2.2e-16Eseguo nuovamente la diagnostica dei residui:
Il test Breusch-Pagan di omoschedasticità:
##
## studentized Breusch-Pagan test
##
## data: mod_final_corr
## BP = 37.686, df = 6, p-value = 1.294e-06Il test Shapiro-Wilk di normalità dei residui:
##
## Shapiro-Wilk normality test
##
## data: residuals(mod_final_corr)
## W = 0.98897, p-value = 5.711e-13
L’osservazione empirica della distribuzione dei residui non sembra evidenziare delle problematiche molto gravi, se non una coda destra più lunga.
Previsioni e risultati
Validato il modello, pur consapevole delle sue potenziali limitazioni, si procede ad una stima inferenziale per il peso di una bambina, nata alla 39esima settimana per una madre che è alla sua terza gravidanza. Il modello che ho scelto prevede anche l’uso delle variabili cranio e lunghezza.
Per procedere, quindi, considererò la media campionaria per entrambe le variabili. Stimerò inoltre un intervallo di previsione che darà una stima del range attorno al valore puntuale emesso dal modello.
## fit lwr upr
## 1 3267.778 2729.244 3806.311Visualizzazioni
Concluderò con delle visualizzazioni in grado di offire un’immagine intuitiva delle principali relazioni tra variabili e risultati del modello. Ad esempio, in questo grafico si può vedere la relazione tra peso reale e peso predetto dal modello:
I punti si concentrano lungo la bisettrice, indicando una discreta qualità del modello, pur con qualche dispersione dovuta alla variabilità naturale del fenomeno. Si nota, forse, una leggera sottostima del modello nella coda sinistra (valori bassi di peso osservato e previsto).
A questo proposito, ripropongo il grafico dei residui e dei valori previsti:
Anche da questo grafico, si nota un pattern: i residui tendono ad essere positivi per neonati con peso predetto basso, suggerendo che il modello sottostima i pesi più bassi. Questo comportamento è coerente con l’esito del test di Breusch-Pagan, che indica eteroschedasticità residua. È, inoltre, possibile che fenomeni non lineari o interazioni non incluse nel modello contribuiscano a questo pattern.