{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import gzip\n", "import matplotlib.pyplot as plt\n", "import numpy\n", "import random\n", "import scipy\n", "import tensorflow as tf\n", "from collections import defaultdict\n", "from fastFM import als\n", "from scipy.spatial import distance" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Data is available at http://cseweb.ucsd.edu/~jmcauley/pml/data/. Download and save to your own directory" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "dataDir = \"/home/jmcauley/pml_data/\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Factorization Machine (fastFM)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Parse the Goodreads comic book data (excluding review text)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "def parseData(fname):\n", " for l in gzip.open(fname):\n", " d = eval(l)\n", " del d['review_text'] # Discard the reviews to save memory\n", " d['year'] = int(d['date_added'][-4:]) # Use this for exercises\n", " yield d" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "data = list(parseData(dataDir + \"goodreads_reviews_comics_graphic.json.gz\"))" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "random.shuffle(data)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For example..." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'book_id': '15799191',\n", " 'date_added': 'Sun Jun 30 06:12:24 -0700 2013',\n", " 'date_updated': 'Thu Jul 04 08:02:07 -0700 2013',\n", " 'n_comments': 0,\n", " 'n_votes': 0,\n", " 'rating': 4,\n", " 'read_at': 'Thu Jul 04 08:02:07 -0700 2013',\n", " 'review_id': '86b43a448463928ac5d7887b364e2fcd',\n", " 'started_at': 'Sun Jun 30 00:00:00 -0700 2013',\n", " 'user_id': '23852e2647c217deb24964aadb26be64',\n", " 'year': 2013}" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "data[0]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Utility data structures. Most importantly, each user and item is mapped to an ID from 1 to nUsers/nItems" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "userIDs,itemIDs = {},{}\n", "\n", "for d in data:\n", " u,i = d['user_id'],d['book_id']\n", " if not u in userIDs: userIDs[u] = len(userIDs)\n", " if not i in itemIDs: itemIDs[i] = len(itemIDs)\n", "\n", "nUsers,nItems = len(userIDs),len(itemIDs)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(59347, 89311)" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nUsers,nItems" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Build the factorization machine design matrix. Note that each instance is a row, and the columns encode both users and items. Other features could straightforwardly be added." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "X = scipy.sparse.lil_matrix((len(data), nUsers + nItems))" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "for i in range(len(data)):\n", " user = userIDs[data[i]['user_id']]\n", " item = itemIDs[data[i]['book_id']]\n", " X[i,user] = 1 # One-hot encoding of user\n", " X[i,nUsers + item] = 1 # One-hot encoding of item" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Target (rating) to predict for each row" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "y = numpy.array([d['rating'] for d in data])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Initialize the factorization machine" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "fm = als.FMRegression(n_iter=1000, init_stdev=0.1, rank=5, l2_reg_w=0.1, l2_reg_V=0.5)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Split data into train and test portions" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "X_train,y_train = X[:400000],y[:400000]\n", "X_test,y_test = X[400000:],y[400000:]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Train the model" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "FMRegression(init_stdev=0.1, l2_reg=0, l2_reg_V=0.5, l2_reg_w=0.1, n_iter=1000,\n", " random_state=123, rank=5)" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fm.fit(X_train, y_train)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Extract predictions on the test set" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "y_pred = fm.predict(X_test)" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([2.23866277, 4.37847526, 3.84866594, 5.03002627, 3.94993096,\n", " 4.36740166, 4.22497778, 4.30029268, 3.28282377, 4.05036905])" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "y_pred[:10]" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([5, 4, 2, 5, 4, 5, 5, 5, 5, 5])" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "y_test[:10]" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [], "source": [ "def MSE(predictions, labels):\n", " differences = [(x-y)**2 for x,y in zip(predictions,labels)]\n", " return sum(differences) / len(differences)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.5940755947213534" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "MSE(y_pred, y_test)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Exercises" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.1" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Simple example, just incorporating a one-hot encoding of the year (see data extraction in examples above)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [], "source": [ "minYear = min([d['year'] for d in data])\n", "maxYear = max([d['year'] for d in data])\n", "nYears = maxYear - minYear + 1" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(2005, 2017, 13)" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "minYear, maxYear, nYears" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [], "source": [ "userIDs,itemIDs = {},{}\n", "\n", "for d in data:\n", " u,i = d['user_id'],d['book_id']\n", " if not u in userIDs: userIDs[u] = len(userIDs)\n", " if not i in itemIDs: itemIDs[i] = len(itemIDs)\n", "\n", "nUsers,nItems = len(userIDs),len(itemIDs)" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [], "source": [ "X = scipy.sparse.lil_matrix((len(data), nUsers + nItems + nYears))" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [], "source": [ "for i in range(len(data)):\n", " user = userIDs[data[i]['user_id']]\n", " item = itemIDs[data[i]['book_id']]\n", " year = data[i]['year'] - minYear\n", " X[i,user] = 1 # One-hot encoding of user\n", " X[i,nUsers + item] = 1 # One-hot encoding of item\n", " X[i,nUsers + nItems + year] = 1 # One-hot encoding of year" ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [], "source": [ "y = numpy.array([d['rating'] for d in data])" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [], "source": [ "fm = als.FMRegression(n_iter=1000, init_stdev=0.1, rank=5, l2_reg_w=0.1, l2_reg_V=0.5)" ] }, { "cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [], "source": [ "X_train,y_train,data_train = X[:400000],y[:400000],data[:400000]\n", "X_test,y_test,data_test = X[400000:],y[400000:],data[400000:]" ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "FMRegression(init_stdev=0.1, l2_reg=0, l2_reg_V=0.5, l2_reg_w=0.1, n_iter=1000,\n", " random_state=123, rank=5)" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fm.fit(X_train, y_train)" ] }, { "cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [], "source": [ "y_pred_with_features = fm.predict(X_test)" ] }, { "cell_type": "code", "execution_count": 30, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.6068283763129323" ] }, "execution_count": 30, "metadata": {}, "output_type": "execute_result" } ], "source": [ "MSE(y_pred_with_features, y_test)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Cold start plots. Count training instances per item (could also measure coldness per user if we had user features)." ] }, { "cell_type": "code", "execution_count": 31, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'book_id': '31423340',\n", " 'date_added': 'Mon Sep 18 19:24:36 -0700 2017',\n", " 'date_updated': 'Thu Sep 21 04:54:01 -0700 2017',\n", " 'n_comments': 6,\n", " 'n_votes': 13,\n", " 'rating': 4,\n", " 'read_at': 'Thu Sep 21 15:03:21 -0700 2017',\n", " 'review_id': '4f39bbb5aeab832c794a3d2e40943949',\n", " 'started_at': 'Tue Sep 19 20:23:13 -0700 2017',\n", " 'user_id': '1d945500234cbc7a6138a4d017dbfe4b',\n", " 'year': 2017}" ] }, "execution_count": 31, "metadata": {}, "output_type": "execute_result" } ], "source": [ "d" ] }, { "cell_type": "code", "execution_count": 32, "metadata": {}, "outputs": [], "source": [ "nTrainPerItem = defaultdict(int)\n", "for d in data_train:\n", " nTrainPerItem[d['book_id']] += 1" ] }, { "cell_type": "code", "execution_count": 33, "metadata": {}, "outputs": [], "source": [ "rmsePerNtrainFeatures = defaultdict(list)\n", "rmsePerNtrain = defaultdict(list)\n", "for d,y,ypf,yp in zip(data_test,y_test,y_pred_with_features,y_pred):\n", " e2_features = (y-ypf)**2\n", " e2 = (y-yp)**2\n", " nt = nTrainPerItem[d['book_id']]\n", " if nt < 5:\n", " rmsePerNtrainFeatures[str(nt)].append(e2_features)\n", " rmsePerNtrain[str(nt)].append(e2)\n", " elif nt < 10:\n", " rmsePerNtrainFeatures['5-9'].append(e2_features)\n", " rmsePerNtrain['5-9'].append(e2)\n", " elif nt < 20:\n", " rmsePerNtrainFeatures['10-19'].append(e2_features)\n", " rmsePerNtrain['10-19'].append(e2)\n", " else:\n", " rmsePerNtrainFeatures['20+'].append(e2_features)\n", " rmsePerNtrain['20+'].append(e2)" ] }, { "cell_type": "code", "execution_count": 34, "metadata": {}, "outputs": [], "source": [ "for r in rmsePerNtrain:\n", " rmsePerNtrainFeatures[r] = sum(rmsePerNtrainFeatures[r]) / len(rmsePerNtrainFeatures[r])\n", " rmsePerNtrain[r] = sum(rmsePerNtrain[r]) / len(rmsePerNtrain[r])" ] }, { "cell_type": "code", "execution_count": 35, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "defaultdict(list,\n", " {'0': 1.1034747614259885,\n", " '1': 1.286213343267511,\n", " '10-19': 1.8579978486101703,\n", " '2': 1.3303440529731148,\n", " '20+': 1.696268806936738,\n", " '3': 1.4912379121450905,\n", " '4': 1.4174955135661096,\n", " '5-9': 1.670982336713272})" ] }, "execution_count": 35, "metadata": {}, "output_type": "execute_result" } ], "source": [ "rmsePerNtrain" ] }, { "cell_type": "code", "execution_count": 36, "metadata": {}, "outputs": [], "source": [ "Xlab = ['0','1','2','3','4','5-9','10-19','20+']\n", "X = [0,1,2,3,4,5.5,7,8.5] # Average of above ranges\n", "YF = [rmsePerNtrainFeatures[x] for x in Xlab]\n", "Y = [rmsePerNtrain[x] for x in Xlab]" ] }, { "cell_type": "code", "execution_count": 37, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEZCAYAAACJjGL9AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd4VFX6wPHvG3oJJTSlhSagCCIgQmgRFVSkl4AIIogV1l1XV6zguurasf0UWESRhYSiQAQiCEzoooCIIiBIURCkhQ4p8/7+mEk2kEkj05K8n+eZh7lz7z3nvcNk3rn33HOOqCrGGGNMZkICHYAxxpjgZonCGGNMlixRGGOMyZIlCmOMMVmyRGGMMSZLliiMKeREpKGIPBjoOEzwskRRSIjIv0TksIgcCHQsJiMReUpEJgao+l+A3iLSw9cVicg9IrLS1/UY77JEEaREZI+InBWRkyLyh4hMEZHSl1lWLeAxoLGqVvdupMYbVPUVVb0fQETCRcQpIn75+1RXZ6rBwDMiUt4fVfqhjgxEpJOILA9E3fmdJYrgpUA3VS0HtABaAc/mthARKQKEA0dU9ehl7l+oiYj4u0pc//9+q1dVj6jqjap6wl91Boj1ML4MliiCmwCo6h/AIuBaABEpJyL/EZEDIvKbiLyY+mXmPrVfJSJvicgRYDmwGKjhPjv52L1dDxH5UUSOicgyEWmcVqnIbhH5h4hsBk6LSBH3a4+LyGYROSUik0SkqogsdJe7OP2vURGZ6T4TOi4iDhG5Jt26KSLyvoh86d53rYjUTbe+ibu8o+4yxrhfFxEZIyI73ZfRokWkgsc3TmSriNyRbrmIiPwpIs3dy21EZLU7vk0i0indtsvdl+pWicgZoK6IDBORXe54d4nIIPe2Y0Xks3T7XnQ2kNl+HuIdKyJT3Yvx7n8T3Pvd6N5muPu4jorIIhGpnW5/p4g8JCI7ROSEiPxTROq5jzHB/V4V9VR3VkRkpLvOk+7PS+r719j9Ph0XkS0i0j3dPuVEZKr7/d4tIs9kUb5TRB5wx31MRN6/ZH1Wx/y2iBxyH+/m1M+YiNwhIj+5Y/5NRB7L7XGbS6iqPYLwAewGOruf1wJ+BMa5l78A/g8oCVQG1gEj3evuAZKAh3H9ECgBdAL2pSu7IXAa6AwUAZ7AdZ26aLq6NwLVgRLpXlvjru9K4BDwHdAMKA4sBZ5LV8cwoDRQDHgL2JRu3RTgMNDSHeM0YLp7XVngAPBXd7llgBvc6x51x3Clu9wPU/fz8P49C0xLt9wN+Mn9vAZwBOjqXr7ZvVzJvbwc2AM0dsdXDjgBNHCvrwZc7X4+Fpiarp5wIMW9X+nM9vMQb1o56cqQdOt7Ajvc/3chwNPA6nTrne7PRRngauA8sMRdVijwEzAkl5/B/sBvQAv3cj1cn8Wi7s/Lk+7nNwEngavc2011x1LaXf924N50n88Vl8Q93x1jLeBPoEt2xwx0Ab4FQt3LjYBq7ucHgAj38/JA80D/Pef3R8ADsEcm/zGuL+aTwDH38/dwfelXdX8JlEi37UBgmfv5PcCeS8q6NFE8C0SnWxbgd6Bjurrv8RDPoHTLs4EP0i2PAj7P5FgquL8QUv+opwAT062/Hdjqfj4I2JBJOVuBm9ItXwkkAiEetq3vfv9KupenAc+6n/8D+PSS7eNSv0hxJYpx6daVdv8/9E4tL9267BKFx/08xOspUYSkW78Q95etezkEOAPUci87gTbp1n8HPJFu+Q3grVx+BuOA0R5ebw8cuOS16cDz7rguAI3Srbv/ks/npYmibbrlGOAf2R0zruS0DbiRdAnVvd0eYGTq580eeX/Ypafg1lNVw1S1rqqOVtULuL5EigF/uE/VjwMf4fqln+q3bMqtDuxNXVDXX9dvuH5pp/rdw36H0j0/52G5LICIhIjIv92XiBJwJRm9JMaD6Z6fTd0XqAnsyiTucOAL93Efw5U4knD9Ur+Iqu5yr+8uIqWAHsB/05UzILUc93vYDrgiXRG/pSvrLBAFPITrfY8VkYaZxJg+Bk/7Ncpuv0yEA++kO/ajuN7T9P9nf6Z7nun/Ty7UwvP/RXUyfsb2umOpjOvzuc/DusykjzP9ZyHTY1bV5cD7wAfAIRH5SERS9+uL6wxyr/vyWJusD9NkxxJFcPPUmPkbrjOKSu4kUlFVK6hqs3TbZNdgdwDXH2F6tbg4OeSl0W8w0B3XpbMKQB1cx5KTxtnfcJ0NeLIPuN193KnHXkZdbTieRAN34bqE8ZOq7k5Xx9RLyglV1dfT7XvR8avqElXtgiuZbAcmuVedwXXmkOrKHO6XFU/v/T7ggUtiLquq63JQ3uXK7P/iAK7PS3q1gf24LuElcfHnK9y97nLqz/SYVfV9VW0FXIPr0tMT7tc3qGovoAowD5h5GXWbdCxR5DOqehBX4/TbIhLqbuCtJyIdc1HMTKCbiNwkIkVF5HFcyWetl8Isi+vyw3ERKQO8Qs4Tz5fAFSLyFxEpLiJlRaS1e90E4OXUBk0RqSJZ3/sfjeta9kO4Lo2kmobrTKOL++ynpLhunfR467C4Gu17iOv25CRc7TtO9+rvgY4iUktcjfljstkvJQfvwWF3+em/pCcAT6drsC0vIv1yUFZe/Ad4XERauOusL65brb8BzorrhoeiIhIJ3AnMUFUnrstHL7n/78KBvwGfea4iSx+RyTGLSCsRae1uoD+H6/PrFJFiInKXiJRT1RTgFDl7z00WLFEEr6y+WIfiaujdiusa+CwuvmySdcGqO4C7cZ26H8Z1mt5dVZOzqPvS17KKbyquX8D7cTXCr8lFbKeBW3FdKjqIqzEz0r36HVy/EBeLyAl3ua09FJNa1kFcya8Nri+v1Nd/x3WW8TSu498LPM7//h4uPbYQXP1QUn8xd8SVfFDVr91l/4CrcTU2J/tl8x6cA14CVrsvu7RW1bnAv4Fo9+W8H4Db0u92aTHZ1ZODOGa745guIidxNVCHqWoSrjPGO3Ad1/u42nd+ce/6F1yXkH4FVuC6qWBKZtVktpzNMZfDdXaW2oZ3BEg9IxwC7Hbvcz+us0qTB+K6PO2jwkUm4/qlceiSSyOp6ysAH+P65XQOGK6qW30WkDHGmFzz9RnFFKBrFuufxnXb5HW47oZ418fxGGOMySWfJgpVXQUcz2KTa4Bl7m23A3VEpIovYzLGGJM7gW6j2Az0AXA3WNbGdXukMcaYIBHoRPFvoKKIbAQeATZhdygYY0xQyfXYL96kqqeA4anLIrIb150SGYiI71rdjTGmAFPVPA0w6Y8zikw7Wrnviy7mfj4SiHffHulRoLuxZ/UYO3ZswGOw2ApffHasdqzZPbzBp2cUIjId1z3wlURkH67xbIrjGjViIq7Byz4VESeuQctG+DIeY4wxuefTRKGqWXZ0UVdX/Msd+8YYY4wfBLoxu8CIjIwMdAiZstguX7DH5012rCYzPu2Z7U0iovklVmNM8Fq1ahXly5enadOmgQ7FL0QEzQeN2T5Vp04dRMQe9sj1o06dOoH++Bo/++STT+jbty+dO3fms88uZ5zCwimgt8d6w969e73Wsm8KFxF/T4VtAkVVefnll/nPf/7DihUrSExMpE+fPnzzzTe89dZbFC9ePNAhBrV8f0ZhjDFZSUlJYdSoUcyePZs1a9bQqFEjmjZtyrfffsu+ffuIjIxk//7LmS6j8LBEYYwpsM6dO0e/fv3Yvn078fHxXHnl/+aVqlChAnPnzqVbt27ccMMNrFixIoCRBjdLFMaYAunYsWPccsstlC5dmoULF1KuXLkM24SEhPDMM88wZcoU+vfvz9tvv22Xsj2wRBEkQkND2bNnT6br69aty7Jly3Jc3rPPPkuVKlWoXt3jpG3GFGh79+6lXbt2tGvXjs8++yzbNoiuXbuybt06PvvsMwYNGsTp05kOEFEoWaIIEqdOnUq7C+fee+/l+eefv+yyfvvtN9566y22bdvGgQMH8hRXfHw8tWpdOj2yMcFr8+bNtGvXjgcffJDXXnuNkBDX11xSUlKGbZOTk9Oe161bl9WrV1OqVCnatGnDL7/8kmH7wsoSRQG0d+9eKleuTKVKlfJclqrm6e6glBQbDNj4z7Jly7j11lt56623ePTRR9Ne//HHH5k4cWKGz+OXX37JwoUL0xJGqVKl+Pjjjxk1ahTt2rVj3rx5fo0/aAV6wKpcDGylnmT2ejCYMmWKdu/ePW25QYMGOmDAgLTlWrVq6ebNm1VVVUR0165dOnHiRC1WrJiWKFFCQ0NDtUePHqqqWqdOHX3jjTe0WbNmWqFCBR04cKBeuHAhQ51ff/21lipVSosUKaKhoaF67733qqrq2rVrNSIiQitUqKDNmzdXh8NxUZxXX321hoaGav369XXChAmqqnrmzJm0ssqWLauhoaH6xx9/6LBhw/S5555L29/hcGjNmjXTluvUqaOvvvqqNmvWTEuWLKkpKSl64MAB7du3r1apUkXr1aun7777btr269ev11atWmm5cuX0iiuu0L///e95et9zKpg/Oyb3ZsyYoVWrVtXly5enveZ0OnX16tX61ltv6cGDBzPsc+7cOY2OjtaJEyfq8ePHL1q3du1arVmzpj799NOanJzs6/B9xv05z9v3b14L8NcjPyaKX3/9VStWrKiqqgcOHNDw8HCtVauWqqru2rVLw8LC0rYNCQnRXbt2qapm+CJWdX353njjjXrw4EE9fvy4Xn311Wlf6JdyOBxp9aiq7t+/XytVqqRxcXGq6komlSpV0iNHjqiq6sKFC3X37t2qqrpixQotXbq0btq0yWNZnuK7dJs6dero9ddfr/v379fz58+r0+nUli1b6r/+9S9NTk7W3bt3a/369XXx4sWqqtq2bVudNm2aqrqS0zfffJPte+sNwfzZMbnz5ptvaq1atfSHH35Iey0lJUUXLlyo//d//6cnTpzIdN/UZPL666/rjh07Llp36NAhjYyM1C5duqT9veQ33kgUheLSkzd68V6OunXrEhoayvfff8+KFSvo2rUr1atXZ8eOHaxYsYIOHTqkbas5uNPi0UcfpVq1alSoUIHu3bvz/fff5yiOadOm0a1bN7p2dU1ffvPNN9OqVSsWLlwIwO23357WPtKhQwe6dOnCypUrc3m0GWOtXr06JUqU4Ntvv+XIkSM888wzFClShDp16nDfffcRHR0NQLFixdi5cydHjx6ldOnStG7dOk91m8LD6XTy2GOPMXnyZFavXp02LIeqMmfOHP7880/uvfdej3c8pRIRIiIiGDBgAF9++SW//vq/KXGqVq3KkiVLaNasGa1atWLDhg0+P6ZglO97ZudETr6EfaVTp04sX76cnTt3EhkZScWKFXE4HKxdu5ZOnTrlqqxq1aqlPS9dujR//PFHjvbbu3cvM2fOJDY2FnC9H8nJyXTu3BmARYsW8c9//pMdO3bgdDo5d+4czZo1y1Vsl6pZ838z2u7du5f9+/cTFhaWVr/T6aRjx44AfPzxxzz33HM0btyYevXq8fzzz9OtW7c81W8KvgsXLjBs2DD279/PqlWrqFixYto6EaFFixaEh4dTtGjOvuZq167NAw88QMmSJS96vWjRorz++uvceOON3Hbbbbz66qsMHz48k1IKpkKRKAKpY8eOxMbGsmfPHp555hnKly/Pf//7X9atW8fo0aM97uPtoSVq1arF0KFDmTBhQoZ1iYmJ9OvXj2nTptGzZ09CQkLo3bt3WnL1FEuZMmU4e/Zs2rKnhJV+v1q1alGvXj22b9/uMb769eszffp0AObMmUO/fv04duwYpUqVyt2BmkLjxIkT9O7dm7CwMBYvXpzhyx1cn6vcKl26dKbr+vXrR5MmTejduzfffPMN7777LiVKlMh1HflRobj0FEipZxTnzp2jevXqdOjQgbi4OI4ePcr111/vcZ9q1apddPqbV3fffTexsbEsXrwYp9PJ+fPniY+P58CBAyQmJpKYmEjlypUJCQlh0aJFLF68+KJYjh49ysmTJ9Nea968OQsXLuT48eMcPHiQd955J8v6W7duTWhoKK+99hrnz58nJSWFn376ie+++w6A//73vxw5cgSA8uXLIyJptzQac6kDBw7QsWNHrrnmGmJiYjwmCV+5+uqrWb9+PYcPH6Zjx4789ttvfqs7kOyv0ceuuuoqQkND0y6zhIaGUr9+fdq3b3/Rr+70z0eMGMFPP/1EWFgYffr0ybA+t2rWrMm8efN4+eWXqVKlCuHh4bzxxhs4nU7Kli3Lu+++S//+/QkLCyM6OpqePXum7duoUSMGDRpEvXr1CAsL4+DBgwwZMoRmzZpRp04dbrvtNgYOHHhRfZfGGhISwpdffsn3339P3bp1qVq1KiNHjkxLPnFxcTRp0oRy5crxt7/9jZiYmELzS83kzs8//0xERAR33XUX7733HkWKFAHg/PnzPqvz1KlTzJ07l3PnzgFQrlw55syZQ58+fWjdunWuOsLmV/l+PgoRCWgbhMm/7LOTv6xatYq+ffvyxhtvMGTIkLTXN27cyJo1a3jooYfSEoc3paSksGTJErZv307//v0vGu1g6dKlDB48mL///e88/vjjQTkisftznqfALFGYQss+O/nHF198wQMPPMC0adPo0qUL4LopwuFwsGXLFgYPHuyVDqZZ2bp1KwsWLCAyMpJWrVqlJYV9+/bRt29fwsPDmTJlCqGhoT6NI7e8kSh8eulJRCaLyCER+SGT9eVEZL6IfC8iW0RkmC/jMcbkPx9++CGjRo0iLi4uLUmkpKQwf/58du7cyfDhw32eJACuueYahg8fzoYNG5g7d27aj4zatWuzcuVKwsLCaN26Ndu2bfN5LP7m0zMKEWkPnAamqmqG+y1F5CmgnKo+JSKVge1ANVVN9rCtnVEYr7LPTnBTVZ599llmzZpFXFwc9erVS3s9OjoaVaVfv35+n3QoKSmJ3bt307BhwwzrJk+ezJgxY/joo4/o27evX+PKjDfOKHx6e6yqrhKR8Kw2AVLP00KBo56ShDGmcElKSmLkyJH8/PPPrF69mipVqqStExHat29PjRo1AnJ3XLFixTwmCXDdiHLdddfRt29fvvnmG15++eUc9+MIZoG+6+l94BoROQBsBh7NZntjTAF3+vRpunfvztGjR1m2bNlFSSJVrVq1gvYW6tQe3Js2baJr164cPnw40CHlWaBTXVdgk6p2FpH6wBIRaaaqHgeDHzduXNrzyMhIIiMj/RKkMcY/Dh06RLdu3bj++uv58MMP89Wv8X379hEaGkrFihWpXLkycXFxPPfcc7Rs2ZLZs2f7bWgah8OBw+HwbqF5HSwquwcQDvyQybovgXbplpcCrTLZNqsBr4zJNfvsBJdffvlF69evr2PHjlWn05n2+qlTpwIYVc5t2LBBX3vtNd22bdtFr3/++edapUoVnTBhwkXH5S/kk0EBxf3wZC9wC4CIVAMaAt7rkmyMyRfWr19Phw4dePLJJxk3blzaradr165lypQp+WJekxYtWjBo0CAWLVrEkiVLcDqdAPTu3ZuVK1fyzjvvcN9996V13MtPfH177HRgDdBQRPaJyL0i8oCI3O/e5F9AhPv22SXAP1T1mC9jClbengrVm/788086duxI+fLleeKJJwISgym4Fi5cyJ133snEiRMZOXIk4LrSERcXx6ZNmxg6dKhPOtL5Qs2aNbn//vs5dOgQn376KadOnQJcIxx88803nDp1ivbt22f5tx6U8npK4q8HhejSU2bzUSxdutTrdX3yySfavn37LLd58cUXtW/fvl6pz9OxBUpB/OzkN5MnT9Zq1arp2rVr015LSkrSmTNn6pQpU/Ts2bMBjO7yOZ1OjY+P119//TXD62+++aZWq1ZNv/rqK7/EghcuPeWfliLjE6rZT3W6d+9errnmGj9FlLWUlJR88+vSZE5Veemll/j444+Jj4+nUaNGaa9HR0dTsmRJ7r777nzVmJ2eiKSN73bp64899hgtWrTgrrvuYtSoUYwZMyZo7+BKk9dM468H+fCMwp9ToU6cOFEbNGiglSpV0p49e+qBAwdUVXXPnj0qIpqSkpK2bWRkpE6ePFl//vlnLVmypBYtWlTLli2bNhtfesOGDdNixYpp8eLFNTQ0VJcuXapOp1NfeeUVrV+/vlauXFmjoqL02LFjafv0799fr7jiCq1QoYJ26tRJt27dmhajp2NLPfb0daaedaROs/rqq6/qFVdcoUOHDlVV1djYWG3evLlWqFBB27Vrd9HMZv/+97+1Ro0aGhoaqo0bN9Zly5Z5/P8J5s9OQZacnKwPPvigNm/ePO1zmt7BgwcD0ujrb7/99pu2adNGe/bsqQkJCT6rB5sKNbj/2P01FerSpUu1cuXK+v3332tiYqKOHj1aO3bsqKquRBESEuIxUai6Lj116NAhy+O4NJ7x48dr27Zt9cCBA5qYmKgPPvigDho0KG39lClT9MyZM5qYmKh/+9vftHnz5pmWdemxX7qNw+HQokWL6lNPPaWJiYl6/vx53bhxo1atWlW//fZbdTqdOnXqVK1Tp44mJibq9u3btVatWmnzI+/duzfD6X+qYP7sFFRnzpzRnj176i233JLl9KQF2cmTJ9Oenz9/Xh9++GG96qqrdMuWLT6pzxuJIn+e1+WSw+EgPj4+w+udOnXy2Bfj0u0z2y476adC3b59O127dmXz5s3s2LGDNWvWXPZUqMBFU6FOnz49rUcowCuvvELFihXZt29frmPOiQkTJvDBBx9w5ZVXAvD8888THh7OtGnTCAkJYdiwYWnbPv/884wfP55Tp05lOlhadsdepEgRXnjhBYoVKwbApEmTePDBB2nVqhUAQ4YM4aWXXmLdunVUr16dxMREfvzxRypVqkTt2rW9cMTGG44ePUr37t2pX78+M2fO9PvQG8Hg3LlzTJgwgQ4dOtC6dWtKlCjBBx98wNSpU7npppt4//33iYqKCnSYGRSKRJHbznne7Mznj6lQDxw4QMuWLdPWlSlThkqVKrF///6LhkT2lr1799K7d++066qqSrFixTh06BDVqlXj6aefZvbs2Rw5ciRtzvEjR45c9qiaVapUSUsSqfVPnTqV9957L63+pKQkDhw4QIcOHRg/fjzjxo1j69atdO3alTfffDMtqZnA2LNnD7fddhu9evXi5ZdfTvvsJCQkUKFChQBH5z+lSpVixIgRzJw5k3379tGjRw9KlCjB0KFDadq0adrQH6+++upFn/lAC/IWlPyvY8eOOBwOVq1aRadOnejYsSPx8fGsWLEi00SR2zHtq1evzt69e9OWz5w5w9GjR6lZsyZlypQBuGjq0oMHD152XeAaLXPRokUcO3aMY8eOcfz4cc6cOcOVV17J9OnTiY2NZdmyZSQkJLBnz570lw891le6dOlM4/O0T61atXjmmWcuqv/06dNpv8QGDhzIypUr096TMWPG5PoYjfd8//33tG/fnkceeYR///vfhISEoKrEx8czbdq0fNFHwpsqVqzIiBEjKFmyJJMmTeLQoUMAXH/99Xz33Xf8/PPP3HLLLRn+DgLJEoWP+WMq1EGDBjFlyhR++OEHLly4wNNPP02bNm2oVasWlStXpkaNGkybNg2n08nHH3/Mrl27Lqrr999/JykpKcf1PfDAAzz99NNpl7YOHz7M/PnzAddsYCVKlKBixYqcOXOGp5566qIvek/Hdv311zN9+nScTidxcXEeLxOmN3LkSD766CPWr18PuBLjwoULOXPmDDt27GD58uUkJiZSvHhxSpUqFfx3lBRgS5cupUuXLrz99ttpc8Q7nU5iY2PZvn07w4YNK5R3sRUtWpTu3bvToUMHPv/887TOeWFhYXz55Zd06tSJVq1asWbNmgBH6mJ/QT7mj6lQb775Zl588UX69OlDjRo12L17N9HR0WnrJ02axGuvvUblypX5+eefadeuXdq6zp0706RJE6644gqqVq3qsfxL63700Ufp2bMnXbp0oXz58kRERKR9aQ8dOpTatWtTo0YNrr32WiIiIi7a19OxjR8/nvnz51OxYkVmzJhB7969M39DgZYtWzJp0iRGjRpFWFgYDRs25NNPPwXgwoULjBkzhipVqlC9enUOHz7MK6+8kmV5xjemT5/OXXfdxaxZs+jfvz8AiYmJREdHc+rUKYYNG0bZsmUDHGVgXXfdddx///0X/ZgpUqQI//znP/noo4/o1asXH3zwQY7aMH3JZrgzhZZ9dnxDVXnzzTd59913WbhwIddee23a65999hkVKlSgW7duhfJMIrd27txJnz59aN68OR999BGlS5fOdRk2FSr2x24un312vM/pdPL3v/+dJUuWsGjRImrVqnXR+mPHjlGxYsWgnFs6mKj+ryPsmTNnuP/++/npp5/4/PPP0yZwyqmgnwrVGFN4XLhwgUGDBrFx40ZWrlyZIUmA6xq8JYnsrVmzhq+++oqUlBTKlCnDtGnTGDFiBG3btmXhwoV+j8cShTEmzxISErjttttwOp189dVXVKxYMdAh5WstWrTg6NGjfPrpp5w8eRIRYfTo0cyZM4eRI0fywgsvpDWA+4MlCmNMnuzfv5+OHTvStGnTtHGawDXqsLk8pUqVYtCgQTRs2JBJkyal3anYvn17vvvuO77++mt69OjB8ePH/RKPJQpjzGX76aefiIiIYMiQIbzzzjsUKVIEVWXJkiXMnj270PWR8KbUucH79u3LvHnz+OGHHwC48sorWbZsGfXr1+eGG25Ie92nseSXxjxrzDbeZp+dvFm1ahV9+/blzTff5O677wYgOTmZefPmceLECQYNGkSpUqUCHGXBcPr0aUQkrQNtqunTp/Poo4/y9ttvp/0fXMobjdkBH+wvpw8yGcAtPDxcAXvYI9eP8PBwj58pk705c+ZolSpVdPHixWmvnTt3TqdMmaIxMTGamJgYwOgKl82bN2uDBg101KhRF40onQryPihgvj+jMMb41wcffMDLL7/Ml19+mTa6gKry6aefcsUVV9ClSxfrDe9nCQkJDBkyhGPHjjFr1qyLxnizfhTGGL9RVZ555hnmzJlDXFwcdevWvWh9ViMEG+9zOp2sWLGCNm3aULJkSZxOJy+99BKY9/HaAAAgAElEQVQffvgh0dHRaaNBWD8KY4xfJCUlMWzYMJYtW8bq1aszJAnAkoSfqSpnz55l4sSJHDx4kJCQEJ577jk+/vhj+vfvz/jx473WBufTMwoRmQzcCRxS1WYe1j8ODMZ1zbgYcDVQWVUTPGxrZxTGBMCpU6fo168fxYsXJyYm5rKGkTC+s2XLFuLi4rj55pu5/vrrERF2795Nnz59uPrqq5kxY0bQn1FMAbpmtlJV31DV61W1BfAU4PCUJIwxgXHo0CFuuukmwsPD+eKLLyhdujSqyu+//x7o0Ixb06ZNGTZsGOvWrWPevHkkJydTt25d1qxZ47XJoXyaKFR1FZDTHiGDgBk+DMcYkwu//PILERER9OjRgwkTJlC0aFGcTicLFixgwYIF1kciiFSpUoX77ruPKlWqpN1IUKpUKaZMmeKV8n3emC0i4UCsp0tP6bYpBfwO1M/sjMIuPRnjP+vXr6dnz568+OKL3HfffYBriPA5c+aQkpJC//79KVGiRICjNDnhjcbsYJkKtTuwKrvLTuPGjUt77s3pSo0x/7NgwQKGDRvGlClTuPPOOwHXCKYzZsygSpUq3HnnnTZEeBBzOBw4HA6vlhksZxSfAzNVNTqLbeyMwhgfmzx5Ms888wzz5s3jxhtvBP7XRyI8PJzIyEgb/TWfyRf9KESkDq5E0TST9eWBX4Gaqnoui3IsURjjI6rKiy++yCeffEJcXBwNGza8aP25c+dsOI58KugvPYnIdCASqCQi+4CxQHFcXconujfrBXyVVZIwxvhOUlISDz/8MBs3bmTNmjVcccUVGbaxJFG4Wc9sYwqxo0eP0q9fP8qWLcv06dOt01wBZD2zjTGX7eeff+bGG2/khhtuYO7cuYSGhqKqaXMfGJMqWO56Msb40VdffcWQIUN47bXXGDZsGOAaInz+/PkcP36c8PBwiha1rwfjYpeejClEVJX333+fl19+mVmzZtG+fXsAzp8/T0xMDCVLlqRPnz4UK1YswJEabwn6xmxjTPBISkpi9OjRrF69mjVr1qQN7JeUlMS0adOoXr06t912mw0RbjKwRGFMIXDs2DH69+9PqVKlWL16NeXKlUtbt2TJEsLCwrj99tutj4TxyC49GVPAbdu2je7du9OrVy/+/e9/Z+hVferUKUqVKmVtEgWUX+56EpHXRKSciBQTkaUiclhEPE/OaowJKosXL6Zjx4489dRTvP766x6H3ggNDbUkYbKUk4uRXVT1JK55JfYADYAnfBmUMSZvUhuthw4dyuzZsxk+fHigQzL5WE5+RqRu0w2Ypaon7DqmMcErKSmJv/zlL6xYsYI1a9ZQr169QIdk8rmcnFF8KSLbgJbAUhGpApz3bVjGmMtx7NgxbrvtNvbt28fatWszJImkpCS++eYbr02RaQqHbBOFqo4BIoBWqpoEnAF6+jowY0zubN++nTZt2tC8eXPmz59/0Z1N4LocNXfuXA4cOBCgCE1+le2lJxEpArQH6ohI+u3f8llUxphcWbJkCXfffTcvv/wyI0aM8LiNw+Hg1KlTDB061G6DNbmSkzaKWFyXmrYATt+GY4zJrQ8++IAXX3yRmTNn0qlTJ4/bbNmyhc2bNzNy5Ei7w8nkWk4+MTWzmnTIGBMYSUlJ/PWvf8XhcGTZaL1//37i4uIYOnQoZcqU8XOUpiDItsOdiLwKLFXVxf4JKdM4rMOdMW7Hjx+nf//+FC9enBkzZlC+fPlMtz1z5gyHDx+mTp06/gvQBA1/DTO+DvhCRM6JyEkROSUiJ/NSqTHm8u3YsYM2bdrQrFkzYmNjs0wSAGXKlLEkYfIkJ4niLaAtUFpVy6lqqKqWy24nY4z3LV26lA4dOvDEE0/w1ltveexpbYy35SRR/Ab8aNd9jAmsDz/8kMGDBxMTE8N9990X6HBMIZKTxuxfAYeILAIupL6oqnZ7rDF+kJyczF//+leWLl3KqlWraNCgQZbb79+/nyuvvNKGCzdek5NP0m5gKVAcCE33yJaITBaRQyLyQxbbRIrIJhH5UUSW56RcYwqL48ePc8cdd7Bz507WrVuXbZL47bffmD59OgkJCX6K0BQGOR5mXERKq+rZXBUu0h44DUz1dIutiJQH1uAaeHC/iFRW1SOZlGVXv0yh8ssvv3DnnXdy++2388Ybb2Tb/yEhIYHJkyfTvXt3GjZs6KcoTbDz1zDjbUVkK7DNvXydiPxfTgpX1VXA8Sw2uQuYo6r73dt7TBLGFDbLli2jffv2/P3vf2f8+PHZJonExESio6OJiIiwJGG8LieXnsYDXYGjAKq6GejopfobAmEislxEvhWRIV4q15h8a8KECdx1111ER0dz//33Z7u9qvL5559TvXp12rRp44cITWGTo778qvrbJWPDpHix/hZAZ6AMsFZE1qrqTk8bjxs3Lu15ZGQkkZGRXgrDmMBLTk7mscceY8mSJTlqtE6VlJRE5cqVuemmm2wMJ4PD4cDhcHi1zJz0zJ6Nqy/F+8CNwKO4RpIdmKMKRMKB2EzaKJ4ESqrqC+7l/wCLVHWOh22tjcIUWAkJCURFRQEQExNDhQoVAhyRKSj81TP7QeARoAawH2gOPJyLOsT98GQe0F5EiohIaVyJ6OdclG1Mvrdz507atm1Lo0aNWLBggSUJE3RycumpkaoOTv+CiLQDVme3o4hMByKBSiKyDxiL6zZbVdWJqrpNRL4CfsB1OWuiqm7N5TEYk28tX76cgQMH8sILL/Dggw8GOhxjPMrJpaeNqtoiu9d8zS49mYJmwoQJPP/888yYMYPOnTvneL/ExERCQkJsuHCTI9649JTpJ01E2uKa2a6KiDyWblU5wAaYMeYyJScn8/jjjxMXF8eqVau46qqrcryv0+lkzpw5hIeHExER4cMojfmfrH6SFAfKurdJ3xP7JNDPl0EZU1CdOHGCqKgonE4n69aty3V7xNdff01iYiI33nijjyI0JqOcXHoKV9W9foonqzjs0pPJ13bt2kX37t25+eabefvtt3N96WjTpk2sWrWKESNGULp0aR9FaQoaX196Gq+qfwXeF5EM39Cq2iMvFRtTmDgcDgYOHMjYsWN56KGHcr3/3r17Wbp0KcOGDbMkYfwuq580n7n/fcMfgRhTUE2aNIlnn32W//73v9xyyy2XVcaWLVvo3bs3lStX9nJ0xmQvx4MCBppdejL5TXJyMk888QQLFy4kNjbWxmAyAeHTS0/GmMt34sQJBg4cSHJyMuvWraNixYqBDsmYy2YzmxjjZbt27aJt27bUq1ePhQsXWpIw+Z4lCmO8KD4+nnbt2vHII4/wwQcfUKxYscsqxy6zmmCSk/koGorIJBFZLCLLUh/+CM6Y/GTy5MkMGDCAzz77jEceeeSyy9mzZw8zZsywZGGCRk7aKGYBHwGT8N7w4sYUGCkpKTzxxBN8+eWXrFixgkaNGl12WceOHWP27Nn06dPHhgw3QSMniSJZVT/0eSTG5EMnT55k0KBBXLhwgXXr1hEWFnbZZZ0/f54ZM2bQqVMn6tWr58UojcmbTC89iUiYiIQBsSLysIhcmfqa+3VjCrVff/2Vtm3bEh4ezqJFi/KUJJxOJ7Nnz6Zu3brccMMNXozSmLzL6oxiA6D8by6JJ9KtU8B+8phCa8WKFQwYMIBnn32WUaNG5bm8LVu2AHDbbbfluSxjvM063BmTSx9//DFjxoxh2rRpdOnSxStlqirJycmXfZeUMZnxS4c7EXkE+K+qJriXKwKDVPX/8lKxMflNSkoKTz75JPPmzWPFihU0btzYa2WLiCUJE7RyMnrs96ra/JLXNqnq9T6NLGMcdkZhAubkyZPcddddnD17ltmzZ+epPcIYf/LXnNlFJN19eiJSBNdcFcYUCrt37yYiIoKaNWvy1VdfWZIwhU5OEkUcECMiN4vIzcAM92vGFHgrV64kIiKCBx54gA8//NArl4dSUlJYsGABp0+f9kKExvheThLFk8By4CH3Yynwj5wULiKTReSQiPyQyfpOIpIgIhvdj2dzGrgxvpSQkMD48ePp27cvn3zyCaNHj/ZaB7i4uDgSEhJsXgmTb2TbmK2qTuBD9yO3pgDvAVOz2GaFTYJkgsGpU6eYP38+MTExxMfH07lzZ+Lj47n66qu9Vsf69evZu3cvw4cPJyTEhloz+UNO7nrajavfxEVUNdt+FKq6SkTCs6siu3KM8ZWzZ8+yYMECoqOj+frrr2nfvj0DBw5k2rRplCtXzqt17dq1i5UrVzJ8+HBKlizp1bKN8aWcDOHRKt3zkkB/wJuteW1EZBNwAHhCVbd6sWxjMjh//jxxcXHExMSwaNEiWrduTVRUFJMmTfJZQ/Xp06f54osv6N+/vw07bvKdy+pwJyIbVLVlDrcNB2JVtZmHdWUBp6qeFZHbgXdU1eM0YCKiY8eOTVuOjIwkMjIy17GbwikxMZGvv/6amJgYYmNjadasGVFRUfTt25eqVav6vH5V5fDhw36pyxRuDocDh8ORtvzCCy/k+fbYnPSjaJFuMQTXGcZDqnpdjirIIlF42HY30FJVj3lYZ/0oTK4kJyfjcDiIjo5m7ty5NGrUiKioKPr160f16tUDHZ4xfuGvqVDfTPc8GdgDDMhFHUIm7RAiUk1VD7mft8aVuDIkCWNyyul0snLlSmJiYpgzZw61a9cmKiqKjRs3Urt27UCHZ0y+lJO7nm663MJFZDoQCVQSkX3AWFyd9VRVJwL9ROQhIAk4B0Rdbl0mfzt8+DAJCQk0aNAg17ehqirr1q0jJiaGWbNmUblyZaKiolizZg3169f3UcTGFB45ufT0KK7bXE/hmryoBTBGVRf7PryL4rBLTwWE0+nMcGvovn37iI2NpXTp0nTu3Jnw8KxvllNVNm7cSHR0NDNnzqR06dIMHDiQqKgor47BdLl2795NaGgolStXDnQoppDzxqWnnCSKzap6nYh0BR4EngU+U9UWWe7oZZYoCoakpCQ+++wzbr755gzJwOl0smXLFhwOB5UrV+amm266qC1BVdmyZQsxMTHExMQAEBUVRVRUFE2bNg2aGeEOHz7MJ598QlRUlF3uMgHnrzaK1AruAKaq6k8SLH+RJl9xOp18/vnnVKhQweMXaEhICNdddx3XXnstGzdu5PPPP2fkyJHs3r07LTmcOXOGqKgoYmJiaNGiRdAkh1Rnz55lxowZ3HrrrZYkTIGRk0SxQUQWA3WBp0QkFHD6NixTEC1evJjz58/Tt2/fLL/gixQpQlhYGCdOnKB169YcOXKE/v37M3nyZNq0aRN0ySFVSkoKM2fO5Oqrr6Z58+bZ72BMPpGTS08hQHPgV1VNEJFKQA1V9Th+k6/Ypaf8be3atWzatCnLXsn79u1j5syZxMTEsG/fPvr27UtUVBTt27enSJEiadupalAmi9jYWM6cOcOAAQNseA4TNPxy6ck91tPGdMtHgaN5qdQULhcuXODHH39k8ODBGZLEH3/8waxZs4iJiWH79u306tWLV155hcjISIoW9fzxjImJoVKlSrRr1y6oBtarV68eDRo0sCRhChybCtX4haezgC+++ILhw4fTo0cPoqKiuOWWWyhePPupTk6ePMmKFSvYunUrrVu3pm3btpQoUcJXoRuTr/nlrqdgYYmiYFm/fj3dunUjLi6Oli1zNBpMBseOHSM+Pp5du3Zx0003XXY5xhRkPk8U7tnsflLVgN+Ybomi4NizZw8RERF89NFH9OiR9xHm//zzTxISEmjY0OMwYcYUaj6fClVVU4DtImL3+ZkcczozvykuISGBO+64gzFjxnglSQBUrVrV70kiOTmZw4cP+7VOYwIlJ61uFYGfRGSpiMxPffg6MJM/OZ1OYmJi2LZtW4Z1iYmJ9O3bl1tuuYW//OUvfoll586dePtMVFVZsGABK1eu9Gq5xgSrnPSjeM7nUZgCQVVZtGgRycnJXHXVVRnWPfjgg5QpU4a3337bL/GcPn2apUuXsnz5cjp37ky9evW8clvt2rVrOXjwIPfee68XojQm+OWoMVtEqgE3uBfXq+qfPo3KcwzWRhHkVq1axY8//si9996b4S6kl156ic8//5z4+HjKli3rt5hUla1bt7J8+XJCQ0Pp3LkztWrVuuzytm/fzoIFCxgxYgTly5f3YqTG+IbP2yjclQwA1uOa2W4A8I2I9MtLpabg2bJlC99++y133XVXhiQxY8YMJk6cSGxsrF+TBLj+SJo0acLDDz9Ms2bNmDNnDrt27bqssg4dOsT8+fMZMGCAJQlTqORoUEDg1tSzCBGpAnyd04mLvMXOKIJXSkoKU6ZMoXv37lSrVu2idatWraJPnz58/fXXNGuW7dxVPpecnExISMhldYrbtm0bycnJXHvttT6IzBjf8NfosVtUtWm65RBgc/rX/MESRXDz1KHul19+oUOHDnz66ad07do1QJEZU7j5a/TYOBH5CpjhXo4CFualUlPwXJokjh49Srdu3fjnP/+ZL5LEd999x+HDh+nQoYPfL48ZE+wyPaMQkRKqesH9vA/Q3r1qpap+4af40sdjZxT5xPnz57n11luJiIjg1VdfDXQ4OXL69GlWrVrFDz/8QMuWLYmIiKBUqVKBDsuYPPPppScR2aiqLUTkM1UdkpdKvMESRfDIavRWVWXw4MEkJSURExOT7wbIO3HiBPHx8Wzfvp1OnTrRunXrQIdkTJ74+tJTcRG5C4hwn1FcRFU/z0vFJn9SVWJjY6lZsyYtWmSc5PD5559n9+7dLFu2LN8lCYDy5cvTo0cPjh49yvr167lw4YINOGgKvaz+kh8EOgAVgO6XPO7MSeEiMllEDolIlnNXiMgNIpLkKSGZ4LJy5Ur++OMPmjRpkmHdJ598wvTp05k3b16+v2xTqVIlbr/9dksSxpCzu55GqOrkyypcpD1wGtcUqh7vjXTfRbUEOAd8nNmZil16CrzNmzezfPlyRowYQWho6EXrli1bxqBBg3A4HFx99dUBitAYcym/dLi73CTh3ncVcDybzUYDswG/9/Y2Offrr7+yZMkSBg8enCFJbN26lYEDBxITE2NJwpgCKKAXkUWkOtBLVT8Egm9uSwO42iXi4+Pp168fVapUuWjdoUOH6NatG6+//jqRkZGBCdAY41M56UfhS+OBJ9MtZ5ksxo0bl/Y8MjLSvpj8RES45557MjROnz17lh49ejB06FDuueeeAEVnjEnP4XDgcDi8WmZWt8ferarT3M/bqerqdOtGqer7OapAJByI9dRGISK/pj4FKgNngPtVNcMw5tZGEVycTif9+/endOnSTJ061SujshpjvM/XbRSPpXv+3iXrhueiDiGTMwVVred+1MXVTvGwpyRhgs+TTz7JkSNH+M9//mNJwpgCLqtLT5LJc0/LngsQmQ5EApVEZB8wFigOqKpOvGRzO10IEqlnbpklgI8++oj58+ezdu1au33UmEIgq0ShmTz3tOy5ANW7chqIqubmLMX4UHx8PEWKFKFDhw4Z1i1atIhx48axatUqwsLCAhCdMcbfskoUjd0d5QSon67TnAD1fB6ZCYhNmzaxefNmRowYkWHd5s2bueeee5g7dy4NGjQIQHTGmEDIKlHYDfGFzM6dO1m6dCnDhg3LMILq/v376d69O++99x4REREBitAYEwiZJgpV3Zt+WUQqAR2Bfaq6wdeBGf86ePAgX3zxBVFRUVSuXPmidadPn+bOO+/koYceIioqKkARGmMCJdO7nkTkSxG51v38SuBHXHc7fSYif/VTfMZPli9fzh133EHt2rUvej05OZmBAwfSsmVLxowZE6DojDGBlFU/ip9UtYn7+dNAY1UdKiKhwOrMxm7yFetH4VtOpzNDhzpVZfTo0Wzfvp2FCxdSrFixAEVnjLlcvh5mPCnd85uBSQCqekpEnHmp1AQfT0OCv/POOzgcDlavXm1JwphCLKtE8ZuIjAZ+B1oAcQAiUgqwb40Cbu7cubz++uusWbOG8uXLBzocY0wAZdUzewTQBBgGRKlqgvv1NsAUH8dlfCyry3jfffcdI0eOZO7cuYSHh/sxKmNMMMp2PopgYW0U3rNhwwaOHj1Kly5dMqzbu3cvERERfPDBB/Tq1SsA0RljvMmnbRQikuWYS6raIy8Vm8D45ZdfWL58OcOHZ+wIf+LECbp168bjjz9uScIYkyarNoq2wG/ADOAbbL6IfO/AgQPMnTuXQYMGZRh+IykpiX79+hEZGclf/2p3Pxtj/ier22OLALcCg4BmwAJghqr+5L/wLorHLj3lQUJCAh9//DG33357hlnoVJWRI0dy8OBB5s6dS9GigZ6mxBjjLT4dZlxVU1Q1TlXvwdWAvRNwiMiovFRoAsPhcNCuXTuPU5W++uqrbNiwgejoaEsSxpgMsmzMFpESQDdcZxV1gPnAx6q63y/RXRyLnVHkQUpKCkWKFMnwekxMDE888QRr166lRo0aAYjMGONL3jijyOrS01TgWmAhEK2qP+aloryyROF9a9asoWfPnnz99ddcd911gQ7HGOMDvk4UTlxTk8LF808IromHyuWl4tyyRJG15ORk9uzZw7Zt2wgPD6dp06ZZbr9r1y7at2+f1m5hjCmYfHp7rKpm1RnPBIELFy7wyy+/sG3bNnbt2kWVKlVo3LgxtWrVynK/o0ePcscdd/D8889bkjDGZMs63OVje/fuZfXq1TRu3JiGDRtmmEPCkwsXLtClSxduuOEG3njjDT9EaYwJJJ9eego2hTlRJCQkUKFChTyXo6oMHTqUs2fPMmvWLI8DARpjChZfjx6bZyIyGbgTOORpWHIR6QG8CDhxjVb7N1Vd7cuY8gNVZf/+/Wzbto1t27aRlJTEqFGj8jyC6wsvvMCOHTtYvny5JQljTI759IxCRNoDp4GpmSSK0qp61v28KTBTVT1OwVpYzigcDgcbNmygVKlSNGrUiMaNG1O9enVE8tYxfurUqYwbN461a9dSrVo1L0VrjAl2QX9GoaqrRCTT4UdTk4RbWVxnFoVa7dq1adq0KZUqVfJamQ6Hg8cffxyHw2FJwhiTawHvhisivYBXgCq4OvcVaKdOnWL79u2UKVPGYy/pevXq5al8VeXEiRMcO3aMY8eO8fvvv/PAAw8wY8YMrrnmmjyVbYwpnAKeKFR1LjDXfZnqX7jGl/Jo3Lhxac8jIyOJjIz0dXheceTIEbZt28b27ds5cuQIV111Fddff32W+zidTk6ePMmxY8c4evRo2hd/dsvHjx+ndOnSVKpUibCwMMLCwnjvvfe4+eab/XS0xphAcjgcOBwOr5bp87ue3JeeYnMyx7aI7AJuUNVjHtblyzaK/fv3M2PGDGrUqEH58uUREY4fP57tF39CQgJlypQhLCzsoi99T8vpX6tYsaJNW2qMSRP0bRRuQiZDlItIfVXd5X7eAijuKUkEM1Xlww8/ZPv27R6/+BMSEihbtmymX/J16tShRYsWGdZXqFDBvvCNMUHB17fHTgcigUoisg8YCxTHNQTIRKCviAwFEoFzwABfxuMLb775Jj/++CN169alefPmVK1a9aIv/QoVKtiIrMaYfM063OXBsmXLeP/997njjjsYOHBgjnpGG2OMP/l0PgqTtX379vHcc8/RunVrBg8ebEnCGFNgWaK4DOfPn2fIkCF06dKFoUOHUqpUqUCHZIwxPmOJIpdUlYcffphatWrRv39/qlevHuiQjDHGp6yVNZc++ugjvv32W9auXWuXm4wxhYI1ZufCmjVr6NWrF2vWrKFBgwYBjcUYY3LCGrP96I8//mDAgAFMmTLFkoQxplCxRJEDiYmJREVFMXLkSLp1K/DDURljzEUsUeTAY489RqtWrejRo0egQzHGGL+zRJGNTz/9lD/++IMmTZpw3XXXBTocY4zxO0sUWdiwYQPvvfceN954I1FRUTYrnDGmULJvvkwcPnyYe+65h169etnwHMaYQs1uj/UgOTmZrl27EhERwR133EHbtm39Uq8xxnibN26PtUThwT/+8Q82bdrEggULKFasWJ7nqzbGmEDJL/NR5CszZ85k1qxZfPfddxQvXjzQ4RhjTMDZGUU6P/74IzfddBOLFy/OdqpSY4zJD6xnthclJCTQu3dv3nrrLUsSxhiTjp1RAE6nk+7du9OkSRNefPFFSpQo4ZN6jDHG3+yMwkv++c9/UrFiRa666iqKFCkS6HCMMSao+DRRiMhkETkkIj9ksv4uEdnsfqwSkaa+jMeT2NhYFi5cSLNmzRgwYIDNb22MMZfw9RnFFKBrFut/BTqq6nXAv4BJPo7nIjt27GD06NH06dOHvn37Ur58eX9Wb4wx+YJPfz6r6ioRCc9i/bp0i+uAGr6MJ71Tp07Rp08f7r//ftq2bUv9+vX9VbUxxuQrwXSd5T5gkT8qUlWGDx9OmzZtuPXWW2nVqpU/qjXGmHwpKBKFiNwE3Au090d9r732Gnv37mXFihWULFnSH1UaY0y+FfBEISLNgInAbap6PKttx40bl/Y8MjKSyMjIXNe3ZMkSxo8fz/r16y1JGGMKHIfDgcPh8GqZPu9HISJ1gFhVzXBHk4jUBpYCQy5pr/BUTp77UezZs4c2bdoQHR19WUnGGGPym6AfFFBEpgORQCXgEDAWKA6oqk4UkUlAH2AvIECSqrbOpKw8JYpz587Rrl077r77bh577LHLLscYY/KToE8U3pSXRKGq3HPPPZQtW5YOHTowaNAgL0dnjDHByXpm59D777/Pvn37qFu3Ll26dAl0OMYYk68U+ESxcuVK3njjDe688066d+9OpUqVAh2SMcbkKwG/68mX9u/fz6BBg3j88cdp0qQJjRs3DnRIxhiT7xTYNooLFy4QGRnJHXfcQa1atRg6dCghIQX+BMoYYy5ijdlZePDBBzl06BBz5sxBVW1UWGNMoWRToWZi8uTJxMfH880339hZhDHG5FGBO6NYv3493bp1Y7p7HKIAAA0RSURBVOXKldYmYYwp9Oz22Ev8+eef9OvXj4kTJ1qSMMYYLykwiSI5OZmoqCiGDh3KjTfeGOhwjDGmwCgwieIf//gHFStWpHz58jZLnTHGeFGBSBTTp09n4cKFdOzYka5du1K1atVAh2SMMQVGvm/M3rx5M7fccgsvvfQSNWrUoFu3bgGIzhhjglOh70dx7NgxbrjhBh577DGKFy/Ovffea5edjDEmnULdjyIlJYXBgwfTo0cPbrnlFqpVq2ZJwhhjfCDffrOOHTuWc+fO8dprr1GsWLFAh2OMMQVWvkwUc+fOZerUqXz33XeWJIwxxsfyXRvFtm3b6NChAwsWLKB1a4+T4RljjHErdD2zT548Sa9evXj55ZctSRhjjJ/kq0Rxzz330LlzZ06dOsXZs2cDHY4xxhQKPk0UIjJZRA6JyA+ZrG8kImtE5LyIPJZdeYcPH+baa6+lY8eOlC5d2vsBG2OMycDXZxRTgK5ZrD8KjAZez0lhDz30ENWrV6dly5beiM2rHA5HoEPIlMV2+YI9Pm+yYzWZ8WmiUNVVwPEs1h9R1Q1Ack7KO3nyJN26dUMkT+0yPhHMHzyL7fIFe3zeZMdqMpOv2igGDBhgt8MaY4yf5atEUalSpUCHYIwxhY7P+1GISDgQq6rNsthmLHBKVd/KYpv80eHDGGOCTH4Y60ncj5xsl6m8HqgxxpjL49MzChGZDkQClYBDwFigOKCqOlFEqgHfAaGAEzgNXKOqp30WlDHGmFzJN0N4GGOMCYx80ZgtIreJyDYR2SEiTwY6nvSy61QYSCJSU0SWichPIrJFRP4S6JhSiUgJEflGRDa5Yxsb6JguJSIhIrJRROYHOhZvEZE9IrLZ/b6vz2SbR93/J0H1mcmMp79BEakoIotFZLuIfCUi5TPZ9xER+UVEUkQkLN3rFUTkc/d7tU5ErvHHseRGZn/fOT323Aj6RCEiIcD7uDruNQEGiUjjwEZ1kew6FQZSMvCYqjYB2gKPBMt7p6oXgJtU9XqgOXC7iATbAF6PAlsDHYSXOYFIVb1eVTO83yLSBBgBtML1/3KniNTzc4y55elvcAzwtao2ApYBT2Wy7yrgZmDvJa8/DWxS1euAe4B3vReu12T2953tsYvIFBHpmNOKgj5RAK2BX1R1r6omAdFAzwDHlCa7ToWBpKoHVfV79/P/b+/8g62qqjj++SKiAiPl7x/IQ1EnLW1EoRTzYalj4lg4aiBNlAg2ToY/pqls1JlqBsjJxJpq0CJStNJS1FBBeuo4ipD8fIE2xjOsMW1GkhdjqLzVH3td3/Z673n34oN78a3PzJ27z/511j7n7LP2r7P2f4H1wKGNlaobMysZ7NqDtLCiacZBJQ0FzgFua7QsvYworvfHAM+Y2VYz2wY8AZy/UyTbTqrUwc8B89w9D/h8lbSrzWwj711McyzpJYuZPQ8Ml7R/rwndC1Sp30Opsez1sCsoikOBl7Ljf9BEL7tdBUnDSS3EZxorSTc+tLMS+Bew2MyWN1qmjB8B36CJlFcvYcAjkpZLmlohvB34lA9fDCQpy8N2qoS9wwFm9gqkFypwQJ3pV+MK0nu6w0gv4aYkq99LgQNrLHvNK0l3yY2LgvqQNBi4B5jeTCvKzKwLOEHS3sB9ko41s4YP9UgaB7xiZqskjaWOCrULMMbMXvbW8WJJ671FDoCZPSdpFrCYtApxJbCtQbL2JvUq/JnAbEkrgLU08XUor98Vvjkzj3cWMMuPW4AxkrYA/zOzk4vOsSv0KP5J0uYlhrpfUAOS+pMeotvNbEGj5amEmW0G2oCzGy2LMwY4T9IG4C7gdEm/brBMvYKZvez//wbuBcb6xPYKSdM8bK6ZnWRmY4H/AH9tmMDbzyu+/B5JBwGvuvthL+ucsvjvermaWaeZXWJmI81sMqlVvmFnCF4PVep3xbKb2SKfmxoJLAAu9eNCJQG7hqJYDhwpqUXSAGAC0GyrUGr9qLAR/BJYZ2azGy1IjqT9SqsxJO0FnAk811ipEmZ2rZkNM7MjSM/bn8zsS42W6/0iaaC3PpE0CDgLWFZ6eZjZHA/b3/+HAeOBOxslcx2U18H7gS+7ezLpxYiZne1lnVaUXtIQSbu7eyrweDP1xjMq1e+KZS+jrvdV0ysKn1D7GrAI+AvwGzNb31ipuvGPCp8Cjpa0UdJXGi1TCUljgEnAp7NWY7O02g8G2iStIs2bPGJmCxss0wedA4EnfV5oKcm0zqIK8X4vqZ30grnce3xNS5U6OBM4U9LzpFVNM6ukvULSS6R5z9VZT+MYoF3SetKKquk7uhz1UlC/Z9Fz2esaiosP7oIgCIJCmr5HEQRBEDSWUBRBEARBIaEogiAIgkJCUQRBEASFhKIIgiAICglFEQRBEBQSiqIPIKlL0o3Z8TWSru+lvOdK2uFG4yRdIGmdpCVl/i2SJm5nnk/WEGfOzra4q2Tme8/s+EE3c9L0SDpR0s3ubpXU41e/QfMTiqJvsBU4P7e33wxI2q2O6FNIJgc+U+Z/OHDx9uRvZqf2dFIzm2Zmvf7FuKSiL2OvBAZmMpzbrB+9lV9jM3vWzK70w7HAKTtdqKDXCUXRN3gbmANcXR5Q3iOQ1On/rZIek3SfpBckzZB0sdJmQ6slHZ5lc6ZbI33ODeqVLMP+wOOvKlkq9XyfkLSA9KV9uTwTJa3x3wz3uw44FfiFG6zLmQGc6l+lTpc0WdIC73k8KmmQpEcl/dnlPq9KWdsk3S1pvaTbszhtkkaW4kv6vpfnqczUxRGSnvb8v1fKt6xcLX595klaCwyV9FNJy5Rt3CTpCuAQ0lfrS9yvQ9I+nsc67+W0K9kt2sPjjPLzr/Drvtb9j/V7sMLlHlFBtk5JN3meiyXtm5XrIb+3j0s6OntmfiZpKekr4DyvVkkPSGoBvgpc6eceo2S25R6X55lSb0PSDZJ+5c9Fh6Txkmb5M7CwzgZFsCMws/h9wH/AZmAw0EHan/wa4HoPmwucn8f1/1bgNZIxtAEk8+43eNjXgZuy9AvdfSTJJPwAYCpwrfsPINnsavF8O4FhFeQ8mLSBzD6kRswS4DwPawNOqJCmFbg/O54MbASG+HE/YLC79yXtbVKprJv8/CKZgzglO+9Id3cB57h7Vla+B4CL3H1ZKd8yOVtICntU5vehTMY24GN+vAH4cBZvg1+TFuBN4Dj3/y1wsbvXAqPdPQNY4+5bgInu7g/sUUG2LmCCu68DbnH3o8AId48GlmT3/P7yfMrvB3ADaWOdUtj87LoeRrJRVIr3hF+H44EtwFke9ofSMxC/xv3CzHgfwZL54XkkmzVv1JhsuZmVrG7+jWRvC9JLaWwW73d+jhc83kdIBueOk3Shx9kbOAp4i2SIbmOF840C2szsNT/nfOA0uo1A1mrIbLGZve7ufsAMpd28uoBDJB1QKlfGMnPLqkr2p4aTFEbOVuu2R/UscIa7T6Z7M607gRupzN/t3XtuTPCeVn/gINJmOe2818Bd7u4ws7WZDMOVjCsONrPS1qZ3AuPc/TTwHaWNmO41sxcqyLUNv4fAHSRbT4NIw0Z3S+8Mk+2epbm7ShmLOAM4JstvsNKeFwAPmVmX94T6WbcNqrWkexE0kFAUfYvZwApSi7DE2/gQpFfgAVnY1szdlR138e5nJzcYJj8WcIWZLc4FkNRKajFWozes8Ob5TwL2I/VGuiR1AHtWSJOXdRuV68ZbVeKUl79HuZQ2mrkGONHMNkuaW0WunuQspal4XjO7y4eIzgUWSppmZo/1cA4jPRObLJmkrkTRPayGgE9Y2qmy2zPpja0ur0nKw8uftaABxBxF30AAZraJ1HKckoW9SNofGVKreHfq50IlRpAml58HHgEuV7KXj6SjstZjNZYBp/l4/G7AROCxHtJ0kobTqjEEeNWVxOmk4ZsS9SqlavGXAhe4e0KN6fcmbQzUqbR3wGezsM0eXpMM3nvaLGlUuQySDjezDjP7Mcka7PEV8twtk38S8KSZdQIdkkr+SKqUtojOsnIsIrPCKunjVdI1q8n+Pksoir5B3uL9IWmsvuR3K9CqZHr6k1RvKRaZGd5Iesn/EbjMzN4k7TW9Dljhwwk/J72QqguZtm38Fkk5rCQNfT3Yw/nXAF1KZpanV4g3HxglaTXwRdK+wj2VyWpw51wFXO1DViOA16vEeye9ma0BVrk8dwD5Ut1bgYfVvRS4FhkuBW5T2pFtYCbDRT5JvRL4KFBpA6YtwGi/T2OB77r/JGCKT4K3A6WFALWanH4AGF+azCbNbZ3kk+7tpPmcSoRJ6yYjzIwHwftE0l5m9oa7v0CaGB6/k2UYZGZb3P1N4CAzu6rGtJ1mVtQrC/o4MfYXBO+fEyX9hDRksgm4pAEyjJP0bVKdfpHuHc5qIVqLQSHRowiCIAgKiTmKIAiCoJBQFEEQBEEhoSiCIAiCQkJRBEEQBIWEogiCIAgKCUURBEEQFPJ/a8eJ5n/E3p4AAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.xlim(0, max(X))\n", "plt.plot(X,YF,color='k',label='with features')\n", "plt.plot(X,Y,color='grey',linestyle='--',label='without features')\n", "plt.xticks(X,Xlab)\n", "plt.xlabel(\"Number of training ratings per item\")\n", "plt.ylabel(\"MSE for such items\")\n", "plt.title(\"Performance versus item `coolness'\")\n", "plt.legend(loc=\"best\")\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.3" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Read social data from epinions" ] }, { "cell_type": "code", "execution_count": 38, "metadata": {}, "outputs": [], "source": [ "userIDs = {}\n", "itemIDs = {}\n", "interactions = []\n", "\n", "socialTrust = defaultdict(set)\n", "\n", "f = open(dataDir + \"epinions_data/epinions.txt\", 'rb')\n", "header = f.readline()\n", "\n", "for l in f:\n", " try:\n", " l = l.decode('utf-8')\n", " l = l.split()\n", " except Exception as e:\n", " continue\n", " i = l[0]\n", " u = l[1]\n", " \n", " if not u in userIDs: userIDs[u] = len(userIDs)\n", " if not i in itemIDs: itemIDs[i] = len(itemIDs)\n", " interactions.append((u,i))\n", "\n", "f.close()\n", "\n", "f = open(dataDir + \"epinions_data/network_trust.txt\", 'r')\n", "for l in f:\n", " try:\n", " u,_,v = l.strip().split()\n", " except Exception as e:\n", " continue\n", " if u in userIDs and v in userIDs:\n", " socialTrust[u].add(v)\n", " \n", "f.close()" ] }, { "cell_type": "code", "execution_count": 39, "metadata": {}, "outputs": [], "source": [ "random.shuffle(interactions)" ] }, { "cell_type": "code", "execution_count": 40, "metadata": {}, "outputs": [], "source": [ "items = list(itemIDs.keys())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "BPR model. First we'll use a regular BPR model just to assess similarity between friends' latent representations. Later we can implement different social sampling assumptions just by passing different samples to the same model." ] }, { "cell_type": "code", "execution_count": 41, "metadata": {}, "outputs": [], "source": [ "class BPRbatch(tf.keras.Model):\n", " def __init__(self, K, lamb):\n", " super(BPRbatch, self).__init__()\n", " # Initialize variables\n", " self.betaI = tf.Variable(tf.random.normal([len(itemIDs)],stddev=0.001))\n", " self.gammaU = tf.Variable(tf.random.normal([len(userIDs),K],stddev=0.001))\n", " self.gammaI = tf.Variable(tf.random.normal([len(itemIDs),K],stddev=0.001))\n", " # Regularization coefficient\n", " self.lamb = lamb\n", "\n", " # Prediction for a single instance\n", " def predict(self, u, i):\n", " p = self.betaI[i] + tf.tensordot(self.gammaU[u], self.gammaI[i], 1)\n", " return p\n", "\n", " # Regularizer\n", " def reg(self):\n", " return self.lamb * (tf.nn.l2_loss(self.betaI) +\\\n", " tf.nn.l2_loss(self.gammaU) +\\\n", " tf.nn.l2_loss(self.gammaI))\n", " \n", " def score(self, sampleU, sampleI):\n", " u = tf.convert_to_tensor(sampleU, dtype=tf.int32)\n", " i = tf.convert_to_tensor(sampleI, dtype=tf.int32)\n", " beta_i = tf.nn.embedding_lookup(self.betaI, i)\n", " gamma_u = tf.nn.embedding_lookup(self.gammaU, u)\n", " gamma_i = tf.nn.embedding_lookup(self.gammaI, i)\n", " x_ui = beta_i + tf.reduce_sum(tf.multiply(gamma_u, gamma_i), 1)\n", " return x_ui\n", "\n", " def call(self, sampleU, sampleI, sampleJ):\n", " x_ui = self.score(sampleU, sampleI)\n", " x_uj = self.score(sampleU, sampleJ)\n", " return -tf.reduce_mean(tf.math.log(tf.math.sigmoid(x_ui - x_uj)))" ] }, { "cell_type": "code", "execution_count": 42, "metadata": {}, "outputs": [], "source": [ "def trainingStepBPR(model, interactions):\n", " Nsamples = 50000\n", " with tf.GradientTape() as tape:\n", " sampleU, sampleI, sampleJ = [], [], []\n", " for _ in range(Nsamples):\n", " u,i = random.choice(interactions) # positive sample\n", " j = random.choice(items) # negative sample\n", " while j in itemsPerUser[u]:\n", " j = random.choice(items)\n", " sampleU.append(userIDs[u])\n", " sampleI.append(itemIDs[i])\n", " sampleJ.append(itemIDs[j])\n", "\n", " loss = model(sampleU,sampleI,sampleJ)\n", " loss += model.reg()\n", " gradients = tape.gradient(loss, model.trainable_variables)\n", " optimizer.apply_gradients((grad, var) for\n", " (grad, var) in zip(gradients, model.trainable_variables)\n", " if grad is not None)\n", " return loss.numpy()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First, train a regular BPR model (no social terms)" ] }, { "cell_type": "code", "execution_count": 43, "metadata": {}, "outputs": [], "source": [ "optimizer = tf.keras.optimizers.Adam(0.1)\n", "modelBPR = BPRbatch(10, 0.00001)" ] }, { "cell_type": "code", "execution_count": 44, "metadata": {}, "outputs": [], "source": [ "nTrain = int(len(interactions) * 0.9)\n", "nTest = len(interactions) - nTrain\n", "interactionsTrain = interactions[:nTrain]\n", "interactionsTest = interactions[nTrain:]" ] }, { "cell_type": "code", "execution_count": 45, "metadata": {}, "outputs": [], "source": [ "itemsPerUser = defaultdict(list)\n", "usersPerItem = defaultdict(list)\n", "for u,i in interactionsTrain:\n", " itemsPerUser[u].append(i)\n", " usersPerItem[i].append(u)" ] }, { "cell_type": "code", "execution_count": 46, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "iteration 10, objective = 0.5586551\n", "iteration 20, objective = 0.5336623\n", "iteration 30, objective = 0.5308751\n", "iteration 40, objective = 0.53656584\n", "iteration 50, objective = 0.5435766\n", "iteration 60, objective = 0.543326\n", "iteration 70, objective = 0.54365164\n", "iteration 80, objective = 0.54245454\n", "iteration 90, objective = 0.5428503\n", "iteration 100, objective = 0.5451498\n" ] } ], "source": [ "for i in range(100):\n", " obj = trainingStepBPR(modelBPR, interactions)\n", " if (i % 10 == 9): print(\"iteration \" + str(i+1) + \", objective = \" + str(obj))" ] }, { "cell_type": "code", "execution_count": 47, "metadata": {}, "outputs": [], "source": [ "interactionsTestPerUser = defaultdict(set)\n", "itemSet = set()\n", "for u,i in interactionsTest:\n", " interactionsTestPerUser[u].add(i)\n", " itemSet.add(i)" ] }, { "cell_type": "code", "execution_count": 48, "metadata": {}, "outputs": [], "source": [ "def AUCu(model, u, N):\n", " win = 0\n", " if N > len(interactionsTestPerUser[u]):\n", " N = len(interactionsTestPerUser[u])\n", " positive = random.sample(interactionsTestPerUser[u],N)\n", " negative = random.sample(itemSet.difference(interactionsTestPerUser[u]),N)\n", " for i,j in zip(positive,negative):\n", " si = model.predict(userIDs[u], itemIDs[i]).numpy()\n", " sj = model.predict(userIDs[u], itemIDs[j]).numpy()\n", " if si > sj:\n", " win += 1\n", " return win/N" ] }, { "cell_type": "code", "execution_count": 49, "metadata": {}, "outputs": [], "source": [ "def AUC(model):\n", " av = []\n", " for u in interactionsTestPerUser:\n", " av.append(AUCu(model, u, 10))\n", " return sum(av) / len(av)" ] }, { "cell_type": "code", "execution_count": 50, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.6800902434544706" ] }, "execution_count": 50, "metadata": {}, "output_type": "execute_result" } ], "source": [ "AUC(modelBPR)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute similarities among friends' latent representations" ] }, { "cell_type": "code", "execution_count": 51, "metadata": {}, "outputs": [], "source": [ "sims = []\n", "simFriends = []\n", "while len(sims) < 10000:\n", " try:\n", " u,i = random.choice(interactions)\n", " v = random.sample(socialTrust[u],1)[0] # trust link\n", " j = random.sample(itemsPerUser[v],1)[0] # friend's item\n", " k = random.choice(items) # random item\n", " except Exception as e:\n", " continue\n", " s1 = 1 - distance.cosine(modelBPR.gammaI[itemIDs[i]],modelBPR.gammaI[itemIDs[k]])\n", " s2 = 1 - distance.cosine(modelBPR.gammaI[itemIDs[i]],modelBPR.gammaI[itemIDs[j]])\n", " if s1 > 1:\n", " print(\"?\")\n", " break\n", " sims.append(s1)\n", " simFriends.append(s2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Similarity between randomly chosen pairs of items" ] }, { "cell_type": "code", "execution_count": 52, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.003534959910223597" ] }, "execution_count": 52, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sum(sims)/len(sims)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Similarity between an item and one consumed by a friend" ] }, { "cell_type": "code", "execution_count": 53, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.04061316501272158" ] }, "execution_count": 53, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sum(simFriends)/len(simFriends)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "(similarity is not particularly high, but still significantly higher than random pairs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 6.4" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Implement the social model. Uses the model above, just with different samples." ] }, { "cell_type": "code", "execution_count": 54, "metadata": {}, "outputs": [], "source": [ "def trainingStepBPRsocial(model, interactions):\n", " Nsamples = 50000\n", " with tf.GradientTape() as tape:\n", " sampleU, sampleI, sampleJ = [], [], []\n", " while len(sampleU) < Nsamples/2:\n", " try:\n", " u,i = random.choice(interactions) # positive sample\n", " v = random.sample(socialTrust[u],1)[0] # trust link\n", " j = random.sample(itemsPerUser[v],1)[0] # friend's item\n", " k = random.choice(items) # negative item\n", " if j in itemsPerUser[u] or k in itemsPerUser[u]:\n", " continue\n", " except Exception as e:\n", " continue\n", " while j in itemsPerUser[u]:\n", " j = random.choice(items)\n", " sampleU.append(userIDs[u])\n", " sampleI.append(itemIDs[i]) # Positive\n", " sampleJ.append(itemIDs[j]) # greater than social\n", " sampleU.append(userIDs[u])\n", " sampleI.append(itemIDs[j]) # Social\n", " sampleJ.append(itemIDs[k]) # greater than negative \n", "\n", " loss = model(sampleU,sampleI,sampleJ)\n", " loss += model.reg()\n", " gradients = tape.gradient(loss, model.trainable_variables)\n", " optimizer.apply_gradients((grad, var) for\n", " (grad, var) in zip(gradients, model.trainable_variables)\n", " if grad is not None)\n", " return loss.numpy()" ] }, { "cell_type": "code", "execution_count": 55, "metadata": {}, "outputs": [], "source": [ "optimizer = tf.keras.optimizers.Adam(0.1)\n", "modelBPRsocial = BPRbatch(10, 0.00001)" ] }, { "cell_type": "code", "execution_count": 56, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "iteration 10, objective = 0.6443011\n", "iteration 20, objective = 0.618395\n", "iteration 30, objective = 0.6267633\n", "iteration 40, objective = 0.6299376\n", "iteration 50, objective = 0.6279181\n", "iteration 60, objective = 0.6186671\n", "iteration 70, objective = 0.6207659\n", "iteration 80, objective = 0.6173624\n", "iteration 90, objective = 0.61636496\n", "iteration 100, objective = 0.60944057\n" ] } ], "source": [ "for i in range(100):\n", " obj = trainingStepBPRsocial(modelBPRsocial, interactions)\n", " if (i % 10 == 9): print(\"iteration \" + str(i+1) + \", objective = \" + str(obj))" ] }, { "cell_type": "code", "execution_count": 57, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.6322810869785714" ] }, "execution_count": 57, "metadata": {}, "output_type": "execute_result" } ], "source": [ "AUC(modelBPRsocial)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.5.2" } }, "nbformat": 4, "nbformat_minor": 2 }