PythonNotebook.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "a20290f5-2a26-4622-9e8b-ce94cc1dcb1a"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "\n",
    "# <center>Introduction to Python</center>\n",
    "<p>\n",
    "<center><i>Yale Center for Research Computing</i></center>\n",
    "<p>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "a8801e71-32d1-4497-876e-4c193189da28"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## What is the Yale Center for Research Computing?\n",
    "\n",
    "\n",
    "- Independent center under the Provost's office\n",
    "- Created to support your research computing needs\n",
    "- Focus is on high performance computing and storage\n",
    "- ~20 staff, including applications specialists and system engineers\n",
    "- Available to consult with and educate users\n",
    "- Manage compute clusters and support users\n",
    "- Located at 160 St. Ronan St, at the corner of Edwards and St. Ronan\n",
    "- [ycrc.yale.edu](https://research.computing.yale.edu/)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "494d6527-1db0-4bca-a3a9-ea7d28eed532"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Why Python?\n",
    "- Free, portable, easy to learn\n",
    "- Wildly popular, huge and growing community\n",
    "- Intuitive, natural syntax\n",
    "- Large ecosystem of libraries/packages\n",
    "- Ideal for rapid prototyping but also for large applications\n",
    "- Very efficient to write, reasonably efficient to run as is\n",
    "- Can be very efficient (numpy, cython, ...)\n",
    "- Huge number of packages (modules)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Installing Python\n",
    "\n",
    "We recommend Anaconda:\n",
    "- easy to install\n",
    "- easy to add additional packages\n",
    "- allows creation of custom environments\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "12ea1d90-4e9e-435e-9531-27bf0426ae6f"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## You can use Python to...\n",
    "- Convert or filter files\n",
    "- Automate repetitive tasks\n",
    "- Compute statistics\n",
    "- Build processing pipelines\n",
    "- Build simple web applications\n",
    "- Perform large numerical computations\n",
    "- ...\n",
    "\n",
    "You can use Python instead of bash, Java, or C\n",
    "\n",
    "Python can be run interactively or as a program"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "f8c8a2ab-9b15-45e6-bb89-f37ae176b4f2"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Different ways to run Python\n",
    "\n",
    "1. Jupyter notebook\n",
    "    ``` bash\n",
    "   jupyter notebook notebook.ipynb\n",
    "   ```  \n",
    "\n",
    "1. Run interpreter interactively\n",
    "\n",
    "   ``` bash\n",
    "   python\n",
    "   ```\n",
    "\n",
    "1. Create a file using editor, then:\n",
    "\n",
    "   ``` bash\n",
    "   python myscript.py\n",
    "   ```\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Follow along using Binder:\n",
    "- https://mybinder.org/v2/gh/ycrc/Python-Bootcamp/master?filepath=PythonNotebook.ipynb\n",
    "\n",
    "or \n",
    "\n",
    "- [`https://tinyurl.com/y9l66yly`](https://tinyurl.com/y9l66yly)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Installing Anaconda Python and getting tutorial notebook (later)\n",
    "\n",
    "1. Install Anaconda Python (includes Jupyter notebook)\n",
    "\n",
    "Go to the anaconda.com download page https://www.anaconda.com/products/individual and get the python 3.8 for your OS.  Follow the instructions to install it (this will depend on your OS).  If you don't want to create an account, you can just use one of these links:\n",
    "\n",
    "```\n",
    "https://repo.anaconda.com/archive/Anaconda3-2021.05-Windows-x86_64.exe\n",
    "https://repo.anaconda.com/archive/Anaconda3-2021.05-MacOSX-x86_64.pkg\n",
    "https://repo.anaconda.com/archive/Anaconda3-2021.05-Linux-x86_64.sh\n",
    "```\n",
    "2. Get tutorial files\n",
    "\n",
    "The files are in a github repository.  If you are familiar with git, you can do:\n",
    "     git clone https://github.com/ycrc/Python-Bootcamp.git\n",
    "\n",
    "Alternatively, you can download the files as a zip, and unzip them:\n",
    "https://github.com/ycrc/Python-Bootcamp/archive/master.zip\n",
    "\n",
    "3.  Run the Jupyter notebook and select  PythonNotebook.ipynb from the tutorial folder.\n",
    "Again, exactly how to do this will depend on your OS.  \n",
    "In windows type anaconda into the search box and open \n",
    "\n",
    "Jupyter should connect to your web browser and open the notebook."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "1ebb4259-2081-4ec5-af7b-c14dedde0872"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Python 2 versus 3\n",
    "\n",
    "- Python 2 no longer supported, but some old projects have not migrated\n",
    "- Any new code should be Python 3\n",
    "- Key changes in Python 3:\n",
    " - print is a function\n",
    " - integer division returns float\n",
    " - functions returning iterators instead of lists (range, dic.keys()...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "%%python2\n",
    "print 7/3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "e675a434-eb22-4d12-a159-d6464ec9d916"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Basic Python types: _integers, floats, strings, booleans_\n",
    "\n",
    "The building blocks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "radius=2\n",
    "pi=3.14\n",
    "diam=radius*2\n",
    "area=pi*(radius**2)\n",
    "title=\"fun with strings\"\n",
    "pi='cherry'\n",
    "radius=2.5\n",
    "delicious=True\n",
    "print(radius)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "cf3cb8c6-45f0-48de-b171-598ebf92904f"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "\n",
    "- Variables do not need to be declared or typed\n",
    "- Integers and floating points can be used together\n",
    "- The same variable can hold different types\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "a304d336-0970-4929-b6ed-13c646c59725"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### _strings_ ...\n",
    "\n",
    "- are defined with `'` or `\"`\n",
    "- cannot be modified\n",
    "- have lots of useful methods, e.g. `strip`, `split`, etc\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "4765d0a2-2a8a-497a-b969-9d6b80486aec"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "s=\"Some String's\"\n",
    "s[2]='3'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "41936fb6-e55b-4679-9285-f603d2fe9732"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "s=\"hello\"\n",
    "s.split('l')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "01c5e397-ce22-4e60-bece-26892ff66c77"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Data structures: _lists_\n",
    "\n",
    "Like arrays in other languages.  \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "32dcb56e-4dae-4bef-821a-48f9334a6f3d"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers=[1,2,3,4,5,6,7,8,9]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "8dc3a7fd-a5ce-41f2-b336-c9df5e112dfd"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers[3]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "fb27453f-b80c-48e1-96d6-d77ec61835cd"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers[5:7]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers[1:6:2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "6fc8d55d-6c8e-4b5b-ac6f-dc4da49e5dc1"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers[2]=3.14\n",
    "numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "4f14437f-6521-4091-bea7-e3c7a4506478"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers.reverse()\n",
    "numbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "a4ee966d-3b8a-4abb-919f-8430acedf470"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Lists are more flexible than arrays\n",
    "\n",
    "You can ...\n",
    "\n",
    "- insert or append new elements\n",
    "- remove elements\n",
    "- nest lists\n",
    "- combine values of different types into lists\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-04-12T11:37:57.973507Z",
     "start_time": "2020-04-12T11:37:57.959003Z"
    },
    "nbpresent": {
     "id": "b4655176-1398-4987-8ff0-6f46bd24c9e9"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers=[1,2,3,4,5,6,7,8,9]\n",
    "numbers[2]=[11,12,13]\n",
    "numbers[2][0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "378a93fb-d87e-46da-80c4-ac69fd10baaf"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers[3:6]=['four to six']\n",
    "numbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "82714e0c-2a74-45d8-b7cc-99a09924f080"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Data structures: _tuples_\n",
    "\n",
    "Like lists, but not modifiable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "e772b731-9081-4946-92c9-992f9c176d99"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "tup=(1,2,3,4,5,6,7,8,9)\n",
    "tup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "3797ae0a-3333-4e11-b5f8-89569c0d37ba"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "tup[4]=4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "f4c0903e-8146-4f12-9b0f-f09553cabfaa"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Data structures: _dictionaries_\n",
    "\n",
    "Dicts are what python calls \"hash tables\"\n",
    "\n",
    "- dicts associate keys with values, which can be of (almost) any type\n",
    "- dicts have length, but are not ordered\n",
    "- looking up values in dicts is very fast, even if the dict is BIG.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "812f1017-5af8-477d-9f6f-ded5379ca997"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "coins={'penny':1, 'nickle':5, 'dime':10, 'quarter':25}\n",
    "coins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "coins['dime']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "coins['half-dollar'] = 50\n",
    "coins"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "138e7c96-d712-484e-8e82-5b21a548b0af"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Control Flow Statements: _if_\n",
    "\n",
    "- `if` statements run a test, then do something based on the result\n",
    "- `else` is optional\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "8b9c7e6e-3770-4bec-ab1d-4219b015f840"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import random\n",
    "number=random.randint(0,100)\n",
    "if number < 50:\n",
    "    if number < 5: \n",
    "        print (\"really small\")\n",
    "    print (\"small\", number)\n",
    "    print (\"another line\")\n",
    "else: \n",
    "    print (\"big\", number)  \n",
    "print (\"after else\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "9e1755ac-fc3a-470f-ad85-81383bcc5548"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Control Flow Statements: _while_\n",
    "\n",
    "- While statements execute one or more statements repeatedly until the\n",
    "test is false"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "1209ed01-100c-4249-9790-7ceaf4b8d25f"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import random\n",
    "count=0\n",
    "while count<100:\n",
    "   count=count+random.randint(0,10)\n",
    "   print (count)\n",
    "print(\"done with loop\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "15c9bfc0-08a7-43e4-bfcc-fcd6ce5810ad"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Control Flow Statements: _for_\n",
    "\n",
    "For statements take some sort of _iterable_ object and loop once for\n",
    "every value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for i in [0,1,2,3]:\n",
    "    print(i)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "2e40b847-7de3-426b-93e2-66fa8507cb37"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for i in range(10):\n",
    "   print(i)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "3289621e-dc24-4f9a-a9ca-c4ed005fb5e7"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for letter in \"this string\":\n",
    "   print(letter)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true,
    "nbpresent": {
     "id": "d4bddfab-61ba-4fb3-91d3-69bc8065c1f5"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Using `for` loops and `dicts`\n",
    "\n",
    "Do something for each key (Use `items()` for keys and values)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "0601b993-3329-4700-a4c1-36e64b485275"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for val in coins:  \n",
    "   print (val)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for val in coins.items():  \n",
    "   print (val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "568d1355-af0f-4d7c-8a72-d8eb8e6a3480"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Control Flow Statements: altering loops\n",
    "While and For loops can skip steps (`continue`) or terminate early (`break`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "aedb70f1-1599-42f5-9ea7-25177e656040"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for i in range(10):\n",
    "   if i%2 != 0: continue # is i odd?\n",
    "   print (i)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "1ad39692-d322-4b4b-8b41-d7705e3ff75c"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for i in range(10):\n",
    "   if i>5: break\n",
    "   print (i)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "81547778-1b73-4649-90a9-c0b3dc386014"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Note on code blocks\n",
    "\n",
    "In the previous example:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "53639b58-c255-4deb-9d31-1db8a85ef2e8"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "for i in range(10):\n",
    "   if i>5: break\n",
    "   print(i)\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "1708f57b-4d71-4035-897b-8d4245fce9a5"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "How did we know that `print(i)` was part of the loop?  What defines a loop?\n",
    "\n",
    "- Many programming languages use `{ }` or Begin End to delineate blocks of code to treat as a single unit.\n",
    "\n",
    "- Python uses white space.  Code indented to the same level is one block.\n",
    "\n",
    "- By convention and for readability, indent a consistent number (many editors will do this for you)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## List comprehensions\n",
    "\n",
    "can replace simple loops, are fast and concise"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers=[1,2,3,4,5,6,7,8,9]\n",
    "new_numbers=[]\n",
    "for number in numbers:\n",
    "    new_numbers.append(number*5) \n",
    "numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "numbers=[1,2,3,4,5,6,7,8,9]\n",
    "new_numbers = [number * 5 for number in numbers]\n",
    "new_numbers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "69d16429-24ce-4d40-885c-8847d95c18b0"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Functions\n",
    "\n",
    "allow you to write code once and use it many times\n",
    "\n",
    "compartmentalize detail so code is more understandable\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "52ea5585-d6b7-4f6d-b565-23d11b1c819a"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "def area(w, h):\n",
    "   return w*h\n",
    "\n",
    "area(3, 10) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "b1618f18-c6fc-4d71-b0dd-af66de188f3e"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Summary of basic elements of Python\n",
    "\n",
    "- Basics: int, float, string, boolean\n",
    "- More complex: list, dict, tuple\n",
    "- Control constructs: if, while, for, list comprehension, def\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Printing & Formatting\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "i=5\n",
    "print(\"Simple \"+\"Stuff \"+str(i))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import math\n",
    "x=16\n",
    "\n",
    "# old-school way\n",
    "print(\"The sqrt of %i is %f\" % (x, math.sqrt(x)))\n",
    "\n",
    "# format function\n",
    "print(\"The sqrt of {} is {}\".format(x, math.sqrt(x)))\n",
    "print(\"{1:.0f} squared is {0}\". format(x, math.sqrt(x)))\n",
    "\n",
    "# f-strings (3.6+)\n",
    "print(f\"The sqrt of {x} is {math.sqrt(x):.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "ad7a2c1e-02a0-4b73-b35e-9beb2ae0c863"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Objects\n",
    "\n",
    "- All python values are objects (lists, strings, dicts, etc.)\n",
    "- Objects combine value(s) and methods (functions)\n",
    "- Advanced users can create their own classes of objects\n",
    "- All the usual OO stuff: inheritance, data hiding, etc.\n",
    "- Use `dir()` to discover an object's methods & attributes\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "89c2153d-7867-4fc5-be84-756c5b5f6266"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "\n",
    "numbers=[1,2,3]\n",
    "numbers.clear()\n",
    "numbers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "eee8e363-d77a-47c9-aec3-3b9d74813686"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "# try: index, upper, isupper, startswith, replace\n",
    "s=\"This is a String\"\n",
    "print(s.partition(\"i\"))\n",
    "print(s.split(\"s\"))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "df32239b-8fce-4f9a-8556-0712e98e6e49"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "dir(s) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "efc2efdd-8beb-49d3-8905-1f73fecb687d"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Example 1: File Reformatter\n",
    "\n",
    "#### Task: given a file of hundreds or thousands of lines\n",
    "\n",
    "```\n",
    "FCID,Lane,Sample_ID,SampleRef,index,Description,Control,Recipe,...\n",
    "160212,1,A1,human,TAAGGCGA-TAGATCGC,None,N,Eland-rna,Mei,Jon_mix10\n",
    "160212,1,A2,human,CGTACTAG-CTCTCTAT,None,N,Eland-rna,Mei,Jon_mix10\n",
    "160212,1,A3,human,AGGCAGAA-TATCCTCT,None,N,Eland-rna,Mei,Jon_mix10\n",
    "...\n",
    "```\n",
    "\n",
    "#### Remove the last 3 letters from the 5th column\n",
    "\n",
    "```\n",
    "FCID,Lane,Sample_ID,SampleRef,index,Description,Control,Recipe,...\n",
    "160212,1,A1,human,TAAGGCGA-TAGAT,None,N,Eland-rna,Mei,Jon_mix10\n",
    "160212,1,A2,human,CGTACTAG-CTCTC,None,N,Eland-rna,Mei,Jon_mix10\n",
    "160212,1,A3,human,AGGCAGAA-TATCC,None,N,Eland-rna,Mei,Jon_mix10\n",
    "...\n",
    "```\n",
    "\n",
    "`TAAGGCGA-TAGATCGC` -> `TAAGGCGA-TAGAT` and so on\n",
    "\n",
    "In this example, we'll show:\n",
    "- reading lines of a file\n",
    "- parsing and modifying the lines\n",
    "- writing them back out\n",
    "- creating a script to do the above and running it\n",
    "- passing the script the file to modify"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "18c0ff90-5bef-4435-996a-1a015a4ced49"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## In pseudocode \n",
    "```\n",
    " open the input file\n",
    " read the first header line, and print it out\n",
    " for each remaining line in the file\n",
    "   read the line\n",
    "   find the value in the 5th column\n",
    "   truncate it by removing the last three letters\n",
    "   put the line back together\n",
    "   print it out\n",
    "```\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "4b900509-1a2e-493a-b91d-6f7733b2eb5c"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 1: open the input file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "9df59d49-02d8-4108-b063-d58b0f01bbf5"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer=open('badfile.txt')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "4fea7e01-5e99-486d-b238-6275ff3659ce"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "e67dde7c-a815-4650-a5db-c68087ea73a4"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "`open()` takes a filename, and returns a file pointer\n",
    "\n",
    "Use that to read from the file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "684785a7-608d-4eae-968d-75d8daebdaed"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 2: read the first header line, and print it out"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "778bfab0-219a-4ebd-8d36-57f0e55d7641"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer=open('badfile.txt')\n",
    "print (file_pointer.readline().rstrip())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "2d8dd505-f5bf-4b59-ba48-b4a2d6cd5f57"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "- Call `readline()` on the file pointer to get a single line from the file\n",
    "(the header line)\n",
    "\n",
    "- `rstrip()` removes the return character at the end of the line\n",
    "\n",
    "- Then print it"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "3f85fedf-db72-44af-b673-e3d1b66c67ad"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 3: for each remaining line in the file, read the line"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "19f8b102-3960-4f9a-af07-efda4fd07d88"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer=open('badfile.txt')\n",
    "print (file_pointer.readline().rstrip())\n",
    "for line in file_pointer:\n",
    "  line=line.rstrip()\n",
    "  print(line)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "cb0419a7-7c7d-4e8a-88c2-72bbbdcddc2f"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "A file pointer is an example of an iterator.\n",
    "\n",
    "Instead of explicitly calling `readline()` for each line, we can just loop on the file\n",
    "pointer, getting one line each time.\n",
    "\n",
    "Since we already read the header, we won't get that line."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "18f07bd0-4807-400f-b98b-6da00a29a45d"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 4: find the value in the 5th column, and remove last 3 letters\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "d0c99a58-e78f-40c6-a2f9-356ba462abb3"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer=open('badfile.txt')\n",
    "print (file_pointer.readline().strip())  \n",
    "for line in file_pointer:\n",
    "    fields=line.strip().split(',')\n",
    "    fields[4]=fields[4][:-3]\n",
    "    print(fields)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "798ae994-7b95-4d34-95db-fde9158ff3ba"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Like before, we strip the return from the line.\n",
    "\n",
    "We split it into\n",
    "individual elements where we find commas.\n",
    "\n",
    "The 5th field is referenced by\n",
    "flds[4], since python starts indexing with 0.  [:-3] takes all characters\n",
    "of the string until the last 3."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Brief detour: string splitting and joining"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "s='this,is,a,string'\n",
    "l=s.split(',')\n",
    "print(l)\n",
    "print(':'.join(l))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "36640e7d-fdbc-4c4c-abc6-fbbd71a87d45"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 5: put the line back together, and print it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "1b380c2b-de9d-4045-a7c8-efb7c1d047b7"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "file_pointer=open(\"badfile.txt\")\n",
    "print (file_pointer.readline().strip())\n",
    "for line in file_pointer:\n",
    "    fields=line.strip().split(',')\n",
    "    fields[4]=fields[4][:-3]\n",
    "    print (','.join(fields))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "065395a1-8d57-470c-bd53-0bfa04850206"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    " \n",
    "Join takes a list of strings, and combines them into one string using the\n",
    "string provided. Then we just print that string."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "5438e79a-1925-4bb6-ac97-069bb89d98f1"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    " \n",
    "We would invoke it like this:\n",
    "```\n",
    "$ python Ex1.py badfile.txt\n",
    "\n",
    "$ python Ex1.py badfile.txt > fixedfile.txt\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "cfab7eb8-e915-4da8-881f-18c7030d53d9"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Example 2: directory walk with file ops\n",
    "\n",
    "Imagine you have a directory tree with many subdirectories.\n",
    "\n",
    "In those directories are files named *.fastq.  You want to:\n",
    "\n",
    "- find them\n",
    "- compress them to fastq.gz using a program\n",
    "- delete them if the conversion was successful\n",
    "\n",
    "In this example, we'll demonstrate:\n",
    "\n",
    "- traversing an entire directory tree\n",
    "- executing a program on files in that tree\n",
    "- testing for successful program execution\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "814659d4-6d87-41ab-a8da-4a89a73eb682"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    " ## In psuedocode\n",
    "```    \n",
    "for each directory\n",
    "   get a list of files in that directory\n",
    "   for each file in that directory\n",
    "     if that file's name ends with .fastq\n",
    "       create a new file name with .gz added\n",
    "       create a command to do the compression\n",
    "       run that command and check for success\n",
    "       if success\n",
    "         delete the original\n",
    "       else\n",
    "         stop\n",
    "```\n",
    "The conversion command is: \n",
    "```gzip -c file.fastq > file.fastq.gz```\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "0be403ec-10ea-4f29-994a-7ca8dd6c4ea8"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 1: directory traversal\n",
    "\n",
    "We need a way to traverse all the files and directories.\n",
    "```os.walk(dir)``` starts at dir and visits every subdirectory below it.\n",
    "It returns a list of files and subdirectories at each subdirectory.\n",
    "\n",
    "For example, imagine we have the following dirs and files:\n",
    "\n",
    "```\n",
    "Ex2dir\n",
    "Ex2dir/d1\n",
    "Ex2dir/d1/d2\n",
    "Ex2dir/d1/d2/f2.fastq\n",
    "Ex2dir/d1/f1.fastq\n",
    "```\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "506ef877-d76b-445c-af55-870d87f0d6e5"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "for path , dirs, files in os.walk('Ex2dir'):\n",
    "   print (path, dirs, files)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "34aeed1f-91e0-425e-94bd-a4d031956de2"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 2: Invoking other programs from python\n",
    "\n",
    "The [subprocess module](https://docs.python.org/3/library/subprocess.html) has a variety of ways to do this. A simple one:\n",
    "\n",
    "```\n",
    "import subprocess\n",
    "\n",
    "ret=subprocess.call(cmd, shell=True)\n",
    "\n",
    "```\n",
    "\n",
    "ret is 0 on success, non-zero error code on failure.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "1a02756f-a031-4887-a1ae-e07fc72ea137"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import subprocess\n",
    "ret_code=subprocess.call('gzip -c myfile.fastq > myfile.fastq.gz', shell=True)\n",
    "ret_code"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "d65d542b-e74e-477a-bd04-21a2e9d39e8c"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Put it all together"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "3bdbdc71-a376-4eb1-b1ab-0d28a818f904"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "import os, sys, subprocess\n",
    "sys.argv=['Ex2.py', 'Ex2dir'] # for Jupyter we'll cheat\n",
    "start=sys.argv[1]\n",
    "for path, subdirs, files in os.walk(start):\n",
    "    for file in files:\n",
    "        if file.endswith('.fastq'):\n",
    "            file_name=f'{path}/{file}'\n",
    "            cmpress_file=file_name.replace('.fastq', '.fastq.gz')\n",
    "            cmd=f'gzip -c {file_name} > {cmpress_file}'\n",
    "            print (f\"running {cmd}\")\n",
    "            ret_code=subprocess.call(cmd, shell=True)\n",
    "            if ret_code==0:\n",
    "                if os.path.exists(cmpress_file):\n",
    "                    os.remove(file_name)\n",
    "            else:\n",
    "                print (\"Failed on \", file_name)\n",
    "                sys.exit(1)\n",
    "print(\"Done\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "4b58e3b7-112a-4f65-b942-e2161e143d19"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "\n",
    "We would invoke it like this:\n",
    "```\n",
    "$ python Ex2.py Ex2dir\n",
    "```\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "8f579d6b-c04b-4b17-87c7-4c4b58ead239"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Example 3: Nested Dictionaries\n",
    "\n",
    "Dictionaries associate names with data, and allow quick retrieval by name.\n",
    "\n",
    "By nesting dictionaries, powerful lookups are fast and easy.\n",
    "\n",
    "In this example, we'll:\n",
    "- create a dict containing objects\n",
    "- load the objects with search data\n",
    "- use the dict to retrieve the appropriate object for a search\n",
    "- perform the search\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "bbe27243-669f-4da2-84c0-ac0226bcf0bc"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "\n",
    "genes.txt describes the locations of genes:\n",
    "\n",
    "(name, chrom, strand, start, end)\n",
    "\n",
    "```\n",
    "uc001aaa.3      chr1    +       11873   14409 \n",
    "uc010nxr.1      chr1    +       11873   14409 \n",
    "uc010nxq.1      chr1    +       11873   14409  \n",
    "uc009vis.3      chr1    -       14361   16765  \n",
    "uc009vit.3      chr1    -       14361   19759  \n",
    "...\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true,
    "nbpresent": {
     "id": "59a6a2fa-0c16-41bd-82b9-7d79e0371db1"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "\n",
    "mappedreads.txt describes mapped dna sequences\n",
    "\n",
    "(name, chrom, position, sequence)\n",
    "\n",
    "```\n",
    "seq1 chr1  674540   ATCTGTGCAGAGGAGAACGCAGCTCCGCCCTCGCGGT\n",
    "seq2 chr19 575000   AGAGGAGAACGCAGCTCCGCCCTCGCGGTGCTCTCCG\n",
    "seq3 chr5  441682   TCTGCATCTGCTCTGGTGTCTTCTGCCATATCACTGC\n",
    "...\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "b4d9919d-acf4-4662-bb5b-bd6c170f255f"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "We'd like to be able to quickly determine the genes overlapped by a dna sequence.\n",
    "\n",
    "First, we need a simple but efficient way to determine if two intervals overlap.\n",
    "\n",
    "intervaltree is a python module that makes that easy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "d9bf9eab-fbbf-4385-9ccf-a9b854a44638"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "from intervaltree import IntervalTree\n",
    "it=IntervalTree()\n",
    "it[4:7]='gene1'\n",
    "it[5:10]='gene2'\n",
    "it[1:11]='gene3'\n",
    "it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "1b4c4b15-ab2d-49a0-9d2d-e684a340d0b2"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "it[3]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "15ffcb66-d804-42f4-8864-59a8cb3154f7"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## General plan\n",
    "\n",
    "- create an interval tree for each chromosome\n",
    "- organize the trees in a dictionary by chromosome\n",
    "- store an interval for each gene the tree for its chromosome\n",
    "\n",
    "```\n",
    "{'chr1': IntervalTree([Interval(1000, 1100, 'GeneA'), \n",
    "                       Interval(2000, 2100, 'GeneB'), ...\n",
    " 'chr2': IntervalTree([Interval(4000, 5100, 'GeneC'), \n",
    "                       Interval(7000, 8100, 'GeneD'), ...\n",
    " 'chr3':\n",
    " ...\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "250a9cd0-3744-4d27-a902-d211a387afb2"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    " # In psuedocode\n",
    " ### Step 1: Set up the lookup table\n",
    " ```\n",
    " create empty dict\n",
    " open the gene file\n",
    " for each line in the file\n",
    "    get gene name, chrom, start, end\n",
    "    initialize an intervaltree for the chrom, if needed, and add to dict\n",
    "    add the interval and gene name to the interval tree\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "efbf1716-a4cc-4be3-af5a-a90170c6ee06"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "import sys\n",
    "from intervaltree import IntervalTree\n",
    "\n",
    "print(\"initializing table\")\n",
    "table={}\n",
    "sys.argv=['Ex3.py', 'genes.txt', 'mappedreads.txt', 'results.txt'] # for Jupyter\n",
    "for line in open(sys.argv[1]):\n",
    "    genename, chrm, strand, start, end = line.split()\n",
    "    if not chrm in table:\n",
    "        table[chrm]=IntervalTree()\n",
    "    table[chrm][int(start):int(end)]=genename\n",
    "print(\"done\")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "fc1e602a-2e38-4de8-b24c-6a879129b14e"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "table['chr1'][770000:780000]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "2e97e704-e457-46e0-9efa-387a21066fff"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Step 2: Use the interval trees to find overlapped genes\n",
    "```\n",
    " open the dna sequence file\n",
    " for each line in the file:\n",
    "   get chrom, mapped position, and dna sequence\n",
    "   look up the interval tree for that chrom in the dict\n",
    "   search the interval tree for overlaps [pos, pos+len]\n",
    "   print out the gene names\n",
    " ```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpresent": {
     "id": "bdd1c433-9194-4c66-9fe5-5fcba100cf2a"
    },
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "print(\"reading sequences\")\n",
    "\n",
    "outfp=open(sys.argv[3], 'w')\n",
    "for line in open(sys.argv[2]):\n",
    "    name, chrm, pos, seq = line.strip().split()\n",
    "    genes=table[chrm][int(pos):int(pos)+len(seq)]\n",
    "    if genes:\n",
    "        print(\"\\t\".join([name, chrm, pos, seq]))\n",
    "        for gene in genes:\n",
    "            print (f'\\t{gene.data}')\n",
    "print(\"done\")\n",
    "outfp.close()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Example 4: Travelling salesman\n",
    "- Given N cities, find shortest route visiting each once\n",
    "- NP-complete -> very expensive to find optimal solution\n",
    "- We'll use a heuristic approach, simulated annealing:\n",
    "\n",
    "```\n",
    "T = initialT\n",
    "curr = best = initialsolution\n",
    "while T > stopT:\n",
    "    candidate = pick two cities, swap and reverse path between them\n",
    "    if cost(candidate) < cost(best):\n",
    "        curr=candidate\n",
    "        best=candidate\n",
    "    elif cost(candidate)-cost(best) < T\n",
    "        curr=candidate\n",
    "    T = T * alpha\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "import random\n",
    "def initial(n):\n",
    "    return [[random.randint(0,1000), random.randint(0,1000)] for i in range(n)]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "initial(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "def plot(path):\n",
    "    xs = []; ys = []\n",
    "    for x,y in path:\n",
    "        xs.append(x)\n",
    "        ys.append(y)\n",
    "    plt.plot(xs, ys, 'co-')\n",
    "    plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "cities=initial(100)\n",
    "plot(cities)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "import math\n",
    "\n",
    "def dist(a,b):\n",
    "    return math.sqrt((a[0]-b[0])**2+(a[1]-b[1])**2)\n",
    "\n",
    "def cost(curr):\n",
    "    total=0\n",
    "    n=len(curr)\n",
    "    for i in range(n):\n",
    "        total+=dist(curr[i], curr[(i+1)%n])\n",
    "    return total\n",
    "\n",
    "cost(cities)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def doswap(curr):\n",
    "    cand=curr[:]\n",
    "    n=len(cand)\n",
    "    i,j = sorted(random.sample(range(n),2))\n",
    "    cand[i:j+1] = reversed(cand[i:j+1])\n",
    "    return cand\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "curr=initial(5)\n",
    "cand=doswap(curr)\n",
    "print(curr)\n",
    "print(cand)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def p_accept(candcost, currcost, T):\n",
    "    p = math.exp(-abs(candcost - currcost) / T)\n",
    "    return p\n",
    "\n",
    "p_accept(20010,20000,5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def main(n):\n",
    "    #random.seed(0)\n",
    "    curr=initial(n)\n",
    "    plot(curr)\n",
    "    currcost=cost(curr)\n",
    "    T=math.sqrt(n)\n",
    "    while T > 1e-10:\n",
    "        cand=doswap(curr)\n",
    "        candcost=cost(cand)\n",
    "        if candcost<currcost:\n",
    "            curr=cand\n",
    "            currcost=candcost\n",
    "            print(\"forward T %f cost %f\" % (T, currcost))\n",
    "        elif p_accept(candcost, currcost, T) > random.random():\n",
    "            curr=cand\n",
    "            currcost=candcost\n",
    "            print(\"backward T %f cost %f\" % (T, currcost))\n",
    "        T=T*.999\n",
    "    plot(curr)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "main(100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "788cac40-7b07-4690-9aa9-5c7fe7c296c1"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Important Packages \n",
    "- numpy: high performance arrays\n",
    "- scipy: stats, linear alg, etc.\n",
    "- multiprocessing: easy parallelism\n",
    "- matplotlib: plotting\n",
    "- pandas: stats with R-like dataframes\n",
    "- flask, web2py: web applications"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpresent": {
     "id": "7598763d-f96e-475c-8824-10d1c42a120b"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Python Resources we like\n",
    "\n",
    "- anaconda python: www.continuum.io\n",
    "- Jupyter notebook\n",
    "- pycharm debugger: www.jetbrains.com\n",
    "- _Introducing Python_, Bill Lubanovic, O'Reilly\n",
    "- _Python in a Nutshell_, Alex Martelli, O'Reilly\n",
    "- _Python Cookbook_, Alex Martelli, O'Reilly\n",
    "- Google's python class: https://www.youtube.com/watch?v=tKTZoB2Vjukxo\n",
    "- https://docs.python.org/3/tutorial\n",
    "- codecademy https://www.codecademy.com/learn/learn-python\n"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "celltoolbar": "Slideshow",
  "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.8.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}