Alignment_Fast3.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def pairwise_NW(hl_list1, hl_list2, score_matched, score_mismatched, score_gap, isret):\n",
    "\n",
    "    isH1Long = True if len(hl_list1) >= len(hl_list2) else False \n",
    "    long_list = hl_list1 if isH1Long else hl_list2\n",
    "    short_list = hl_list2 if isH1Long else hl_list1\n",
    "    \n",
    "    m = len(long_list)\n",
    "    n = len(short_list)\n",
    "    \n",
    "    # create a  score matrix and initial first vertical row\n",
    "    # scoreMatrix element = (score, directionsList())\n",
    "    # directionsList contains direction\n",
    "    # (because there may be multiple directions in one score cell)\n",
    "    #\n",
    "    # direction: 0-matched, 1-mismatched, 2-insert verticalGap, 3-insert horizontalGap\n",
    "    scoreMatrix = [[(score_gap*i,list()) if j==0 else (0,list()) for j in range(m+1)] for i in range(n+1)]\n",
    "    scoreMatrix[0] = [(0+score_gap*j,list()) for j in range(m+1)] # initial first horizontal column values\n",
    "    \n",
    "    # run Needleman-Wunsch algorithm\n",
    "    for i in range(1, n+1): # vertical(short)\n",
    "#         s = \"\"\n",
    "        for j in range(1, m+1): # horizontal(long)\n",
    "            # count gap scores\n",
    "            verticalGapScore = scoreMatrix[i][j-1][0] + score_gap\n",
    "            horizontalGapScore = scoreMatrix[i-1][j][0] + score_gap\n",
    "            \n",
    "            matchScore = scoreMatrix[i-1][j-1][0] # initial\n",
    "            isMatchOrNot = 0\n",
    "            \n",
    "            # count match/mismatch score\n",
    "            if compareSequence( long_list[j-1], short_list[i-1] ):\n",
    "                matchScore = matchScore + score_matched\n",
    "                isMatchOrNot = 1\n",
    "            else:\n",
    "                matchScore = matchScore + score_mismatched\n",
    "                isMatchOrNot = -1\n",
    "            \n",
    "            maxScore = max(verticalGapScore, horizontalGapScore, matchScore) # find maximum score\n",
    "            \n",
    "            directions = list()\n",
    "            # add (direction to each cell of scoreMatrix) into directions\n",
    "            if(maxScore == matchScore):\n",
    "                if isMatchOrNot==1: # matched\n",
    "                    directions.append(0)\n",
    "                else:               # mismatched\n",
    "                    directions.append(1)\n",
    "            if(maxScore == verticalGapScore): # vertical add gap\n",
    "                directions.append(2)\n",
    "            if(maxScore == horizontalGapScore): # horizontal add gap\n",
    "                directions.append(3)\n",
    "            \n",
    "            scoreMatrix[i][j] = (maxScore,directions)\n",
    "\n",
    "    output_longList = []\n",
    "    output_shortList = []\n",
    "\n",
    "    # walk sequence path(from end cell to first cell)\n",
    "\n",
    "    \n",
    "    while(m>0 and n>0) : # if a direction touch margin than finish walking\n",
    "        \n",
    "        ## only pick up one path\n",
    "        ## if want to get all paths, go through all directions in 'scoreMatrix[n][m][1]'\n",
    "        direction = scoreMatrix[n][m][1][0]\n",
    "\n",
    "        # the condition of direction check\n",
    "        if direction == 0:\n",
    "            output_longList.append((m, long_list[m-1]))\n",
    "            output_shortList.append((n, short_list[n-1]))\n",
    "            m=m-1\n",
    "            n=n-1\n",
    "\n",
    "        elif direction == 1:\n",
    "            output_longList.append((m, long_list[m-1]))\n",
    "            output_shortList.append((n, short_list[n-1]))\n",
    "            m=m-1\n",
    "            n=n-1\n",
    "\n",
    "        elif direction == 2:\n",
    "            output_shortList.append((-1, \"=\"))\n",
    "            output_longList.append((m, long_list[m-1]))\n",
    "            m=m-1\n",
    "\n",
    "        else:\n",
    "            output_longList.append((-1, \"=\"))\n",
    "            output_shortList.append((n, short_list[n-1]))\n",
    "            n=n-1\n",
    "    \n",
    "    # two below loop are for walking path which didn't stop at (0,0)\n",
    "    # if not stop at (0,0), should insert gap to push back until (0,0)\n",
    "    while(n>0):\n",
    "        output_longList.append((-1, \"=\"))\n",
    "        output_shortList.append((n, short_list[n-1]))\n",
    "        n=n-1\n",
    "        \n",
    "    while(m>0):\n",
    "        output_shortList.append((-1, \"=\"))\n",
    "        output_longList.append((m, long_list[m-1]))\n",
    "        m=m-1\n",
    "                \n",
    "    alignment_result = list(zip(output_longList, output_shortList)) if isH1Long else list(zip(output_shortList, output_longList))\n",
    "    alignment_result = alignment_result[::-1]\n",
    "    \n",
    "    return None, None, alignment_result # MIKE: 20170811, for the compatibility"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def compareSequence(seq1, seq2):\n",
    "    return True if seq1 == seq2 else False"
   ]
  }
 ],
 "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.4.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}