Patrizia Favaron
Introduzione
“Media”. Cosa c’è, di più semplice?
Ed infatti è cosa semplice: basta prendere i numeri che ci interessano, sommarli tra loro, e dividere per quanti erano.
Basta una calcolatrice.
Sino a che i numeri sono, diciamo, una decina.
Da dieci a cento, un foglio Excel comincia a non fare schifo.
Oltre i cento… 🙄
In questo articolo vedremo come si può calcolare una media-monstre di 263089 misure di velocità del vento (quel numero lì di misure, e non un altro, perché vogliamo calcolare la media complessiva delle velocità mediate su 10 minuti).
Quanto ai dati, sono le velocità del vento raccolte dalla stazione meteorologica ARPA Lombardia di Bormio Eliporto nel periodo dal 01. 01. 2020 al 01. 01. 2025, per un totale di 5 anni complessivi: non abbastanza per fare della climatologia, ma più che sufficienti per indicarci qualche domanda e, se possibile, ottenere delle risposte.
I dati – intanto
Acci. Duecentosessantatremilaottantanove numeri…
Mica pochi.
E tutti corredati da una “marca temporale”, che indica quando sono stati raccolti.
Tempo, più dati, uguale “serie temporale”: una raccolta di misure, più o meno lunga, che in linea di principio ci può permettere di capire qualcosa (se c’è qualcosa da capire) sull’andamento del vento a Bormio.
Serie temporali come questa sono oggetti concettualmente diversi dagli insiemi di dati che capita spesso di dover comprendere a colpi di statistica. Per vederlo, consideriamo un insieme di dati più “normale”, come ad esempio le lunghezze dei lingotti di acciaio che arrivano al laminatoio. In teoria queste lunghezze dovrebbero essere tutte esattamente uguali, ma nella realtà tra le misure dell’una e dell’altra qualche differenza c’è sempre – o potrebbe esserci.
In questo caso, se desideriamo quantificare la produzione giornaliera di lamiera dovremo conoscere qualcosa delle lunghezze (ed anche altezze, spessori, pesi, presenza o no di bolle e cricche, …) dei lingotti che entrano nel processo. E per farlo, uno dei modi migliori è di calcolare la media aritmetica delle lunghezze (altezze, spessori, pesi, …) e, per ogni lingotto, lo scostamento in più od in meno rispetto alle medie. Così, tanto per farsi un’idea: lingotto piccolo, o grande? Pieno di fratture che magari denotano una qualità del materiale sotto standard? E via così.
In un certo senso, però, ogni lingotto è un caso a sé. E tutti sono stati prodotti da uno stesso processo, con effetti sperabilmente identici.
Ho detto “lingotti”, ma avrei potuto parlare di risultati di misure di un esperimento di laboratorio, o di tempo passato dalle persone in un supermercato, oppure di prezzo di esemplari diversi di uno stesso bene. O di un’infinità di altre cose. Nel caso, la media ci da un’idea piuttosto precisa di cosa potremmo aspettarci, arrivassero nuovi esemplari dell’oggetto, o della misura. E la media del valore assoluto, o del quadrato, degli scostamenti ci direbbe quanto sono “dispersi” i dati.
Ma qui siamo in una situazione molto diversa: ogni dato, etichettato da un tempo, è sì diverso dagli altri, ma in un certo senso di più. In una serie temporale, i valori misurati non sono necessariamente il risultato di misure tra loro indipendenti: nel caso delle velocità del vento, per esempio, se in corrispondenza di una certa marca temporale abbiamo una lettura di, facciamo finta, 10 m/s allora sarebbe piuttosto sorprendente scoprire che il dato immediatamente successivo è di 1 m/s o 50 m/s: è molto più facile che troviamo un valore non troppo distante da 10 m/s.
Più facile, certo, e non “sicuro”: comunque, siamo di fronte ad un fenomeno caratterizzato da un certo grado di variabilità interna, ed incertezza, ma anche da una qualche struttura di dipendenza tra un valore e quelli che lo seguono nel tempo.
Di più: le serie temporali sono, spesso, la manifestazione visibile di processi fisici ben precisi, che in qualche modo ci interessano.
Poi, la natura di questo interesse può essere la più varia: scientifica, oppure economica (che so, quanta energia sarebbe estraibile dal vento se piazzassimo una turbina eolica a fianco della stazione), oppure, perché no, di semplice curiosità personale. Anche questo conta, e lo vedremo: ma viene dopo le misure, ed il processo che si nasconde dietro ad esse.
Comunque, il fatto che una serie temporale non sia la stessa cosa di un insieme di misure di laboratorio non vuol dire che sui suoi valori non possiamo determinare la media, ed altri indicatori statistici.
Popolazione, e campioni
Supponiamo, così, di voler determinare il valor medio di tutte le misure che abbiamo a disposizione. Questo insieme, l’insieme su cui desideriamo compiere qualche affermazione “forte” basata sui dati, sarà la nostra popolazione.
I dati li ho scaricati dal sito di ARPA Lombardia e, almeno nominalmente, ci sono tutti.
Potrei caricarli dentro Excel, e fare tutto lì. Ma personalmente preferisco usare il programma statistico “R”, che da alcuni vantaggi:
- Intanto, a differenza di Excel, è gratuito.
- Poi (soprattutto) è fatto apposta per la statistica, ed ha un’infinità di funzioni che permettono di fare davvero “tutto”.
Se non lo conosci, ti invito a prenderne visione, e magari scaricarlo dal sito ufficiale e provare a usarlo, con l’aiuto di uno dei tanti bellissimi tutorial che si trovano in rete. Se già non conosci “R” questa esplorazione potrebbe comportare un certo sforzo da parte tua. Ma, assicuro, ne vale la pena.
Prima di calcolare medie e quant’altro, però, dobbiamo leggere i dati. In R, disponendo dei dati (che puoi scaricare qui – dopo averli scaricati, per poterli adoperare questi vanno “decompressi”), la lettura ed altre operazioni si fanno usando dei “comandi”, cioè delle righe di testo con le istruzioni per il calcolatore, che una volta scritte si possono fare eseguire premendo il tasto “Invio”. Nel caso, per leggere i dati (che ho sistemato in un file di nome “Vel.csv”) il comando da usare è
d <- read.csv(“Vel.csv”)
Il senso di questo comando è di leggere i dati in formato CSV ( read(“Vel.csv”) ), e riversare il contenuto nella “variabile” ‘d’. A comando eseguito, ‘d’ conterrà tutti i dati presenti nel file. Una volta che i dati sono finiti in ‘d’, possiamo farci tutto quello che vogliamo. Per esempio, con l’altro comando
names(d)
farci dire come si chiamano le colonne della tabella di dati (che sono “Id.Sensore”, “Data.Ora” e “Medio”). Un altro comando ancora,
str(d)
permette di vedere oltre ai nomi delle colonne anche qualche valore – cosa che qui non faremo: lascio a voi la prova.
Adesso, con
mean(d$Medio)
calcoliamo il valore della velocità media sui nostri cinque anni di dati. La risposta, istantanea, è
-110.9604
Come? meno centoundici metri al secondo? Ma le velocità non dovrebbero essere positive, o al peggio nulle? 🧐
E sì, in effetti sì: dovrebbero. Ma evidentemente non solo.
Quando ARPA Lombardia spedisce i dati, insieme ad essi invia sempre anche un a “legenda”, un breve file di testo con tutte le istruzioni per interpretare i dati.
Nel nostro caso il contenuto della legenda è:
======<Inizio>======
Legenda
-999 Valore mancante o invalido
Legenda dei valori speciali
Valore, Descrizione
777,calma
7777,calma
888,variabile
8888,variabile
Nota: i gradi di Direzione del Vento sono riferiti al Nord
Riepilogo estrazione
Id_Stazione,Nome_Stazione,Id_Sensore,Nome_Sensore,Unita_Misura,Periodo_Dal,Periodo_Al,Operatore
1650,Bormio eliporto,19395,Velocità Vento,m/s,2020/01/01 00:00,2025/01/01 00:00,Medio
I dati sono estratti in ORA SOLARE (UTC + 1)
=====<Fine>=====
Vediamo subito una cosa: esiste una codifica di dati invalidi, di vari tipi. Mancanti, per esempio (-999). Oppure corrispondenti ad una “calma” (777 o 7777), o ad un “vento variabile” (888 o 8888). Calme e venti variabili sono artefatti strumentali: non dipendono dal vento in quanto tale, ma dall’anemometro che stiamo usando. A volte la velocità del vento non è sufficiente a vincere i pur piccoli attriti interni, ed altre è in grado di eccitare il sensore di direzione con una risonanza meccanica. Per noi, dati invalidi propriamente detti, o calme, o venti variabili, sono tutti da escludere dalla media.
Si, ma come?
Secondo me è in cose come queste che si tocca con mano la potenza e la flessibilità di R: per prima cosa “ripuliamo” i dati, dichiarando non disponibili tutti i valori negativi, o che superano una soglia ragionevole (60 m/s). Così:
d$Medio[d$Medio < 0 | d$Medio > 60] <- NA
done “NA” è uno speciale valore numerico facile da distinguere dai valori, invece, disponibili. E questa volta calcoliamo la media escludendo i valori non disponibili:
mean(d$Medio, na.rm=TRUE)
(“na.rm=TRUE” chiede ad R proprio di escludere i dati invalidi).
Il risultato, questa volta, sembra molto più sensato:
1.858702
Quasi due metri al secondo: non sembra moltissimo, ma, ad esperienza, “ci sta”. Se la stazione si fosse trovata più vicina alla pianura mi sarei aspettata un valor medio più basso.
Un momento, però. Per calcolare la media ho dovuto escludere alcuni dati. In quei casi non è che il vento non era misurabile: c’era, eccome, ma la stazione per una ragione o per l’altra non è stata capace di misurarli.
Ora, io avevo detto che la cosa per me interessante era la media su tutti i 5 anni. Il numero che abbiamo ottenuto, però, è la media di un sottoinsieme. E cioè di un campione.
Grande quanto? Vale la pena saperlo, perché non è esattamente la stessa cosa dire che “la media vale tot” se i dati su cui l’ho calcolata sono, per esempio, solo due o tre, piuttosto che “tutti”.
Per capire, chiediamo ad R il numero di dati validi (nel caso, “non invalidi”: il punto esclamativo significa “not”):
sum(!is.na(d$Medio))
La risposta:
233433
Su quanti dati complessivi?
length(d$Medio)
Risposta:
263089
Insomma: i dati validi ammontano all’88.73% del totale. Gli invalidi, quindi, ammontano all’11.27%.
L’insieme dei dati validi, così, abbiamo veduto, è un campione: un sottoinsieme su cui calcoliamo per davvero la media. In questo caso molto specifico il campione non lo abbiamo scelto noi, ma, in pratica, il caso. Ed è, comunque, grandissimo.
In molti casi pratici, il campione invece siamo noi a sceglierlo. I modi che possiamo adoperare sono, lo vediamo tra poco, molti. Ma il motivo è sempre lo stesso: raccogliere dati è un atto di forza. I dati non crescono sugli alberi, e vanno “acquisiti” per mezzo di un processo di acquisizione che, in tutti i casi, comporta un costo.
Costo, che può essere quello del lavoro delle persone che si occupano della sua selezione ed elaborazione. Ma che include, anche, i costi vivi della manutenzione della stazione (impressionanti), della corrente elettrica, delle pratiche burocratiche di asservimento dei terreni, e di un’infinità di altri fattori.
Una volta acquisiti, i dati costituiscono un patrimonio preziosissimo. Ma, appunto, costoso.
A Bormio Eliporto abbiamo la fortuna della presenza di una stazione di misura fissa, mantenuta in vita da ARPA Lombardia, automatica, e programmata per raccogliere tutto quello che può nelle condizioni di tempo più svariate.
Ma nella pratica delle misure, questa è più spesso che no l’eccezione: la regola è che la raccolta dei dati avviene per “campagne di misura”, di durata più o meno lunga (di solito non molto).
Ed anche questo è un caso fortunato. Se invece che duecentomila e rotti misure di velocità del vento avessimo dovuto determinare la media del diametro di altrettante rondelle, allora con ogni probabilità le misure le avremmo dovute compiere a mano. Nel caso, la dimensione del campione sarebbe stata relativamente piccola, magari cento esemplari.
Campioni casuali, e non
Nel bellissimo libro “Mentire con le statistiche“, l’Autore, Darrell Huff, ci mostra perché, nel caso di statistiche generiche, è il caso di diffidare di campioni scelti in modo non-casuale. E in effetti, uno dei cardini del design degli esperimenti è proprio la scelta accurata, casuale appunto, del campione (vedi Montgomery, D. C. 2001, Progettazione e analisi degli esperimenti, McGraw-Hill).
I campioni casuali sono, infatti, particolarmente “efficienti” – se lo scopo che ci prefiggiamo è di stimare il valor medio di una popolazione.
Torniamo, allora, a noi. Ed eleggiamo (lo so, è una scorrettezza) a “popolazione” l’insieme di tutti i dati disponibili – per prenderlo, basta che usiamo il comando
population <- d$Medio[!is.na(d$Medio)]
(con questa operazione andiamo a mettere nella “variabile” ‘population’ tutti e soli i dati validi, cioè diversi da NA). La media che avevamo calcolato prima è ancora valida:
mean(population)
1.858702
Ora, facciamo finta di voler stimare questa media (che potremmo non conoscere, se non avessimo a disposizione una stazione automatica), usando campioni casuali di 100 elementi. In R possiamo farlo usando il comando
mean(sample(population,100))
La risposta:
1.601
Rifacciamolo:
mean(sample(population,100))
1.772
mean(sample(population,100))
1.884
Immaginiamo di rifarlo molte volte, per esempio 100 volte. Ogni campione è lungo 100, e così 100 campioni sono formati dal 10000 numeri: una frazione molto piccola delle dimensioni della popolazione.
E questi, i valori:
1.632 1.804 1.915 1.885 1.502 1.856 2.649 1.856 1.828 1.758 2.290 2.086 2.171 1.835 1.914 2.355 1.603 1.999 1.941 1.972 2.016 1.598
1.807 1.808 1.974 2.071 1.437 1.727 1.804 1.962 1.814 1.905 1.891 1.855 2.024 1.992 1.657 1.934 1.685 1.807 1.650 1.717 1.571 1.580
1.707 1.915 1.740 1.741 1.934 1.676 2.048 1.735 1.618 1.745 1.998 1.620 1.871 2.035 1.973 1.464 1.753 2.115 1.828 1.970 1.769 1.948
2.166 1.677 1.782 1.844 1.895 1.686 1.847 1.808 1.587 1.983 1.431 1.895 1.772 2.060 1.797 1.820 1.695 1.828 2.023 1.951 1.802 1.655
1.474 1.689 1.791 2.000 1.831 2.283 1.750 1.836 2.061 1.867 1.886 1.497
La loro media?
1.84309
Dai, non è lontanissima dalla media di popolazione!
Certo, le singole stime, prese una per una, variano parecchio. E per vederlo possiamo calcolare la deviazione standard, cioè la radice quadrata della media dei quadrati delle differenze delle singole medie campionarie dalla loro media complessiva. Il risultato è
0.2003141
che rapportato ad 1.84, la media delle medie campionarie, non è proprio poco: il 10.87%.
Potevamo, anche, prendere un singolo campione di 10000 elementi: la sua media è:
mean(sample(population,10000))
1.87717
Bello, vero? Con un lavoro tutto sommato “piccolo” (vabbuò…) siamo riuscite e riusciti a stimare un’approssimazione dell’intera media di popolazione!
Quasi non sembra vero.
Del modo più appropriato per scegliere i campioni si è occupato, letteralmente, un mare di persone. E quando le cose stanno così, è facile trovare libri di testo (invece che paper accademici) già fatti, e magari molto antichi. Ve ne dico uno (per la sola ragione che ce l’ho – onestamente, comunque, è facile da capire e molto esteso, nonostante le dimensioni): Thompson, S. K. (1992) Sampling, Wiley.
Peccato, però, che le campagne di misura non garantiscono campioni casuali!
Prima di chiudere questa sezione, osserviamo una cosa. Ricordi le famose 100 medie campionarie calcolate ciascuna su un campione casuale di 100 elementi? Bene: immaginiamo di calcolare non 100, ma 10000 medie campionarie. Lo so: è uno spreco enorme, visto che 100 per 10000 vale 1000000, circa 5 volte di più della popolazione. Ma possiamo determinare la loro “distribuzione statistica empirica”, usando un diagramma chiamato “istogramma”.
Una “curva a campana”… Per adesso non dico nulla, salvo farti notare che questa forma è ricorrente nella pratica. E che è quasi sempre visibile, quando si adoperino campioni casuali.
Ma torniamo alle campagne di misura. Nel caso, queste iniziano ad una certa data, e durano un tot (prestabilito, in base di solito ai vincoli di budget).
Facciamo finta che le campagne durino 100 osservazioni – la stessa dimensione dei campioni casuali. Per generarne una, scegliamo a caso il punto di inizio, e determiniamo la media dei valori da lì in poi. Questi i risultati:
1.873 1.192 1.294 1.312 2.833 3.029 3.247 1.651 1.702 1.023 2.272 0.496 3.190 1.193 2.369 1.266 0.981 1.622 0.824 2.009 2.473 2.977
2.855 0.468 2.637 4.995 2.124 2.282 2.014 1.137 0.728 2.946 2.813 2.311 3.062 2.939 0.528 2.506 2.868 0.935 2.639 2.691 1.823 2.634
0.905 0.432 2.764 1.901 2.086 3.701 2.888 0.662 1.229 2.275 1.225 2.223 2.931 5.472 0.440 1.241 1.860 3.512 0.557 1.199 3.236 4.711
1.873 1.721 1.066 3.225 1.862 0.945 2.077 2.693 3.670 2.015 3.213 2.729 0.949 1.739 1.308 1.248 0.459 1.939 1.056 0.967 2.765 4.439
0.698 0.499 0.977 0.980 2.893 3.023 2.568 2.143 1.033 1.066 1.749 1.506
A prima vista non sembrano così diversi da quelli ottenuti dai campioni casuali: un’accozzaglia di 100 valori, apparentemente messi a casaccio, non diversamente da prima.
Solo che la loro media adesso vale
2.01306
(questa volta non ci siamo proprio). Mentre la deviazione standard è
1.06511
cioè, un’enormità rispetto a prima! E la distribuzione? È ancora una curva a campana? Vediamola, sempre con 10000 campioni:
Propri no. Niente curva a campana.
Magari non sembra, ma questa apparente piccolezza ha conseguenze profondissime sul modo di raccogliere e trattare i dati.
Sì, e allora?
Mica finisce qui.
Forse, abbiamo appena cominciato.
Nelle puntate successive vedremo altri tipi di medie, ma per il momento ci fermiamo qui, con una domanda provocatoria.
Facciamo finta di voler costruire un registratore automatico di dati (“data logger“) di velocità del vento, e di collegarlo ad un anemometro a dinamo. Per costruzione, questo eroga una tensione funzione della velocità angolare delle coppe, che per definizione varia con continuità.
Quindi possiamo compiere letture al rateo che vogliamo.
Poi, calcolare la media di queste letture (facciamo finta: ogni minuto).
Cosa puoi dire del risultato di questa operazione di media su un minuto? In particolare, è una media di popolazione, o campionaria?
E: puoi dirmi qualcosa su quanto è rappresentativa?
(Non è un esame: solo, un invito a riflettere su quanto abbiamo veduto).