3. Ontwikkeling model voor het vo

3.1 Inleiding

In dit hoofdstuk bespreken we de totstandkoming van het model waarmee voor het vo op leerlingniveau de kans op voortijdig schoolverlaten (vsv) kan worden geschat volgens de aanpak zoals beschreven in paragraaf 2.7. De eerste stap betreft de selectie van achtergrondkenmerken op basis van beschrijvende en bivariate analyses. In de tweede stap wordt een stepwise selectieprocedure toegepast op het cohort 2018/’19 en op het validatiecohort 2017/’18. Vervolgens bespreken we in de derde stap de resultaten van de kruisvalidaties. Bovendien onderzoeken we in de vierde stap of het model met multilevel-component toegevoegde waarde heeft. Tot slot presenteren we in de laatste paragraaf een vergelijking van de modelschattingen op basis van een DUO-cohort en een eerste beeld van de resultaten op het niveau van RMC-regio. 

3.2 Stap 1: Voorselectie kenmerken 

Om het effect van de achtergrondkenmerken, zoals beschreven in paragraaf 2.5, te onderzoeken is er een voorselectie gemaakt op basis van beschrijvende statistieken en bivariate analyses. Op basis van deze analyses zijn er keuzes gemaakt in de codering van de variabelen, zijn referentiecategorieën bepaald en is een keuze gemaakt bij (inhoudelijk) vergelijkbare kenmerken. 

De belangrijkste wijzigingen zijn hier uitgelicht:

• Ongeoorloofd verzuim: vsv kwam zeer beperkt voor in de groep leerlingen die doorverwezen was naar Halt wegens overtreding van de leerplicht. Dit is niet wenselijk bij het uitvoeren van regressieanalyses, omdat dit kan leiden tot onbetrouwbare en moeilijk te interpreteren resultaten. Het kenmerk is daarom niet meegenomen in de verdere analyses.
• Huishoudinkomen onder de lage inkomensgrens: Dit kenmerk overlapt sterk met de welvaartspercentielen, waardoor we, op advies van een inhoudelijk expert van het CBS, uiteindelijk ervoor hebben gekozen om de welvaartspercentielen mee te nemen in de verdere analyses in plaats van deze variabele. Dit kenmerk bevat daarnaast alleen inkomen, terwijl de welvaartspercentielen ook rekening houden met het vermogen. Hetzelfde geldt voor het kenmerk op buurtniveau. 
• Lage welvaart in het huishouden: Omdat de welvaartspercentielen in vijf categorieën een vollediger beeld geven dan alleen het laagste percentiel, is er voor gekozen de lage welvaart niet mee te nemen in de verdere analyses. 
• Aantal gewerkte uren: Het hebben van een bijbaan kan – theoretisch gezien – de kans op vsv zowel vergroten als verkleinen. Dit kenmerk is daarom lastig te interpreteren en werkelijke effecten kunnen daardoor vertekend zijn. Daarom is besloten dit kenmerk niet mee te nemen in verdere analyses. 
• Problematische schulden: Omdat vo-leerlingen meestal jonger zijn dan 18 jaar kunnen ze geen persoonlijke schulden hebben binnen de definitie zoals beschreven in bijlage 3. Daarom nemen we enkel de schulden op huishoudniveau mee voor het vo. Dit betekent dat de indeling van dit kenmerk is aangepast voor het vo naar twee categorieën: geen problematische schulden (0) en wel problematische schulden (1) in het huishouden. 

3.3 Stap 2: Stepwise analyse

