multiple-regression-gradient-descent.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multiple Regression (gradient descent)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fire up graphlab create"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import graphlab"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load in house sales data\n",
    "\n",
    "Dataset is from house sales in King County, the region where the city of Seattle, WA is located."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "[INFO] graphlab.cython.cy_server: GraphLab Create v2.1 started. Logging: /tmp/graphlab_server_1529355203.log\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "This non-commercial license of GraphLab Create for academic use is assigned to nanlee_89@yahoo.com and will expire on December 07, 2018.\n"
     ]
    }
   ],
   "source": [
    "sales = graphlab.SFrame('kc_house_data.gl/')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Convert SFrame to Numpy Array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np # note this allows us to refer to numpy as np instead "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def get_numpy_data(data_sframe, features, output):\n",
    "    data_sframe['constant'] = 1 # this is how you add a constant column to an SFrame\n",
    "    # add the column 'constant' to the front of the features list so that we can extract it along with the others:\n",
    "    features = ['constant'] + features # this is how you combine two lists\n",
    "    # select the columns of data_SFrame given by the features list into the SFrame features_sframe (now including constant):\n",
    "    features_sframe = data_sframe[features]\n",
    "    \n",
    "    # the following line will convert the features_SFrame into a numpy matrix:\n",
    "    feature_matrix = features_sframe.to_numpy()\n",
    "    # assign the column of data_sframe associated with the output to the SArray output_sarray\n",
    "    output_sarray = data_sframe[output]\n",
    "    \n",
    "    # the following will convert the SArray into a numpy array by first converting it to a list\n",
    "    output_array = output_sarray.to_numpy()\n",
    "    return(feature_matrix, output_array)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Test the function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.00e+00 1.18e+03]\n",
      "221900.0\n"
     ]
    }
   ],
   "source": [
    "(example_features, example_output) = get_numpy_data(sales, ['sqft_living'], 'price') # the [] around 'sqft_living' makes it a list\n",
    "print example_features[0,:] # this accesses the first row of the data the ':' indicates 'all columns'\n",
    "print example_output[0] # and the corresponding output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Predicting output given regression weights"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Suppose we had the weights [1.0, 1.0] and the features [1.0, 1180.0] and we wanted to compute the predicted output 1.0\\*1.0 + 1.0\\*1180.0 = 1181.0 this is the dot product between these two arrays. If they're numpy arrayws we can use np.dot() to compute this:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1181.0\n"
     ]
    }
   ],
   "source": [
    "my_weights = np.array([1., 1.]) # the example weights\n",
    "my_features = example_features[0,] # we'll use the first data point\n",
    "predicted_value = np.dot(my_features, my_weights)\n",
    "print predicted_value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Function to predict output given feature matrix and weight vector"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def predict_output(feature_matrix, weights):\n",
    "    # assume feature_matrix is a numpy matrix containing the features as columns and weights is a corresponding numpy array\n",
    "    # create the predictions vector by using np.dot()\n",
    "    predictions =  np.dot(feature_matrix, weights)\n",
    "    return(predictions)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Test the prediction function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1181.0\n",
      "2571.0\n"
     ]
    }
   ],
   "source": [
    "test_predictions = predict_output(example_features, my_weights)\n",
    "print test_predictions[0] # should be 1181.0\n",
    "print test_predictions[1] # should be 2571.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Computing the Derivative"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the derivative of a sum is the sum of the derivatives we can compute the derivative for a single data point and then sum over data points. We can write the squared difference between the observed output and predicted output for a single point as follows:\n",
    "\n",
    "(w[0]\\*[CONSTANT] + w[1]\\*[feature_1] + ... + w[i] \\*[feature_i] + ... +  w[k]\\*[feature_k] - output)^2\n",
    "\n",
    "Where we have k features and a constant. So the derivative with respect to weight w[i] by the chain rule is:\n",
    "\n",
    "2\\*(w[0]\\*[CONSTANT] + w[1]\\*[feature_1] + ... + w[i] \\*[feature_i] + ... +  w[k]\\*[feature_k] - output)\\* [feature_i]\n",
    "\n",
    "The term inside the paranethesis is just the error (difference between prediction and output). So we can re-write this as:\n",
    "\n",
    "2\\*error\\*[feature_i]\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def feature_derivative(errors, feature):\n",
    "    # Assume that errors and feature are both numpy arrays of the same length (number of data points)\n",
    "    # compute twice the dot product of these vectors as 'derivative' and return the value\n",
    "    derivative = 2 * np.dot(errors, feature)\n",
    "    return(derivative)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### test feature derivartive function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-23345850022.0\n",
      "-23345850022.0\n"
     ]
    }
   ],
   "source": [
    "(example_features, example_output) = get_numpy_data(sales, ['sqft_living'], 'price') \n",
    "my_weights = np.array([0., 0.]) # this makes all the predictions 0\n",
    "test_predictions = predict_output(example_features, my_weights) \n",
    "# just like SFrames 2 numpy arrays can be elementwise subtracted with '-': \n",
    "errors = test_predictions - example_output # prediction errors in this case is just the -example_output\n",
    "feature = example_features[:,0] # let's compute the derivative with respect to 'constant', the \":\" indicates \"all rows\"\n",
    "derivative = feature_derivative(errors, feature)\n",
    "print derivative\n",
    "print -np.sum(example_output)*2 # should be the same as derivative"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Gradient Descent"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here is a function that performs a gradient descent. Given a starting point we update the current weights by moving in the negative gradient direction. The gradient is the direction of *increase* and therefore the negative gradient is the direction of *decrease* and we're trying to *minimize* a cost function. \n",
    "\n",
    "The amount by which we move in the negative gradient *direction*  is called the 'step size'. We stop when we are 'sufficiently close' to the optimum. We define this by requiring that the magnitude (length) of the gradient vector to be smaller than a fixed 'tolerance'. For each step in the gradient descent we update the weight for each feature befofe computing our stopping criteria."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from math import sqrt # recall that the magnitude/length of a vector [g[0], g[1], g[2]] is sqrt(g[0]^2 + g[1]^2 + g[2]^2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def regression_gradient_descent(feature_matrix, output, initial_weights, step_size, tolerance):\n",
    "    converged = False \n",
    "    weights = np.array(initial_weights) # make sure it's a numpy array\n",
    "    while not converged:\n",
    "        # compute the predictions based on feature_matrix and weights using your predict_output() function\n",
    "        predictions = predict_output(feature_matrix, weights)\n",
    "        \n",
    "        # compute the errors as predictions - output\n",
    "        errors = predictions - output\n",
    "\n",
    "        gradient_sum_squares = 0 # initialize the gradient sum of squares\n",
    "        # while we haven't reached the tolerance yet, update each feature's weight\n",
    "        for i in range(len(weights)): # loop over each weight\n",
    "            # Recall that feature_matrix[:, i] is the feature column associated with weights[i]\n",
    "            # compute the derivative for weight[i]:\n",
    "            derivative = feature_derivative(errors, feature_matrix[:, i])\n",
    "\n",
    "            # add the squared value of the derivative to the gradient sum of squares (for assessing convergence)\n",
    "            squared_derivative = derivative * derivative \n",
    "            gradient_sum_squares = gradient_sum_squares + squared_derivative\n",
    "            \n",
    "            # subtract the step size times the derivative from the current weight\n",
    "            weights[i] = weights[i] - step_size * derivative\n",
    "            \n",
    "        # compute the square-root of the gradient sum of squares to get the gradient magnitude:\n",
    "        gradient_magnitude = sqrt(gradient_sum_squares)\n",
    "        if gradient_magnitude < tolerance:\n",
    "            converged = True\n",
    "    return(weights)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the gradient is a sum over all the data points and involves a product of an error and a feature the gradient itself will be very large since the features are large (squarefeet) and the output is large (prices). So while you might expect \"tolerance\" to be small, small is only relative to the size of the features. \n",
    "\n",
    "For similar reasons the step size will be much smaller than you might expect but this is because the gradient has such large values."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Running the Gradient Descent as Simple Regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Split the data into training and test data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "train_data,test_data = sales.random_split(.8,seed=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Although the gradient descent is designed for multiple regression since the constant is now a feature we can use the gradient descent function to estimat the parameters in the simple regression on squarefeet. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Set up the feature_matrix, output, initial weights and step size for the first model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# test out the gradient descent\n",
    "simple_features = ['sqft_living']\n",
    "my_output = 'price'\n",
    "(simple_feature_matrix, output) = get_numpy_data(train_data, simple_features, my_output)\n",
    "initial_weights = np.array([-47000., 1.])\n",
    "step_size = 7e-12\n",
    "tolerance = 2.5e7"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Run gradient descent with the above parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "weights = regression_gradient_descent(simple_feature_matrix, output, initial_weights, step_size, tolerance)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The value of the weight for sqft_living:  [-46999.88716555    281.91211912]\n"
     ]
    }
   ],
   "source": [
    "print \"The value of the weight for sqft_living: \", weights"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use your newly estimated weights and your predict_output() function to compute the predictions on all the TEST data (you will need to create a numpy array of the test feature_matrix and test output first:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "(test_simple_feature_matrix, test_output) = get_numpy_data(test_data, simple_features, my_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compute predictions using test_simple_feature_matrix and weights from above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "predictions_test = predict_output(test_simple_feature_matrix, weights) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The predicted price for the 1st house in the TEST data set for model 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "356134.4431709297"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "predictions_test[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compute the RSS on the test data set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def get_residual_sum_of_squares(predictions, output):\n",
    "\n",
    "    residuals = output - predictions\n",
    "\n",
    "    # Then square and add them up\n",
    "    squared_residuals = residuals*residuals\n",
    "    \n",
    "    RSS = squared_residuals.sum()\n",
    "\n",
    "    return(RSS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "275400047593155.94"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "get_residual_sum_of_squares(predictions_test, test_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Running a multiple regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we will use more than one actual feature. Use the following code to produce the weights for a second model with the following parameters:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Produce the weights for a second model with the following parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model_features = ['sqft_living', 'sqft_living15'] # sqft_living15 is the average squarefeet for the nearest 15 neighbors. \n",
    "my_output = 'price'\n",
    "(feature_matrix, output) = get_numpy_data(train_data, model_features, my_output)\n",
    "initial_weights = np.array([-100000., 1., 1.])\n",
    "step_size = 4e-12\n",
    "tolerance = 1e9"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Use the above parameters to estimate the model weights"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "model_weights = regression_gradient_descent(feature_matrix, output, initial_weights, step_size, tolerance)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Use newly estimated weights and the predict_output function to compute the predictions on the TEST data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "(test_model_feature_matrix, test_output) = get_numpy_data(test_data, model_features, my_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The predicted price for the 1st house in the TEST data set for model 2 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "predictions_test_model_2 = predict_output(test_model_feature_matrix, model_weights) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "366651.4120365591"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "predictions_test_model_2[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The actual price for the 1st house in the test data set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "310000.0"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "test_data['price'][0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compute the RSS for model 2 on TEST data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "270263446465244.06"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "get_residual_sum_of_squares(predictions_test_model_2, test_output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compare the performance of the two model on 1st house on the Test data set"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "** Which estimate was closer to the true price for the 1st house on the TEST data set, model 1 or model 2?**\n",
    "\n",
    "The simple model or model 1 is more closer to the actual price.¶\n",
    "\n",
    "** Which model (1 or 2) has lowest RSS on all of the TEST data? **\n",
    "\n",
    "The multiple regression model has lower RSS than Simple model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [default]",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}