<jupyter_start><jupyter_text>Information Retrieval I  Assignment 2: retrieval models [100 points + 10 bonus points] **TA**: <NAME> (<EMAIL>; C3.258B, Science Park 904)**Secondary TAs**: <NAME>, <NAME> In this assignment you will get familiar with basic information retrieval concepts. You will implement and evaluate different information retrieval ranking models and evaluate their performance.We provide you with a VirtualBox image that comes pre-loaded with an index and a Python installation. To query the index, you'll use a Python package ([pyndri](https://github.com/cvangysel/pyndri)) that allows easy access to the underlying document statistics.For evaluation you'll use the [TREC Eval](https://github.com/usnistgov/trec_eval) utility, provided by the National Institute of Standards and Technology of the United States. TREC Eval is the de facto standard way to compute Information Retrieval measures and is frequently referenced in scientific papers.This is a **groups-of-two assignment**, the deadline is **23:59 - 25 January, 2017**. Code quality, informative comments and convincing analysis of the results will be considered when grading. Submission should be done through blackboard, questions can be asked on the course [Piazza](https://piazza.com/class/ixoz63p156g1ts). Technicalities (must-read!) This assignment comes pre-loaded on a VirtualBox running Ubuntu. We have configured the indexing software and Python environment such that it works out of the box. You are allowed to extract the files from the VirtualBox and set-up your own non-virtualized environment. However, in this case you are on your own w.r.t. software support.The assignment directory is organized as follows:   * `./assignment.ipynb` (this file): the description of the assignment.   * `./index/`: the index we prepared for you.   * `./ap_88_90/`: directory with ground-truth and evaluation sets:      * `qrel_test`: test query relevance collection (**test set**).      * `qrel_validation`: validation query relevance collection (**validation set**).      * `topics_title`: semicolon-separated file with query identifiers and terms.      `Python + Jupyter`, `Indri`, `Gensim` and `Pyndri` come pre-installed (see `$HOME/.local`). TREC Eval can be found in `$HOME/Downloads/trec_eval.9.0`. The password of the `student` account on the VirtualBox is `<PASSWORD>`.  TREC Eval primer The TREC Eval utility can be downloaded and compiled as follows:    git clone https://github.com/usnistgov/trec_eval.git    cd trec_eval    makeTREC Eval computes evaluation scores given two files: ground-truth information regarding relevant documents, named *query relevance* or *qrel*, and a ranking of documents for a set of queries, referred to as a *run*. The *qrel* will be supplied by us and should not be changed. For every retrieval model (or combinations thereof) you will generate a run of the top-1000 documents for every query. The format of the *run* file is as follows:    $query_identifier Q0 $document_identifier $rank_of_document_for_query $query_document_similarity $run_identifier    where   * `$query_identifier` is the unique identifier corresponding to a query (usually this follows a sequential numbering).   * `Q0` is a legacy field that you can ignore.   * `$document_identifier` corresponds to the unique identifier of a document (e.g., APXXXXXXX where AP denotes the collection and the Xs correspond to a unique numerical identifier).   * `$rank_of_document_for_query` denotes the rank of the document for the particular query. This field is ignored by TREC Eval and is only maintained for legacy support. The ranks are computed by TREC Eval itself using the `$query_document_similarity` field (see next). However, it remains good practice to correctly compute this field.   * `$query_document_similarity` is a score indicating the similarity between query and document where a higher score denotes greater similarity.   * `$run_identifier` is an identifier of the run. This field is for your own convenience and has no purpose beyond bookkeeping.   For example, say we have two queries: `Q1` and `Q2` and we rank three documents (`DOC1`, `DOC2`, `DOC3`). For query `Q1`, we find the following similarity scores `score(Q1, DOC1) = 1.0`, `score(Q1, DOC2) = 0.5`, `score(Q1, DOC3) = 0.75`; and for `Q2`: `score(Q2, DOC1) = -0.1`, `score(Q2, DOC2) = 1.25`, `score(Q1, DOC3) = 0.0`. We can generate run using the following snippet:<jupyter_code>import logging
import sys

def write_run(model_name, data, out_f,
              max_objects_per_query=sys.maxsize,
              skip_sorting=False):
    """
    Write a run to an output file.
    Parameters:
        - model_name: identifier of run.
        - data: dictionary mapping topic_id to object_assesments;
            object_assesments is an iterable (list or tuple) of
            (relevance, object_id) pairs.
            The object_assesments iterable is sorted by decreasing order.
        - out_f: output file stream.
        - max_objects_per_query: cut-off for number of objects per query.
    """
    for subject_id, object_assesments in data.items():
        if not object_assesments:
            logging.warning('Received empty ranking for %s; ignoring.',
                            subject_id)

            continue

        # Probe types, to make sure everything goes alright.
        # assert isinstance(object_assesments[0][0], float) or \
        #     isinstance(object_assesments[0][0], np.float32)
        assert isinstance(object_assesments[0][1], str) or \
            isinstance(object_assesments[0][1], bytes)

        if not skip_sorting:
            object_assesments = sorted(object_assesments, reverse=True)

        if max_objects_per_query < sys.maxsize:
            object_assesments = object_assesments[:max_objects_per_query]

        if isinstance(subject_id, bytes):
            subject_id = subject_id.decode('utf8')

        for rank, (relevance, object_id) in enumerate(object_assesments):
            if isinstance(object_id, bytes):
                object_id = object_id.decode('utf8')

            out_f.write(
                '{subject} Q0 {object} {rank} {relevance} '
                '{model_name}\n'.format(
                    subject=subject_id,
                    object=object_id,
                    rank=rank + 1,
                    relevance=relevance,
                    model_name=model_name))
            
# The following writes the run to standard output.
# In your code, you should write the runs to local
# storage in order to pass them to trec_eval.
with open('./example.run','w') as f:
    write_run(
        model_name='example',
        data={
            'Q1': ((1.0, 'DOC1'), (0.5, 'DOC2'), (0.75, 'DOC3')),
            'Q2': ((-0.1, 'DOC1'), (1.25, 'DOC2'), (0.0, 'DOC3')),
        },
        out_f=f,
        max_objects_per_query=1000)<jupyter_output><empty_output><jupyter_text>Now, imagine that we know that `DOC1` is relevant and `DOC3` is non-relevant for `Q1`. In addition, for `Q2` we only know of the relevance of `DOC3`. The query relevance file looks like:    Q1 0 DOC1 1    Q1 0 DOC3 0    Q2 0 DOC3 1    We store the run and qrel in files `example.run` and `example.qrel` respectively on disk. We can now use TREC Eval to compute evaluation measures. In this example, we're only interested in Mean Average Precision and we'll only show this below for brevity. However, TREC Eval outputs much more information such as NDCG, recall, precision, etc.    $ trec_eval -m all_trec -q example.qrel example.run | grep -E "^map\s"    > map                   	Q1	1.0000    > map                   	Q2	0.5000    > map                   	all	0.7500    Now that we've discussed the output format of rankings and how you can compute evaluation measures from these rankings, we'll now proceed with an overview of the indexing framework you'll use.  Pyndri primer For this assignment you will use [Pyndri](https://github.com/cvangysel/pyndri) [[1](https://arxiv.org/abs/1701.00749)], a python interface for [Indri](https://www.lemurproject.org/indri.php). We have indexed the document collection and you can query the index using Pyndri. We will start by giving you some examples of what Pyndri can do:First we read the document collection index with Pyndri:<jupyter_code>import pyndri

index = pyndri.Index('index/')<jupyter_output><empty_output><jupyter_text>The loaded index can be used to access a collection of documents in an easy manner. We'll give you some examples to get some idea of what it can do, it is up to you to figure out how to use it for the remainder of the assignment.First let's look at the number of documents, since Pyndri indexes the documents using incremental identifiers we can simply take the lowest index and the maximum document and consider the difference:<jupyter_code>print("There are %d documents in this collection." % (index.maximum_document() - index.document_base()))<jupyter_output><empty_output><jupyter_text>Let's take the first document out of the collection and take a look at it:<jupyter_code>example_document = index.document(index.document_base())
print(example_document, 'lol')<jupyter_output><empty_output><jupyter_text>Here we see a document consists of two things, a string representing the external document identifier and an integer list representing the identifiers of words that make up the document. Pyndri uses integer representations for words or terms, thus a token_id is an integer that represents a word whereas the token is the actual text of the word/term. Every id has a unique token and vice versa with the exception of stop words: words so common that there are uninformative, all of these receive the zero id.To see what some ids and their matching tokens we take a look at the dictionary of the index:<jupyter_code>token2id, id2token, id2df = index.get_dictionary()
print(list(id2token.items())[:15])<jupyter_output><empty_output><jupyter_text>Using this dictionary we can see the tokens for the (non-stop) words in our example document:<jupyter_code>print([id2token[word_id] for word_id in example_document[1] if word_id > 0])
print(max(id2token.keys()))<jupyter_output><empty_output><jupyter_text>The reverse can also be done, say we want to look for news about the "University of Massachusetts", the tokens of that query can be converted to ids using the reverse dictionary:<jupyter_code>query_tokens = index.tokenize("University of Massachusetts")
print("Query by tokens:", query_tokens)
query_id_tokens = [token2id.get(query_token,0) for query_token in query_tokens]
print("Query by ids with stopwords:", query_id_tokens)
query_id_tokens = [word_id for word_id in query_id_tokens if word_id > 0]
print("Query by ids without stopwords:", query_id_tokens)<jupyter_output><empty_output><jupyter_text>Naturally we can now match the document and query in the id space, let's see how often a word from the query occurs in our example document: While this is certainly not everything Pyndri can do, it should give you an idea of how to use it. Please take a look at the [examples](https://github.com/cvangysel/pyndri) as it will help you a lot with this assignment.**CAUTION**: Avoid printing out the whole index in this Notebook as it will generate a lot of output and is likely to corrupt the Notebook.<jupyter_code>matching_words = sum([True for word_id in example_document[1] if word_id in query_id_tokens])
print("Document %s has %d word matches with query: \"%s\"." % (example_document[0], matching_words, ' '.join(query_tokens)))
print("Document %s and query \"%s\" have a %.01f%% overlap." % (example_document[0], ' '.join(query_tokens),matching_words/float(len(example_document[1]))*100))<jupyter_output><empty_output><jupyter_text>Parsing the query fileYou can parse the query file (`ap_88_89/topics_title`) using the following snippet:<jupyter_code>import collections
import io
import logging
import sys

def parse_topics(file_or_files,
                 max_topics=sys.maxsize, delimiter=';'):
    assert max_topics >= 0 or max_topics is None

    topics = collections.OrderedDict()

    if not isinstance(file_or_files, list) and \
            not isinstance(file_or_files, tuple):
        if hasattr(file_or_files, '__iter__'):
            file_or_files = list(file_or_files)
        else:
            file_or_files = [file_or_files]

    for f in file_or_files:
        assert isinstance(f, io.IOBase)

        for line in f:
            assert(isinstance(line, str))

            line = line.strip()

            if not line:
                continue

            topic_id, terms = line.split(delimiter, 1)

            if topic_id in topics and (topics[topic_id] != terms):
                    logging.error('Duplicate topic "%s" (%s vs. %s).',
                                  topic_id,
                                  topics[topic_id],
                                  terms)

            topics[topic_id] = terms

            if max_topics > 0 and len(topics) >= max_topics:
                break

    return topics

with open('./ap_88_89/topics_title', 'r') as f_topics:
    queries = parse_topics([f_topics])
    print('Num of Queries:', len(queries))<jupyter_output><empty_output><jupyter_text>Task 1: Implement and compare lexical IR methods [45 points]  In this task you will implement a number of lexical methods for IR using the **Pyndri** framework. Then you will evaluate these methods on the dataset we have provided using **TREC Eval**.Use the **Pyndri** framework to get statistics of the documents (term frequency, document frequency, collection frequency; **you are not allowed to use the query functionality of Pyndri**) and implement the following scoring methods in **Python**:- [TF-IDF](http://nlp.stanford.edu/IR-book/html/htmledition/tf-idf-weighting-1.html). **[5 points]**- [BM25](http://nlp.stanford.edu/IR-book/html/htmledition/okapi-bm25-a-non-binary-model-1.html) with k1=1.2 and b=0.75. **[5 points]**- Language models ([survey](https://drive.google.com/file/d/0B-zklbckv9CHc0c3b245UW90NE0/view))    - Jelinek-Mercer (explore different values of 𝛌 in the range [0.1, 0.2, ..., 0.9]). **[5 points]**    - Dirichlet Prior (explore different values of 𝛍 [500, 1000, ..., 2000]). **[5 points]**    - Absolute discounting (explore different values of 𝛅 in the range [0.1, 0.2, ..., 0.9]). **[5 points]**    - [Positional Language Models](http://sifaka.cs.uiuc.edu/~ylv2/pub/sigir09-plm.pdf) define a language model for each position of a document, and score a document based on the scores of its PLMs. The PLM is estimated based on propagated counts of words within a document through a proximity-based density function, which both captures proximity heuristics and achieves an effect of “soft” passage retrieval. Implement the PLM, all five kernels, but only the Best position strategy to score documents. Use 𝛔 equal to 50, and Dirichlet smoothing with 𝛍 optimized on the validation set (decide how to optimize this value yourself and motivate your decision in the report). **[10 points]**    Implement the above methods and report evaluation measures (on the test set) using the hyper parameter values you optimized on the validation set (also report the values of the hyper parameters). Use TREC Eval to obtain the results and report on `NDCG@10`, Mean Average Precision (`MAP@1000`), `Precision@5` and `Recall@1000`.For the language models, create plots showing `NDCG@10` with varying values of the parameters. You can do this by chaining small scripts using shell scripting (preferred) or execute trec_eval using Python's `subprocess`.Compute significance of the results using a [two-tailed paired Student t-test](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.ttest_rel.html) **[10 points]**. Be wary of false rejection of the null hypothesis caused by the [multiple comparisons problem](https://en.wikipedia.org/wiki/Multiple_comparisons_problem). There are multiple ways to mitigate this problem and it is up to you to choose one.Analyse the results by identifying specific queries where different methods succeed or fail and discuss possible reasons that cause these differences.**NOTE**: Don’t forget to use log computations in your calculations to avoid underflows.<jupyter_code>%pylab inline
#%load_ext memory_profiler
import collections
import scipy.spatial
import scipy.stats as st
import bottleneck
from functools import reduce
from subprocess import Popen, PIPE
import pyndri
import pickle
import copy
import gensim
import logging
import pyndri.compat
import sys
import pandas as pd
import seaborn as sns
sns.set(style="white", color_codes=True, font_scale = 1.3)
logger = logging.getLogger()
logger.setLevel(logging.DEBUG)
## Sign/Binomial test
# Performs a sign/binomial test on the proportion of impressions where E outperformed P across all ranking pairs.
# Inputs:
# - Proportions: a vector of statistics representing a difference or ratio of sample A and sample B
# - Threshold: a value t which is used to assign the label +1 if the proportions are above t, or -1 if below
# The sign test does not take into account examples where both algorithms are tied (a sample value is equal to the threshold)
# Returns: p-value
def sign_binomial_test(proportions, threshold = 0.5, verbose=True, allow_ties=False, alternative = 'two-sided', sig = 0.05):
  k = len(np.where(proportions > threshold)[0])
  
  if alternative == 'two-sided':
    sig /= 2
    
  if allow_ties:
    N = len(proportions)
  else:
    num_ties = len(np.where(proportions == threshold)[0])
    N = len(proportions) - num_ties
    
  if verbose:
    print('Running %s Binomial test on k = %s, n = %s' % (alternative, k, N))
    print('Percentage of examples where E>P: %s' % (k * 1.0 / N))

  pval = st.binom_test(k, N, 0.5, alternative)

  if verbose:
    if pval < sig:
      print('Should reject H0 given p-val = %s ' % (pval))
    else:
      print('Not enough evidence to reject H0 given p-val = %s \n' % (pval))
  
  return pval
## Load cache data into memory
# To improve performance and reducing IO overhead, at the expensive of high-memory usage, 
# we load a lot of data into dictionaries for O(1) reading operations

# Pyndri
index = pyndri.Index('index/')
token2id, id2token, id2df = index.get_dictionary() # id2df: id->num of docs containing term
id2tf = index.get_term_frequencies() # id->num of occurences in whole collection

doc_idx = range(index.document_base(), index.maximum_document()) # range of document indices

# Vocabulary size
vocab_size = len(id2tf) 

# A dictionary of document_index:bow pairs. Each bow is a dictionary of token_id:tf(tokenid, document). We remove stop words
docs_idx = {_d:{k:v for k,v in collections.Counter(index.document(_d)[1]).items() if k != 0} for _d in doc_idx}

# A dictionary of document_index:document_length pairs
doc_lengths = {_d:float(sum([docs_idx[_d][_t] for _t in docs_idx[_d] if not _t == 0])) for _d in doc_idx}

# Average document length
avg_doc_length = mean(list(doc_lengths.values()))

# Size of collection excluding stop words
collection_size_nostop = sum(list(doc_lengths.values()))

# number of documents in the collection
num_docs = len(docs_idx)

del doc_idx

# Load validation and test queries into two separate OrderedDicts

with open("./ap_88_89/qrel_test", "r") as ins:
    test_ids = set()
    for line in ins:
        test_ids.add(line.split(' ')[0])
    
    test_queries = collections.OrderedDict()
    for qid in test_ids:
      test_queries[qid] = queries[qid]


with open("./ap_88_89/qrel_validation", "r") as ins:
    val_ids = set()
    for line in ins:
        val_ids.add(line.split(' ')[0])
    
    val_queries = collections.OrderedDict()
    for qid in val_ids:
      val_queries[qid] = queries[qid]
# Normalizes and tokenizes a string of text
# Removes unusual characters, tokenizes and then stems documents
# Returns a list of normalized tokens
def normalize(text,index):
  try:
    return [pyndri.stem(x) for x in index.tokenize(pyndri.escape(text)) if not x =='' ]
  except(OSError):
    print('NORMALIZE ERROR: ', text)
  return []
# Returns the document in position pos in collection index
def get_doc(pos, index):
  return index.document(pos + index.document_base())

# Prints the contents of a document
def print_doc(doc):
  print([id2token[id] for id in doc[1] if id>0])

# Return number of documents in which term t appears
def df(t):
  return id2df.get(t,0)

# Returns idf for term t in a collection of size n documents
# Adds 1 to both the numerator and denominator avoid division by zero
def idf(t, n):
  return np.log(float(n + 1)/(df(t) + 1))

# Counts the occurences of term t in document d
def tf(t,d):
  return np.count_nonzero(np.array(d[1])==t)

# Returns the number of times term t occurs in document d (given as a dictionary)
def fast_tf(t,d):
  return d[t] if t in d else 0

# Calculates the tf-idf score of a term t that appears in document d
def tfidf(t,d,n,tf_fn = fast_tf):
  return np.log(1 + tf_fn(t,d)) * idf(t,n)

## Perform sorting on the top k elements of an array
# Uses partitioning to create an ordering on left and right sides of the kth element of input array
# Then sorts the i=[0,k-1] results in descending order
# Inputs:
# - values: array to be sorted
# - k: cut-off value. Entires before k are sorted
# Returns:
# - length k array of sorted idx corresponding to input array
def top_k_sort(values, k):
  
  assert k <= len(values), "Partition index larger than array length"
  
  partition_idx = bottleneck.argpartition(values, len(values) - k)[-k:]
  partition_score = values[partition_idx]
  
  sorted_idx = np.argsort(partition_score)[::-1]
  
  return partition_idx[sorted_idx]

# Receives a list of scores for the documents and returns the indices that give the top rank_size ranking 
# by calling function top_k_sort
def rank_results(scores, index, rank_size=1000):
  sorted_idx = top_k_sort(scores, rank_size) 
  return [x+index.document_base() for x in sorted_idx], scores[sorted_idx] # indexing of pyndri docs begin at index.document_base(), so the sorted_idx is adjusted accordingly


# Calls TREC EVAL and reports results
def run_trec_eval(query_results, model_name, metrics, groundtruth_qrel='./qrel_validation', qrun='./tmp_run', max_results=1000):
  
  # Create run file
  with open(qrun,'w') as f:
      write_run(
          model_name=model_name,
          data=query_results,
          out_f=f,
          max_objects_per_query=max_results, skip_sorting=False)
  
  results = np.zeros(len(metrics))
  
  status_ok = 1
  
  for m in range(len(metrics)):
    p = Popen(['trec_eval -m %s %s %s' % (metrics[m], groundtruth_qrel, qrun)],cwd='/home/student/assignment', stdin=PIPE, stdout=PIPE, stderr=PIPE, shell=True, executable='/bin/bash')
    out, err = p.communicate()

    if(err):
      print('ERROR: ', err.decode())
      results[m]= -666
      status_ok = 0

    results[m] = float(out.decode().split('\n')[0].split('\t')[-1])
    
  return status_ok, results<jupyter_output><empty_output><jupyter_text>TFIDF<jupyter_code># TF-IDF scoring model
# Inputs:
# - q: query tuple (id, query string)
# - index: Pyndri index
# - rank_size: number of results to return. Sorting done efficiently here.
# - docs_idx: document:bag of words dictionary
# Returns: 
# - sorted_idx: the idx of documents in the index object, sorted by their score
# - scores: td-idf scores of each ranking
def tfidf_scoring(query, index, docs_idx, doc_lengths, token2id, rank_size = 1000, sort_results = True):
  
  n = len(docs_idx) # num of docs
  
  unique_tokens = list(set(normalize(query, index))) # unique tokens in query
  token_ids = [token2id[token] for token in unique_tokens if not token == '' and token in token2id] # unique token ids
  
  scores = np.zeros(len(docs_idx))  # score for each document
  i = 0
  
  for d_idx, doc in docs_idx.items():
    
    common_tokens = [t_id for t_id in token_ids if t_id in doc] # find common tokens in query and document
    
    # If document is empty or no lexical matches
    if doc_lengths[d_idx] == 0 or len(common_tokens) == 0:
      scores[i] = -Inf
      i += 1
      continue
  
    for t_id in token_ids:
      if t_id in doc:
        scores[i] += tfidf(t_id, doc, n) # for each common token, accumulate tfidf score
     
    scores[i] *= 1
    
    i+=1
  
  if not sort_results:
    return scores
    
  return rank_results(scores, index, rank_size) # sorts results

## Runs TFIDF on test queries

metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']

q_db = test_queries
bn_name = 'test'
tfidf_test_matrix = np.zeros((len(q_db),len(metrics)))
results_dic = {}

for qnum, qid in enumerate(q_db):
  results_dic.clear()
  print('Q ' + str(qnum) + '/' + str(len(q_db)))
  idx, scores = tfidf_scoring(q_db[qid], index, docs_idx, doc_lengths, token2id, 1000)
  result_ids = [index.document(i)[0] for i in idx]
  results_dic[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx))) 
  code, res = run_trec_eval(results_dic, 'TFIDF', metrics, groundtruth_qrel='./ap_88_89/qrel_'+bn_name, qrun='./tmp_run', max_results=1000)
  tfidf_test_matrix[qnum,:] = np.array(res)
   
