Εισαγωγή και Παρουσίαση του Dataset

Για την παρούσα εργασία χρησιμοποιείται το dataset Life Expectancy Data (του Παγκόσμιου Οργανισμού Υγείας). Στόχος μας είναι να δημιουργήσουμε ένα μοντέλο Λογιστικής Παλινδρόμησης για να προβλέψουμε αν μια χώρα υψηλό ή χαμηλό προσδόκιμο ζωής.

Επιλέξαμε τις εξής ανεξάρτητες μεταβλητές:

  • Schooling: Έτη εκπαίδευσης.

  • Adult.Mortality: Ποσοστό θνησιμότητας ενηλίκων.

  • Income.composition.of.resources: Δείκτης σύνθεσης εισοδήματος.

  • HIV.AIDS: Θάνατοι από HIV/AIDS.

Επειδή η Λογιστική Παλινδρόμηση απαιτεί δίτιμη εξαρτημένη μεταβλητή, μετατρέπουμε το Life.expectancy σε μια νέα μεταβλητή High_Life, η οποία παίρνει την τιμή 1 αν το προσδόκιμο είναι \(\ge 70\) έτη, και 0 διαφορετικά.

Προετοιμασία Δεδομένων & Βιβλιοθήκες

Αρχικά, φορτώνουμε τις απαραίτητες βιβλιοθήκες, εισάγουμε τα δεδομένα και κρατάμε μόνο τις στήλες που μας ενδιαφέρουν.

# Φόρτωση βιβλιοθηκών
library(readr)
library(ggplot2)
library(dplyr)
## 
## 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
library(caTools)
## Warning: package 'caTools' was built under R version 4.5.3
library(ROCR)
## Warning: package 'ROCR' was built under R version 4.5.3
# Εισαγωγή δεδομένων
Life_Data <- read.csv("Life_Expectancy_Data.csv")

str(Life_Data)
## '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 ...
summary(Life_Data)
##    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

Διαχωρισμός σε Train και Test Sets

Διαχωρίζουμε τα καθαρά δεδομένα τυχαία, δίνοντας το 65% στο εκπαιδευτικό σύνολο. Ορίζουμε το seed στο 940 . Ονομάζουμε τα sets train και test.

set.seed(940)

split <- sample.split(df_subset$ High_Life,SplitRatio=0.65)
#Δημιουρφία των 2 νέων sets(train, test)
train = subset(df_subset,split==TRUE)
test = subset(df_subset,split==FALSE)

nrow(train) 
## [1] 1903
nrow(test) 
## [1] 1035
View(train)
View(test)  

Σχολιασμός Set

Για το 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
summary(LifeExpectancyLog)
## 
## 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, μειώνεται η πιθανότητα για υψηλό προσδόκιμο ζωής.

Προβλέψεις (Predictions) στο Test Set

Χρησιμοποιούμε το μοντέλο για να κάνουμε προβλέψεις πάνω σε νέα δεδομένα (στο test).

predictTest <- predict(LifeExpectancyLog, type = "response", newdata = test)

summary(predictTest)
##    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.

Confusion Matrix και Αξιολόγηση Μοντέλου

Ορίζουμε ως κατώφλι (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
cat("Ευαισθησία (Sensitivity):", round(sensitivity, 4), "\n")
## Ευαισθησία (Sensitivity): 0.9413
cat("Ειδικότητα (Specificity):", round(specificity, 4), "\n")
## Ειδικότητα (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ότε το μοντέλο μας όντως μαθαίνει χρήσιμες πληροφορίες από τις ανεξάρτητες μεταβλητές.

Δημιουργία train2 και test2 με na.omit

Σύμφωνα με την οδηγία, δημιουργούμε νέα σύνολα χωρίς κενές τιμές για την ανάλυση ROCR.

train2 <- na.omit(train)
test2 <- na.omit(test)

cat("Καταχωρήσεις στο train2:", nrow(train2), "\n")
## Καταχωρήσεις στο train2: 1793
cat("Καταχωρήσεις στο test2:", nrow(test2), "\n")
## Καταχωρήσεις στο test2: 975

Σχολιασμός Set

Για το set train2 έχουμε 1793 καταχωρήσεις και για το set test2 έχουμε 975 καταχωρήσεις.

Καμπύλη ROC και AUC

Τέλος, χρησιμοποιούμε τη βιβλιοθήκη 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 , επιβεβαιώνει την εξαιρετική ικανότητα διάκρισης του μοντέλου μεταξύ χωρών με χαμηλό και υψηλό προσδόκιμο.