AJO_1415_Ejemplo_regresión_dos_rectas

857 days ago by sevillad1415

Ejemplo de regresión lineal

 

A continuación calculamos la recta de regresión según la definición clásica (valores de $x$ exactos y de $y$ aproximados) y según la definición geométrica (valores de $x$ e $y$ igualmente inexactos).

puntos = [[0,0],[1,1],[2,1],[4,4]] scatter_plot(puntos,figsize=4,aspect_ratio=1) 
       

La siguiente función es la suma de cuadrados de residuos, que minimizaremos.Las variables son los coeficientes de la recta: $y=ax+b$.

Primero consideramos las coordenadas $x$ exactas, luego los residuos son los desplazamientos verticales (corregimos un punto llevándolo verticalmente hasta la recta).

var('a,b') f = (b-0)^2 + (a+b-1)^2 + (a*2+b-1)^2 + (a*4+b-4)^2 plot3d(f,(a,-2,4),(b,-4,3),figsize=4) 
       
Sleeping...
If no image appears re-execute the cell. 3-D viewer has been updated.

Como vemos la función a minimizar es un paraboloide, luego su mínimo será su punto crítico.

solve([f.diff(a),f.diff(b)],[a,b]) 
       
[[a == (34/35), b == (-1/5)]]
[[a == (34/35), b == (-1/5)]]

Dibujamos la recta de regresión. Esta es la recta que minimiza las desviaciones verticales.

r1 = 34/35*x-1/5 scatter_plot(puntos,figsize=4,aspect_ratio=1) + plot(r1,(x,0,4)) 
       

Ahora calculamos otra recta de regresión (en la literatura la recta de regresión es solamente la anterior).

Esta vez los residuos son las distancias geométricas (perpendiculares) de cada punto a la recta buscada. La fórmula de la distancia punto-recta está por ejemplo en la introducción de http://es.wikipedia.org/wiki/Distancia_de_un_punto_a_una_recta. Ojo, tenemos que elevar esas distancias al cuadrado en nuestra función.

var('a,b') f(a,b) = (a*0-0+b)^2/(a^2+1) + (a*1-1+b)^2/(a^2+1) + (a*2-1+b)^2/(a^2+1) + (a*4-4+b)^2/(a^2+1) 
       

Ahora la función no es un paraboloide (no es un polinomio de grado dos). Al calcular los puntos críticos nos salen dos, puedes verlos a continuación (el formato es distinto que antes).

s = solve([f.diff(a),f.diff(b)],[a,b],solution_dict=True) for ss in s: print ss 
       
{b: 35/272*sqrt(37)*sqrt(5) + 401/272, a: -5/68*sqrt(185) + 1/68}
{b: -35/272*sqrt(37)*sqrt(5) + 401/272, a: 5/68*sqrt(185) + 1/68}
{b: 35/272*sqrt(37)*sqrt(5) + 401/272, a: -5/68*sqrt(185) + 1/68}
{b: -35/272*sqrt(37)*sqrt(5) + 401/272, a: 5/68*sqrt(185) + 1/68}

Para ver cuál es el mínimo de los dos evaluamos $f$, resultando ser el segundo. En realidad habría que asegurarse razonando de que el mínimo tiene que estar en un punto crítico, pero supongámoslo por simplicidad.

for ss in s: f.subs(ss).n() 
       
17.3759190679597
0.374080932040348
17.3759190679597
0.374080932040348

En el gráfico siguiente tienes $f$, el punto mínimo (negro) y el otro (rojo).

plot3d(f,(a,-5,5),(b,-5,5)) + point3d([qq.subs(s[0]).n() for qq in [a,b,f]],color="red") + point3d([qq.subs(s[1]).n() for qq in [a,b,f]],color="black") 
       
Sleeping...
If no image appears re-execute the cell. 3-D viewer has been updated.

Tomamos esos valores encontrados, guardados en $s[1]$, para la nueva recta. Fíjate en que a continuación no nos hace falta copiar y pegar los números a mano; eso es debido a cómo se ha usado el  comando solve más arriba para el formato de las soluciones.


r2 = s[1][a]*x+s[1][b] scatter_plot(puntos,figsize=4,aspect_ratio=1) + plot(r2,(x,0,4),color="black") 
       

Para terminar comparamos las dos rectas.

scatter_plot(puntos,figsize=6,aspect_ratio=1) + plot(r1,(x,0,4)) + plot(r2,(x,0,4),color="black")