# Store results
with open('./tfidf_test_matrix', 'wb') as f:
  pickle.dump(tfidf_test_matrix, f)

print('Result shape: ' + str(tfidf_test_matrix.shape))
print(np.mean(tfidf_test_matrix, axis=0))

del tfidf_test_matrix<jupyter_output><empty_output><jupyter_text>BM25<jupyter_code># Performs BM25 scoring on a query string
# Default paramters are as instructured in task 1
# Returns sorted document indexes and scores
def bm25(query, index, docs_idx, doc_lengths, token2id, avg_doc_length, k1=1.2, b=0.75, rank_size=1000, sort_results = True):
  
    n = len(docs_idx) # num of docs

    unique_tokens = list(set(normalize(query, index)))
    token_ids = [token2id[token] for token in unique_tokens if not token == '' and token in token2id]

    scores = np.zeros(len(docs_idx)) 
    i = 0
    for d_idx, doc in docs_idx.items():

        if doc_lengths[d_idx] == 0:
            scores[i] = -Inf
            i += 1
            continue

        for t_id in token_ids:
            if t_id in doc:  # if term not in document, tf(t,d) = 0, then bm25 score of that term is 0
                tf = doc[t_id]
                scores[i] += ((k1+1) * tf) / (k1*((1-b) + b*(doc_lengths[d_idx]/avg_doc_length)) + tf) * idf(t_id, n)
        i+=1

    if not sort_results:
        return scores
  
    return rank_results(scores, index, rank_size)
