Eșantionarea în jos pentru modelarea predictivă | R-BLOGGERS

Exemplu de lucru

library(medley)
data('pisa', package="medley")
data('pisa_variables', package="medley")

Programul de evaluare internațională a studenților (PISA) este un studiu internațional realizat de Organizația pentru Cooperare și Dezvoltare Economică (OCDE) la fiecare trei ani. Evaluează elevii în vârstă de 15 ani în matematică, știință și lectură în timp ce colectează informații despre elevi și școlile lor. pisa set de date incluse în medley Pachetul provine de la administrația din 2009 și este utilizat pentru a demonstra prezicerea prezenței în școală privată sau publică. Există 5.233 de observații pe 44 de variabile cu 93,4% studenți ai școlii publice și 6,6% studenți ai școlii private.

Pentru început, vom împărți datele într -un set de instruire și validare folosind splitstackshape::stratified() Funcție pentru a se asigura că raportul dintre elevii de școală public-privat este același în ambele seturi de date.

pisa_formu <- Public ~ .
names(pisa) <- pisa_variables(names(pisa))
pisa_splits <- splitstackshape::stratified(
    pisa, group = "Public", size = 0.75, bothSets = TRUE)
pisa_train <- pisa_splits((1)) |> as.data.frame()
pisa_valid <- pisa_splits((2)) |> as.data.frame()

table(pisa$Public, useNA = 'ifany') |> print() |> prop.table()

     FALSE       TRUE 
0.06592777 0.93407223 

table(pisa_train$Public, useNA = 'ifany') |> print() |> prop.table()

     FALSE       TRUE 
0.06598726 0.93401274 

table(pisa_valid$Public, useNA = 'ifany') |> print() |> prop.table()

     FALSE       TRUE 
0.06574924 0.93425076 

Putem estima un model de regresie logistică și obținem probabilitățile prezise pentru setul de date de validare.

pisa_lr_out <- glm(pisa_formu, data = pisa_train, family = binomial(link = 'logit'))
pisa_predictions <- predict(pisa_lr_out, newdata = pisa_valid, type="response")

Figura de mai jos arată distribuția probabilităților prevăzute pentru setul de date de validare. Există o oarecare separare între elevii de școală publică și privați, dar densitățile sunt clar centrate în partea dreaptă a gamei.

ggplot(data.frame(Public = pisa_valid$Public, 
                  Prediction = pisa_predictions), 
       aes(x = Prediction, color = Public)) +
  geom_density()

Figura de mai jos oferă o curbă caracteristică a operatorului de receptor (ROC) împreună cu o diagramă a preciziei, sensibilității și specificității.

calculate_roc(predictions = pisa_predictions, 
              observed = pisa_valid$Public) |> plot()

Matricea de confuzie de mai jos, împărțind la 0,5, indică faptul că acest model nu este mai bun decât modelul nul (adică procentul de studenți ai școlii publice este de 93,4%). Desigur, am putea ajusta această valoare redusă pentru a optimiza fie specificitate sau sensibilitate.

confusion_matrix(observed = pisa_valid$Public, 
                 predicted = pisa_predictions > 0.5)

           predicted              
  observed     FALSE          TRUE
     FALSE 1 (0.08%)    85 (6.50%)
      TRUE 4 (0.31%) 1218 (93.12%)
Accuracy: 93.2%
Sensitivity: 1.16%
Specificity: 99.67%

Micșorarea valorilor montate

Se dovedește că gama de valori montate din regresia logistică se va micșora pe măsură ce cantitatea de dezechilibru în variabila dependentă crește. Am întâlnit pentru prima dată această problemă când am estimat scoruri de propensiune pentru disertația mea într -un studiu al Charter versus elevii tradiționali ai școlii publice. În studiul respectiv folosind evaluarea națională a progresului educațional (NAEP), aproximativ 3% dintre studenți au participat la o școală charter. În acel studiu, gama de scoruri de propensiune a fost grav restricționată. Pentru a explora de ce, funcția multilevel :: PSRANGE () a fost dezvoltată rezultatul acestei funcții este figura de mai jos. Începând cu partea de jos, 345 de studenți ai școlii publice au fost selectați la întâmplare, astfel încât regresia logistică ar putea fi estimată acolo unde există un echilibru perfect în variabila dependentă. Pe măsură ce urcăm în sus, creștem raportul de la 1: 1 la 1:13. Pentru fiecare raport, sunt desenate 20 de probe aleatorii, modelul de regresie logistică estimat și se înregistrează valorile minime și maxime montate (adică probabilități prezise) (sunt reprezentate de punctele negre și barele verzi). Distribuțiile pe toate modelele sunt de asemenea incluse.

