Gli Obbiettivi

L’obbiettivo principale di questo articolo è la ricerca di cambiamenti riguardanti i ragazzi che entrano per la prima volta nel mondo universitario, di conseguenza la variabile di interesse che abbiamo scelto per l’analisi è il numero di immatricolati, definito come:

• Sono considerati immatricolati solo gli studenti che si sono iscritti per la prima volta nella loro vita ad un corso di laurea. Questo numero rappresenta effettivamente quanti nuovi studenti entrano nel mondo universitario dal momento che non considera evenutali cambi di corso o di ateneo.

L’analisi partirà da una visione del fenomeno a livello nazionale, per poi stringersi fino a porre l’attenzione sulle singole regioni, inoltre nel corso dell’articolo cercheremo di rispondere ad alcune domande che ci siamo posti:

1. Qual è stata l’evoluzione delle immatricolazioni negli ultimi 50 anni? E’ possibile individuare delle variabili che possano spiegare questa dinamica?
2. Qual è stato il reale impatto del Covid19 sulle immatricolazioni a livello nazionale?
3. E’ possibile costruire un modello per prevedere il numero di immatricolati?
4. Come sono cambiate le preferenze degli studenti con il Covid19? Quali sono le regioni che attirano più studenti? E quali quelle da cui ne vanno via di piu?
5. Il Covid19 ha disincentivato gli spostamenti fuori regione per ragioni di studio?

Analisi Esplorativa delle Immatricolazioni

Per iniziare, abbiamo deciso di analizzare l’andamento della serie storica delle immatricolazioni in Italia in modo da monitorarne l’evoluzione negli ultimi anni. Inizialmente avevamo sfruttato solo i dati reperibili dal sito del Miur, tuttavia essendo dati riferiti agli anni che vanno dal 1998 al 2020 non riuscivamo a cogliere una visione generale del fenomeno. Di conseguenza per estendere la nostra comprensione, abbiamo deciso di ricostruire delle serie storiche aggiungendo anche i dati disponibili sul sito dell’ Istat, ottenendo cosi delle informazioni più dettagliate che spaziano dal 1963 fino ad oggi.

spot_points = tibble(anno = c(1969, 1981, 1993, 1999, 2003, 2013, 2020), imm_totali = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$imm_totali, label = c("1969", 
    "1981", "1993", "1999", "2003", "2013", "2020"), nudge_y = c(-10, -10, 10, -10, 
    10, -10, 10), nudge_x = c(1.5, 0, 0, 0, 0, 0, 0))

