LAM.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# <u>Technical support : LAM-HTGTS</u>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1) Introduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;1.1) Background"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://www.ncbi.nlm.nih.gov/pubmed/27031497"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;1.2) LAM-HTGTS"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://www.ncbi.nlm.nih.gov/pubmed/27031497"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;1.3) Documentation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This documentation contains the requierements and installation process to set up LAM-HTGTS process, the preprocessing of your data, following by the LAM-HTGTS program itself. Last but not least, a solid postprocessing of this data in order to dig throught results from LAM-HTGTS technique.<br>\n",
    "\n",
    "This is dedicated to beginner linux users comfortable with computer science concepts.<br>\n",
    "\n",
    "All scripts and documentations are accessible on github repositories :<br>\n",
    "\n",
    "- LAM-HTGTS technique : https://github.com/robinmeyers/transloc_pipeline\n",
    "- My postprocessing of LAM-HTGTS technique : https://github.com/Baboon61/LAM\n",
    "\n",
    "While reading, you will find 2 kinds of block instructions : <br>\n",
    "\n",
    "- (### CODE ###) you need to write the line(s) below in a terminal. Modify the command line at your convenance. To open a terminal, right click on desktop or search it in the applications tool bar. To copy/paste in a ternimal Ctrl+Shift+C/Ctrl+Shift+V (You will be able to copy/paste after the CD Guest addition (2.2.2))\n",
    "- (### INFORMATION ###) meaning there are information about how a program works, like all the options a scripts can handle. It is possible to use an option in the short way, using one letter (-m) or in the long way, using the complete word for the option (--metadata). Note the unique hypehn and the double hypen. Options can be followed by arguments which can be integer, float or string. Non essential options are in brakets. This block should not be coopy paste in a terminal, this is not a command line.<br>\n",
    "\n",
    "All command lines invinting the user to input parameters have output/log file where you can find the information about the running program, like parameters set up, libraries checked, errors and warnings.\n",
    "\n",
    "To create a new file, open a terminal in the location you want to save your file then :<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "gedit name_of_file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "If you got error about permissions on script programs run the following :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/LAM/Scripts\n",
    "chmod 755 *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;1.4) Good practises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "All scripts contains simple comments about what is going on major parts. If you want more information about Python or R scripts, open them in text editor and read the commented lines"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Please, try to avoid special characters (space, hyphen, star, hash...) in all your inputs (parameter arguments, parameter files...)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The minimum hardware config is 30Go of hard drive space for mm10 with one reference genome and 6Go of RAM"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2) Installations & Requierements"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In installation steps, make sure you have no error at the end of each command lines. If you do, please contact bastienherve61340@orange.fr"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.1) Install a Virtual Machine tool"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You can install VirtualBox from this link : https://www.virtualbox.org/wiki/Downloads"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You can use any virtual machine tool you want."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.1.1) Download Linux iso"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "An iso is a file that looks like a CD for the machine. You can download an iso to mount a virtual machine in a virtual box.<br >\n",
    "This pipeline works with Ubuntu 18.04, could work with an other distribution or lower Ubuntu version (not tested)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Get the iso here : https://www.ubuntu.com/download/desktop (take amd64)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.1.2) Create your Linux Virtual Machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Run your VirtualBox application.<br >\n",
    "- Go to new, I named mine \"LAM\", choose type -> Linux, version Ubuntu 64 bits, memory space 8000Mio (or less if you don't have), hard drive -> create a virtual hard drive now<br >\n",
    "- File size -> 30Go (If you have 1 genome), VDI, allocate dynamicaly<br >\n",
    "- And create!\n",
    "- Restart your machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.2) Install Linux"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Run your linux virtual machine\n",
    "- First run, it will ask you the iso file (click on the folder and look for the iso file)\n",
    "- Install it with all the updates (you can erase disk and install ubuntu, it is a virtual machine anyway)\n",
    "- Name your machine whatever you want but remember it, it is your user name (mine is : lam)\n",
    "- Wait...\n",
    "- Then, open a terminal, write the line below"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "#To update the machine\n",
    "sudo apt-get update\n",
    "sudo apt-get -y upgrade"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Shutdown the linux virtual machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.2.1) Install shared folder"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "No need to copy/paste all your data on the virtual machine, just create a shared folder"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Create a LAM directory wherever you want, in your Documents for example"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- In machine manager -> share folder\n",
    "- Add new share directory, pathway -> C:\\Users\\your_name\\Documents\\LAM\n",
    "- Select auto mount !\n",
    "- In machine manager, go to Settings -> General -> Advanced\n",
    "- Shared paper press to bidirectionnal. Same for copy/paste"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Run linux virtual machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.2.2) Install CD Guest addition"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- To install the guest addiction. Go in the Virtual Machine menu, Devices -> Insert CD Guest Addition<br>\n",
    "- If you only got a kernel error you can go to the next step. If you got a \"Building the main Guest Additions module... fail\", that means the install failed.\n",
    "If you have trouble to install the CD Guest Addition, you can try to get a more recent version of the CD Guest Addition ( https://download.virtualbox.org/virtualbox/5.2.18/VBoxGuestAdditions_5.2.18.iso ), or in general ( https://download.virtualbox.org/virtualbox/ )<br>\n",
    "Then go the Virtual Machine menu, Devices -> optical drives -> Choose disk image, and look for the file you just downloaded\n",
    "- You can eject VBox disk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Change lam by your machine name\n",
    "sudo usermod -aG vboxsf lam"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Restart linux virtual machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "From here, you should be able to access the shared directory as an hard drive (with \"sf_\" in front of the share folder name).<br>\n",
    "It will also be accessible at path : /media/sf_LAM<br>\n",
    "The path will be very usefull in all command line."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.3) Install GIT"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install git"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4) Install Python (2.7)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Minimal requierement : python2.7"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "python -V\n",
    "###If you don't have it (Python 2.7)\n",
    "sudo apt-get install python"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.1) Install Pip"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install python-pip"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.2) Install BioPython"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install biopython"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.3) Install Pandas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install pandas"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.4) Install Pyensembl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install pyensembl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.5) Install SciPy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install scipy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.6) Install IPython"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install ipython==5.0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.4.7) Install ipywidgets"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This import is not needed if you don't want to use CSReport"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo pip install ipywidgets"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.3.3) Install python-tk"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This import is not needed if you don't want to use CSReport"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install python-tk"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.5) Install Blast"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install ncbi-blast+"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.6) Install Bedtools"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install bedtools"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.7) Import data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Fastq you got from experiments need to be see by the virtual machine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create a Data folder in your LAM shared folder"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###On your computer copy/paste your fastq.gz in Data (inside your share folder) (C:\\Users\\your_name\\Documents\\LAM\\Data\\name_of_your_data)\n",
    "###Create a new folder per study (Example : Puce1, like :\\Users\\your_name\\Documents\\LAM\\Data\\Puce1)\n",
    "###You can access the directory created at this path /media/sf_LAM/Data/Puce1/\n",
    "ls /media/sf_LAM/Data/Puce1/\n",
    "###You should be able to see your file here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.8) Import Scripts"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Postprocessing scripts are store on github"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remember the first \"lam\" in the next command line is you user name, change it if needed\n",
    "mkdir /home/lam/Documents/LAM\n",
    "cd /home/lam/Documents/LAM\n",
    "mkdir Scripts\n",
    "cd Scripts\n",
    "git clone https://github.com/Baboon61/LAM.git ."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9) Install Transloc_pipeline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.1) Install application"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/\n",
    "mkdir Software\n",
    "cd Software/\n",
    "mkdir robinmeyers\n",
    "cd robinmeyers\n",
    "git clone https://github.com/robinmeyers/transloc_pipeline.git"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "**Apply modifications till github patch by the author**<br>\n",
    "**Not needed if Perl version between 5.14 and 5.22**<br>\n",
    "With the Perl new version (>=5.23) a feature is depreciated, some modification are needed in 2 files"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sed -i \"s/values \\$R1_alns_ref/values %{\\$R1_alns_ref}/g\" /home/lam/Documents/Software/robinmeyers/transloc_pipeline/lib/TranslocHelper.pl\n",
    "sed -i \"s/values \\$R2_alns_ref/values %{\\$R2_alns_ref}/g\" /home/lam/Documents/Software/robinmeyers/transloc_pipeline/lib/TranslocHelper.pl\n",
    "sed -i \"s/keys \\$R1_alns/keys %{\\$R1_alns}/g\" /home/lam/Documents/Software/robinmeyers/transloc_pipeline/lib/TranslocFilters.pl\n",
    "sed -i \"s/keys \\$R2_alns/keys %{\\$R2_alns}/g\" /home/lam/Documents/Software/robinmeyers/transloc_pipeline/lib/TranslocFilters.pl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.2) Install source paths"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "gedit ~/.bashrc\n",
    "### in ~/.bashrc (end of file), add these\n",
    "export PATH=$PATH:/home/lam/Documents/Software/robinmeyers/transloc_pipeline/bin\n",
    "export PATH=$PATH:/home/lam/Documents/Software/robinmeyers/transloc_pipeline/R\n",
    "export GENOME_DB=/home/lam/Documents/Software/robinmeyers/genomes\n",
    "export BOWTIE2_INDEXES=/home/lam/Documents/Software/robinmeyers/genomes/bowtie2_indexes\n",
    "### out ~/.bashrc\n",
    "source ~/.bashrc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.3) Install Bowtie2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install bowtie2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.4) Set up a genome"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The LAM-HTGTS technique needs a reference genome to map your experimental reads. If you used multiple genomes or/and modified genomes, there are all require."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.4.1) Classic Genome"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/Software/robinmeyers\n",
    "mkdir genomes\n",
    "mkdir -p $BOWTIE2_INDEXES\n",
    "cd $BOWTIE2_INDEXES\n",
    "###For MOUSE, update the link if necessary (may 2018)\n",
    "###If you want HUMAN (ftp://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_28/GRCh38.p12.genome.fa.gz)\n",
    "###If Gencode is down or you want an other plateform, you can go to GRC or Ensembl (becareful about fasta header names)\n",
    "wget ftp://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_mouse/release_M17/GRCm38.p6.genome.fa.gz\n",
    "gunzip GRCm38.p6.genome.fa.gz\n",
    "mv GRCm38.p6.genome.fa mm10.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.4.2) Modified Genome"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "If you have sequenced a modified animal, a modified genome will be needed. You will modify the reference genome.\n",
    "First, download the reference genome you want to modify (2.6.4.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Then, the **lookmod_genome.py** script will help you to modifiy this genome with 3 fonctions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "lookmod_genome.py -m <mode : list, look, modify> -g <genome type> -r <reference fasta file> -c <construction fasta file> -p <position> [ -s <surroundings length> -i <info file>]\n",
    "    -h or --help : display this help\n",
    "    -m or --mode : list, look, modify\n",
    "    -g or --genome : genome type (if mode is modify)\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -c or --file_construction : construction fasta file (if mode is modify)\n",
    "    -p or --position : position chrmX:start:end (if mode is look)\n",
    "    -s or --surrondings : surroundings length (default : 10) (if mode is look)\n",
    "    -i or --info : output in file the terminal ouput"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "**How it works :**<br>\n",
    "A script allowing to modify a reference genome according to your experimental genome.<br>\n",
    "3 functions :<br>\n",
    "* **List** : Will allow you to have a look on chromosomes present in your reference genome file<br>\n",
    "* **Look** : Will allow you to have a look around the area you want to modify. Before any modification you can have a look at the position you want to insert or delete in the reference genome. If you want to insert something in chr12 between 113 000 000 and 114 000 000, you will have a look before and after this position.<br> Ex : If you want to look at the surroundings spot of 9 bases let's say between 12 and 22 on chr1, where you want to insert a sequence, you can write<br> python lookmod_genome.py -m look -f file.fasta -p chr1:12:22 -s 10<br> You will get the following area 12-9, so 3, to 22+9, so 31<br> from3to12\\**from22to31<br>\n",
    "* **Modify** : Will allow you to modify a reference genome. This method needs a contruction file (option -c). The trick is in this file.<br> In the file, information is carried in fasta format, with headers and sequences.The header contains the type of modification you want to apply, then \":\", following by the location of this modification (chr:start:end). 4 Types of modification are supported :<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**deletion :** Delete a part of the reference genome (Ex: chr4 from 40 to 50, will delete base from 41 to 49), no sequence needed<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;>deletion:chr4:40:50<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**insertion :** Insertion of sequence in the reference genome (Ex: chr1 at position 25), don't forget your sequence below<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;>insertion:chr1:25:26<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Info : You can also replace with insertion type, just enlarge your location, like :<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;>insertion:chr1:25:35<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**add :** Add a new chromosome to the reference genome, don't forget your sequence below<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;>add:your_chromosome_name<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;**remove :** Remove a chromosome from the reference genome, no sequence needed<br>\n",
    "&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;>remove:your_chromosome_name<br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Info :<br>\n",
    "* You had to combine DELETION, INSERTION, ADD and REMOVE in the same fasta file<br>\n",
    "* All positions will be calculated on reference genome, this way a modification will not impact the genome position for the next modification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "I saved my construction fasta file under /media/sf_LAM/Data/mm10_Jk_construction.fasta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "###### List "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "python /home/lam/Documents/LAM/Scripts/lookmod_genome.py -m list -r mm10.fa -i output_list_chromosomes_mm10.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "###### Look"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "python /home/lam/Documents/LAM/Scripts/lookmod_genome.py -m look -r mm10.fa -p chr12:113000000:114000000 -s 100 -i output_look_chr12_113000000_114000000_mm10.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "###### Modify"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "python /home/lam/Documents/LAM/Scripts/lookmod_genome.py -m modify -g mm10_Jk -r mm10.fa -c /media/sf_LAM/Data/mm10_Jk_construction.fasta -i output_modify_mm10_Jk_construction_mm10.txt\n",
    "###Result : mm10_Jk.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.4.3) Check modification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You can check your insertion length, doing quick math, using the look_genome.py mode list to get the length of the non-modified chromosome and modified chromosome."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.5) Filtering out no \"conventionnal\" chromosomes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To keep the most junctions as possible only \"conventional chromosomes\" (readable with Circos plot) will be kept. Using to the GRC documentation (https://www.ncbi.nlm.nih.gov/grc/help/definitions ), chromosomes and contigs outside the conventional assemblies will be deleted from the reference genome.<br>\n",
    "\n",
    "Reads will not have the possibility to be align on fixed or patched chromosomes. This action will \"force\" reads to be align on the conventional chromosomes (to be correct : will not have the possibility to be align on alternative chromosomes). So in case a read is aligned on a secondary assembly but also have a good match on conventional chromosome, this solution will save a read that we would have lost downstream in the analysis (circos_plot, karyotype_plot, locus detection... all focus on conventional chromosomes)\n",
    "\n",
    "This is a manipulation to beware of. This works in that particular case because we want to know were is the read on our genome, not the intrinsic read composition. Typicaly, in SNPs discovery the step is very dangerous and could lead to false positive discovery.<br>\n",
    "\n",
    "Furthermore, chrM known as mitochondrial chromosome has been remove due to the almost nul probability of junctions with this chromosome. Mitochondrial DNA is very difficult to access due to the mitochondrial membrane"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The **filter_reference_genome.py** will just need the number of chromosome you have in your reference genome and if you want to keep or not sexual chromosomes in your analysis."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "filter_reference_genome.py -r <reference fasta file> -a <autosome number> [-s -i <info file>]\n",
    "    -h or --help : display this help\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -a or --autosomes : number of automome in the reference genome\n",
    "    -s or --sexual : remove sexual chromosomes\n",
    "    -i or --info : output in file the terminal ouput"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "cd $BOWTIE2_INDEXES\n",
    "python /home/lam/Documents/LAM/Scripts/filter_reference_genome.py -r mm10.fa -a 19 -i output_chromosomes_mm10_filtered.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "cd $BOWTIE2_INDEXES\n",
    "python /home/lam/Documents/LAM/Scripts/filter_reference_genome.py -r mm10_Jk.fa -a 19 -i output_chromosomes_mm10_Jk_filtered.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.6) Create bowtie2 index"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###If you work on hg19, replace mm10 by hg19\n",
    "cd $BOWTIE2_INDEXES\n",
    "bowtie2-build mm10_filtered.fa mm10 > output_build_mm10.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###If you work on hg19, replace mm10 by hg19\n",
    "cd $BOWTIE2_INDEXES\n",
    "bowtie2-build mm10_Jk_filtered.fa mm10_Jk > output_build_mm10_Jk.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.7) Create genome folder"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###If you work on hg19, replace mm10 by hg19\n",
    "cd $BOWTIE2_INDEXES\n",
    "mkdir -p $GENOME_DB/mm10\n",
    "bowtie2-inspect mm10 > $GENOME_DB/mm10/mm10.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###If you work on hg19, replace mm10 by hg19\n",
    "cd $BOWTIE2_INDEXES\n",
    "mkdir -p $GENOME_DB/mm10_Jk\n",
    "bowtie2-inspect mm10_Jk > $GENOME_DB/mm10_Jk/mm10_Jk.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.8) Install Samtools"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/Software/robinmeyers\n",
    "wget https://sourceforge.net/projects/samtools/files/samtools/0.1.17/samtools-0.1.17.tar.bz2\n",
    "tar xfvj samtools-0.1.17.tar.bz2\n",
    "cd samtools-0.1.17\n",
    "gedit Makefile\n",
    "### in Makefile\n",
    "###replace\n",
    "CFLAGS=\t\t-g -Wall -O2 #-m64 #-arch ppc\n",
    "###by\n",
    "CFLAGS=\t\t-g -Wall -O2 -fPIC #-m64 #-arch ppc\n",
    "### out Makefile\n",
    "sudo apt-get install libz-dev\n",
    "sudo apt-get install libncurses5-dev #libncursesw5-dev\n",
    "make\n",
    "make razip\n",
    "sudo cp samtools /usr/bin/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.9) Install Perl >= 5.16"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install libgd-perl\n",
    "sudo apt-get install libxml2-dev\n",
    "sudo apt-get install graphviz\n",
    "sudo cpan YAML\n",
    "#sudo cpan inc::latest\n",
    "sudo cpan Module::Build\n",
    "sudo cpan Bio::DB::Sam\n",
    "###Please enter the location of the bam.h and compiled libbam.a files:\n",
    "/home/lam/Documents/Software/robinmeyers/samtools-0.1.17\n",
    "sudo cpan File::Which Getopt::Long Text::CSV List::MoreUtils Data::GUID Interpolation IPC::System::Simple Storable Switch::Plain\n",
    "sudo cpan CJFIELDS/BioPerl-1.6.924.tar.gz\n",
    "###accept default options\n",
    "###[../apache_x.x/src] --> Ctrl+c (until you get the prompt line)\n",
    "###Once again\n",
    "sudo cpan CJFIELDS/BioPerl-1.6.924.tar.gz\n",
    "###Do not worry if you get an XML::DOM::XPath error. You got this error because of the Perl version. This module will be not used.\n",
    "###make: *** [test-dynamic] Error 255\n",
    "###   MIROD/XML-DOM-CPath-0.14.tar.gz\n",
    "###   ..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.10) Install ea-utils"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install libgsl-dev\n",
    "cd /home/lam/Documents/Software\n",
    "mkdir ea-utils\n",
    "cd ea-utils\n",
    "wget https://storage.googleapis.com/google-code-archive-downloads/v2/code.google.com/ea-utils/ea-utils.1.1.2-537.tar.gz\n",
    "tar xzvf ea-utils.1.1.2-537.tar.gz\n",
    "cd ea-utils.1.1.2-537\n",
    "make\n",
    "gedit ~/.bashrc\n",
    "### in ~/.bashrc\n",
    "export PATH=$PATH:/home/lam/Documents/Software/ea-utils/ea-utils.1.1.2-537\n",
    "### out ~/.bashrc\n",
    "source ~/.bashrc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.11) Install SeqPrep"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/Software\n",
    "git clone https://github.com/jstjohn/SeqPrep.git\n",
    "cd SeqPrep\n",
    "make\n",
    "gedit ~/.bashrc\n",
    "### in ~/.bashrc\n",
    "export PATH=$PATH:/home/lam/Documents/Software/SeqPrep\n",
    "### out ~/.bashrc\n",
    "source ~/.bashrc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2.9.12) Install R >= 3.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo apt-get install r-base\n",
    "sudo apt-get -y install libxml2-dev libssl-dev libcurl4-openssl-dev libmariadb-client-lgpl-dev libcairo2-dev libxt-dev\n",
    "sudo su\n",
    "echo \"deb http://www.stats.bris.ac.uk/R/bin/linux/ubuntu trusty/\" >> /etc/apt/sources.list\n",
    "apt-key adv --keyserver keyserver.ubuntu.com --recv-keys E084DAB9\n",
    "apt-get update\n",
    "apt-get upgrade\n",
    "exit\n",
    "R\n",
    "###In R\n",
    ">install.packages(c(\"curl\", \"readr\", \"argparser\", \"magrittr\", \"stringr\", \"plyr\", \"dplyr\", \"RCircos\", \"Cairo\", \"getopt\", \"seqinr\"))\n",
    "###twice yes\n",
    "###If you had to, select a mirror (France Lyon)\n",
    "> source(\"https://bioconductor.org/biocLite.R\")\n",
    "> biocLite(c(\"GenomicRanges\", \"BSgenome\", \"karyoploteR\", \"biomaRt\", \"regioneR\"))\n",
    "###enter a, then y to answer the question\n",
    "###Out R (quit())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "sudo add-apt-repository ppa:marutter/rrutter3.5\n",
    "sudo apt-get update\n",
    "sudo apt install r-api-3.5\n",
    "sudo apt-get -y install libxml2-dev libssl-dev libcurl4-openssl-dev libmariadb-client-lgpl-dev libcairo2-dev libxt-dev\n",
    "sudo su\n",
    "echo \"deb http://www.stats.bris.ac.uk/R/bin/linux/ubuntu trusty/\" >> /etc/apt/sources.list\n",
    "apt-key adv --keyserver keyserver.ubuntu.com --recv-keys E084DAB9\n",
    "apt-get update\n",
    "apt-get upgrade\n",
    "exit\n",
    "R\n",
    "###In R\n",
    ">install.packages(c(\"curl\", \"readr\", \"argparser\", \"magrittr\", \"stringr\", \"plyr\", \"dplyr\", \"RCircos\", \"Cairo\", \"getopt\", \"seqinr\"))\n",
    "###twice yes, if needed\n",
    "###If you had to, select a mirror (France Lyon)\n",
    "> source(\"https://bioconductor.org/biocLite.R\")\n",
    "> biocLite(\"BiocUpgrade\")\n",
    "> install.packages(\"BiocManager\")\n",
    "> BiocManager::install(c(\"GenomicRanges\", \"BSgenome\", \"karyoploteR\", \"biomaRt\", \"regioneR\", \"chromPlot\"))\n",
    "###Out R (quit())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**From here I will take the Puce1 example<br>\n",
    "All files are txt files and can be opened with a text editor (care to the file size)<br>**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3) Demultiplexing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you need to merge fastq.gz files (test ship and normal ship, for example)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1/\n",
    "cat normal_1.fastq.gz nano_1.fastq.gz > merged_1.fastq.gz\n",
    "cat normal_2.fastq.gz nano_2.fastq.gz > merged_2.fastq.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Your result files will be in merged.fastaq.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;3.1) Create metadata file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "To use transloc_pipeline, according to the LAM-HTGTS pipeline documentation (https://robinmeyers.github.io/transloc_pipeline/thedocs.html), a description file of your libraries is required (metadata.txt).<br>\n",
    "<b>Every element of this file is separated by a tabulation</b><br>\n",
    "So in this file, the first line (the header of the file) is always : <br><br>\n",
    "Library&nbsp;&nbsp;&nbsp;&nbsp;Assembly&nbsp;&nbsp;&nbsp;&nbsp;Chr&nbsp;&nbsp;&nbsp;&nbsp;Start&nbsp;&nbsp;&nbsp;&nbsp;End&nbsp;&nbsp;&nbsp;&nbsp;Strand&nbsp;&nbsp;&nbsp;&nbsp;MID&nbsp;&nbsp;&nbsp;&nbsp;Primer&nbsp;&nbsp;&nbsp;&nbsp;Adapter&nbsp;&nbsp;&nbsp;&nbsp;Description<br>\n",
    "\n",
    "- Library : The name of your library, (do not use space or special characters)\n",
    "- Assembly : The reference genome used to create this library (mm10, mm10_Jk, hg19...)\n",
    "- Chr : The chromosome sequenced\n",
    "- Start : First base of your primer on your reference genome (+1 using lookmod_genome.py, mode look)\n",
    "- End : Last base of your primer on your reference genome (+1 using lookmod_genome.py, mode look)\n",
    "- Strand : The strand of your sequence\n",
    "- MID : The MID sequence\n",
    "- Primer : The primer sequence\n",
    "- Adaptater : The adaptater sequence\n",
    "- Description : The description you want<br>\n",
    "\n",
    "Then, every libraries you made will follow line by line<br>\n",
    "For example : <br><br>\n",
    "<hr>\n",
    "Library&nbsp;&nbsp;&nbsp;&nbsp;Assembly&nbsp;&nbsp;&nbsp;&nbsp;Chr&nbsp;&nbsp;&nbsp;&nbsp;Start&nbsp;&nbsp;&nbsp;&nbsp;End&nbsp;&nbsp;&nbsp;&nbsp;Strand&nbsp;&nbsp;&nbsp;&nbsp;MID&nbsp;&nbsp;&nbsp;&nbsp;Primer&nbsp;&nbsp;&nbsp;&nbsp;Adapter&nbsp;&nbsp;&nbsp;&nbsp;Description<br>\n",
    "WT1-Smu&nbsp;&nbsp;&nbsp;&nbsp;mm10&nbsp;&nbsp;&nbsp;&nbsp;chr12&nbsp;&nbsp;&nbsp;&nbsp;113426679&nbsp;&nbsp;&nbsp;&nbsp;113426699&nbsp;&nbsp;&nbsp;&nbsp;-&nbsp;&nbsp;&nbsp;&nbsp;TTATCC&nbsp;&nbsp;&nbsp;&nbsp;cacacaaagactctggacct&nbsp;&nbsp;&nbsp;&nbsp;ACCACGCGTGCCCTATAGT&nbsp;&nbsp;&nbsp;&nbsp;WT-Smu replicat 1<br>\n",
    "KIKS1-Smu&nbsp;&nbsp;&nbsp;&nbsp;mm10&nbsp;&nbsp;&nbsp;&nbsp;chr12&nbsp;&nbsp;&nbsp;&nbsp;113426679&nbsp;&nbsp;&nbsp;&nbsp;113426699&nbsp;&nbsp;&nbsp;&nbsp;-&nbsp;&nbsp;&nbsp;&nbsp;TTAAGG&nbsp;&nbsp;&nbsp;&nbsp;cacacaaagactctggacct&nbsp;&nbsp;&nbsp;&nbsp;ACCACGCGTGCCCTATAGT&nbsp;&nbsp;&nbsp;&nbsp;KIKS-Smu replicat 1<br>\n",
    "<hr>\n",
    "\n",
    "If you failed at your start/end positions or at your library names, the preprocessing will fail.<br>\n",
    "\n",
    "Info : You can also create a csv file with LibreOffice or Excel, do not forget to save as csv file using tabulation delimiter\n",
    "\n",
    "Save the metadata file in /media/sf_LAM/Data/Puce1 directory (metadata.txt)\n",
    "\n",
    "**The metadata file will be the most important option in quite all scripts along this workflow. We will manage and filter files in different ways according to the reference genome set up in metadata file**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;3.2) Run demultiplexing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "1h for 13 000 000 reads (16Go RAM)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "TranslocPreprocess.pl metadata.txt preprocess --read1 merged_1.fastq.gz --read2 merged_2.fastq.gz > output_preprocess.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4) Quality check"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Reads quality check is a vital step to ensure that libraries do not have abnormal metrics as reads length, GC percent... Final results would be accessible as html file, openable with a browser"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;4.1) Install FastQC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/Software\n",
    "wget https://launchpad.net/ubuntu/+archive/primary/+files/fastqc_0.11.5+dfsg-3_all.deb\n",
    "sudo dpkg -i fastqc_0.11.5+dfsg-3_all.deb\n",
    "###If you got an error about bad installation package, you can try\n",
    "sudo apt --fix-broken install\n",
    "###Then redo the above command\n",
    "sudo apt-get install -f"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;4.2) Install MultiQC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/Software\n",
    "git clone https://github.com/ewels/MultiQC.git\n",
    "cd MultiQC\n",
    "sudo python setup.py install"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;4.3) Run FastQC"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The **run_FastQC.py** will launch FastQC for every libraries in the metadata file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "run_FastQC.py -m <metadata file> -d <preprocess directory> -o <output directory>\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -p or --dir_prep : preprocess directory\n",
    "    -o or --dir_output : output directory"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "mkdir FastQC\n",
    "python /home/lam/Documents/LAM/Scripts/run_FastQC.py -m metadata.txt -p preprocess -o /media/sf_LAM/Data/Puce1/FastQC/ > output_FastQC.txt 2>&1 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;4.4) Run MultiQC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "multiqc -o FastQC/Overall FastQC/*zip > output_MultiQC.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;4.5) Quality overview"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can have a global overview of your quality data in /media/sf_LAM/Data/Puce1/FastQC/Overall/ opening the multiqc_report.html or either library by library in /media/sf_LAM/Data/Puce1/FastQC/ opening the library_name_fastqc.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5) Junctions detection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;5.1) Run junction detection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "4h for 13 000 000 reads (16Go RAM)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "TranslocWrapper.pl metadata.txt preprocess results > output_results.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;5.2) Remove _result.tlx file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "This could sounds weird to remove these files (_result.tlx, _result_filter.txt, _result_filter_stats.txt) but we will not need them.\n",
    "Those were create with default parameters of TranslocPipeline and we want to apply our own filters, so remove them."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1/results\n",
    "find . -type f -name '*result*' -delete"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;5.3) Alignment check"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The purpose of the alignment check is to have a look at the localisation of aligned reads over the reference genome."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "###For R1 sam file\n",
    "awk -F\"\\t\" '{print $1}' metadata.txt | sed '1d' | while read lib; do samtools view -h -F 4 -S results/$lib/$lib'_R1.sam' | samtools view -F 256 -S - | awk ' $3 !~ /@/ {print $3}' - | sort | uniq -c > results/$lib/output_check_alignment_R1.txt; done;\n",
    "###For R2 sam file\n",
    "awk -F\"\\t\" '{print $1}' metadata.txt | sed '1d' | while read lib; do samtools view -h -F 4 -S results/$lib/$lib'_R2.sam' | samtools view -F 256 -S - | awk ' $3 !~ /@/ {print $3}' - | sort | uniq -c > results/$lib/output_check_alignment_R2.txt; done;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If a \"_Inline\" folder still remains in \"/media/sf_LAM/Data/Puce1\" after this command, you can delete this directory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<hr size=\"5\">\n",
    "<b>\n",
    "As filtering and spreading junctions is inconditionnaly related to user choice and can not be automated, it is hard to predict input and output files for the next step. To try to counter this, at each filtering and spreading step, the user will be invited to put a mark on the output file, let's say \"Dup\", for duplicate (mark that will be present on the output file name). The user can now chose to filter the raw file or to filter his new \"Dup\" filtered file. In order to do that, the user will be invited to input markers already applied to a file to guide the program towards the good file to filter.</b><br>\n",
    "<b>\n",
    "All manipulation junctions and visualization programs (part 6 and part 7) will have an option (-o) to apply a new mark and an option (-i) to list marks already used on this file (Please, do not change file names, that will trigger the algorith).<b><br>\n",
    "    \n",
    "- -o : You can add every word you want, as long as you want (something useful to remember what you did). Do not use space, underscore, hyphen... To glue multiple words in one use this strategy : UseMoreThanOneWord\n",
    "- -i : On this option the user will have to remember all marks applied to the file to filter.\n",
    "\n",
    "<b>\n",
    "All marks will be separated by comma like this : Dup,UseMoreThanOneWord. All marks have to be set in the good order from oldest to lastest. To filter the raw file, just delete the option (-i) from your command line\n",
    "</b><br><br>\n",
    "\n",
    "<b>\n",
    "Also, in some scripts the construction file is an option. The file is what you have modified from the reference genome. If you add this construction file in the command line, the script will modify junction positions according to the construction file, if you don't it will not.\n",
    "</b>\n",
    "<br>\n",
    "\n",
    "<hr size=\"5\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6) Junctions manipulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1) Junctions filtering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After junctions detection from TranslocWrapper.pl, the script output a file which contains all junctions (name_of_your_library.tlx).<br>\n",
    "This file needs to be filtered to keep only worthwhile junctions.<br>\n",
    "\n",
    "The python script **filter_transloc_results.py** is an options manager to use the TranslocFilter.pl script from https://robinmeyers.github.io/transloc_pipeline/thedocs.html documentation.<br>\n",
    "\n",
    "These options could have a binary value as 0 and 1 (0 the read has not this status, 1 the read has this status) or non binary value (for example, X number of base aligned after the bait).<br>\n",
    "\n",
    "Remember that the script will filter out reads that matched with the option you set up.<br>\n",
    "For exemple, if you do not want duplicate, set up option duplicate to 1 will remove reads with the status duplicate.<br>\n",
    "For options with non binary value, with X the threshold you want to set up, you can filter out the read with value bigger than X with letter G before X, like G10, or you can filter out the read with value smaller than X with letter L before X, like L10.<br>\n",
    "Last thing you can do with the non binary value, if you want to filter out a specific value, just write your X without G or L.\n",
    "\n",
    "To not use a filter, do not write it in the command line.<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "filter_transloc_results.py -m <metadata file> -g <genome type> -o <output mark> -t <results directory> [-i <input mark> --unaligned <value> --baitonly <value> --uncut <value> --misprimed <value> --freqcut <value> --largegap <value> --mapqual <value> --breaksite <value> --sequential <value> --repeatseq <value> --duplicate <value>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter libraries results with this genome\n",
    "    -o or --output_mark : mark added to the output file name\n",
    "    -t or --dir_results : results directory\n",
    "    -i or --input_mark : marks from input file\n",
    "    --unaligned : No OCS alignments (0 or 1)\n",
    "    --baitonly : Bait alignment is either the only alignent in the OCS or only followed by adapter alignment (0 or 1)\n",
    "    --uncut : Bait alignment runs greater than some number of bases past the cutsite (L or G, then your value)\n",
    "    --misprimed : Bait alignment runs fewer than some number of bases past the primer (L or G, then your value)\n",
    "    --freqcut : Restriction enzyme site within some number of bases from the junction (0 or 1)\n",
    "    --largegap : More than some number of bases between the bait and prey alignments (L or G, then your value)\n",
    "    --mapqual : OCS had a competing prey junction (0 or 1)\n",
    "    --breaksite : Prey alignment maps into non-endogenous breaksite cassette (0 or 1)\n",
    "    --sequential : Junction occurs downstream on read from first bait-prey junction (0 or 1)\n",
    "    --repeatseq : Repeated sequences (0 or 1)\n",
    "    --duplicate : Duplicate junctions (0 or 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b>In our study, we copy a part of chr12 to chr6 in the reference genome called mm10_Jk<br>\n",
    "This lead to commun genome parts between chr12 and ch6 on some position which results in competing prey junctions.<br>\n",
    "To avoid this issue, we do not filter out mapqual for libraries with mm10_Jk reference genome, but we will further manually filter out mapqual that does not apply on commun region in chr12 and chr6 (part 6.1.3).</b>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b>We will keep both, duplicate junctions and no duplicate junctions in separated files</b>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.1) Filter junctions, keep duplicate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### mm10 genome"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_transloc_results.py -m metadata.txt -g mm10 -o Dup -t results --unaligned 1 --baitonly 1 --misprimed L10 --freqcut 1 --largegap G30 --mapqual 1 --breaksite 1 --sequential 1 --repeatseq 1 > output_mm10_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### mm10_Jk genome"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_transloc_results.py -m metadata.txt -g mm10_Jk -o Dup -t results --unaligned 1 --baitonly 1 --misprimed L10 --freqcut 1 --largegap G30 --breaksite 1 --sequential 1 --repeatseq 1 > output_mm10_Jk_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.2) Filter junctions, remove duplicate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### mm10 genome"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_transloc_results.py -m metadata.txt -g mm10 -o NoDup -t results --unaligned 1 --baitonly 1 --misprimed L10 --freqcut 1 --largegap G30 --mapqual 1 --breaksite 1 --sequential 1 --repeatseq 1 --duplicate 1 > output_mm10_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### mm10_Jk genome"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_transloc_results.py -m metadata.txt -g mm10_Jk -o NoDup -t results --unaligned 1 --baitonly 1 --misprimed L10 --freqcut 1 --largegap G30 --breaksite 1 --sequential 1 --repeatseq 1 --duplicate 1 > output_mm10_Jk_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.3) Modified genome :  filter out non-essential mapqual junctions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we have a modified genome and duplicated parts in chr12 and chr6, we removed the mapqual filter from the filtering part for libraries mm10_Jk. But competing prey junctions can only append (in theory) between those two chromosomes, to be precise between duplicated parts only.<br>\n",
    "The following script will filter out mapqual junctions outside duplicated sequences boundaries."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.3.1) Duplicate locus file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You will have to create a file with all the duplicated sequences you have created in your modified genome (and only the one you have created). Junctions falling into donor areas will be translated to acceptor area<br>\n",
    "**The duplicate file do not need a header file** but here are the column of this file :<br><br>\n",
    "Chr_donor&nbsp;&nbsp;&nbsp;&nbsp;Start_donor&nbsp;&nbsp;&nbsp;&nbsp;End_donor&nbsp;&nbsp;&nbsp;&nbsp;Strand_donor&nbsp;&nbsp;&nbsp;&nbsp;Chr_acceptor&nbsp;&nbsp;&nbsp;&nbsp;Start_acceptor&nbsp;&nbsp;&nbsp;&nbsp;End_acceptor&nbsp;&nbsp;&nbsp;&nbsp;Strand_acceptor&nbsp;&nbsp;&nbsp;&nbsp;Area_name<br>\n",
    "\n",
    "- Chr_donor : From chromosome\n",
    "- Start_donor : The start position\n",
    "- End_donor : The end position\n",
    "- Strand_donor : The strand of the donor area\n",
    "- Chr_acceptor : To chromosome\n",
    "- Start_acceptor : The start position\n",
    "- End_acceptor : The end position\n",
    "- Strand_acceptor : The strand of the acceptor area\n",
    "- Area_name : The name fo the area\n",
    "\n",
    "Every donor/acceptor couple will follow line by line, items separated by tabulation<br>\n",
    "For example : <br><br>\n",
    "<hr size=\"1\">\n",
    "chr12&nbsp;&nbsp;&nbsp;&nbsp;113423000&nbsp;&nbsp;&nbsp;&nbsp;113425000&nbsp;&nbsp;&nbsp;&nbsp;-&nbsp;&nbsp;&nbsp;&nbsp;chr6&nbsp;&nbsp;&nbsp;&nbsp;70724000&nbsp;&nbsp;&nbsp;&nbsp;70726000&nbsp;&nbsp;&nbsp;&nbsp;+&nbsp;&nbsp;&nbsp;&nbsp;area1<br>\n",
    "chr6&nbsp;&nbsp;&nbsp;&nbsp;69378000&nbsp;&nbsp;&nbsp;&nbsp;69379000&nbsp;&nbsp;&nbsp;&nbsp;+&nbsp;&nbsp;&nbsp;&nbsp;chr6&nbsp;&nbsp;&nbsp;&nbsp;70727000&nbsp;&nbsp;&nbsp;&nbsp;70728000&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;+&nbsp;&nbsp;&nbsp;area2<br>\n",
    "<hr size=\"1\">\n",
    "\n",
    "Save your duplicate file in /media/sf_LAM/Data (mm10_Jk_duplicate_locus.txt)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Becareful, positions you write in this file have to be relative to your modified genome not Ensembl or UCSC browser.<br>\n",
    "Example, if you have inserted 700 bases chr6 position 70 000 000, for each locus chr6 after 70 000 000 you will have to add 700 to all positions. If you have inserted multiple sequences, you will have to do the calculation.\n",
    "\n",
    "Tips : Blast all your locus against your modified genome to find locations and refine positions, if necessary with lookmod_genome.py, mode look.\n",
    "You can do a nucleotide blast in command line or using the blast API on NCBI ( https://blast.ncbi.nlm.nih.gov/Blast.cgi?PROGRAM=blastn&PAGE_TYPE=BlastSearch&LINK_LOC=blasthome ).<br>\n",
    "Unfortunatly, the blast API can not hundle a full modified genome, it is too heavy. The solution will be to crop your modified genome. Let's say you know that all duplications take place chr6 from 70 720 000 to 70 740 000 and chr12 from 113 000 000 to 114 000 000. Extract these positions from your modified genome using bedtools, you will get your sub-modified genome with chr6 : 20 000 bases and chr12 : 1 000 000 bases. Then blast all your locus against your sub-modified genome on blast API.<br>\n",
    "\n",
    "To create this subreference genome, make a new file (sub_modified_genome.bed) in $BOWTIE2_INDEXES path and write the sequences parts where duplicates take place, separated by tabulation:<br><br>\n",
    "<hr>\n",
    "chr6&nbsp;&nbsp;&nbsp;&nbsp;70720000&nbsp;&nbsp;&nbsp;&nbsp;70740000<br>\n",
    "chr12&nbsp;&nbsp;&nbsp;&nbsp;113000000&nbsp;&nbsp;&nbsp;&nbsp;114000000\n",
    "<hr>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "bedtools getfasta -fi mm10_Jk_filtered.fa -bed sub_modified_genome.bed -fo sub_1M_mm10_Jk_filtered.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Then on the blast website, check \"Align two or more sequences\". Enter your locus as query sequence and as subect sequences select your sub-genome."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.3.2) Filter out mapqual junctions outside duplicate locus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Once you have created the mm10_Jk_duplicate_contruction file (/media/sf_LAM/Data), then process the script **filter_non_essential_mapqual.py**<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "filter_non_essential_mapqual.py -m <metadata file> -g <genome type> -o <output mark> -t <results directory> -n <duplicate file> [-i <input mark>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file name\n",
    "    -t or --dir_results : results directory\n",
    "    -n or --file_duplicate : duplicate position from the experimental construction (chr,start,end)\n",
    "    -i or --input_mark : marks from input file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.3.2.1) File with Duplicate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_non_essential_mapqual.py -m metadata.txt -g mm10_Jk -o FilOutMapqual -t results -n ../mm10_Jk_duplicate_locus.txt -i Dup > output_mm10_Jk_FilOutMapqual_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.1.3.2.2) File with NoDuplicate"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/filter_non_essential_mapqual.py -m metadata.txt -g mm10_Jk -o FilOutMapqual -t results -n ../mm10_Jk_duplicate_locus.txt -i NoDup > output_mm10_Jk_FilOutMapqual_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2) Junctions spreading"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Junctions files are now clean from unessessary junctions. The next step is to spread these jonctions to different groups. 3 groups of junctions have been made :\n",
    "- Legitimate junctions : junctions filling your construction, what you were waiting for\n",
    "- Illegitimate junctions : junctions outside your construction but still biologicaly coherent\n",
    "- Trash junctions : junction that can not occur biologicaly, or junctions created by the construction\n",
    "\n",
    "Using these 3 groups of junctions, 4 files will be output :\n",
    "- The legitimate file : contains all the legitimate junctions with locus information for each junctions (if locus option is added to the command line). The script will affect to each junction the closest locus and the junction length. All legitimate junctions have to fall into legitimate locus file positions. If duplicate file is add, junctions falling in the duplicate areas will be moved accordingg to the file ; strand could also be change according to this file.\n",
    "- The illegitimate file : contains illegitimate junctions not in the trash file. If a RepeatMasker file is given each junction will be confront to the RepeatMasker file to find a matching event to the junction.\n",
    "- The percent illegitimate file : contains an aggregate of junctions close to X base, called megajunction (X given by the user), plus the percentage of junctions embeddied in the megajunction. Then, with the species, the release number of the reference genome and Ensembl database, each megajunction will be confront to Ensembl database to try to match gene coordinates. Gene anmes will be reported if any.\n",
    "- The trash file : contains others junctions coming from the construction itself or impossible junctions. By example, junctions below the min_gap parameters\n",
    "\n",
    "According to your construction, you could have, in your reference genome, some duplicated locus from a chromosome (donor) to an other one (acceptor). It is possible to input a duplicate file, to modify junction falling into the donor position to acceptor position, see further. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1) Requirements"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.1) RepeatMasker"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "RepeatMasker database contains pool of information about low complexity and repeat events in multiple species. Database files are text file with chromosome information, start and end position and the event engaged in the locus."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- Go to http://www.repeatmasker.org/genomicDatasets/RMGenomicDatasets.html\n",
    "- Select your species, here mouse (if you can't find it, look for an other database or process without it)\n",
    "- Right click on the last update (or the one you want) of your_species.fa.out.gz (here, mm10.fa.out.gz)\n",
    "- copy link adress (http://www.repeatmasker.org/genomes/mm10/RepeatMasker-rm405-db20140131/mm10.fa.out.gz)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.1.1) Download"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data\n",
    "wget http://www.repeatmasker.org/genomes/mm10/RepeatMasker-rm405-db20140131/mm10.fa.out.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.1.2) Modification"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this part I did not try to modify something else than mm10 and hg38 files. Maybe you will have a different file with an other species so the following command will maybe not work.<br>\n",
    "At the end we need a file with 4 columns (sequence, begin, end and class/family as header)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "# Change mm10 by your species name\n",
    "cd /media/sf_LAM/Data\n",
    "gzip -d mm10.fa.out.gz\n",
    "mv mm10.fa.out repeatMasker_mm10.fa.out\n",
    "sed 1d repeatMasker_mm10.fa.out | sed -e '1s/(left)\\s\\s\\s\\srepeat/(left)    bad_row    repeat/g' - | sed 's/^[ \\t]*//' - | sed -e \"s/[[:space:]]\\+/\\t/g\" - | sed '/^\\s*$/d' - | awk -F \"\\t\" '{print $5\"\\t\"$6\"\\t\"$7\"\\t\"$11}' - > repeatMasker_Filtered_mm10.csv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.2) Locus file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The locus file is one of the most important file in this workflow. A locus file keeps information about all loci of interest.<br>\n",
    "**The locus file do not need a header file** but here are the column of this file :<br><br>\n",
    "Chr&nbsp;&nbsp;&nbsp;&nbsp;Start&nbsp;&nbsp;&nbsp;&nbsp;End&nbsp;&nbsp;&nbsp;&nbsp;Strand&nbsp;&nbsp;&nbsp;&nbsp;Locus_name<br>\n",
    "\n",
    "- Chr : The locus chromosome\n",
    "- Start : The start position of the locus\n",
    "- End : The end position of the locus\n",
    "- Strand : The strand of the locus\n",
    "- Locus_name : The name of the locus\n",
    "\n",
    "Every locus will follow line by line, items separated by tabulation<br>\n",
    "For example : <br><br>\n",
    "<hr size=\"1\">\n",
    "chr12&nbsp;&nbsp;&nbsp;&nbsp;113335000&nbsp;&nbsp;&nbsp;&nbsp;113336000&nbsp;&nbsp;&nbsp;&nbsp;-&nbsp;&nbsp;&nbsp;&nbsp;locus1<br>\n",
    "chr6&nbsp;&nbsp;&nbsp;&nbsp;70100000&nbsp;&nbsp;&nbsp;&nbsp;70200000&nbsp;&nbsp;&nbsp;&nbsp;+&nbsp;&nbsp;&nbsp;&nbsp;locus2<br>\n",
    "<hr size=\"1\">\n",
    "\n",
    "Note that, as the locus name, for the karyo plot representation (7.3.2) you will be able to use a \"greek\" option to convert greek word to greek letter (E.g \"gamma\" will become \"γ\")\n",
    "\n",
    "Save your locus file in /media/sf_LAM/Data (mm10_locus.txt)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.3) Legitimate locus file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The legitimate locus file contains informations about the legitimate areas. Junctions inside these areas will be send to legitimate file.<br>\n",
    "**The legitimate locus file do not need a header file** but here are the column of this file :<br><br>\n",
    "Chr_bait&nbsp;&nbsp;&nbsp;&nbsp;Chr_prey&nbsp;&nbsp;&nbsp;&nbsp;Start&nbsp;&nbsp;&nbsp;&nbsp;End&nbsp;&nbsp;&nbsp;&nbsp;Area_name<br>\n",
    "\n",
    "- Chr_bait : From chromosome juntion\n",
    "- Chr_prey : To chromosome juntion\n",
    "- Start : The start of the prey chromosome\n",
    "- End : The end of the prey chromosome\n",
    "- Area_name : The name of the area\n",
    "\n",
    "Every locus will follow line by line, items separated by tabulation<br>\n",
    "For example : <br><br>\n",
    "<hr size=\"1\">\n",
    "chr6&nbsp;&nbsp;&nbsp;&nbsp;chr12&nbsp;&nbsp;&nbsp;&nbsp;113000000&nbsp;&nbsp;&nbsp;&nbsp;116000000&nbsp;&nbsp;&nbsp;&nbsp;area1<br>\n",
    "chr6&nbsp;&nbsp;&nbsp;&nbsp;chr6&nbsp;&nbsp;&nbsp;&nbsp;67000000&nbsp;&nbsp;&nbsp;&nbsp;71000000&nbsp;&nbsp;&nbsp;&nbsp;area2<br>\n",
    "<hr size=\"1\">\n",
    "\n",
    "Save your legitimate file in /media/sf_LAM/Data (mm10_legitimate.txt)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.1.4) Install Ensembl python module"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "This python module will bring genes information to the percent illegitimate file. You have to choose your species (latin name), then, choose your release https://media.readthedocs.org/pdf/pyensembl/latest/pyensembl.pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "pyensembl install --release 87 --species mus_musculus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.2) Launch spreading"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "mkdir postprocess"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Once you got all the requiered files, you can got to the first postprocessing step, with **split_junctions.py**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "split_junctions.py -m <metadata file> -g <genome type> -o <output mark> -p <postprocess directory> -t <results directory> -q <legitimate file> [-i <input mark> -r <reference fasta file> -c <contruction fasta file> -s <pool size> -l <locus file> -v <species name> -a <release number> -b <gap minimum> -u <repeatMasker file> -n <duplicate file>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -p or --dir_post : postprocess directory\n",
    "    -t or --dir_results : results directory\n",
    "    -q or --file_legitimate : legitimate positions (bait, prey, start, end, flag)\n",
    "    -i or --input_mark : marks from input file\n",
    "    -r or --file_reference : reference fasta file (Mandatory is modified reference genome)\n",
    "    -c or --file_construction : construction fasta file for modified genome (in $BOWTIE2_INDEXES), this option will ajust genes position to modified genome (Mandatory is modified reference genome)\n",
    "    -s or --size_pool : the number of bases between two junctions to pool them for illegitimate junctions (Default : 100)\n",
    "    -l or --file_locus : locus bed file (locus, start, end, strand, flag)\n",
    "    -v or --species : species in latin known in Ensembl database (homo_sapiens, mus_musculus...)\n",
    "    -a or --release : release number of your species known in Ensembl database (87,86...)\n",
    "    -b or --min_gap : minimum of bases between bait and prey to call junction (Default : 15)\n",
    "    -u or --file_repeat : repeatMasker library (.csv)\n",
    "    -n or --file_duplicate : duplicate area file for modified genome, this option will translate junctions from donor site to acceptor site"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.2.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/split_junctions.py -m metadata.txt -g mm10 -o Pool100Gap15 -p postprocess -t results -q ../mm10_legitimate.txt -i Dup -s 100 -l ../mm10_locus.txt -v mus_musculus -a 87 -b 15 -u ../repeatMasker_Filtered_mm10.csv > output_mm10_Pool100Gap15_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/split_junctions.py -m metadata.txt -g mm10_Jk -o Pool100Gap15 -p postprocess -t results -q ../mm10_Jk_legitimate.txt -i Dup,FilOutMapqual -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -c ../mm10_Jk_construction.fasta -s 100 -l ../mm10_Jk_locus.txt -v mus_musculus -a 87 -b 15 -u ../repeatMasker_Filtered_mm10.csv -n ../mm10_Jk_duplicate_locus.txt > output_mm10_Jk_Pool100Gap15_FilOutMapqual_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;6.2.2.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/split_junctions.py -m metadata.txt -g mm10 -o Pool100Gap15 -p postprocess -t results -q ../mm10_legitimate.txt -i NoDup -s 100 -l ../mm10_locus.txt -v mus_musculus -a 87 -b 15 -u ../repeatMasker_Filtered_mm10.csv > output_mm10_Pool100Gap15_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/split_junctions.py -m metadata.txt -g mm10_Jk -o Pool100Gap15 -p postprocess -t results -q ../mm10_Jk_legitimate.txt -i NoDup,FilOutMapqual -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -c ../mm10_Jk_construction.fasta -s 100 -l ../mm10_Jk_locus.txt -v mus_musculus -a 87 -b 15 -u ../repeatMasker_Filtered_mm10.csv -n ../mm10_Jk_duplicate_locus.txt > output_mm10_Jk_Pool100Gap15_FilOutMapqual_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7) Visualization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.1) Visualization files build up"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some visualizations need specific data files to create complete representation. To create histograms and links visualization for circos plot but also to visualize links in karyo plot and histogram in distribution plot. In order to get full information from plot the script **generate_visualization_files.py** will be use to generate these data files."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "generate_visualization_files.py -m <metadata file> -g <genome type> -o <output mark> -r <reference fasta file> -p <postprocess directory> [-i <input mark> -s <pool size> -t <percent translocation> -b <bin size>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -p or --dir_post : postprocess directory\n",
    "    -q or --file_legitimate : legitimate positions (bait, prey, start, end, flag)\n",
    "    -i or --input_mark : marks from input file\n",
    "    -s or --size_pool : the number of bases between two junctions to pool them (Default : 100)\n",
    "    -t or --percent_transloc_circos : percentage of translocation inside a pool of junctions to be display as link on Circos plot (Default : 2.0)\n",
    "    -c or --percent_transloc_karyo : percentage of translocation inside a pool of junctions to be display as link on Karyo plot (Default : 2.0)\n",
    "    -b or --bin_size : size of bins to create histogram on Circos plot (Default : 5000000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "5 kinds of output files :\n",
    "\n",
    "- **Library_Histo...** :\n",
    "\n",
    "Used in the Circos representation. For each chromosome, a panel of bins is created. The last bin start from the length of the chromosome - \"--bin_size\" + 1 to the length of the chromosome etc... The first bin could be smaller than \"--bin_size\", starting from 1. For each bin, the count of the legitimate and illegitimate (acceptable) junctions falling into the bin is counted as 1. With the total number of acceptable junctions, for each bin, the junction count is transformed to a percentage of acceptable junctions.\n",
    "\n",
    "- **Library_Link...** :\n",
    "\n",
    "Used in the Circos representation. The aim of the file is to highlight representative groups of junction. For each junction the script will try to reach an other junction far from \"--size_pool\". If it reach a junction a megajunction is created and then the megajunction look for an other junction in the \"--size_pool\" range, if a junction is reached, the megajunction grow. When the megajunction can not be enlarge anymore (no junction found in \"--size_pool\" range), the megajunction position is saved with the number of junction involved. With the total number of acceptable junctions (legitimate and illegitimate), for each megajunction, the junction count is transformed to a percentage of acceptable junctions. Only megajunctions with more than \"--percent_transloc_circos\" percent is reported in the output file.\n",
    "        \n",
    "- **Library_Karyoplot_freq...** :\n",
    "\n",
    "Used in the Karyoplot frequency representation. Only used in the zoom_in mode of karyo_plot.R script. For each legitimate region given in the legitimate file, bins are created starting from the legitimate start position to legitimate start position plus \"--size_pool\"-1 etc... The last bin could be smaller than \"--size_pool\", starting from the end of the last bin +1 to the legitimate end position. For each bin, the count of the legitimate junctions falling into the bin is counted as 1. With the total number of legitimate junctions, for each bin, the junction count is transformed to a percentage of legitimate junctions. Each bin is created twice, one for the forward prey junction, another one for the reverse prey junction. Junction are spread over bins taking position and strandness into account. Percentage is calculated by strand ; meaning that, to get the frequency in bin A plus strand, the number of plus strand junctions falling in bin A is divided by the total plus strand junction in legitimate file.\n",
    "\n",
    "- **Library_Karyoplot_link...** :\n",
    "\n",
    "Used in the Karyoplot link representation. The aim of the file is to highlight representative groups of junction and its intensity. For each junction the script will try to reach an other junction far from \"--size_pool\". If it reach a junction a megajunction is created and then the megajunction look for an other junction in the \"--size_pool\" range, if a junction is reached, the megajunction grow. When the megajunction can not be enlarge anymore (no junction found in \"--size_pool\" range), the megajunction position is saved with the number of junction involved. With the total number of acceptable junctions (legitimate and illegitimate), for each megajunction, the junction count is transformed to a percentage of acceptable junctions. Only megajunctions with more than \"--percent_transloc_karyo\" percent is reported in the output file. This file is the same as the \"Library_Link...\", except that the percent value is output in the result file allowing to display a junction intensity.\n",
    "\n",
    "- **Library_Distribution...** :\n",
    "\n",
    "Used in the Distribution representation. A simple output of legitimate and illegitimate junctions with the prey name, start and end of a junction (start and end are similar) and the strand of the prey junction."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.1.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/generate_visualization_files.py -m metadata.txt -g mm10 -o Pool100Perc2Bin5M -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -q ../mm10_legitimate.txt -i Dup,Pool100Gap15 -s 100 -t 2 -c 2 -b 5000000 > output_mm10_Pool100Circos2Karyo2Bin5M_Pool100Gap15_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/generate_visualization_files.py -m metadata.txt -g mm10_Jk -o Pool100Perc2Bin5M -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -q ../mm10_Jk_legitimate.txt -i Dup,FilOutMapqual,Pool100Gap15 -s 100 -t 2 -c 2 -b 5000000 > output_mm10_Jk_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_Dup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.1.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/generate_visualization_files.py -m metadata.txt -g mm10 -o Pool100Circos2Karyo2Bin5M -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -q ../mm10_legitimate.txt -i NoDup,Pool100Gap15 -s 100 -t 2 -c 2 -b 5000000 > output_mm10_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/generate_visualization_files.py -m metadata.txt -g mm10_Jk -o Pool100Circos2Karyo2Bin5M -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -q ../mm10_Jk_legitimate.txt -i NoDup,FilOutMapqual,Pool100Gap15 -s 100 -t 2 -c 2 -b 5000000 > output_mm10_Jk_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2) Circos plot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Circos representation is a way to display data and information in a circular layout. This is ideal for exploring connections and relationships between objects and positions. The whole genome of the species studied will be display in a circle, chromosome by chromosome, and links will be made between positions to reflect double strand breaks.\n",
    "\n",
    "Histograms are created with \"Library_Histo...\" file (7.1). The maximum is set to 100%\n",
    "\n",
    "Links are created with \"Library_Link...\" file (7.1). Red links means intra-chromosomic megajunctions whereas blue links means inter-chromosomic megajunction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.1) Cytoband information"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- This script use data from UCSC (http://genome.ucsc.edu)\n",
    "- Go to http://hgdownload.cse.ucsc.edu/downloads.html and look for your species\n",
    "- Click on \"Annotation database\" of the last release and copy the URL of the web page (http://hgdownload.cse.ucsc.edu/goldenPath/mm10/database/)\n",
    "- Run to command below"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "<b>Example for mm10</b>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "wget ftp://hgdownload.cse.ucsc.edu/goldenPath/mm10/database/cytoBandIdeo.txt.gz\n",
    "gunzip cytoBandIdeo.txt.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "As for 2.6.5) we will filter out \"no conventional\" chromosome out from cytoband file.<br>\n",
    "If you have a modified genome, you will have to transform your cytoband file to match your modified reference genome."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "For both normal reference and modified reference, the script **modiy_cytoband.py** will manage the cytoband file to fit your data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "modify_cytoband.py -g <genome type> -r <reference fasta file> -y <cytoband file> [-c <contruction fasta file>]\n",
    "    -h or --help : display this help\n",
    "    -g or --genome : will add the genome to the cytoband file\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -y or --file_cytoband : cytoband file\n",
    "    -c or --file_construction : construction fasta file for modified genome (in $BOWTIE2_INDEXES)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You should, at the end, get a .txt file openable under a sheet table. This table contains 5 columns (Chromosome, Start, End, Band and Stain).<br>\n",
    "For some species, Band and/or Stain are not available. I suggest you to add blanck columns to the cytoband file."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.1.1) Filter normal genome (mm10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "cd /media/sf_LAM/Data/\n",
    "python /home/lam/Documents/LAM/Scripts/modify_cytoband.py -g mm10 -r $BOWTIE2_INDEXES/mm10_filtered.fa -y cytoBandIdeo.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.1.2) Modify cytoband information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "##Use the genomic construction to create new cytoband\n",
    "cd /media/sf_LAM/Data/\n",
    "python /home/lam/Documents/LAM/Scripts/modify_cytoband.py -g mm10_Jk -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -y cytoBandIdeo.txt -c mm10_Jk_construction.fasta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.2) Circos plot creation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "All requirements are completed, it is now possible to create the circos plot representation for each library, using **circos_plot.R**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "circos_plot.R -m <metadata file> -g <genome type> -o <output mark> -p <postprocess directory> -y <cytoband file> [-i <input mark>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -p or --dir_post : postprocess directory \n",
    "    -y or --file_cytoband : cytoband file\n",
    "    -i or --input_mark : marks from input file\n",
    "    -q or --max_histo : histogram maximum value (Default : 100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.2.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/circos_plot.R -m metadata.txt -g mm10 -o Circos -p postprocess -y ../mm10_cytoBandIdeo.txt -i Dup,Pool100Gap5,Pool100Circos2Karyo2Bin5M > output_mm10_Circos_Pool100Circos2Karyo2Bin5M_Pool100Gap5_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/circos_plot.R -m metadata.txt -g mm10_Jk -o Circos -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -i Dup,FilOutMapqual,Pool100Gap5,Pool100Circos2Karyo2Bin5M > output_mm10_Jk_Circos_Pool100Circos2Karyo2Bin5M_Pool100Gap5_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.2.2.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/circos_plot.R -m metadata.txt -g mm10 -o Circos -p postprocess -y ../mm10_cytoBandIdeo.txt -i NoDup,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Circos_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/circos_plot.R -m metadata.txt -g mm10_Jk -o Circos -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -i NoDup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Jk_Circos_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3) Karyo plot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Karyo plot representation is a way to display data and information in a chromosome line layout. This is ideal for exploring connections and relationships between chromosomes."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As said further option -v or --visualization allows to visualize a karyoplot under three differents point of view. \n",
    "\n",
    "- The \"**full_genome**\" option will allow you to generate a karyoplot with links over the full genome (using the \"Library_Karyoplot_link...\" file (7.1)). Green markers are bait start and red markers are prey start. Links heigth are not informative but intensity related to the frequency of acceptable junction, as said in 7.1. Scale is added automatically at chromosomes length. Genes falling in the megajunctions can be display using -z option.\n",
    "\n",
    "- The \"**selected_chromosomes**\" option will allow you to generate a karyoplot with links over the selected chromosomes set in option -a and -b (using the \"Library_Karyoplot_link...\" file (7.1)). Green markers are bait start and red markers are prey start. Links heigth are not informative but intensity related to the frequency of acceptable junction, as said in 7.1. Scale is only available for single chromosome display and added automatically to the chromosome length. Genes falling in the megajunctions can be display using -z option.\n",
    "\n",
    "- The \"**zoom_in**\" option will allow you to generate more precise visualizations. As the two others it will display a karyoplot with links over the selected chromosomes parts set in option -a and -b (using the \"Library_Karyoplot_link...\" file (7.1)). Green markers are bait start and red markers are prey start. Scale is only available for single chromosome display and added automatically to the chromosome length. Links heigth are not informative but intensity related to the frequency of acceptable junction, as said in 7.1. Locus are generated by the input locus file. Scale ticks are only generated for single chromosome visualization and are set according to the chromosome length to display. But also display the frequency karyoplot using the \"Library_Karyoplot_freq...\" file (7.1)). Locus are generated by the input locus file. Scale ticks are only generated for single chromosome visualization and are set according to the chromosome length to display. Green histograms are plus strand frequency junctions, red histograms are minus strand frequency junctions. Histograms which could be straddling on display and not display are actuallt not display. Maximum of both histogram panels are set to the maximum integer value of plus and minus frequency junctions values plus 1. Scale on the left side are set automatically according to the histogram panel maximum. The last display is the multi frequency karyoplot display. Every frequency karyoplot from the libraries set in option -w will be merge in a single frequency plot (multi frequency karyoplot). As the histogram panel maximum can not be set, the option -q allows the user to chose this maximum."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This representation highly depends on chromosome lengths. To be able to correctly display this figure, the chromosomes length will be retrieve from the reference genome."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3.1) Genome length creation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "samtools faidx $BOWTIE2_INDEXES/mm10_filtered.fa && echo -e \"chr\\tstart\\tend\\n$(awk -F \"\\t\" '{print $1\"\\t1\\t\"$2}' $BOWTIE2_INDEXES/mm10_filtered.fa.fai)\\n\" > mm10_length.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "samtools faidx $BOWTIE2_INDEXES/mm10_Jk_filtered.fa && echo -e \"chr\\tstart\\tend\\n$(awk -F \"\\t\" '{print $1\"\\t1\\t\"$2}' $BOWTIE2_INDEXES/mm10_Jk_filtered.fa.fai)\\n\" > mm10_Jk_length.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3.2) Rename chromosome"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "To get a very detailed karyo plot, it is possible to rename chromosomes you want to display, or partial chromosomes. In order to do the modification a rename file is requiered.<br>\n",
    "**The rename file do not need a header file** but here are the column of this file :<br><br>\n",
    "chr_part_name&nbsp;&nbsp;&nbsp;&nbsp;new__name<br>\n",
    "\n",
    "- chr_part_name : This has to be one of the chr_bait or chr_prey input you set up to filter junctions, in -a (--chr_bait) or/and -b (--chr_prey), if you set up a full chromosome filtering or all chromosomes (no filter) just enter the chromosome name to rename\n",
    "- new_name : The new name of the chromosome or partial chromosome. For once, spaces are accepted !\n",
    "\n",
    "Every rename couples will follow line by line, items separated by tabulation<br>\n",
    "For example : <br><br>\n",
    "<hr size=\"1\">\n",
    "chr12:113200000:113500000&nbsp;&nbsp;&nbsp;&nbsp;chr12 IgH<br>\n",
    "chr6:70723910:70735540&nbsp;&nbsp;&nbsp;&nbsp;KIKS c3'RR<br>\n",
    "<hr size=\"1\">\n",
    "\n",
    "Save your rename file in /media/sf_LAM/Data (mm10_rename_igh.txt)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3.3) Karyo plot creation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "All requirements are completed, it is now possible to create the karyo plot representation for each library, using **karyo_plot.R**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The option -k or --greek refers to part 6.2.1.2 about locus file. With this option the script will change \"gamma\" to \"γ\" or \"Gamma\" to \"Γ\" in the karyo plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "karyo_plot.R -m <metadata file> -g <genome type> -o <output mark> -r <reference fasta file> -p <postprocess directory> -y <cytoband file> -n <chromosome length file> -v <visualization type> -a <chromosome bait> -b <chromosome prey> [-i <input mark> -z <UCSC genes> -c <construction file> -l <locus file> -k <greek locus names> -u <unlink> -t <threshold> -x <rename file>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -p or --dir_post : postprocess directory \n",
    "    -y or --file_cytoband : cytoband file\n",
    "    -n or --file_chrom_length : the chromosome length file corresponding to your genome\n",
    "    -v or --visualization : visualization type for your karyo plot (full_genome, selected_chromosomes or zoom_in)\n",
    "    -a or --chr_bait : chromosomes bait to display, separated by comma, write 'all' to display all chromosomes\n",
    "    -b or --chr_prey : chromosomes prey to display, separated by comma, write 'all' to display all chromosomes\n",
    "    -i or --input_mark : marks from input file\n",
    "    -z or --UCSC : add UCSC gene loci\n",
    "    -c or --file_construction : construction fasta file for modified genome (in $BOWTIE2_INDEXES)\n",
    "    -l or --file_locus : file to set up some locus labels\n",
    "    -k or --greek : (only if -l option), transform locus name to greek name (alpha, Alpha, beta...)\n",
    "    -u or --unlink : unlink the bait1/prey1, bait2/prey2 process and allows to search for bait1/prey2 junctions\n",
    "    -x or --file_rename : file to rename chromosome to a better display\n",
    "    -w or --multi : library to display in the multi visualization, separated by comma (in the input order) (Default : 'all') (zoom_in mode only)\n",
    "    -q or --max_y : Y axis maximum for multi visualization only (Default : 100) (zoom_in mode only)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "As in request the Ensembl database in this script, please if this error occurs, check that Ensembl server and retry later :\n",
    "\n",
    "- The BioMart web service you're accessing may be down."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### **Additionnal explanation** : <br>\n",
    "- The visualization option allow to display either the complete genome (full_genome), somes chromosomes (selected_chromosomes) or a small part of chromosomes (zoom_in)<br>\n",
    "- Multiples chr_bait and chr_prey are allowed, separated by comma<br>\n",
    "- To set a part of a chromosome, set it like this : chr:start:end\n",
    "- chr_bait and chr_prey options will filter junctions according to the input. The display will be set on these input, except when larger visualization is set (full_genome or selected_chromosomes).<br>\n",
    "- A rename file is mentatory when using zoom_in option with different bait/prey location on the same chromosome. In the zoom_in option every bait/prey will generate a chromosome visualization, having multiples times the same chromosome is forbiddden.\n",
    "- Also in the zoom_in option, please do not visualize a bait inside a prey or reverse, it is not working for now.\n",
    "- Option -z may take a while if you have a lot of junctions\n",
    "\n",
    "**Exemple** : <br>\n",
    "- --chr_bait chr6, --chr_prey chr12 and --visualization full_genome, will filter junctions from chromosome 6 to chromosome 12 but will display the complete genome.<br>\n",
    "- --chr_bait chr6:70000000:71000000, --chr_prey chr12:113000000:114000000 and --visualization selected_chromosomes, will filter junctions from chromosome 6 from 70000000 to 71000000, to chromosome 12 from 113000000 to 114000000, but will display complete chromosomes 6 and 12.<br>\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3.3.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###ALL GENOME\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_All -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v full_genome -a all -b all -i Dup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -u > output_mm10_Karyo_All_Pool100Circos2Karyo2Bin5M_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr12\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_Selected_Chromosomes -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v selected_chromosomes -a chr12 -b chr12 -i Dup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -u > output_mm10_Karyo_Selected_Chromosomes_Pool100Circos2Karyo2Bin5M_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr12:113200000:113500000\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_Zoom_In -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v zoom_in -a chr12:113200000:113500000 -b chr12:113200000:113500000 -i Dup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -l ../mm10_locus.txt -k -u -x ../mm10_rename_igh.txt -w all -q 100 > output_mm10_Karyo_Zoom_In_Pool100Circos2Karyo2Bin5M_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###ALL GENOME\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_All -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v full_genome -a all -b all -i Dup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -u > output_mm10_Jk_Karyo_All_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr6-chr12\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_Selected_Chromosomes -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v selected_chromosomes -a chr6 -b chr6,chr12 -i Dup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -u > output_mm10_Jk_Karyo_Selected_Chromosomes_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr6:70723910:70735540,chr12:113200000:113500000\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_Zoom_In -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v zoom_in -a chr6:70723910:70735540 -b chr6:70723910:70735540,chr12:113200000:113500000 -i Dup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -l ../mm10_Jk_locus.txt -k -u -x ../mm10_Jk_rename_igh.txt -w all -q 100 > output_mm10_Jk_Karyo_Zoom_In_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.3.3.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###ALL GENOME\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_All -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v full_genome -a all -b all -i NoDup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -u > output_mm10_Karyo_All_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr12\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_Selected_Chromosomes -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v selected_chromosomes -a chr12 -b chr12 -i NoDup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -u > output_mm10_Karyo_Selected_Chromosomes_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr12:113200000:113500000\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10 -o Karyo_Zoom_In -r $BOWTIE2_INDEXES/mm10_filtered.fa -p postprocess -y ../mm10_cytoBandIdeo.txt -n ../mm10_length.txt -v zoom_in -a chr12:113200000:113500000 -b chr12:113200000:113500000 -i NoDup,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -l ../mm10_locus.txt -k -u -x ../mm10_rename_igh.txt -w all -q 100 > output_mm10_Karyo_Zoom_In_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###ALL GENOME\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_All -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v full_genome -a all -b all -i NoDup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -u > output_mm10_Jk_Karyo_All_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr6-chr12\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_Selected_Chromosomes -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v selected_chromosomes -a chr6 -b chr6,chr12 -i NoDup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -u > output_mm10_Jk_Karyo_Selected_Chromosomes_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###chr6:70723910:70735540,chr12:113200000:113500000\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/karyo_plot.R -m metadata.txt -g mm10_Jk -o Karyo_Zoom_In -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -n ../mm10_Jk_length.txt -v zoom_in -a chr6:70723910:70735540 -b chr6:70723910:70735540,chr12:113200000:113500000 -i NoDup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M -z -c ../mm10_Jk_construction.fasta -l ../mm10_Jk_locus.txt -k -u -x ../mm10_Jk_rename_igh.txt -w all -q 100 > output_mm10_Jk_Karyo_Zoom_In_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4) Distribution plot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The distribution representation allows to detect areas of splicing events along a genome according to the junction strandness using \"Library_Distribution...\" file. Chromosomes are display verticaly, red histograms on the left of the chromosome histograms represent minus strand junctions while on the right the side green histograms represent plus strand junctions. Scale and cytobands are set automatically with the input files."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This representation highly depends on chromosome gap informations. To be able to correctly display this figure, the genome gaps will be retrieve from the UCSC server and filtered according to your construction."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.1) Gap information"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "- This script use data from UCSC (http://genome.ucsc.edu)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "<b>Example for mm10</b>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "As for 2.6.5) we will filter out \"no conventional\" chromosome out from gap file.<br>\n",
    "If you have a modified genome, you will have to transform your cytoband file to match your modified reference genome."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "wget ftp://hgdownload.cse.ucsc.edu/goldenPath/mm10/database/gap.txt.gz\n",
    "gunzip gap.txt.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "For both normal reference and modified reference, the script **modiy_gap.py** will manage the gap file to fit your data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "modify_gap.py -g <genome type> -r <reference fasta file> -y <gap file> [-c <contruction fasta file>]\n",
    "    -h or --help : display this help\n",
    "    -g or --genome : will add the genome to the cytoband file\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -w or --file_gap : gap file\n",
    "    -c or --file_construction : construction fasta file for modified genome (in $BOWTIE2_INDEXES)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "You should, at the end, get a .txt file openable under a sheet table. This table contains 4 columns (Chromosome, Start, End and Name).<br>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.1.1) Filter normal genome (mm10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "cd /media/sf_LAM/Data/\n",
    "python /home/lam/Documents/LAM/Scripts/modify_gap.py -g mm10 -r $BOWTIE2_INDEXES/mm10_filtered.fa -w gap.txt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.1.2) Modify cytoband information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###Remove secondary assemblies from reference genome\n",
    "##Use the genomic construction to create new gap\n",
    "cd /media/sf_LAM/Data/\n",
    "python /home/lam/Documents/LAM/Scripts/modify_gap.py -g mm10_Jk -r $BOWTIE2_INDEXES/mm10_Jk_filtered.fa -w gap.txt -c mm10_Jk_construction.fasta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.2) Distribution plot creation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This distribution representation is generated with <b>distribution_plot.R</b>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "distribution_plot.R -m <metadata file> -g <genome type> -o <output mark> -p <postprocess directory> -y <cytoband file> -w <gap file> [-i <input mark>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -p or --dir_post : postprocess directory \n",
    "    -y or --file_cytoband : cytoband file\n",
    "    -w or --file_gap : gap file\n",
    "    -i or --input_mark : marks from input file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.2.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/distribution_plot.R -m metadata.txt -g mm10 -o Distribution -p postprocess -y ../mm10_cytoBandIdeo.txt -w ../mm10_gap.txt -i Dup,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Distribution_Pool100Circos2Karyo2Bin5M_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/distribution_plot.R -m metadata.txt -g mm10_Jk -o Distribution -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -w ../mm10_Jk_gap.txt -i Dup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Jk_Distribution_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.4.2.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/distribution_plot.R -m metadata.txt -g mm10 -o Distribution -p postprocess -y ../mm10_cytoBandIdeo.txt -w ../mm10_gap.txt -i NoDup,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Distribution_Pool100Circos2Karyo2Bin5M_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "/home/lam/Documents/LAM/Scripts/distribution_plot.R -m metadata.txt -g mm10_Jk -o Distribution -p postprocess -y ../mm10_Jk_cytoBandIdeo.txt -w ../mm10_Jk_gap.txt -i NoDup,FilOutMapqual,Pool100Gap15,Pool100Circos2Karyo2Bin5M > output_mm10_Jk_Distribution_Pool100Circos2Karyo2Bin5M_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5) CSReport"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An automated tool able to analize large sata sets of CSR junction sequences produced by high-throughput sequencing (HTS), to support analysis of CSR recombination junctions ( https://www.ncbi.nlm.nih.gov/pubmed/28416601 ).\n",
    "This tool output some interesting plots as the proportion of micro homology, blunt, insert junctions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.1) Extract locus sequences"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CSReport needs junction sequences as input,so in order to investigate junctions locus by locus, sequences will be dispatched by loci.<br>\n",
    "In order to control the dispatching, the script **extract_loci_sequences.py** will help"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b>In this CSReport study, for mm10_Jk we do not want to analyze information about chr12. So, only chr6 locis will be kept in mm10_Jk_locus.txt. The new file will be called mm10_Jk_locus_CSReport.txt</b>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "extract_loci_sequences.py -m <metadata file> -g <genome type> -p <postprocess directory> -l <locus file> [-i <mark input>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -p or --dir_post : postprocess directory\n",
    "    -l or --file_locus : locus bed file (locus, start, end, strand, flag)\n",
    "    -i or --input_mark : mark from input file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.1.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/extract_loci_sequences.py -m metadata.txt -g mm10 -p postprocess -l ../mm10_locus.txt -i Dup,Pool100Gap15 > output_mm10_Sequences_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/extract_loci_sequences.py -m metadata.txt -g mm10_Jk -p postprocess -l ../mm10_Jk_locus_CSReport.txt -i Dup,FilOutMapqual,Pool100Gap15 > output_mm10_Jk_Sequences_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.1.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/extract_loci_sequences.py -m metadata.txt -g mm10 -p postprocess -l ../mm10_locus.txt -i NoDup,Pool100Gap15 > output_mm10_Sequences_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/extract_loci_sequences.py -m metadata.txt -g mm10_Jk -p postprocess -l ../mm10_Jk_locus_CSReport.txt -i NoDup,FilOutMapqual,Pool100Gap15 > output_mm10_Jk_Sequences_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.2) Get CSReport"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Please, take contact with persons from CSReport publication in order to get CSReport"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/LAM/Scripts/\n",
    "###Copy/Paste CSReport repository in this folder"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Comment all plt.show() line"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /home/lam/Documents/LAM/Scripts/CSReport\n",
    "sed -i s/plt.show\\(\\)/\\plt.close\\(\\)/g CSReport.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.3) Set up minimum and maximum locus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CSReport needs a small reference genome file to map on. We are stuying junctions falling into locus file positions. So, the reference genome file will be split using the minimum start position and the maximum end position of the locus file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "###mm10\n",
    "###chr12    113225832    113255223    3'RR\n",
    "###chr12    113426697    116009954    5'bait\n",
    "min_mm10=$(awk -F $'\\t' '{if($1==\"chr12\"){print}}' /media/sf_LAM/Data/mm10_locus.txt | sort -k2,2n | awk 'NR == 1 {print $2}')\n",
    "max_mm10=$(awk -F $'\\t' '{if($1==\"chr12\"){print}}' /media/sf_LAM/Data/mm10_locus.txt | sort -k3,3rn | awk 'NR == 1 {print $3}')\n",
    "###mm10_Jk\n",
    "###chr6    67550000    70723951    5'bait\n",
    "###chr6    70733460    70735539    mEik/Ck2\n",
    "min_mm10_Jk=$(awk -F $'\\t' '{if($1==\"chr6\"){print}}' /media/sf_LAM/Data/mm10_Jk_locus_CSReport.txt | sort -k2,2n | awk 'NR == 1 {print $2}')\n",
    "max_mm10_Jk=$(awk -F $'\\t' '{if($1==\"chr6\"){print}}' /media/sf_LAM/Data/mm10_Jk_locus_CSReport.txt | sort -k3,3rn | awk 'NR == 1 {print $3}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Commands part 7.5.4 and 7.5.5 have to be done in the **same** terminal. Otherwise minimum and maximum variables would be lost."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.4) Splice reference genome into IgH locus reference file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "samtools faidx mm10_filtered.fa chr12:$min_mm10-$max_mm10 > /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.fasta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd $BOWTIE2_INDEXES\n",
    "samtools faidx mm10_Jk_filtered.fa chr6:$min_mm10_Jk-$max_mm10_Jk > /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.fasta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.5) Adapt annotation file to IgH locus reference file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "awk -v min_mm10=$min_mm10 -F $'\\t' '{if($1==\"chr12\"){$1=$4=\"\"; print $5 \"\\t\" $2 - min_mm10+1 \"\\t\" $3 - min_mm10+1}}' mm10_locus.txt > /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.annot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "awk -v min_mm10_Jk=$min_mm10_Jk -F $'\\t' '{if($1==\"chr6\"){$1=$4=\"\"; print $5 \"\\t\" $2 - min_mm10_Jk+1 \"\\t\" $3 - min_mm10_Jk+1}}' mm10_Jk_locus_CSReport.txt > /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.annot"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.6) Run CSReport"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this part, you will be invited to input loci names. To not struggle **run_CSReport.py** script, please take care to special characters as \" ' \", use \" \\\\' \" or \" \", use \"\\\\ \", otherwise, the command line will fail."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Also, it is possible to input multiples from_annot and to_annot, just separate entries by comma"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**For now, I blocked CSReport to is minimal functionality, because with small data sets CSReport can fail.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "run_CSReport.py -m <metadata file> -g <genome type> -r <reference fasta file> -p <postprocess directory> -f <annotation file> -a <from annotation> -b <to annotation> [-i <mark input> -u <unlink> -s <cluster size> -z <reverse CSReport>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -r or --file_reference : reference fasta file\n",
    "    -p or --dir_post : postprocess directory\n",
    "    -f or --file_annot : annotation bed file (locus, start, end)\n",
    "    -a or --from_annot : search junction from loci, separated by comma, write 'all' to confront from all loci\n",
    "    -b or --to_annot : search junction to loci, separated by comma, write 'all' to confront to all loci\n",
    "    -i or --input_mark : mark from input file\n",
    "    -u or --unlink : unlink the from_annot1/to_annot1, from_annot2/to_annot2 process and allows to search for from_annot1/to_annot2 junctions\n",
    "    -s or --cluster_size : count number of clustered junctions (Default : 2)\n",
    "    -z or --reverse_CSR : do normal + reverse CSReport (change from_annot and to_annot)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As I have troubble with file names, I have translated \"/\" to \"-\", \" \" to \"+\" and \"'\" to \"!\"<br>\n",
    "\n",
    "If CSReport find one or more junctions in the sequences file, a folder is generated, even if no junction is finally output. If no junctions are found, no folder will be generated."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.6.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/run_CSReport.py -m metadata.txt -g mm10 -r /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.fasta -p postprocess -f /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.annot -a muIgh -b all -i Dup,Pool100Gap15 -u -s 1 -z > output_mm10_CSReport_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/run_CSReport.py -m metadata.txt -g mm10_Jk -r /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.fasta -p postprocess -f /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.annot -a Intron\\ IgK -b all -i Dup,FilOutMapqual,Pool100Gap15 -u -s 1 -z > output_mm10_Jk_CSReport_Pool100Gap15_FilOutMapqual_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.5.6.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/run_CSReport.py -m metadata.txt -g mm10 -r /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.fasta -p postprocess -f /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_chr12.annot -a Smu -b all -i NoDup,Pool100Gap15 -u -s 1 -z > output_mm10_CSReport_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "python /home/lam/Documents/LAM/Scripts/run_CSReport.py -m metadata.txt -g mm10_Jk -r /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.fasta -p postprocess -f /home/lam/Documents/LAM/Scripts/CSReport/Ref_mm10_Jk_chr6.annot -a Intron\\ IgK -b all -i NoDup,FilOutMapqual,Pool100Gap15 -u -s 1 -z > output_mm10_Jk_CSReport_Pool100Gap15_FilOutMapqual_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.6) Break positions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The break position visualization, generated with **break_positions_plot.py** script is a representation of the total of junctions break positions occuring in libraries. This display allows to know where breaks take place in each libraries. From bait location to prey location, at every location, the number of breaks is calculated. The color intensity is related to the number of breaks according to normal counts or log counts. Red makers are related to illegitimate breaks instead black markers are related to legitimate breaks.\n",
    "\n",
    "The number of breakpoint per base is divided by the number of breakpoint of the base which received the higher number of breakpoint for a given library (legitimate + illegitimate)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### INFORMATION ###\n",
    "break_positions_plot.py -m <metadata file> -g <genome type> -o <output_mark> -p <postprocess directory> -e <method type> [-i <mark input>]\n",
    "    -h or --help : display this help\n",
    "    -m or --file_metadata : metadata file\n",
    "    -g or --genome : only filter librairies results with this genome\n",
    "    -o or --output_mark : mark added to the output file\n",
    "    -p or --dir_post : postprocess directory\n",
    "    -e (0,1) or --method (0,1) : 0 = output picture using counts, 1 = output picture using log2(counts)\n",
    "    -i or --input_mark : mark from input file\n",
    "    -w or --libraries : library to display in the break poitions visualization, separated by comma (in the input order) (Default : all)\n",
    "    -t or --type : only display one type of junction : legitimate or illegitimate (Default : all)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.6.1) Files with duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "###Normal counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10 -o BreakPosition_counts -p postprocess -e 0 -i Dup,Pool100Gap15 -w all -t all > output_mm10_breakPositions_counts_Pool100Gap15_Dup.txt 2>&1\n",
    "###Log2 counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10 -o BreakPosition_log2counts -p postprocess -e 1 -i Dup,Pool100Gap15 -w all -t all > output_mm10_breakPositions_logCounts_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "###Normal counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10_Jk -o BreakPosition_counts -p postprocess -e 0 -i Dup,FilOutMapqual,Pool100Gap15 -w all -t all > output_mm10_Jk_breakPositions_counts_Pool100Gap15_Dup.txt 2>&1\n",
    "###Log2 counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10_Jk -o BreakPosition_log2counts -p postprocess -e 1 -i Dup,FilOutMapqual,Pool100Gap15 -w all -t all > output_mm10_Jk_breakPositions_logCounts_Pool100Gap15_Dup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;7.6.2) Files without duplicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "###Normal counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10 -o BreakPosition_counts -p postprocess -e 0 -i NoDup,Pool100Gap15 -w all -t all > output_mm10_breakPositions_counts_Pool100Gap15_NoDup.txt 2>&1\n",
    "###Log2 counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10 -o BreakPosition_log2counts -p postprocess -e 1 -i NoDup,Pool100Gap15 -w all -t all > output_mm10_breakPositions_logCounts_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### mm10_Jk"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/Puce1\n",
    "###Normal counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10_Jk -o BreakPosition_counts -p postprocess -e 0 -i NoDup,FilOutMapqual,Pool100Gap15 -w all -t all > output_mm10_Jk_breakPositions_counts_Pool100Gap15_NoDup.txt 2>&1\n",
    "###Log2 counts\n",
    "python /home/lam/Documents/LAM/Scripts/break_positions_plot.py -m metadata.txt -g mm10_Jk -o BreakPosition_log2counts -p postprocess -e 1 -i NoDup,FilOutMapqual,Pool100Gap15 -w all -t all > output_mm10_Jk_breakPositions_logCounts_Pool100Gap15_NoDup.txt 2>&1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8) Merge libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you want to inspect biological replicate libraries as one, you will need to merge the libraries.<br>\n",
    "Create a new puce (PuceM) in /media/sf_LAM/Data/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/\n",
    "mkdir PuceM"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Demultiplexing and junctions discovery can be done separatly. So, the starting point will be the junction filtering.<br>\n",
    "Create a results folder inside the new puce folder and a library folder inside the results folder"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/PuceM\n",
    "mkdir results\n",
    "cd results\n",
    "mkdir new_library"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;8.1) Merge junctions files"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then copy/paste in this new folder, the output junctions file (*.tlx file) from LAM-HTGTS technique for each library you want to merge"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "### CODE ###\n",
    "cd /media/sf_LAM/Data/PuceM\n",
    "for library in results/*; do cat $library/*.tlx > $library/$(basename $library)\".tlx\" ; sed -i '1!{/^Qname/d;}' $library/$(basename $library)\".tlx\"; done;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;8.2) Create a new metadata file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Last but not least, create a new metadata.txt file as shown 3.1) corresponding to your new Puce<br>\n",
    "Only Library and Assembly are important after the junctions detection step, so MID, primer etc can be reuse from one or another library.<br>\n",
    "\n",
    "Save your metadata file in /media/sf_LAM/Data/PuceM (metadata.txt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, start over from point 6)"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Bash",
   "language": "bash",
   "name": "bash"
  },
  "language_info": {
   "codemirror_mode": "shell",
   "file_extension": ".sh",
   "mimetype": "text/x-sh",
   "name": "bash"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}