{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# This is a basic Starter Kernel for the New York City Taxi Fare Prediction Playground Competition \nHere we'll use a simple linear model based on the travel vector from the taxi's pickup location to dropoff location which predicts the `fare_amount` of each ride.\n\nThis kernel uses some `pandas` and mostly `numpy` for the critical work.  There are many higher-level libraries you could use instead, for example `sklearn` or `statsmodels`.  ","metadata":{"_uuid":"b4578d48b219735043a4d2102119fb307d2fc83f"}},{"cell_type":"code","source":"# Initial Python environment setup...\nimport numpy as np # linear algebra\nimport pandas as pd # CSV file I/O (e.g. pd.read_csv)\nimport os # reading the input files we have access to\n\nprint(os.listdir('../input'))","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-08-11T23:04:39.169785Z","iopub.execute_input":"2022-08-11T23:04:39.170066Z","iopub.status.idle":"2022-08-11T23:04:39.185531Z","shell.execute_reply.started":"2022-08-11T23:04:39.169998Z","shell.execute_reply":"2022-08-11T23:04:39.185088Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Setup training data\nFirst let's read in our training data.  Kernels do not yet support enough memory to load the whole dataset at once, at least using `pd.read_csv`.  The entire dataset is about 55M rows, so we're skipping a good portion of the data, but it's certainly possible to build a model using all the data.","metadata":{"_uuid":"fb969a26e52931bcaced3cbb7a36d8d8b1b04556"}},{"cell_type":"code","source":"train_df =  pd.read_csv('../input/train.csv', nrows = 10_000_000)\ntrain_df.dtypes","metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","execution":{"iopub.status.busy":"2022-08-11T23:06:11.543191Z","iopub.execute_input":"2022-08-11T23:06:11.543444Z","iopub.status.idle":"2022-08-11T23:06:29.150490Z","shell.execute_reply.started":"2022-08-11T23:06:11.543402Z","shell.execute_reply":"2022-08-11T23:06:29.149878Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"len(train_df.columns)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:11:35.906705Z","iopub.execute_input":"2022-08-11T23:11:35.907019Z","iopub.status.idle":"2022-08-11T23:11:35.911626Z","shell.execute_reply.started":"2022-08-11T23:11:35.906959Z","shell.execute_reply":"2022-08-11T23:11:35.911001Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's create two new features in our training set representing the \"travel vector\" between the start and end points of the taxi ride, in both longitude and latitude coordinates.  We'll take the absolute value since we're only interested in distance traveled. Use a helper function since we'll want to do the same thing for the test set later.","metadata":{"_uuid":"25df18156ed90f583efdbbc028c58a9d2bdfdc7b"}},{"cell_type":"code","source":"# Given a dataframe, add two new features 'abs_diff_longitude' and\n# 'abs_diff_latitude' reprensenting the \"Manhattan vector\" from\n# the pickup location to the dropoff location.\ndef add_travel_vector_features(df):\n    df['abs_diff_longitude'] = (df.dropoff_longitude - df.pickup_longitude).abs()\n    df['abs_diff_latitude'] = (df.dropoff_latitude - df.pickup_latitude).abs()\n\nadd_travel_vector_features(train_df)","metadata":{"_uuid":"59f0595db44dd60044cfd0404824651a7c2bee87","execution":{"iopub.status.busy":"2022-08-11T23:09:27.449268Z","iopub.execute_input":"2022-08-11T23:09:27.449539Z","iopub.status.idle":"2022-08-11T23:09:27.620158Z","shell.execute_reply.started":"2022-08-11T23:09:27.449491Z","shell.execute_reply":"2022-08-11T23:09:27.619545Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Explore and prune outliers\nFirst let's see if there are any `NaN`s in the dataset.","metadata":{"_uuid":"b1dbc7610bd467f1dfaf9042b5ec638eb2014aaf"}},{"cell_type":"code","source":"print(train_df.isnull())","metadata":{"_uuid":"e808c7e75338b45ca30f9f261dfbc90845700624","execution":{"iopub.status.busy":"2022-08-11T23:13:22.833484Z","iopub.execute_input":"2022-08-11T23:13:22.833967Z","iopub.status.idle":"2022-08-11T23:13:23.626004Z","shell.execute_reply.started":"2022-08-11T23:13:22.833903Z","shell.execute_reply":"2022-08-11T23:13:23.625357Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"False+False+True","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:14:39.397618Z","iopub.execute_input":"2022-08-11T23:14:39.397856Z","iopub.status.idle":"2022-08-11T23:14:39.402408Z","shell.execute_reply.started":"2022-08-11T23:14:39.397824Z","shell.execute_reply":"2022-08-11T23:14:39.401459Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:15:30.512542Z","iopub.execute_input":"2022-08-11T23:15:30.512792Z","iopub.status.idle":"2022-08-11T23:15:32.451957Z","shell.execute_reply.started":"2022-08-11T23:15:30.512744Z","shell.execute_reply":"2022-08-11T23:15:32.451330Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df['dropoff_longitude'].isnull()","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:17:20.743044Z","iopub.execute_input":"2022-08-11T23:17:20.743323Z","iopub.status.idle":"2022-08-11T23:17:20.760876Z","shell.execute_reply.started":"2022-08-11T23:17:20.743271Z","shell.execute_reply":"2022-08-11T23:17:20.760375Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"There are a small amount, so let's remove them from the dataset.","metadata":{"_uuid":"29bc86f2fa8baa37f0c4eb4300f77a8cb69f12aa"}},{"cell_type":"code","source":"train_df[train_df['dropoff_longitude'].isnull()]","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:18:11.402016Z","iopub.execute_input":"2022-08-11T23:18:11.402307Z","iopub.status.idle":"2022-08-11T23:18:11.479060Z","shell.execute_reply.started":"2022-08-11T23:18:11.402246Z","shell.execute_reply":"2022-08-11T23:18:11.478444Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('Old size: %d' % len(train_df))\ntrain_df = train_df.dropna(how = 'any', axis = 'rows')\nprint('New size: %d' % len(train_df))","metadata":{"_uuid":"9d8f28e24f3d4ca55ad93692329680774c341376","execution":{"iopub.status.busy":"2022-08-11T23:18:56.037214Z","iopub.execute_input":"2022-08-11T23:18:56.037494Z","iopub.status.idle":"2022-08-11T23:19:00.292062Z","shell.execute_reply.started":"2022-08-11T23:18:56.037450Z","shell.execute_reply":"2022-08-11T23:19:00.291393Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"help(train_df.dropna)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:19:32.893173Z","iopub.execute_input":"2022-08-11T23:19:32.893539Z","iopub.status.idle":"2022-08-11T23:19:32.898335Z","shell.execute_reply.started":"2022-08-11T23:19:32.893504Z","shell.execute_reply":"2022-08-11T23:19:32.897581Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now let's quickly plot a subset of our travel vector features to see its distribution.","metadata":{"_uuid":"6a045ef14c636ec726a5e8c349ca7e5fbb3a87c1"}},{"cell_type":"code","source":"plot = train_df.iloc[:2000].plot.scatter('abs_diff_longitude', 'abs_diff_latitude')","metadata":{"_uuid":"97d0aa1deab1c6cf0c97a4a3a12ba7007aada6c5","execution":{"iopub.status.busy":"2022-08-11T23:21:33.534166Z","iopub.execute_input":"2022-08-11T23:21:33.534666Z","iopub.status.idle":"2022-08-11T23:21:33.700885Z","shell.execute_reply.started":"2022-08-11T23:21:33.534626Z","shell.execute_reply":"2022-08-11T23:21:33.700277Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We expect most of these values to be very small (likely between 0 and 1) since it should all be differences between GPS coordinates within one city.  For reference, one degree of latitude is about 69 miles.  However, we can see the dataset has extreme values which do not make sense.  Let's remove those values from our training set. Based on the scatterplot, it looks like we can safely exclude values above 5 (though remember the scatterplot is only showing the first 2000 rows...)","metadata":{"_uuid":"22277d77f75e3177a5acaec9b820e0de6e869663"}},{"cell_type":"code","source":"print('Old size: %d' % len(train_df))\ntrain_df = train_df[(train_df['abs_diff_longitude'] < 5.0) & (train_df.abs_diff_latitude < 5.0)]\nprint('New size: %d' % len(train_df))","metadata":{"_uuid":"9703895e6c7e67b32c504f843b5ef19be2023964","execution":{"iopub.status.busy":"2022-08-11T23:24:27.928456Z","iopub.execute_input":"2022-08-11T23:24:27.928726Z","iopub.status.idle":"2022-08-11T23:24:28.761595Z","shell.execute_reply.started":"2022-08-11T23:24:27.928680Z","shell.execute_reply":"2022-08-11T23:24:28.760865Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Train our model\nOur model will take the form $X \\cdot w = y$ where $X$ is a matrix of input features, and $y$ is a column of the target variable, `fare_amount`, for each row. The weight column $w$ is what we will \"learn\".\n\nFirst let's setup our input matrix $X$ and target column $y$ from our training set.  The matrix $X$ should consist of the two GPS coordinate differences, plus a third term of 1 to allow the model to learn a constant bias term.  The column $y$ should consist of the target `fare_amount` values.","metadata":{"_uuid":"2151480a168d291bc2f4fd014fdac4ab7b5f6560"}},{"cell_type":"code","source":"len(train_df)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:26:34.132661Z","iopub.execute_input":"2022-08-11T23:26:34.132997Z","iopub.status.idle":"2022-08-11T23:26:34.137674Z","shell.execute_reply.started":"2022-08-11T23:26:34.132940Z","shell.execute_reply":"2022-08-11T23:26:34.137080Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Construct and return an Nx3 input matrix for our linear model\n# using the travel vector, plus a 1.0 for a constant bias term.\ndef get_input_matrix(df):\n    return np.column_stack((df.abs_diff_longitude, df.abs_diff_latitude, np.ones(len(df))))\n\ntrain_X = get_input_matrix(train_df)\ntrain_y = np.array(train_df['fare_amount'])\n\nprint(train_X.shape)\nprint(train_y.shape)","metadata":{"_uuid":"fb752441a1c1ce3e01d78452389ec48c95d52dc6","execution":{"iopub.status.busy":"2022-08-11T23:28:06.583114Z","iopub.execute_input":"2022-08-11T23:28:06.583610Z","iopub.status.idle":"2022-08-11T23:28:06.894890Z","shell.execute_reply.started":"2022-08-11T23:28:06.583573Z","shell.execute_reply":"2022-08-11T23:28:06.894334Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_X","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:28:39.213424Z","iopub.execute_input":"2022-08-11T23:28:39.213881Z","iopub.status.idle":"2022-08-11T23:28:39.218596Z","shell.execute_reply.started":"2022-08-11T23:28:39.213838Z","shell.execute_reply":"2022-08-11T23:28:39.217818Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now let's use `numpy`'s `lstsq` library function to find the optimal weight column $w$.","metadata":{"_uuid":"0dd0d8fe3f6478c050079df934381317a65f7d2b"}},{"cell_type":"code","source":"dir(np.linalg)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:47:45.549415Z","iopub.execute_input":"2022-08-11T23:47:45.549637Z","iopub.status.idle":"2022-08-11T23:47:45.553883Z","shell.execute_reply.started":"2022-08-11T23:47:45.549608Z","shell.execute_reply":"2022-08-11T23:47:45.553447Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"help(np.linalg.lstsq)","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:48:18.599852Z","iopub.execute_input":"2022-08-11T23:48:18.600124Z","iopub.status.idle":"2022-08-11T23:48:18.604202Z","shell.execute_reply.started":"2022-08-11T23:48:18.600093Z","shell.execute_reply":"2022-08-11T23:48:18.603614Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# The lstsq function returns several things, and we only care about the actual weight vector w.\n(w, _, _, _) = np.linalg.lstsq(train_X, train_y, rcond = None)\nprint(w)","metadata":{"_uuid":"85abbb09a27d2e1e2a15b261264b3c7cbdde39e4","execution":{"iopub.status.busy":"2022-08-11T23:53:23.063195Z","iopub.execute_input":"2022-08-11T23:53:23.063433Z","iopub.status.idle":"2022-08-11T23:53:24.085792Z","shell.execute_reply.started":"2022-08-11T23:53:23.063402Z","shell.execute_reply":"2022-08-11T23:53:24.085126Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"fare_amount = 180.02* (abs_diff_longitude) + 22.91*(abs_diff_latitide) + 6.45\n\nis the \"formula\" that the least squares method gives","metadata":{}},{"cell_type":"code","source":" np.linalg.lstsq(train_X, train_y, rcond = None)\n","metadata":{"execution":{"iopub.status.busy":"2022-08-11T23:56:30.242627Z","iopub.execute_input":"2022-08-11T23:56:30.243059Z","iopub.status.idle":"2022-08-11T23:56:30.960493Z","shell.execute_reply.started":"2022-08-11T23:56:30.243017Z","shell.execute_reply":"2022-08-11T23:56:30.959872Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"These weights pass a quick sanity check, since we'd expect the first two values -- the weights for the absolute longitude and latitude differences -- to be positive, as more distance should imply a higher fare, and we'd expect the bias term to loosely represent the cost of a very short ride.\n\nSidenote:  we can actually calculate the weight column $w$ directly using the [Ordinary Least Squares](https://en.wikipedia.org/wiki/Ordinary_least_squares) method:\n$w = (X^T \\cdot X)^{-1} \\cdot X^T \\cdot y$","metadata":{"_uuid":"4c11c9993467cd31c6be525f864eae24b0da364d"}},{"cell_type":"code","source":"train_df['abs_diff_longitude']","metadata":{"execution":{"iopub.status.busy":"2022-08-12T00:04:53.785313Z","iopub.execute_input":"2022-08-12T00:04:53.785671Z","iopub.status.idle":"2022-08-12T00:04:53.792737Z","shell.execute_reply.started":"2022-08-12T00:04:53.785629Z","shell.execute_reply":"2022-08-12T00:04:53.792198Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"w_OLS = np.matmul(np.matmul(np.linalg.inv(np.matmul(train_X.T, train_X)), train_X.T), train_y)\nprint(w_OLS)","metadata":{"_uuid":"4a629cdacdddd48a7ba9e8492b0e748cde819829","execution":{"iopub.status.busy":"2022-08-12T00:09:37.114432Z","iopub.execute_input":"2022-08-12T00:09:37.114883Z","iopub.status.idle":"2022-08-12T00:09:37.350354Z","shell.execute_reply.started":"2022-08-12T00:09:37.114838Z","shell.execute_reply":"2022-08-12T00:09:37.349708Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"w","metadata":{"execution":{"iopub.status.busy":"2022-08-12T00:08:20.144677Z","iopub.execute_input":"2022-08-12T00:08:20.144930Z","iopub.status.idle":"2022-08-12T00:08:20.149957Z","shell.execute_reply.started":"2022-08-12T00:08:20.144876Z","shell.execute_reply":"2022-08-12T00:08:20.149466Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Make predictions on the test set\nNow let's load up our test inputs and predict the `fare_amount`s for them using our learned weights!","metadata":{"_uuid":"a70ed21b43d720282bbae70e934b1188be2bc382"}},{"cell_type":"code","source":"test_df = pd.read_csv('../input/test.csv')\ntest_df.dtypes","metadata":{"_uuid":"3cbf4836cf8c71dfb67d13a9621b18a8d487197e","execution":{"iopub.status.busy":"2022-08-12T00:10:20.284679Z","iopub.execute_input":"2022-08-12T00:10:20.285015Z","iopub.status.idle":"2022-08-12T00:10:20.338895Z","shell.execute_reply.started":"2022-08-12T00:10:20.284956Z","shell.execute_reply":"2022-08-12T00:10:20.338264Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df.dtypes","metadata":{"execution":{"iopub.status.busy":"2022-08-12T00:13:02.445573Z","iopub.execute_input":"2022-08-12T00:13:02.445857Z","iopub.status.idle":"2022-08-12T00:13:02.451808Z","shell.execute_reply.started":"2022-08-12T00:13:02.445791Z","shell.execute_reply":"2022-08-12T00:13:02.451184Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-12T00:13:22.826136Z","iopub.execute_input":"2022-08-12T00:13:22.826411Z","iopub.status.idle":"2022-08-12T00:13:22.851146Z","shell.execute_reply.started":"2022-08-12T00:13:22.826360Z","shell.execute_reply":"2022-08-12T00:13:22.850568Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_df.head()","metadata":{"execution":{"iopub.status.busy":"2022-08-12T00:13:55.803643Z","iopub.execute_input":"2022-08-12T00:13:55.803906Z","iopub.status.idle":"2022-08-12T00:13:55.826714Z","shell.execute_reply.started":"2022-08-12T00:13:55.803864Z","shell.execute_reply":"2022-08-12T00:13:55.826240Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Reuse the above helper functions to add our features and generate the input matrix.\nadd_travel_vector_features(test_df)\ntest_X = get_input_matrix(test_df)\n# Predict fare_amount on the test set using our model (w) trained on the training set.\ntest_y_predictions = np.matmul(test_X, w).round(decimals = 2)\n\n# Write the predictions to a CSV file which we can submit to the competition.\nsubmission = pd.DataFrame(\n    {'key': test_df.key, 'fare_amount': test_y_predictions},\n    columns = ['key', 'fare_amount'])\nsubmission.to_csv('submission.csv', index = False)\n\nprint(os.listdir('.'))","metadata":{"_uuid":"ddba4a856ff617411a641dfdf7635e47f969dff8","execution":{"iopub.status.busy":"2022-08-12T00:13:25.268203Z","iopub.execute_input":"2022-08-12T00:13:25.268645Z","iopub.status.idle":"2022-08-12T00:13:25.487212Z","shell.execute_reply.started":"2022-08-12T00:13:25.268607Z","shell.execute_reply":"2022-08-12T00:13:25.485851Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Ideas for Improvement\nThe output here will score an RMSE of $5.74, but you can do better than that!  Here are some suggestions:\n\n* Use more columns from the input data.  Here we're only using the start/end GPS points from columns `[pickup|dropoff]_[latitude|longitude]`.  Try to see if the other columns -- `pickup_datetime` and `passenger_count` -- can help improve your results.\n* Use absolute location data rather than relative.  Here we're only looking at the difference between the start and end points, but maybe the actual values -- indicating where in NYC the taxi is traveling -- would be useful.\n* Use a non-linear model to capture more intricacies within the data.\n* Try to find more outliers to prune, or construct useful feature crosses.\n* Use the entire dataset -- here we're only using about 20% of the training data!","metadata":{"_uuid":"80ed89470e25d75c0b99008b9c88861be9739da3"}},{"cell_type":"markdown","source":"Special thanks to Dan Becker, Will Cukierski, and Julia Elliot for reviewing this Kernel and providing suggestions!","metadata":{"_uuid":"8fd559ff5ca72a73091d5dfd5b7032522832e999"}}]}