ggplot(df_model) + geom_line(aes(anno, imm_totali/1000, group = 1), col = "black", 
    size = 1) + geom_point(data = spot_points, aes(anno, imm_totali/1000, group = 1), 
    col = "red", size = 1.5) + geom_text(data = spot_points, aes(anno + nudge_x, 
    imm_totali/1000 + nudge_y, label = label)) + ggthemes::theme_solarized() + scale_fill_manual(values = c("#F58772", 
    "#B7F28D")) + theme(axis.text.x = element_text(angle = 0, face = "bold"), axis.text.y = element_text(face = "bold"), 
    axis.title = element_text(face = "bold"), plot.title = element_text(face = "bold"), 
    plot.subtitle = element_text(face = "italic"), plot.caption = element_text(face = "italic"), 
    panel.grid.major.x = element_blank(), panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(x = "Anno Accademico", 
    y = "Immatricolati (migliaia)", title = "Immatricolati in Italia", subtitle = "Anni Accademici 1963/1964 - 2020/2021", 
    caption = "Fonte: Istat, Miur Open Data") + scale_x_continuous(breaks = seq(1965, 
    2020, 5)) + annotate("rect", ymin = c(230), ymax = c(380), xmin = c(1991), xmax = c(2022), 
    alpha = 0.05, color = "#F58772", fill = "#F58779")

Prime considerazioni

• Il punto di partenza è il 1969 anno di emanazione della legge Codignola che ha costituito un punto di svolta: fino ad allora infatti l’accesso all’università era consentito solo agli studenti provenienti dai licei; da quel momento in poi, venendo a mancare questo sbarramento, l’accesso all’università era permesso, almeno teoricamente, a tutti gli studenti. Il grafico mostra senz’altro un’andamento crescente nel periodo successivo, ma non elevato come ci si potrebbe ragionevolmente aspettare.

• Periodo di stagnazione dal 1975-85: nel decennio tra gli anni settanta fino alla prima metà degli anni ottanta sembra esserci stato un periodo di stallo: uno dei fatti storici caratterizzanti questo periodo è la Crisi Energetica e i disagi alla mobilità (e non solo) che ne conseguirono. Si potrebbe pensare che in questo periodo ci siano stati ulteriori fattori esterni che abbiano portato la situazione a essere leggeremente decrescente o stabile, registrando il minimo di immatricolati nel 1981.

Negli ultimi 30 anni l’andamento potrebbe sembrare ciclico, tuttavia:

• Possiamo evidenziare due massimi: uno nel 1993 e l’altro nel 2003. Il massimo del 2003 potrebbe essere ricondotto agli effetti della legge Berlinguer dei primi anni duemila, in cui venne operata una distinzione tra lauree triennali, lauree magistrali e magistrali a ciclo unico. L’effetto sembrerebbe essere stato particolarmente positivo, ma di carattere transitorio avendo interessato solo il periodo dal 1999 al 2003.

• Minimo del 2013: questo valore risulta interessante nella nostra analisi, considerando che dal 2003 al 2013 gli immatricolati sono in calo costante. In questi anni sono state emanate diverse leggi con lo scopo di aumentare l’indipendenza economica ed amministrativa degli atenei. Un altro dato rilevate è stata la crisi economica del 2007, la quale ha provocato una situazione di particolare disagio in Italia (toccando vari ambiti) fino a raggiungere il proprio apice proprio nel 2012-13, con l’aumento dello spread, le dimissioni dell’allora Presidente del Consiglio Berlusconi e la costituzione del governo tecnico Monti. Senza dubbio in quel periodo in Italia vigeva un clima di profonda incertezza sul futuro e, probabilmente, si potrebbe ipotizzare anche sull’utilità di conseguire una laurea.

La Demografia dei ragazzi in Italia

Dopo aver esaminato l’andamento generale della nostra variabile di interesse, vogliamo cercare di capire se effettivamente in italia il numero medio di studenti che decide di intraprendere gli studi sia aumentato.

Dare il giusto rilievo a questo aspetto risulta particolarmente difficoltoso, infatti se volessimo rapportare il numero di immatricolati al numero di ragazzi tra i 19 e i 25 anni, potremmo giungere a deduzioni errate, in quanto il risultato sarebbe influenzato dalla (de)crescita demografica, che in Italia ha seguito un calo repentino negli ultimi venti anni.

spot_points = tibble(anno = c(1969, 1981, 1993, 1999, 2003, 2013, 2020), pop1925 = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$pop1925, pop1419 = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$pop1419, label = c("1969", 
    "1981", "1993", "1999", "2003", "2013", "2020"), nudge_y = c(-50, 30, 2, 50, 
    170, 100, 80), nudge_x = c(1.5, 2, -2, 1.5, 1, 0, 0))
ggplot(df_model) + geom_line(aes(anno, pop1925/1000, group = 1), col = "black", size = 1) + 
    geom_point(data = spot_points, aes(anno, pop1925/1000, group = 1), col = "red", 
        size = 1.5) + geom_text(data = spot_points, aes(anno + nudge_x, pop1925/1000 + 
    nudge_y, label = label)) + xlab("Anno Accademico") + ylab("Valore (in migliaia)") + 
    ggtitle("Popolazione di ragazzi tra i 19 e i 25 anni", subtitle = "Italia, Anni 1963 - 2021") + 
    ggthemes::theme_solarized() + scale_fill_manual(values = c("#F58772", "#B7F28D")) + 
    theme(axis.text.x = element_text(angle = 0, face = "bold"), axis.text.y = element_text(face = "bold"), 
        axis.title = element_text(face = "bold"), plot.title = element_text(face = "bold"), 
        plot.subtitle = element_text(face = "italic"), plot.caption = element_text(face = "italic"), 
        panel.grid.major.x = element_blank(), panel.grid.minor.x = element_blank(), 
        panel.grid.major.y = element_line(colour = "grey60", linetype = "dotted"), 
        legend.position = "none") + labs(caption = "Fonte: Istat") + scale_x_continuous(breaks = c(1965, 
    1970, 1975, 1980, 1985, 1990, 1995, 2000, 2005, 2010, 2015, 2020))

Osservando il grafico, a partire dagli anni ottanta è possibile notare una drastica dimiuzione della popolazione tra i 19 e i 25 anni fino a raggiungere il minimo nel 2003.

Negli ultimi venti anni la situazione, seppur con dei livelli nettamente inferiori rispetto al passato, sembrerebbe essere caratterizzata da una certa stabilità. Pur non volendo entrare nel merito di previsioni demografiche, dato che sarebbe un argomento molto vasto, è comunque possibile esprimere delle preoccupazioni rigurardanti il futuro delle università italiane, che potranno contare su un numero sempre minore di giovani di nazionalità italiana.

Un indicatore per le Immatricolazioni

Dopo aver fornito un quadro generale, abbiamo deciso di costruire un indicatore che potesse riuscire a rappresentare un aumento di immatricolazioni a parità di altre condizioni: il più esplicativo che abbiamo individuato è stato il rapporto tra il numero di immatricolati e il numero di ragazzi tra i 19-25 anni, (indicato con il nome di Pop1925) al quale sono stati sottratti i ragazzi che sono iscritti all’università per rappresentare i giovani che non studiano. L’indicatore, che abbiamo chiamato semplicemente con il nome di Rapporto Immatricolati, (RI) è stato calcolato come segue:

\[ RI = \frac{Immatricolati}{Pop1925 - Iscritti}*100\]

La finalità di questo indicatore è quella di evidenziare la percentuale di immatricolati in rapporto alle persone attualmente non iscritte all’iniversità.

spot_points = tibble(anno = c(1969, 1981, 1993, 1999, 2003, 2013, 2020), imm_totali = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$imm_totali, isc_totali = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$isc_totali, pop1925 = filter(df_model, 
    anno %in% c(1969, 1981, 1993, 1999, 2003, 2013, 2020))$pop1925, label = c("1969", 
    "1981", "1993", "1999", "2003", "2013", "2020"), nudge_y = c(1.5, 1.1, 0, -1, 
    1, -1, 1), nudge_x = c(0, 0, -1.5, 0, 0, 0, 0))
ggplot(df_model) + geom_line(aes(anno, imm_totali/(pop1925 - isc_totali) * 100, group = 1), 
    col = "black", size = 1) + geom_point(data = spot_points, aes(anno, imm_totali/(pop1925 - 
    isc_totali) * 100, group = 1), col = "red", size = 1.5) + geom_text(data = spot_points, 
    aes(anno + nudge_x, imm_totali/(pop1925 - isc_totali) * 100 + nudge_y, label = label)) + 
    xlab("Anno Accademico") + ylab("Rapp (%)") + ggtitle("Rapporto Immatricolati*", 
    subtitle = "Popolazione di 19-25 anni, Anni Accademici 1963/1964 - 2020/2021") + 
    ggthemes::theme_solarized() + scale_fill_manual(values = c("#F58772", "#B7F28D")) + 
    theme(axis.text.x = element_text(angle = 0, face = "bold"), axis.text.y = element_text(face = "bold"), 
        axis.title = element_text(face = "bold"), plot.title = element_text(face = "bold"), 
        plot.subtitle = element_text(face = "italic"), plot.caption = element_text(face = "italic"), 
        panel.grid.major.x = element_blank(), panel.grid.minor.x = element_blank(), 
        panel.grid.major.y = element_line(colour = "grey60", linetype = "dotted"), 
        legend.position = "none") + labs(caption = "Fonte: (*) Rielaborazione Dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = c(1965, 1970, 1975, 1980, 1985, 1990, 1995, 2000, 
        2005, 2010, 2015, 2020))

Confrontando il grafico delle immatricolazioni totali con il rapporto calcolato vengono evidenziati tre aspetti interessanti:

1. Tra il 1969 e il 1981, nonostante le immatricolazioni sembravano essere aumentate in valore assoluto, in realtà non è stato così: una volta rapporto il valore al numero di giovani potenziali sembrerebbe quasi stagnante, ad indicare che quell’aumento di immatricolati in realtà potrebbe essere rincoducibile a delle variazioni demografiche.

2. L’ aumento di immatricolazioni alla fine degli anni ottanta e quello dei primi anni duemila hanno portato ad un aumento in percentuale dei ragazzi che entrano all’università, dato estremamente positivo, ma transitorio.

3. L’ipotesi di un’andamento ciclico degli immatricolati negli ultimi 20 anni sembra venire rafforzata dalla similitudine dei due grafici in quel perido. Tuttavia dato che il ciclo riguarda anche il rapporto calcolato possiamo fare alcune considerazioni riguardo la sua natura. Potremmo ipotizzare che la diminuzione dal 2003 al 2013, possa essere ricondotta a fattori esterni, ad esempio legati alla Crisi Economica, come ad esempio la perdita di speranza nel futuro dal mmento che a una diminuzione di immatricolazioni corrisponde ad una diminuzione della percentuale di giovani che decidono di studiare. Dal 2013, la situazione sembra ritornare, anche se molto lentamente, ai livelli precedenti.

