In [5]:
# Configuracion para recargar módulos y librerías
%reload_ext autoreload
%autoreload 2
# Si quieren graficos dentro del ipython notebook, sacar comentario a linea siguiente.
#%matplotlib inline
from IPython.core.display import HTML
HTML(open("style/mat281.css", "r").read())
Out[5]:
In [6]:
from mat281_code.lab import greetings
alumno_1 = ("Sebastian Flores", "2004001-7")
alumno_2 = ("Maria Jose Vargas", "2004007-8")
HTML(greetings(alumno_1, alumno_2))
Out[6]:
Este laboratorio utiliza la librería sklearn (oficialmente llamada scikit learn), puesto que buscamos aplicar la técnica de regresión lineal a datos tal como se haría en una aplicación real.
El laboratorio consiste de 2 secciones:
Generemos un conjunt de datos sintéticos para aplicar regresión lineal, utilizando numpy y sklearn. El modelo será inicialmente muy simple, del tipo
$$ y = \theta_0 + \theta_1 \ x_1 + \theta_2 x_2 + \varepsilon $$con $x_i$, $\theta_i$ e $y$ valores reales, y $\varepsilon$ el error de medición.
In [7]:
# Generando los datos
import numpy as np
# Model
m = 1000 # Number of data examples
true_coefficients = np.array([27., 5.])
true_intercept = 1982.
n = len(true_coefficients) + 1
# Data
nx1, nx2 = 21, 26
x1, x2 = np.meshgrid(np.linspace(-1.,1.,nx1), np.linspace(80.,120.,nx2))
X = np.array([x1.flatten(), x2.flatten()]).T
Y = true_intercept + np.dot(X, true_coefficients)
# Measurement contamination
measurement_error = np.random.normal(size=Y.shape)
# Mesured data
Y = Y + 200 * measurement_error/measurement_error.max()
Visualicemos los datos que hemos generado.
In [8]:
# Visualización de los datos
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# Plot of data
fig = plt.figure(figsize=(16,8))
ax = fig.add_subplot(111, projection='3d')
ax.scatter(X[:,0], X[:,1], Y, 'rs')
plt.xlabel("X[:,0]")
plt.ylabel("X[:,1]")
plt.title("y")
plt.show()
Apliquemos un modelo de regresión lineal, obteniendo el intercepto y coeficientes. Utilizando sklearn podemos obtener además la predicción bajo el modelo lineal de manera extremadamente sencilla.
In [9]:
# Aplicando regresion lineal con sklearn a los datos
from sklearn import linear_model
regr = linear_model.LinearRegression(fit_intercept=True, normalize=False)
regr.fit(X, Y) # No es necesario agregar la columna de 1s, el algooritmo lo hace por nosotros.
print "Intercepto hallado: ", regr.intercept_
print "Coeficientes hallados: ", regr.coef_
Y_pred = regr.predict(X) # Puedo verificar que valores daría el modelo en los puntos conocidos.
Visualicemos simultáneamente los datos y la superficie generada con la regresión lineal.
In [12]:
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# Plot of data
fig = plt.figure(figsize=(16,8))
ax = fig.add_subplot(111, projection='3d')
ax.scatter(X[:,0], X[:,1], Y, 'rs')
ax.plot_trisurf(X[:,0], X[:,1], Y_pred, color="r", alpha=0.5, lw=0)
plt.xlabel("X[:,0]")
plt.ylabel("X[:,1]")
plt.title("y")
plt.show()
Hmm, se ve bastante bien. ¿Qué tal los residuos? ¿Responden a nuestra hipótesis de distribución normal? Realicemos una exploración visual.
In [13]:
# Analisando gráfico de los residuos
res = Y_pred - Y
plt.figure(figsize=(16,8))
plt.hist(res)
plt.show()
El gráfico anterior parece suficientemente normal. Utilicemos un test de normalidad para verificar que efectivamente los residuos tengan una distribución normal.
In [14]:
# Analisis estadístico de residuos
import scipy.stats as stats
z, pval = stats.normaltest(res)
if pval < 0.05:
print "Not normal distribution"
print "p-value =", pval
Calculemos el error promedio en varias normas
In [19]:
print "Error promedio en norma 1:\t", np.linalg.norm(res,1)/len(res)
print "Error promedio en norma 2:\t", np.linalg.norm(res,2)/len(res)
print "Error promedio en norma inf:\t", np.linalg.norm(res,np.inf)
Utilice el código anterior para generar los mismos datos, pero esta vez no contamine los datos con errores gausianos, sino con elija utilizar un error de tipo gamma, uniforme o poisson (o bien, alguna otra de las posibles). Utilice alguno de los siguientes códigos:
gamma_measurement_error = np.random.gamma(1, size=Y.shape)
uniform_measurement_error = np.random.uniform(size=Y.shape)-0.5
weibull_measurement_error = np.random.weibull(1, size=Y.shape)
In [1]:
# 1.1 Desafio
# Codigo para generar datos con otro tipo de error de medicion
# y para realizar regresión lineal y análisis de los datos.
FIX ME: COMENTARIO AQUI
FIX ME: COMENTARIO AQUI
FIX ME: COMENTARIO AQUI
FIX ME: COMENTARIO AQUI
Un conjunto de datos bastante interesante es el Motor Dataset, que contiene características de 392 automobiles, entre los años 1970 y 1982.
Este conjunto de datos tiene los siguientes valores:
Como ya es tradición, exploraremos los datos utilizando head (o tail, o cat).
In [21]:
%%bash
head data/auto-mpg.data.txt
Realice una regresión lineal a los datos, buscando predecir el rendimiento (en mpg) utilizando como predictores: el desplazamiento, los caballos de fuerza, el peso y la aceleración, correspondientes a las columnas 2, 3, 4 y 5. Utilice el código provisto anteriormente para ajustar un modelo de regresión lineal, obtener los coeficientes del modelo, realizar predicciones y analice los residuos.
No es necesario incluir gráficos (pero tampoco está prohibido).
In [5]:
import numpy as np
# Lectura de datos
Y = np.loadtxt("data/auto-mpg.data.txt", delimiter=" ", usecols=(0,))
X = np.loadtxt("data/auto-mpg.data.txt", delimiter=" ", usecols=(2,3,4,5))
# Procesamiento de datos
FIX ME: COMENTARIO AQUI
FIX ME: COMENTARIO AQUI
FIX ME: COMENTARIO AQUI
Calculo el rendimiento en mpg tendría el Delorean, si en base a cierta información sabemos que posee:
Consejo: Utilice el metodo .predict() apropiadamente.
In [2]:
# Calcular el mpg aqui
mpg = 0
print mpg
FIX ME: COMENTARIO AQUI
