The goal of this first notebook is to explore multiple regression and feature engineering with existing graphlab functions.
In this notebook you will use data on house sales in King County to predict prices using multiple regression. You will:
In [3]:
import graphlab
Dataset is from house sales in King County, the region where the city of Seattle, WA is located. this dataset you can download from this Link data and keep both unzip file and notebook in same place
In [4]:
sales = graphlab.SFrame('kc_house_data.gl')
sales.head(2) ## view first 2 rows of the data set
Out[4]:
In [5]:
train_data,test_data = sales.random_split(.8,seed=0)
Recall we can use the following code to learn a multiple regression model predicting 'price' based on the following features: example_features = ['sqft_living', 'bedrooms', 'bathrooms'] on training data with the following code:
(Aside: We set validation_set = None to ensure that the results are always the same)
In [6]:
example_features = ['sqft_living', 'bedrooms', 'bathrooms']
example_model = graphlab.linear_regression.create(train_data, target = 'price', features = example_features,
validation_set = None)
Now that we have fitted the model we can extract the regression weights (coefficients) as an SFrame as follows:
In [7]:
example_weight_summary = example_model.get("coefficients")
print example_weight_summary
In the gradient descent notebook we use numpy to do our regression. In this book we will use existing graphlab create functions to analyze multiple regressions.
Recall that once a model is built we can use the .predict() function to find the predicted values for data we pass. For example using the example model above:
In [8]:
example_predictions = example_model.predict(train_data)
print example_predictions[0] # should be 271789.505878
Now that we can make predictions given the model, let's write a function to compute the RSS of the model. Complete the function below to calculate RSS given the model, data, and the outcome.
Residual sum of squares (RSS) is an error metric for regression. GraphLab Create uses root mean squared error (RMSE). These are two common measures of error regression, and RMSE is simply the square root of the RSS
In [9]:
import numpy as np
def get_residual_sum_of_squares(model, data, outcome):
# First get the predictions
predictions = model.predict(data) # this the y heat and real output is ___outcome___
# Then compute the residuals/errors
# error = model.evaluate(data)
# rmse = graphlab.evaluation.rmse(outcome, predictions) ## root mean square error
diff = np.subtract(outcome,predictions)
# square the residuals and add them up
RSS = np.vdot(diff,diff)
#RSS = rmse
# Then square and add them up
return(RSS)
Test your function by computing the RSS on TEST data for the example model:
In [10]:
rss_example_train = get_residual_sum_of_squares(example_model, test_data, test_data['price'])
print rss_example_train # should be 2.7376153833e+14
Although we often think of multiple regression as including multiple different features (e.g. # of bedrooms, squarefeet, and # of bathrooms) but we can also consider transformations of existing features e.g. the log of the squarefeet or even "interaction" features such as the product of bedrooms and bathrooms.
You will use the logarithm function to create a new feature. so first you should import it from the math library.
In [11]:
from math import log
Next create the following 4 new features as column in both TEST and TRAIN data:
In [12]:
train_data['bedrooms_squared'] = train_data['bedrooms'].apply(lambda x: x**2)
test_data['bedrooms_squared'] = test_data['bedrooms'].apply(lambda x: x**2)
In [13]:
# create the remaining 3 features in both TEST and TRAIN data
train_data['bed_bath_rooms'] = train_data['bedrooms'] * train_data['bathrooms']
test_data['bed_bath_rooms'] = test_data['bedrooms'] * test_data['bathrooms']
train_data['log_sqft_living'] = train_data['sqft_living'].apply(lambda x: log(x))
test_data['log_sqft_living'] = test_data['sqft_living'].apply(lambda x: log(x))
train_data['lat_plus_long'] = train_data['lat'] + train_data['long']
test_data['lat_plus_long'] = test_data['lat'] + test_data['long']
What is the mean (arithmetic average) value of your 4 new features on TEST data? (round to 2 digits)
In [15]:
print sum(test_data['bedrooms_squared'])/len(test_data['bedrooms_squared'])
print sum(test_data['bed_bath_rooms'])/len(test_data['bed_bath_rooms'])
print sum(test_data['log_sqft_living'])/len(test_data['log_sqft_living'])
print sum(test_data['lat_plus_long'])/len(test_data['lat_plus_long'])
Now we will learn the weights for three (nested) models for predicting house prices. The first model will have the fewest features the second model will add one more feature and the third will add a few more:
In [16]:
model_1_features = ['sqft_living', 'bedrooms', 'bathrooms', 'lat', 'long']
model_2_features = model_1_features + ['bed_bath_rooms']
model_3_features = model_2_features + ['bedrooms_squared', 'log_sqft_living', 'lat_plus_long']
Now that you have the features, learn the weights for the three different models for predicting target = 'price' using graphlab.linear_regression.create() and look at the value of the weights/coefficients:
In [17]:
# Learn the three models: (don't forget to set validation_set = None)
model1 = graphlab.linear_regression.create(train_data, target = 'price', features = model_1_features,
validation_set = None)
model2 = graphlab.linear_regression.create(train_data, target = 'price', features = model_2_features,
validation_set = None)
model3 = graphlab.linear_regression.create(train_data, target = 'price', features = model_3_features,
validation_set = None)
In [18]:
# Examine/extract each model's coefficients:
model1_coefficients = model1.get("coefficients")
print model1_coefficients
model2_coefficients = model2.get("coefficients")
print model2_coefficients
model3_coefficients = model3.get("coefficients")
print model3_coefficients
What is the sign (positive or negative) for the coefficient/weight for 'bathrooms' in model 1?
What is the sign (positive or negative) for the coefficient/weight for 'bathrooms' in model 2?
Think about what this means.
First use your functions from earlier to compute the RSS on TEST Data for each of the three models.
In [19]:
# Compute the RSS on TESTING data for each of the three models and record the values:
rss_model1_test = get_residual_sum_of_squares(model1, test_data, test_data['price'])
print rss_model1_test
rss_model2_test = get_residual_sum_of_squares(model2, test_data, test_data['price'])
print rss_model2_test
rss_model3_test = get_residual_sum_of_squares(model3, test_data, test_data['price'])
print rss_model3_test
Which model (1, 2 or 3) has lowest RSS on TRAINING Data? Is this what you expected?
Now compute the RSS on on T data for each of the three models.
In [20]:
# Compute the RSS on TRAINING data for each of the three models and record the values:
rss_model1_train = get_residual_sum_of_squares(model1, train_data, train_data['price'])
print rss_model1_train
rss_model2_train = get_residual_sum_of_squares(model2, train_data, train_data['price'])
print rss_model2_train
rss_model3_train = get_residual_sum_of_squares(model3, train_data, train_data['price'])
print rss_model3_train
Quiz Question: Which model (1, 2 or 3) has lowest RSS on TESTING Data? Is this what you expected?Think about the features that were added to each model from the previous.