Na de voorselectie van kenmerken is er een stepwise procedure toegepast om de kenmerken te selecteren die het beste model vormen. Dit is gedaan met de forward en backward search en een combinatie van beide. Uiteindelijk zijn de modelschattingen van de methoden vergeleken. Het doel is om een compact model over te houden, met voldoende verklaringskracht. Waar de forward methode kenmerken stapsgewijs toevoegt, verwijdert de backward methode deze stapsgewijs. De forward methode stopt als er geen modelverbetering meer optreedt op basis van de BIC-waarde. De backward methode doet hetzelfde bij het verwijderen van de kenmerken. De forward methode leverde in de analyses een compactere selectie kenmerken op dan de backward methode, waardoor deze stepwise methode als uitgangspunt is genomen voor de verdere selectie van kenmerken. Bovendien lagen zowel de BIC als de pseudo-R²-maten voor beide methoden dicht bij elkaar (forward: BIC = 61313, McKelveyZavoina R² = 0,36; backward: BIC = 61291, McKelveyZavoina R² = 0,37). 

De forward analyse resulteerde in de volgende selectie van kenmerken, in volgorde van toegevoegde waarde voor het model:

1. Leeftijd van de leerling
2. Onderwijssoort
3. Vertraging
4. Ouderlijke structuur
5. Welvaart
6. Verdachte van een misdrijf
7. Problematische schulden in huishouden
8. Geslacht 
9. Langdurige gezondheidsproblemen 
10. Psychosociale problemen moeder 

Na de selectie van het tiende kenmerk stopte de stepwise procedure, omdat er volgens het model geen extra verklaringskracht meer werd toegevoegd. 

De forward methode voegt telkens één kenmerk toe aan het model. Per stap in deze methode is de BIC uitgerekend om te bepalen in hoeverre er nog modelverbetering optreedt, zie figuur 3.3.1. Hierbij geldt dat een lagere BIC-waarde een betere modelkwaliteit betekent. 

De figuur laat zien dat hoe meer kenmerken er worden opgenomen, hoe lager de BIC wordt en des te beter het totale model de kans op vsv dus kan schatten. De grootste afname van de BIC ligt bij het model met de eerste twee kenmerken, wat te zien is aan de ‘knik’ in de figuur. De volgende knik in de daling ligt rond het model met vijf kenmerken, gevolgd door het model met zeven kenmerken. Daarna vlakt de lijn van de BIC af; de modelkwaliteit neemt bij het toevoegen van de laatste kenmerken nauwelijks nog toe.

Daarnaast is de stabiliteit van het model onderzocht door de forward stepwise procedure toe te passen op het validatiecohort (2017/’18). Er is dus opnieuw een stepwise procedure toegepast, waarbij opnieuw bepaald werd welke modelkenmerken relevant zijn voor dat cohort. Vervolgens konden de geselecteerde kenmerken en hun volgorde vergeleken worden tussen de twee cohorten. De resultaten worden weergegeven in tabel 3.3.2.

De selectie en volgorde van de eerste zeven kenmerken zijn in beide cohorten identiek. Daarna zijn er verschillen te zien tussen beide modellen. Dit komt overeen met het beeld in figuur 3.3.1. De verandering in de BIC-waarden was bij de laatste modellen zo klein dat het toevoegen van extra kenmerken aan het model willekeuriger wordt.   

3.4 Stap 3: Kruisvalidaties 

Na de stepwise analyse zijn er als derde stap kruisvalidaties uitgevoerd met de volgorde van kenmerken zoals beschreven in paragraaf 3.3. Tijdens deze analyse werd er eerst een leeg model geschat, om te onderzoeken wat de modelkwaliteit was zonder verklarende kenmerken. Daarna is er herhaaldelijk een nieuw model geschat waarbij telkens een extra kenmerk is toegevoegd op basis van de eerder vastgestelde volgorde. Uiteindelijk resulterend in het complete model met de tien kenmerken in het laatste model. 

De kruisvalidaties zijn geëvalueerd met behulp van de fitmaten zoals beschreven in de bijlages 4.1.2 en 4.1.3 en worden weergegeven in tabel 3.4.1. Om de recall-, precision- en F1-waarde te kunnen berekenen zijn leerlingen ingedeeld in twee categorieën: geen vsv (0), en wel vsv (1). Dit is gedaan met een grenswaarde, zoals beschreven in bijlage 4.1.3. Bij de daadwerkelijke toepassing van het model zullen we niet gaan werken met een classificatie van 0 of 1, maar met de daadwerkelijke kansen per leerling om vsv’er te worden (waarde tussen 0 en 1). Deze fitmaten geven dus vooral een globaal beeld van de modelkwaliteit en dienen gebruikt te worden voor onderlinge modelvergelijkingen. Dit geldt niet voor de (relatieve) entropie en gemiddelde R².