## BM25 on test queries

metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']

q_db = test_queries
bn_name = 'test'

bm25_test_matrix = np.zeros((len(q_db),len(metrics)))


results_dic = {}

for qnum, qid in enumerate(q_db):  # for each query
    
    results_dic.clear()
    print('Q ' + str(qnum) + '/' + str(len(q_db)))
    idx, scores = bm25(q_db[qid], index, docs_idx, doc_lengths, token2id, avg_doc_length, k1=1.2, b=0.75, rank_size=1000) 
    result_ids = [index.document(i)[0] for i in idx]
    results_dic[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx))) 
    code, res = run_trec_eval(results_dic, 'BM25', metrics, groundtruth_qrel='./ap_88_89/qrel_'+bn_name, qrun='./tmp_run', max_results=1000)
    bm25_test_matrix[qnum,:] = np.array(res)
   

with open('./bm25_test_matrix', 'wb') as f:
    pickle.dump(bm25_test_matrix, f)

print('Result shape: ' + str(bm25_test_matrix.shape))
print(np.mean(bm25_test_matrix, axis=0))

del bm25_test_matrix<jupyter_output><empty_output><jupyter_text>Jelinek-Mercer<jupyter_code># Jelinek-Mercer with query likelihood
# Inputs:
# - query: query string
# - Parameter _lambda 
# Outputs:
# - sorted document index
# - sorted corresponding scores
def jelinek_mercer(query, index, docs_idx, doc_lengths, collection_size_nostop, term_counts, _lambda, rank_size=1000, sort_results = True):
  
  n = len(docs_idx) # num of docs
  
  tokens = normalize(query, index)
  token_ids = [token2id[token] for token in tokens if not token == '' and token in token2id]
  
  scores = np.zeros(len(docs_idx))  # score for each doc
  i = 0
  for d_idx, doc in docs_idx.items():
    
    if doc_lengths[d_idx] == 0: # ignore empty docs
      scores[i] = -Inf
      i += 1
      continue
    
    for t_id in token_ids: # for each token
        tf = doc[t_id] if t_id in doc else 0 # get tf
        scores[i] += np.log(_lambda * tf / doc_lengths[d_idx] + (1 - _lambda) * term_counts[t_id] / collection_size_nostop)  
    i+=1
  
  if not sort_results:
    return scores
  
  return rank_results(scores, index, rank_size)
## JM - optimization on validation queries

eval_metric = ['ndcg_cut.10']
_lambdas = [0.1,0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9]
_lambda_scores = np.zeros(len(_lambdas))


for i_lambda, _lambda in enumerate(_lambdas):
  print('Evaluating lambda=',_lambda)
  tmp_res = {}
  
  for qid, query in val_queries.items(): # run validation set queries
    idx, scores = jelinek_mercer(query, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, _lambda, rank_size=1000) 
    result_ids = [index.document(i)[0] for i in idx]
    tmp_res[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx)))
  
  print('Running treceval lambda=',_lambda)
  code, _lambda_scores[i_lambda] = run_trec_eval(tmp_res, 'JM', eval_metric, groundtruth_qrel='./ap_88_89/qrel_validation', qrun='./tmp_run', max_results=1000)
    
best_lambda = _lambdas[np.argmax(_lambda_scores)]
print('Best lambda:', best_lambda )
res = [(_lambdas[i],_lambda_scores[i]) for i in range(len(_lambdas))]
print(res)
  
with open('./jm_val_ndcg', 'wb') as f:
  pickle.dump(res, f)
  
plot(_lambdas, _lambda_scores)
plt.ylabel('NDCG@10')
plt.xlabel('$\lambda$')
plt.title('Performance of Jelinek-Mercer on NDCG@10 for several $\lambda$ values')
plt.savefig('./jm_optim.png')
## JM on test queries

metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
best_lambda = 0.4
q_db = test_queries
bn_name = 'test'

jm_test_matrix = np.zeros((len(q_db),len(metrics)))

results_dic = {}
for qnum, qid in enumerate(q_db):
  results_dic.clear()
  print('Q ' + str(qnum) + '/' + str(len(q_db)))
  idx, scores = jelinek_mercer(q_db[qid], index, docs_idx, doc_lengths, collection_size_nostop, id2tf, best_lambda, rank_size=1000) 
  result_ids = [index.document(i)[0] for i in idx]
  results_dic[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx))) 
  code, res = run_trec_eval(results_dic, 'JM', metrics, groundtruth_qrel='./ap_88_89/qrel_'+bn_name, qrun='./tmp_run', max_results=1000)
  jm_test_matrix[qnum,:] = np.array(res)
   
with open('./jm_test_matrix', 'wb') as f:
  pickle.dump(jm_test_matrix, f)

print('Result shape: ' + str(jm_test_matrix.shape))
print(np.mean(jm_test_matrix, axis=0))

del jm_test_matrix<jupyter_output><empty_output><jupyter_text>Absolute Discounting<jupyter_code># Absolute Discounting with query likelihood 
# Inputs:
# - query: query string
# - Parameter delta 
# Outputs:
# - sorted document index
# - sorted corresponding scores
def abs_discounting(query, index, docs_idx, doc_lengths, collection_size_nostop, term_counts, delta, rank_size=1000, sort_results = True):
  
  n = len(docs_idx) # num of docs
  
  tokens = normalize(query, index)
  token_ids = [token2id[token] for token in tokens if not token == '' and token in token2id]
  
  scores = np.zeros(len(docs_idx)) 
  i = 0
  for d_idx, doc in docs_idx.items():
    
    if doc_lengths[d_idx] == 0:
      scores[i] = -Inf
      i += 1
      continue
      
    for t_id in token_ids:
        tf = doc[t_id] if t_id in doc else 0
        scores[i] += np.log( (max( tf - delta, 0)  + delta * len(doc) * (term_counts[t_id] / collection_size_nostop) ) / doc_lengths[d_idx] )  
    i+=1
  
  if not sort_results:
    return scores
  
  return rank_results(scores, index, rank_size)

