w207_keith_lomurray_project2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Project 2: Topic Classification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this project, you'll work with text data from newsgroup postings on a variety of topics. You'll train classifiers to distinguish between the topics based on the text of the posts. Whereas with digit classification, the input is relatively dense: a 28x28 matrix of pixels, many of which are non-zero, here we'll represent each document with a \"bag-of-words\" model. As you'll see, this makes the feature representation quite sparse -- only a few words of the total vocabulary are active in any given document. The bag-of-words assumption here is that the label depends only on the words; their order is not important.\n",
    "\n",
    "The SK-learn documentation on feature extraction will prove useful:\n",
    "http://scikit-learn.org/stable/modules/feature_extraction.html\n",
    "\n",
    "Each problem can be addressed succinctly with the included packages -- please don't add any more. Grading will be based on writing clean, commented code, along with a few short answers.\n",
    "\n",
    "As always, you're welcome to work on the project in groups and discuss ideas on the course wall, but please prepare your own write-up and write your own code."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/klomurray/anaconda3/lib/python3.6/site-packages/sklearn/cross_validation.py:41: DeprecationWarning: This module was deprecated in version 0.18 in favor of the model_selection module into which all the refactored classes and functions are moved. Also note that the interface of the new CV iterators are different from that of this module. This module will be removed in 0.20.\n",
      "  \"This module will be removed in 0.20.\", DeprecationWarning)\n",
      "/Users/klomurray/anaconda3/lib/python3.6/site-packages/sklearn/grid_search.py:42: DeprecationWarning: This module was deprecated in version 0.18 in favor of the model_selection module into which all the refactored classes and functions are moved. This module will be removed in 0.20.\n",
      "  DeprecationWarning)\n"
     ]
    }
   ],
   "source": [
    "# This tells matplotlib not to try opening a new window for each plot.\n",
    "%matplotlib inline\n",
    "\n",
    "# General libraries.\n",
    "import re\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# SK-learn libraries for learning.\n",
    "from sklearn.pipeline import Pipeline\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.naive_bayes import BernoulliNB\n",
    "from sklearn.naive_bayes import MultinomialNB\n",
    "from sklearn.grid_search import GridSearchCV\n",
    "\n",
    "# SK-learn libraries for evaluation.\n",
    "from sklearn.metrics import confusion_matrix\n",
    "from sklearn.metrics import f1_score\n",
    "from sklearn import metrics\n",
    "from sklearn.metrics import classification_report\n",
    "\n",
    "# SK-learn library for importing the newsgroup data.\n",
    "from sklearn.datasets import fetch_20newsgroups\n",
    "\n",
    "# SK-learn libraries for feature extraction from text.\n",
    "from sklearn.feature_extraction.text import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load the data, stripping out metadata so that we learn classifiers that only use textual features. By default, newsgroups data is split into train and test sets. We further split the test so we have a dev set. Note that we specify 4 categories to use for this project. If you remove the categories argument from the fetch function, you'll get all 20 categories."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "training label shape: (2034,)\n",
      "test label shape: (677,)\n",
      "dev label shape: (676,)\n",
      "labels names: ['alt.atheism', 'comp.graphics', 'sci.space', 'talk.religion.misc']\n"
     ]
    }
   ],
   "source": [
    "categories = ['alt.atheism', 'talk.religion.misc', 'comp.graphics', 'sci.space']\n",
    "newsgroups_train = fetch_20newsgroups(subset='train',\n",
    "                                      remove=('headers', 'footers', 'quotes'),\n",
    "                                      categories=categories)\n",
    "newsgroups_test = fetch_20newsgroups(subset='test',\n",
    "                                     remove=('headers', 'footers', 'quotes'),\n",
    "                                     categories=categories)\n",
    "\n",
    "num_test = int(len(newsgroups_test.target)/2)\n",
    "test_data, test_labels = newsgroups_test.data[num_test:], newsgroups_test.target[num_test:]\n",
    "dev_data, dev_labels = newsgroups_test.data[:num_test], newsgroups_test.target[:num_test]\n",
    "train_data, train_labels = newsgroups_train.data, newsgroups_train.target\n",
    "\n",
    "print('training label shape:', train_labels.shape)\n",
    "print('test label shape:', test_labels.shape)\n",
    "print('dev label shape:', dev_labels.shape)\n",
    "print('labels names:', newsgroups_train.target_names)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(1) For each of the first 5 training examples, print the text of the message along with the label."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Example: 1\n",
      "\n",
      "Label: comp.graphics\n",
      "\n",
      "Text of Message:\n",
      "\n",
      " Hi,\n",
      "\n",
      "I've noticed that if you only save a model (with all your mapping planes\n",
      "positioned carefully) to a .3DS file that when you reload it after restarting\n",
      "3DS, they are given a default position and orientation.  But if you save\n",
      "to a .PRJ file their positions/orientation are preserved.  Does anyone\n",
      "know why this information is not stored in the .3DS file?  Nothing is\n",
      "explicitly said in the manual about saving texture rules in the .PRJ file. \n",
      "I'd like to be able to read the texture rule information, does anyone have \n",
      "the format for the .PRJ file?\n",
      "\n",
      "Is the .CEL file format available from somewhere?\n",
      "\n",
      "Rych\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example: 2\n",
      "\n",
      "Label: talk.religion.misc\n",
      "\n",
      "Text of Message:\n",
      "\n",
      " \n",
      "\n",
      "Seems to be, barring evidence to the contrary, that Koresh was simply\n",
      "another deranged fanatic who thought it neccessary to take a whole bunch of\n",
      "folks with him, children and all, to satisfy his delusional mania. Jim\n",
      "Jones, circa 1993.\n",
      "\n",
      "\n",
      "Nope - fruitcakes like Koresh have been demonstrating such evil corruption\n",
      "for centuries.\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example: 3\n",
      "\n",
      "Label: sci.space\n",
      "\n",
      "Text of Message:\n",
      "\n",
      " \n",
      " >In article <1993Apr19.020359.26996@sq.sq.com>, msb@sq.sq.com (Mark Brader) \n",
      "\n",
      "MB>                                                             So the\n",
      "MB> 1970 figure seems unlikely to actually be anything but a perijove.\n",
      "\n",
      "JG>Sorry, _perijoves_...I'm not used to talking this language.\n",
      "\n",
      "Couldn't we just say periapsis or apoapsis?\n",
      "\n",
      " \n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example: 4\n",
      "\n",
      "Label: alt.atheism\n",
      "\n",
      "Text of Message:\n",
      "\n",
      " I have a request for those who would like to see Charley Wingate\n",
      "respond to the \"Charley Challenges\" (and judging from my e-mail, there\n",
      "appear to be quite a few of you.)  \n",
      "\n",
      "It is clear that Mr. Wingate intends to continue to post tangential or\n",
      "unrelated articles while ingoring the Challenges themselves.  Between\n",
      "the last two re-postings of the Challenges, I noted perhaps a dozen or\n",
      "more posts by Mr. Wingate, none of which answered a single Challenge.  \n",
      "\n",
      "It seems unmistakable to me that Mr. Wingate hopes that the questions\n",
      "will just go away, and he is doing his level best to change the\n",
      "subject.  Given that this seems a rather common net.theist tactic, I\n",
      "would like to suggest that we impress upon him our desire for answers,\n",
      "in the following manner:\n",
      "\n",
      "1. Ignore any future articles by Mr. Wingate that do not address the\n",
      "Challenges, until he answers them or explictly announces that he\n",
      "refuses to do so.\n",
      "\n",
      "--or--\n",
      "\n",
      "2. If you must respond to one of his articles, include within it\n",
      "something similar to the following:\n",
      "\n",
      "    \"Please answer the questions posed to you in the Charley Challenges.\"\n",
      "\n",
      "Really, I'm not looking to humiliate anyone here, I just want some\n",
      "honest answers.  You wouldn't think that honesty would be too much to\n",
      "ask from a devout Christian, would you?  \n",
      "\n",
      "Nevermind, that was a rhetorical question.\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example: 5\n",
      "\n",
      "Label: sci.space\n",
      "\n",
      "Text of Message:\n",
      "\n",
      " AW&ST  had a brief blurb on a Manned Lunar Exploration confernce\n",
      "May 7th  at Crystal City Virginia, under the auspices of AIAA.\n",
      "\n",
      "Does anyone know more about this?  How much, to attend????\n",
      "\n",
      "Anyone want to go?\n",
      "\n",
      " ************************************************************ \n",
      "\n"
     ]
    }
   ],
   "source": [
    "def P1(num_examples=5):\n",
    "### STUDENT START ###\n",
    "    \n",
    "    for i in range(num_examples):\n",
    "        print(\"Example:\", i+1)\n",
    "        print(\"\\nLabel:\", newsgroups_train.target_names[train_labels[i]])\n",
    "        print(\"\\nText of Message:\\n\\n\", train_data[i])\n",
    "        print('\\n', '*'*60, '\\n')\n",
    "\n",
    "### STUDENT END ###\n",
    "P1()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(2) Use CountVectorizer to turn the raw training text into feature vectors. You should use the fit_transform function, which makes 2 passes through the data: first it computes the vocabulary (\"fit\"), second it converts the raw text into feature vectors using the vocabulary (\"transform\").\n",
    "\n",
    "The vectorizer has a lot of options. To get familiar with some of them, write code to answer these questions:\n",
    "\n",
    "a. The output of the transform (also of fit_transform) is a sparse matrix: http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.sparse.csr_matrix.html. What is the size of the vocabulary? What is the average number of non-zero features per example? What fraction of the entries in the matrix are non-zero? Hint: use \"nnz\" and \"shape\" attributes.\n",
    "\n",
    "b. What are the 0th and last feature strings (in alphabetical order)? Hint: use the vectorizer's get_feature_names function.\n",
    "\n",
    "c. Specify your own vocabulary with 4 words: [\"atheism\", \"graphics\", \"space\", \"religion\"]. Confirm the training vectors are appropriately shaped. Now what's the average number of non-zero features per example?\n",
    "\n",
    "d. Instead of extracting unigram word features, use \"analyzer\" and \"ngram_range\" to extract bigram and trigram character features. What size vocabulary does this yield?\n",
    "\n",
    "e. Use the \"min_df\" argument to prune words that appear in fewer than 10 documents. What size vocabulary does this yield?\n",
    "\n",
    "f. Using the standard CountVectorizer, what fraction of the words in the dev data are missing from the vocabulary? Hint: build a vocabulary for both train and dev and look at the size of the difference."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a.\n",
      "Size of vocabulary:  26879\n",
      "Non-zero entries:  96.70599803343165\n",
      "Fraction non-zero entries:  0.0035978272269590263\n",
      "\n",
      "b.\n",
      "First Feature:  00\n",
      "Last Feature:  zyxel\n",
      "\n",
      "c. \n",
      "Shape:  4\n",
      "Average number of non-zero features per example:  0.26843657817109146\n",
      "\n",
      "d.\n",
      "Vocaculary size:  35478\n",
      "\n",
      "e. \n",
      "With pruning, size of vocabulary: 3064\n",
      "\n",
      "f.\n",
      "Fraction of the words in the dev data missing not in training dataset:  0.24787640034470024\n"
     ]
    }
   ],
   "source": [
    "def P2():\n",
    "### STUDENT START ###\n",
    "    cv = CountVectorizer()\n",
    "    cv_feature=cv.fit_transform(train_data)\n",
    "    \n",
    "    print(\"a.\")\n",
    "    print(\"Size of vocabulary: \", cv_feature.shape[1])\n",
    "    print(\"Non-zero entries: \", cv_feature.nnz/len(train_data))\n",
    "    print(\"Fraction non-zero entries: \", float(cv_feature.nnz) / (cv_feature.shape[0] * cv_feature.shape[1]))\n",
    "\n",
    "    \n",
    "    print(\"\\nb.\")\n",
    "    print(\"First Feature: \", cv.get_feature_names()[0])\n",
    "    print(\"Last Feature: \", cv.get_feature_names()[-1])\n",
    "    \n",
    "    print(\"\\nc. \")\n",
    "    cv2 = CountVectorizer(vocabulary=[\"atheism\",\"graphics\",\"space\",\"religion\"])\n",
    "    cv_feature2 = cv2.fit_transform(train_data)\n",
    "    print(\"Shape: \", cv_feature2.shape[1])\n",
    "    print('Average number of non-zero features per example: ', cv_feature2.nnz/float(cv_feature2.shape[0]))\n",
    "    \n",
    "    print(\"\\nd.\")\n",
    "    cv3 = CountVectorizer(ngram_range =(2,3), analyzer='char')\n",
    "    bi_tri_feature3 = cv3.fit_transform(train_data)  \n",
    "    print(\"Vocaculary size: \", bi_tri_feature3.shape[1])\n",
    "    \n",
    "    print(\"\\ne. \")\n",
    "    cv4=CountVectorizer(min_df=10)\n",
    "    cv_feature_min_df=cv4.fit_transform(train_data)\n",
    "    print(\"With pruning, size of vocabulary:\",cv_feature_min_df.shape[1])\n",
    "          \n",
    "    print(\"\\nf.\") \n",
    "    cv = CountVectorizer()\n",
    "    cv_feature=cv.fit_transform(train_data)\n",
    "    train_vocab=cv.vocabulary_\n",
    "    dev_feature=cv.fit_transform(dev_data)\n",
    "    dev_vocab=cv.vocabulary_\n",
    "    print(\"Fraction of the words in the dev data missing not in training dataset: \", \n",
    "        sum(list(map(lambda x: x not in train_vocab,list(dev_vocab.keys()))))/len(dev_vocab))\n",
    "\n",
    "    ### STUDENT END ###\n",
    "P2()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(3) Use the default CountVectorizer options and report the f1 score (use metrics.f1_score) for a k nearest neighbors classifier; find the optimal value for k. Also fit a Multinomial Naive Bayes model and find the optimal value for alpha. Finally, fit a logistic regression model and find the optimal value for the regularization strength C using l2 regularization. A few questions:\n",
    "\n",
    "a. Why doesn't nearest neighbors work well for this problem?\n",
    "\n",
    "b. Any ideas why logistic regression doesn't work as well as Naive Bayes?\n",
    "\n",
    "c. Logistic regression estimates a weight vector for each class, which you can access with the coef\\_ attribute. Output the sum of the squared weight values for each class for each setting of the C parameter. Briefly explain the relationship between the sum and the value of C."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal k: {'n_neighbors': 151}\n",
      "KNN f1 score for k = 151 is 0.4028149344405508\n",
      "\n",
      "Optimal alpha for MNB:  {'alpha': 0.01}\n",
      "MNB f1-score for optimal alpha is: 0.7552894198013128\n",
      "\n",
      "Optimal C for LR:  {'C': 0.5}\n",
      "Logistic Regression f1 Score for C = 0.5 is 0.6920151915112704\n",
      "C:  0.0001 ;\tSum of Squared weights:\t [0.00770175 0.0119412  0.00943508 0.00910284]\n",
      "C:  0.001 ;\tSum of Squared weights:\t [0.16509345 0.20095275 0.18067094 0.18724278]\n",
      "C:  0.01 ;\tSum of Squared weights:\t [2.54148768 2.93970937 2.8624691  2.25002919]\n",
      "C:  0.1 ;\tSum of Squared weights:\t [27.13623511 24.65956732 27.45872348 23.02600455]\n",
      "C:  0.5 ;\tSum of Squared weights:\t [102.51345719  83.10215817  98.9828885   89.00084388]\n",
      "C:  1.0 ;\tSum of Squared weights:\t [166.94186105 130.89897162 157.91550429 145.74696578]\n",
      "C:  2.0 ;\tSum of Squared weights:\t [257.43567163 197.83380949 239.94818001 226.51650405]\n",
      "C:  10.0 ;\tSum of Squared weights:\t [586.54870411 448.43848855 538.92744264 530.92736584]\n"
     ]
    }
   ],
   "source": [
    "def P3():\n",
    "### STUDENT START ###\n",
    "\n",
    "    cv=CountVectorizer()\n",
    "    train_cv=cv.fit_transform(train_data)\n",
    "    dev_cv=cv.transform(dev_data)\n",
    " \n",
    "    train_data_matrix = cv.transform(train_data)\n",
    "    dev_data_matrix = cv.transform(dev_data)\n",
    "    \n",
    "    #KNN\n",
    "    #Use GridSearch to find optimal KNN\n",
    "    num = np.arange(1,(len(dev_data) + 1))\n",
    "    n_range = {'n_neighbors': num}\n",
    "    knn = KNeighborsClassifier()\n",
    "    knn_clf = GridSearchCV(knn, n_range)\n",
    "    knn_clf.fit(train_data_matrix, train_labels)\n",
    "    print('Optimal k:', knn_clf.best_params_)\n",
    "    \n",
    "    #Identify F-Score for opimal KNN\n",
    "    knn_clf1 = KNeighborsClassifier(n_neighbors = knn_clf.best_params_['n_neighbors'])\n",
    "    knn_clf1.fit(train_data_matrix, train_labels)\n",
    "    knn_clf_preds = knn_clf1.predict(dev_data_matrix)\n",
    "    print('KNN f1 score for k = {0} is {1}'.format(knn_clf.best_params_['n_neighbors'],\n",
    "        metrics.f1_score(y_true = dev_labels, y_pred = knn_clf_preds, average=\"macro\")))\n",
    "\n",
    "    # Multinomial Naive Bayes\n",
    "    alphas = {'alpha': [0.0001, 0.001, 0.01, 0.1, 0.5, 1.0, 2.0, 10.0]}\n",
    "    mnb_clf = MultinomialNB()\n",
    "    clf = GridSearchCV(mnb_clf, alphas)\n",
    "    clf.fit(train_data_matrix, train_labels)\n",
    "    print('\\nOptimal alpha for MNB: ', clf.best_params_)\n",
    "    \n",
    "    mnb_clf = MultinomialNB(alpha = clf.best_params_['alpha'])\n",
    "    mnb_clf.fit(train_data_matrix, train_labels)\n",
    "    mnb_clf_preds = mnb_clf.predict(dev_data_matrix)\n",
    "    print('MNB f1-score for optimal alpha is:', \n",
    "          metrics.f1_score(y_true = dev_labels, y_pred = mnb_clf_preds, average=\"macro\"))\n",
    "\n",
    "    #Logistic Regression\n",
    "    C = {'C': [0.0001, 0.001, 0.01, 0.1, 0.5, 1.0, 2.0, 10.0]}\n",
    "    \n",
    "    log_clf = LogisticRegression(penalty = 'l2')\n",
    "    clf = GridSearchCV(log_clf, C)\n",
    "    clf.fit(train_data_matrix, train_labels)\n",
    "    log_clf_preds = clf.predict(dev_data_matrix)\n",
    "    print('\\nOptimal C for LR: ', clf.best_params_)\n",
    "    print('Logistic Regression f1 Score for C = {0} is {1}'.format(clf.best_params_['C'],\n",
    "        metrics.f1_score(y_true = dev_labels, y_pred = log_clf_preds, average=\"macro\")))\n",
    "    \n",
    "    #Print sum of the squared weight values for each class \n",
    "    for c in C['C']:\n",
    "        log_clf = LogisticRegression(penalty = 'l2', C = c)\n",
    "        log_clf.fit(train_data_matrix, train_labels)\n",
    "        log_clf_preds = log_clf.predict(dev_data_matrix)\n",
    "        print(\"C: \", c, \";\\tSum of Squared weights:\\t\", np.sum(log_clf.coef_**2, axis=1))\n",
    "    \n",
    "    \n",
    "### STUDENT END ###\n",
    "P3()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ANSWER:\n",
    "\n",
    "a) In this case we have a sparse matrix. KNN doesn't perform well in this space where the KNN uses the value of features to calculate the distance between them. In this case, see P2, we have almost 25% of features in the dev data missing from the training dataset. Due to the nature of the problem, having the dev data not be represented in the training data, results in a poor fit for the KNN model. \n",
    "\n",
    "b) In this case the number of observations is too low for logistic regression to out perform naive bayes. Due to the smallish number of observations, logistic regression overfits the training data. Since Naive Bayes assumes indepedence between features it performs better. Potentially, if we had a larger sample logistic regression would likely improve. \n",
    "\n",
    "c) C is the inverse of the regularization strength. Regularization penalizes large coefficients to reduce overfitting. Small C values will increase the regularization strength, which creates simpler models that under fit the data. As C increases, the squared error also increases. Since smaller values of C correspond to smaller maginitudes of error we end up with smaller squared weights. As C increases we have larger error due to more complicated models, which increases the squared weights.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(4) Train a logistic regression model. Find the 5 features with the largest weights for each label -- 20 features in total. Create a table with 20 rows and 4 columns that shows the weight for each of these features for each of the labels. Create the table again with bigram features. Any surprising features in this table?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Top Features: alt.atheism\n",
      "['atheists' 'bobby' 'religion' 'atheism' 'deletion']\n",
      "\n",
      "Top Features: comp.graphics\n",
      "['graphics' 'image' 'file' '3d' 'computer']\n",
      "\n",
      "Top Features: sci.space\n",
      "['space' 'orbit' 'nasa' 'launch' 'spacecraft']\n",
      "\n",
      "Top Features: talk.religion.misc\n",
      "['christians' 'christian' 'blood' 'fbi' 'order']\n",
      "\n",
      "            alt.atheism  comp.graphics  sci.space  talk.religion.misc\n",
      "atheists       0.832119      -0.093456  -0.267107           -0.637606\n",
      "bobby          0.823690      -0.191303  -0.286562           -0.387156\n",
      "religion       0.798272      -0.508368  -0.658891           -0.050604\n",
      "atheism        0.793995      -0.350952  -0.367167           -0.377007\n",
      "deletion       0.723416      -0.249358  -0.272419           -0.270464\n",
      "graphics      -0.641205       1.625936  -1.103125           -0.626308\n",
      "image         -0.468683       1.117015  -0.671884           -0.381255\n",
      "file          -0.275453       1.059703  -0.682746           -0.506153\n",
      "3d            -0.303641       0.941824  -0.571095           -0.317591\n",
      "computer       0.070543       0.844575  -0.569522           -0.396109\n",
      "space         -1.061110      -1.116107   1.863778           -0.972125\n",
      "orbit         -0.350376      -0.559592   1.012597           -0.498962\n",
      "nasa          -0.464907      -0.410400   0.851649           -0.401573\n",
      "launch        -0.376461      -0.396786   0.783752           -0.280916\n",
      "spacecraft    -0.297418      -0.322026   0.737135           -0.291824\n",
      "christians    -0.607123      -0.315351  -0.404374            0.929423\n",
      "christian     -0.483107      -0.339053  -0.247399            0.925525\n",
      "blood         -0.428940      -0.096922  -0.230526            0.846944\n",
      "fbi           -0.244036      -0.217507  -0.374844            0.754721\n",
      "order         -0.649120      -0.065593  -0.124128            0.742824\n",
      "\n",
      "Using bigram features\n",
      "Top Features: alt.atheism\n",
      "['claim that' 'cheers kent' 'was just' 'you are' 'are you']\n",
      "\n",
      "Top Features: comp.graphics\n",
      "['looking for' 'in advance' 'comp graphics' 'out there' 'is there']\n",
      "\n",
      "Top Features: sci.space\n",
      "['the space' 'the moon' 'sci space' 'and such' 'it was']\n",
      "\n",
      "Top Features: talk.religion.misc\n",
      "['the fbi' 'cheers kent' 'ignorance is' 'but he' 'of jesus']\n",
      "\n",
      "               alt.atheism  comp.graphics  sci.space  talk.religion.misc\n",
      "claim that        0.605549      -0.199038  -0.274388           -0.140375\n",
      "cheers kent       0.555719      -0.697909  -0.663787            0.534813\n",
      "was just          0.482050      -0.131420  -0.128857           -0.227463\n",
      "you are           0.472741      -0.279912  -0.481446            0.028322\n",
      "are you           0.446954      -0.248255  -0.097167           -0.305621\n",
      "looking for      -0.630341       1.108368  -0.499929           -0.571839\n",
      "in advance       -0.459351       0.832570  -0.438491           -0.418467\n",
      "comp graphics    -0.292166       0.801203  -0.370902           -0.285188\n",
      "out there        -0.274802       0.758680  -0.479018           -0.277106\n",
      "is there         -0.340882       0.754981  -0.468253           -0.257073\n",
      "the space        -0.268187      -0.530085   0.871149           -0.273826\n",
      "the moon         -0.350563      -0.492717   0.831217           -0.213777\n",
      "sci space        -0.258008      -0.329148   0.621124           -0.221244\n",
      "and such         -0.204269      -0.337525   0.590562           -0.218188\n",
      "it was           -0.203030      -0.309666   0.525331           -0.313575\n",
      "the fbi          -0.131974      -0.211345  -0.295326            0.552001\n",
      "cheers kent       0.555719      -0.697909  -0.663787            0.534813\n",
      "ignorance is     -0.158011      -0.171066  -0.139583            0.504306\n",
      "but he           -0.190019      -0.216871  -0.137369            0.491867\n",
      "of jesus         -0.091910      -0.173757  -0.212350            0.423749\n"
     ]
    }
   ],
   "source": [
    "def P4():\n",
    "### STUDENT START ###\n",
    "\n",
    "    import pandas as pd\n",
    "    \n",
    "    def gram_feature(cv):\n",
    "        cv\n",
    "        feature_train = cv.fit_transform(train_data)\n",
    "        feature_dev = cv.transform(dev_data)\n",
    "\n",
    "        #C = 0.5 from P3\n",
    "        log_reg = LogisticRegression(C = 0.5)\n",
    "        log_reg.fit(feature_train, train_labels)\n",
    "        vocab=np.array(cv.get_feature_names())\n",
    "\n",
    "        inx=np.zeros(20,dtype=int)\n",
    "        r_num = range(len(newsgroups_train.target_names))\n",
    "\n",
    "        for i in r_num:\n",
    "            log_reg_coef=log_reg.coef_[i,]\n",
    "            ind1=np.array(log_reg_coef).argsort()[-5:][::-1]\n",
    "            inx[i*5:(i+1)*5]=ind1\n",
    "\n",
    "            print(\"Top Features: %s\"%newsgroups_train.target_names[i])\n",
    "            print(vocab[ind1])\n",
    "            print(\"\")\n",
    "\n",
    "        dictionary={}\n",
    "        for i in r_num:\n",
    "            dictionary[newsgroups_train.target_names[i]]=log_reg.coef_[i,inx]\n",
    "        table=pd.DataFrame(data=dictionary,index=vocab[inx])\n",
    "        print( table)\n",
    "\n",
    "    cv = CountVectorizer()\n",
    "    gram_feature(cv)\n",
    "    \n",
    "    print(\"\\nUsing bigram features\")\n",
    "    cv = CountVectorizer(ngram_range=(2, 2))\n",
    "    gram_feature(cv)\n",
    "    \n",
    "\n",
    "### STUDENT END ###\n",
    "P4()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Example of Cheer Kent 1 \n",
      "\n",
      "\n",
      "Seems he didn't understand anything about realities, liar, lunatic\n",
      "or the real thing is a very narrow view of the possibilities of Jesus\n",
      "message.\n",
      "\n",
      "Sigh, it seems religion makes your mind/brain filter out anything\n",
      "that does not fit into your personal scheme. \n",
      "\n",
      "So anyone that thinks the possibilities with Jesus is bound to the\n",
      "classical Lewis notion of 'liar, lunatic or saint' is indeed bound\n",
      "to become a Christian.\n",
      "\n",
      "Cheers,\n",
      "Kent\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example of Cheer Kent 2 \n",
      "\n",
      "\n",
      "Most likely the tragic situation in Bosnia is a combination of ethnical\n",
      "and religious motives, where religion is just one attribute that separates\n",
      "the groups from each other.\n",
      "\n",
      "But I must agree that the sad saga in Bosnia is a terrible example\n",
      "of a case where religion is not helping, instead it is used as a weapon\n",
      "against other humans. And my sympathies are mostly on the Bosnian side,\n",
      "it looks like the Serbs are the oppressors, willing to use even\n",
      "Christianity as a weapon against their former friends.\n",
      "\n",
      "Cheers,\n",
      "Kent\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example of Cheer Kent 3 \n",
      "\n",
      "\n",
      "Glad to hear this, just a note, Osiris, Mithras and many other\n",
      "cult gods resurrected as well, so there's a good chance for all of\n",
      "us to maybe end up in a virtual reality simulator, and live forever,\n",
      "hurrah!\n",
      "\n",
      "Sorry, this was a joke, some sort of one anyway. I'm the first\n",
      "that connected Osiris with a virtual reality personality database.\n",
      "Time to write a book.\n",
      "\n",
      "\n",
      "Cheers,\n",
      "Kent\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example of Cheer Kent 4 \n",
      "\n",
      "Sorry, the San Jose based Rosicrucian order is called A.M.O.R.C, \n",
      "I don't remember for the time being what the A.M. stand for\n",
      "but O.R.C is Ordo Rosae Crucis, in other words latin for\n",
      "Order of the Rose Cross. Sigh, seems l'm loosing more and more\n",
      "of my long term memory.\n",
      "\n",
      "Otherwise their headquarters in San Jose has a pretty decent\n",
      "metaphysical bookstore, if any of you are interested in such books.\n",
      "And my son loves to run around in their Egyptian museum.\n",
      "\n",
      "Cheers,\n",
      "Kent\n",
      "\n",
      " ************************************************************ \n",
      "\n",
      "Example of Cheer Kent 5 \n",
      "\n",
      "\n",
      "Well, looking at the videos it seems that this fire started in various\n",
      "places at the same time, which would indicate that this was a planned\n",
      "action. I'm sure FBI and BATF didn't *deliberately* start a possible\n",
      "fire, having a sniper kill Korresh would have been a far easier \n",
      "method. Looking at the careful operation, and use of tear gas\n",
      "that as I know don't start fires, it is less likely that this \n",
      "was the case.\n",
      "\n",
      "Sorry, but my bets are on fanatical people keen to start\n",
      "Armageddon -- theirs.\n",
      "\n",
      "Cheers,\n",
      "Kent\n",
      "\n",
      " ************************************************************ \n",
      "\n"
     ]
    }
   ],
   "source": [
    "vectorizer_bigram = CountVectorizer(ngram_range=(2,2)) \n",
    "vectorized_bigram_train = vectorizer_bigram.fit_transform(train_data)\n",
    "cheer_kent = list(vectorized_bigram_train[:,37174].nonzero()[0])\n",
    "\n",
    "\n",
    "\n",
    "for i in range(0,5):\n",
    "    print(\"Example of Cheer Kent\", i+1, '\\n')\n",
    "    print(train_data[cheer_kent[i]])\n",
    "    print('\\n', '*'*60, '\\n')\n",
    "    \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ANSWER: In the bigram, \"cheers kent\" appears in both religion and atheism. Looking at some examples, it appears a person named Kent might been active in both the religion and atheism. It appears he signed off with \"Cheers, Kent\" ever time. The model learned that \"Cheers Kent\" is a feature with the largest weight for that category. \n",
    "\n",
    "Overall the bigram has a lot of 'filler' phrases such as \"is there\", \"are you\", etc. that it assigned to high weights on one category. A possible consideration is that certain people with distinct writing styles were very active in a specific topic. As a result of the smallish amount of data, the bigram model is picking these up and assigning high weights. If we had a large training set, logistic regression might improve. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(5) Try to improve the logistic regression classifier by passing a custom preprocessor to CountVectorizer. The preprocessing function runs on the raw text, before it is split into words by the tokenizer. Your preprocessor should try to normalize the input in various ways to improve generalization. For example, try lowercasing everything, replacing sequences of numbers with a single token, removing various other non-letter characters, and shortening long words. If you're not already familiar with regular expressions for manipulating strings, see https://docs.python.org/2/library/re.html, and re.sub() in particular. With your new preprocessor, how much did you reduce the size of the dictionary?\n",
    "\n",
    "For reference, I was able to improve dev F1 by 2 points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Empty Pre-processor\n",
      "Vocab Size:  33291\n",
      "F1-score 0.7023340087555402\n",
      "\n",
      "Better Pre-Processor\n",
      "Vocab Size:  14965\n",
      "F1-score 0.7196013153193458\n"
     ]
    }
   ],
   "source": [
    "def empty_preprocessor(s):\n",
    "    return s\n",
    "\n",
    "def better_preprocessor(s):\n",
    "### STUDENT START ###\n",
    "\n",
    "    #Removes lower case, numbers, special characters, \n",
    "    s = s.lower()\n",
    "    s = re.sub(r'\\d{3,}', \" numbers \", s)\n",
    "    s = re.sub(r'\\W', \" \", s)\n",
    "    s = re.sub(r'_', \" \", s)\n",
    "\n",
    "    #Trim words larger than 6 characters\n",
    "    s1 = \"\"\n",
    "    for w in s.split():\n",
    "        if len(w) > 6:\n",
    "            w = w[:6]\n",
    "        s1 = s1 + \" \" + w\n",
    "    \n",
    "    return s1\n",
    "\n",
    "\n",
    "### STUDENT END ###\n",
    "\n",
    "def P5():\n",
    "### STUDENT START ###\n",
    "    \n",
    "    def preprocess(ppf):\n",
    "        cv = CountVectorizer(preprocessor = ppf)\n",
    "        feature_train = cv.fit_transform(train_data)\n",
    "        feature_dev = cv.transform(dev_data)\n",
    "        \n",
    "        log_reg = LogisticRegression()\n",
    "        log_reg.fit(feature_train, train_labels)\n",
    "        log_reg_pred = log_reg.predict(feature_dev)\n",
    "        f1 = f1_score(dev_labels, log_reg_pred, average=\"weighted\")\n",
    "        print(\"Vocab Size: \", feature_train.shape[1])\n",
    "        print(\"F1-score\", f1)\n",
    "\n",
    "    print(\"Empty Pre-processor\")\n",
    "    preprocess(empty_preprocessor)\n",
    "    \n",
    "    print(\"\\nBetter Pre-Processor\")\n",
    "    preprocess(better_preprocessor)\n",
    "### STUDENT END ###\n",
    "P5()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Answer: \n",
    "Vocab size decreased by 15,149\n",
    "\n",
    "F1-score increased by 0.017 or ~1.78 points.\n",
    "\n",
    "By removing special characters, numbers, case sensitvity, and long words, we were able to improve the performance of the model. Through testing, I found that 6 was the optimal trimming value that maximized performance. For example, a trimming value of 5 characters reduced the model performance by 0.06.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(6) The idea of regularization is to avoid learning very large weights (which are likely to fit the training data, but not generalize well) by adding a penalty to the total size of the learned weights. That is, logistic regression seeks the set of weights that minimizes errors in the training data AND has a small size. The default regularization, L2, computes this size as the sum of the squared weights (see P3, above). L1 regularization computes this size as the sum of the absolute values of the weights. The result is that whereas L2 regularization makes all the weights relatively small, L1 regularization drives lots of the weights to 0, effectively removing unimportant features.\n",
    "\n",
    "Train a logistic regression model using a \"l1\" penalty. Output the number of learned weights that are not equal to zero. How does this compare to the number of non-zero weights you get with \"l2\"? Now, reduce the size of the vocabulary by keeping only those features that have at least one non-zero weight and retrain a model using \"l2\".\n",
    "\n",
    "Make a plot showing accuracy of the re-trained model vs. the vocabulary size you get when pruning unused features by adjusting the C parameter.\n",
    "\n",
    "Note: The gradient descent code that trains the logistic regression model sometimes has trouble converging with extreme settings of the C parameter. Relax the convergence criteria by setting tol=.01 (the default is .0001)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "C = 0.01\n",
      "\tNumber non-zero weights w/ l1 penalty: 35\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.02\n",
      "\tNumber non-zero weights w/ l1 penalty: 79\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.03\n",
      "\tNumber non-zero weights w/ l1 penalty: 119\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.04\n",
      "\tNumber non-zero weights w/ l1 penalty: 171\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.05\n",
      "\tNumber non-zero weights w/ l1 penalty: 209\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.06\n",
      "\tNumber non-zero weights w/ l1 penalty: 227\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.07\n",
      "\tNumber non-zero weights w/ l1 penalty: 258\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.08\n",
      "\tNumber non-zero weights w/ l1 penalty: 299\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.09\n",
      "\tNumber non-zero weights w/ l1 penalty: 328\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.1\n",
      "\tNumber non-zero weights w/ l1 penalty: 359\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.25\n",
      "\tNumber non-zero weights w/ l1 penalty: 749\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 0.5\n",
      "\tNumber non-zero weights w/ l1 penalty: 1231\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 1\n",
      "\tNumber non-zero weights w/ l1 penalty: 1871\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 1.25\n",
      "\tNumber non-zero weights w/ l1 penalty: 2131\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 1.5\n",
      "\tNumber non-zero weights w/ l1 penalty: 2254\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 1.75\n",
      "\tNumber non-zero weights w/ l1 penalty: 2390\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 2\n",
      "\tNumber non-zero weights w/ l1 penalty: 2818\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 2.5\n",
      "\tNumber non-zero weights w/ l1 penalty: 2859\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 3\n",
      "\tNumber non-zero weights w/ l1 penalty: 3072\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 5\n",
      "\tNumber non-zero weights w/ l1 penalty: 4318\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 7\n",
      "\tNumber non-zero weights w/ l1 penalty: 4870\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n",
      "C = 10\n",
      "\tNumber non-zero weights w/ l1 penalty: 5463\n",
      "\tNumber non-zero weights w/ l2 penalty: 107516\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEWCAYAAACXGLsWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3Xl8nWWd///XJ0uTtmnTJWnTpvtG\nS1tooUBlVwotqICMo6A4oAijM7iMjCMoo/x0HHVmvs44MzgOIogjCshmVaAUkEWg0FYo3aELbdI2\nS5tmabOf8/n9cd8Jp+lJTtrk9CQn7+fjcR657+u+73Nf1+np/TnXcl+3uTsiIiJdyUh1BkREpO9T\nsBARkYQULEREJCEFCxERSUjBQkREElKwEBGRhBQsBAAz+ycz229mZeH6R8ysxMwOmdnCFOarT+Qj\nlpndYWa/7KX3+rqZ3d0b7yWSTAoWA4SZvWtmDeFFt+313+G2icAtwMnuXhQe8m/Aze6e5+5v9OC8\nbmYzepD1TvNhZlvM7DNxzvklM1vTg3P2KjO7wszeNLPaMCA/a2ZTANz9n939s6nNYec6+d6MD7fd\nZWZbzSxqZtcneJ8JZvZIWP4aM1uf6BjpW7JSnQE5oT7s7s/ESZ8MHHD3ig5pG09MtrrUVT7uA/4K\nuKdD+qfCbSkXBspfAFcBzwF5wCVANJX5OkadfW/WAQ8CP+jGe/xfuP9koAmYDxR1ecQxMrMsd2/t\nzfeU96hmMcCZ2RJgJTA+/NX4azM7BGQC68xse7jf+PCXYaWZ7TSzL8a8R2bYnLLdzOrMbK2ZTTSz\nF8Nd1oXv/fE4588ws9vNbJeZVZjZL8ws38xy4uWjg/8DzjWzyTHvNwc4Bfh1TL6Xm1mVmW0zsxsT\n5Tvc9qOw+as2TD+vw7lzzezB8Lg/m9mpnXzEC4Cd7v6sB+rc/RF33x2ep71Jy8z+u8Mv+FYzuyPR\n5x/nM80PP8fK8HO93cwywm3Xm9mfzOzfzOxg+F6XdvZeXXH3O939WaCxG7ufAfzc3Q+7e6u7v+Hu\nT8bk+Vwze8XMqsPP/fpuluVlM/t3M6sC2j6rz5jZ5rB8K2K/H9ID7q7XAHgB7wJLOtl2IVDaIc2B\nGeFyBrAW+CYwCJgG7ACWhtu/CqwHTgIMOBUY3fF9Ojn3Z4Bt4XvmAY8C/xcvH50cvxK4PWb9e8Dj\nMesvAD8Gcgku3JXARd3I97XAaILa9y1AGZAbbrsDaAE+CmQDfw/sBLLj5G8awcX034H3A3kdtt8B\n/DLOcW15XZjo849z7C+A3wLDgCnA28AN4bbrw7zfSBCIPw/sBexYvzcx+/wJuD7BPs8ALwNXA5M6\nbJsE1AHXhJ/naGBBN8vSCnwh/HcaDFwZfp/mhGm3A6+k+v9fOrxSngG9TtA/dPCf/hBQHfO6Mdx2\nIV0Hi7OA3R223wbcGy5vBa7o5LyJLvbPAn8Ts35SeDHL6ubx1wJbw+UMYDfwkXB9IhABhsXs/z2C\nX7hd5jvOeQ4Cp4bLdwCrYrZlAPuA8zo5djHwEMHFvxH4OWHQIE6wAArDf6+ru/P5d0jPJGjmOTkm\n7a+B58Pl64FtMduGhJ9xUTe/N4/H2ac7wWIk8H2CJsUI8CZwRkxZHjvOsnT8XJ4kDCYx/zb1wOQT\n/X8u3V5qhhpYrnT3ETGvn3bzuMkEzVTVbS/g68DYcPtEIF4zUXeMB3bFrO8i+EU4Nv7uR3kUGGdm\niwmC3hDgDzHvXeXudR3evzhc7jTfZnZL2JRRE5Y3HyiI2aWkbcHdo0BpeL6juPsqd/+YuxcC5wHn\nA9/o5LzZwMPAr9z9gTA50ecfq4Cg9tHxMy2OWS+LyVt9uJgXLz+h2O/NlV3s1yl3P+jut7r73DDf\nbwKPm5nR+b9Dd8pSwpEmAz+K+ZyqCGqNxUiPqINbuqOEoN19ZhfbpwMbjuO99xL8B28ziaBpobw7\nB7t7vZk9TNDRPRh4wN2bY957lJkNiwkYk4A9XeU77J/4GnARsNHdo2Z2kOCi02ZizP4ZwITwfIny\nu9rMHgXmdbLLfxE0ydwek5bo84+1n6BmNhnYFKbFljnl3H2/mf0bcB0wiqB8Z8bZtTtl6Thtdgnw\nXXe/v1czLapZSLe8DtSa2dfMbHDYMTzPzM4It98NfMfMZlrgFDMbHW4rJ2hj78yvgb8zs6lmlgf8\nM/CgH9uolvuAjwN/QcwoKHcvAV4BvmdmuWZ2CnAD0HYh6SzfwwgCViWQZWbfBIZ3OOfpZnaVmWUB\nXyZoLlnVMWNhx+2NZjYmXJ8NXN7Jvn8NXAB8IqyttEn0+bdz9whBk9d3zWxY2Ln7FaBX7gvpkN9B\nZpZLEESzw8847jXFzH4Q5jnLzIYR9JVsc/cDBP8eS8zsY+H20Wa24DjL8hPgNjObG54338z+svdK\nPXApWAwsv+sw2uax7hwU/qf9MOHIHoJffHcTNM0A/JDgP/XTQC3wM4Jf+RC0yd8XNgt8LM7b30Mw\nqunF8L0bCTosj8WLQA2wx91Xd9h2DUHH6F7gMeBb7r4yQb5XELR9v03Q7NHI0c0dvyUIUAcJhupe\n5e4tcfJWTRAc1lswuuupMB//EmffawgC696Yf6Ovd+Pz7+gLwGGCTvA/Ab/i6OHFveFpoAE4G7gr\nXD6/k32HEJS7OszXZILPBQ9Ghl1GMJCgiqCJqm102TGVxd0fIxjK+4CZ1RLUGo9rtJccycJOIBER\nkU6pZiEiIgkpWIiISEIKFiIikpCChYiIJJQ291kUFBT4lClTUp0NEZF+Ze3atfvDG0a7lDbBYsqU\nKaxZ02dmpRYR6RfMbFfivdQMJSIi3aBgISIiCSlYiIhIQgoWIiKSkIKFiIgkpGAhIiIJKViIiEhC\naXOfhYgMPCVV9fxxawUFeTlceFIhQwbpkpYs+mRFpN9wd7aU1bFiYxkrNpazeV9t+7bc7AwumFXI\nsnlFXDRnLMNzs1OY0/SjYCEifVo06vx590FWbCzj6U3l7DpQjxmcPmkk37hsDktOHktZTSNPbdjH\nU2EQyc40zplRwKXzirj45CJGDR2U6mL0e2nz8KNFixa5pvsQSQ/NrVFe3XGAFRvLWLmpnMq6JrIz\njbOnF7B0bhFLTh7DmGG5Rx0XjTpvlFSzYmMZT27YR0lVA5kZxllTR3HpvCIumVvE2OFHHzeQmdla\nd1+UcD8FCxHpCw43tfLC25Ws2FjGc1sqqGtsZcigTC48qZClc4t4/+wxx9S05O5s3FvLUxuCwLG9\n8jBmcNqkkVw6r4ilc4uYOGpIEkvUPyhYiEifV3W4mWc2l/P0xjJefGc/za1RRg7JZsmcsSydW8S5\nMwvIzc7slXO9U17HkxvKeHJDWXtfx/zifJbNK+LSeUVMK8zrlfP0NwoWItIn7alu4OmNZazYWMbr\nO6uIOozPz+WSucGv/TOmjCQrM7mj+ncdOBzWOMp4s6QagFlj81g2bxyXzitidtEwzCypeegrFCxE\npM/YVlHHUxuCzuf1e2oAmDkmj6VhgJhXPDxlF+e91Q1hH0cZq9+twh2mjB7SHjhOmZCf1oFDwUIk\njtZIlG88toEDh5u561Onk5FxYi8C5bWN/Oq13Vx1WjGTRw89oec+kaJR5609NeEQ1zJ2VB4GYMHE\nEWGAGNsnm30q65p4elMZT20o49XtB2iNOuPzc1k6r4hL543j9MkjyTzB35lk626w0NBZGTBaI1H+\n7qF1/G7dXgB+99ZerlhQfELOXd/cyv++sIO7XtxBQ0uExpYIt10254Sc+0RpiUR5fWdVMMR1Yzll\ntY1kZRiLp43m02dP4eKTiyjK79sjkQqH5fDJsybzybMmU13fzDObK3hqwz7uf2039778LgV5OSyd\nO5ZL543jrGmjyE5yc1lfopqFDAiRqPOVh97kt2/u5R+WncTv1+2jpqGFZ2+5oNc6UDs778NrS/h/\nT79NRV0TH5w/jnWl1cwuGsbd152RtPOeKA3NEV58JxjB9OzmCmoaWtpvjls6t4iLZo8lf0j/vznu\nUFMrz22pYMWGYKRWQ0uEEWFH/MJJI5g4cgiTRg1h/IjBDMrqXwFENQuRUCTqfPU36/jtm3v56tKT\n+JsLZ7Bgwgg+cfdr/PyVd/ncBdOTct4X367kn5/YzJayOhZOGsH/XHsap08exed/uZYtZXVJOeeJ\nUFPfwrNbylmxsYwX3q6ksSXK8NysYATTvCLOn1nI4EHJC8CpkJeTxeWnjufyU8fT2BLhhbcrwz6Y\nMh5eW9q+X4ZB0fBcJo4awsRRQQCZOGpwezApHJbTb/s/FCwkrUWjztceeYtH39jDLRfP4m/fPwOA\ns2cU8IHZY7jzuW385ekTGJ2X02vn3FpWx3ef2MyLb1cycdRg7vzEaVw2v6j9IjFjTB5PbyqnqTVC\nTlb/uKiW1zaGI5jKWbUjaMsfOzyHvzx9IkvnFg2oJpnc7Mz2jvlI1CmvbaSkqp7dVfWUHGygpKqe\nkqp6XnqnkvLapiOOzcnKCALJyMHtwWTCyPeCyrA+PEWJgoWkrWjUufXRt3h4bSlfXjKTL1w084jt\nX79sNkv/4yV+9Ow7fPuKeT0+X0VtIz9c+TYPrSkhLyeL2z84h0+9b/JRAWF6YR6RqLPrQD2zxg7r\n8XmTZef+w6zYGHT2tg0vnVYwlM+eN42lc8dy6oQRJ3yAQF+TmWGMHzGY8SMGc9a00Udtb2yJUHqw\ngZKD9e1BpKSqgd1V9azZdZC6xtYj9h8xJDsIHCOHhLWTwX2miUvBQtJSNOp8/bH1PLSmlC9+YAZf\nXjLrqH1mjBnGNWdO5P7XdnPd2VOYfpyjc+qbW/npizv53xe30xKJcv3ZU/niRTMYMST+fEQzxgTn\n2V5xqE8Fi7Y7nttGML1dfggIblz7+0tmsXRuETPG5PXbZpRUyM3OZMaYvPZ/845q6lvCGklYMwlr\nJ5v21bJyUznNkWj7vl02cY0eEnf6k96kYCFppak1wlMbyvjFq7tYu+sgf/v+6fzdxUcHijZfXjKL\nx9/Yy/ee2MLd1yXs4ztCJOo8sraUf3t6KxV1TVw6r4ivLZvNlIKuh8ROKwy2b6s4dEznS4ZI1Fn9\n7nsjmPZUN5BhcObUUXzrwydzydwiikcMTnU201b+kGzmD8ln/oT8o7ZFo055XSO7DxzZvFVysJ4/\nvbOfstrG9n3nFQ/n9184L6l5VbDow2oaWnhy/T6e2VzBvOLhXLt4MgW92LaeTnYfqOf+13fxmzWl\nVB1uZvLoIXznynlce9akLn8JF+Tl8PkLp/OvK7by6vYDvG/60U0J8bz0TiXf/UPQeb1g4gh+/MnT\nWDRlVLeOHTIoi+IRg9lemZpg0dgS4eVt+1mxsYxnNldQdbiZQVkZnD+zgC8tmcmSOWM1S2sfkJFh\njMsfzLj8wZwVZ3tjS4Q91UGT1omo6ylY9DHNrVH+uLWCx9/Yw7ObK2iORBmXn8szm8v58fPbuWph\nMZ85d2qfar5IldZIlOe2VPDL13bz4tuVZGYYS+aM4drFkzlnekG329NvOHcq96/axXef2MTyvz23\ny+PeLq/jn5/YzPNbK5kwcjD/dc1CPnTKuGNumplWOJTt4Y1qJ0JdYwvPbang6Y3lPL+1gsPNEYbl\nZPGBOWNYOreIC2YVMjRHl4P+JDc7k+mFecfdfHqs9O3oA9ydtbsO8tgbe/jD+n1U17dQkDeITy6e\nxEcWFjO/OJ/tlYe59+WdPPLnUh5YXcL5swr57LlTOW9mwYBrQy6vbeSB10t4YPVu9tU0UjQ8ly8v\nmcnVZ0w6rpu+crMz+eqyk/i7B9fx+Jt7uOq0CUftU1HXyL+vfIcHV+9maE4W37hsDn919tGd1901\nvTCPh9aUEI160jqJK+uaWLkpGOL6yvb9tEScgrwcLl9QzNK5Yzl7ekG/uydAUkc35aXQ9spD/PaN\nPTz25h5KqhrIzc5g6dwirlxYzHkzCuJOplZ1uJlfvbaL+17dRWVdE7PG5nHDuVO5YkFxUm8uS7Vo\n1Hll+wHuf20XT28qJxJ1zptZwLWLJ3PR7DE9nnguGnWu/PHLVNY18dwtF7bfJ9DQHOGnL+3gJy9s\np7k1yqfeN5kvfmAmI3vYTPPLVbu4/fENvHLrBxjfi30CJVX17R3Ua3YdxB0mjRrC0rnBLK4LJ6Xf\ndBXSM5obqo/af6iJ363by+Nv7GFdaQ0ZBufMKODKBcUsnVdEXjebAppaI/x+3T5+9qedbNpXy+ih\ng7h28WSuXTyZwmHp069x8HAzj/y5lPtf283O/YcZOSSbjy2ayDVnTkrYkXysVu04wNV3reKrS0/i\ncxdM59E/B53X5bVNLJtbxNcunc3UXjrnq9sPcM1PV/F/N5zJeTMLj/t9OnvM6Jxxw9sDxECaQVWO\nnYJFH9LQHOHpTWU8/sYeXnxnP5Goc/K44Vx1WjEfPnV8j57c5e68uuMAP3tpJ89uqWBQVgZXLhjP\nDedO46Si/tWv0dgSYWtZHW/tqWF9aTXr99Tydnkdkahz+uSRXLt4EpfOG5fUGtSNv1jDK9v2M3n0\nUDbtq+XUiSO4/YNzOKObndfdVVnXxBnffYZvffhkPn3O1GM6NvYxoys2lgcdnAaLJo9k6dwiLjm5\niEmj9VAf6R5N95Fikajz6vYDPPbGHp7asI/DzRHG5+dy0/nTuHJBca9dyM2CR02ePb2A7ZWHuPfl\nnTy8tpSH1pRy3swCbjh3KhfMKuxzvyybW6NhYKhmw54a3iqtYWtZHa3R4MfLyCHZzJ8wgotmj+GD\np4xjzrjhJyRft106m2X/8RI1DS385zUL+dD8cUnpUyjIG8Tw3Cw276uloq6R5tYoza1RmsK/zZG2\n9Uh7elNrlDd2V7NyUzn7DzW1P2f68xdOZ8mcsWlVo5S+RzWLXuTubNpXy+Nv7OG3b+6loq6JYTlZ\nXDZ/HFcuLOasqaNOyB2vBw8386vXd3PfK+9SUdfEzDF5fObcqXxkYWr6NVoiQWBYv6cmeIWBoe2G\no/zB2ZwyIZ95xfmcUhyMOS8eMThlAW5fTQMjhwxK+mf1F//zCmt3HTymY4YOyuTC2cEIpvefVNin\np4eQ/kHNUCnwmZ+v5rktFWRlGBeeNIarTivmA7PHpKzjubk1yu/f2svdLwX9GqPCfo1PJbFfozUS\n5Z2KQ6wvreGtPdWsL61hc1kdza1BYBiWm8X8MCCcUjyC+cX5TByVusCQShv31rBqRxU5WRkMysoI\n/ma2LWcyKEx/Ly2DMcNz+s18UtI/9IlgYWbLgB8BmcDd7v79OPt8DLgDcGCdu38iTI8A68Pddrv7\n5V2dK9XBork1yux/fJIPnTKeOy6f26duanJ3Vu2o4md/2sGzWyrIzsjgigXjueG8qcwuOv7mndZI\nlO2Vh3mrtLq91rBpby1NYWDIy8liXvFwTpkwor3WMGnUkAE/n5BIX5LyPgszywTuBC4GSoHVZrbc\n3TfF7DMTuA04x90PmtmYmLdocPcFycpfb9tb3UDU4fxZhX0qUEDQr/G+6aN53/TR7Kg8xL0vv8vD\na0v5zdpSzp1RwA3nTeWCmYVdXsQjUWdH5SHeKq1pDwwb99bQ2BIEhqGDMplbnM+nFk9m/oR85hfn\nM2X0UAUGkTSRzA7uM4Ft7r4DwMweAK4ANsXscyNwp7sfBHD3iiTmJ6l2V9UDwZj2vmxaYR7fuXIe\nt1wyq71f49P3rmZ64VBuOHcaV51WzKDMDHbsP9ze8bx+TzUb99ZS3xwBYHB2JvOKh3PNmZM4ZUI+\n84tHMLVgqMbvi6SxZAaLYqAkZr0UjpriZBaAmb1M0FR1h7s/FW7LNbM1QCvwfXd/PIl57bH+Eiza\njBgyiL+5cAafPXcaf1i/l5/9aSdff2w9P3hqC5Goc6gpmDo5NzuDk8cN52OLJrb3NUwvzFNgEBlg\nkhks4l1NOnaQZAEzgQuBCcBLZjbP3auBSe6+18ymAc+Z2Xp3337ECcxuAm4CmDRpUm/n/5jsrqpn\nUFYGY/rZ8MVBWRl8ZOEErlxQzGs7q3hodQlDc7KCDugJ+cwozOvx3dEi0v8lM1iUAhNj1icAe+Ps\ns8rdW4CdZraVIHisdve9AO6+w8yeBxYCRwQLd78LuAuCDu5kFKK7dh+o79edt2bG4mmjWRznAS4i\nIsn8ybgamGlmU81sEHA1sLzDPo8D7wcwswKCZqkdZjbSzHJi0s/hyL6OPmd3VX2/aYISETlWSQsW\n7t4K3AysADYDD7n7RjP7tpm1DYNdARwws03AH4GvuvsBYA6wxszWhenfjx1F1de4OyUKFiKSxpI6\n3Ye7PwE80SHtmzHLDnwlfMXu8wowP5l5603V9S3UNbUyUcFCRNKUei57QX8bCSUicqwULHpBW7CY\nrJk+RSRNKVj0grZgMXGkgoWIpCcFi16w+0A9BXk57U9XExFJNwoWvaDkYD2TRvXeozFFRPoaBYte\nUHKwXiOhRCStKVj0UGskyt7qRvVXiEhaU7DooX01jUSizkQ1Q4lIGlOw6KGSg8FIqAmqWYhIGlOw\n6KHSgw2Ahs2KSHpTsOih0qp6MgzGjchNdVZERJJGwaKHSg42MC5/MNl65oOIpDFd4XqopKpendsi\nkvYULHqo5GC9+itEJO0pWPRAY0uE8tomjYQSkbSnYNEDe6vDkVBqhhKRNKdg0QNtw2ZVsxCRdKdg\n0QNltY0AjMvXsFkRSW8KFj1QWdcEQOGwnBTnREQkuRQseqCitpH8wdnkZus5FiKS3hQseqC8tokx\nqlWIyACgYNEDFXWNjBmuYCEi6U/BogfKa5sYO0yd2yKS/hQsjpO7U1nXRKFqFiIyAChYHKfq+haa\nI1HVLERkQFCwOE4V4bBZ9VmIyECgYHGcysMb8sYOV81CRNKfgsVxaq9ZaOisiAwAChbHqaIuqFmM\nUZ+FiAwAChbHqaK2iWG5WQwepLu3RST9KVgcp4q6RjVBiciAkdRgYWbLzGyrmW0zs1s72edjZrbJ\nzDaa2a9i0q8zs3fC13XJzOfxKK9tUue2iAwYWcl6YzPLBO4ELgZKgdVmttzdN8XsMxO4DTjH3Q+a\n2ZgwfRTwLWAR4MDa8NiDycrvsaqoa+T0SSNTnQ0RkRMimTWLM4Ft7r7D3ZuBB4ArOuxzI3BnWxBw\n94owfSmw0t2rwm0rgWVJzOsxcfdgEkHVLERkgEhmsCgGSmLWS8O0WLOAWWb2spmtMrNlx3AsZnaT\nma0xszWVlZW9mPX43J1tFYeobWiluTWqPgsRGTCSGSwsTpp3WM8CZgIXAtcAd5vZiG4ei7vf5e6L\n3H1RYWFhD7Ob2C9e3cWSH77A797aC6CahYgMGMkMFqXAxJj1CcDeOPv81t1b3H0nsJUgeHTn2BOq\ntrGFHz37DgD/9Vzwd6xqFiIyQCQzWKwGZprZVDMbBFwNLO+wz+PA+wHMrICgWWoHsAK4xMxGmtlI\n4JIwLWWe21xB1eFmzpw6ivLatnmhVLMQkYEhacHC3VuBmwku8puBh9x9o5l928wuD3dbARwws03A\nH4GvuvsBd68CvkMQcFYD3w7TUmZH5SEyDP7pynntaeqzEJGBImlDZwHc/QngiQ5p34xZduAr4avj\nsfcA9yQzf8di+/7DTBw1hFljh3HezALWlVQzNCepH5+ISJ+hq1037aw8zNSCoQB876r57K6qT3GO\nREROHAWLbnB3du4/zFnTRgEwYeQQJowckuJciYicOJobqhvKahtpaIkwrTAv1VkREUkJBYtu2Fl5\nGIBpYTOUiMhAo2DRDTv2B8FiqoKFiAxQChbdsKPyMIOzMynSfRUiMkApWHTDzv2HmFIwlIyMeLOQ\niIikv4TBwsyGmNk/mtlPw/WZZvah5Get79i5/7D6K0RkQOtOzeJeoAl4X7heCvxT0nLUxzS3Rik5\n2MC0QgULERm4uhMsprv7vwAtAO7eQPxZYdPS7qp6IlFX57aIDGjdCRbNZjaYcIpwM5tOUNMYEHZq\nJJSISLfu4P4W8BQw0czuB84Brk9mpvqS7ZWHAJhWoBvyRGTg6jJYmJkBW4CrgMUEzU9fcvf9JyBv\nfcLWsjqKhueSPyQ71VkREUmZLoOFu7uZPe7upwN/OEF56lM276tl9rhhqc6GiEhKdafPYpWZnZH0\nnPRBza1RtlceYnbR8FRnRUQkpbrTZ/F+4K/NbBdwmKApyt39lKTmrA/Ysf8QLRFnjmoWIjLAdSdY\nXJr0XPRRW/bVAahmISIDXsJmKHffBYwAPhy+RoRpaW9zWS2DMjN0Q56IDHjdme7jS8D9wJjw9Usz\n+0KyM9YXbNlXx4wxeWRnagotERnYutMMdQNwlrsfBjCzHwCvAv+VzIz1BVvKajlnRkGqsyEiknLd\n+clsQCRmPcIAmO6j6nAz5bVNzFF/hYhIt2oW9wKvmdlj4fqVwM+Sl6W+YUtZLYDusRARoRvBwt1/\naGbPA+cS1Cg+7e5vJDtjqaaRUCIi70kYLMxsMbDR3f8crg8zs7Pc/bWk5y6FNu+rpSBvEIXDclKd\nFRGRlOtOn8X/AIdi1g+HaWltS1mdahUiIqFudXC7u7etuHuU7vV19FutkShvl9cxu0j9FSIi0L1g\nscPMvmhm2eHrS8COZGcsld49UE9Ta5TZ41SzEBGB7gWLzwFnA3vC11nATcnMVKptLWvr3FbNQkQE\nujcaqgK4+gTkpc+oqGsEYPyIwSnOiYhI39BpzcLMbjSzmeGymdk9ZlZjZm+Z2WknLosnXnV9CwDD\nc9O6a0ZEpNu6aob6EvBuuHwNcCowDfgK8KPkZiu1ahpaGJabRZbmhBIRAboOFq3u3hIufwj4hbsf\ncPdngG5Nw2pmy8xsq5ltM7Nb42y/3swqzezN8PXZmG2RmPTlx1Konqqub2aEHqMqItKuq3aWqJmN\nAw4CFwHfjdmWsDHfzDKBO4GLgVJgtZktd/dNHXZ90N1vjvMWDe6+INF5kqG6oYWRQwal4tQiIn1S\nVzWLbwJrCJqilrv7RgAzu4DuDZ09E9jm7jvcvRl4ALiiZ9k9MarrW8gfrJqFiEibToOFu/8emAzM\ncfcbYzatAT7ejfcuBkpi1kvDtI7+Iuw0f9jMJsak55rZGjNbZWZXxjuBmd0U7rOmsrKyG1nqnqAZ\nSjULEZE2Xfbgunurux/skHYgGYu7AAAPl0lEQVTY3Q91dkyMeNOYe4f13wFTwud5PwPcF7Ntkrsv\nAj4B/IeZTY+Tv7vcfZG7LyosLOxGlrqnuqGFEapZiIi0S+Zwn1IgtqYwAdgbu0PYYd4Urv4UOD1m\n297w7w7geWBhEvPaLhp1ahpa1MEtIhIjmcFiNTDTzKaa2SCCG/uOGNUUdqC3uRzYHKaPNLOccLkA\nOAfo2DGeFHWNrbijZigRkRjHddeZmc129y1d7ePurWZ2M7ACyATucfeNZvZtYI27Lwe+aGaXA61A\nFXB9ePgc4H/NLEoQ0L4fZxRVUlQ3NAOoGUpEJMbx3qL8NDAp0U7u/gTwRIe0b8Ys3wbcFue4V4D5\nx5m3Hmm7e1vNUCIi7+k0WJjZf3a2CRiRnOyk3sH6sGahYCEi0q6rmsWngVuApjjbrklOdlKvpiGo\nWeQPVp+FiEibroLFamBD2CR0BDO7I2k5SrG2ZqiRqlmIiLTrKlh8FGiMt8HdpyYnO6nXFix0B7eI\nyHu6Gjqb5+71JywnfUR1QzPDcjTjrIhIrK6uiI+3LZjZIycgL31CTX0L+WqCEhE5QlfBIna6jmnJ\nzkhfUa27t0VEjtJVsPBOltPawfpmTU8uItJBVx3cp5pZLUENY3C4TLju7j486blLgZr6For17G0R\nkSN0GizcPfNEZqSvUDOUiMjRNOQnRjTqwbMsdEOeiMgRFCxiHGpuJeqa6kNEpCMFixg17ZMIqmYh\nIhJLwSJG+ySCuntbROQIChYxND25iEh8ChYxqhsULERE4lGwiFETNkNpenIRkSMpWMRQM5SISHwK\nFjGqG1rIy8kiWzPOiogcQVfFGDUNLQzPPd7HkouIpC8FixgNLREGDxqQs5yIiHRJwSJGY3OE3GwF\nCxGRjhQsYjS2RhisYCEichQFixgNzWqGEhGJR8EiRmNLlJwsBQsRkY4ULGI0qoNbRCQuBYsYDS0R\ncrP0kYiIdKQrYwzVLERE4lOwiNHQotFQIiLxKFiE3D3o4FawEBE5ioJFqKk1CqCahYhIHEkNFma2\nzMy2mtk2M7s1zvbrzazSzN4MX5+N2Xadmb0Tvq5LZj4huMcCYHC24qeISEdJmzXPzDKBO4GLgVJg\ntZktd/dNHXZ90N1v7nDsKOBbwCLAgbXhsQeTld+GliBYaLoPEZGjJfNn9JnANnff4e7NwAPAFd08\ndimw0t2rwgCxEliWpHwCwUgoQKOhRETiSGawKAZKYtZLw7SO/sLM3jKzh81s4rEca2Y3mdkaM1tT\nWVnZo8y21Sx0B7eIyNGSGSwsTpp3WP8dMMXdTwGeAe47hmNx97vcfZG7LyosLOxRZhtbwg5u1SxE\nRI6SzGBRCkyMWZ8A7I3dwd0PuHtTuPpT4PTuHtvb2puh1GchInKUZAaL1cBMM5tqZoOAq4HlsTuY\n2biY1cuBzeHyCuASMxtpZiOBS8K0pGkbDZWr0VAiIkdJ2mgod281s5sJLvKZwD3uvtHMvg2scffl\nwBfN7HKgFagCrg+PrTKz7xAEHIBvu3tVsvIKwbMsQDULEZF4kvrAaXd/AniiQ9o3Y5ZvA27r5Nh7\ngHuSmb9Y79UsFCxERDpSm0uoUfdZiIh0SsEipNFQIiKdU7AItd/BredZiIgcRVfGUENLhOxMIytT\nH4mISEe6MoYaWyLqrxAR6YSCRUjBQkSkcwoWocaWqO6xEBHphIJFqKFZj1QVEemMgkWooSWiqT5E\nRDqhq2NIfRYiIp1TsAg1tkR0Q56ISCcULEINLRFy9eAjEZG4FCxCjS1R1SxERDqhYBFSB7eISOd0\ndQw1NquDW0SkMwoWocZW3WchItIZBQugJRKlJeKqWYiIdELBgvcefKSahYhIfAoWvPfgo1yNhhIR\niUvBgphHqurBRyIicenqSEwzlGoWIiJxKVgQ+0hVBQsRkXgULAimJwfVLEREOqNgATS2hh3cGg0l\nIhKXggXv1Sw03YeISHy6OqL7LEREElGwQKOhREQSUbBAo6FERBJRsOC9YKGahYhIfAoWvDfdR47u\n4BYRiSupV0czW2ZmW81sm5nd2sV+HzUzN7NF4foUM2swszfD10+Smc/G8MFHZpbM04iI9FtZyXpj\nM8sE7gQuBkqB1Wa23N03ddhvGPBF4LUOb7Hd3RckK3+xGlv0LAsRka4ks2ZxJrDN3Xe4ezPwAHBF\nnP2+A/wL0JjEvHSpQU/JExHpUjKDRTFQErNeGqa1M7OFwER3/32c46ea2Rtm9oKZnRfvBGZ2k5mt\nMbM1lZWVx53RBtUsRES6lMxgEa8DwNs3mmUA/w7cEme/fcAkd18IfAX4lZkNP+rN3O9y90Xuvqiw\nsPC4M9rYElXNQkSkC8kMFqXAxJj1CcDemPVhwDzgeTN7F1gMLDezRe7e5O4HANx9LbAdmJWsjLZ1\ncIuISHzJvEKuBmaa2VQzGwRcDSxv2+juNe5e4O5T3H0KsAq43N3XmFlh2EGOmU0DZgI7kpXRhpaI\n7rEQEelC0oKFu7cCNwMrgM3AQ+6+0cy+bWaXJzj8fOAtM1sHPAx8zt2rkpXXxpaI7t4WEelC0obO\nArj7E8ATHdK+2cm+F8YsPwI8ksy8xWpoiej52yIiXVBDPdDYrNFQIiJdUbAgePiROrhFRDqnKyTB\nTXmqWYiIdG7ABwt3p7FVwUJEpCsDPlg0tUZxhxwFCxGRTg34YKFHqoqIJDbgg4WZ8cFTxjF9TF6q\nsyIi0mcl9T6L/iB/cDZ3fuK0VGdDRKRPG/A1CxERSUzBQkREElKwEBGRhBQsREQkIQULERFJSMFC\nREQSUrAQEZGEFCxERCQhc/dU56FXmFklsOs4Di0A9vdydvoSla9/U/n6t/5QvsnuXphop7QJFsfL\nzNa4+6JU5yNZVL7+TeXr39KpfGqGEhGRhBQsREQkIQULuCvVGUgyla9/U/n6t7Qp34DvsxARkcRU\nsxARkYQULEREJKEBGyzMbJmZbTWzbWZ2a6rz011mdo+ZVZjZhpi0UWa20szeCf+ODNPNzP4zLONb\nZnZazDHXhfu/Y2bXpaIs8ZjZRDP7o5ltNrONZvalMD0tymhmuWb2upmtC8v3/4XpU83stTCvD5rZ\noDA9J1zfFm6fEvNet4XpW81saWpKFJ+ZZZrZG2b2+3A9bcpnZu+a2Xoze9PM1oRpafH97JK7D7gX\nkAlsB6YBg4B1wMmpzlc3834+cBqwISbtX4Bbw+VbgR+Ey5cBTwIGLAZeC9NHATvCvyPD5ZGpLluY\nt3HAaeHyMOBt4OR0KWOYz7xwORt4Lcz3Q8DVYfpPgM+Hy38D/CRcvhp4MFw+Ofze5gBTw+9zZqrL\nF1POrwC/An4frqdN+YB3gYIOaWnx/ezqNVBrFmcC29x9h7s3Aw8AV6Q4T93i7i8CVR2SrwDuC5fv\nA66MSf+FB1YBI8xsHLAUWOnuVe5+EFgJLEt+7hNz933u/udwuQ7YDBSTJmUM83koXM0OXw58AHg4\nTO9YvrZyPwxcZGYWpj/g7k3uvhPYRvC9TjkzmwB8ELg7XDfSqHydSIvvZ1cGarAoBkpi1kvDtP5q\nrLvvg+BiC4wJ0zsrZ78of9gksZDg13falDFsonkTqCC4SGwHqt29NdwlNq/t5Qi31wCj6cPlA/4D\n+AcgGq6PJr3K58DTZrbWzG4K09Lm+9mZrFRnIEUsTlo6jiHurJx9vvxmlgc8AnzZ3WuDH5vxd42T\n1qfL6O4RYIGZjQAeA+bE2y3826/KZ2YfAircfa2ZXdiWHGfXflm+0DnuvtfMxgArzWxLF/v2x/LF\nNVBrFqXAxJj1CcDeFOWlN5SHVVvCvxVhemfl7NPlN7NsgkBxv7s/GianVRkB3L0aeJ6gLXuEmbX9\neIvNa3s5wu35BM2QfbV85wCXm9m7BM27HyCoaaRL+XD3veHfCoJgfyZp+P3saKAGi9XAzHCExiCC\njrXlKc5TTywH2kZTXAf8Nib9r8IRGYuBmrCKvAK4xMxGhqM2LgnTUi5sr/4ZsNndfxizKS3KaGaF\nYY0CMxsMLCHol/kj8NFwt47layv3R4HnPOghXQ5cHY4mmgrMBF4/MaXonLvf5u4T3H0Kwf+r59z9\nk6RJ+cxsqJkNa1sm+F5tIE2+n11KdQ97ql4EoxTeJmgv/kaq83MM+f41sA9oIfh1cgNBG++zwDvh\n31HhvgbcGZZxPbAo5n0+Q9BpuA34dKrLFZOvcwmq428Bb4avy9KljMApwBth+TYA3wzTpxFcDLcB\nvwFywvTccH1buH1azHt9Iyz3VuDSVJctTlkv5L3RUGlRvrAc68LXxrZrR7p8P7t6aboPERFJaKA2\nQ4mIyDFQsBARkYQULEREJCEFCxERSUjBQkREElKwkAHLzJ7vOJupmX3ZzH7ci+f4uZl9NME+GeHM\npBvC2UxXh/cWYGZPtN2XIZJKA3W6DxEI7lm5miNvhroa+OoJzsfHgfHAKe4eDSfiOwzg7ped4LyI\nxKWahQxkDwMfMrMcaJ+4cDzwp/CO23+N+bX/8baDzOwfwrR1Zvb9MO3GsEawzsweMbMhMedZYmYv\nmdnb4dxJHY0D9rl7FMDdSz2YibTt2QkFZvY5C56f8KaZ7TSzP4bbLzGzV83sz2b2m3BOLZFep2Ah\nA5a7HyC4a7htaui25yk4cBWwADiVYEqOfzWzcWZ2KcH002e5+6kEzzEAeNTdzwjTNhPcWd9mCnAB\nwbTdPzGz3A5ZeQj4cBgI/p+ZLYyT15+4+wLgDII7939oZgXA7cASdz8NWEPwHAmRXqdgIQNdW1MU\n4d9fh8vnAr9294i7lwMvEFyolwD3uns9gLu3PVtkXlh7WA98Epgbc46H3D3q7u8QPORmdmwG3L0U\nOAm4jWBa72fN7KJO8vsjgvmTfkcwAeHJwMvhlOfXAZOP50MQSUR9FjLQPU7wK/00YLCHD14i/hTS\nbenx5sj5OXClu68zs+sJ5kVq03H/o4539yaCJ6o9aWblBLWXZ484cfC+k4GbY/Ky0t2v6SSvIr1G\nNQsZ0Dx4at3zwD28V6sAeBH4ePigokKCx9m+DjwNfKatT8LMRoX7DwP2hdOrf7LDaf4yHPE0nWAi\nuq2xG83sNDMbHy5nEEw2uKvDPqcDfw9c29a3AawCzjGzGeE+Q8xs1vF9EiJdU81CJAgSj/JecxQE\nzyl4H8Hsog78g7uXAU+Z2QJgjZk1A08AXwf+keCJfrsIZhcdFvNeWwmascYCn3P3xg7nHwP8tK2j\nnSAo/XeHfW4meF7zH4NZ3Fnj7p8Naxu/jjn2doLZlEV6lWadFRGRhNQMJSIiCSlYiIhIQgoWIiKS\nkIKFiIgkpGAhIiIJKViIiEhCChYiIpLQ/w+fqGX3BVJ69AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def P6():\n",
    "    # Keep this random seed here to make comparison easier.\n",
    "    np.random.seed(0)\n",
    "\n",
    "    ### STUDENT START ###\n",
    "    cv=CountVectorizer()\n",
    "    feature_train=cv.fit_transform(train_data)\n",
    "    feature_dev=cv.transform(dev_data)\n",
    "    feature_names = cv.get_feature_names()\n",
    "    C_list=[0.01, 0.02,0.03,0.04, .05,0.06,0.07,0.08,0.09, 0.1, 0.25, 0.5, 1, 1.25, 1.5, 1.75, 2, 2.5, 3, 5, 7, 10]\n",
    "    \n",
    "    \n",
    "    vocab_size = []\n",
    "    f1_scores = []  \n",
    "    \n",
    "    for c in C_list:\n",
    "        print('C =', c)\n",
    "        \n",
    "        def log_reg(lnum):\n",
    "            \n",
    "            log_reg = LogisticRegression(penalty=lnum, C=c, tol=0.01)\n",
    "            log_reg.fit(feature_train, train_labels)\n",
    "            nonzero_weights = np.argwhere(log_reg.coef_!=0)\n",
    "            print('\\tNumber non-zero weights w/', lnum, 'penalty:', nonzero_weights.shape[0])\n",
    "        \n",
    "            if 'l1' == lnum:\n",
    "                vocab_size.append(nonzero_weights.shape[0])\n",
    "                vocab_idx = np.unique(nonzero_weights[:,-1])\n",
    "                \n",
    "                new_vocab = []\n",
    "                for i in vocab_idx:\n",
    "                    new_vocab.append(feature_names[i])\n",
    "                    \n",
    "                cv2 = CountVectorizer(vocabulary=new_vocab)\n",
    "                feature_train2=cv2.fit_transform(train_data)\n",
    "                feature_dev2=cv2.transform(dev_data)\n",
    "                log_reg_l2 = LogisticRegression(penalty='l2', C=c, tol=0.01)\n",
    "                log_reg_l2.fit(feature_train2, train_labels)\n",
    "                pred = log_reg_l2.predict(feature_dev2)\n",
    "                f1_scores.append(metrics.f1_score(pred, dev_labels, average='macro'))\n",
    "\n",
    "                \n",
    "        log_reg('l1')\n",
    "        log_reg('l2')\n",
    "\n",
    "    plt.plot(vocab_size, f1_scores)\n",
    "    plt.title('Effect of Vocab Size on F1 Score')  \n",
    "    plt.ylabel('F1 Score') \n",
    "    plt.xlabel('Vocab Size') \n",
    "    plt.show()\n",
    "\n",
    "        \n",
    "    ### STUDENT END ###\n",
    "P6()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ANSWER: \n",
    "\n",
    "L1 results in significntly fewer non-zero weights compared to L2. The plot shown above highlights that the performance based on the F1 increases significantly as the vocb size increases from about 0 to ~500. Once the vocab size reaches that threshold the performancr don't progress much.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(7) Use the TfidfVectorizer -- how is this different from the CountVectorizer? Train a logistic regression model with C=100.\n",
    "\n",
    "Make predictions on the dev data and show the top 3 documents where the ratio R is largest, where R is:\n",
    "\n",
    "maximum predicted probability / predicted probability of the correct label\n",
    "\n",
    "What kinds of mistakes is the model making? Suggest a way to address one particular issue that you see."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "F1 Score:  0.7473913101645842 \n",
      "\n",
      "************************************************************ \n",
      "\n",
      "Ratio R = 287.1791925998367 \n",
      "\n",
      "Predicted Category:  talk.religion.misc\n",
      "Correct Category:  alt.atheism\n",
      "\n",
      "Document Text: \n",
      "\n",
      "The 24 children were, of course, killed by a lone gunman in a second story\n",
      "window, who fired eight bullets in the space of two seconds...\n",
      " \n",
      "\n",
      "************************************************************ \n",
      "\n",
      "Ratio R = 325.00406840436005 \n",
      "\n",
      "Predicted Category:  comp.graphics\n",
      "Correct Category:  talk.religion.misc\n",
      "\n",
      "Document Text: \n",
      "Can anyone provide me a ftp site where I can obtain a online version\n",
      "of the Book of Mormon. Please email the internet address if possible. \n",
      "\n",
      "************************************************************ \n",
      "\n",
      "Ratio R = 929.3580586193167 \n",
      "\n",
      "Predicted Category:  comp.graphics\n",
      "Correct Category:  talk.religion.misc\n",
      "\n",
      "Document Text: \n",
      "I am pleased to announce that a *revised version* of _The Easy-to-Read Book\n",
      "of Mormon_ (former title: _Mormon's Book_) by Lynn Matthews Anderson is now\n",
      "available through anonymous ftp (see information below). In addition to the\n",
      "change in title, the revised ETR BOM has been shortened by several pages\n",
      "(eliminating many extraneous \"that's\" and \"of's\"), and many (minor) errors\n",
      "have been corrected. This release includes a simplified Joseph Smith Story,\n",
      "testimonies of the three and eight witnesses, and a \"Words-to-Know\"\n",
      "glossary.\n",
      "\n",
      "As with the previous announcement, readers are reminded that this is a\n",
      "not-for-profit endeavor. This is a copyrighted work, but people are welcome\n",
      "to make *verbatim* copies for personal use. People can recuperate the\n",
      "actual costs of printing (paper, copy center charges), but may not charge\n",
      "anything for their time in making copies, or in any way realize a profit\n",
      "from the use of this book. See the permissions notice in the book itself\n",
      "for the precise terms.\n",
      "\n",
      "Negotiations are currently underway with a Mormon publisher vis-a-vis the\n",
      "printing and distribution of bound books. (Sorry, I'm out of the wire-bound\n",
      "\"first editions.\") I will make another announcement about the availability\n",
      "of printed copies once everything has been worked out.\n",
      "\n",
      "FTP information: connect via anonymous ftp to carnot.itc.cmu.edu, then \"cd\n",
      "pub\" (you won't see anything at all until you do).\n",
      "\n",
      "\"The Easy-to-Read Book of Mormon\" is currently available in postscript and\n",
      "RTF (rich text format). (ASCII, LaTeX, and other versions can be made\n",
      "available; contact dba@andrew.cmu.edu for details.) You should be able to\n",
      "print the postscript file on any postscript printer (such as an Apple\n",
      "Laserwriter); let dba know if you have any difficulties. (The postscript in\n",
      "the last release had problems on some printers; this time it should work\n",
      "better.) RTF is a standard document interchange format that can be read in\n",
      "by a number of word processors, including Microsoft Word for both the\n",
      "Macintosh and Windows. If you don't have a postscript printer, you may be\n",
      "able to use the RTF file to print out a copy of the book.\n",
      "\n",
      "-r--r--r--  1 dba                   1984742 Apr 27 13:12 etrbom.ps\n",
      "-r--r--r--  1 dba                   1209071 Apr 27 13:13 etrbom.rtf\n",
      "\n",
      "For more information about how this project came about, please refer to my\n",
      "article in the current issue of _Sunstone_, entitled \"Delighting in\n",
      "Plainness: Issues Surrounding a Simple Modern English Book of Mormon.\"\n",
      "\n",
      "Send all inquiries and comments to:\n",
      "\n",
      "    Lynn Matthews Anderson\n",
      "    5806 Hampton Street\n",
      "    Pittsburgh, PA 15206 \n",
      "\n",
      "************************************************************ \n",
      "\n",
      "Confusion Matrix\n",
      "['alt.atheism', 'comp.graphics', 'sci.space', 'talk.religion.misc']\n",
      "[[102   5  10  31]\n",
      " [ 11 168  20  10]\n",
      " [ 19  10 165   5]\n",
      " [ 33   2   4  81]]\n"
     ]
    }
   ],
   "source": [
    "def P7():\n",
    "### STUDENT START ###\n",
    "\n",
    "    cv = TfidfVectorizer()\n",
    "    features_train = cv.fit_transform(train_data)\n",
    "    features_dev =  cv.transform(dev_data)\n",
    "    log_reg = LogisticRegression(C=100)\n",
    "    log_reg.fit(features_train, train_labels)\n",
    "    log_reg_pred = log_reg.predict(features_dev)\n",
    "    \n",
    "    print(\"F1 Score: \", metrics.f1_score(log_reg_pred, dev_labels, average = 'macro'), '\\n')\n",
    "\n",
    "    prob = log_reg.predict_proba(features_dev)\n",
    "    max_prob = np.amax(prob, axis=1)\n",
    "  \n",
    "    R = []\n",
    "    for i in range(prob.shape[0]):\n",
    "        R.append(prob[i, dev_labels[i]])\n",
    "     \n",
    "    r_ratio = max_prob / R\n",
    "    \n",
    "    top_3_r_idx = np.argpartition(r_ratio, -3)[-3:]\n",
    "    top_3_r = list(np.take(r_ratio, top_3_r_idx))\n",
    "    \n",
    "    for r in top_3_r_idx:\n",
    "        print('*'*60, '\\n')\n",
    "        print(\"Ratio R =\", top_3_r.pop(0), '\\n')\n",
    "        print(\"Predicted Category: \", newsgroups_train.target_names[log_reg_pred[r]])\n",
    "        print(\"Correct Category: \", newsgroups_train.target_names[dev_labels[r]])\n",
    "        print(\"\\nDocument Text: \")\n",
    "        print(dev_data[r], '\\n')\n",
    "\n",
    "        \n",
    "\n",
    "    cm = confusion_matrix(log_reg_pred, dev_labels)\n",
    "    print('*'*60, '\\n')\n",
    "    print(\"Confusion Matrix\")\n",
    "    print(newsgroups_train.target_names)\n",
    "    print(cm)\n",
    "## STUDENT END ###\n",
    "P7()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ANSWER: TF-IDF goes beyond taking a simple count of the frequency of each word in the document. It also takes into account the frequency of the word in the corpus of documents. So it will weight the word higher if it appears frequently in the document, but infrequently in all other documents. \n",
    "\n",
    "According to the confusion matrix, the most common misclassifications were between religion and atheism. This makes sense considering the similar types of words between those two categories. The two largest R ratios were associated with religious topics that got confused with comp.graphics. Checking the documents, but mention computer related words such as FTP. Since we're using TF_IDF this makes sense that its inaccurately predicting these documents. Terms like FTP and ASCII is infrequently used in general documents, but much more likely in documents related to computers. The text that was mis-categorized, was religious content, but that was discussing stuff related to computers such as FTP. \n",
    "\n",
    "There are a few options to improve the model. \n",
    "1. Remove features that have a high true positive, but also a high false positive rate. These are features that are predictive of the category, but are also used enough in other categories that it drives false positives. \n",
    "2. Remove stop words\n",
    "3. For the misclassification of religion as comp.graphics, its likely due to terms like FTP having a higher weight then terms like Book of Mormom. One option would be to identify those proper nouns and weight them higher. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(8) EXTRA CREDIT\n",
    "\n",
    "Try implementing one of your ideas based on your error analysis. Use logistic regression as your underlying model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "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.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}