In de tabel zien we dat de entropie-waarde afneemt, naarmate het model uitgebreider wordt. Dit betekent dat hoe uitgebreider het model, des te beter het model wordt in het schatten van de kans op vsv. De relatieve entropie geeft de relatieve verbetering ten opzichte van het lege model weer. De toename blijft oplopen, maar vlakt af rond het model met zeven kenmerken. 

De recall-waarde ligt tussen de 0,73 en 0,84. Dit betekent dat het model de leerlingen die werkelijk vsv’er worden vaak als zodanig classificeert. De precision ligt echter tussen de 0,03 en 0,06 in, wat een relatief lage waarde is,  maar die wel in lijn is met het lage aandeel vsv’ers in de populatie (0,8% op het vo, zie ook paragraaf 2.4). De precision laat zien dat de modellen met achtergrondkenmerken het een stuk beter doen dan een leeg model. Stel de achtergrondkenmerken worden niet meegenomen en iedereen zou als vsv’er geclassificeerd worden, dan zou de precision gelijk zijn aan het aandeel vsv’ers in de vo populatie, dus 0,008. Dan is een precision van 0,03 tot 0,06 weer een relatieve verbetering. De F1-waarde vat bovenstaande resultaten van de recall- en precision-waarde samen.

De gemiddelde McKelveyZavoina R² over de vijf kruisvalidaties ligt tussen de 0,32 en 0,41, maar ligt voor de meeste modelvarianten dicht bij elkaar. Deze R² variant moet met voorzichtigheid worden geïnterpreteerd en de grootte van het effect is daarbij ook context-afhankelijk. We gebruiken de R² in de kruisvalidaties dan ook voornamelijk om modelvergelijkingen te maken. Daaraan zien we dat model twee de hoogste waarde heeft. Voor de overige modellen ligt de waarde wat lager en loopt de waarde op naarmate het model uitgebreider wordt. Omdat een model met twee kenmerken wel erg beperkt is en weinig inzicht geeft in de achtergrondkenmerken die bijdragen aan het risico op vsv, gaat de voorkeur uit om naar de uitgebreidere varianten te kijken.

3.5 Conclusie modelselectie 

Op basis van bovenstaande resultaten is een definitief voorkeursmodel voor het vo geselecteerd. Het is daarbij van belang om een model te selecteren op basis van de modelkwaliteit en stabiliteit. Daarnaast wil je een zo informatief mogelijk model, dat toch eenvoudig, transparant en goed uit te leggen blijft. Op basis van deze criteria is het model met zeven kenmerken gekozen. Ten eerste zien we in figuur 3.3.1 dat de BIC-waarde niet veel verder afneemt na zeven kenmerken. Ten tweede zien we dat de selectie en volgorde van de eerste zeven kenmerken stabiel blijven tussen de twee onderzochte cohorten. Ook blijkt uit stap 3 dat de modelkwaliteit bij zeven kenmerken vergelijkbaar blijft met modellen die uitgebreider zijn. Het verder uitbreiden van het model heeft dan ook weinig toegevoegde waarde. Indien een verdere versimpeling van het model door OCW gewenst is en in fase 2 vergelijkbare resultaten oplevert met het voorkeursmodel, zou eventueel gekozen kunnen worden voor het model met vijf kenmerken (ook wel het “back-upmodel” genoemd). De coëfficiënten en Odds Ratio’s (OR) behorende tot het model met vijf en zeven kenmerken worden weergegeven in bijlage 5. 

