3_dsl_showcase.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Domain Specific Language (DSL) Showcase\n",
    "## Contents\n",
    "1. [Catalyst.jl](https://github.com/SciML/Catalyst.jl/tree/master) (Existing Biochemical Reaction DSL) \n",
    "2. Custom DSL for Project"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Catalyst.jl\n",
    "We will now explore Catalyst.jl for specifiying biochemical reactions in human-readable format for conversion into an Ordinary Differential Equation (ODE) system."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import Catalyst Package, Latexify (For viewing reactions)\n",
    "using Catalyst, Latexify"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Review of original transcription and translation reactions:\n",
    "If you recall our simple framework of reactions for transcription and translation are:\n",
    "\n",
    "*Transcription*:\n",
    "\\begin{align*}\n",
    "    g_i + \\pi &\\xrightleftharpoons[k_{1;i}^{-}]{k_{1;i}^{+}} X_i \\\\\n",
    "    X_i &\\xrightarrow{k^{\\theta}_{1;i}} g_i + \\pi + m_i \\\\\n",
    "    m_i &\\xrightarrow{\\delta_i} \\varnothing  \n",
    "\\end{align*}\n",
    "\n",
    "*Translation*:\n",
    "\\begin{align*}\n",
    "    m_i + \\rho &\\xrightleftharpoons[k_{2;i}^{-}]{k_{2;i}^{+}} \\Xi_i \\\\\n",
    "    \\Xi_i &\\xrightarrow{k^{\\theta}_{2;i}} m_i + \\rho + x_i \\\\\n",
    "    x^{\\omega}_i &\\xrightarrow{\\lambda} \\varnothing  \n",
    "\\end{align*}\n",
    "\n",
    "Let's look at what it looks like with the Catalyst.jl DSL!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "node_1 = @reaction_network begin\n",
    "    (k1⁺₁, k1⁻₁), g₁ + p <--> X₁ # Gene 1 binds with RNA polymerase to form the intermediate complex in a reversible reaction\n",
    "    (k1θ₁), Xᵢ --> g₁ + p + m₁ # The intermediate complex reforms gene 1 and RNAP after transcription and mRNA 1 is formed\n",
    "    (δ₁), m₁ --> ∅ # mRNA 1 decays at a rate δ₁\n",
    "    (k2⁺₁, k2⁻₁), m₁ + r <--> Ξ₁ # mRNA 1 binds with a ribosome to form the intermediate complex in a reversible reaction\n",
    "    (k2θ₁), Ξ₁ --> m₁ + r + x₁ # The intermediate complex reforms mRNA 1 and the ribosome after translation and protein 1 is formed\n",
    "    (λ), x₁ --> ∅ # Protein 1 decays at a rate λ\n",
    "end\n",
    "\n",
    "latexify(node_1)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Look how close it looks to the original reactions! But if we want to apply our resource conservation framework, we'll have to implement new code to generate the ODE system and code a new DSL."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Custom DSL For Project\n",
    "To implement the mathematical model in this project, a whole back-end \"engine\" was developed to generate the ODE system. A DSL front-end was also made to allow you to easily enter the original reactions. Let's start with the front-end."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "@TFinactivate (macro with 1 method)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Import necessary files\n",
    "include(\"DSL_backend/dsl.jl\")\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's create a node! You can see that a different macro (@createnode) is used compared to Catalyst.jl. Also, another feature is that you can specify initial conditions very easily when you create a node unlike Catalyst! The reactions specified in the DSL are processed into a custom data structure (Node type) in the back-end."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Node(:A, 4.0, 1.0, 5.0, 4.0, 1.0, 5.0, 0.04, 0.5)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Create Node\n",
    "node_1 = @createnode begin\n",
    "    A ~ g = 0.5 # Node A with initial gene concentration of 0.5\n",
    "    [k1⁺₁ = 4, k1⁻₁ = 1], g₁ + p <--> C₁ # gene (g₁) binds with polymerase (p) to create the intermediate complex C₁\n",
    "    [k1θ₁ = 5], C --> g₁ + p + m₁ # the intermediate complex reforms the gene and polymerase while giving mRNA (m₁) after transcription\n",
    "    [k2⁺₁ = 4, k2⁻₁ = 1], m₁ + r <--> X₁ # m₁ binds with a ribosome (r) to create the intermediate complex X₁\n",
    "    [k2θ₁ = 5], X₁ --> m₁ + r + x₁ # The intermediate complex gives the protein (x₁) while returning the gene and ribosome\n",
    "    [d₁ = 0.04], m₁ --> ∅ # mRNA decays at a rate of d = 0.04, protein dilution rate is dealt with later\n",
    "end"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now to specify the edges. There are a few types: activation (@activate), repression (@repress), a combination of the two (@combi), TF activation (@TFactivate), TF inactivation (@TFinactivate) and annihilation (@annihilate). \n",
    "\n",
    "When creating these, you can conveniently specify which the nodes that these edges connect. (e.g. A -> B means that node A activates node B)\n",
    "\n",
    "A key feature is to the notation to specify the relations between nodes with an arrow notation to conveniently understand your system after you create it! (\"->\" is activation, \"-!\" is repression/inactivation and \"-X-\" is annihilation)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "ProtoIndEdge(2, :I₂, :A, 1000.0, 1.0, 2, 30)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Activation\n",
    "edge_1 = @activate begin\n",
    "    A -> B # Node A activates Node B\n",
    "    [k'⁺₂ = 100, k'⁻₂ = 1], g₂ + 3x₁ <--> c₂ # gene (g₂) + 3x₁ (TF) gives activated complex (c₂ )\n",
    "end\n",
    "\n",
    "# Repression\n",
    "edge_2 = @repress begin\n",
    "    B -! A # Node B represses Node A\n",
    "    [k'⁺₂ = 100, k'⁻₂ = 1], g₁ + 3x₂ <--> c⁰₁ # gene (g₁) + 3x₂ (TF) gives repressed complex (c⁰₁) \n",
    "end\n",
    "\n",
    "# Combinatorial promoter (Simplifies notation for specifically specifying activation and repression)\n",
    "edges_3 = @combi begin\n",
    "    A -> C, B -! C  # A activates C and B represses C\n",
    "    [k⁺ₐ₃ = 100, k⁻ₐ₃ = 1], g₃ + x₁ <--> c₃ # Activation reaction\n",
    "    [k⁺ᵣ₃ = 100, k⁻ᵣ₃ = 1], g₃ + 3x₂ <--> c⁰₃ # Repression reaction\n",
    "end\n",
    "# Note that this produces 2 edges, one for activation and one for repression\n",
    "\n",
    "# Annihilation\n",
    "edge_4 = @annihilate begin\n",
    "    A -X- B # A and B annihilate each other \n",
    "    [k⁺₁₂ = 1000, k⁻₁₂ = 0.1], m₁ + m₂ --> m∅ # In this case it is mRNA annihilation (you must use m in the variable)\n",
    "end\n",
    "# Using any other symbol without \"m\" is detected as protein annihilation\n",
    "\n",
    "# TF activation\n",
    "edge_5 = @TFactivate begin\n",
    "    I₁ -> A ~ I₁ = 30 # The TF produced by A is activated by external inducer I₁. The initial value of I₁ is 30.\n",
    "    [κ⁺₁ = 1000, κ⁺₂ = 1], x⁰₁ + 2I₁ <--> x⁺₁ # Inactive TF is activated by I₁ with a cooperativity of 2\n",
    "end\n",
    "\n",
    "# TF inactivation\n",
    "edge_5 = @TFinactivate begin\n",
    "    I₂ -! A ~ I₂ = 30 # The TF produced by A is inactivated by inducer I₂\n",
    "    [κ⁺₁ = 1000, κ⁺₂ = 1], x⁺₁ + 2I₂ <--> x⁰₁ # Active TF is inactivated by I₂ with a cooperatvity of 2\n",
    "end"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### DSL Showcase (Actual Use)\n",
    "Now that we've seen the showcase of all the macros, we can now see how you would actually use it. \n",
    "\n",
    "One thing to note is that the symbols for nodes are specified in the edges by they do not mean anything yet.\n",
    "\n",
    "To link the symbols to nodes, we must create the nodes first. Now let's create a system of 3 nodes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# DSL Showcase (Actual Use)\n",
    "node_1 = @createnode begin\n",
    "    A ~ g = 0.5\n",
    "    [k1⁺₁ = 4, k1⁻₁ = 1], g₁ + p <--> C₁\n",
    "    [k1θ₁ = 5], C --> g₁ + p + m₁\n",
    "    [k2⁺₁ = 4, k2⁻₁ = 1], m₁ + r <--> X₁\n",
    "    [k2θ₁ = 5], X₁ --> m₁ + r + x₁\n",
    "    [d₁ = 0.04], m₁ --> ∅\n",
    "end\n",
    "\n",
    "node_2 = @createnode begin\n",
    "    B ~ g = 0.5\n",
    "    [k1⁺₂ = 4, k1⁻₂ = 1], g₂ + p <--> C₂\n",
    "    [k1θ₂ = 5], C₂ --> g₂ + p + m₂\n",
    "    [k2⁺₂ = 4, k2⁻₂ = 1], m₂ + r <--> X₂\n",
    "    [k2θ₂ = 5], X₂ --> m₂ + r + x₂\n",
    "    [d₂ = 0.04], m₂ --> ∅\n",
    "end\n",
    "\n",
    "node_3 = @createnode begin\n",
    "    C ~ g = 0.5\n",
    "    [k1⁺₃ = 4, k1⁻₃ = 1], g₃ + p <--> C₃\n",
    "    [k1θ₃ = 5], C₃ --> g₃ + p + m₃\n",
    "    [k2⁺₃ = 4, k2⁻₃ = 1], m₃ + r <--> X₃\n",
    "    [k2θ₃ = 5], X₃ --> m₃ + r + x₃\n",
    "    [d₃ = 0.04], m₃ --> ∅\n",
    "end\n",
    "\n",
    "# You must use the compilenodes() function to compile all the nodes into a variable\n",
    "nodes = compilenodes(node_1, node_2, node_3)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, you have to first create the nodes that have symbols as names (A, B and C) in this case. Nodes then have to be combined using the \"compilenodes()\" function. \n",
    "\n",
    "Now when we create edges, we need to use the \"edge(nodes, macro)\" or \"edges(nodes, macro)\" functions to connect the edge to the nodes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Note the use of the edge function which contains both the vector of nodes specified and the macros\n",
    "edge1 = edge(nodes, @repress begin \n",
    "    A -! B\n",
    "    [k'⁺₂ = 100, k'⁻₂ = 1], g₂ + 3x₁ <--> c₂\n",
    "end)\n",
    "\n",
    "edge2_3 =  edges(nodes, @combi begin\n",
    "    A -> C, B -! C\n",
    "    [k⁺ₐ₃ = 100, k⁻ₐ₃ = 1], g₃ + x₁ <--> c₃\n",
    "    [k⁺ᵣ₃ = 99, k⁻ᵣ₃ = 2], g₃ + 3x₂ <--> c⁰₃\n",
    "end)\n",
    "\n",
    "edge4 = edge(nodes, @annhilate begin\n",
    "    A -X- B\n",
    "    [k⁺₁₂ = 1000, k⁻₁₂ = 0.1], m₁ + m₂ --> m∅\n",
    "end)\n",
    "\n",
    "edge5 = edge(nodes, @TFinactivate begin\n",
    "    I₁ -! A ~ I₁ = 30\n",
    "    [κ⁺₁ = 1000, κ⁺₂ = 1], x⁰₁ + 2I₁ <--> x⁺₁\n",
    "end)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And with this, we can create some systems!"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.8.2",
   "language": "julia",
   "name": "julia-1.8"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.8.2"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}