Ricerca delle variabili

Per verificare statisticamente se ci sia stato o meno un effetto Covid19, abbiamo deciso di costrurire un modello econometrico in grado di spiegare l’andamento degli immatricolati nel tempo, e di ottenere così delle informazioni sull’evoluzione del fenomeno, in modo da poter fare affermazioni su un possibile effetto Covid19.

Per costruire un modello in grado di spiegare il nostro fenomeno occorre ricercare delle variabili correlate con il fenomeno oggetto di interesse. Dal momento che abbiamo deciso di creare un modello occorre fare particolare attenzione ai concetto di correlazione e causalità di cui non riporteremo definizioni matematiche o statistiche, ma per spiegarlo cercheremo di usare questa foto:

Correlazione e Causalità

• La correlazione indica semplicemente che esiste un legame tra due fenomeni, può essere misurata matematicamente e ci può aiutare a formulare modelli matematici e a fare previsioni, tuttavia non possiamo fare affermazioni circa la causalità tra i due fenomeni.

• La causalità: indica che un fenomeno causa un altro, e siamo sicuri che questa influenza vada da una direzione all’altra.

Gli Immatricolati e i Diplomati

Una variabile correlata con gli immatricolati è il numero di diplomati: è ragionevole pensare che una variazione di diplomati quest’anno si rifletta sulle immatricolazioni dell’anno prossimo dato che sono presenti piu giovani potenziali che hanno accesso al mondo universitario.

Tuttavia, dopo le considerazioni fatte, vogliamo ribadire che anche se il modello dovesse funzionare, non potremmo mai essere certi che effettivamente una variazione nei diplomati abbia un effetto causale (ovvero che ne sia la causa) sulla variazione di immatricolati.

Osserviamo ora l’andamento dei diplomati dal 1963 ad oggi:

spot_points = tibble(anno = c(1986, 1992, 1999), diplomati_totali = c(331861, 490871, 
    444367), label = c("1986", "1992", "1999"), nudge_y = c(-20, 20, -20))

ggplot(df_model) + geom_line(aes(anno, diplomati_totali/1000, group = 1), col = "black", 
    size = 1) + geom_point(data = spot_points, aes(anno, diplomati_totali/1000, group = 1), 
    col = "red", size = 1.5) + geom_text(data = spot_points, aes(anno, diplomati_totali/1000 + 
    nudge_y, label = label)) + xlab("Anno Accademico") + ylab("Diplomati (migliaia)") + 
    ggtitle("Diplomati in Italia", subtitle = "Anni Accademici 1963/1964 - 2020/2021") + 
    ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, face = "bold"), 
    axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Istat") + 
    scale_x_continuous(breaks = c(1965, 1970, 1975, 1980, 1985, 1990, 1995, 2000, 
        2005, 2010, 2015, 2020))

• L’andamento di diplomati sembra in crescita costante fino al 1992, il che sembrerebbe un dato positivo.

• Negli ultimi trent’anni possiamo osservare un numero pressochè costante di diplomati, probailmente frutto del fatto che l’istruzione è diventata obbligatoria almeno fino ai 16 anni di età.

ggplot(df_model[-c(1, 2), ], aes(lag1_var_diplomati, var_immatricolati)) + geom_line(size = 0.9) + 
    geom_smooth(method = "lm", se = FALSE, col = "red", formula = y ~ x) + xlab("Diplomati (%)") + 
    ylab("Immatricolati (%)") + ggtitle("Relazione tra Diplomati e Immatricolati", 
    subtitle = "Relazione tra la variazione di immatricolati nell'anno corrente (y) e la variazione \n di immatricolati nell'anno passato (x). Il periodo considerato va dal 1963 al 2020.") + 
    ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, face = "bold"), 
    axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Rielaborazione dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = c(-5, 0, 5, 10, 15)) + scale_y_continuous(breaks = c(-10, 
    0, 10, 20))

Il grafico conferma in parte il nostro ragionamento: sembrerebbe sussistere una relazione lineare positiva tra le due variabili, seppur caratterizzata da un’elevata variabilità. Ora andremo a rispondere alle seguenti domande:

• Una variazione nel numero di diplomati provoca una variazione nel numero di immatricolati?
• Esiste un’effetto Covid19?

Il passo successivo sarà verificare statisticamente se ci siano stati degli effetti particolari durante il 2021. Iniziamo con il dire che la previsione riportata nell’articolo era completamente errata: a livello nazionale si è verificato un aumento delle immatricolazioni di 13876 studenti, ovvero il 4,43% in più rispetto al 2019/20.

Formulazione del modello

Nella costruzione del modello vogliamo tenere conto dell’effetto covid, per farlo utilizzeremo quella che in econometria viene definita “variabile Dummy”. Questa assume valore 1 in determinate condizioni e zero in tutti gli altri casi.

Le variabili dummy vengono utilizzate per “catturare” degli effetti che avvengono in determinate condizioni, nel nostro caso l’effetto della pandemia sulle immatricolazioni. La variabile Dummy che abbiamo creato sarà indicata con il nome di Covid19 e assumerà valore uno nell’anno 2020/2021 e zero in tutti gli altri anni.

Il modello costruito è lineare e molto semplice dal momento che la finalità della sua costruzione non è previsiva, ma essenzialmente descrittiva. La forma funzionale che abbiamo ipotizzato è la seguente:

\[ \Delta IMM_t = \alpha_0 + \alpha_1 \Delta DIP_{t-1} + \alpha_2 Covid19\]

-\(\Delta IMM_t\): variazione percentuale degli immatricolati nell’anno t.

-\(\Delta DIP_t\): variazione percentuale dei diplomati nell’anno t.

Di seguito riportiamo le stime ottenute per i coefficenti \(\alpha_0\), \(\alpha_1\), \(\alpha_2\), che rappresentano rispettivamente:

• \(\alpha_0\): intercetta del modello, ovvero il livello minimo del fenomeno.

• \(\alpha_1\): coefficente per l’effetto diplomati e misura l’influenza della variazione dei diplomati nell’anno precedente sulle immatricolazioni nell’anno corrente.

• \(\alpha_2\): coefficente stimati per l’effetto Covid19, ci dice, in media, quanto la variazione di immatricolazioni sia dovuta a fattori presenti nell’anno 2020 (ovvero la Pandemia).

modello_ridotto = lm(var_immatricolati ~ lag1_var_diplomati + covid19, data = df_model)