## Abs Discounting - optimization on validation queries
# Plots results
eval_metric = ['ndcg_cut.10']
validation_jm = {}
deltas = [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9]
delta_scores = np.zeros(len(deltas))


for i_delta, delta in enumerate(deltas):
  print('Evaluating delta=',delta)
  tmp_res = {}
  
  for qid, query in val_queries.items(): # run validation set queries
    idx, scores = abs_discounting(query, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, delta, rank_size=1000)
    result_ids = [index.document(i)[0] for i in idx]
    tmp_res[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx)))
  
  print('Running treceval delta=',delta)
  code, delta_scores[i_delta] = run_trec_eval(tmp_res, 'AD', eval_metric, groundtruth_qrel='./ap_88_89/qrel_validation', qrun='./tmp_run', max_results=1000)
    
best_delta = deltas[np.argmax(delta_scores)]
print('Best delta:', best_delta )
res = [(deltas[i],delta_scores[i]) for i in range(len(deltas))]
print(res)
  
with open('./ad_val_ndcg', 'wb') as f:
  pickle.dump(res, f)
  
plot(deltas, delta_scores)
plt.ylabel('NDCG@10')
plt.xlabel('$\delta$')
plt.title('Performance of Abs. Discounting on NDCG@10 for several $\delta$ values')
plt.savefig('./ad_optim.png')
## Abs Discounting on test queries

metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
q_db = test_queries
bn_name = 'test'
best_delta = 0.8

ad_test_matrix = np.zeros((len(q_db),len(metrics)))

results_dic = {}
for qnum, qid in enumerate(q_db):
  results_dic.clear()
  print('Q ' + str(qnum) + '/' + str(len(q_db)))
  idx, scores = abs_discounting(q_db[qid], index, docs_idx, doc_lengths, collection_size_nostop, id2tf, best_delta, rank_size=1000)
  result_ids = [index.document(i)[0] for i in idx]
  results_dic[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx))) 
  code, res = run_trec_eval(results_dic, 'AD', metrics, groundtruth_qrel='./ap_88_89/qrel_'+bn_name, qrun='./tmp_run', max_results=1000)
  ad_test_matrix[qnum,:] = np.array(res)
   
with open('./ad_test_matrix', 'wb') as f:
  pickle.dump(ad_test_matrix, f)

print('Result shape: ' + str(ad_test_matrix.shape))
print(np.mean(ad_test_matrix, axis=0))

del ad_test_matrix<jupyter_output><empty_output><jupyter_text>Dirichlet Smoothing<jupyter_code># Dirichlet smoothing with query likelihood
# Inputs:
# - query: query string
# - Parameter _lambda 
# Outputs:
# - sorted document index
# - sorted corresponding scores
def dirichlet_smoothing(query, index, docs_idx, doc_lengths, collection_size_nostop, term_counts, mu, rank_size=1000, sort_results = True):
  
  tokens = normalize(query, index)
  token_ids = [token2id[token] for token in tokens if not token == '' and token in token2id]
  
  scores = np.zeros(len(docs_idx))
  i=0
  
  for d_idx, doc in docs_idx.items():
    
    if doc_lengths[d_idx] == 0:
      scores[i] = -Inf
      i += 1
      continue
    
    for t_id in token_ids:
      tf = doc[t_id] if t_id in doc else 0
      p = (tf + mu * (term_counts[t_id]/collection_size_nostop)) / (doc_lengths[d_idx] + mu) 
      scores[i] += np.log(p)
      
    i+=1
  
  if not sort_results:
    return scores
  
  return rank_results(scores, index, rank_size)
## Dirichlet smoothing - optimization on validation set
# Plots results
eval_metric = ['ndcg_cut.10']
validation_dirichlet = {}
mus = [500,1000,1500,2000,3500,5000]
mu_scores = np.zeros(len(mus))

for i_mu, mu in enumerate(mus):
  print('Evaluating mu=',mu)
  tmp_res = {}
  
  for qid, query in val_queries.items(): # run validation set queries
    idx, scores = dirichlet_smoothing(query, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, mu, 1000)
    result_ids = [index.document(i)[0] for i in idx]
    tmp_res[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx)))
  
  print('Running treceval mu=',mu)
  code, mu_scores[i_mu] = run_trec_eval(tmp_res, 'Dirichlet', eval_metric, groundtruth_qrel='./ap_88_89/qrel_validation', qrun='./tmp_run', max_results=1000)
  
best_mu = mus[np.argmax(mu_scores)]
print('Best mu:', best_mu )
res = [(mus[i],mu_scores[i]) for i in range(len(mus))]
print(res)
  
with open('./dir_val_ndcg', 'wb') as f:
  pickle.dump(res, f)

plot(mus, mu_scores)
plt.ylabel('NDCG@10')
plt.xlabel('$\mu$')
plt.title('Performance of Dirichlet Smoothing on NDCG@10 for several $\mu$ values')
plt.savefig('./dir_optim.png')
## Dirichlet Discounting on test queries

metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
q_db = test_queries
bn_name = 'test'
best_mu = 2000

dir_test_matrix = np.zeros((len(q_db),len(metrics)))

results_dic = {}
for qnum, qid in enumerate(q_db):
  results_dic.clear()
  print('Q ' + str(qnum) + '/' + str(len(q_db)))
  idx, scores = dirichlet_smoothing(test_queries[qid], index, docs_idx, doc_lengths, collection_size_nostop, id2tf, best_mu, 1000)
  result_ids = [index.document(i)[0] for i in idx]
  results_dic[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx))) 
  code, res = run_trec_eval(results_dic, 'DIR', metrics, groundtruth_qrel='./ap_88_89/qrel_'+bn_name, qrun='./tmp_run', max_results=1000)
  dir_test_matrix[qnum,:] = np.array(res)
   
with open('./dir_test_matrix', 'wb') as f:
  pickle.dump(dir_test_matrix, f)

print('Result shape: ' + str(dir_test_matrix.shape))
print(np.mean(dir_test_matrix, axis=0))

del dir_test_matrix<jupyter_output><empty_output><jupyter_text>Positional Language Model<jupyter_code># Heavily based on PLM's paper author <NAME>'s code
# http://sifaka.cs.uiuc.edu/~ylv2/pub/plm/PLMRetEval.cpp


def passage_kernel(dis, sigma):
  return float(np.abs(dis) <= sigma)


def gaussian_kernel(dis, sigma):
  return np.exp(-dis**2 / (2 * sigma**2))


def triangle_kernel(dis, sigma):
  if np.abs(dis) <= sigma:
    return 1 - np.abs(dis)/sigma
  else:
    return 0

def cosine_kernel(dis, sigma):
  if np.abs(dis) <= sigma:
    return 0.5 * (1 + np.cos(pi*np.abs(dis)/sigma))
  else:
    return 0

  
def circle_kernel(dis, sigma):
  if np.abs(dis) <= sigma:
    return np.sqrt(1 - (dis/sigma)**2)
  else:
    return 0
  
def gaussian_cdf(x, mean, sigma):
  x = (x - mean) / sigma
  if x==0:
    res = 0.5
  else:
    t = 1 / (1 + 0.2316419 * np.abs(x))
    t *= t * (1 / np.sqrt(2 * pi)) * np.exp(-0.5 * x * x) * (0.31938153 + t * (-0.356563782 + t * (1.781477937 + t * (-1.821255978 + t * 1.330274429))))
    if x >= 0:
      res = 1 - t
    else:
      res = t
  return sqrt(2 * pi) * sigma * res

def triangle_cdf(x, mean, sigma):
  x = (x - mean) / sigma
  if (x >= 1):
    res = sigma
  elif(x < -1):
    res = 0
  elif(x < 0):
    res = sigma * (1 - np.abs(x)) * (1 - np.abs(x)) / 2.0
  else:
    res = sigma - sigma * (1 - x) * (1 - x) / 2.0
  return res

def cosine_cdf(x, mean, sigma):
  x = (x - mean) / sigma
  if (x >= 1):
    res = sigma
  elif(x < -1):
    res = 0
  elif(x < 0):
    res = sigma * (1 + x - np.sin(pi * x) / pi) / 2.0
  else:
    res = sigma - sigma * (1 - x + np.sin(pi * x) / pi) / 2.0
  return res

def circle_cdf(x, mean, sigma):
  x = (x - mean) / sigma
  if (x >= 1):
    res = (pi - 2.0) * sigma
  elif(x < -1):
    res = 0
  elif (x < 0):
    res = sigma * (np.arcsin(x) + pi / 2.0 - np.sqrt(1 - x * x))
  else:
    res = (pi - 2.0) * sigma - sigma * (np.arcsin(-x) + pi / 2.0 - np.sqrt(1 - x * x))
  return res

def passage_cdf(x, mean, sigma):
  res = 0
  if sigma > x- mean:
    res += x - mean
  else:
    res += sigma

  if sigma > mean:
    res += mean
  else:
    res += sigma
  return res