Met behulp van de gegevens in figuur 3.5.1 kan een voorbeeld gegeven worden van de toepassing van de odds ratio’s voor een fictieve leerling. Deze leerling heeft bepaalde kenmerken, zoals een leeftijd van 19 jaar of ouder en een lage welvaart. Elke categorie waarin een leerling valt, heeft een odds ratio ten opzichte van de referentiecategorie. Door vervolgens deze met elkaar te vermenigvuldigen, komen we uit op de odds op vsv van de betreffende leerling. Zie paragraaf 2.7 voor een uitgebreidere uitleg van odds en odds ratio’s. Volgens de figuur is de kans op vsv voor die leerling 1,303 keer zo groot als de kans op geen vsv. Bij deze kansverhouding hoort een geschatte kans op vsv van 0,566. De kans dat een leerling met deze combinatie van achtergrondkenmerken vsv’er wordt is dus 56,6%. 

3.6 Stap 4: Uitbreiding met multilevel-component 

Leerlingen zijn geclusterd binnen scholen en scholen zijn weer geclusterd binnen RMC-regio’s. Zoals eerder is opgemerkt in paragraaf 2.7 zou het kunnen dat twee jongeren die op dezelfde school zitten of binnen dezelfde regio naar school gaan relatief vaker allebei wel of allebei geen vsv’er worden dan twee jongeren op verschillende scholen of uit verschillende regio’s. Om te onderzoeken hoe sterk dit clustereffect is en of hier in de analyses rekening mee gehouden dient te worden, hebben we in stap 4 multilevel modellen geschat. 

Allereerst is gekeken naar clustering op het hoogste niveau, te weten RMC-regio. In een model waarin alleen een random intercept op RMC-regioniveau is opgenomen was de Median Odds Ratio (MOR) 1,30 (95% betrouwbaarheidsinterval (BI): 1,23-1,41). Dat wil zeggen dat wanneer een leerling verhuist van een RMC-regio met een lagere odds op vsv naar een RMC-regio met een hogere odds op vsv, de mediane odds op vsv 1,3 keer zo groot zijn. Aangezien de MOR een odds ratio is, kan hij ook direct vergeleken worden met de andere odds ratio’s van de variabelen in het model. In verhouding is dit effect van RMC-regio dusdanig klein, dat wij hebben besloten hier in de modelontwikkeling geen rekening mee te houden. Aanvullend is wel nog onderzocht of deze MOR nog kleiner werd na het toevoegen van de zeven geselecteerde verklarende variabelen, wat inderdaad het geval was. 

Vervolgens is de clustering binnen scholen onderzocht. In het model met alleen een random intercept op schoolniveau was de MOR 2,91 (95% BI:2,72-3,13). Het schooleffect is dus aanzienlijk groter dan het effect van RMC-regio, d.w.z. er is een veel sterkere clustering binnen scholen dan binnen RMC-regio’s. Vervolgens zijn aan het model met het random intercept op schoolniveau de zeven geselecteerde verklarende variabelen toegevoegd. Na toevoeging van deze variabelen was de MOR nog maar 1,45 (95% BI:1,40-1,51). De variabelen konden de clustering binnen de scholen dus al voor een groot deel verklaren. Verder hebben we de geschatte coëfficiënten van de verklarende variabelen vergeleken tussen het model met een random intercept op schoolniveau en het model zonder multilevel-component. Er waren geen duidelijke verschillen zichtbaar tussen deze geschatte coëfficiënten, wat betekent dat de geschatte kansen op vsv uit beide modellen in de praktijk dicht bij elkaar zouden liggen. Daarom is besloten dat in de modelontwikkeling geen rekening hoeft te worden gehouden met de clustering binnen scholen. De conclusie is dat een multilevel component niet nodig is voor het vo model. 

Bij het validatiecohort vonden wij vergelijkbare resultaten. In het model met alleen een random intercept op het niveau van RMC-regio was de MOR 1,31 (95% BI:1,23-1,42). In het model met alleen een random intercept op schoolniveau was de MOR 2,96 (95% BI:2,77-3,19) en na toevoeging van de zeven verklarende variabelen nog maar 1,47 (95% BI:1,41-1,53). 