broom::tidy(modello_ridotto) %>% mutate(estimate = round(estimate, 3)) %>% mutate_if(is.numeric, 
    round, 2) %>% select(alpha = "term", stima = "estimate", errore = "std.error", 
    p.value) %>% mutate(alpha = c("Intercetta", "Effetto Diplomati (t-1)", "Effetto Covid19 (2020)")) %>% 
    knitr::kable(caption = "Modello con Diplomati") %>% kableExtra::kable_classic() %>% 
    kable_styling() %>% row_spec(c(1, 2, 3), font_size = 14) %>% row_spec(c(2, 3), 
    bold = T, color = "black", background = "yellow") %>% column_spec(c(3, 4), background = "white", 
    color = "black")

Un modello di questo tipo potrebbe essere utilizzato in due modi: descrivere e prevedere. In questo caso a noi interessano entrambe, soprattutto la descrizione. Dal momento che abbiamo stimato una relazione lineare tra le due variabili, possiamo facilmente interpretare i coefficenti del modello:

• Ad una variazione dell’1 % dei diplomati nell’anno precedente corrisponde, in media, una variazione del 0,88 % del numero di Immatricolati nell’anno corrente. L’impatto dei diplomati sembra esistere ed essere positivo. Per capirne l’utilità pratica, proviamo a fare un esempio: immaginiamo di sapere che quest’anno i diplomati sono aumentati del 5% rispetto allo scorso anno, in base al modello potremmo dire che l’impatto di questa variazione sugli immatricolati l’anno prossimo sarà del 4,4% in più, quindi prevediamo un aumento. Al contrario, immaginiamo che i diplomati siano diminuiti del 5% rispetto allo scorso anno, l’impatto sarà sempre 4,4%, ma in meno.

• Il coefficente stimato per l’effetto Covid-19 sta ad indicare che la variazione di immatricolati per il 2020 è stata in media, del 2 % maggiore rispetto agli altri anni. Tuttavia il coefficente stimato non risulta significativo, il che in statistica equivale a dire che la probabilità che il coefficente stimato sia in realtà uguale a zero è alta. Di conseguenza, potremmo affermare che: non sembrano esserci delle evidenze statistiche circa un effetto particolare (attribuibile alla pandemia) sulle immatricolazioni a livello nazionale nel 2020.

tibble(anno_acc = df_model$anno[-c(1, 2)], empiric = df_model$imm_totali[-c(1, 2)], 
    prediction = (1 + predict(modello_ridotto)/100) * df_model$lag1_imm_totali[-c(1, 
        2)], lower_bound = as.vector(confint(modello_ridotto, level = 0.99)[, 1] %*% 
        t(model.matrix(modello_ridotto))), upper_bound = as.vector(confint(modello_ridotto, 
        level = 0.99)[, 2] %*% t(model.matrix(modello_ridotto)))) %>% ggplot() + 
    geom_line(aes(anno_acc, empiric/1000, group = 1), size = 0.8, alpha = 0.5) + 
    geom_line(aes(anno_acc, prediction/1000, group = 1), col = "red", size = 1, alpha = 0.8) + 
    # geom_line(aes(anno_acc, lower_bound, group = 1), col = 'blue', alpha = 0.5) +
# geom_line(aes(anno_acc, upper_bound, group = 1), col = 'blue', alpha = 0.5) +
xlab("Anno Accademico") + ylab("Immatricolati (migliaia)") + ggtitle("Previsione Immatricolati: Valori Reali vs Previsti", 
    subtitle = "Modello con i Diplomati") + ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, 
    face = "bold"), axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Rielaborazione dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = c(1965, 1970, 1975, 1980, 1985, 1990, 1995, 2000, 
        2005, 2010, 2015, 2020))

Considerazioni

• Come possiamo vedere dal grafico il modello si adatta abbastanza bene alla realtà, tuttavia tende a sovrastimare in alcuni periodi con grandi variazioni, soprattutto per il 1993 il valore previsto si discosta molto dalla realtà.

• Utilizzando questo modello non sono state riscontrate evidenze statistiche su un eventuale effetto covid19.

• Se volessimo utilizzare questo modello per prevedere le immatricolazioni future rimarremmo vincolati ai dati sui diplomati. Una previsione è utile solo se viene fatta con un certo anticipo, di conseguenza questa non sarebbe la scelta più adeguata.

Previsione per il 2022

Il modello prevede la variazione percentuale di immatricolati, quindi per ricavare il valore assoluto occorre articolare il procedimemento in due parti: in primo luogo una previsione delle variazioni percentuali e in secondo luogo un’estrapolazione dei valori assoluti, il dettaglio della procedura è riportato in appendice.

df_test = df_model[56:58,]

var2019 = predict(modello_ar1, df_test[2,]) 
pred2019 = (1+var2019/100) * df_test[1,]$imm_totali

var2020 = predict(modello_ar1, df_test[3,])
pred2020 = (1+var2020/100) * df_test[2,]$imm_totali

var2021 = coef(modello_ar1)[1]+coef(modello_ar1)[2]*var2020
pred2021 = (1+var2021/100) * pred2020

low.bound = confint(modello_ar1, level = 0.99)[1,1]+confint(modello_ar1, level = 0.9)[2,1]*var2020
upp.bound = confint(modello_ar1, level = 0.9)[1,2]+confint(modello_ar1, level = 0.7)[2,2]*var2020
low.pred = ((1+low.bound/100) * pred2020) %>% as.integer()
upp.pred = ((1+upp.bound/100) * pred2020) %>% as.integer()


df_pred = 
bind_rows(
tibble(variabile  = c("variazione 2019/20", "variazione 2020/21", "variazione 2021/22"), 
       previsione = c(var2019, var2020, var2021),
       realta = c(df_test$var_immatricolati[2:3], NA),
       errore = realta - previsione)%>%
  mutate_if(is.numeric, round, 2),

tibble(variabile  = c("immatricolati 2019/20", "immatricolati 2020/21", "immatricolati 2021/22"), 
       previsione = c(pred2019, pred2020, pred2021),
       realta = c(df_test$imm_totali[2:3], NA),
       errore = realta - previsione) %>%
  mutate_if(is.numeric, round, 1) 
)  %>%
  mutate(realta =as.character(realta)) %>%
  mutate(errore = as.character(errore))

df_pred[3,3] = "?"
df_pred[6,3] = "?"
        
df_pred[3,4] = paste(round(low.bound,2), "-", round(upp.bound,2))
df_pred[6,4] = paste(low.pred,"-", upp.pred)

df_pred%>%
knitr::kable(caption = "Previsione Immatricolati con Modello AR(1)") %>%
kableExtra::kable_classic() %>% 
kable_styling() %>%
  row_spec(c(5,6), bold = T, color = "black", background = "yellow", font_size = 13)  

