tutorials/sat_mut.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Precursors!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os, subprocess\n",
    "\n",
    "if not os.path.isdir('models/heart'):\n",
    "    os.mkdir('models/heart')\n",
    "if not os.path.isfile('models/heart/model_best.tf.meta'):\n",
    "    subprocess.call('curl -o models/heart/model_best.tf.index https://storage.googleapis.com/basenji_tutorial_data/model_best.tf.index', shell=True)\n",
    "    subprocess.call('curl -o models/heart/model_best.tf.meta https://storage.googleapis.com/basenji_tutorial_data/model_best.tf.meta', shell=True)\n",
    "    subprocess.call('curl -o models/heart/model_best.tf.data-00000-of-00001 https://storage.googleapis.com/basenji_tutorial_data/model_best.tf.data-00000-of-00001', shell=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Saturation mutagenesis is a powerful tool both for dissecting a specific sequence of interest and understanding what the model learned. [basenji_sat.py](https://github.com/calico/basenji/blob/master/bin/basenji_sat.py) enables this analysis from a test set of data. [basenji_sat_vcf.py](https://github.com/calico/basenji/blob/master/bin/basenji_sat_vcf.py) lets you provide a VCF file for variant-centered mutagenesis.\n",
    "\n",
    "To do this, you'll need\n",
    " * Trained model\n",
    " * Input file (FASTA or HDF5 with test_in/test_out)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First, you can either train your own model in the [Train/test tutorial](https://github.com/calico/basenji/blob/master/tutorials/train_test.ipynb) or use one that I pre-trained from the models subdirectory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll bash the GATA4 promoter to see what motifs drive its expression. I placed a 131 kb FASTA file surrounding the GATA4 TSS in data/gata4.fa, so we'll use [basenji_sat.py](https://github.com/calico/basenji/blob/master/bin/basenji_sat.py).\n",
    "\n",
    "The most relevant options are:\n",
    "\n",
    "| Option/Argument | Value | Note |\n",
    "|:---|:---|:---|\n",
    "| -g | | Plot the nucleotides proportional to the gain score, too. |\n",
    "| -f | 20 | Figure width, that I usually scale to 10x the saturation mutageneis region |\n",
    "| -l | 200 | Saturation mutagenesis region in the center of the given sequence(s) |\n",
    "| -o | gata4_sat | Outplot plot directory. |\n",
    "| --rc | | Predict forward and reverse complement versions and average the results. |\n",
    "| -t | 0,1,2 | Target indexes to analyze. |\n",
    "| params_file | models/params_small.txt | Table of parameters to setup the model architecture and optimization parameters. |\n",
    "| model_file | models/heart/model_best.tf | Trained saved model prefix. |\n",
    "| input_file | data/gata4.fa | Either FASTA or HDF5 with test_in/test_out keys. |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/Users/davidkelley/anaconda3/lib/python3.5/importlib/_bootstrap.py:222: RuntimeWarning: compiletime version 3.6 of module 'tensorflow.python.framework.fast_tensor_util' does not match runtime version 3.5\n",
      "  return f(*args, **kwds)\n",
      "/Users/davidkelley/anaconda3/lib/python3.5/site-packages/matplotlib/__init__.py:1405: UserWarning: \n",
      "This call to matplotlib.use() has no effect because the backend has already\n",
      "been chosen; matplotlib.use() must be called *before* pylab, matplotlib.pyplot,\n",
      "or matplotlib.backends is imported for the first time.\n",
      "\n",
      "  warnings.warn(_use_error_msg)\n",
      "{'cnn_dilation': [1, 1, 1, 1, 1, 2, 4, 8, 16, 32, 64, 1], 'link': 'softplus', 'loss': 'poisson', 'cnn_dense': [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0], 'cnn_dropout': 0.1, 'cnn_filter_sizes': [20, 7, 7, 7, 3, 3, 3, 3, 3, 3, 3, 1], 'optimizer': 'adam', 'target_pool': 128, 'adam_beta1': 0.97, 'adam_beta2': 0.98, 'cnn_pool': [2, 4, 4, 4, 1, 0, 0, 0, 0, 0, 0, 0], 'learning_rate': 0.002, 'num_targets': 3, 'cnn_filters': [128, 128, 192, 256, 256, 32, 32, 32, 32, 32, 32, 384], 'batch_buffer': 4096, 'batch_size': 4}\n",
      "Targets pooled by 128 to length 1024\n",
      "Convolution w/ 3 384x1 filters to final targets\n",
      "Model building time 14.755834\n",
      "2018-05-16 17:54:20.283767: I tensorflow/core/platform/cpu_feature_guard.cc:140] Your CPU supports instructions that this TensorFlow binary was not compiled to use: AVX2 FMA\n",
      "Mutating sequence 1 / 1\n"
     ]
    }
   ],
   "source": [
    "! basenji_sat.py -g -f 20 -l 200 -o output/gata4_sat --rc -t 0,1,2 models/params_small.txt models/heart/model_best.tf data/gata4.fa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The saturated mutagenesis heatmaps go into output/gata4_sat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "        <iframe\n",
       "            width=\"1200\"\n",
       "            height=\"400\"\n",
       "            src=\"output/pim1_sat/seq0_t0.pdf\"\n",
       "            frameborder=\"0\"\n",
       "            allowfullscreen\n",
       "        ></iframe>\n",
       "        "
      ],
      "text/plain": [
       "<IPython.lib.display.IFrame at 0x10eb37c50>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.display import IFrame\n",
    "IFrame('output/pim1_sat/seq0_t0.pdf', width=1200, height=400)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [default]",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}