Transmiterea doar a intervalelor împreună cu media valorilor montate pentru elevii de școală publici (albaștri) și privați (verzi) arată că, odată ce raportul este mai mare de 3-la-1, media valorilor montate pentru clasa zero (școlile private din acest exemplu) este mai mare de 0,5.

Eșantionare în jos

După cum s -a discutat mai sus, unul dintre dezavantajele cheie ale eșantionării în jos este că, în situațiile în care există un dezechilibru semnificativ, excludem o mulțime de date din analiză. downsample() Funcția va determina mai întâi câte modele trebuie estimate astfel încât fiecare observație din clasa mai mare este utilizată exact o dată. Pentru acest exemplu, utilizăm un raport al studenților public-privat de 2-la-1, astfel încât pentru fiecare model estimat că există 259 de observații private și 518 de studenți publici. Având în vedere că există 3925 de observații în setul nostru de antrenament, dowmsample() Funcția va estima 7 modele.

pisa_ds_out <- downsample(
  formu = pisa_formu,
  data = pisa_train,
  model_fun = glm,
  ratio = 2,
  family = binomial(link = 'logit'))
length(pisa_ds_out)

Putem folosi predict() Funcție pentru a obține un cadru de date de predicții. Fiecare coloană corespunde valorii prevăzute pentru fiecare din cele 7 modele.

pisa_predictions_ds <- predict(pisa_ds_out,
                               newdata = pisa_valid, 
                               type="response")
head(pisa_predictions_ds)

     model1    model2    model3    model4    model5    model6    model7
1 0.8511828 0.7437341 0.8921605 0.8424369 0.7347052 0.8697531 0.6928875
2 0.7393116 0.6822714 0.9466815 0.7959953 0.8118642 0.9441840 0.9580830
3 0.4944206 0.3813138 0.5575741 0.3586561 0.5023435 0.5805062 0.6281852
4 0.8525691 0.8268514 0.8293386 0.8372777 0.9464037 0.8843848 0.9016058
5 0.1823382 0.3670335 0.4556063 0.1408078 0.1899378 0.3578418 0.2657968
6 0.9216160 0.8192096 0.9040353 0.9213184 0.8080822 0.9076342 0.8768295

Putem media predicțiile pentru a obține un singur vector.

pisa_predictions_ds2 <- pisa_predictions_ds |> apply(1, mean)

Distribuțiile de densitate sunt prezentate mai jos. Aceste distribuții seamănă mai mult cu distribuțiile pe care le așteptăm atunci când avem date echilibrate, chiar dacă am folosit toate observațiile pentru a obține aceste probabilități prezise.

ggplot(data.frame(Public = pisa_valid$Public, 
                  Prediction = pisa_predictions_ds2), 
       aes(x = Prediction, color = Public)) +
  geom_density()

Deși downsample() Funcția pare să abordeze problema valorilor montate de micșorare și oprire centrată, valorile de performanță ale modelului furnizate mai jos sugerează că nu a îmbunătățit performanța generală a predicțiilor modelului.

roc <- calculate_roc(predictions = pisa_predictions_ds2, 
                     observed = pisa_valid$Public)
plot(roc)

confusion_matrix(observed = pisa_valid$Public, 
                 predicted = pisa_predictions_ds2 > 0.5)

              predicted              
  observed        FALSE          TRUE
     FALSE   45 (3.44%)    41 (3.13%)
      TRUE 204 (15.60%) 1018 (77.83%)
Accuracy: 81.27%
Sensitivity: 52.33%
Specificity: 83.31%

Apendicele: rezumate ale modelului

Mai sus am înregistrat o medie a valorilor prevăzute pe toate modelele pentru a obține o singură predicție pentru fiecare observație din setul nostru de date de validare. Cu toate acestea, este posibil să grupați modele folosind mice::pool() Funcție pentru a obține un singur set de coeficienți de regresie. Tabelul de mai jos oferă coeficienții de regresie combinați de la downsample Funcționați împreună cu coeficienții din modelul de regresie logistică folosind toate datele.