Studiul datelor evaluării naţionale din 2016 (folosind R)
[1] Bacalaureat 2015 - de la forma microsoftizată, la R
[2] Studiul datelor examenului de bacalaureat din 2016 (folosind R)
Am lansat LibreOffice pe adresa furnizată de data.gov.ro pentru descărcarea setului de date conţinând "rezultatele anonimizate ale Evaluării Naționale 2016":
vb@Home:~/16_evn$ url_res="http://data.gov.ro/dataset/\
> 5eac43d1-b364-468f-bdff-18e668abfb20/resource/8d31a801-a430-47e3-96de-04bc4cc3ce55/\
> download/evnat20162.ods"
vb@Home:~/16_evn$ libreoffice $url_res
Interpretorul de comenzi bash permite separarea pe mai multe rânduri a unui text (încadrat de ghilimele), folosind '\' (caracterul "backslash") la capătul fiecărui rând; valoarea variabilei se referă apoi folosind '$'.
Am descărcat fişierul în format "ODS" (care - salvat pe disc şi apoi dezarhivat - are aproape 400 MB) nu pentru a-l folosi în vreun fel, ci doar pentru a obţine de sub LibreOffice fişierul în format "CSV" standard (cu acesta se va putea lucra).
Este drept că puteam descărca direct "evnat20162.csv" (în loc să folosesc intermediar, fişierul ".ods") - dar acesta are un format nestandard şi necesita transformări (v. [1] §1); de exemplu, foloseşte ',' şi nu '.' ca separator de zecimale (ceea ce poate că este convenabil pentru afişarea datelor, dar nu convine pentru a opera cu numere în programe); în plus, foloseşte codul specific Windows (doi octeţi, 0x0D0A) pentru separarea liniilor - ori de sub LibreOffice am obţinut direct codul specific Linux (un singur octet, 0x0A).
Inspectăm întâi fişierul CSV obţinut astfel (folosind vizualizatorul încorporat în GNOME Commander):
Trebuind să ne referim în programe la câmpurile respective, vom înlocui scrierea ALL CAPS şi vom simplifica denumirile (evitând în orice caz să folosim caracterul "spaţiu", pentru a nu fi obligaţi ulterior să încadrăm între ghilimele); vom renunţa la unele câmpuri (de exemplu, nu are de ce să ne intereseze codurile "SIIIR" - încât le vom înlocui cu denumirile judeţelor).
Lansăm R şi constituim un obiect de clasă 'data.frame' pentru setul de date din fişierul nostru:
vb@Home:~/16_evn$ R -q # --quiet (omite mesajul introductiv) > evn <- read.csv("evnat20162.csv") > str(evn) # inspectăm structura de date creată 'data.frame': 153673 obs. of 21 variables: # am adăugat aici şi indecşii câmpurilor [1] $ COD.UNIC.CANDIDAT : int 158870 136940 106543 152070 121592 106178 ... [2] $ SEX : Factor w/ 2 levels "F","M": 2 2 1 2 1 1 1 1 2 ... [3] $ MEDIU : Factor w/ 2 levels "RURAL","URBAN": 1 2 1 1 1 ... [4] $ COD.SIIIR : num 2.56e+09 3.26e+09 2.56e+09 2.26e+09 ... [5] $ STATUS.ROMANA : Factor w/ 3 levels "ABSENT","ELIMINAT",..: 3 3 ... [6] $ STATUS.LIMBA.MATERNA : Factor w/ 3 levels "-","ABSENT","PREZENT": 1 1 ... [7] $ STATUS.MATEMATICA : Factor w/ 3 levels "ABSENT","ELIMINAT",..: 3 3 ... [8] $ NOTA.ROMANA : num 3.2 5.1 6.05 6.4 7.4 9.5 4.1 9.9 8.9 8.65 ... [9] $ NOTA.LIMBA.MATERNA : num NA NA NA NA NA NA NA NA NA NA ... [10] $ NOTA.MATEMATICA : num 2.2 6.5 4.15 6.8 3.85 7.95 4.1 8 9.85 7.6 ... [11] $ CONTESTATIE.ROMANA : Factor w/ 2 levels "DA","NU": 2 2 2 2 2 2 2 2 1 ... [12] $ NOTA.CONTESTATIE.ROMANA : num NA NA NA NA NA NA NA NA 8.85 NA ... [13] $ CONTESTATIE.LIMBA.MATERNA : Factor w/ 2 levels "DA","NU": 2 2 2 2 2 2 2 2 2 ... [14] $ NOTA.CONTESTATIE.LB.MATERNA: num NA NA NA NA NA NA NA NA NA NA ... [15] $ CONTESTATIE.MATEMATICA : Factor w/ 2 levels "DA","NU": 2 2 2 2 2 2 2 2 2 ... [16] $ NOTA.CONTESTATIE.MATEMATICA: num NA NA NA NA NA NA NA NA NA NA ... [17] $ NOTA.FINALA.ROMANA : num 3.2 5.1 6.05 6.4 7.4 9.5 4.1 9.9 8.9 8.65 ... [18] $ NOTA.FINALA.LB.MATERNA : num NA NA NA NA NA NA NA NA NA NA ... [19] $ NOTA.FINALA.MATEMATICA : num 2.2 6.5 4.15 6.8 3.85 7.95 4.1 8 9.85 7.6 ... [20] $ MEDIA : num 2.7 5.8 5.1 6.6 5.62 8.72 4.1 8.95 9.37 ... [21] $ MEDIA.V.VIII : num 6.94 8.69 9.21 8 8.84 9.58 7.9 9.76 8.66 ...
Pecizăm că mai toate comenzile pe care le folosim aici au fost deja descrise în [1] şi [2]; în redarea secvenţelor de comenzi vom omite de obicei, promptul consolei R (caracterele iniţiale "> "); de altfel, secvenţele mai lungi le scriem într-un fişier "prog.R" separat (deci, fără prompt) şi folosim mecanismul asigurat de comanda source("prog.R").
Faţă de 2015 - v. [3] - numărul candidaţilor s-a micşorat cu aproape 10000. Ca şi în 2015, avem 21 de câmpuri - dar acum lipseşte câmpul "DATA NASTERE" şi în schimb, avem câmpul (cu siguranţă mai interesant) "MEDIA V-VIII".
Ne întrebăm dacă n-am putea renunţa la câmpurile "STATUS" (de indecşi 5-7):
summary(evn[5:7]) STATUS.ROMANA STATUS.LIMBA.MATERNA STATUS.MATEMATICA ABSENT : 4543 - :144528 ABSENT : 5020 ELIMINAT: 2 ABSENT : 238 ELIMINAT: 5 PREZENT :149128 PREZENT: 8907 PREZENT :148648
Au fost eliminaţi numai 7 elevi, iar proporţia totală a absenţilor este 9563 / 153673 ≈ 6.22%; nu zicem că ar fi mică, dar situaţia de absent (sau eliminat) este reflectată pe câmpul "MEDIA": în aceste cazuri, media lipseşte - fiind înlocuită de către R cu valoarea specifică 'NA'. Prin urmare, putem colecta indicii acestor câmpuri, pentru eliminare:
col_to_remove <- c(5:7) # colectează indicii coloanelor pe care le vom elimina în final
De asemenea, putem elimina cele trei câmpuri "CONTESTATIE" (cu valori "DA/NU"): pentru "DA" avem înscrisă nota corespunzătoare contestaţiei în câmpul asociat "NOTA CONTESTATIE" (iar altfel, acest câmp conţine valoarea 'NA'):
append(col_to_remove, c(11, 13, 15)) -> col_to_remove
Câmpul "COD SIIIR" conţine pentru fiecare candidat, codul SIIIR corespunzător acestuia; extragem codul judeţului (prima sau primele două cifre, după cum SIIIR are 9 sau 10 cifre) - folosind funcţia jud_from_siiir() din [1] - şi înfiinţăm coloana 'jud' (de clasă factor) conţinând denumirile judeţelor - folosind pentru aceasta datele din [1]:
evn$jud <- jud_from_siiir(evn[[4]]) # codul "NN" al judeţului candidatului (din COD.SIIIR) evn$jud <- as.factor(evn$jud) load("../15_bac/bac5.RData") # recuperează 'bac5', pentru nivelurile factorului 'jud' levels(evn$jud) <- levels(bac5$jud) rm(bac5) # elimină 'bac5' din sesiunea de lucru curentă (rm() = "remove")
Noul câmp 'jud' a fost adăugat la sfârşit, încât nu afectează indicii redaţi la început ai câmpurilor existente deja. Ştergem acum, coloanele colectate în acest scop (incluzând şi 'COD.SIIIR'):
evn[, c(col_to_remove, 4)] <- NULL names(evn) # numele coloanelor rămase [1] "COD.UNIC.CANDIDAT" "SEX" "MEDIU" [4] "NOTA.ROMANA" "NOTA.LIMBA.MATERNA" "NOTA.MATEMATICA" [7] "NOTA.CONTESTATIE.ROMANA" "NOTA.CONTESTATIE.LB.MATERNA" "NOTA.CONTESTATIE.MATEMATICA" [10] "NOTA.FINALA.ROMANA" "NOTA.FINALA.LB.MATERNA" "NOTA.FINALA.MATEMATICA" [13] "MEDIA" "MEDIA.V.VIII" "jud"
Ne bazăm de-acum pe noii indici de coloane, dar modificăm denumirile astfel:
names(evn) <- c("id", "Sex", "Mediu", "Română1", "Maternă1", "Matematică1", # notele iniţiale "Română2", "Maternă2", "Matematică2", # note contestaţii "Română", "Maternă", "Matematică", # note finale "Media", "Media8", "Jud")
Mai departe - rămâne să vedem ce spun datele colectate în setul de date "evn"; dacă ar fi să ne luăm numai după ce ştim din [1]-[3] - fără alte investigaţii - am zice de exemplu, că rezultatele examenului sunt mai bune pentru fete decât pentru băieţi şi la fel, pentru mediul "urban" faţă de cel "rural" (şi la fel pentru judeţul Cluj faţă de judeţul Ilfov, ş.a.m.d.).
1. Repartiţia candidaţilor după sex, mediu, sau judeţ
Faţă de 2015 - v. [3] - proporţia "F" = 51%, "M" = 49% se menţine, dar a scăzut cu câte 1% la "RURAL" şi a crescut cu câte 1% la "URBAN":
mediu_sex <- addmargins(table(evn[c('Mediu', 'Sex')])) nrcd <- mediu_sex[3,3] # total: 153673 candidaţi mediu_sex <- 100*round(mediu_sex / nrcd, 2) mediu_sex # 2016 Sex # valori procentuale (aproximative) Mediu F M Sum RURAL 23 22 45 # cu câte 1% mai puţin ca în 2015 URBAN 28 27 55 # cu câte 1% mai mult ca în 2015 Sum 51 49 100
Graficul redat mai sus l-am obţinut (v. şi [2] §5) prin:
require(lattice) barchart(mediu_sex, aspect=1, stack=FALSE, xlab="Procente", auto.key = list(columns=3, cex=0.8, size=1.5, points=FALSE, rectangles=TRUE), panel = function(...) { panel.barchart(...) panel.abline(v=seq(0,100,10), lwd=0.6, col="lightgrey") } )
Să vedem şi repartiţia pe judeţe, a candidaţilor:
nrcd <- nrow(evn) rep_jud <- as.data.frame(sort(table(evn$Jud) / nrcd)) require(ggplot2) require(scales) # scales::percent() (formatare procentuală) ggplot(rep_jud, aes(x = Var1, y = Freq)) + geom_bar(stat='identity', col="white", fill="lightgreen", width=0.85) + labs(x="", y="", title = paste0("Repartiţia pe judeţe (", nrcd, " candidaţi, 2016)")) + scale_y_continuous(labels = percent) + coord_flip()
Candidaţii din "M.Bucureşti" sunt cu aproape 4% până la 7% mai mulţi, decât în celelalte judeţe. Următoarea secvenţă produce diferenţele faţă de 2015 (dar de data aceasta nu ne convine calculul procentual):
load("evna.RData") # recuperăm datele din 2015 (163418 candidaţi) - v. [3] rep_jud5 <- as.data.frame(table(evna$Jud)) rep_jud6 <- as.data.frame(table(evn$Jud)) rep_jud6$difer <- rep_jud6$Freq - rep_jud5$Freq rm(list = c("evna", "rep_jud5")) # elimină obiectele de care nu mai avem nevoie print(data.frame(rep_jud6[1:14, ], rep_jud6[15:28, ], rep_jud6[29:42, ])) Var1 Freq difer Var1.1 Freq.1 difer.1 Var1.2 Freq.2 difer.2 Alba 2588 -256 Dâmboviţa 4290 -234 Prahova 5279 -284 Arad 3172 -191 Dolj 4571 -765 Satu-Mare 2565 -194 Argeş 4930 -149 Galaţi 3917 -595 Sălaj 1997 -136 Bacău 5226 -326 Gorj 3169 -214 Sibiu 2808 -142 Bihor 4704 -234 Harghita 2774 7 Suceava 6678 -466 Bistriţa-Năsăud 2557 -244 Hunedoara 3097 -155 Teleorman 2530 -327 Botoşani 3986 -121 Ialomiţa 2049 -193 Timiş 4477 -170 Braşov 3622 -88 Iaşi 6332 -414 Tulcea 1594 -26 Brăila 2023 -53 Ilfov 2257 -59 Vaslui 3593 -407 Buzău 3380 -302 Maramureş 4053 -231 Vâlcea 2970 -209 Caraş-Severin 2040 -242 Mehedinţi 2030 -80 Vrancea 2817 -144 Cluj 4102 -291 Mureş 4298 -190 M.Bucureşti 11983 -69 Constanţa 5555 -161 Neamţ 4372 -486 Călăraşi 2085 -288 Covasna 1532 -54 Olt 3280 -524 Giurgiu 2391 -38
Exceptând numai Harghita (unde a crescut cu 7), în toate judeţele numărul candidaţilor a scăzut faţă de 2015 (cel mai mult la Dolj (-765), Galaţi, Neamţ, Olt, Iaşi, Vaslui (-407)) - scăderea totală (obţinută prin sum(rep_jud6$difer)) fiind -9745.
Dacă am fi lucrat procentual - raportând numărul candidaţilor din fiecare judeţ la numărul total de candidaţi din 2016 şi respectiv, la cel din 2015 - atunci diferenţele procentelor respective ar avea suma 0 (reducându-se la diferenţa dintre suma procentelor judeţene din 2016 şi aceea din 2015, deci la 100% - 100%) şi ar fi fost greu de imaginat cum să le interpretăm corect.
2. Distribuţia mediilor de examen
Să verificăm întâi, notele şi mediile (amintindu-ne surpriza din [1] §9.5, când în coloana de note am găsit şi valori ca "-1" şi "-2" - însemnând calitatea de "absent", sau "eliminat"); aplicăm funcţia summary() fiecăreia dintre coloanele respective şi apoi transpunem matricea rezultată:
sapply(evn[, 4:14], FUN = summary) -> brev t(brev) -> brev cbind(brev, brev[, 3] - brev[, 4]) # adaugă coloana diferenţelor Median-Mean Min. 1st Qu. Median Mean 3rd Qu. Max. NA's Median - Mean Română1 1.00 5.40 7.25 6.870 8.60 10 4545 0.380 # note iniţiale Maternă1 1.00 6.35 7.80 7.415 8.80 10 144766 0.385 Matematică1 1.00 4.30 6.10 6.056 7.95 10 5025 0.044 Română2 1.20 7.05 8.15 7.803 8.85 10 141739 0.347 # contestaţii Maternă2 3.45 7.10 8.10 7.821 8.70 10 153281 0.279 Matematică2 1.00 5.90 7.30 6.971 8.20 10 144774 0.329 Română 1.00 5.40 7.25 6.884 8.65 10 4545 0.366 # note finale Maternă 1.00 6.40 7.80 7.432 8.80 10 144766 0.368 Matematică 1.00 4.30 6.10 6.058 7.95 10 5025 0.042 Media 1.00 5.00 6.70 6.509 8.22 10 5025 0.191 # medii finale Media8 5.13 7.89 8.79 8.628 9.48 10 5 0.162 # medii V-VIII
De data aceasta, toate notele sunt ceea ce trebuie să fie - valori din gama [1, 10]. Observăm că pentru toate cele 11 coloane, valoarea medie 'Mean' este mai mică decât valoarea mediană Median; aceasta înseamnă de regulă, că valorile respective au o distribuţie asimetrică, având o plajă spre stânga medianei ceva mai întinsă faţă de coada (mai abruptă) spre mediile mai mari ca mediana (vom vedea acest aspect pe graficele distribuţiilor, mai jos).
Am adăugat (prin cbind()) şi o coloană pe care am calculat diferenţele dintre valorile mediane şi cele medii - mult mai mari la română decât la matematică (0.38 faţă de 0.04), încât probabil că asimetria intuită mai sus este mai pronunţată la română, faţă de matematică.
Cei care nu au notă la una sau alta dintre probe (fiind absenţi sau eliminaţi) au fost contorizaţi în coloana 'NA's şi constatăm că numărul acestora la oricare probă este sub 10000. Pe de altă parte, observăm că valorile medianei şi mediei pentru notele iniţiale aproape coincid, cu cele pentru notele finale (în cazul fiecăreia dintre cele trei probe). Prin urmare - proporţia de contestaţii pe fiecare probă fiind relativ mică şi în plus, neavând în urma rezolvării contestaţiilor, modificări sensibile ale valorilor centrale - putem să ignorăm de-acum "contestaţiile" (considerându-le ca nerelevante din punct de vedere statistic).
Valorile medianei şi mediei pentru coloana 'Media' sunt cu 2 puncte mai mici decât cele pentru 'Media8' - dar nu prea putem interpreta acum, această diferenţă: valorile respective încep într-un caz de la 1.00, iar în celălalt de la 5.13 (iar pe de altă parte, valorile din 'Medie8' acoperă toate obiectele şcolare, nu numai pe cele două sau trei din care au rezultat valorile din 'Medie').
Există obiceiul de a categorisi notele pe intervale: "sub 5", [5, 6), [6, 7), ..., [9, 10]; dar dacă vrem să evidenţiem distribuţia lor, cel mai convenabil este să le considerăm ca atare, lăsând gruparea acestora în seama funcţiilor care produc histograme sau curbe de densitate.
Programul următor (scris într-un fişier "Media.R" şi lansat din consola R prin comanda source("Media.R")) trasează o histogramă a mediilor finale (grupându-le pe intervale de câte o zecime), curba densităţii estimate pentru acestea şi curba distribuţiei normale având ca parametri valoarea medie şi dispersia acestora:
opar <- par(mar=c(2,2,0,0)+0.1, cex.axis=0.9, xpd=FALSE) # histograma: hist(evn$Media, breaks = "FD", freq = FALSE, right = FALSE, col="cornsilk", border="lightgray", ylim=c(0, 0.21), main="") grid(lwd = 1.1) # estimarea densităţii: dst <- density(evn$Media, from = 1, to = 10, adjust = 0.75, na.rm=TRUE) lines(dst, col="blue", lwd=1.5) # densitatea distribuţiei normale N(mev, sev): mev <- mean(evn$Media, na.rm=TRUE) # media notelor (= 6.509183) sev <- sd(evn$Media, na.rm=TRUE) # dispersia (= 2.044166) lines(dst$x, dnorm(dst$x, mean = mev, sd = sev), lwd=1.5, col="magenta") # zona mev ± sev cuprinde 68.27% dintre valori (în distribuţia normală): msev <- c(mev-sev, mev, mev+sev) rmsev <- round(msev, 2) abline(v = msev, lwd=1, lty="dotted", col="red") # legende şi etichetări: legend("topleft", legend=c("densitatea estimată", expression("distribuţia normală N("*mu*","~sigma*")"), "", as.expression(bquote(mu ~"="~ .(rmsev[2]))), bquote(sigma ~"="~ .(rmsev[2]-rmsev[1]))), text.col=c("blue", "magenta", "white", "black", "black"), bty="n", cex=0.9, inset=c(0, 0.15)) text(c(mev-sev, mev, mev+sev), 0, labels = c(as.expression(bquote(mu ~ "-" ~ sigma ~ "=" ~ .(rmsev[1]))), bquote(mu == .(rmsev[2])), bquote(mu ~ "+" ~sigma == .(rmsev[3]))), cex=0.9, col="red", pos=1, offset=0.22) legend("top", legend="Densitatea distribuţiei mediilor, Evaluare Naţională 2016", bty="n", text.font=3, cex=1.2, xjust=0.5, yjust=0.5) par(opar) # reconstituie parametrii grafici standard
În maniera în care este folosită mai sus, funcţia expression() îşi interpretează argumentul ca "expresie matematică" (cu sensul şi sintaxa din TeX; numele 'mu', 'sigma' produc la afişarea expresiei literele greceşti astfel denumite, iar operatorii '*' şi '~' asigură concatenarea şi spaţierea).
Funcţia bquote() permite (prin construcţia .( )) combinarea unei expresii cu valori ale unor variabile. Dar există o problemă la folosirea ei într-un vector de expresii construit cu funcţia c(), cum avem mai sus în funcţia text() (anume, se returnează un obiect de clasă 'list', în loc de 'expression'); la stackoverflow am găsit ca soluţie convertirea prin as.expression() a fiecăreia dintre aplicările bquote() existente în construcţia c() respectivă (folosind eventual mecanismul sapply()) - dar în secvenţa de mai sus am demonstrat deja că este suficient să facem aceasta doar pentru una singură dintre componente.
Graficul obţinut mai sus arată că distribuţia mediilor finale nu poate fi considerată ca rezultând din cea normală. Avem o confirmare a acestui fapt şi printr-un calcul mai simplu: într-o distribuţie N(μ, σ) între limitele μ-σ şi μ+σ avem 68.27% dintre valorile respective, pe când în distribuţia mediilor finale avem numai 62.74%:
medii <- subset(evn, Media >= msev[1] & Media <= msev[3]) nrow(medii) / (nrow(evn) - 5025) # minus 5025 valori 'NA' (v. tabelul 'brev' mai sus) [1] 0.6273815 # 62.74% medii în intervalul [4.47, 8.55] ("normal" ar fi 68.27%)
Calcule similare celui redat mai sus ne arată că avem 18.45% medii sub 4.47 şi 18.81% medii peste 8.55.
3. 'F' versus 'M' (distribuţia mediilor condiţionată de sex)
Programul următor determină şi plotează densitatea mediilor pentru fiecare valoare a variabilei 'Sex':
dens <- density(evn$Media[evn$Sex == "F"], from = 1, to = 10, na.rm = TRUE) dens1 <- density(evn$Media[evn$Sex == "M"], from = 1, to = 10, na.rm = TRUE) plot(dens, col = "firebrick2", lwd = 1.5, cex.axis = 0.9, xlab = "", ylab = "", main = "Densitatea mediilor după sex", cex.main = 1, font.main = 3) lines(dens1, col="blue", lwd=1.5) legend("topleft", bty = "n", inset = c(0, 0.05), horiz=TRUE, legend=c("F", "M"), title="Sex", title.col="black", text.col= c("firebrick2", "blue")) grid()
În privinţa mediilor finale, avem un decalaj important între 'M' şi 'F'; proporţia mediilor sub 6 este sensibil mai mare, iar a acelora peste 6 este sensibil mai mică, la 'M' faţă de 'F'. Nu avem cum explica acest decalaj; observăm doar că distribuţiile redate mai sus seamănă (iarăşi inexplicabil) cu cele pe care le obţinem mai jos pentru notele de la probele de bază.
4. Română versus Matematică
Tabelul "brev" constituit mai sus indică deosebiri importante, între categoriile de note finale (dar precizăm că ignorăm 'Maternă' - proba suplimentară, susţinută de doar 8907 candidaţi, reprezentând abia 5.8% din total). Quartilele întâia şi a treia (în spatele cărora sunt 25% şi respectiv 75% dintre valori) sunt cu 1.1 şi respectiv 0.7 mai mici pentru 'Matematică', decât pentru 'Română' iar diferenţa 'Median - Mean' este mult mai mare la 'Română' (0.366 faţă de 0.042).
În secvenţa următoare preluăm quartilele celor două probe din tabelul "brev" (pentru a le marca ulterior pe graficul densităţilor), folosim mecanismul lapply() pentru a obţine "într-un singur pas" (spre deosebire de programul redat în §3) cele două densităţi (pentru coloanele 2 şi respectiv 10 din "evn", corespunzătoare notelor finale pentru cele două probe) şi apoi plotăm curbele densităţilor:
q13 <- c(brev[9, 2], brev[9, 5], brev[7, 2], brev[7, 5]) # quartilele 1 (25%) şi 3 (75%) cmr <- c("red", "red", "blue", "blue") # pentru linii şi etichete, la cele două probe dst <- lapply(evn[c(10, 12)], FUN = density, # estimarea densităţilor, pe cele două probe from = 1, to = 10, adjust = 1.75, na.rm = TRUE) plot(dst[[1]], col = "blue", lwd = 1.5, cex.axis = 0.9, xlab = "", ylab = "", main = "Română vs. Matematică (Ev. Nat. 2016)", cex.main = 1, font.main = 3) lines(dst[[2]], col = "red", lwd = 1.5) abline(v = q13, lty="dotted", col = cmr) # marchează quartilele fiecărei probe mtext(c("q1", "q3", "q1", "q3"), side=1, line=0, cex=0.8, at = q13, col = cmr) legend("topleft", bty = "n", cex = 0.9, inset = c(0, 0.05), legend=c("Română", "Matematică", "", "q1 = 25%, q3 = 75%"), text.col= c("blue", "red", "white", "black")) grid()
Asimetria intuită când am constituit tabelul 'brev' este într-adevăr, vizibil mai pronunţată la curba albastră (pentru 'Română') decât la cealaltă. Cele două curbe au un singur punct de intersecţie, pe la abscisa 6.25; observând zonele subântinse de cele două curbe în stânga şi în dreapta acestei abscise, găsim că proporţia notelor sub 6.25 este de aproape două ori mai mare la "Matematică" decât la "Română", iar pentru notele mai mari ca 6.25 raportul tocmai menţionat se inversează.
Dar trebuie menţionat că modificând parametrul "adjust" al funcţiei density() (mai sus am ales 'adjust=1.75'), putem obţine curbe "mai dese"; de exemplu, pentru valoarea 0.3 obţinem curbele:
Graficul rezultat astfel nu contrazice observaţiile anterioare (exceptând faptul că acum, pe lângă intersecţia de abscisă 6.25, apar şi alte puncte de intersecţie). "Zimţii" ca de pieptene, care au apărut pe cele două curbe reflectă faptul că notele respectă un anumit şablon - de exemplu, apar note de 5.50, 5.55, 5.60, 5.65 etc., dar nu şi note "intermediare" (ca 5.51, 5.52, etc.); ne putem edifica asupra acestui aspect folosind funcţia unique(), prin care obţinem reprezentanţii unici ai notelor:
sort(unique(evn$Matematică)) # curba colorată cu roşu, pe graficul de mai sus [1] 1.00 1.50 1.65 1.70 1.90 1.95 2.00 2.05 2.10 2.15 2.20 2.25 [13] 2.30 2.35 2.40 2.45 2.50 2.55 2.60 2.65 2.70 2.75 2.80 2.85 # ... ... ... [169] 9.95 10.00 # numai 4.57, 6.66, 7.88 sunt în afara şablonului "din 5 în 5 sutimi"
În cazul redat mai sus, exceptând intervalul [1, 1.50] şi încă trei valori (4.57, 6.66, 7.88) - secvenţa notelor respectă şablonul "din 5 în 5 sutimi" (conducând la forma zimţată a curbei).
5. Densităţile mediilor finale dintr-un interval dat, pe judeţe
Vrem să vizualizăm curbele care estimează densitatea distribuţiei mediilor finale cuprinse într-un interval dat, pe fiecare judeţ, în ordinea valorilor mediei judeţene pe acel interval. Selectăm din "evn" numai datele care ne interesează acum (valorile "Media" diferite de 'NA' şi câmpul "Jud") şi definim o funcţie în care mai întâi, restrângem datele corespunzător intervalului indicat şi apoi procedăm ca în [2] §9: creem un "factor" ordonat după mediile judeţene (pentru datele restrânse), ordonăm crescător vectorul din atributul 'scores' al acestuia şi instituim ca etichete pentru afişări alipirile de nume de judeţe şi respectiv, medii judeţene:
medii <- subset(evn, !is.na(Media), select = c("id", "Media", "Jud"))
dens_gap <- function(from = 1, to = 10) { medii <- subset(medii, Media >= from & Media <= to) rejud <- reorder(medii$Jud, medii$Media, FUN = mean, order = TRUE) medii$Jud <- factor(rejud, levels = levels(rejud), labels = paste(levels(rejud), round(sort(attr(rejud,"scores")), 2))) medii <- medii[order(medii$Jud), ] require(lattice) densityplot( ~ Media | Jud, # densitatea mediilor finale, pe fiecare judeţ în parte data = medii, from = from, to = to, n = 128, # pentru calculul densităţii type = "g", aspect = 0.5, # grid, înălţime/lăţime=0.5 par.strip.text = list(cex=0.8), # reduce mărimea de scriere a etichetelor layout = c(7, 6), # 6 linii a câte 7 panouri (42 judeţe) as.table = TRUE, # de la stânga la dreapta, de sus în jos main = list(label = bquote("Medii ["*.(from)*","~.(to)*"] în judeţe)"), font = 3, col = "cornsilk4") ) -> dst_list # listă cu toate datele şi setările necesare trasării graficului plot(dst_list) }
Apelul dens_gap() produce densităţile respective pentru intervalul implicit [1, 10]:
Pe măsură ce media judeţeană creşte - adică, vizând panourile de la stânga spre dreapta de sus în jos (sau mai simplu, parcurgând oricare coloană de panouri de sus în jos) - curba densităţii se aplatizează spre stânga (scade proporţia de medii mici, vizibil în special între gradaţiile 1 şi 4) şi se înalţă spre dreapta (culminând la "M.Bucureşti 7.41", unde înălţimea depăşeşte puţin 25%, la o abscisă apropiată de media 9). Observând înălţimea curbelor în punctul de abscisă 10, deducem că proporţia mediilor 10 în judeţele respective este cel mult 5% (în câteva cazuri fiind totuşi puţin peste 5%) şi doar la "M.Bucureşti 7.41" se înalţă până pe la 10%.
Să mai observăm că pe intervalul vizat [1, 10], mediile judeţene variază de-a lungul a două puncte (de la 'Giurgiu 5.39' la 'M.Bucureşti 7.41').
Apelul dens_gap(5) produce densităţile respective pentru intervalul [5, 10]:
Pentru intervalul vizat acum [5, 10], mediile judeţene variază de-a lungul unui punct - de la 'Harghita 6.95', la 'M.Bucureşti 7.93'.
Cu dens_gap(1, 5) obţinem densităţile respective pentru intervalul [1, 5], cu dens_gap(7.5, 9) pe cele corespunzătoare intervalului [7.50, 9], ş.a.m.d.
Într-un jurnal de ştiri am găsit acest aspect interesant: intervalul care cuprinde cele mai multe note este [8, 8.49]; pentru a releva acest aspect putem folosi funcţia cut(), clasificând mediile din jumătate în jumătate de punct:
medii$gap <- cut(medii$Media, breaks = c(seq(1, 9.5, by=0.5), 10.01), right=FALSE) as.data.frame(sort(table(medii$gap))) -> tst print(data.frame(tst[1:6, ], tst[7:12, ], tst[13:18, ]), row.names=FALSE) Var1 Freq Var1.1 Freq.1 Var1.2 Freq.2 [1,1.5) 823 [3.5,4) 6684 [6,6.5) 11083 [1.5,2) 1276 [4,4.5) 7968 [6.5,7) 11441 [2,2.5) 2276 [4.5,5) 9278 [7,7.5) 11835 [2.5,3) 3624 [5,5.5) 10302 [7.5,8) 12541 [3,3.5) 5121 [9,9.5) 10482 [8.5,9) 13290 [9.5,10] 6344 [5.5,6) 10758 [8,8.5) 13522
Constatăm că într-adevăr, intervalul [8, 8.49] conţine cel mai mare număr de note (anume, 13522 note), dintre toate intervalele "din 0.5 în 0.5 puncte" (urmat de [8.50, 9) cu 13290 note şi de [7.50, 8) cu 12541 note). Cu dens_gap(8, 8.49) putem vizualiza curbele de densitate corespunzătoare pe judeţe; menţionăm doar că în acest caz, plaja mediilor judeţene (calculate pentru acest interval de note) măsoară numai 5 sutimi - de la 'Covasna 8.21' la 'Sălaj 8.26'.
Putem folosi dens_gap() şi pentru mediile V-VIII: redefinim 'medii', selectând 'Media8' în loc de 'Media' (şi păstrăm filtrarea valorilor 'NA', fiindcă în tabelul 'brev' din §2 vedem că există 5 candidaţi la care lipseşte valoarea mediei V-VIII); renumim 'medii$Media8' cu 'medii$Media' (evitând necesitatea unor modificări în corpul funcţiei dens_gap()) şi apoi putem apela funcţia de mai sus:
medii <- subset(evn, !is.na(Media8), select = c("id", "Media8", "Jud")) names(medii)[2] <- "Media" # redenumeşte "Media8", pentru a invoca direct dens_gap() dens_gap(6)
Am optat pentru intervalul de note [6, 10] fiindcă dacă există medii generale V-VIII mai mici de 6, cu siguranţă ele sunt în număr foarte mic. Redăm doar primul si ultimul panou dintre cele 42 ale graficului (acestea fiind obţinute după includerea în funcţia densityplot() a parametrului subset = as.numeric(Jud) %in% c(1, 42), şi modificarea layout = c(2, 1) (pentru rând grafic cu două panouri)):
Vedem astfel că mediile judeţene V-VIII sunt cuprinse într-un interval de 60 de sutimi, de la 'Tulcea 8.33', la 'M.Bucureşti 8.94'. În general (observând toate cele 42 de panouri grafice), mediile judeţene de examen sunt cam cu 1-1.50 puncte mai mici decât mediile judeţene V-VIII.
vezi Cărţile mele (de programare)