# -----------------------------------------------------------------------------
# LIBRERIE
# -----------------------------------------------------------------------------
library(tidyverse)
library(car)
library(MASS)
library(broom)
library(lmtest)
library(caret)
library(GGally)
# -----------------------------------------------------------------------------
# 1. CARICAMENTO E PREPARAZIONE DEI DATI
# -----------------------------------------------------------------------------
df <- read.csv("neonati.csv")
# Conversione variabili qualitative in factor
df$Fumatrici <- factor(df$Fumatrici, levels = c(0, 1),
labels = c("Non Fumatrice", "Fumatrice"))
df$Tipo.parto <- as.factor(df$Tipo.parto)
df$Ospedale <- as.factor(df$Ospedale)
df$Sesso <- as.factor(df$Sesso)
# Ispezione dataset
str(df)## '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 : Factor w/ 2 levels "Non Fumatrice",..: 1 1 1 1 1 1 1 1 1 1 ...
## $ 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 ...
## Anni.madre N.gravidanze Fumatrici Gestazione
## Min. : 0.00 Min. : 0.0000 Non Fumatrice:2396 Min. :25.00
## 1st Qu.:25.00 1st Qu.: 0.0000 Fumatrice : 104 1st Qu.:38.00
## Median :28.00 Median : 1.0000 Median :39.00
## Mean :28.16 Mean : 0.9812 Mean :38.98
## 3rd Qu.:32.00 3rd Qu.: 1.0000 3rd Qu.:40.00
## Max. :46.00 Max. :12.0000 Max. :43.00
## Peso Lunghezza Cranio Tipo.parto Ospedale Sesso
## Min. : 830 Min. :310.0 Min. :235 Ces: 728 osp1:816 F:1256
## 1st Qu.:2990 1st Qu.:480.0 1st Qu.:330 Nat:1772 osp2:849 M:1244
## Median :3300 Median :500.0 Median :340 osp3:835
## Mean :3284 Mean :494.7 Mean :340
## 3rd Qu.:3620 3rd Qu.:510.0 3rd Qu.:350
## Max. :4930 Max. :565.0 Max. :390
Considerazioni:
# --- 2.1 Distribuzione del peso e ricerca outlier ---
ggplot(df, aes(y = Peso)) +
geom_boxplot(fill = "lightblue") +
theme_minimal() +
labs(title = "Distribuzione del Peso Neonatale", y = "Peso (grammi)")Possiamo notare che ci sono più outlier nella parte sinistra della distribuzione che nella parte destra.
# --- 2.2 Matrice di correlazione tra variabili numeriche ---
# Permette di anticipare problemi di multicollinearità
df |>
dplyr::select(Peso, Lunghezza, Cranio, Gestazione, Anni.madre, N.gravidanze) |>
ggpairs(
title = "Matrice di Correlazione - Variabili Numeriche",
upper = list(continuous = wrap("cor", size = 3)),
lower = list(continuous = wrap("points", alpha = 0.2, size = 0.5)),
diag = list(continuous = wrap("densityDiag", fill = "lightblue"))
)Nella matrice di correlazione possiamo notare che: Il peso, ovvero la nostra variabile target, ha una correlazione elevata con la lunghezza e il diametro del cranio del neonato e la gestazione. Un’altra correlazione interessante è quella tra la gestazione e la lunghezza del cranio. Da tenerne conto per un’eventuale multicollinearità.
# --- 2.3 Distribuzione del peso per variabili categoriali ---
ggplot(df, aes(x = Sesso, y = Peso, fill = Sesso)) +
geom_boxplot() +
theme_minimal() +
labs(title = "Peso per Sesso")ggplot(df, aes(x = Fumatrici, y = Peso, fill = Fumatrici)) +
geom_boxplot() +
scale_fill_manual(values = c("Non Fumatrice" = "#2c3e50", "Fumatrice" = "#e74c3c")) +
theme_minimal() +
labs(title = "Peso per Abitudine al Fumo")ggplot(df, aes(x = Ospedale, y = Peso, fill = Ospedale)) +
geom_boxplot() +
theme_minimal() +
labs(title = "Peso per Ospedale")# --- 3.1 Chi-quadro: parti cesarei diversi tra ospedali? ---
tabella_parti <- table(df$Ospedale, df$Tipo.parto)
print(tabella_parti)##
## Ces Nat
## osp1 242 574
## osp2 254 594
## osp3 232 602
# verifica delle frequenze attese affinché siano >= 5
# (assunzione del test chi-quadro)
cat("\n--- Frequenze attese (devono essere tutte >= 5) ---\n")##
## --- Frequenze attese (devono essere tutte >= 5) ---
##
## Ces Nat
## osp1 237.8094 578.1906
## osp2 247.1353 600.8647
## osp3 243.0552 590.9448
# Se tutte >= 5, possiamo usare il chi-quadro; altrimenti usiamo il test di Fisher
risultato_chi <- chisq.test(tabella_parti)
print(risultato_chi)##
## Pearson's Chi-squared test
##
## data: tabella_parti
## X-squared = 1.083, df = 2, p-value = 0.5819
COMMENTO: Abbiamo p-value > 0.05 che ci porta a non rifiutare l’ipotesi nulla H0. I tre ospedali hanno proporzioni di cesarei statisticamente equivalenti.
# Visualizzazione proporzioni per ospedale
df |>
count(Ospedale, Tipo.parto) |>
group_by(Ospedale) |>
mutate(prop = n / sum(n)) |>
ggplot(aes(x = Ospedale, y = prop, fill = Tipo.parto)) +
geom_col(position = "dodge") +
scale_y_continuous(labels = scales::percent) +
theme_minimal() +
labs(title = "Proporzione Tipi di Parto per Ospedale",
y = "Proporzione", fill = "Tipo di Parto")# --- 3.2 T-test a un campione: confronto con valori di riferimento ---
# Valori di riferimento clinici letteratura scientifica: Peso = 3300g, Lunghezza = 500mm (50cm)
# H0: la media campionaria è uguale al valore di riferimento
cat("\n--- T-test: Peso (mu = 3300g) ---\n")##
## --- T-test: Peso (mu = 3300g) ---
##
## One Sample t-test
##
## data: df$Peso
## t = -1.505, df = 2497, p-value = 0.1324
## alternative hypothesis: true mean is not equal to 3300
## 95 percent confidence interval:
## 3263.577 3304.791
## sample estimates:
## mean of x
## 3284.184
COMMENTO: Non significativo. La media campionaria di 3284g è compatibile con il riferimento clinico di 3300g. Il campione è rappresentativo della popolazione di riferimento per questa misura.
##
## --- T-test: Lunghezza (mu = 500mm) ---
##
## One Sample t-test
##
## data: df$Lunghezza
## t = -10.069, df = 2497, p-value < 2.2e-16
## alternative hypothesis: true mean is not equal to 500
## 95 percent confidence interval:
## 493.6628 495.7287
## sample estimates:
## mean of x
## 494.6958
COMMENTO: Il test in questo caso è significativo, in quanto non abbiamo elementi sufficienti per accettare l’ipotesi H0. La media campionaria è 494.7mm contro il riferimento di 500mm.
# --- 3.3 T-test a due campioni: differenze antropometriche per sesso ---
# H0: le medie sono uguali tra maschi e femmine
cat("\n--- T-test: Peso per Sesso ---\n")##
## --- T-test: Peso per Sesso ---
##
## Welch Two Sample t-test
##
## data: Peso by Sesso
## t = -12.115, df = 2488.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.4841 -207.3844
## sample estimates:
## mean in group F mean in group M
## 3161.061 3408.496
##
## --- T-test: Lunghezza per Sesso ---
##
## Welch Two Sample t-test
##
## data: Lunghezza by Sesso
## t = -9.5823, df = 2457.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.939001 -7.882672
## sample estimates:
## mean in group F mean in group M
## 489.7641 499.6750
##
## --- T-test: Cranio per Sesso ---
##
## Welch Two Sample t-test
##
## data: Cranio by Sesso
## t = -7.4366, df = 2489.4, p-value = 1.414e-13
## alternative hypothesis: true difference in means between group F and group M is not equal to 0
## 95 percent confidence interval:
## -6.110504 -3.560417
## sample estimates:
## mean in group F mean in group M
## 337.6231 342.4586
COMMENTO: I test sono tutti significativi. I maschi sono sistematicamente più grandi in tutte le misure. La differenza di peso di 247g è clinicamente rilevante oltre che statisticamente significativa.
# --- 4.1 Modello completo con interazione teoricamente motivata ---
# L'interazione Fumatrici*Gestazione testa se l'effetto delle settimane
# di gestazione sul peso cambia in funzione dello status di fumo
modello_full <- lm(Peso ~ Anni.madre + N.gravidanze + Fumatrici * Gestazione +
Lunghezza + Cranio + Tipo.parto + Ospedale + Sesso,
data = df)
summary(modello_full)##
## Call:
## lm(formula = Peso ~ Anni.madre + N.gravidanze + Fumatrici * Gestazione +
## Lunghezza + Cranio + Tipo.parto + Ospedale + Sesso, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1122.52 -181.96 -14.95 162.19 2611.87
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6753.4971 142.3219 -47.452 < 2e-16 ***
## Anni.madre 0.7781 1.1468 0.678 0.4976
## N.gravidanze 11.4666 4.6689 2.456 0.0141 *
## FumatriciFumatrice 830.2103 756.1258 1.098 0.2723
## Gestazione 33.1481 3.8533 8.603 < 2e-16 ***
## Lunghezza 10.2864 0.3009 34.185 < 2e-16 ***
## Cranio 10.4670 0.4263 24.555 < 2e-16 ***
## Tipo.partoNat 29.8321 12.0910 2.467 0.0137 *
## Ospedaleosp2 -10.7584 13.4494 -0.800 0.4238
## Ospedaleosp3 28.5942 13.5080 2.117 0.0344 *
## SessoM 78.1833 11.1987 6.981 3.73e-12 ***
## FumatriciFumatrice:Gestazione -21.9184 19.2473 -1.139 0.2549
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274 on 2486 degrees of freedom
## Multiple R-squared: 0.729, Adjusted R-squared: 0.7278
## F-statistic: 608.1 on 11 and 2486 DF, p-value: < 2.2e-16
# --- 4.2 Test esplicito dell'interazione con F-test ---
# Confronto diretto tra modello con e senza interazione
modello_no_interazione <- lm(Peso ~ Anni.madre + N.gravidanze + Fumatrici + Gestazione +
Lunghezza + Cranio + Tipo.parto + Ospedale + Sesso,
data = df)
cat("\n--- F-test: il termine di interazione Fumatrici:Gestazione è significativo? ---\n")##
## --- F-test: il termine di interazione Fumatrici:Gestazione è significativo? ---
## Analysis of Variance Table
##
## Model 1: Peso ~ Anni.madre + N.gravidanze + Fumatrici + Gestazione + Lunghezza +
## Cranio + Tipo.parto + Ospedale + Sesso
## Model 2: Peso ~ Anni.madre + N.gravidanze + Fumatrici * Gestazione + Lunghezza +
## Cranio + Tipo.parto + Ospedale + Sesso
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 2487 186743194
## 2 2486 186645831 1 97363 1.2968 0.2549
COMMENTO: Effettuando un F-test considerando il modello senza interazione e quello con l’interazione, possiamo giungere alla conclusione che l’interazione risulta essere non significativa.
# --- 4.3 Controllo multicollinearità (VIF) ---
# VIF > 5: potenziale problema; VIF > 10: problema serio
cat("\n--- Variance Inflation Factor (VIF) ---\n")##
## --- Variance Inflation Factor (VIF) ---
## GVIF Df GVIF^(1/(2*Df))
## Anni.madre 1.190241 1 1.090982
## N.gravidanze 1.189278 1 1.090540
## Fumatrici 1.007426 1 1.003706
## Gestazione 1.695675 1 1.302181
## Lunghezza 2.086879 1 1.444604
## Cranio 1.631049 1 1.277125
## Tipo.parto 1.004227 1 1.002111
## Ospedale 1.004267 2 1.001065
## Sesso 1.040743 1 1.020168
COMMENTO: Non è presente un problema di multicollinearità.
# --- 4.4 Selezione del modello ottimale con stepAIC ---
modello_aic <- stepAIC(modello_full, direction = "both", trace = FALSE)
cat("\n--- Modello ottimale (AIC) ---\n")##
## --- Modello ottimale (AIC) ---
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio +
## Tipo.parto + Ospedale + Sesso, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1113.07 -181.71 -16.66 161.08 2619.57
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6707.9252 136.0257 -49.314 < 2e-16 ***
## N.gravidanze 12.3360 4.3344 2.846 0.00446 **
## Gestazione 32.0386 3.7925 8.448 < 2e-16 ***
## Lunghezza 10.3059 0.3006 34.286 < 2e-16 ***
## Cranio 10.4920 0.4257 24.648 < 2e-16 ***
## Tipo.partoNat 29.4080 12.0875 2.433 0.01505 *
## Ospedaleosp2 -10.8939 13.4447 -0.810 0.41786
## Ospedaleosp3 28.7917 13.4969 2.133 0.03301 *
## SessoM 77.4657 11.1842 6.926 5.48e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 274 on 2489 degrees of freedom
## Multiple R-squared: 0.7287, Adjusted R-squared: 0.7278
## F-statistic: 835.7 on 8 and 2489 DF, p-value: < 2.2e-16
# Confronto anche con selezione BIC
modello_bic <- stepAIC(modello_full, direction = "both",
k = log(nrow(df)), trace = FALSE)
cat("\n--- Modello ottimale (BIC) ---\n")##
## --- Modello ottimale (BIC) ---
##
## Call:
## lm(formula = Peso ~ N.gravidanze + Gestazione + Lunghezza + Cranio +
## Sesso, data = df)
##
## 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
# Confronto AIC e BIC tra i due modelli selezionati
cat("\n--- Confronto AIC/BIC tra modello AIC-ottimale e BIC-ottimale ---\n")##
## --- Confronto AIC/BIC tra modello AIC-ottimale e BIC-ottimale ---
## df AIC
## modello_aic 10 35143.29
## modello_bic 7 35152.89
## df BIC
## modello_aic 10 35201.52
## modello_bic 7 35193.65
COMMENTO: I due modelli concordano sul nucleo (le 4 variabili biologiche + sesso) ma divergono sui predittori ospedalieri. La differenza di R^2 risulta essere minima. Inoltre risulta essere interessante il fatto che la variabile Fumatrici è stata eliminata, considerata l’assunzione teorica.
# --- 5.1 RMSE in-sample
rmse_insample <- sqrt(mean(modello_ottimale$residuals^2))
cat("\nRMSE in-sample:", round(rmse_insample, 2), "grammi\n")##
## RMSE in-sample: 273.51 grammi
# --- 5.2 RMSE con 10-fold cross-validation
set.seed(123)
ctrl <- trainControl(method = "cv", number = 10)
modello_cv <- train(formula(modello_ottimale),
data = df,
method = "lm",
trControl = ctrl)
cat("\nRMSE 10-fold CV:", round(modello_cv$results$RMSE, 2), "grammi\n")##
## RMSE 10-fold CV: 274.24 grammi
## R² 10-fold CV: 0.7266
COMMENTO: Considerati i risultati ottenuti dal modello in-sample e utilizzando la cross validation per testarlo out-of-sample, possiamo dire che la differenza minima è trascurabile e il modello non fa overfitting.
# --- 5.3 Test formali sulle assunzioni della regressione ---
# Normalità dei residui (Shapiro-Wilk; affidabile per n <= 5000)
cat("\n--- Test di Shapiro-Wilk (normalità dei residui) ---\n")##
## --- Test di Shapiro-Wilk (normalità dei residui) ---
##
## Shapiro-Wilk normality test
##
## data: residuals(modello_ottimale)
## W = 0.97414, p-value < 2.2e-16
# Omoschedasticità (Breusch-Pagan)
# H0: varianza dei residui costante (omoschedasticità)
cat("\n--- Test di Breusch-Pagan (omoschedasticità) ---\n")##
## --- Test di Breusch-Pagan (omoschedasticità) ---
##
## studentized Breusch-Pagan test
##
## data: modello_ottimale
## BP = 91.829, df = 8, p-value < 2.2e-16
COMMENTO NORMALITÀ: Considerando il p-value del test BP andiamo a “rifiutare” l’ipotesi H0 di normalità. Dal grafico Q-Q dei residui possiamo notare una leggera coda destra. Possiamo concludere che si tratta di una deviazione tollerabile.
COMMENTO OMOSCHEDASTICITÀ: Considerando il p-value del test BP andiamo a “rifiutare” l’ipotesi H0 di omoschedasticità. Dal grafico della Scale-Location emerge un aumento della variabilità all’aumentare del peso dei neonati.
# --- 5.5 Valori influenti: distanza di Cook ---
soglia_cook <- 4 / nrow(df)
plot(cooks.distance(modello_ottimale), type = "h",
main = "Distanza di Cook (Valori Influenti)",
ylab = "Distanza di Cook", xlab = "Osservazione")
abline(h = soglia_cook, col = "red", lty = 2)
legend("topright", legend = paste("Soglia =", round(soglia_cook, 5)),
col = "red", lty = 2)# Ispezione dei valori influenti
influenti_idx <- which(cooks.distance(modello_ottimale) > soglia_cook)
cat("\nNumero di osservazioni influenti (Cook > soglia):", length(influenti_idx), "\n")##
## Numero di osservazioni influenti (Cook > soglia): 116
## --- Dettaglio osservazioni influenti ---
## Anni.madre N.gravidanze Fumatrici Gestazione Peso Lunghezza Cranio
## 13 36 5 Non Fumatrice 38 3060 455 325
## 34 27 0 Non Fumatrice 39 3150 480 382
## 119 31 0 Non Fumatrice 40 3410 550 372
## 130 30 2 Non Fumatrice 39 4240 485 352
## 134 38 6 Non Fumatrice 37 3950 500 350
## 140 29 1 Fumatrice 41 4420 530 362
## 146 24 1 Non Fumatrice 40 3820 500 320
## 155 30 0 Non Fumatrice 36 3610 410 330
## 161 35 9 Non Fumatrice 42 3760 540 348
## 220 23 1 Non Fumatrice 40 3520 445 363
## 274 28 0 Non Fumatrice 39 4320 530 350
## 295 18 0 Non Fumatrice 40 1850 460 305
## 310 40 3 Non Fumatrice 28 1560 420 379
## 329 25 1 Non Fumatrice 40 4560 540 340
## 375 33 4 Non Fumatrice 38 4270 510 353
## 377 31 1 Non Fumatrice 39 3400 460 325
## 390 38 0 Non Fumatrice 40 3700 470 320
## 455 24 0 Non Fumatrice 38 2920 505 351
## 472 25 2 Non Fumatrice 41 3990 495 335
## 478 32 3 Non Fumatrice 42 4310 505 365
## 516 40 8 Non Fumatrice 38 3520 470 341
## 582 30 7 Non Fumatrice 35 2220 470 316
## 615 27 1 Non Fumatrice 35 2100 440 345
## 616 29 0 Non Fumatrice 42 3540 540 368
## 632 21 0 Non Fumatrice 37 2750 510 333
## 633 30 0 Non Fumatrice 40 3860 480 338
## 648 35 1 Non Fumatrice 38 2850 500 360
## 656 38 3 Non Fumatrice 41 2320 450 330
## 657 32 3 Non Fumatrice 41 2970 520 318
## 684 30 1 Non Fumatrice 39 3000 475 390
## 791 30 1 Non Fumatrice 41 4440 510 335
## 805 30 2 Non Fumatrice 29 1190 360 272
## 828 32 6 Non Fumatrice 40 4200 510 350
## 890 33 2 Non Fumatrice 39 3660 490 318
## 908 35 2 Non Fumatrice 39 3300 460 320
## 928 25 0 Non Fumatrice 28 830 310 254
## 950 25 0 Non Fumatrice 41 2980 500 356
## 1008 32 5 Non Fumatrice 36 2400 470 333
## 1014 17 0 Non Fumatrice 37 2050 390 295
## 1036 26 0 Non Fumatrice 40 4330 500 355
## 1130 33 11 Non Fumatrice 43 3400 475 360
## 1132 21 0 Non Fumatrice 39 3950 500 330
## 1137 26 0 Non Fumatrice 40 3000 520 340
## 1181 30 1 Non Fumatrice 36 4070 500 373
## 1188 21 0 Non Fumatrice 40 4140 550 320
## 1192 24 2 Non Fumatrice 40 2560 490 323
## 1194 41 1 Non Fumatrice 39 2160 450 325
## 1212 26 0 Non Fumatrice 34 2050 460 325
## 1230 29 1 Non Fumatrice 41 4010 470 353
## 1253 35 4 Non Fumatrice 39 3450 523 365
## 1267 24 0 Non Fumatrice 39 4260 545 343
## 1268 28 1 Non Fumatrice 40 3790 460 332
## 1273 32 1 Non Fumatrice 33 2040 480 307
## 1287 35 3 Non Fumatrice 40 2660 485 335
## 1293 30 3 Non Fumatrice 38 4600 485 380
## 1306 23 0 Non Fumatrice 41 4900 510 352
## 1317 28 2 Non Fumatrice 39 3650 500 315
## 1341 32 3 Non Fumatrice 36 3780 480 349
## 1395 30 2 Non Fumatrice 39 3790 505 304
## 1399 42 2 Non Fumatrice 38 2560 525 349
## 1402 35 1 Non Fumatrice 39 2660 460 364
## 1426 32 4 Fumatrice 39 2250 460 319
## 1429 24 4 Non Fumatrice 29 1280 390 355
## 1433 31 4 Non Fumatrice 41 4810 530 364
## 1450 36 8 Non Fumatrice 41 3730 480 335
## 1472 37 4 Non Fumatrice 37 3720 480 355
## 1499 30 1 Non Fumatrice 39 2380 495 325
## 1505 30 8 Non Fumatrice 39 2860 490 337
## 1541 30 0 Non Fumatrice 38 4540 530 343
## 1551 35 1 Non Fumatrice 38 4370 315 374
## 1553 30 4 Non Fumatrice 35 4520 520 360
## 1556 37 0 Non Fumatrice 41 2420 490 300
## 1588 22 0 Non Fumatrice 42 3100 510 361
## 1593 41 3 Non Fumatrice 35 1500 420 304
## 1619 31 0 Non Fumatrice 31 990 340 278
## 1635 32 2 Non Fumatrice 39 3430 445 322
## 1639 39 3 Non Fumatrice 40 4760 550 365
## 1694 23 1 Non Fumatrice 36 3850 460 334
## 1702 21 0 Non Fumatrice 41 3000 495 357
## 1712 28 0 Non Fumatrice 39 3800 520 300
## 1718 34 4 Non Fumatrice 42 2660 500 320
## 1780 25 2 Non Fumatrice 25 900 325 253
## 1838 34 3 Non Fumatrice 40 4580 515 360
## 1856 33 0 Non Fumatrice 36 3110 465 300
## 1868 29 0 Non Fumatrice 40 3470 525 390
## 1893 22 1 Non Fumatrice 34 3030 470 312
## 1915 27 1 Non Fumatrice 38 2480 500 315
## 1920 26 0 Fumatrice 39 4930 550 350
## 1926 29 0 Non Fumatrice 41 3840 485 358
## 1937 36 1 Non Fumatrice 41 2940 500 358
## 1962 38 4 Non Fumatrice 38 4370 530 340
## 1963 27 0 Non Fumatrice 42 4700 540 362
## 2023 27 1 Non Fumatrice 39 4650 510 354
## 2040 27 1 Non Fumatrice 38 3240 410 359
## 2076 23 0 Non Fumatrice 40 4720 540 360
## 2086 26 8 Non Fumatrice 40 3250 500 355
## 2114 36 0 Non Fumatrice 31 1180 355 270
## 2115 35 1 Non Fumatrice 32 1890 500 309
## 2120 32 0 Non Fumatrice 27 1140 370 267
## 2123 32 0 Non Fumatrice 39 3950 520 330
## 2135 33 3 Non Fumatrice 39 3950 480 345
## 2175 37 8 Non Fumatrice 28 930 355 235
## 2195 38 1 Non Fumatrice 40 3980 480 335
## 2204 33 0 Non Fumatrice 39 4000 490 360
## 2219 37 1 Non Fumatrice 39 2500 490 352
## 2221 35 10 Non Fumatrice 39 2950 495 335
## 2225 27 0 Non Fumatrice 35 3140 465 290
## 2287 33 3 Non Fumatrice 40 3640 470 337
## 2315 24 0 Non Fumatrice 42 2800 520 340
## 2317 25 6 Non Fumatrice 38 2530 460 340
## 2392 28 1 Non Fumatrice 40 4720 540 355
## 2421 20 0 Non Fumatrice 41 4160 500 358
## 2422 33 10 Non Fumatrice 40 3090 485 353
## 2437 28 1 Non Fumatrice 27 980 320 265
## 2452 28 0 Non Fumatrice 26 930 345 245
## 2471 34 10 Non Fumatrice 38 2880 470 345
## Tipo.parto Ospedale Sesso
## 13 Ces osp1 F
## 34 Nat osp1 F
## 119 Nat osp2 M
## 130 Nat osp2 M
## 134 Nat osp3 M
## 140 Ces osp2 F
## 146 Ces osp2 F
## 155 Nat osp1 M
## 161 Nat osp2 F
## 220 Nat osp1 F
## 274 Ces osp3 M
## 295 Nat osp3 F
## 310 Nat osp3 F
## 329 Nat osp1 M
## 375 Nat osp2 M
## 377 Nat osp1 M
## 390 Nat osp1 M
## 455 Ces osp3 M
## 472 Nat osp3 M
## 478 Nat osp1 M
## 516 Nat osp3 M
## 582 Nat osp3 M
## 615 Ces osp2 F
## 616 Nat osp1 M
## 632 Nat osp1 M
## 633 Nat osp1 M
## 648 Nat osp1 F
## 656 Nat osp2 F
## 657 Ces osp2 M
## 684 Ces osp2 F
## 791 Nat osp3 M
## 805 Nat osp2 F
## 828 Nat osp1 M
## 890 Nat osp2 F
## 908 Ces osp1 M
## 928 Nat osp1 F
## 950 Nat osp3 M
## 1008 Nat osp1 F
## 1014 Nat osp2 F
## 1036 Nat osp3 F
## 1130 Nat osp1 M
## 1132 Nat osp1 F
## 1137 Ces osp3 M
## 1181 Nat osp2 M
## 1188 Ces osp1 M
## 1192 Nat osp1 M
## 1194 Ces osp3 M
## 1212 Nat osp1 F
## 1230 Nat osp2 M
## 1253 Nat osp1 M
## 1267 Ces osp2 F
## 1268 Nat osp2 F
## 1273 Ces osp1 F
## 1287 Ces osp3 M
## 1293 Nat osp1 M
## 1306 Nat osp2 F
## 1317 Nat osp1 M
## 1341 Nat osp3 M
## 1395 Ces osp3 M
## 1399 Ces osp2 M
## 1402 Ces osp1 F
## 1426 Nat osp2 F
## 1429 Nat osp1 F
## 1433 Ces osp3 M
## 1450 Nat osp3 M
## 1472 Ces osp1 F
## 1499 Nat osp1 F
## 1505 Ces osp2 F
## 1541 Ces osp3 M
## 1551 Nat osp3 F
## 1553 Nat osp2 F
## 1556 Ces osp1 M
## 1588 Nat osp2 F
## 1593 Nat osp1 M
## 1619 Ces osp2 F
## 1635 Ces osp1 F
## 1639 Nat osp2 F
## 1694 Ces osp3 F
## 1702 Ces osp1 M
## 1712 Nat osp3 F
## 1718 Nat osp2 F
## 1780 Nat osp3 F
## 1838 Nat osp2 M
## 1856 Ces osp3 M
## 1868 Nat osp1 M
## 1893 Nat osp2 F
## 1915 Ces osp2 M
## 1920 Ces osp2 F
## 1926 Ces osp2 F
## 1937 Nat osp3 M
## 1962 Nat osp3 F
## 1963 Nat osp3 F
## 2023 Nat osp2 M
## 2040 Ces osp1 F
## 2076 Nat osp3 F
## 2086 Nat osp2 M
## 2114 Nat osp3 F
## 2115 Nat osp2 F
## 2120 Nat osp3 F
## 2123 Ces osp1 F
## 2135 Nat osp3 M
## 2175 Nat osp1 F
## 2195 Nat osp1 F
## 2204 Ces osp3 F
## 2219 Ces osp1 M
## 2221 Nat osp1 F
## 2225 Nat osp2 F
## 2287 Nat osp3 F
## 2315 Ces osp2 M
## 2317 Nat osp1 F
## 2392 Nat osp3 M
## 2421 Ces osp1 M
## 2422 Nat osp3 M
## 2437 Nat osp1 M
## 2452 Ces osp3 F
## 2471 Ces osp2 M
COMMENTO: 116 osservazioni (circa il 4.6% del campione) superano la soglia di Cook. Tra le prime si trovano casi come l’osservazione 140 (madre fumatrice, neonato di 4420 g a 41 settimane) e la 130 (4240 g con lunghezza pari a 485 mm), che presentano combinazioni di valori piuttosto insolite. Si tratta di punti che il modello riesce a rappresentare con più difficoltà e che possono avere un certo peso nella stima complessiva.
# --- 6.1 Impatto gestazione e fumo sul peso ---
ggplot(df, aes(x = Gestazione, y = Peso, color = Fumatrici)) +
geom_point(alpha = 0.4) +
geom_smooth(method = "lm", se = TRUE) +
theme_minimal() +
scale_color_manual(values = c("Non Fumatrice" = "#2c3e50",
"Fumatrice" = "#e74c3c")) +
labs(
title = "Impatto della Gestazione e del Fumo sul Peso Neonatale",
subtitle = "Pendenze diverse indicano un effetto di interazione",
x = "Settimane di Gestazione",
y = "Peso del Neonato (g)",
color = "Status Madre"
)COMMENTO: A parità di settimane di gestazione, i neonati di madri fumatrici tendono a pesare meno. In entrambi i gruppi all’aumentare della gestazione -> aumenta il peso.
# --- 6.3 Coefficienti del modello con intervalli di confidenza ---
tidy(modello_ottimale, conf.int = TRUE) |>
filter(term != "(Intercept)") |>
ggplot(aes(x = reorder(term, estimate), y = estimate,
ymin = conf.low, ymax = conf.high)) +
geom_pointrange() +
geom_hline(yintercept = 0, linetype = "dashed", color = "red") +
coord_flip() +
theme_minimal() +
labs(title = "Coefficienti del Modello con IC al 95%",
x = NULL, y = "Stima (grammi)")# Neonata femmina, madre alla 3a gravidanza, 39 settimane di gestazione.
nuovo_caso <- data.frame(
Anni.madre = 30,
N.gravidanze = 3,
Fumatrici = factor("Non Fumatrice",
levels = levels(df$Fumatrici)),
Gestazione = 39,
Lunghezza = round(mean(df$Lunghezza)), # media campionaria
Cranio = round(mean(df$Cranio)), # media campionaria
Tipo.parto = factor("Nat",
levels = levels(df$Tipo.parto)),
Ospedale = factor("osp1",
levels = levels(df$Ospedale)),
Sesso = factor("F",
levels = levels(df$Sesso))
)
previsione <- predict(modello_ottimale,
newdata = nuovo_caso,
interval = "prediction",
level = 0.95)
cat("\n=== PREVISIONE PESO NEONATO ===\n")##
## === PREVISIONE PESO NEONATO ===
## Stima puntuale: 3277 grammi
cat("Intervallo di previsione al 95%: [",
round(previsione[2], 0), ",",
round(previsione[3], 0), "] grammi\n")## Intervallo di previsione al 95%: [ 2739 , 3815 ] grammi
cat("(L'intervallo di PREVISIONE, più ampio dell'IC sulla media,",
"riflette l'incertezza sulla singola nascita futura)\n")## (L'intervallo di PREVISIONE, più ampio dell'IC sulla media, riflette l'incertezza sulla singola nascita futura)