# Positional language model with KL-divergence scoring
# See report for details into this algorithm
# Due to performance considerations, we only re-rank other documents (from tfidf)
# Inputs:
# - query: query string
# - filtered_idx: tfidf top K document indices
# - Parameters: mu, kernel, kernel_cdf, sigma
# Returns:
# - Ranked list of document indices
# - Ranked list of corresponding scores
def plm(query, index, docs_idx, filtered_idx, doc_lengths, collection_size_nostop, term_counts, mu, kernel, kernel_cumul, sigma, rank_size=1000, sort_results = True):
  
  tokens = normalize(query, index)# query tokens
  token_ids = [token2id[token] for token in tokens if not token == '' and token in token2id]
  num_unique_tokens= len(set(token_ids)) 
  scores = np.zeros(len(filtered_idx)) # score for each doc
  
  # For each doc
  for i, d_idx in enumerate(filtered_idx):
    
    doc = docs_idx[d_idx]
    
    common_tokens = [t_id for t_id in token_ids if t_id in doc]
    
    q_len = len(common_tokens) # number of non-unique common tokens
   
    if doc_lengths[d_idx] == 0 or q_len == 0:
      scores[i] = -Inf
      continue
    
    unique_common_tokens = dict(collections.Counter(common_tokens))
    
    # Normalized counts of query words in document
    #q = np.array([unique_common_tokens[t] for t in unique_common_tokens]) # See proof for simplification on notepad
    
    # Query vector: uniform prior (division done later)
    q = np.ones(len(unique_common_tokens))
      
    # construct inverted index of positions of common terms in a document
    seq = index.document(d_idx)[1]
    inv_idx = {t_id:set() for t_id in unique_common_tokens}
    for pos,x in enumerate(seq):
      if x in inv_idx:
        inv_idx[x].add(pos+1)      
    
    # All positions where commons words occur. As an approximation, we only consider these positions
    qwpos = set()
    for key in inv_idx:
      qwpos = qwpos.union(inv_idx[key])
    
    p_mu = np.zeros(len(unique_common_tokens)) # probability vector for each word in query
    
    best_k = -1
    best_score = -Inf
    
    # k is the current position in the document
    
    for k in qwpos: # positions where query words occur
        
      Z_k = kernel_cumul(doc_lengths[d_idx],k,sigma) - kernel_cumul(0, k, sigma) # normalizing term
                   
      for _pos, tok in enumerate(unique_common_tokens):
        
        sum_kernels = 0
        for j in inv_idx[tok]:
          sum_kernels += kernel(k - j, sigma) 
        
        p_mu[_pos] =  (sum_kernels + mu * (term_counts[tok]/collection_size_nostop)) / (Z_k + mu)
      
      # Larger dot product is preferred
      kl = np.dot(q, np.log2(p_mu)) 
      
      if kl > best_score:
        best_k = k
        best_score = kl
        
    #print('Best Score: ' + str(best_score))
    scores[i] = best_score/ len(unique_common_tokens) + np.log2(len(token_ids))  
  
  if not sort_results:
    return scores
  
  sorted_idx, sorted_scores = rank_results(scores, index, rank_size)
  return sorted_idx, sorted_scores
## PLM on validation set
import time
results_plm = {}
mu = 2000 # best value from dirichlet smoothing
sigma = 50
rank_size=100 # for performance considerations
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.100']

kernels = [('circle', circle_kernel, circle_cdf), ('triangle', triangle_kernel, triangle_cdf), 
           ('cosine', cosine_kernel, cosine_cdf), ('passage',passage_kernel, passage_cdf), ('gaussian',gaussian_kernel, gaussian_cdf)]

tmp = {}
query_set = val_queries #val_queries

# Initialize a dictionary to store results 
for kernel_name, _, _ in kernels:
  results_plm[kernel_name] = np.zeros((len(query_set), len(metrics)))

# For each query
for i, qid in enumerate(query_set):
  
  # Rerank BM25 scores
  tfidf_idx, tfidf_scores = bm25(query_set[qid], index, docs_idx, doc_lengths, token2id, avg_doc_length, k1=1.2, b=0.75, rank_size=rank_size)
  #tfidf_idx, tfidf_scores = tfidf_scoring(query_set[qid], index, docs_idx, doc_lengths, token2id, rank_size)
  
    # for each kernel
  for kernel_name, kernel, kernel_cdf in kernels:

    start_time = time.time()
    
    idx, scores = plm(query_set[qid], index, docs_idx, tfidf_idx, doc_lengths, collection_size_nostop, id2tf, mu, kernel, kernel_cdf, sigma, rank_size)
    
    idx = [tfidf_idx[dum-1] for dum in idx] # index correction. ignore this!
    
    # trec eval
    result_ids = [index.document(dum)[0] for dum in idx]
    tmp[qid] = tuple((scores[dum], result_ids[dum]) for dum in range(len(idx)))
    code, res = run_trec_eval(tmp, 'PLM-'+kernel_name, metrics, groundtruth_qrel='./ap_88_89/qrel_validation', qrun='./tmp_run', max_results=rank_size)
    
    if code == 0:
      print('CODE ZERO',res)
      
    results_plm[kernel_name][i,:] = np.array(res)
    
    tmp.clear()
    
    if i % 5 == 0:
      print('- Query %d of %d: Kernel %s, elapsed time: %f' % (i+1,len(query_set), kernel_name, time.time() - start_time))
  #break
  
with open('./plm_perkernel_val_results', 'wb') as f:
  pickle.dump(f, results_plm)
# Plots kernel metrics

kn_list = []
ndcg_list = []

plm_res = np.zeros((len(query_set),len(kernels)))

print('')

i=0
for kname, _, _ in kernels:
  kn_list.append(kname)
  plm_res[:,i] = results_plm[kname][:,0]
  ndcg_list.append(np.mean(results_plm[kname], axis=0)[0])
  i+=1

plm_val = pd.DataFrame(plm_res, columns=kn_list)

sns.barplot(data = plm_val)
plt.ylabel('NDCG@10')
#plt.title('Performance of each kernel on NDCG@10 with 95% confidence intervals')
plt.savefig('./kernel_selection_validation_ndcg.png')

_ = sign_binomial_test(plm_res[:,2] - plm_res[:,1], 0, True, alternative='greater');
## PLM with Triangle Kernel on test set. 
import time

mu = 2000
sigma = 50
rank_size=1000
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']

kernels = [('triangle', triangle_kernel, triangle_cdf)]
tmp = {}

query_set = test_queries
results_plm = np.zeros((len(query_set), len(metrics)))

for i, qid in enumerate(query_set):
  
  # Rerank TFIDF scores
  tfidf_idx, tfidf_scores = tfidf_scoring(query_set[qid], index, docs_idx, doc_lengths, token2id, rank_size = 1000)
  
  for kernel_name, kernel, kernel_cdf in kernels:

    start_time = time.time()
    
    idx, scores = plm(query_set[qid], index, docs_idx, tfidf_idx, doc_lengths, collection_size_nostop, id2tf, mu, kernel, kernel_cdf, sigma, rank_size)
    
    idx = [tfidf_idx[dum-1] for dum in idx] # index correction
  
    result_ids = [index.document(dum)[0] for dum in idx]
     
    tmp[qid] = tuple((scores[dum], result_ids[dum]) for dum in range(len(idx)))

    code, res = run_trec_eval(tmp, 'PLM-'+kernel_name, metrics, groundtruth_qrel='./ap_88_89/qrel_test', qrun='./tmp_run', max_results=rank_size)
    
    if code == 0:
      print('CODE ZERO',res)
      
    results_plm[i,:] = np.array(res)
    
    tmp.clear()
    
    if i % 20 == 0:
      print('- Query %d of %d: Kernel %s, elapsed time: %f' % (i+1,len(query_set), kernel_name, time.time() - start_time))

with open('./plm_test_matrix', 'wb') as f:
  pickle.dump(results_plm, f)<jupyter_output><empty_output><jupyter_text>Plots<jupyter_code>metrics = ['NDCG@10', 'P@5', 'MAP@10', 'R@1000']

