Για την παρούσα εργασία χρησιμοποιείται το dataset Life Expectancy Data (του Παγκόσμιου Οργανισμού Υγείας). Στόχος μας είναι να δημιουργήσουμε ένα μοντέλο Λογιστικής Παλινδρόμησης για να προβλέψουμε αν μια χώρα υψηλό ή χαμηλό προσδόκιμο ζωής.
Επιλέξαμε τις εξής ανεξάρτητες μεταβλητές:
Schooling: Έτη εκπαίδευσης.
Adult.Mortality: Ποσοστό θνησιμότητας ενηλίκων.
Income.composition.of.resources: Δείκτης σύνθεσης εισοδήματος.
HIV.AIDS: Θάνατοι από HIV/AIDS.
Επειδή η Λογιστική Παλινδρόμηση απαιτεί δίτιμη εξαρτημένη μεταβλητή, μετατρέπουμε το Life.expectancy σε μια νέα μεταβλητή High_Life, η οποία παίρνει την τιμή 1 αν το προσδόκιμο είναι \(\ge 70\) έτη, και 0 διαφορετικά.
Προετοιμασία Δεδομένων & Βιβλιοθήκες
Αρχικά, φορτώνουμε τις απαραίτητες βιβλιοθήκες, εισάγουμε τα δεδομένα και κρατάμε μόνο τις στήλες που μας ενδιαφέρουν.
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
## Warning: package 'caTools' was built under R version 4.5.3
## Warning: package 'ROCR' was built under R version 4.5.3
## 'data.frame': 2938 obs. of 22 variables:
## $ Country : chr "Afghanistan" "Afghanistan" "Afghanistan" "Afghanistan" ...
## $ Year : int 2015 2014 2013 2012 2011 2010 2009 2008 2007 2006 ...
## $ Status : chr "Developing" "Developing" "Developing" "Developing" ...
## $ Life.expectancy : num 65 59.9 59.9 59.5 59.2 58.8 58.6 58.1 57.5 57.3 ...
## $ Adult.Mortality : int 263 271 268 272 275 279 281 287 295 295 ...
## $ infant.deaths : int 62 64 66 69 71 74 77 80 82 84 ...
## $ Alcohol : num 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.03 0.02 0.03 ...
## $ percentage.expenditure : num 71.3 73.5 73.2 78.2 7.1 ...
## $ Hepatitis.B : int 65 62 64 67 68 66 63 64 63 64 ...
## $ Measles : int 1154 492 430 2787 3013 1989 2861 1599 1141 1990 ...
## $ BMI : num 19.1 18.6 18.1 17.6 17.2 16.7 16.2 15.7 15.2 14.7 ...
## $ under.five.deaths : int 83 86 89 93 97 102 106 110 113 116 ...
## $ Polio : int 6 58 62 67 68 66 63 64 63 58 ...
## $ Total.expenditure : num 8.16 8.18 8.13 8.52 7.87 9.2 9.42 8.33 6.73 7.43 ...
## $ Diphtheria : int 65 62 64 67 68 66 63 64 63 58 ...
## $ HIV.AIDS : num 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 ...
## $ GDP : num 584.3 612.7 631.7 670 63.5 ...
## $ Population : num 33736494 327582 31731688 3696958 2978599 ...
## $ thinness..1.19.years : num 17.2 17.5 17.7 17.9 18.2 18.4 18.6 18.8 19 19.2 ...
## $ thinness.5.9.years : num 17.3 17.5 17.7 18 18.2 18.4 18.7 18.9 19.1 19.3 ...
## $ Income.composition.of.resources: num 0.479 0.476 0.47 0.463 0.454 0.448 0.434 0.433 0.415 0.405 ...
## $ Schooling : num 10.1 10 9.9 9.8 9.5 9.2 8.9 8.7 8.4 8.1 ...
## Country Year Status Life.expectancy
## Length:2938 Min. :2000 Length:2938 Min. :36.30
## Class :character 1st Qu.:2004 Class :character 1st Qu.:63.10
## Mode :character Median :2008 Mode :character Median :72.10
## Mean :2008 Mean :69.22
## 3rd Qu.:2012 3rd Qu.:75.70
## Max. :2015 Max. :89.00
## NA's :10
## Adult.Mortality infant.deaths Alcohol percentage.expenditure
## Min. : 1.0 Min. : 0.0 Min. : 0.0100 Min. : 0.000
## 1st Qu.: 74.0 1st Qu.: 0.0 1st Qu.: 0.8775 1st Qu.: 4.685
## Median :144.0 Median : 3.0 Median : 3.7550 Median : 64.913
## Mean :164.8 Mean : 30.3 Mean : 4.6029 Mean : 738.251
## 3rd Qu.:228.0 3rd Qu.: 22.0 3rd Qu.: 7.7025 3rd Qu.: 441.534
## Max. :723.0 Max. :1800.0 Max. :17.8700 Max. :19479.912
## NA's :10 NA's :194
## Hepatitis.B Measles BMI under.five.deaths
## Min. : 1.00 Min. : 0.0 Min. : 1.00 Min. : 0.00
## 1st Qu.:77.00 1st Qu.: 0.0 1st Qu.:19.30 1st Qu.: 0.00
## Median :92.00 Median : 17.0 Median :43.50 Median : 4.00
## Mean :80.94 Mean : 2419.6 Mean :38.32 Mean : 42.04
## 3rd Qu.:97.00 3rd Qu.: 360.2 3rd Qu.:56.20 3rd Qu.: 28.00
## Max. :99.00 Max. :212183.0 Max. :87.30 Max. :2500.00
## NA's :553 NA's :34
## Polio Total.expenditure Diphtheria HIV.AIDS
## Min. : 3.00 Min. : 0.370 Min. : 2.00 Min. : 0.100
## 1st Qu.:78.00 1st Qu.: 4.260 1st Qu.:78.00 1st Qu.: 0.100
## Median :93.00 Median : 5.755 Median :93.00 Median : 0.100
## Mean :82.55 Mean : 5.938 Mean :82.32 Mean : 1.742
## 3rd Qu.:97.00 3rd Qu.: 7.492 3rd Qu.:97.00 3rd Qu.: 0.800
## Max. :99.00 Max. :17.600 Max. :99.00 Max. :50.600
## NA's :19 NA's :226 NA's :19
## GDP Population thinness..1.19.years
## Min. :1.681e+00 Min. :3.400e+01 Min. : 0.10
## 1st Qu.:4.639e+02 1st Qu.:1.958e+05 1st Qu.: 1.60
## Median :1.767e+03 Median :1.387e+06 Median : 3.30
## Mean :7.483e+03 Mean :1.275e+07 Mean : 4.84
## 3rd Qu.:5.911e+03 3rd Qu.:7.420e+06 3rd Qu.: 7.20
## Max. :1.192e+05 Max. :1.294e+09 Max. :27.70
## NA's :448 NA's :652 NA's :34
## thinness.5.9.years Income.composition.of.resources Schooling
## Min. : 0.10 Min. :0.0000 Min. : 0.00
## 1st Qu.: 1.50 1st Qu.:0.4930 1st Qu.:10.10
## Median : 3.30 Median :0.6770 Median :12.30
## Mean : 4.87 Mean :0.6276 Mean :11.99
## 3rd Qu.: 7.20 3rd Qu.:0.7790 3rd Qu.:14.30
## Max. :28.60 Max. :0.9480 Max. :20.70
## NA's :34 NA's :167 NA's :163
# Επιλογή μεταβλητών που μας ζητήθηκαν
df_subset <- Life_Data[, c("Life.expectancy", "Schooling", "Adult.Mortality",
"Income.composition.of.resources", "HIV.AIDS")]
# Δημιουργία δίτιμης μεταβλητής
df_subset$High_Life <- ifelse(df_subset$Life.expectancy >= 70, 1, 0)
# Αφαιρούμε την αρχική συνεχή μεταβλητή
df_subset$Life.expectancy <- NULLΔιαχωρίζουμε τα καθαρά δεδομένα τυχαία, δίνοντας το 65% στο εκπαιδευτικό σύνολο. Ορίζουμε το seed στο 940 . Ονομάζουμε τα sets train και test.
#Δημιουρφία των 2 νέων sets(train, test)
train = subset(df_subset,split==TRUE)
test = subset(df_subset,split==FALSE)
nrow(train) ## [1] 1903
## [1] 1035
Για το set train έχουμε 1909 καταχωρήσεις και για το set test έχουμε 1029 καταχωρήσεις.
Δημιουργούμε το μοντέλο μας χρησιμοποιώντας τα δεδομένα του train.
LifeExpectancyLog <- glm(High_Life ~ Schooling + Adult.Mortality + Income.composition.of.resources + HIV.AIDS,
data = train, family = binomial)## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
##
## Call:
## glm(formula = High_Life ~ Schooling + Adult.Mortality + Income.composition.of.resources +
## HIV.AIDS, family = binomial, data = train)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.737977 0.554829 -8.540 < 2e-16 ***
## Schooling 0.521814 0.049185 10.609 < 2e-16 ***
## Adult.Mortality -0.013638 0.001224 -11.141 < 2e-16 ***
## Income.composition.of.resources 2.101216 0.551047 3.813 0.000137 ***
## HIV.AIDS -0.935123 0.179638 -5.206 1.93e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 2456.30 on 1792 degrees of freedom
## Residual deviance: 995.53 on 1788 degrees of freedom
## (110 observations deleted due to missingness)
## AIC: 1005.5
##
## Number of Fisher Scoring iterations: 9
Από την περίληψη (summary) του μοντέλου, κοιτάμε τη στήλη Pr(>|z|) (p-value).
Οι μεταβλητές με p-value < 0.05 θεωρούνται στατιστικά σημαντικές (δηλώνονται με αστεράκια ***). Όπου σε αυτή την περίπτωση όλες οι ανεξάρτητες μεταβλητές είναι στατιστικά σημαντικές
Το Schooling (0.521814) και το Income.composition.of.resources (2.101216) έχουν θετικό συντελεστή. Αυτό σημαίνει ότι όσο αυξάνονται τα έτη εκπαίδευσης και το εισόδημα, αυξάνεται σημαντικά και η πιθανότητα το προσδόκιμο ζωής να είναι υψηλό (δηλαδή \(\ge\) 70).
Το Adult.Mortality (-0.013638) και το HIV.AIDS (-0.935123) έχουν αρνητικό συντελεστή. Αυτό σημαίνει ότι όσο αυξάνεται η θνησιμότητα στους ενήλικες και τα κρούσματα HIV, μειώνεται η πιθανότητα για υψηλό προσδόκιμο ζωής.
Χρησιμοποιούμε το μοντέλο για να κάνουμε προβλέψεις πάνω σε νέα δεδομένα (στο test).
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.00000 0.06357 0.70291 0.55687 0.94502 0.99894 60
Τι ακριβώς μας δείχνει η predict; Η
predict με την παράμετρο type="response" δεν
μας δίνει 0 ή 1, αλλά πιθανότητες (από 0 έως 1).
Δηλαδή, για κάθε γραμμή (χώρα) του test set, μας δείχνει
ποια είναι η στατιστική πιθανότητα η συγκεκριμένη χώρα να έχει
Προσδόκιμο Ζωής \(\ge 70\) χρόνια,
σύμφωνα με αυτά που έμαθε το μοντέλο μας.
Δηλαδή στην περίπτωση αυτη η μέση στατιστική πίθανότητα η συγκεκριμένη χώρα να έχει Προσδόκιμο Ζωής \(\ge 70\) χρόνια είναι 0.55687.
Ορίζουμε ως κατώφλι (threshold) το 0.5. Αν η πιθανότητα > 0.5 προβλέπουμε 1, αλλιώς 0.
#table(test$High_Life, predictTest>0.5)
# Confusion Matrix
conf_matrix <- table(test$High_Life, predictTest > 0.5)
print(conf_matrix)##
## FALSE TRUE
## 0 363 67
## 1 32 513
# Υπολογισμός Μετρικών (προσαρμογή στα ονόματα του πίνακα)
# Τα κελιά είναι:
# conf_matrix[1,1] -> True Negatives (TN)
# conf_matrix[1,2] -> False Positives (FP)
# conf_matrix[2,1] -> False Negatives (FN)
# conf_matrix[2,2] -> True Positives (TP)
TN <- conf_matrix[1,1]
FP <- conf_matrix[1,2]
FN <- conf_matrix[2,1]
TP <- conf_matrix[2,2]
accuracy <- (TP + TN) / (TP + TN + FP + FN)
sensitivity <- TP / (TP + FN) # True Positive Rate
specificity <- TN / (TN + FP) # True Negative Rate
cat("\nΑκρίβεια (Accuracy):", round(accuracy, 4), "\n")##
## Ακρίβεια (Accuracy): 0.8985
## Ευαισθησία (Sensitivity): 0.9413
## Ειδικότητα (Specificity): 0.8442
# Baseline Model: Πρόβλεψη πάντα της πιο συχνής κατηγορίας στο test set
baseline_accuracy <- max(table(test$High_Life)) / nrow(test)
cat("\nΑκρίβεια Baseline Μοντέλου:", round(baseline_accuracy, 4), "\n")##
## Ακρίβεια Baseline Μοντέλου: 0.5507
Παρατηρούμε το Accuracy του μοντέλου μας (0.8985) είναι σημαντικά μεγαλύτερο από το Baseline Accuracy(0.5507) (δηλαδή την ακρίβεια που θα είχαμε αν μαντεύαμε τυφλά την πιο συχνή κατηγορία). Tότε το μοντέλο μας όντως μαθαίνει χρήσιμες πληροφορίες από τις ανεξάρτητες μεταβλητές.
Σύμφωνα με την οδηγία, δημιουργούμε νέα σύνολα χωρίς κενές τιμές για την ανάλυση ROCR.
## Καταχωρήσεις στο train2: 1793
## Καταχωρήσεις στο test2: 975
Για το set train2 έχουμε 1793 καταχωρήσεις και για το set test2 έχουμε 975 καταχωρήσεις.
Τέλος, χρησιμοποιούμε τη βιβλιοθήκη ROCR για να οπτικοποιήσουμε την απόδοση του μοντέλου. Χρησιμοποιούμε το test2 για να παράξουμε την καμπύλη ROC, καθώς η βιβλιοθήκη ROCR απαιτεί δεδομένα χωρίς NAs για να λειτουργήσει σωστά.
# Προβλέψεις στο καθαρό test2
predictTest2 <- predict(LifeExpectancyLog, type = "response", newdata = test2)
# Δημιουργία αντικειμένου ROCRpred
ROCRpred <- prediction(predictTest2, test2$High_Life)
# Καμπύλη ROC
ROCRperf <- performance(ROCRpred, "tpr", "fpr")
plot(ROCRperf, colorize=TRUE, print.cutoffs.at=seq(0,1,by=0.1))# Υπολογισμός AUC
auc_val <- performance(ROCRpred, measure = "auc")@y.values[[1]]
cat("Τιμή AUC:", auc_val, "\n")## Τιμή AUC: 0.9535481
Το μοντέλο Λογιστικής Παλινδρόμησης απέδειξε ότι δείκτες όπως τα έτη εκπαίδευσης, η θνησιμότητα ενηλίκων, το εισόδημα και η εξάπλωση του HIV παίζουν καθοριστικό ρόλο στην πρόβλεψη ενός υψηλού προσδόκιμου ζωής. Το μοντέλο μας πέτυχε ακρίβεια που ξεπερνά αισθητά το baseline model, ενώ η τιμή της περιοχής κάτω από την καμπύλη (AUC), που βρίσκεται κοντά στο 1.0 , επιβεβαιώνει την εξαιρετική ικανότητα διάκρισης του μοντέλου μεταξύ χωρών με χαμηλό και υψηλό προσδόκιμο.