Risulta molto soddisfacente vedere che pur non considerando nessun effetto Covid19, e non utilizzando nelle stime alcun dato nè sul 2019 nè sul 2020, nel 2018 avremmo potuto fare una stima abbastanza precisa della situazione attuale. La previsione del 2019 risulta minore e soggetta ad un errore più elevato. Abbiamo voluto aggiungere anche una previsione ottenuta con il nostro modello per l’anno accademico 2021/2022: prevediamo tra i 324mila e i 339mila immatricolati con un intervallo di confidenza pari al 90%. Possiamo affermare che in ambito previsivo il secondo modello si rivela quello piu utile essendo in grado di fornire una stima approssimativa conoscendo solo i dati sugli immatricolati.

Per valutare la bontà di adattamento del modello autoregressivo nel suo complesso possiamo fare un confronto tra i valori previsti (in rosso) e i valori reali:

tibble(anno_acc = df_model$anno[-c(1, 2)], empiric = df_model$imm_totali[-c(1, 2)], 
    prediction = (1 + predict(modello_ar1, newdata = df_model)[-c(1, 2)]/100) * df_model$lag1_immatricolati[-c(1, 
        2)], lower_bound = (1 + as.vector(confint(modello_ar1, level = 0.99)[, 1] %*% 
        t(model.matrix(var_immatricolati ~ lag1_var_immatricolati, df_model)))/100) * 
        df_model$lag1_immatricolati[-c(1, 2)], upper_bound = (1 + as.vector(confint(modello_ar1, 
        level = 0.99)[, 2] %*% t(model.matrix(var_immatricolati ~ lag1_var_immatricolati, 
        df_model)))/100) * df_model$lag1_immatricolati[-c(1, 2)]) %>% filter(anno_acc > 
    1965) %>% ggplot() + geom_line(aes(anno_acc, empiric/1000, group = 1), size = 0.8, 
    alpha = 0.5) + geom_line(aes(anno_acc, prediction/1000, group = 1), col = "red", 
    size = 1, alpha = 0.8) + xlab("Anno Accademico") + ylab("Immatricolati (migliaia)") + 
    ggtitle("Previsione Immatricolati: Valori Reali vs Previsti", subtitle = "Modello Autoregressivo") + 
    ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, face = "bold"), 
    axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Rielaborazione dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = c(1965, 1970, 1975, 1980, 1985, 1990, 1995, 2000, 
        2005, 2010, 2015, 2020))

Nel complesso l’adattamento del modello autoregressivo è migliorata rispetto al modello precedente che utilizzava i diplomati, seppur la stima è ancora soggetta ad un errore elevato.

• Nel 1993 quest modello è molto più preciso di quello con i diplomati, avvicinandosi di molto alla realà.

• Dal 2010 e nel 2018, il modello prevede quasi perfettamente il numero di immmatricolati.

• Tuttavia le inversioni di trend si presentano con un leggero ritardo.

Un modello con le Riforme Universitarie

Dopo aver delineato il quadro legislativo degli ultimi 20 anni, ci siamo chiesti se fosse possibile misurare le variazioni di immatricolati in concomitanza delle riforme. Per cercare di rispondere alla nostra domanda seguiremo l’approccio utilizzato nella formulazione del modello con i diplomati, includendo delle variabili dummy, che avranno valore uno in determinati anni e zero in altri:

• Riforma Berlinguer: effetto stimato dal 2000 al 2003.
• Riforma Moratti: effetto stimato dal 2005 al 2006.
• Riforma Brunetta: effetto stimato dal 2008 al 2009.
• Riforma Gelmini: effetto stimato dal 2010 al 2013.
• CrisiIT: è una variabile di controllo per catturare gli effetti della crisi italiana nel 2012.

Dal momento che per ottenere deglie errori minori è stato stimato un modello più complesso che è in grado cogliere anche degli effetti non lineari, abbiamo deciso di riportare solo le stime delle variabili di interesse. I dettagli completi sulla stima con alcuni commenti sono riportati nell’appendice.

modello_completo = lm(var_immatricolati ~ lag1_var_diplomati + I(lag1_var_diplomati^2) + 
    lag1_var_immatricolati + RiformaBerlinguer + RiformaMoratti + RiformaBrunetta + 
    CrisiIT + RiformaGelmini, data = df_model)

broom::tidy(modello_completo)[-c(1:4), ] %>% mutate(estimate = round(estimate, 3)) %>% 
    mutate_if(is.numeric, round, 2) %>% select(coefficente = "term", stima = "estimate", 
    errore = "std.error", p.value) %>% knitr::kable(caption = "Modello Completo: gli effetti delle Riforme") %>% 
    kableExtra::kable_classic() %>% kable_styling() %>% row_spec(c(1, 3), bold = TRUE, 
    color = "black", background = "yellow", font_size = 13) %>% row_spec(c(2, 4, 
    5), italic = TRUE, color = "black", background = "white", font_size = 13)

Interpretazione dei risultati:

Considerando che nella realtà l’impatto di queste riforme è influenzato da moltissimi fattori esterni sconosciuti, proveremo a dare delle interpretazioni sui coefficenti stimati, che, tuttavia, devono essere considerate con molta cautela. Inoltre è bene precisare che gli effetti stimati e le interpretazioni devono essere considereate in relazione agli anni di riferimento, volendo fare un esempio con la riforma Berlinguer potremmo dire che: tra gli anni 2000 e 2003 (periodo in cui era stata appena emanata la riforma Berlinguer) la variazione di immatricolati è aumentata, a parità di altre condizioni, del 4,93%.

• La Riforma Berlinguer sembra aver avuto un impatto positivo sulle immatricolazioni, portando, in media, per gli anni che vanno dal 2000 al 2003, ad un aumento del 4,93% di immatricolazioni. L’ipotesi che la divisione in lauree triennali e magistrali possa aver avuto nel complesso un impatto significativo viene rafforzata. Inoltre la probabilità che il coefficente stimato sia errato è molto bassa (5%) e ci concede un certo grado di sicurezza riguardo le affermazioni fatte precedentemente.

• La Riforma Moratti: l’impatto della Riforma Moratti potrebbe sembrare a primo impatto negativo, dal momento che si stima una variazione negativa del -2,46% tra il 2005 e il 2006. Tuttavia, come anticipato nella descrizione della legge, questa riforma non ha influenzato direttamente il mondo degli studenti universitari, bensi quello dei docenti. Se consideriamo inoltre che la probabilità di errore è del 48%, possiamo affermare che: non ci sono evidenze statistiche circa l’influenza della Rifoma Moratti sulle immatricolazioni.

• La Riforma Brunetta: le novità introdotte con questa legge sono state molteplici e hanno sicuramente influenzato il comportamento degli atenei italiani, tuttavia il periodo riguardante la riforma, coincide con quello della Crisi del Debito Pubblico, iniziata nel 2007. Di conseguenza, anche se la probabilità che il coefficente sia errato è del 22%, il che non esclude completamente un’effetto della legge, è molto rischioso attribuirne il calo stimato del -4,25% unicamente alla riforma Brunetta.

