{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Source of the materials**: Biopython cookbook (adapted)\n",
    "<font color='red'>Status: Draft</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<font color='red'>Consider looking at the Substitution Matrices example in [Chapter 19 - Cookbook](19 - Cookbook - Cool things to do with it.ipynb#Substitution-Matrices)</font>\n",
    "\n",
    "# Advanced\n",
    "\n",
    "[Parser Design](#Parser-Design)\n",
    "\n",
    "[Substitution Matrices](#Substitution-Matrices)\n",
    "\n",
    "[FreqTable](#FreqTable)\n",
    "\n",
    "## Parser Design\n",
    "\n",
    "Many of the older Biopython parsers were built around an event-oriented\n",
    "design that includes Scanner and Consumer objects.\n",
    "\n",
    "Scanners take input from a data source and analyze it line by line,\n",
    "sending off an event whenever it recognizes some information in the\n",
    "data. For example, if the data includes information about an organism\n",
    "name, the scanner may generate an `organism_name` event whenever it\n",
    "encounters a line containing the name.\n",
    "\n",
    "Consumers are objects that receive the events generated by Scanners.\n",
    "Following the previous example, the consumer receives the\n",
    "`organism_name` event, and the processes it in whatever manner necessary\n",
    "in the current application.\n",
    "\n",
    "This is a very flexible framework, which is advantageous if you want to\n",
    "be able to parse a file format into more than one representation. For\n",
    "example, the `Bio.GenBank` module uses this to construct either\n",
    "`SeqRecord` objects or file-format-specific record objects.\n",
    "\n",
    "More recently, many of the parsers added for `Bio.SeqIO` and\n",
    "`Bio.AlignIO` take a much simpler approach, but only generate a single\n",
    "object representation (`SeqRecord` and `MultipleSeqAlignment` objects\n",
    "respectively). In some cases the `Bio.SeqIO` parsers actually wrap\n",
    "another Biopython parser - for example, the `Bio.SwissProt` parser\n",
    "produces SwissProt format specific record objects, which get converted\n",
    "into `SeqRecord` objects.\n",
    "\n",
    "## Substitution Matrices\n",
    "\n",
    "\n",
    "### SubsMat\n",
    "\n",
    "This module provides a class and a few routines for generating\n",
    "substitution matrices, similar to BLOSUM or PAM matrices, but based on\n",
    "user-provided data. Additionally, you may select a matrix from\n",
    "MatrixInfo.py, a collection of established substitution matrices. The\n",
    "`SeqMat` class derives from a dictionary:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "class SeqMat(dict)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "The dictionary is of the form `{(i1,j1):n1, (i1,j2):n2,...,(ik,jk):nk}`\n",
    "where i, j are alphabet letters, and n is a value.\n",
    "\n",
    "1.  Attributes\n",
    "\n",
    "    1.  `self.alphabet`: a class as defined in Bio.Alphabet\n",
    "\n",
    "    2.  `self.ab_list`: a list of the alphabet’s letters, sorted. Needed\n",
    "        mainly for internal purposes\n",
    "\n",
    "2.  Methods\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "__init__(self,data=None,alphabet=None, mat_name='', build_later=0):\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        1.  `data`: can be either a dictionary, or another\n",
    "            SeqMat instance.\n",
    "\n",
    "        2.  `alphabet`: a Bio.Alphabet instance. If not provided,\n",
    "            construct an alphabet from data.\n",
    "\n",
    "        3.  `mat_name`: matrix name, such as “BLOSUM62” or “PAM250”\n",
    "\n",
    "        4.  `build_later`: default false. If true, user may supply only\n",
    "            alphabet and empty dictionary, if intending to build the\n",
    "            matrix later. this skips the sanity check of alphabet\n",
    "            size vs. matrix size.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "entropy(self,obs_freq_mat)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        1.  `obs_freq_mat`: an observed frequency matrix. Returns the\n",
    "            matrix’s entropy, based on the frequency in `obs_freq_mat`.\n",
    "            The matrix instance should be LO or SUBS.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "sum(self)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        Calculates the sum of values for each letter in the matrix’s\n",
    "        alphabet, and returns it as a dictionary of the form\n",
    "        `{i1: s1, i2: s2,...,in:sn}`, where:\n",
    "\n",
    "        -   i: an alphabet letter;\n",
    "\n",
    "        -   s: sum of all values in a half-matrix for that letter;\n",
    "\n",
    "        -   n: number of letters in alphabet.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "print_mat(self,f,format=\"%4d\",bottomformat=\"%4s\",alphabet=None)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        prints the matrix to file handle f. `format` is the format field\n",
    "        for the matrix values; `bottomformat` is the format field for\n",
    "        the bottom row, containing matrix letters. Example output for a\n",
    "        3-letter alphabet matrix:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "A 23\n",
    "        B 12 34\n",
    "        C 7  22  27\n",
    "          A   B   C\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        The `alphabet` optional argument is a string of all characters\n",
    "        in the alphabet. If supplied, the order of letters along the\n",
    "        axes is taken from the string, rather than by\n",
    "        alphabetical order.\n",
    "\n",
    "3.  Usage\n",
    "\n",
    "    The following section is laid out in the order by which most people\n",
    "    wish to generate a log-odds matrix. Of course, interim matrices can\n",
    "    be generated and investigated. Most people just want a log-odds\n",
    "    matrix, that’s all.\n",
    "\n",
    "    1.  Generating an Accepted Replacement Matrix\n",
    "\n",
    "        Initially, you should generate an accepted replacement\n",
    "        matrix (ARM) from your data. The values in ARM are the counted\n",
    "        number of replacements according to your data. The data could be\n",
    "        a set of pairs or multiple alignments. So for instance if\n",
    "        Alanine was replaced by Cysteine 10 times, and Cysteine by\n",
    "        Alanine 12 times, the corresponding ARM entries would be:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "('A','C'): 10, ('C','A'): 12\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        as order doesn’t matter, user can already provide only one\n",
    "        entry:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "('A','C'): 22\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        A SeqMat instance may be initialized with either a full (first\n",
    "        method of counting: 10, 12) or half (the latter method, 22)\n",
    "        matrices. A full protein alphabet matrix would be of the size\n",
    "        20x20 = 400. A half matrix of that alphabet would be 20x20/2 +\n",
    "        20/2 = 210. That is because same-letter entries don’t change.\n",
    "        (The matrix diagonal). Given an alphabet size of N:\n",
    "\n",
    "        1.  Full matrix size: N\\*N\n",
    "\n",
    "        2.  Half matrix size: N(N+1)/2\n",
    "\n",
    "        The SeqMat constructor automatically generates a half-matrix, if\n",
    "        a full matrix is passed. If a half matrix is passed, letters in\n",
    "        the key should be provided in alphabetical order: (’A’,’C’) and\n",
    "        not (’C’,A’).\n",
    "\n",
    "        At this point, if all you wish to do is generate a log-odds\n",
    "        matrix, please go to the section titled Example of Use. The\n",
    "        following text describes the nitty-gritty of internal functions,\n",
    "        to be used by people who wish to investigate their\n",
    "        nucleotide/amino-acid frequency data more thoroughly.\n",
    "\n",
    "    2.  Generating the observed frequency matrix (OFM)\n",
    "\n",
    "        Use:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```python\n",
    "OFM = SubsMat._build_obs_freq_mat(ARM)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        The OFM is generated from the ARM, only instead of replacement\n",
    "        counts, it contains replacement frequencies.\n",
    "\n",
    "    3.  Generating an expected frequency matrix (EFM)\n",
    "\n",
    "        Use:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "EFM = SubsMat._build_exp_freq_mat(OFM,exp_freq_table)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        1.  `exp_freq_table`: should be a FreqTable instance. See\n",
    "            section \\[sec:freq\\_table\\] for detailed information\n",
    "            on FreqTable. Briefly, the expected frequency table has the\n",
    "            frequencies of appearance for each member of the alphabet.\n",
    "            It is implemented as a dictionary with the alphabet letters\n",
    "            as keys, and each letter’s frequency as a value. Values sum\n",
    "            to 1.\n",
    "\n",
    "        The expected frequency table can (and generally should) be\n",
    "        generated from the observed frequency matrix. So in most cases\n",
    "        you will generate `exp_freq_table` using:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "```\n",
    "from Bio import SubsMat\n",
    "from Bio.SubsMat import _build_obs_freq_mat\n",
    "OFM = _build_obs_freq_mat(ARM)\n",
    "exp_freq_table = SubsMat._exp_freq_table_from_obs_freq(OFM)\n",
    "EFM = SubsMat._build_exp_freq_mat(OFM, exp_freq_table)\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        But you can supply your own `exp_freq_table`, if you wish\n",
    "\n",
    "    4.  Generating a substitution frequency matrix (SFM)\n",
    "\n",
    "        Use:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "SFM = SubsMat._build_subs_mat(OFM,EFM)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        Accepts an OFM, EFM. Provides the division product of the\n",
    "        corresponding values.\n",
    "\n",
    "    5.  Generating a log-odds matrix (LOM)\n",
    "\n",
    "        Use:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "LOM=SubsMat._build_log_odds_mat(SFM[,logbase=10,factor=10.0,round_digit=1])\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "        1.  Accepts an SFM.\n",
    "\n",
    "        2.  `logbase`: base of the logarithm used to generate the\n",
    "            log-odds values.\n",
    "\n",
    "        3.  `factor`: factor used to multiply the log-odds values. Each\n",
    "            entry is generated by log(LOM\\[key\\])\\*factor And rounded to\n",
    "            the `round_digit` place after the decimal point,\n",
    "            if required.\n",
    "\n",
    "4.  Example of use\n",
    "\n",
    "    As most people would want to generate a log-odds matrix, with\n",
    "    minimum hassle, SubsMat provides one function which does it all:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "make_log_odds_matrix(acc_rep_mat,exp_freq_table=None,logbase=10,\n",
    "                          factor=10.0,round_digit=0):\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "    1.  `acc_rep_mat`: user provided accepted replacements matrix\n",
    "\n",
    "    2.  `exp_freq_table`: expected frequencies table. Used if provided,\n",
    "        if not, generated from the `acc_rep_mat`.\n",
    "\n",
    "    3.  `logbase`: base of logarithm for the log-odds matrix. Default\n",
    "        base 10.\n",
    "\n",
    "    4.  `round_digit`: number after decimal digit to which result should\n",
    "        be rounded. Default zero.\n",
    "\n",
    "## FreqTable\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "FreqTable.FreqTable(UserDict.UserDict)\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "1.  Attributes:\n",
    "\n",
    "    1.  `alphabet`: A Bio.Alphabet instance.\n",
    "\n",
    "    2.  `data`: frequency dictionary\n",
    "\n",
    "    3.  `count`: count dictionary (in case counts are provided).\n",
    "\n",
    "2.  Functions:\n",
    "\n",
    "    1.  `read_count(f)`: read a count file from stream f. Then convert\n",
    "        to frequencies.\n",
    "\n",
    "    2.  `read_freq(f)`: read a frequency data file from stream f. Of\n",
    "        course, we then don’t have the counts, but it is usually the\n",
    "        letter frequencies which are interesting.\n",
    "\n",
    "3.  Example of use: The expected count of the residues in the database\n",
    "    is sitting in a file, whitespace delimited, in the following format\n",
    "    (example given for a 3-letter alphabet):\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "A   35\n",
    "    B   65\n",
    "    C   100\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "    And will be read using the\n",
    "    `FreqTable.read_count(file_handle)` function.\n",
    "\n",
    "    An equivalent frequency file:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "A  0.175\n",
    "    B  0.325\n",
    "    C  0.5\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "    Conversely, the residue frequencies or counts can be passed as\n",
    "    a dictionary. Example of a count dictionary (3-letter alphabet):\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "{'A': 35, 'B': 65, 'C': 100}\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "    Which means that an expected data count would give a 0.5 frequency\n",
    "    for ’C’, a 0.325 probability of ’B’ and a 0.175 probability of ’A’\n",
    "    out of 200 total, sum of A, B and C)\n",
    "\n",
    "    A frequency dictionary for the same data would be:\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "{'A': 0.175, 'B': 0.325, 'C': 0.5}\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "    Summing up to 1.\n",
    "\n",
    "    When passing a dictionary as an argument, you should indicate\n",
    "    whether it is a count or a frequency dictionary. Therefore the\n",
    "    FreqTable class constructor requires two arguments: the dictionary\n",
    "    itself, and FreqTable.COUNT or FreqTable.FREQ indicating counts or\n",
    "    frequencies, respectively.\n",
    "\n",
    "    Read expected counts. readCount will already generate the\n",
    "    frequencies Any one of the following may be done to geerate the\n",
    "    frequency table (ftab):\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "```python\n",
    "from Bio.SubsMat import *\n",
    "ftab = FreqTable.FreqTable(my_frequency_dictionary, FreqTable.FREQ)\n",
    "ftab = FreqTable.FreqTable(my_count_dictionary, FreqTable.COUNT)\n",
    "ftab = FreqTable.read_count(open('myCountFile'))\n",
    "ftab = FreqTable.read_frequency(open('myFrequencyFile'))\n",
    "```"
   ]
  },
  {
   "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.5.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}