{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "# Building a Model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "This simple example demonstrates how to create a model, create a reaction, and then add the reaction to the model.\n",
    "\n",
    "We'll use the '3OAS140' reaction from the STM_1.0 model:\n",
    "\n",
    "1.0 malACP[c] + 1.0 h[c] + 1.0 ddcaACP[c] $\\rightarrow$ 1.0 co2[c] + 1.0 ACP[c] + 1.0 3omrsACP[c]\n",
    "\n",
    "First, create the model and reaction."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "from __future__ import print_function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "from cobra import Model, Reaction, Metabolite\n",
    "# Best practise: SBML compliant IDs\n",
    "model = Model('example_model')\n",
    "\n",
    "reaction = Reaction('3OAS140')\n",
    "reaction.name = '3 oxoacyl acyl carrier protein synthase n C140 '\n",
    "reaction.subsystem = 'Cell Envelope Biosynthesis'\n",
    "reaction.lower_bound = 0.  # This is the default\n",
    "reaction.upper_bound = 1000.  # This is the default"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "We need to create metabolites as well. If we were using an existing model, we could use `Model.get_by_id` to get the appropriate Metabolite objects instead."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "ACP_c = Metabolite(\n",
    "    'ACP_c',\n",
    "    formula='C11H21N2O7PRS',\n",
    "    name='acyl-carrier-protein',\n",
    "    compartment='c')\n",
    "omrsACP_c = Metabolite(\n",
    "    '3omrsACP_c',\n",
    "    formula='C25H45N2O9PRS',\n",
    "    name='3-Oxotetradecanoyl-acyl-carrier-protein',\n",
    "    compartment='c')\n",
    "co2_c = Metabolite('co2_c', formula='CO2', name='CO2', compartment='c')\n",
    "malACP_c = Metabolite(\n",
    "    'malACP_c',\n",
    "    formula='C14H22N2O10PRS',\n",
    "    name='Malonyl-acyl-carrier-protein',\n",
    "    compartment='c')\n",
    "h_c = Metabolite('h_c', formula='H', name='H', compartment='c')\n",
    "ddcaACP_c = Metabolite(\n",
    "    'ddcaACP_c',\n",
    "    formula='C23H43N2O8PRS',\n",
    "    name='Dodecanoyl-ACP-n-C120ACP',\n",
    "    compartment='c')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Adding metabolites to a reaction requires using a dictionary of the metabolites and their stoichiometric coefficients. A group of metabolites can be added all at once, or they can be added one at a time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'ddcaACP_c + h_c + malACP_c --> 3omrsACP_c + ACP_c + co2_c'"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reaction.add_metabolites({\n",
    "    malACP_c: -1.0,\n",
    "    h_c: -1.0,\n",
    "    ddcaACP_c: -1.0,\n",
    "    co2_c: 1.0,\n",
    "    ACP_c: 1.0,\n",
    "    omrsACP_c: 1.0\n",
    "})\n",
    "\n",
    "reaction.reaction  # This gives a string representation of the reaction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "The gene_reaction_rule is a boolean representation of the gene requirements for this reaction to be active as described in [Schellenberger et al 2011 Nature Protocols 6(9):1290-307](http://dx.doi.org/doi:10.1038/nprot.2011.308). We will assign the gene reaction rule string, which will automatically create the corresponding gene objects."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "frozenset({<Gene STM1197 at 0x7f2d85786898>, <Gene STM2378 at 0x7f2dc45437f0>})"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reaction.gene_reaction_rule = '( STM2378 or STM1197 )'\n",
    "reaction.genes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "At this point in time, the model is still empty"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0 reactions initially\n",
      "0 metabolites initially\n",
      "0 genes initially\n"
     ]
    }
   ],
   "source": [
    "print('%i reactions initially' % len(model.reactions))\n",
    "print('%i metabolites initially' % len(model.metabolites))\n",
    "print('%i genes initially' % len(model.genes))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "We will add the reaction to the model, which will also add all associated metabolites and genes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 reaction\n",
      "6 metabolites\n",
      "2 genes\n"
     ]
    }
   ],
   "source": [
    "model.add_reactions([reaction])\n",
    "\n",
    "# Now there are things in the model\n",
    "print('%i reaction' % len(model.reactions))\n",
    "print('%i metabolites' % len(model.metabolites))\n",
    "print('%i genes' % len(model.genes))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "We can iterate through the model objects to observe the contents"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Reactions\n",
      "---------\n",
      "3OAS140 : ddcaACP_c + h_c + malACP_c --> 3omrsACP_c + ACP_c + co2_c\n",
      "\n",
      "Metabolites\n",
      "-----------\n",
      "    co2_c : CO2\n",
      " malACP_c : C14H22N2O10PRS\n",
      "      h_c : H\n",
      "3omrsACP_c : C25H45N2O9PRS\n",
      "ddcaACP_c : C23H43N2O8PRS\n",
      "    ACP_c : C11H21N2O7PRS\n",
      "\n",
      "Genes\n",
      "-----\n",
      "STM1197 is associated with reactions: {3OAS140}\n",
      "STM2378 is associated with reactions: {3OAS140}\n"
     ]
    }
   ],
   "source": [
    "# Iterate through the the objects in the model\n",
    "print(\"Reactions\")\n",
    "print(\"---------\")\n",
    "for x in model.reactions:\n",
    "    print(\"%s : %s\" % (x.id, x.reaction))\n",
    "\n",
    "print(\"\")\n",
    "print(\"Metabolites\")\n",
    "print(\"-----------\")\n",
    "for x in model.metabolites:\n",
    "    print('%9s : %s' % (x.id, x.formula))\n",
    "\n",
    "print(\"\")\n",
    "print(\"Genes\")\n",
    "print(\"-----\")\n",
    "for x in model.genes:\n",
    "    associated_ids = (i.id for i in x.reactions)\n",
    "    print(\"%s is associated with reactions: %s\" %\n",
    "          (x.id, \"{\" + \", \".join(associated_ids) + \"}\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "Last we need to set the objective of the model. Here, we just want this to be the maximization of the flux in the single reaction we added and we do this by assigning the reaction's identifier to the `objective` property of the model."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "model.objective = '3OAS140'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "The created objective is a symbolic algebraic expression and we can examine it by printing it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-1.0*3OAS140_reverse_65ddc + 1.0*3OAS140\n",
      "max\n"
     ]
    }
   ],
   "source": [
    "print(model.objective.expression)\n",
    "print(model.objective.direction)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "which here shows that the solver will maximize the flux in the forward direction."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}