PHP: calcolare la retta di regressione lineare con i minimi quadrati

Uno strumento statistico molto usato in tutti i campi è la regressione lineare con il metodo dei minimi quadrati. Questa procedura serve a trovare una curva che interpoli al meglio una serie di dati campionari, che siano in una certa relazione fra di loro. Esaminiamo il caso più semplice, quello in cui la relazione sia lineare e quindi la migliore curva interpolante sia una retta.
Il metodo dei minimi quadrati si basa sul principio per il quale la migliore curva interpolante di un dato insieme di punti è la curva che ha la proprietà di avere minima la somma degli scarti quadratici, ovvero le differenze elevate al quadrato delle singole distanze fra i punti dati e i punti corrispondenti della retta interpolante.

Il procedimento per determinare i punti della nostra retta viene calcolato risolvendo:

sapendo che il termine noto si ricava da:

sostituendo q e sviluppando il quadrato si ottiene:

… a vederla così sembra bruttina, ma con le opportune sostituzioni rispetto a m:

otteniamo:

… adesso è più leggibile, è una equazione di una parabola. Trovare il valore minimo equivale a trovare il vertice. Qui bisognerebbe recuperare qualche reminiscenza sul calcolo delle derivate, comunque la soluzione è:

m, guarda caso… è anche il coefficiente angolare della nostra retta: quello che ci mancava!
Ora abbiamo in mano tutti i dati per costruire l’algoritmo.
Ecco quindi il codice della funzione, sviluppato in linguaggio PHP:

function regressione($X,$Y)
{
  if (!is_array($X) && !is_array($Y)) return false;
  if (count($X) <> count($Y)) return false;
  if (empty($X) || empty($Y)) return false;

  $regres = array();
  $n = count($X);
  $mx = array_sum($X)/$n; // media delle x
  $my = array_sum($Y)/$n; // media delle y
  $sxy = 0;
  $sxsqr = 0;

  for ($i=0;$i<$n;$i++){
    $sxy += ($X[$i] - $mx) * ($Y[$i] - $my);
    $sxsqr += pow(($X[$i] - $mx),2); // somma degli scarti quadratici medi
  }

  $m = $sxy / $sxsqr; // coefficiente angolare
  $q = $my - $m * $mx; // termine noto

  for ($i=0;$i<$n;$i++){
    $regres[$i] = $m * $X[$i] + $q;
  }

  return $regres;
}

La funzione restituisce l’array delle ordinate della retta di regressione. I dati campione si riferiscono al solito rapporto peso / altezza di alcune persone. Per una più chiara rappresentazione dei dati, è necessario realizzare un grafico con il modello “dispersione” per i dati campionari e il modello “retta” per i dati della regressione.

Ho utilizzato il plugin di JQuery: Flot, facendo generare server-side lo script necessario alla rappresentazione del grafico. L’esempio di applicazione di questa funzione e del grafico correlato è visibile in questa demo.

Per completezza occorre dire che bisognerebbe calcolare anche il coefficiente R² che fornisce indicazioni sulla qualità della correlazione rispetto ai dati. Il coefficiente varia tra 0 e 1 e, tanto più si avvicina a 1, tanto più i dati sono ben correlati.

Conclusioni:

Un breve ripasso di statistica e matematica delle superiori si è reso necessario per sviluppare un algoritmo molto usato, ma non sempre con cognizione di causa!

Riferimenti ed approfondimenti:

• Giorgio Meini – consulenza e formazione
• L’algoritmo dei minimi quadrati – Fortran – Guida di SuperEva
• Wikipedia: Regressione lineare
• Wikipedia: Metodo dei minimi quadrati
• Elena Loli Piccolomini: Approssimazione dati: minimi quadrati (pdf)
• Wikipedia: Derivata