{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os, io, random\n",
    "import string\n",
    "import numpy as np\n",
    "\n",
    "from Bio.Seq import Seq\n",
    "from Bio.Align import MultipleSeqAlignment\n",
    "from Bio import AlignIO, SeqIO\n",
    "\n",
    "#import panel as pn\n",
    "#import panel.widgets as pnw\n",
    "#pn.extension()\n",
    "\n",
    "from bokeh.plotting import figure, show\n",
    "from bokeh.models import ColumnDataSource, Plot, Grid, Range1d\n",
    "from bokeh.models.glyphs import Text, Rect\n",
    "from bokeh.layouts import gridplot\n",
    "\n",
    "from bokeh.io import export_svgs\n",
    "import selenium\n",
    "\n",
    "# THIS LAST BIT IS TO NAVIAGATE TO WHERE MY CUSTOM MODULES ARE LOCATED\n",
    "if os.getcwd()[-4:] != 'AIMS':\n",
    "    default_path = os.getcwd()[:-10]\n",
    "    os.chdir(default_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# All of the Code in this script is adapted from this source:\n",
    "https://dmnfarrell.github.io/bioinformatics/bokeh-sequence-aligner"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def view_alignment(aln, fontsize=\"9pt\", plot_width=800):\n",
    "    \"\"\"Bokeh sequence alignment view\"\"\"\n",
    "\n",
    "    #make sequence and id lists from the aln object\n",
    "    seqs = [rec.seq for rec in (aln)]\n",
    "    ids = [rec.id for rec in aln]    \n",
    "    text = [i for s in list(seqs) for i in s]\n",
    "    colors = get_colors(seqs)    \n",
    "    N = len(seqs[0])\n",
    "    S = len(seqs)    \n",
    "    width = .4\n",
    "\n",
    "    x = np.arange(1,N+1)\n",
    "    y = np.arange(0,S,1)\n",
    "    #creates a 2D grid of coords from the 1D arrays\n",
    "    xx, yy = np.meshgrid(x, y)\n",
    "    #flattens the arrays\n",
    "    gx = xx.ravel()\n",
    "    gy = yy.flatten()\n",
    "    #use recty for rect coords with an offset\n",
    "    recty = gy+.5\n",
    "    h= 1/S\n",
    "    #now we can create the ColumnDataSource with all the arrays\n",
    "    source = ColumnDataSource(dict(x=gx, y=gy, recty=recty, text=text, colors=colors))\n",
    "    plot_height = len(seqs)*15+50\n",
    "    x_range = Range1d(0,N+1, bounds='auto')\n",
    "    if N>100:\n",
    "        viewlen=100\n",
    "    else:\n",
    "        viewlen=N\n",
    "    #view_range is for the close up view\n",
    "    view_range = (0,viewlen)\n",
    "    tools=\"xpan, xwheel_zoom, reset, save\"\n",
    "\n",
    "    #entire sequence view (with zoom)\n",
    "    p = figure(title=None, plot_width= plot_width, plot_height=150,\n",
    "               x_range=x_range, y_range=(0,S), tools=tools,\n",
    "               min_border=0, toolbar_location='below')\n",
    "    glyph = Text(x=\"x\", y=\"y\", text=\"text\", text_align='center',text_color=\"black\",\n",
    "                text_font=\"monospace\",text_font_size=fontsize)\n",
    "    rects = Rect(x=\"x\", y=\"recty\",  width=1, height=1, fill_color=\"colors\",\n",
    "                 line_color=None, fill_alpha=0.6)\n",
    "    # This add_glyph line is where we add in letters, and this wasn't in the original code.\n",
    "    p.add_glyph(source, glyph)\n",
    "    p.add_glyph(source, rects)\n",
    "    p.yaxis.visible = False\n",
    "    p.grid.visible = False  \n",
    "    p.output_backend='svg'\n",
    "    \n",
    "    ################################################################\n",
    "    ## THIS SECTION IS FOR PLOTTING AN INTERACTIVE FIGURE, SO FOR ##\n",
    "    ##                   NOW LEAVE OUT OVERALL                    ##\n",
    "    ################################################################\n",
    "    #sequence text view with ability to scroll along x axis\n",
    "    p1 = figure(title=None, plot_width=plot_width, plot_height=plot_height,\n",
    "                x_range=view_range, y_range=ids, tools=\"xpan,reset\",\n",
    "                min_border=0, toolbar_location='below')#, lod_factor=1)          \n",
    "    glyph = Text(x=\"x\", y=\"y\", text=\"text\", text_align='center',text_color=\"black\",\n",
    "                text_font=\"monospace\",text_font_size=fontsize)\n",
    "    rects = Rect(x=\"x\", y=\"recty\",  width=1, height=1, fill_color=\"colors\",\n",
    "                line_color=None, fill_alpha=0.4)\n",
    "    p1.add_glyph(source, glyph)\n",
    "    p1.add_glyph(source, rects)\n",
    "    \n",
    "    p1.grid.visible = False\n",
    "    p1.xaxis.major_label_text_font_style = \"bold\"\n",
    "    p1.yaxis.minor_tick_line_width = 0\n",
    "    p1.yaxis.major_tick_line_width = 0\n",
    "\n",
    "    p = gridplot([[p],[p1]], toolbar_location='below')\n",
    "    return p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def make_seq(length=40):    \n",
    "    return ''.join([random.choice(['A','C','T','G']) for i in range(length)])\n",
    "\n",
    "def mutate_seq(seq):\n",
    "    \"\"\"mutate a sequence randomly\"\"\"\n",
    "    seq = list(seq)\n",
    "    pos = np.random.randint(1,len(seq),6)    \n",
    "    for i in pos:\n",
    "        seq[i] = random.choice(['A','C','T','G'])\n",
    "    return ''.join(seq)\n",
    "\n",
    "def get_colors(seqs):\n",
    "    \"\"\"make colors for bases in sequence\"\"\"\n",
    "    text = [i for s in list(seqs) for i in s]\n",
    "    clrs = {'A':'white','G':'white','L':'white','M':'white','F':'white','W':'white',\n",
    "            'K':'blue','Q':'gray','E':'red','S':'grey',\n",
    "            'P':'white','V':'white','I':'white','C':'yellow','Y':'grey','H':'blue',\n",
    "            'R':'blue','N':'grey','D':'red','T':'grey','.':'white'}\n",
    "    #clrs =  {'A':'red','T':'green','G':'orange','C':'blue','-':'white'}\n",
    "    colors = [clrs[i] for i in text]\n",
    "    return colors\n",
    "\n",
    "def muscle_alignment(seqs):\n",
    "    \"\"\"Align 2 sequences with muscle\"\"\"\n",
    "    filename = 'temp.faa'\n",
    "    SeqIO.write(seqs, filename, \"fasta\")\n",
    "    name = os.path.splitext(filename)[0]\n",
    "    from Bio.Align.Applications import MuscleCommandline\n",
    "    cline = MuscleCommandline(input=filename, out=name+'.txt')\n",
    "    stdout, stderr = cline()\n",
    "    align = AlignIO.read(name+'.txt', 'fasta')\n",
    "    return align"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# These are notes for determining what is the \"Prevalent Standard Gene\" in Mouse Data\n",
    "So we need to make these same alignments for mouse antibodies. All of these are *01*\n",
    "\n",
    "Top poly Heavy: HV3-6, HV5-6, HV3-8, HV5-17, HV7-3, HV6-3\n",
    "Top mono Heavy: HV5-17, HV5-4, HV3-6, HV5-6, HV6-3, HV3-5\n",
    "\n",
    "MASSIVE Connection between HV3-6 and KV8-30 in poly that doesn't exist in mono.\n",
    "HV3-6 def enriched in poly. \n",
    "HV5-6 probably close to equal occurence. HV6-3 probably equal\n",
    "HV5-17 probably enriched in mono.\n",
    "\n",
    "Top poly Light: KV8-30, KV10-96, KV5-48, KV19-93, KV17-127\n",
    "Top mono Light: KV10-96, KV17-127, KV8-30, KV4-59, KV14-111\n",
    "\n",
    "So overall light chain values are VERY different"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# NOTE: Saving as an SVG function for some reason doesn't work well\n",
    "Instead, just \"print\" the resultant web page that pops up to a pdf and use Adobe illustrator to take out the bottom part"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Note, in this alignments folder are all of the fasta files used for alignments in the paper\n",
    "aln = AlignIO.read('app_data/alignments/all_genes_human.aln','fasta')\n",
    "p = view_alignment(aln, plot_width=900)\n",
    "show(p)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python3 (CDRpro)",
   "language": "python",
   "name": "cdrpro"
  },
  "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}