• CrisiIT e Riforma Gelmini: abbiamo inserito la variabile CrisiIT, per cercare di ridurre l’errore sulla stima per la riforma Gelmini, dato che è avvenuta in un periodo di profonda crisi economica. La riforma sembrerebbe aver avuto un impatto negativo -2.37%, tuttavia la probabilità di errore è molto alta (43%), di conseguenza possiamo affermare che nemmeno per la riforma Gelmini ci sono evidenze statistiche circa un effetto sulle immatricolazioni.

Effetto NON-lineare dei Diplomati

Il modello completo, costruito in questo paragrafo, ha messo in risalto un fattore molto interessante che vale la pena commentare: non solo ha confermato l’impatto che ha una variazione di diplomati sulle immatricolazioni, ma ha anche messo in risalto che la relazione reale tra le due variabili non è lineare.

x = seq(-20, 20, 0.01)

ggplot() + geom_line(data = tibble(x = x, y = coef(modello_completo)[2] * x + coef(modello_completo)[3] * 
    x^2), aes(x, y), size = 1, alpha = 1, color = "red") + geom_line(data = tibble(x = x, 
    y = coef(modello_ridotto)[2] * x), aes(x, y), size = 1, alpha = 1, color = "black") + 
    xlab("Variazione Diplomati anno passato (%)") + ylab("Variazione Immatricolati anno corrente (%)") + 
    ggtitle("Effetto non-lineare dei Diplomati", subtitle = "Effetto della variazione di diplomati: non-lineare nel Modello Completo (rosso) \n e lineare nel Modello Diplomati (nero)") + 
    ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, face = "bold"), 
    axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Rielaborazione dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = seq(-20, 20, 5))

La relazione lineare tra le due variabili sembrerebbe valida solo intorno allo zero, quindi per piccole variazioni, mentre nel caso di variazioni più consistenti la relazione diventa non-lineare. Per provare ad interpretarne le implicazioni osserviamo la curvatura della linea rossa, possiamo affermare che:

• Un aumento elevato di diplomati produce un aumento meno che proporzionale nella variazione di immatricolazioni. Volendo ottenere una prova empirica di questo fatto basta osservare le stime ottenute per il periodo che va dal 1990 al 1993 con il modello lineare del paragrafo 2. In quel periodo il numero di diplomati crebbe molto, tuttavia le immatricolazioni crebbero, ma non proporzionalmente difatti il modello sbaglia enormemente le previsioni. Il modello completo, utilizzando una relazione non lineare, riesce invece a catturare più precisamente il reale effetto dei diplomati.

• Una diminuzione drastica dei diplomati produce una diminuzione più che proporzionale nella variazione di immatricolazioni.

I risultati sono molto interessanti, soprattutto se letti in chiave demografica: le previsioni per i prossimi anni indicano chiaramente una diminuzione del numero di giovani, il che si tradurrà inevitabilmente in una diminuzione del numero di diplomati. Di conseguenza, almeno stando al modello, si avrebbe un calo drastico delle immatricolazioni: una diminuzione del 10% nel numero di diplomati provocherebbe una dimuzione del 15% circa di immatricolazioni.

Conclusioni

1. La divisione delle lauree in triennali e magistrali sembrerebbe aver avuto un impatto consistente, ma transitorio e circorscitto al periodo 2000-2003.

2. Non sono state riscontrate evidenze circa gli effetti delle altre riforme. Tuttavia, considerando il contesto storico-politico in cui vennero emanate, possiamo presupporre che in quegli stessi anni siano presenti ulteriori fattori che possano aver influenzato la scelta dei nuovi studenti.

3. L’effetto dei diplomati sugli immatricolati è positivo e non lineare, la relazione stimata tra le due variabili può essere considerata abbastanza stabile e affidabile. Le implicazioni nella realtà come abbiamo accennato sono molteplici e destano particolari preoccupazioni riguardo al futuro.

4. Non abbiamo riscontrato evidenze su effetto del Covid19 sulle immatricolazioni a livello nazionale. Se volessimo considerare la stima ottenuta per il coefficente come valida potremmo affermare che l’impatto globale sembra essere stato positivo. Una possibile spiegazione che ci siamo dati su quest’aumento è la seguente: persone che, in una situazione normale, avrebbero preferito iniziare a lavorare una volta terminata la scuola superiore, essendosi trovate in un periodo di Lockdown, con offerte di lavoro scarse potrebbero aver deciso di iscriversi all’università. Non è da escludere che, una volta ripristinata una situazione normale queste persone abbandonino poi l’università, per verificare quindi se l’aumento sia stato effettivo, occorre aspettare l’anno prossimo e osservare l’evoluzione.

R2 = glance(modello_completo)$adj.r.squared %>% round(2)
sd_mod = glance(modello_completo)$sigma %>% round(2)

