Codon Optimization using Nussinov + CAI

Background

Importance of Codon Optimization:

Codon optimization is essentially identifying the optimal RNA sequence that will encode for the creation of a protein. Proteins drive every biological function in our bodies. For example, the enzymes that break down molecules, the antibodies that fight infections, and receptors that send signals between cells, are all proteins.

The cases where producing a protein would help correct bodily functions are:

• A harmful protein is being produced (a cancer cell's surface protein)
• A benficial protein is no longer being produced (insulin in diabetes, reduced production causing aging)
• An unknown protein is in your immune system (from virus or infection)

In these cases we can either increase protein production, produce a protein to combat a negative protein or prepare your immune system for foreign proteins.

Once we identify our desired protein, we need a way to produce the protein in the desired area. We do this by injecting the strand of mRNA into the body, and letting the ribosomes synthesize the protein inside the cell. Injecting mRNA is prefered over injecting the protein itself, as it is much cheaper and it is more easily modified

This graphic shows the process of synthesizing proteins from mRNA

But there are nearly infinite ways for a mRNA sequence to encode a large protein, so how do we determine the optimal encoding. This is where codon optimization comes in, it looks to optimize:

• Translation efficiency: different organisms have different concentrations of codons, so choosing more common codons will lead to faster translation
• Structural stability: the RNA sequence determines the structure it will fold into, the more stable it is, the less likely it will degrade in the cell.

This table shows which codons encode which amino acids and we can see the degeneracy here

Algorithms

CDSfold

This was the first algorithm to optimize the RNA sequence given a target protein sequence. CDSfold builds off the Zuker algorithm to determine the rna encoding sequence that minimizes minimum free energy. This sequence is our most durable RNA structure.

The Zuker algorithm is a very complicated method that determines the structure that minimizes the free energy. Base pairs contribute differently to free energy, stacking these pairs also contribute differently, and all loops contribute differently: (a) Hairpin loop, (b) bulge loop, (c) internal loops and (d) multi-loop

CDSfold needs to maintain two DP tables to find the optimal structure, $F^{n_i, n_j}(i,j)$ stores the minimum free energy (MFE) of the sub-sequence from $i$ to $j$, where $n_i$ and $n_j$ are the boundary nucleotides:

$$F^{n_i, n_j}(i,j) = \min \begin{cases} \displaystyle\min_{n_{i+1} \,\in\, N_{i+1|n_i}} F^{n_{i+1},\, n_j}(i+1,\, j) & \text{left unpaired} \\\ \displaystyle\min_{n_{j-1} \,\in\, N_{j-1 \wedge n_j}} F^{n_i,\, n_{j-1}}(i,\, j-1) & \text{right unpaired} \\\ C^{n_i,\, n_j}(i,\, j) & \text{force pair at } (i,j) \\\ \displaystyle\min_{\substack{i < k < j \\\ n_k \in N_k,\; n_{k+1} \in N_{k+1|n_k}}} \Bigl\{F^{n_i,\, n_k}(i,\, k) + F^{n_{k+1},\, n_j}(k+1,\, j)\Bigr\} & \text{bifurcation} \end{cases}$$

Base cases:

$$F^{n_i, n_i}(i,\, i) = 0 \qquad F^{n_i, n_i}(i,\, i-1) = 0$$

$C^{n_i, n_j}(i,j)$ stores the MFE given that $n_i$ and $n_j$ form a base pair:

$$C^{n_i, n_j}(i,j) =\begin{cases}\displaystyle\min_{n_{i+1},\, n_{j-1}}\left\{ C^{n_{i+1}, n_{j-1}}(i+1, j-1) + E_{\text{stack}}\langle n_i, n_j, n_{i+1}, n_{j-1}\rangle \right\} & \text{ stacking}\\\ \displaystyle E_{\text{hairpin}}\langle n_i, n_j, n_{i+1}, n_{j-1}, v_h\rangle & \text{hairpin }\\\ \displaystyle \min_{p,\, n_p,\, n_{j-1}} \left\{ C^{n_p, n_{j-1}}(p, j-1) + E_{\text{bulge}}\langle n_i, n_p, n_{j-1}, n_j, v_l\rangle \right\} & \text{bulge}\\\ \displaystyle \min_{p,\, q,\, n_p,\, n_q,\, \psi} \left\{ C^{n_p, n_q}(p, q) + E_{\text{intloop}}\langle \psi, v_l, v_r\rangle \right\} & \text{internal loop} \\\ \displaystyle \min_{n_{i+1},\, n_{j-1}} \left\{ M^{n_{i+1}, n_{j-1}}(i+1, j-1) + \text{ML}_{\text{close}}\langle n_i, n_j\rangle \right\} & \text{multibranch}\end{cases}$$

However, CDSfold does not account for the Translation efficiency of the gerated protein sequence.

---

DERNA

DERNA expands on CDSfold by including translational efficiency to the problem including the Codon Adaptation Index. This way DERNA is optimizing both aspects of the codon optimization problem.

The codon adaptation index (CAI) uses a frequency table of the organism's codon usage, and the CAI weighs the possible codons by their relative frequencies, for a codon $x$ encoding amino acid $w_a$, the log-CAI term is:

$$\gamma(a, x) = \log \frac{g(x)}{\max_{x' \in S(w_a)} g(x')}$$

where:

• $g(x)$ is the relative usage frequency of codon $x$ in the target organism
• $S(w_a)$ is the set of all codons that encode amino acid $w_a$
• $\max_{x' \in S(w_a)} g(x')$ is the frequency of the most common codon for that amino acid

Since $g(x) \leq \max_{x'} g(x')$ always, we have $\gamma(a, x) \leq 0$ for all codons. A value of $\gamma = 0$ means codon $x$ is the most frequently used codon for that amino acid

DERNA requires 4 DP tables with an extentive amount of conditions leading to lengthy and confusing equations, but in essense DERNA is finding a secondary structure that minimizes a combination of MFE and CAI, where we can change the tradeoff parameter:

$$\lambda \cdot \text{MFE} - (1-\lambda) \cdot \text{CAI}$$

When $\lambda = 1$: pure MFE minimization.
When $\lambda = 0$: pure CAI maximization.

---

My approach (Nussinov + CAI)

To simplify the DERNA idea into an approachable mini project, my idea is to optimize the RNA sequence by using Nussinov plus the CAI term as defined earlier, where my overall function will be:

$$\lambda \cdot \text{Nussinov} + (1-\lambda) \cdot \text{CAI}$$

MFE adds the majority of the complexity, and with Nussinov I will only have to manage 1 DP table.

---

DP Setup

$N[i,\, c_i;\, j,\, c_j]$ stores the maximum score over the nucleotide range $[i, j]$, where $c_i \in S(\alpha(i))$ and $c_j \in S(\alpha(j))$ are the codons at the left and right boundaries.

Notation

• $i, j$: 0-indexed nucleotide positions ($0$ to $3m-1$)
• $\alpha(i) = \lfloor i/3 \rfloor$: amino acid index for nucleotide $i$
• $S(a)$: set of synonymous codons for amino acid $a$
• $\text{cai}(a, c) = (1-\lambda)\,\gamma(a, c)$: weighted CAI contribution of codon $c$ for amino acid $a$
• $c[p]$: nucleotide at 0-indexed position $p$ within codon $c$
• $\mathcal{G}$: set of valid Watson-Crick/wobble base pairs
• $\ell_{\min}$: minimum hairpin loop length

Base Cases

$$N[i,\, c;\, i,\, c] \;=\; N[i,\, c;\, i{-}1,\, c] \;=\; \text{cai}(\alpha(i),\, c)$$

Recurrence

$N[i, c_i; j, c_j]$ is the maximum of the four cases below.

Skip left — leave $i$ unpaired, advance the left boundary:

$$\text{skip-L} = \begin{cases} N[i{+}1,\, c_i;\, j,\, c_j] & \alpha(i{+}1) = \alpha(i) \\[6pt] \displaystyle\max_{\substack{c' \in S(\alpha(i+1)) \\ c' = c_j \;\text{if}\; \alpha(i+1) = \alpha(j)}} N[i{+}1,\, c';\, j,\, c_j] + \text{cai}(\alpha(i), c_i) & \alpha(i{+}1) \neq \alpha(i) \end{cases}$$

Skip right — leave $j$ unpaired, retract the right boundary:

$$\text{skip-R} = \begin{cases} N[i,\, c_i;\, j{-}1,\, c_j] & \alpha(j{-}1) = \alpha(j) \\[6pt] \displaystyle\max_{\substack{c' \in S(\alpha(j-1)) \\ c' = c_i \;\text{if}\; \alpha(j-1) = \alpha(i)}} N[i,\, c_i;\, j{-}1,\, c'] + \text{cai}(\alpha(j), c_j) & \alpha(j{-}1) \neq \alpha(j) \end{cases}$$

Pair — form base pair $(i, j)$; only when $(c_i[i \bmod 3],\, c_j[j \bmod 3]) \in \mathcal{G}$ and $j - i - 1 \geq \ell_{\min}$. Let $\Delta_L = [\alpha(i{+}1) \neq \alpha(i)]$ and $\Delta_R = [\alpha(j{-}1) \neq \alpha(j)]$:

$$\text{pair} = \begin{cases} N[i{+}1,\, c_i;\, j{-}1,\, c_j] + \lambda & \lnot\Delta_L,\; \lnot\Delta_R \\[6pt] \displaystyle\max_{\substack{c' \in S(\alpha(j-1)) \\ c' = c_i \;\text{if}\; \alpha(j-1) = \alpha(i)}} N[i{+}1,\, c_i;\, j{-}1,\, c'] + \lambda + \text{cai}(\alpha(j), c_j) & \lnot\Delta_L,\; \Delta_R \\[6pt] \displaystyle\max_{\substack{c' \in S(\alpha(i+1)) \\ c' = c_j \;\text{if}\; \alpha(i+1) = \alpha(j)}} N[i{+}1,\, c';\, j{-}1,\, c_j] + \lambda + \text{cai}(\alpha(i), c_i) & \Delta_L,\; \lnot\Delta_R \\[6pt] \displaystyle\max_{\substack{c'' \in S(\alpha(i+1)),\; c' \in S(\alpha(j-1)) \\ c'' = c_j \;\text{if}\; \alpha(i+1) = \alpha(j) \\ c'\; = c_i \;\text{if}\; \alpha(j-1) = \alpha(i)}} N[i{+}1,\, c'';\, j{-}1,\, c'] + \lambda + \text{cai}(\alpha(i), c_i) + \text{cai}(\alpha(j), c_j) & \Delta_L,\; \Delta_R \end{cases}$$

Bifurcation — split at cut point $k$, $i \leq k &lt; j$. Let $\beta = \alpha(k)$, $\beta' = \alpha(k{+}1)$:

$$\text{bif} = \begin{cases} \displaystyle\max_{\substack{i \leq k < j,\;\; \beta \neq \beta' \\ c_k \in S(\beta),\;\; c_{k+1} \in S(\beta') \\ c_k = c_i \;\text{if}\; \beta = \alpha(i) \\ c_{k+1} = c_j \;\text{if}\; \beta' = \alpha(j)}} N[i,\, c_i;\, k,\, c_k] + N[k{+}1,\, c_{k+1};\, j,\, c_j] & \beta \neq \beta' \\[6pt] \displaystyle\max_{\substack{i \leq k < j,\;\; \beta = \beta' \\ c_k \in S(\beta) \\ c_k = c_i \;\text{if}\; \beta = \alpha(i) \\ c_k = c_j \;\text{if}\; \beta = \alpha(j)}} N[i,\, c_i;\, k,\, c_k] + N[k{+}1,\, c_k;\, j,\, c_j] - \text{cai}(\beta, c_k) & \beta = \beta' \end{cases}$$

The pair cases cover all four combinations of whether each pointer crosses a codon boundary when advancing inward; the inner recursion, $\lambda$ reward, and CAI credits are the same — only which codons must be re-chosen differs.

The codon-forcing constraints (e.g. $c' = c_j$ if $\alpha(i{+}1) = \alpha(j)$) arise when the advancing pointer lands inside the opposite boundary's codon, fixing that codon rather than allowing a free choice.

During bifurication if k and k+1 are pointing to the same amino acid, we subtract one $\text{cai}$ to prevent double-counting: the shared codon at the split appears as the right boundary of the left sub-problem and the left boundary of the right sub-problem simultaneously.

---

Optimal Solution

$$\text{score} = \max_{c_0 \in S(0),\; c_L \in S(m-1)} N[0,\, c_0;\, 3m{-}1,\, c_L]$$

---

Comparison

	CDSfold	DERNA	My Approach
Input	Protein sequence	Protein sequence	Protein sequence
Structural objective	Min MFE (Zuker)	Min MFE (Zuker)	Max base pairs
Optimizes CAI	✗	✔	✔
Pareto front	✗	✔	✔
DP tables	2	4	1
Complexity	$O(n^3)$	$O(m^3 L K^8)$	$O(m^3 K^4)$
Thermodynamically accurate	✔	✔	✗

Dynamic Programming Table Visualized

The green path traces the optimal codon choices and folding structure