Ajuste de curvas 
mediante regresión lineal

Una exposición más extensa del tema de regresión lineal puede ser consultada en el libro Algebra Lineal, de Carlos Arce, William Castillo y Jorge González, tercera edición, Editorial de la Universidad de Costa Rica.
Por otra parte, antes de ejecutar los procedimientos que acompañan la siguiente exposición, por lo menos en lo relativo al manejo de applets, se debe consultar la última sección de este documento, donde se explican los requerimientos exigidos.

Regresión lineal: el caso de dos variables

El caso más simple de regresión lineal,  supone que salvo por un cierto error ϵ, la variable y puede ser explicada por la variable x,  mediante la relación:
          y  =  (Subscript[a, 0]) + (Subscript[a, 1]) x + ϵ .
 El interés se centra en estimar los parámetros (Subscript[a, 0]) y (Subscript[a, 1]), de manera que el error ϵ sea mínimo según un cierto criterio de calidad y a partir de un conjunto de valores observados para las variables  x y y :
                              ({{Valores (Subscript[x, i]):, (Subscript[x, 1]), (Subscript[x, 2]), ⋯, (Subscript[x, n])}, {Valores (Subscript[y, i]):, (Subscript[y, 1]), (Subscript[y, 2]), ⋯, (Subscript[x, n])}})

El siguiente gráfico ilustra posibles puntos ((Subscript[x, i]), (Subscript[y, i]))que se suponen conocidos y la recta de mejor ajuste y  =  (Subscript[a, 0]) + (Subscript[a, 1]) x, que se busca.

La solución más tradicional al problema usa resultados del cálculo diferencial para determinar los valores (Subscript[a, 0]) y (Subscript[a, 1]) que minizan la suma de errores: (∑_(i = 1)^n) ((((Subscript[y, i]) -  (Subscript[a, 0])  -  (Subscript[a, 1]) (Subscript[x, i]) ))^2).  Sin embargo ese no es el mejor enfoque que puede darse al tema.

Representación geométrica del problema en (ℝ^n)

Según la relación supuesta entre las variables  x y y,  los valores observados de éstas, en n distintas realizaciones o experimentos, se pueden escribir vectorialmente como : 
        (({{(Subscript[y, 1])}, {(Subscript[y, 2])}, {⋮}, {(Subscript[y, n])}}))=(Subscript[a, 0])(({{1}, {1}, {⋮}, {1}}))+(Subscript[a, 1])(({{(Subscript[x, 1])}, {(Subscript[x, 2])}, {⋮}, {(Subscript[x, n])}}))+(({{(Subscript[e, 1])}, {(Subscript[e, 2])}, {⋮}, {(Subscript[e, n])}})),
o equivalente,  si  y, (Subscript[1, n]), x, e son los vectores en  (ℝ^n) de esta relación:
                      y = X a + e    donde  X = ((Subscript[1, n])  x )  y  a = (({{(Subscript[a, 0])}, {(Subscript[a, 1])}})) . 
Y aunque el planteo sea en el espacio (ℝ^n), se puede proponer una ilustración gráfica al problema:

El objetivo, como ya se planteó para el espacio (ℝ^2), es estimar ((((Subscript[a, 0]), (Subscript[a, 1])))^t) de manera que ||e|| sea mínima, es decir, minimizando:      \(TraditionalForm\`\(\(||\)\(e\)\( || \^2\)\) = \(e·
            e\  = \ \(\((e\_1, e\_2, …, 
                e\_n)\)·\ \((e\_1, e\_2, …, 
                e\_n)\) = \ \(∑\+\(i = 1\)\%n e\_i\%2\  = \ ∑\+\(i = 1\)\%n\ \((y\_i - \ a\_0\  - \ a\_1\ x\_i\ )\)\^2\)\)\)\) .      
Consecuentemente  Xa  debe ser el vector en  W = 𝒞ℓ{(Subscript[1, n]), x} —el espacio generado por los vectores (Subscript[1, n]) y x—, más cercano a y, como se observa en el el gráfico, es decir: 
                    \(TraditionalForm\`\(\(y\)\(\ \)\)\&^ = \ \(Xa\  = \ Proy\_W\ y\)\),      donde \(TraditionalForm\`\(\(\(y\)\(\ \)\)\&^\)\)denota el vector de valores estimados para y .  
Las siguientes órdenes producen un applet que permite visualizar mejor el anterior gráfico y obtener la conclusión propuesta.

AppletViewer["Live", {"INPUT=" <> 
            ToString[
              InputForm[
                g1]], "WIDTH=500", "HEIGHT=500", 
"\
MAGNIFICATION=1.0"}];

Solución al problema

En el applet o gráfico se observa que como el vector de errores e = (y - X a) debe ser ortogonal a los vectores (Subscript[1, n]), (Subscript[x, 1]) que generan a W,  entonces:  
                         (X^t)(y-X a) = 0 ⇒ (X^t) y - (X^t)X a = 0 ⇒ a = ((((X^t) X))^(-1)) (X^t) y.   
 Donde a =  ((((Subscript[a, 0]), (Subscript[a, 1])))^t) y siempre que   (X^t) X sea invertible, lo cual siempre es así,  a menos que las componentes del vector x sean todas iguales.

El vector obtenido:   a = ((((X^t) X))^(-1)) (X^t) y  corresponde a un valor estimado de los parámetros del modelo y suele denotarse como \(\(a\&^\)\), sin embargo, esta notación  ha sido evitada para lograr una exposición con un simbolismo menos recargado. Similarmente, la notación no hace diferencia entre la variable x y el vector de valores observados para x, con el mismo propósito.

Indice de calidad

La calidad del ajuste se mide mediante el coeficiente de correlación lineal entre y y \(y\&^\) = X a, el cual, en términos de los vectores x, y y \(y\&^\) centrados,
          \(TraditionalForm\`x\_\(\(c\)\(\ \)\) = \ \((x\  - x\&_\ 1\_n)\)\),   \(TraditionalForm\`y\_c = \((y\  - y\&_\ 1\_n)\)\),  \(y\&^\_c\  = \ \((y\&^ - \(y\&^\)\&_\ 1\_n)\)\),
se expresa como:
                  \(r(y, y\&^)\) = \(\(\(||\)\(y\&^\_c\)\(||\)\)\/\(\(||\)\(y\_c\)\(||\)\)\) =((||(Subscript[a, 1])(\(x\  - x\&_\ 1\_n\))||)/(\(\(||\)\(y - \(y\&_\) 1\_n\)\(||\)\))) = cos(θ).

La ilustración anterior se justifica  porque el vector de errores observados e =y-X a = y - (Subscript[a, 0])(Subscript[1, n]) - (Subscript[a, 1])x debe  ser ortogonal al vector (Subscript[1, n]),  o sea,
       (Subscript[1, n])· (y- (Subscript[a, 0])(Subscript[1, n]) - (Subscript[a, 1])x)  = 0, 
 luego
          (Subscript[1, n])·y- (Subscript[a, 0])(Subscript[1, n])·(Subscript[1, n]) - (Subscript[a, 1])(Subscript[1, n])·x  = 0  ⇒  (∑_(i = 1)^n) (Subscript[y, i])  - 
                (Subscript[a, 0]) n  -  (Subscript[a, 1]) (∑_(i = 1)^n) (Subscript[x, i])  =  
            0,
de lo cual se deduce que:  \(y\&_\) = (Subscript[a, 0]) +  (Subscript[a, 1]) \(x\&_\).  
Por otra parte, 
       \(y\&^\) = Xa = (Subscript[a, 0])(Subscript[1, n]) + (Subscript[a, 1])x  =     ( \(y\&_\)  -  (Subscript[a, 1]) \(x\&_\) )(Subscript[1, n])   +   (Subscript[a, 1])x =   \(y\&_\ \ 1\_n\)   +   (Subscript[a, 1])(\(x - \(x\&_\) 1\_n\)) ,
y como \(\(y\&^\)\&_\) = \(TraditionalForm\`y\&_\)  (se verifica fácilmente)  entonces   
                        \(y\&^\_c\)  =   \(y\&^\  - \ \(y\&^\)\&_\ \ 1\_n\) = (Subscript[a, 1])\((x\  - \(x\&_\) 1\_n)\)

Ilustración gráfica del proceso de centraje de datos

El siguiente gráfico,  que es mejor apreciado usando el applet, permite visualizar  el gráfico anterior,  con los vectores 
    \(TraditionalForm\`x\_\(\(c\)\(\ \)\) = \ \((x\  - x\&_\ 1\_n)\)\),   \(TraditionalForm\`y\_c = \((y\  - y\&_\ 1\_n)\)\),  \(y\&^\_c\  = \ \((y\&^ - \(y\&^\)\&_\ 1\_n)\)\)
en conjunto con los vectores  x, y, \(y\&^ = X\ a\), no centrados.

AppletViewer["Live", {"INPUT=" <> 
            ToString[
              InputForm[
                g2]], "WIDTH=500", "HEIGHT=500", 
"\
MAGNIFICATION=1.4"}];

Observe por ejemplo que el vector \(y\&^\) = X a cuando se centra, es justamente la proyección del vector (Subscript[y, c]) sobre (Subscript[x, c]).  Lo cual se complementa con la siguiente solución alternativa al problema inicial de estimar los parámetros (Subscript[a, 0])y (Subscript[a, 1]).

Otra aproximación a la solución

Anteriormente se obtuvo que:  \(y\&_\) = (Subscript[a, 0]) +  (Subscript[a, 1]) \(x\&_\),  luego  restando los términos izquierdos y derechos de las siguientes igualdades entre vectores: 
              y = (Subscript[a, 0])(Subscript[1, n]) + (Subscript[a, 1])x + e   
       y    \(\(y\&_\) 1\_n\) = (Subscript[a, 0])(Subscript[1, n]) +  (Subscript[a, 1]) \(x\&_\) (Subscript[1, n])
se obtiene:     
          \(y - \(y\&_\) 1\_n\) = (Subscript[a, 1])(\(x\  - x\&_\ 1\_n\))+ e  
y se reduce que para determinar (Subscript[a, 1])de manera que se minimice  ||e||,  el problema es determinar el vector  (Subscript[a, 1])(\(x\  - x\&_\ 1\_n\)) más próximo   \(TraditionalForm\`y - \(y\&_\) 1\_n\).

Como el vector  (Subscript[a, 1])(\(x\  - x\&_\ 1\_n\)) más próximo a  \(TraditionalForm\`y - \(y\&_\) 1\_n\) es dado por \(TraditionalForm\`Proy\_\(x\  - x\&_\ 1\_n\)\ \((y - \(y\&_\) 
              1\_n)\)\), entonces:
                       (Subscript[a, 1]) = \(\(\((y - \(y\&_\) 1\_n)\)·\((x\  - 
                    x\&_\ 1\_n)\)\)\/\(\((x\  - 
                    x\&_\ 1\_n)\)·\((x\  - x\&_\ 1\_n)\)\)\),   
 por otra parte,   \(y\&_\) = (Subscript[a, 0]) +  (Subscript[a, 1]) \(x\&_\)   ⇒    (Subscript[a, 0]) = \(y\&_\)  -  (Subscript[a, 1]) \(x\&_\) .

Ejemplo 1

Determinar la recta, y = a + bx,    que mejor ajusta  los puntos (x,y)  en  la   tabla siguiente:

                   x    3   5   8  10 12 15  17  20  22  25 
                   y  12 10 15 20 26 18  28  32  24  26

Ilustración gráfica del problema

\(\(Pts = {{3, 12}, {5, 10}, {8, 15}, {10, 20}, {12, 26}, {15, 18}, {17, 
            28}, {20, 32}, 
		{22, 24}, {25, 26}};\)\)

\(\(graf1 = 
        ListPlot[Pts, Axes -> False, Frame -> True, 
          PlotRange -> {{2, 26}, {8, 
                34}}, 
		\ \ \ \ FrameTicks -> {{3, 5, 8, 10, 12, 15, 
                17, 20, 22, 25}, 
              Automatic}, 
		PlotStyle -> {AbsolutePointSize[3], 
              RGBColor[1, 0, 0]}];\)\)

Estimación de los valores a y b

El resultado del modelo de regresión lineal visto, para estimar los valores a y b es:
               (({{a}, {b}})) = ((((X^t) X))^(-1)) (X^t) y,   donde X = ((Subscript[1, n]) x).
Aunque en la anterior notación, se supone que (Subscript[1, n]) y x son vectores columnas, a continuación se trabajan como filas, pero para definir la matriz X transpuesta, que se denota Xt.

x  = {3, 5, 8, 10, 12, 15, 17, 20, 22, 25};
y  = {12, 10, 15, 20, 26, 18, 28, 32, 24, 26};

(Subscript[1, n])  =  Table[1, {10}]

Xt = {(Subscript[1, n]), x}

Naturalmente X =((((X^t)))^t) , o sea:

X = Transpose[Xt]

A continuación se hacen más cálculos de los estrictamente requeridos,  sólo para explorar los resultados intermedios.

Xt . X

Inverse[Xt . X]

{(Subscript[a, 0]), (Subscript[a, 1])} = N[Inverse[Xt . X] . Xt . y]

Una de las formas para calcular el índice de calidad, es primero calcular y centrar el vector \(TraditionalForm\`y\&^\) :

Yestimado  =  X . {(Subscript[a, 0]), (Subscript[a, 1])}

\(Observe\ que\ la\ media\ aritmética\ \ de\ \(\(y\)\(\ \)\)\&^\) = X(({{(Subscript[a, 0])}, {(Subscript[a, 1])}})) = (Subscript[a, 0]) (Subscript[1, n]) +(Subscript[a, 1]) x    es la media de y 
 o sea: \(y\&_\) =  (Subscript[a, 0]) + (Subscript[a, 1]) \(x\&_\) :

\(Yestimado\&_\  = \ Media[Yestimado]\)

(Subscript[a, 0])  +  (Subscript[a, 1]) Media[x]

\(Yestimado\_c = \ Yestimado\  - \ Yestimado\&_\ \ 1\_n\)

(Subscript[y, c])  =  y  - Media[y] (Subscript[1, n])

R = N[Norma[(Subscript[Yestimado, c])]/Norma[(Subscript[y, c])]]

Cálculos alternativos para (Subscript[a, 0]) y (Subscript[a, 1])

x  = {3, 5, 8, 10, 12, 15, 17, 20, 22, 25};

y  = {12, 10, 15, 20, 26, 18, 28, 32, 24, 26};

(Subscript[1, n])  =  Table[1, {10}]

\(x\&_\  = \ \ Media[x] // N\)

\(y\&_\  = \ Media[y] // N\)

\(Xc\  = \ x - \ x\&_\ 1\_n\)

\(Yc\  = \ y - \ y\&_\ 1\_n\)

(Subscript[a, 1]) = Yc . (Xc/Xc) . Xc

\(a\_0 = y\&_ - a\_1\ x\&_\)

Calidad del ajuste

R = N[Norma[(Subscript[a, 1]) Xc]/Norma[Yc]]

Representación gráfica

\(\(Pts = {{3, 12}, {5, 10}, {8, 15}, {10, 20}, {12, 26}, {15, 18}, {17, 
            28}, {20, 32}, 
		{22, 24}, {25, 26}};\)\)

(Subscript[a, 0])  +  (Subscript[a, 1]) t

\(\(PlotRL[Pts, a\_0 + a\_1\ t, {t, 2, 26}, 
        PlotRange -> {{2, 26}, {8, 34}}, 
	
        FrameTicks -> {{3, 5, 8, 10, 12, 15, 17, 20, 22, 25}, 
            Automatic}, 
	
        PlotStyle -> {AbsolutePointSize[3], 
            RGBColor[0, 0, 1]}];\)\)

Regresión lineal múltiple

Para el caso de que el modelo de regresión lineal contemple más de una variable explicativa, es decir, que salvo por un cierto error ϵ, p variables (Subscript[x, 1]), (Subscript[x, 2]), ..., (Subscript[x, p]) explican una variable y, mediante la relación:
            y  =  
        (Subscript[a, 0]) + (Subscript[a, 1]) (Subscript[x, 1]) + (Subscript[a, 2]) (Subscript[x, 2]) + ⋯ + (Subscript[a, p]) 
            (Subscript[x, p]) + ϵ
 se obtiene un resultado totalmente similar al ya visto.

Basados en n valores observados de las variables: (({{(Subscript[y, 1])}, {(Subscript[y, 2])}, {⋮}, {(Subscript[y, n])}}))=(Subscript[a, 0])(({{1}, {1}, {⋮}, {1}}))+(Subscript[a, 1])(({{(Subscript[x, 11])}, {(Subscript[x, 12])}, {⋮}, {(Subscript[x, 1  n])}}))+⋯+(Subscript[a, p])(({{(Subscript[x, p1])}, {(Subscript[x, p2])}, {⋮}, {(Subscript[x, pn])}}))+(({{(Subscript[e, 1])}, {(Subscript[e, 2])}, {⋮}, {(Subscript[e, n])}}));
que también se escriben en la forma:
                      y = X a + e
 donde a = ((((Subscript[a, 0]), (Subscript[a, 1]), …, (Subscript[a, p])))^t) es el vector de parámetros a estimar y X = ((Subscript[1, n]) (Subscript[x, 1]) ⋯ (Subscript[x, p])).

El objetivo es estimar ((((Subscript[a, 0]), (Subscript[a, 1]), …, (Subscript[a, p])))^t) de manera que ||e|| sea mínima:
      \(TraditionalForm\`\(\(||\)\(e\)\( || \^2\)\) = \(e·
            e\  = \ \(\((e\_1, e\_2, …, 
                e\_n)\)·\ \((e\_1, e\_2, …, 
                e\_n)\) = \ ∑\+\(i = 1\)\%n e\_i\%2\)\)\).      
Esto significa que  Xa  debe ser el vector en  W = 𝒞ℓ{(Subscript[1, n]), (Subscript[x, 1]), …, (Subscript[x, p])} más cercano a y, o sea: 
                    \(TraditionalForm\`\(\(y\)\(\ \)\)\&^ = \ \(Xa\  = \ Proy\_W\ y\)\).

Lo que permite determinar el vector de parámetros a en la siguiente forma: e = (y - X a) es ortogonal a los vectores (Subscript[1, n]), (Subscript[x, 1]), …, (Subscript[x, p]) que generan W, luego se tiene que:  
                         (X^t)(y-X a) = 0 ⇒ (X^t) y - (X^t)X a = 0 ⇒ a = ((((X^t) X))^(-1)) (X^t) y.

Indice de calidad

El coeficiente de correlación lineal entre y y \(y\&^\) = X a, mide la calidad con que el valor estimado de la variable y : \(y\&^\)  =(Subscript[a, 0])+(Subscript[a, 1])(Subscript[x, 1])+(Subscript[a, 2])(Subscript[x, 2])+⋯+(Subscript[a, p])(Subscript[x, p])aproxima el valor observado de y.
               \(TraditionalForm\`r(y, 
          y\&^)\  = \ \(\(||\)\(y\&^\_c\)\(||\)\)\/\(\(||\)\(y\_c\)\(||\)\)\).

AppletViewer["Live", {"INPUT=" <> 
            ToString[
              InputForm[
                g3]], "WIDTH=500", "HEIGHT=500", 
"\
MAGNIFICATION=1.4"}];

Ejemplo 2

Tomado del  libro Algebra Lineal,  de Arce, Castillo y  González:  ejercicio 5 del capítulo 9, en la tercera edición.
Los datos de la tabla siguiente representan el caudal promedio de Río Macho, para cada mes, medidos  en la Estación de Belén, en un estudio que abarcó 10 años.  (Datos suministrados por el ICE).

                            (m^3)/s                            (m^3)/s
      May.             1.88              Nov.        3.89 
      Jun.              3.5                 Dic.        3.38
      Jul.               3.84               Ene.       2.02
      Ago.             4.61               Feb.       1.35
      Set.               5.2                 Mar.       1.02
      Oct.               4.89              Abr.        1.14

  Determinar la  función cuadrática del mes (como enteros del 1 al 12) que mejor  explica el caudal del río.

Planteamiento del problema:

Aunque la tabla de datos dada sugiere que el mes 1 sea mayo, se va a elegir octubre como el mes 1, con lo que se obtiene un modelo con mejor índice de calidad.   Observe que como los datos corresponden a caudales promedios en varios años, los información dada es un año promedio por lo que el punto inicial de observación puede ser cualquiera de los meses.

x = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
y = {4.89, 3.89, 3.38, 2.02, 1.35, 1.02, 1.14, 1.88, 3.5, 3.84, 4.61, 
          5.2};

Puntos = Transpose[{x, y}]

\(\(\(
\)\(ListPlot[Puntos, 
        PlotStyle -> {AbsolutePointSize[
              3], 
			Red}];\)\)\)

\(El\ \ gráfico\ sugiere\ \ un\ polinomio\ de\ segundo\ grado\ , 
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ y = 
          a\_0 + \(a\_1\) x + a\_2\ x\^2 + e, \ 
        como\ modelo\ para\ explicar\ la\ evolución\ de\ \ la\ variable\ y\ \((caudal\ promedio\ del\ Rio\ Macho)\)\ en\ el\ tiempo\ \ \((variable\ x)\) . 

Usando\ el\ cambio\ de\ \(variables : \ x\_1\)\  = \ x\ \ , \ \ x\_2 = \ x\^2, \ 
        se\ obtiene\ un\ modelo\ \(lineal : 
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ y\)\  = a\_0 + a\_1\ x\_1 + a\_2\ x\_2 + e\)

En términos de  los valores observados de las variables:
   (({{(Subscript[y, 1])}, {(Subscript[y, 2])}, {⋮}, {(Subscript[y, 12])}}))=(Subscript[a, 0])(({{1}, {1}, {⋮}, {1}}))+(Subscript[a, 1])(({{1}, {2}, {⋮}, {12}}))+(Subscript[a, 2])(({{1}, {4}, {⋮}, {144}}))+(({{(Subscript[e, 1])}, {(Subscript[e, 2])}, {⋮}, {(Subscript[e, n])}})),  
donde   (Subscript[y, 1]) = 4.89,   (Subscript[y, 2]) = 3.89,  ⋯ , (Subscript[y, 12]) = 5.2

En resumen :   y=X a + e    con   a=((((Subscript[a, 0]), (Subscript[a, 1]), (Subscript[a, 2])))^t) el vector de parámetros a estimar y X=((Subscript[1, n]) (Subscript[x, 1]) (Subscript[x, 2])), \(y\ como\ resultado\ visto\ se\ tiene\ \(que : 
\ \ \ \ \ \ \ \ \ \ \ \ \((a\_0, a\_1, a\_2)\)\^t\)\  = \ \(\((\(X\^t\) X)\)\^\(-1\)\) \(X\^t\) \(\(y\)\(.\)\)\)

Estimación de los valores a y b

Xt = {{1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}, x, (x^2)};
X = Transpose[Xt]

a = {a0, a1, a2} = N[Inverse[Xt . X] . Xt . y]

Modelo obtenido:

"y"  == a . {1, t, (t^2)}

Calidad del ajuste:

Ye  =  X . a

R  =  N[Norma[Ye - Media[Ye]]/Norma[y - Media[y]]]

Ilustración gráfica:

Puntos = Transpose[{x, y}];
\(\(PlotRL[Puntos, a . {1, t, t\^2}, {t, 0, 13}, 
        PlotStyle -> {AbsolutePointSize[3], 
			RGBColor[1, 0, 
              0]}];\)\)

Otro modelo posible:

Si el mes 1 se elige como de  Junio o Julio, los meses de transición de verano a invierno, los datos parecen explicarse mejor mediante un polinomio de tercer grado:
                 y = (Subscript[a, 0]) + (Subscript[a, 1]) x + (Subscript[a, 2]) (x^2) + (Subscript[a, 3]) (x^3) + e.

x = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
y = {3.84, 4.61, 5.2, 4.89, 3.89, 3.38, 2.02, 1.35, 1.02, 1.14, 1.88, 
          3.5};

Puntos = Transpose[{x, y}]

Puntos = Transpose[{x, y}];
\(\(ListPlot[Puntos, 
        PlotStyle -> {AbsolutePointSize[3], 
			Red}];\)\)

\(TextForm\`Se\ propone\ como\ ejercicio, \ 
    plantear\ el\ modelo\ de\ regresión\ lineal\ múltiple\ requerido, \ 
    estimar\ los\ parámetros\ respectivos\ e\ ilustrar\ gráficamente\ el\ ajuste\ de\ los\ datos\ por\ la\ curva\ \(\(obtenida\)\(.\)\)\)

Primitiva Fit

Mathematica aporta recursos para estimar los parámetros de modelos como los anteriores y más. Se sugiere consultar la información proporcionada en el "Help" sobre el procedimiento Fit.

Por ejemplo para resolver el problema planteado en el ejemplo 2, basta dar el siguiente comando:

x = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
y = {4.89, 3.89, 3.38, 2.02, 1.35, 1.02, 1.14, 1.88, 3.5, 3.84, 4.61, 
          5.2};

Puntos = Transpose[{x, y}]

Fit[Puntos, {1, t, (t^2)}, t]

Recursos para el manejo de gráficos

A continuación se presenta una serie de recursos que han sido usados en las exposiciones sobre Regresión Lineal y los cálculos requeridos,  sin embargo, no son de interés para el proceso de enseñanza y aprendizaje del tema y podrían convertirse en un distractor.  Por ello se han agrupado en esta sección y sus  comandos se han declarado en celdas de inicialización, para que se ejecuten al inicio de la sesión  sin que sea necesario ocuparse ellas.

Las siguientes dos bibliotecas son requeridas, para la construcción de gráficos que involucran el uso de colores apartir de sus nombres en inglés, y flechas en el plano.
Ambas bibliotecas pertenecen al conjunto de bibliotecas estándar de Mathematica.

Needs["Graphics`Colors`"];
Needs["Graphics`Arrow`"]

La siguiente orden requiere haber creado la carpeta AlgebraLineal en Archivos de programa/Wolfram Research/Mathematica/4.2/AddOns/Applications y copiar el archivo Arrow3D.m en dicha carpeta.  Esto habilita para usar los comandos Arrow3D, Vector3D y DescribeVector3D, para constuir flechas en el espacio.

Needs["AlgebraLineal`Arrow3D`"]

Para preparar el uso del applet live.jar deben darse las siguientes órdenes, cambiando  la ruta del directorio de trabajo a su situación particular, el cual además debe contener el archivo live.jar.

Needs["JLink`"];
InstallJava[];
AddToClassPath["C:\Cursos\MA1004\Ciclo41"];

live.jar es un applet escrito por Martin Kraus, http://wwwvis.informatik.uni-stuttgart.de/~kraus/LiveGraphics3D/, la habilidad de Mathematica para correr programas Java usando J/Link, hace posible ver los gráficos de Mathematica en ventanas interactivas de Java que permiten rotar las imágenes con simples movimientos del mouse.

Otros procedimientos, que se usan en la construcción de gráficos son:

Vector[v_, P_: {0, 0}]  :=  
        Arrow[P, P + v];

    Norma[u_]  :=  (Sqrt[u . u]);     

    Unitario[u_]  :=  u/Norma[u];

    Orto2D[{x_, y_}]  :=  {y, -x};

Y los siguientes procedimientos se emplean sólo en este capítulo de regresión lineal.

Media[x_] := Apply[Plus, x]/Length[x];
Centra[x_] := x - Media[x];
\(PlotRL[Pts_, curva_, {t_, a_, b_}, Opc___] := 
      Module[{g1, g2, ejes, marco, 
          estilo}, 
		{ejes, marco, 
            estilo} = \({Axes, Frame, PlotStyle} /. {Opc}\) /. {Axes -> 
                False, Frame -> 
                True, 
					PlotStyle -> {AbsolutePointSize[2], 
                  RGBColor[1, 0, 0]}}; 
		g2 = 
          Plot[curva, {t, a, b}, DisplayFunction -> Identity, 
            Axes -> ejes]; 
		g1 = 
          ListPlot[Pts, DisplayFunction -> Identity, Axes -> ejes, 
            Frame -> marco, 
				PlotStyle -> estilo, 
            Opc]; 
		Show[g1, g2, 
          DisplayFunction -> $DisplayFunction]]\)

Construcción de gráficos ilustrativos

\(\(g1\  = \ 
        Plot[a + b\ t, {t, 0, 
            26}, 
PlotRange -> {0, 
              33}, 
Ticks -> None, 
	
          PlotStyle -> {RGBColor[1, 0, 
                1]}, 
AxesOrigin -> {0, 
              0}, 
Epilog -> {Text[\*""y = (Subscript[a, 0])+(Subscript[a, 1])x"", {18, 21}]}];\)\)

\(\(g2\  = \ 
        ListPlot[Pts, 
	
          PlotStyle -> {AbsolutePointSize[3], RGBColor[0, 0, 1]}, 
          Frame -> False, 
          Axes -> True
];\)\)

Show[g1, g2, AxesLabel -> {"x", "y"}];

Module[{v1, x, y, u1, x1n, u2, u3, u3p, 
          er}, 
v1  =  {2, 2, 0};  
        x  =  {3, -1, 1};  
        y =  {5, 1, 4}; 
u1  =  
          N[v1/(Sqrt[v1 . v1])]; 
x1n  =  (x . 
                u1) u1; 
u2  =  
          Unitario[x - x1n]; 
u3p  =  
          y . u1 u1  +  y . u2 u2; 
er = 
          y - u3p;  
g1  = 
          Show[Graphics3D[{AbsolutePointSize[
                  3], 
{Vector3D[x], Vector3D[y], 
                  Text["x", 0.7  x, {0, 1}], 
                  Text["y", 
                    0.8  y, {1, -1}]}, 
{Magenta, 
                  Vector3D[u1], 
Text["1n", 
                    0.8  u1, {0, -1}]}, 
{Green, 
                  Vector3D[u3p], 
                  Vector3D[er, u3p], 
Text["Xa", 
                    0.8  u3p, {0, 1}], 
Text["e", 
                    u3p + 0.8  er, {-1, 0}]}, 
{Cyan, 
                  Line[{-0.3 u1 - 0.4  u2, -0.3 u1 + 6  u2, 
                      6  u1 + 6  u2, 
                      6  u1 - 0.4  u2, -0.3 u1 - 
                        0.4  u2}], 
Text["W", 
                    1.3  u1 + 4.7  u2]}
}, 
              Boxed -> False]
]];

\(\(Show[
        Graphics[{{Blue, Thickness[0.01], 
              Vector[{3, 0}], 
Vector[{4, 
                  4/3}]}, 
{Black, 
              Thickness[0.005], 
Vector[{4, 0}], 
              Vector[{0, 4/3}, {4, 
                  0}]}, 
Text[\*""y-\!\(y\&_\)(Subscript[1, n])"", {3.6, 
                1.2}, {1, \(-1\)}], 
Text[\*""x-\!\(x\&_\)(Subscript[1, n])"", {2.6, 0}, {1, 1}], 
Text["e", {4, 
                0.9}, {\(-1\), 
                0}], 
Text[\*""(Subscript[a, 1])(x-\!\(x\&_\)(Subscript[1, n]))"", {3.9, 0}, {1, 1}]
}, 
          AspectRatio -> Automatic, 
          PlotRange -> All]];\)\)

\(\(Show[
        Graphics[{{Blue, Thickness[0.01], 
              Vector[{3, 0}], 
Vector[{4, 
                  4/3}]}, 
{Black, 
              Thickness[0.005], 
Vector[{4, 0}], 
              Vector[{0, 4/3}, {4, 
                  0}]}, 
Text[\*""y-\!\(y\&_\)(Subscript[1, n])"", {3.6, 
                1.2}, {1, \(-1\)}], 
Text[\*""x-\!\(x\&_\)(Subscript[1, n])"", {2.6, 0}, {1, 1}], 
Text["e", {4, 
                0.9}, {\(-1\), 
                0}], 
Text[\*""(Subscript[a, 1])(x-\!\(x\&_\)(Subscript[1, n]))"", {3.9, 0}, {1, 1}], 
{Red, 
              Circle[{0, 0}, 
                1, {0, ArcTan[3, 
                    1]}], 
Text["θ", {1.1, 
                  0.2}]}
}, AspectRatio -> Automatic, 
          PlotRange -> All]];\)\)

Module[{v1, x, y, u1, xc, yc, x1n, y1n, u2, u3, u3p, 
          er}, 
v1  =  {2, 2, 0};  
        x  =  {3, -1, 1};  
        y =  {5, 1, 4}; 
u1  =  
          N[v1/(Sqrt[v1 . v1])]; 
x1n  =  (x . u1) u1; 
        xc  =  x  - x1n; 
y1n  =  (y . u1) u1;  
        yc  =  y - y1n; 
u2  =  
          Unitario[x - x1n]; 
u3p  =  
          y . u1 u1  +  y . u2 u2; 
er = 
          y - u3p;  
uu  =  
          yc . u2 u2; 
g2  = 
          Show[Graphics3D[{AbsolutePointSize[
                  3], 
{Vector3D[x], Vector3D[y], 
                  Text["x", 0.7  x, {0, 1}], 
                  Text["y", 
                    0.8  y, {1, -1}]}, 
{Blue, 
                  Vector3D[xc], 
Vector3D[
                    yc], 
Text["Xc", 
                    0.8  xc, {1, 0}], 
Text["Yc", 
                    0.8  yc, {1, 0}]}, 
{Magenta, 
                  Vector3D[u1], 
                  Vector3D[x - x1n, x1n], 
Vector3D[x1n], 
                  Vector3D[y1n], 
Vector3D[yc, 
                    y1n], 
Text["1n", 
                    0.8  u1, {0, -1}], 
Text["Xc", 
                    x1n + 0.8  xc, {-1, 0}]}, 
{Green, 
                  Vector3D[u3p], Vector3D[er, u3p], 
                  Vector3D[uu], 
Vector3D[er, 
                    uu], 
Text["Xa", 
                    0.8  u3p, {0, -1}], 
Text["e", 
                    u3p + 0.8  er, {-1, 0}]}, 
{Green, 
                  Dashing[{0.01, 0.01}], 
                  Line[{uu, u3p, y1n}]}, 
{Cyan, 
                  Line[{-0.3 u1 - 0.4  u2, -0.3 u1 + 6  u2, 
                      5  u1 + 6  u2, 
                      5  u1 - 0.4  u2, -0.3 u1 - 
                        0.4  u2}], 
Text["W", 
                    1.3  u1 + 4.7  u2]}
}, 
              Boxed -> False]
]];

Module[{v1, x, y, x1n, u1, u2, u3p, 
          er}, 
v1  =  {2, 2, 0};  
        x  =  {3, -1, 1};  
        y =  {6, -1, 4.5}; 
u1  =  
          N[v1/(Sqrt[v1 . v1])]; 
x1n  =  (x . 
                u1) u1; 
u2  =  
          Unitario[x - x1n]; 
u3p  =  
          y . u1 u1  +  y . u2 u2; 
er = 
          y - u3p;  
g1  = 
          Show[Graphics3D[{AbsolutePointSize[
                  3], 
{Vector3D[y], 
                  Text["y", 
                    0.8  y, {1, -1}]}, 
{Green, 
                  Vector3D[u3p], 
                  Vector3D[er, u3p], 
Text["e", 
                    u3p + 0.8  er, {-1, 
                      0}], 
Text["Xa", 
                    0.8  u3p, {0, 1}]}, 
{Magenta, 
                  Line[{-0.3 u1 - 0.6  u2, -0.3 u1 + 7  u2, 
                      5  u1 + 7  u2, 
                      5  u1 - 0.6  u2, -0.3 u1 - 
                        0.6  u2}], 
Text["W", 
                    0.6  u1 + 5.5  u2]}
}, Boxed -> False, 
              PlotRange -> All]
]];

Module[{v1, x, y, x1n, y1n, xc, yc, u1, u2, u3p, uu, 
          er}, 
v1  =  {2, 2, 0};  
        x  =  {3, -1, 1};  
        y =  {6, -1, 4.5}; 
u1  =  
          N[v1/(Sqrt[v1 . v1])]; 
x1n  =  (x . u1) u1; 
        xc  =  x  - x1n; 
y1n  =  (y . u1) u1;  
        yc  =  y - y1n; 
u2  =  
          Unitario[x - x1n]; 
u3p  =  
          y . u1 u1  +  y . u2 u2; 
er = 
          y - u3p;  
uu  =  
          yc . u2 u2; 
g3  = 
          Show[Graphics3D[{AbsolutePointSize[
                  3], 
{Vector3D[y], 
                  Text["y", 
                    0.8  y, {1, -1}]}, 
{Green, 
                  Vector3D[u3p], 
                  Vector3D[er, u3p], 
Text["e", 
                    u3p + 0.8  er, {-1, 
                      0}], 
Text["Xa", 
                    0.8  u3p, {0, 1}], 
Vector3D[
                    uu], 
Vector3D[er, 
                    uu]}, 
{Cyan, 
                  Line[{-0.7 u1 - 0.6  u2, -0.7 u1 + 8  u2, 
                      5  u1 + 8  u2, 
                      5  u1 - 0.6  u2, -0.7 u1 - 
                        0.6  u2}], 
Text["W", 
                    0.8  u1 + 5.5  u2]}, 
{Magenta, 
                  Vector3D[u1], 
Vector3D[
                    y1n], 
Vector3D[yc, 
                    y1n], 
Line[{y1n, 
                      u3p}], 
Text["1n", 
                    0.8  u1, {0, -1}]}, 
{Blue, 
                  Vector3D[2.5  xc], 
Vector3D[
                    yc], 
Text["Wc", 
                    2.4  xc, {1, 0}], 
Text["Yc", 
                    0.9  yc, {1, 0}]}
}, Boxed -> False, 
              PlotRange -> All]
]];