pkl_file = open('bm25_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

bm25_test = pd.DataFrame(mtrx, columns= ['BM25 ' + m for m in metrics])

pkl_file = open('ad_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

ad_test = pd.DataFrame(mtrx, columns=['AD ' + m for m in metrics])


pkl_file = open('dir_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

dir_test = pd.DataFrame(mtrx, columns=['DIR ' + m for m in metrics])


pkl_file = open('jm_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

jm_test = pd.DataFrame(mtrx, columns=['JM ' + m for m in metrics])


pkl_file = open('tfidf_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

tfidf_test = pd.DataFrame(mtrx, columns=['TFIDF ' + m for m in metrics])

pkl_file = open('plm_test_matrix', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

plm_test = pd.DataFrame(mtrx, columns=['PLM ' + m for m in metrics])

for i in range(len(metrics)):
  col = metrics[i]
  df2 = pd.concat([plm_test['PLM ' + col], jm_test['JM ' + col], ad_test['AD ' + col], bm25_test['BM25 ' + col], dir_test['DIR ' + col], tfidf_test['TFIDF ' + col]] , axis=1)
  figure(figsize=(11, 3))
  sns.barplot(data=df2)
  plt.show()<jupyter_output><empty_output><jupyter_text>Statistical Tests<jupyter_code># Statistical test on dirichlet and plm. For different values, load different pickled files (saved during each test run)
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']

pkl_file = open('plm_test_matrix', 'rb')
table1 = pickle.load(pkl_file)
pkl_file.close()

pkl_file = open('dir_test_matrix', 'rb')
table2 = pickle.load(pkl_file)
pkl_file.close()

# Computes the significance level applying the Šidák correction for having a final 5% alpha after a len(metrics) multiple comparisons test. https://en.wikipedia.org/wiki/%C5%A0id%C3%A1k_correction
#corr_alpha = 1-(1-0.05)**(1/len(metrics))

# Bonferroni correction
corr_alpha = 0.05/len(metrics)


for col in range(len(metrics)):
  print('Measure: ' + metrics[col])
  _ = sign_binomial_test(table2[:,col] - table1[:,col], 0, True, alternative='greater', sig = corr_alpha)
  print('')<jupyter_output><empty_output><jupyter_text>Task 2: Latent Semantic Models (LSMs) [25 points + 10 bonus points] In this task you will experiment with applying a distributional semantics methods ([word2vec](http://arxiv.org/abs/1411.2738)  **[5 points]**, [LSI](http://lsa3.colorado.edu/papers/JASIS.lsi.90.pdf) **[5 points]**, [LDA](https://www.cs.princeton.edu/~blei/papers/BleiNgJordan2003.pdf) **[5 points]** and [doc2vec](https://cs.stanford.edu/~quocle/paragraph_vector.pdf) **[5 points]**) for retrieval.You do not need to implement word2vec, LSI, LDA and doc2vec on your own. Instead, you can use [gensim](http://radimrehurek.com/gensim/index.html) (pre-loaded on the VirtualBox). An example on how to integrate Pyndri with Gensim for word2vec can be found [here](https://github.com/cvangysel/pyndri/blob/master/examples/word2vec.py). For the remaining latent vector space models, you will need to implement connector classes (such as `IndriSentences`) by yourself.In order to use a latent semantic model for retrieval, you need to:   * build a representation of the query **q**,   * build a representation of the document **d**,   * calculate the similarity between **q** and **d** (e.g., cosine similarity, KL-divergence).     The exact implementation here depends on the latent semantic model you are using. For example, in the case of word2vec, you only have vectors for individual words and not for documents or phrases. Try one of the following methods for producing these representations:   * Average or sum the word vectors.   * Cluster words in the document using [k-means](http://scikit-learn.org/stable/modules/generated/sklearn.cluster.KMeans.html) and use the centroid of the most important cluster. Experiment with different values of K for k-means.   * Using the [bag-of-word-embeddings representation](https://ciir-publications.cs.umass.edu/pub/web/getpdf.php?id=1248). **[10 bonus points]**   Each of these LSMs come with various hyperparameters to tune. Make a choice on the parameters, and explicitly mention the reasons that led you to these decisions. You can use the validation set to optimize hyper parameters you see fit; motivate your decisions. In addition, mention clearly how the query/document representations were constructed for each LSM and explain your choices.In this experiment, you will first obtain an initial top-1000 ranking for each query using TF-IDF in **Task 1**, and then re-rank the documents using the LSMs. Use TREC Eval to obtain the results and report on `NDCG@10`, Mean Average Precision (`MAP@1000`), `Precision@5` and `Recall@1000`.Perform significance testing **[5 points]** (similar as in Task 1) in the class of semantic matching methods.<jupyter_code>#Connector classes
import gensim, logging, bz2

# For Word2vec. Provides a token stream from a document
class IndriSentences(gensim.interfaces.CorpusABC):
    """Integrates an Index with Gensim's word2vec implementation."""

    def __init__(self, index, dictionary, max_documents=None):
        assert isinstance(index, pyndri.Index)

        self.index = index
        self.dictionary = dictionary

        self.max_documents = max_documents

    def _maximum_document(self):
        if self.max_documents is None:
            return self.index.maximum_document()
        else:
            return min(
                self.max_documents + self.index.document_base(),
                self.index.maximum_document())

    def __iter__(self):
        for int_doc_id in range(self.index.document_base(),
                                self._maximum_document()):
            ext_doc_id, tokens = self.index.document(int_doc_id)

            yield tuple(
                self.dictionary[token_id]
                for token_id in tokens
                if token_id > 0 and token_id in self.dictionary)

    def __len__(self):
        return self._maximum_document() - self.index.document_base()
      
# For Doc2Vec. Provides a token stream and document_id tuple
class IndriDocs(gensim.interfaces.CorpusABC):
    """Integrates an Index with Gensim's word2vec implementation."""

    def __init__(self, index, dictionary, max_documents=None):
        assert isinstance(index, pyndri.Index)

        self.index = index
        self.dictionary = dictionary

        self.max_documents = max_documents

    def _maximum_document(self):
        if self.max_documents is None:
            return self.index.maximum_document()
        else:
            return min(
                self.max_documents + self.index.document_base(),
                self.index.maximum_document())

    def __iter__(self):
        for int_doc_id in range(self.index.document_base(),
                                self._maximum_document()):
            ext_doc_id, tokens = self.index.document(int_doc_id)

            yield ext_doc_id, tuple(
                self.dictionary[token_id]
                for token_id in tokens
                if token_id > 0 and token_id in self.dictionary)

    def __len__(self):
        return self._maximum_document() - self.index.document_base()      
      
# For LDA/LSI: provides a BOW dictionary constructed as {token_id:count} for each document  
class IndriBow(gensim.interfaces.CorpusABC):
    """Integrates an Index with Gensim's word2vec implementation."""

    def __init__(self, index, docs_idx, max_documents=None):
        assert isinstance(index, pyndri.Index)

        self.index = index
        self.docs_bow = docs_idx
        self.max_documents = max_documents

    def _maximum_document(self):
        if self.max_documents is None:
            return self.index.maximum_document()
        else:
            return min(
                self.max_documents + self.index.document_base(),
                self.index.maximum_document())

    def __iter__(self):
        for int_doc_id in self.docs_bow:
            yield sorted(self.docs_bow[int_doc_id].items())

    def __len__(self):
        return self._maximum_document() - self.index.document_base()<jupyter_output><empty_output><jupyter_text>LDA/LSI<jupyter_code># Utility classes

# Gets document vectors from LDA/LSI
# Inputs: 
# - model: LSA/LSI
# - bow_dict: bag of words dictionary
def get_vector_for_bows(bow_dict, model):
  bow = [[k,v] for k, v in bow_dict.items()]
  topics_scores = model[bow]
  topics_idx = [topic_id for topic_id, score in topics_scores]
  topics_scores = [score for topic_id, score in topics_scores]
  
  vec = np.zeros(model.num_topics)
  vec[topics_idx] = topics_scores
  
  return vec

# Waterfall
def get_vector_for_tokenids(token_ids, model):
  return get_vector_for_bows(collections.Counter(token_ids), model)

# Tokenizes a string and returns vector representation
def get_vector_for_query(query, model,token2id):
  tokens = normalize(query, index)
  token_ids = [token2id[token] for token in tokens if not token == '' and token in token2id]
  return get_vector_for_tokenids(token_ids, model)


# returns negative cosine similarity (to reduce sorting overhead in numpy by avoiding inversion)
def get_query_doc_score(query, document, model, token2id):
  query_vec = get_vector_for_query(query, model, token2id)
  doc_vec = get_vector_for_bows(document, model)
  return - 1 + scipy.spatial.distance.cosine(query_vec, doc_vec)

# LDA/LSI scoring
# Inputs: 
# - query: query string
# - filtered_idx: list of document_ids
# - model: LSA or LSI
# Returns results ranked by cosine similarity scores
# - Ranked list of document ids
# - Corresponding positive cosine similarity scores
def lda_lsi_scoring(query, index, model, filtered_idx, docs_idx, token2id, docs_lengths):
  
  scores = np.zeros(len(filtered_idx))
  
  for i, d_idx in enumerate(filtered_idx):
    doc_bow = docs_idx[d_idx]
    
    if docs_lengths[d_idx] == 0:
      scores[i] = Inf
      continue
    
    scores[i] = get_query_doc_score(query, doc_bow, model, token2id)
    
  sorted_filtered_idx = np.argsort(scores)
  orig_sorted_idx = [filtered_idx[x] for x in sorted_filtered_idx]
  orig_sorted_scores = [-scores[x] for x in sorted_filtered_idx]
  
  return orig_sorted_idx, orig_sorted_scores

# Calls TREC Eval
def lda_lsi_query_vals (model, queries, qrel_fname, metrics, index, docs_idx, doc_lengths, token2id, out_fname):
  
  results = np.zeros((len(queries),len(metrics)))
  tmp = {}
  
  for i, qid in enumerate(queries):
    tmp.clear()
    print('Query %d out of %d: %s' % (i+1, len(queries), queries[qid]))
  
    # TFIDF reranking
    bm25_idx, bm25_scores = tfidf_scoring(queries[qid], index, docs_idx, doc_lengths, token2id, rank_size = 1000)
    idx, scores = lda_lsi_scoring(queries[qid], index, model, bm25_idx, docs_idx, token2id, doc_lengths)
    
    result_ids = [index.document(i)[0] for i in idx]
    tmp[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx)))

    code, res = run_trec_eval(tmp, str(model), metrics, groundtruth_qrel=qrel_fname, qrun='./tmp_run', max_results=1000)
    if code == 0:
      print('CODE ZERO',res)
    results[i,:] = np.array(res)
  
  with open(out_fname,'wb') as f:
    pickle.dump(results,f)
  
  return results<jupyter_output><empty_output><jupyter_text>Training<jupyter_code>## Train LDA
num_topics = 10
eta = 1.0 / num_topics # uniform priors
_alpha = 'symmetric'
lda = gensim.models.ldamulticore.LdaMulticore(
    corpus=IndriBow(index,docs_idx), num_topics=num_topics, 
    id2word=id2token, workers=3, chunksize=5000, passes=1, 
    batch=False, alpha=_alpha, eta=eta, eval_every=1000, iterations=5, minimum_probability=0.1)

lda.save('./lda')
## Train LSI
num_topics = 100
lsi = gensim.models.lsimodel.LsiModel(corpus=IndriBow(index,docs_idx), id2word=id2token, num_topics=num_topics)
lsi.save('./lsi')<jupyter_output><empty_output><jupyter_text>Reported results on test set<jupyter_code>## LSI test reuslts
lsi = gensim.models.lsimodel.LsiModel.load('./lsi')
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
results = lda_lsi_query_vals(lsi, test_queries, './ap_88_89/qrel_test', metrics, index, docs_idx, doc_lengths, token2id, './lsi_test_results')
print(results)
## LDA test reuslts
lda = gensim.models.LdaModel.load('./lda')
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
results = lda_lsi_query_vals(lda, test_queries, './ap_88_89/qrel_test', metrics, index, docs_idx, doc_lengths, token2id, './lda_test_results')
print(results)<jupyter_output><empty_output><jupyter_text>Word2Vec / Doc2Vec<jupyter_code>#  Gets a vector embedding for a document
# Inputs:
# - model: Word2Vec, Doc2Vec
# - doc_bow: bag-of-words dictionary of document {token_id: count}
def get_docbow_embedding(doc_bow, model, doc_len, doc_idx=None ):
  
  if isinstance(model, gensim.models.Doc2Vec) and doc_idx is not None:
    return model.docvecs[index.document(doc_idx)[0]]
  
  elif isinstance(model, gensim.models.Word2Vec) or  doc_idx is None: # Average Word2Vec vectors
    res = np.zeros(model.vector_size)

    for key, val in doc_bow.items():
      if id2token[key] in model:
        res += val * model[id2token[key]]

    return res/doc_len

# Word2Vec / Doc2vec scoring
# Inputs: 
# - query: query string
# - filtered_idx: list of document_ids
# - model: Word2Vec or Doc2vec
# Returns results ranked by cosine similarity scores
# - Ranked list of document ids
# - Corresponding positive cosine similarity scores
def w2v_d2v_scoring(query, index, model, filtered_idx, docs_idx, doc_lengths):
  
  tokens = [token2id[t] for t in normalize(query, index) if not t == '' and t in token2id]
  query_bow = collections.Counter(tokens)
  query_vec = get_docbow_embedding(query_bow, model, len(tokens))
  
  scores = np.zeros(len(filtered_idx))   
  
  for i, d_idx in enumerate(filtered_idx):
    
    doc_bow = docs_idx[d_idx]    
    
    if doc_lengths[d_idx] == 0:
      scores[i] = Inf
      continue
      
    doc_vec = get_docbow_embedding(doc_bow, model, doc_lengths[d_idx], None)
    scores[i] = -1 + scipy.spatial.distance.cosine(query_vec, doc_vec)  #inverted sign for not reversing list after argsort

  
  sorted_filtered_idx = np.argsort(scores)
  orig_sorted_idx = [filtered_idx[x] for x in sorted_filtered_idx]
  orig_sorted_scores = [-scores[x] for x in sorted_filtered_idx]
  
  return orig_sorted_idx, orig_sorted_scores

# Calls trec eval
def w2v_d2v_eval (model, queries, qrel_fname, metrics, index, docs_idx, doc_lengths, token2id, out_fname):
  
  results = np.zeros((len(queries),len(metrics)))
  tmp = {}
  
  for i, qid in enumerate(queries):
    tmp.clear()
    print('Query %d out of %d: %s' % (i+1, len(queries), queries[qid]))
  
    #TFIDF reranking 
    bm25_idx, bm25_scores = tfidf_scoring(queries[qid], index, docs_idx, doc_lengths, token2id, rank_size = 1000)
    idx, scores = w2v_d2v_scoring(queries[qid], index, model, bm25_idx, docs_idx, doc_lengths)
    
    result_ids = [index.document(i)[0] for i in idx]
    tmp[qid] = tuple((scores[i], result_ids[i]) for i in range(len(idx)))

    code, res = run_trec_eval(tmp, str(model), metrics, groundtruth_qrel=qrel_fname, qrun='./tmp_run', max_results=1000)
    if code == 0:
      print('CODE ZERO',res)
    results[i,:] = np.array(res)
  
  with open(out_fname,'wb') as f:
    pickle.dump(results,f)
  
  return results<jupyter_output><empty_output><jupyter_text>Training<jupyter_code># Train word2vec
print('Loading vocabulary.')
sentences = IndriSentences(index, id2token)
print('Initializing word2vec.')

w2v_model = gensim.models.Word2Vec(sentences,
    size=16,  # Embedding size
    window=5,  # One-sided window size
    sg=True,  # Skip-gram.
    min_count=5,  # Minimum word frequency.
    sample=1e-3,  # Sub-sample threshold.
    hs=False,  # Hierarchical softmax.
    negative=10,  # Number of negative examples.
    iter=1,  # Number of iterations.
    workers=3,  # Number of workers.
)
#Train Doc2Vec
from gensim.models.doc2vec import TaggedDocument

py_docs = IndriDocs(index,id2token)

doc_db = [TaggedDocument(words=text, tags=[docid]) for docid, text in py_docs]

doc2vec_model = gensim.models.Doc2Vec(doc_db, size=16, window=5, min_count=1, workers=3)<jupyter_output><empty_output><jupyter_text>Experiments with pre-trained word/document vectors<jupyter_code>## Word2Vec test results
pretrained_word2vec = gensim.models.Word2Vec.load('./word2vec.bin')
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
results = w2v_d2v_eval(pretrained_word2vec, test_queries, './ap_88_89/qrel_test', metrics, index, docs_idx, doc_lengths, token2id, './w2v_test_results')
print(results)
## Doc2Vec test reuslts
pretrained_doc2vec = gensim.models.Doc2Vec.load('./doc2vec.bin')
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.1000']
results = w2v_d2v_eval(pretrained_doc2vec, test_queries, './ap_88_89/qrel_test', metrics, index, docs_idx, doc_lengths, token2id, './d2v_test_results')
print(results)<jupyter_output><empty_output><jupyter_text>Plots (LDA, LSI, W2V, D2V)<jupyter_code>## Plotting

metrics = ['NDCG@10', 'P@5', 'MAP@10', 'R@1000']

pkl_file = open('lda_test_results', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

lda_test = pd.DataFrame(mtrx, columns= ['LDA ' + m for m in metrics])

pkl_file = open('lsi_test_results', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

lsi_test = pd.DataFrame(mtrx, columns=['LSI ' + m for m in metrics])


pkl_file = open('w2v_test_results', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

w2v_test = pd.DataFrame(mtrx, columns=['W2V ' + m for m in metrics])


pkl_file = open('d2v_test_results', 'rb')
mtrx = pickle.load(pkl_file)
pkl_file.close()

d2v_test = pd.DataFrame(mtrx, columns=['D2V ' + m for m in metrics])

for i in range(len(metrics)):
  col = metrics[i]
  df2 = pd.concat([lda_test['LDA ' + col], lsi_test['LSI ' + col], w2v_test['W2V ' + col], d2v_test['D2V ' + col]] , axis=1)
  figure(figsize=(8, 4))
  sns.barplot(data=df2)
  plt.savefig('./semantic_%s.png' % col)
  plt.show()<jupyter_output><empty_output><jupyter_text>Statistical Tests (LDA, LSI, W2V, D2V)<jupyter_code>#%load_ext memory_profiler
#%memit
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10'] # , 'recall.1000'

pkl_file = open('lda_test_results', 'rb')
lda_data = pickle.load(pkl_file)
pkl_file.close()

pkl_file = open('lsi_test_results', 'rb')
lsi_data = pickle.load(pkl_file)
pkl_file.close()

pkl_file = open('w2v_test_results', 'rb')
w2v_data = pickle.load(pkl_file)
pkl_file.close()

pkl_file = open('d2v_test_results', 'rb')
d2v_data = pickle.load(pkl_file)
pkl_file.close()

# Computes the significance level applying the Šidák correction for having a final 5% alpha after a len(metrics) multiple comparisons test. https://en.wikipedia.org/wiki/%C5%A0id%C3%A1k_correction
corr_alpha = 1-(1-0.05)**(1/len(metrics))
print(corr_alpha)
corr_alpha = 0.05/len(metrics)
print(corr_alpha)

def compute_sb_test(table1, table2, corr_alpha):
  print(table1.shape, table2.shape)
  for col in range(len(metrics)):
    print('Measure: ' + metrics[col])
    _ = sign_binomial_test(table1[:,col] - table2[:,col], 0, True, alternative='greater', sig = corr_alpha)
  print('')

print('LDA > LSI')
compute_sb_test(lda_data, lsi_data, corr_alpha)
print('LDA > W2V')
compute_sb_test(lda_data, w2v_data, corr_alpha)
print('LDA > D2V')
compute_sb_test(lda_data, d2v_data, corr_alpha)

print('LSI > LDA')
compute_sb_test(lsi_data, lda_data, corr_alpha)
print('LSI > W2V')
compute_sb_test(lsi_data, w2v_data, corr_alpha)
print('LSI > D2V')
compute_sb_test(lsi_data, d2v_data, corr_alpha)


print('W2V > LDA')
compute_sb_test(w2v_data, lda_data, corr_alpha)
print('W2V > LSI')
compute_sb_test(w2v_data, lsi_data, corr_alpha)
print('W2V > D2V')
compute_sb_test(w2v_data, d2v_data, corr_alpha)


print('D2V > LDA')
compute_sb_test(d2v_data, lda_data, corr_alpha)
print('D2V > W2V')
compute_sb_test(d2v_data, w2v_data, corr_alpha)
print('D2V > LSI')
compute_sb_test(d2v_data, lsi_data, corr_alpha)<jupyter_output><empty_output><jupyter_text>Task 3: Learning to rank (LTR) [10 points] In this task you will get an introduction into learning to rank for information retrieval. You will experiment with a pointwise learning to rank method, logistic regression, implemented in [scikit-learn](http://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LogisticRegression.html).**NOTE**: you can only perform this task if you have completely finished Task 1 and Task 2.In this experiment, you will use the retrieval methods you implemented in Task 1 and Task 2 as features for the learning to rank model. Train your LTR model using 10-fold cross validation on the test set. For every query, first create a document candidate set using the top-1000 documents using TF-IDF. Secondly, compute query-document values using the retrieval models above and use them as features. Note that the feature values of different retrieval methods are likely to be distributed differently.Your approach will definitely not be as good as the state-of-the-art since you are taking a pointwise approach, but we do not ask you to try pair- or listwise methods because they will be the main topic of the next assignment.<jupyter_code># Retrieves features (ranking scores) for a query and a single document
# Returns a vector of size 10
# Inputs:
# - qid: query id
# - filtered_doc_idx: id of document
# - query_set: dictionary of queries
def get_features_for_query(qid, filtered_doc_idx, query_set, index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id):
  
    # Contains the features of the current (d,q)
    features = []
    
    # Adds query identifier for filtering later
    features.append(int(qid))
    
    # Calculates score for 
    score =  tfidf_scoring(query_set[qid], index, mini_docs_idx, doc_lengths, token2id, rank_size = 1, sort_results = False)[0]
    features.append(score)
    
    score = abs_discounting(query_set[qid], index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, 0.8, rank_size = 1, sort_results = False)[0]
    features.append(score)
    
    score = jelinek_mercer(query_set[qid], index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, 0.4, rank_size = 1, sort_results = False)[0]
    features.append(score)
    
    score = bm25(query_set[qid], index, mini_docs_idx, doc_lengths, token2id, avg_doc_length, k1=1.2, b=0.75, rank_size = 1, sort_results = False)[0]
    features.append(score)
    
    score = dirichlet_smoothing(query_set[qid], index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, 2000, rank_size = 1, sort_results = False)[0]
    features.append(score)
    
    score = plm(query_set[qid], index, mini_docs_idx, [filtered_doc_idx], doc_lengths, collection_size_nostop, id2tf, 2000, triangle_kernel, triangle_cdf,\
                50, rank_size=1, sort_results = False)[0]
    features.append(score)
   
    _, score = lda_lsi_scoring(query_set[qid], index, lda, [filtered_doc_idx], mini_docs_idx, token2id, doc_lengths)
    features.append(score[0])
   
    _, score = lda_lsi_scoring(query_set[qid], index, lsi, [filtered_doc_idx], mini_docs_idx, token2id, doc_lengths)
    features.append(score[0])
   
    _, score = w2v_d2v_scoring(query_set[qid], index, pretrained_word2vec, [filtered_doc_idx], mini_docs_idx, doc_lengths)
    features.append(score[0])
    
    _, score = w2v_d2v_scoring(query_set[qid], index, pretrained_doc2vec, [filtered_doc_idx], mini_docs_idx, doc_lengths)
    features.append(score[0])
    
    return features  
# Creates features for each query-document pair in the test set
def create_features(query_set, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id):
  
  # List of lists of features
  res_feat = []
  
  # Lists of targets (1:R, 0:NR)
  res_tgt = []
  
  i = 0

  f = open('./ap_88_89/qrel_test', 'r')
  for line in f:
    
    sp_line = line.split()
    
    qid = sp_line[0]
    
    _idx_tpl = index.document_ids([sp_line[2]])
    
    # Some documents dont have idx, so the previous tuple is empty
    if len(_idx_tpl) == 0:
      continue
      
    # Takes the document id from the tuple
    idx = _idx_tpl[0][1]
    
    # Creates a mini docs_idx with just the document in the current line
    mini_docs_idx = {idx: docs_idx[idx]}
    
    # After constructing the features list append it to res_feat
    res_feat.append(get_features_for_query(qid, idx, query_set, index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id))
    
    rel_label = float(sp_line[-1])
    res_tgt.append(rel_label)
    
    if (i+1) % 1000 == 0:
      print('Q ' + str(i))
      
    i += 1
    
  res_feat = np.array(res_feat)
  res_tgt = np.array(res_tgt)
  
  return res_feat, res_tgt

# Generate and store training data for LTR
res_feat, res_tgt = create_features(test_queries, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id)
with open('./ltr_training','wb') as f:
  pickle.dump(res_feat, f)
  pickle.dump(res_tgt, f)
from sklearn.linear_model import LogisticRegression
from sklearn import preprocessing
from sklearn.model_selection import KFold
from sklearn.neural_network import MLPClassifier
from sklearn.decomposition import PCA
import sklearn

# Trains a binary classifier (logistic, linear or MLP) with k-fold cross validation
# Inputs:
# - X_raw: unscaled dataframe of data, first columns is query id
# - t: targets (binary)
# - cv_k: number of folds (k) in cross validation 
# - model: 
#    - log: logistic regression
#    - lin: linear regression
#    - nn: multi layer perceptron
# - scale: flag for normalizing data
# Returns
# - List of models, their train/test queries, and their scalers
def train_ltr_models_cv(X_raw, t, cv_k = 10, model = 'log', scale=False):
  
  # List of unique query ids
  query_ids = X_raw['query'].unique()
  data_cols = list(X_raw.drop(['query'], axis=1).columns) 
  ltr_models = [] 
  
  # Do PCA and save transformations
  X = pd.DataFrame.copy(X_raw)
  scaler = None
  if scale:
    scaler = preprocessing.MinMaxScaler().fit(X_raw.drop(['query'], axis=1)) # exclude query ids columns in dataframe
    X[data_cols] = scaler.transform(X_raw[data_cols])
  
  
  # Do Kfold split
  kf = KFold(n_splits=cv_k)
  i = 0
  for train, test in kf.split(query_ids): # perform split on query ids
    
    q_tr = [query_ids[x] for x in train] # training queries
    q_te = [query_ids[x] for x in test]  # testing queries, used for reporting test accuracy and storing for bookkeeping
    
    print('Training model %d out of %d'% (i+1, cv_k))
    i+=1
    
    # Slice up data for k-fold train/test split
    tr_idx = X['query'].isin(q_tr) # filter by queries
    X_train = X[tr_idx].drop(['query'], axis=1)
    t_train = t[tr_idx]
    
    te_idx = X['query'].isin(q_te) # filter by queries
    X_test = X[te_idx].drop(['query'], axis=1)
    t_test = t[te_idx]
    
    
    # Initialize model
    if model == 'log':
      ltr = LogisticRegression(tol=0.00001, fit_intercept=True, solver='sag',max_iter=5000, n_jobs=3, random_state=1)
    elif model == 'lin':
      ltr = sklearn.linear_model.LinearRegression()
    else:
      ltr = MLPClassifier(activation = 'logistic', solver='adam', alpha=1e-5, max_iter=500,  tol=1e-6, nesterovs_momentum = True
                          , early_stopping=True, verbose=True, random_state=1, hidden_layer_sizes=(30))
    
    # Fit training data
    ltr.fit(X_train.as_matrix(), t_train.as_matrix().ravel())
    
    # Repotr test accuracy
    print('Model %d '%(i),ltr.score(X_test, t_test))

    # Save
    ltr_models.append({'model': ltr, 'tr_query': q_tr, 'te_query':q_te})
    
  return ltr_models, scaler
# Train models

# Load up training data
with open('./ltr_training','rb') as f:
  res_feat = pickle.load( f)
  res_tgt = pickle.load( f)

# Store in dataframes    
df_X = pd.DataFrame(res_feat, columns =['query','tfidf','ad','jm','bm25','dirichlet','plm','w2v', 'd2v'])
df_t = pd.DataFrame(res_tgt,columns=['target'])
print(df_X.shape, df_t.shape)
del res_feat
del res_tgt


# Clean data by replacing infs (due to design)
df_X['tfidf'].replace(-inf, -200, inplace=True)
df_X['plm'].replace(-inf, -100, inplace=True)

ltr_models, scaler = train_ltr_models_cv(df_X, df_t, cv_k = 10, model = 'log', scale=False)

with open('kfold_cv_results_log','wb') as f:
  pickle.dump(ltr_models,f)
  pickle.dump(scaler,f)

# Run Trec Eval
def run_ltr_evals(ltr_models, scaler, queries, qrel_fname, metrics, index, docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id, out_fname, rank_size):
  
  tmp = {}
  
  for k, model_settings in enumerate(ltr_models):
    
    print('Evaluating LTR model %d' % (k+1))
    
    model = model_settings['model']
    test_queries_cv = model_settings['te_query']
    test_queries_cv = [str(int(a)) for a in test_queries_cv]
    #qid_ls = [ls_test_query_ids[int(_idx)] for _idx in test_queries_cv]
    #print(qid_ls)
    
    for o,qid in enumerate(test_queries_cv):
      
      print('Query %d out of %d'% (o+1,len(test_queries_cv)))
      
    # TFIDF re-ranking
      tfidf_idx, tfidf_scores = tfidf_scoring(queries[qid], index, docs_idx, doc_lengths, token2id, rank_size)
      
      features = []
      
        # features for each document
      for _dox in tfidf_idx:
        mini_docs_idx = {_dox: docs_idx[_dox]}
        f = get_features_for_query(qid, _dox, queries, index, mini_docs_idx, doc_lengths, collection_size_nostop, id2tf, token2id)        
        features.append(f)
      
      np_arr = numpy.array(features)
      np_arr[:,1][np_arr[:,1] <= -float('Inf')] = -200
      np_arr[:,6][np_arr[:,6] <= -float('Inf')] = -100
      
        #scale if necessary
      if scaler:
        np_arr[:,1:] = scaler.transform(np_arr[:,1:])
      
        # predict
      preds = model.predict_proba(np_arr[:,1:])
      
        # sort predictions
      sort_idx = np.argsort(preds[:,0])[::-1]
      
        # store sorted results for trec eval
      tmp[qid] = tuple((preds[i,0], index.document(tfidf_idx[i])[0]) for i in sort_idx)
      #print(qid, np_arr.shape)
      
      del np_arr
      del features
    
    # run trev eval
  code, res = run_trec_eval(tmp, 'LTR', metrics, groundtruth_qrel=qrel_fname, qrun='./tmp_run', max_results=rank_size)
  print(res)
  
  if code == 0:
      print('CODE ZERO',res)
      
  with open(out_fname,'wb') as f:
    pickle.dump(tmp,f)
    pickle.dump(res,f)
# LTR performance on test data
with open('kfold_cv_results_log_nolsilda','rb') as f:
  ltr_models = pickle.load(f)
  scaler = pickle.load(f)

# we only consider re-rank top 100 from tfidf for performance considerations
metrics = ['ndcg_cut.10', 'P.5', 'map_cut.10', 'recall.100']  

# load models into memory
lda = gensim.models.LdaMulticore.load('./lda')
lsi = gensim.models.LsiModel.load('./lsi')
pretrained_word2vec = gensim.models.Word2Vec.load('./word2vec.bin')
pretrained_doc2vec = gensim.models.Doc2Vec.load('./doc2vec.bin')


run_ltr_evals(ltr_models, scaler, queries, './ap_88_89/qrel_test', metrics, index, docs_idx, \
              doc_lengths, collection_size_nostop, id2tf, token2id, './ltr_treceval_results_log_nolsilda', 100)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11687256596084905, 'fraction_numerical': 0.019265984061073427, 'mean_word_length': 3.969050365244137, 'pattern_counts': {'":': 0, '<': 134, '<?xml version=': 0, '>': 153, 'https://': 17, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '28', 'dir': 'jupyter-structured-clean-dedup', 'id': '564178', 'license': 'MIT', 'n_tokens_mistral': 28187, 'n_tokens_neox': 26722, 'n_words': 7522, 'path': 'Lab2/11391014-11390689-hw2.ipynb', 'repo_name': 'danakianfar/information_retrieval', 'size': '77595'}