tibble(anno_acc = df_model$anno[-c(1, 2)], empiric = df_model$imm_totali[-c(1, 2)], 
    prediction = (1 + predict(modello_completo)/100) * df_model$lag1_imm_totali[-c(1, 
        2)]) %>% ggplot() + geom_line(aes(anno_acc, empiric/1000, group = 1), size = 0.8, 
    alpha = 0.5) + geom_line(aes(anno_acc, prediction/1000, group = 1), col = "red", 
    size = 1, alpha = 0.8) + xlab("Anno Accademico") + ylab("Immatricolati (migliaia)") + 
    ggtitle("Valori Reali vs Valori Previsti", subtitle = "Modello Completo con Riforme") + 
    ggthemes::theme_solarized() + theme(axis.text.x = element_text(angle = 0, face = "bold"), 
    axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Rielaborazione dati Istat, Miur Open Data") + 
    scale_x_continuous(breaks = c(1965, 1970, 1975, 1980, 1985, 1990, 1995, 2000, 
        2005, 2010, 2015, 2020))

Indicatori Regionali

Per completare il quadro a livello regionale occore caratterizzare le regioni non solo in base agli studenti che acquisiscono, ma anche a quelli che perdono. Partendo dalle variazioni globali di immatricolati per regione, abbiamo deciso di calcolare una serie di indicatori che ci possano aiutare ad effettuare dei confronti a livello regionale e fornire delle informazioni in grado di spiegare le variazioni di quest’anno. Con i confronti vorremmo cercare di comprenderne due aspetti fondamentali:

1. Quali sono le regioni che hanno un bilancio tra studenti in entrata e in uscita negativo nel 2021? Quali invece, positivo?

2. Quali sono le regioni con un’elevata emigrazione?

Gli indicatori calcolati a livello regionale sonoi seguenti:

• Variazione Studenti Fuorisede: rappresenta l’aumento/diminuzione percentuale di studenti che si sono immatricolati nella regione avendo la residenza in una regione differente.

• Variazione Studenti Emigrati: rappresenta l’aumento/diminuzione percentuale di studenti che essendo residenti nella regione, si sono immatricolati in altre regioni.

• Fuorisede/Totale: è il rapporto tra il numero di fuorisede immatricolati e il totale di immatricolati nella regione.

• Emigrati/Residenti: è il rapporto tra il numero studenti residenti che si immatricola altrove e il numero di studenti residenti che si immatricola nella regione. Il caso in cui questo indicatore è maggiore 100%, indica che l’emigrazione dalla regione è nettamente maggiore del numero di giovani che restano. Se lo utilizziamo insieme al rapporto dei fuorisede sul totale possiamo subito capire quale sia la dinamica regionale, ovvero se quanti studenti lasciano la regione, quanti restano e quanti arrivano da fuori rispetto al totale.

• Immigrati/Residenti: è il rapporto tra il numero studenti fuorisede che si immatricola nella regione e il numero di studenti residenti che si immatricola nella regione. Nel caso in cui nella regione si immatricolino molti fuorisede, e allo stesso tempo si immatricolano meno residenti, l’indice assumerà un valore maggiore di 100.

• Bilancio 2021: è dato dalla differenza (in termini assoluti) tra il numero di studenti fuorisede e il numero di studenti emigrati.

Per facilitare la visione complessiva degli indicatori, abbiamo deciso di costruire una mappa interattiva, in cui ciascun indice può essere visualizzato cliccando sulla regione di interesse. Il colore delle regioni rispecchia le variazioni globali di immatricolazioni nel 2020/21.

#green.scale = c("#E4F2DA", "#D6F2C2", "#C8F2AA", "#BAF291", "#ACF279", "#9DF261", "#8FF249", "#81F230", "#73F218", "#65F200")
#red.scale   = c("#F5E0DC", "#F5CCC4", "#F5B8AB", "#F5A393", "#F58F7A", "#F57A62", "#F56649", "#F55231", "#F53D18", "#F52900")

bins <- c(-10, -5, -2.5, -1,  -0.5, 0, 5, 10, 20, 40, 60)
mypalette = c("#E60400", "#E6312E", "#E65E5C","#E6807E", "#E6918A","#C8F2AA", "#ACF279","#88E645", "#6DE617","#6AFF00")

pal = colorBin(  palette=mypalette, domain = df_emigrazioni$var_totale, bins = bins)

labels = paste0("<b>", df_emigrazioni$regione, 
                 "</b>", "<br/>", "Variazione Immatricolati tra il 2019/20 e il 2020/21 ", "<b>",  round(df_emigrazioni$var_totale, 2), " %<b/>",
                 "</b>", "<br/>", "Variazione Fuorisede: ", "<b>",  round(df_emigrazioni$var_fuorisede, 2), " %<b/>",
                 "</b>", "<br/>", "Fuorisede/Totale: ", "<b>",  round(df_emigrazioni$immigrati_2021/df_emigrazioni$totale_2021*100, 2), " %<b/>",
                 "</b>", "<br/>", "Variazione Studenti Emigrati: ", "<b>",  round(df_emigrazioni$var_emigrati, 2), " %<b/>", 
                 "</b>", "<br/>"," Emigrati/Residenti: ", "<b>",  round(df_emigrazioni$emigrati_2021/df_emigrazioni$residenti_2021*100, 2), " %<b/>",
                 "</b>", "<br/>","Immigrati/Residenti: ", "<b>",  round(df_emigrazioni$immigrati_2021/df_emigrazioni$residenti_2021*100, 2), " %<b/>",
                 "</b>", "<br/>","Bilancio 2021: ", "<b>",  round(df_emigrazioni$bilancio2021, 2), " studenti. <b/>"
                 ) %>% lapply( htmltools::HTML)

map1 = base_map %>%
  addPolygons(
  popup = labels,
  fillColor = ~pal(df_emigrazioni$var_totale),
  weight = 2,
  opacity = 1,
  color = "white",
  dashArray = "3",
  fillOpacity = 0.7,
  labelOptions = labelOptions(
      style = list("font-weight" = "normal", padding = "3px 8px"),
      textsize = "10px",
      direction = "auto"),
  
  highlight = highlightOptions(
    weight = 5,
    color = ~"#E8A025",
    dashArray = "",
    fillOpacity = 0.7,
    bringToFront = TRUE))

map1 %>% 
    addLegend(pal = pal, values = ~df_emigrazioni$var_totale, opacity = 0.7, title = "Var. Immatricolati 2019-21",
                                       position = "topright") %>% 
  addMeasure(
    position = "bottomleft",
    primaryLengthUnit = "kilometers",
    primaryAreaUnit = "sqmeters",
    activeColor = "#3D535D",
    completedColor = "#7D4479")

I Bilanci del 2021

1. Quali sono le regioni che hanno un bilancio tra studenti in entrata e in uscita negativo nel 2021? Il grafico evidenzia come la maggior parte delle regioni con bilancio negativo si trovi nel CENTRO-SUD, con alcune eccezioni come la Campania e il Lazio. Il Veneto è l’unica regione del nord d’Italia a presentare un bilancio fortemente negativo, ad indicare una forte emigrazione di giovani dalla regione. L’Emilia Romagna è la regione che presenta il bilancio maggiore, come possiamo notare dagli indicatori calcolati nel 2021 la percentuale di residenti emigrati è stata circa del 20%, inoltre rispetto al 2019, sono emigrati più studenti, circa il 25% in più rispetto all’anno precedente. L’Emilia è tradizionalmente una regione fortemente dipendente dalle immatricolazioni di fuorisede, nel 2021 il 40% degli studenti iscritti in Emilia Romagna proviene da altre regioni.

df_emigrazioni %>% mutate(Bilancio = bilancio2021, var = ifelse(Bilancio > 0, "Pos", 
    "Neg")) %>% mutate(regione = case_when(regione == "TRENTINO ALTO ADIGE" ~ "TRENTINO", 
    regione == "FRIULI VENEZIA GIULIA" ~ "FRIULI", TRUE ~ regione)) %>% ggplot(aes(Bilancio, 
    reorder(regione, Bilancio), fill = var)) + geom_bar(stat = "identity", width = 0.8, 
    col = "black") + xlab("Bilancio (valore assoluto)") + ylab("") + ggtitle("Il Bilancio delle immatricolazioni nel 2021", 
    subtitle = "Il bilancio è stato calcolato come differenza tra \n gli immatricolati residenti e gli immatricolati fuorisede.") + 
    ggthemes::theme_solarized() + scale_x_continuous(breaks = seq(-6000, 10000, 2000)) + 
    scale_fill_manual(values = c("#F58772", "#B7F28D")) + theme(axis.text.x = element_text(angle = 0, 
    face = "bold"), axis.text.y = element_text(face = "bold"), axis.title = element_text(face = "bold"), 
    plot.title = element_text(face = "bold"), plot.subtitle = element_text(face = "italic"), 
    plot.caption = element_text(face = "italic"), panel.grid.major.x = element_blank(), 
    panel.grid.minor.x = element_blank(), panel.grid.major.y = element_line(colour = "grey60", 
        linetype = "dotted"), legend.position = "none") + labs(caption = "Fonte: Miur Open Data")

Test sul Modello Ridotto

mod = modello_ridotto

bind_rows(lmtest::bgtest(mod) %>% broom::tidy(), lmtest::dwtest(mod) %>% broom::tidy(), 
    lmtest::bptest(mod) %>% broom::tidy(), lmtest::gqtest(mod) %>% broom::tidy(), 
    lmtest::resettest(mod) %>% broom::tidy())[, c(1, 2, 4)] %>% mutate_if(is.numeric, 
    round, 2) %>% knitr::kable(caption = "Testo Sul Modello Diplomati") %>% kableExtra::kable_paper()

I primi quattro test servono ad evidenziare possibili errori di specificazione del modello, in questo caso evidenziano sia la presenza di autocorrelazione, sia la presenza di eteroschedasticità nei residui. Nella realtà questa tipologia di errore si riflette in una stima più elevata degli errori per i coefficenti stimati, tuttavia il valore dei coefficenti rimane stabile, il che significa che l’interpretazione rimarrebbe la stessa, il che è più che sufficente per le nostre finalità. L’ultimo test è quello piu temibile in quanto la sua significatività implicherebbe lo star sbagliano tutto, il che fortunatamente non è il nostro caso.

Previsione per il Modello Autoregressivo

La previsione della variazione è data dalla somma dell’intercetta (1.06) e la moltiplicazione del coefficente stimato (0.53) per la variazione di immatricolati nell’anno precedente. Il modello con i coefficenti stimati diventa:

\[\hat{\Delta IMM_t} = 1.06 + 0.53 * \Delta IMM_{t-1}\]

Gli immatricolati previsti sono dati dal prodotto degli immatricolati al tempo precedente per la variazione prevista.

\[\hat{IMM_t} = (1+\frac{\hat{\Delta IMM_t}}{100}) * IMM_{t-1}\]

Test Modello Autoregressivo

# Test per l'autocorrelazione dei Residui (rifiuto H0 -> Autocorrelazione) Test
# per la presenza di eteroschedasticita (non rifiutiamo H0 -> No evidenza di
# Eteroschedasticità)

mod = modello_ar1

# Test per l'autocorrelazione dei Residui (rifiuto H0 -> Autocorrelazione) Test
# per la presenza di eteroschedasticita (non rifiutiamo H0 -> No evidenza di
# Eteroschedasticità)
bind_rows(lmtest::bgtest(mod) %>% broom::tidy(), lmtest::dwtest(mod) %>% broom::tidy(), 
    lmtest::bptest(mod) %>% broom::tidy(), lmtest::gqtest(mod) %>% broom::tidy(), 
    lmtest::resettest(mod) %>% broom::tidy())[, c(1, 2, 4)] %>% mutate_if(is.numeric, 
    round, 2) %>% knitr::kable(caption = "Testo Sul Modello AR(1)") %>% kableExtra::kable_paper()

Tutti i test sul modello autoregressivo non evidenziano errori particolari, il che è un dato positivo in quanto l’utilità di questo modello è prevalentemente previsiva: un eventuale stima errata negli errori dei coefficenti si ripercuote su una previsione più ampia e con un valore più incerto.

Test sul Modello Completo

Cercare di ottenere dei coefficenti significativi per gli effetti delle Riforme (con i p.value minimi) è stata messa in atto la stategia di utilizzare una regressione polinomiale di secondo grado nelle variabili che ipotizziamo essere la causa delle immatricolazioni. Stimeremo quindi un unico modello in cui verranno inserite sia le variazioni di diplomati nell’anno passato, sia le variazioni degli immatrcolati. Utilizzando una regressione con i quadrati delle variabili siamo in grado di cogliere degli effetti non lineari, tuttavia essendo la regressione non lineare, i coefficenti non sono più facilmente interpretabili. La forma funzionale stimata è la seguente:

\[\Delta IMM_t = \alpha_0 + \alpha_1 \Delta DIP_{t-1} + \alpha_2 \Delta DIP_{t-1}^2 + \alpha_3 \Delta IMM_{t-1} + \\ + \alpha_4 Berlinguer + \alpha_5 Moratti + \alpha_6 Brunetta + \\ + \alpha_7 CrisiEu + \alpha_8 Gelmini + \epsilon_t \]

res = broom::tidy(modello_completo)
res[1, 1] = "Alpha_0 (Intercetta)"
res[2, 1] = "DIP (t-1)"
res[3, 1] = "DIP^2 (t-1)"
res[4, 1] = "IMM (t-1)"
res %>% mutate(estimate = round(estimate, 3)) %>% mutate_if(is.numeric, round, 2) %>% 
    select(alpha = "term", stima = "estimate", errore = "std.error", p.value) %>% 
    knitr::kable(caption = "Stime: Modello Completo") %>% kableExtra::kable_classic() %>% 
    kable_styling() %>% row_spec(c(2, 3, 5), bold = T, color = "black", background = "yellow", 
    font_size = 13)

# Test per l'autocorrelazione dei Residui (rifiuto H0 -> Autocorrelazione) Test
# per la presenza di eteroschedasticita (non rifiutiamo H0 -> No evidenza di
# Eteroschedasticità)

mod = modello_completo

# Test per l'autocorrelazione dei Residui (rifiuto H0 -> Autocorrelazione) Test
# per la presenza di eteroschedasticita (non rifiutiamo H0 -> No evidenza di
# Eteroschedasticità)
bind_rows(lmtest::bgtest(mod) %>% broom::tidy(), lmtest::dwtest(mod) %>% broom::tidy(), 
    lmtest::bptest(mod) %>% broom::tidy(), lmtest::gqtest(mod) %>% broom::tidy(), 
    lmtest::resettest(mod) %>% broom::tidy())[, c(1, 2, 4)] %>% mutate_if(is.numeric, 
    round, 2) %>% knitr::kable(caption = "Testo Sul Modello Completo") %>% kableExtra::kable_paper()

Dati Europei

• European Union Open Data
• European University Institute
• Portale Europeo dei Dati

Altri Open Data

• Datahub
• Data World
• Kaggle
• World Bank