Hoewel de gevonden clusteringseffecten klein zijn, zijn ze wel statistisch significant. Hierbij moet worden bedacht dat de onderzoekspopulatie een groot aantal waarnemingen bevat waardoor de kans op statistisch significante resultaten wordt vergroot.

3.7 Extra analyses

Tot slot zijn er twee aanvullende analyses uitgevoerd. Ten eerste hebben we, om de resultaten van het model te valideren, ook een vergelijkbare analyse uitgevoerd op basis van DUO-data. Zoals beschreven in de inleiding van dit rapport, hanteert DUO een andere definitie voor zowel vsv als de populatie dan het CBS. We willen deze data daarom vooral gebruiken om te zien of een model op basis van DUO-data tot een vergelijkbare selectie en volgorde van kenmerken komt als met CBS-data. Om de resultaten te valideren, is de stepwise procedure daarom opnieuw toegepast. 

In tabel 3.7.1 worden de resultaten vergeleken. De stepwise procedure resulteert voor de DUO-data in een selectie van twaalf kenmerken, in vergelijking met tien op basis van de CBS-data. De top negen van kenmerken is vergelijkbaar tussen beide modellen. Wel verschuift de volgorde onderling wat. Daarnaast komt het tiende kenmerk uit de CBS-data wel terug op plaats twaalf in de DUO-data en zijn de psychosociale problemen van de leerling en of ouders wanbetaler van de ZVW-premie er als kenmerk bij gekomen. 

Zoals in paragraaf 3.5 besproken is het model met zeven kenmerken ons voorkeursmodel. De eerste zeven kenmerken komen overeen tussen beide datasets; er is slechts lichte verschuiving in de volgorde van de kenmerken. Het is daarbij belangrijker dat dezelfde kenmerken geselecteerd worden dan dat dit in exact dezelfde volgorde gebeurt. Op basis van deze analyse concluderen wij dan ook dat beide datasets tot dezelfde selectie van kenmerken leiden. 

Ten tweede hebben we op basis van het eindmodel met zeven kenmerken voor de totale populatie op het vo een simpele aggregatie uitgevoerd van de geschatte kansen op vsv. Dit betekent dat de geschatte kansen van alle leerlingen in een bepaalde RMC-regio bij elkaar zijn opgeteld.  Hiermee krijgen we al een eerste indicatie van de verschillen tussen het werkelijk en geschatte aantal vsv’ers en hoe dit tussen regio’s verschilt. In fase 2 van dit onderzoek zal deze aggregatie naar RMC-niveau uitvoerig onderzocht worden. Hierbij zullen ook verschillende beleidskeuzes door het ministerie van OCW gemaakt moeten worden. 

Om een eerste indicatie te geven van de samenhang tussen het werkelijke en geschatte aantal vsv’ers per RMC-regio op het vo, hebben we de Pearson correlatiecoëfficiënt uitgerekend. Deze geeft een sterke samenhang aan tussen het werkelijke en geschatte aantal vsv’ers met \( \rho  \)= 0,974. Omdat aantallen afhankelijk zijn van het aantal leerlingen op RMC-regio, hebben we in figuur 3.7.2 de aantallen uitgedrukt als percentage van het totaal aantal leerlingen per regio. Het totaal aantal leerlingen is gedefinieerd zoals in dit gehele rapport: leerlingen zonder startkwalificatie. Daarbij is te zien dat er wat minder spreiding is tussen regio’s in de geschatte waarden t.o.v. in de werkelijke waarden. Ook is te zien dat er vaker een overschatting is van de werkelijke waarde (onder de blauwe lijn) dan een onderschatting (boven de blauwe lijn). Wat de mogelijke oorzaken hiervan zijn, wordt in fase 2 van het onderzoekstraject onderzocht.