# Dynamic Programming II

###### tags: `Dynamic Programming`

>Prerequisite: Recursion

## Content
[ToC]

## Longest Common Subsequence

### Problem Statement

Given two sequences X = (x1, x2,..., xm) and Y = (y1, y2,...yn), find the longest common subsequence (LCS) Z of X and Y, where a subsequence is a sequence that appears in the same relative order, but not necessarily contiguous.

### Example

Input:
    X = "ABCBDAB"
    Y = "BDCABA"
Output:
    Z = "BCBA"
    
Input:
    X = “AGGTAB”
    Y = “GXTXAYB”
Output:
    Z = "GTAB”

### Optimal Substructure

Let X = (x1,x2,...xm) and Y =(y1,y2,..., yn), Z = (z1, z2,...zk) is the LCS of X and Y.

Consider three conditions:
    1.If xm = yn, then zk = xm = yn and Zk-1 is LCS of Xm-1 and Yn-1
    2.If xm ≠ yn, then zk ≠ xm implies Zk is the LCS of Xm-1 and Yn.
    3.If xm ≠ yn, then zk ≠ yn implies Zk is the LCS of Xm and Yn-1.

In previous conditions 2 and 3, if xi ≠ yj, we will need to select a maximum (LCS of Xm-1 and Yn-1, LCS of Xm and Yn-1) . Let c[i,j] represent the LCS of subsequence of Xi and Yj. By the previous optimal substructure, we can write the recursion as follow:

c[i,j] = 0 if i = 0 or j = 0  
c[i,j] = c[i-1,j-1]+1 if xi = yj
c[i,j] = max(c[i-1,j], c[i,j-1]) if xi ≠ yj

### Pseudo Code

function LCS(X[], Y[]):
    m = X.lenght;
    n = Y.length;
    for i = 0 to m:
        c[i,0] = 0;
    for j = 0 to n:
        c[0,j] = 0;
    for i = 1 to m:
        for j = 1 to n:
            if xi == yj:
                c[i,j] = c[i-1,j-1]+1;
                b[i,j] = "left-up"
            else if c[i,j-1] ≥ c[i-1,j]:
                c[i,j] = c[i,j-1];
                b[i,j] = "left"
            else if c[i,j-1] ≤ c[i-1,j]:
                c[i,j] = c[i-1,j];
                b[i,j] = "up";

As usual, we use table b[i,j] to store the information in order to backtrack the entire LCS. Pseudo code for backtracking:

function Backtracking(i, j):
    if i == 0 || j == 0:
        return;
    else:
        if B[i][j] == "left-up":
            Backtrack(i, j);
            print(Xi);
        else if B[i][j] = "up":
            Backtracking(i-1, j);
        else:
            Backtracking(i, j-1);

### Complexity

Complexity of dynamic programming will be the size of table times the operation complexity to fill a single cell of the table. Therefore, the overall complexity will be $O(mn)$.

### C++ Code

#include <iostream>
#include <string>
#include <stdio.h>
using namespace std;

//C[i][j] = 0 when i=0;j=0 no pointer
//C[i][j]++ when X[i] == Y[j] point to \(g)
//C[i][j] = max{C[i-1][j], C[i][j-1]} point to <-(l) or |(u)

string s;
char B[100][100];
void backtrack(int m, int n, string x, string y){
    if(m == 0 || n == 0) return;
    else{
        if(B[m][n] == 'g'){
            backtrack(m-1, n-1, x, y);
            cout << x[m-1];
        }
        else if(B[m][n] == 'u'){
            backtrack(m-1, n, x, y);
        }
        else if(B[m][n] == 'l'){
            backtrack(m, n-1, x, y);
        }
    }
}

int LCS(string x, string y){
    unsigned long m = x.size();
    unsigned long n = y.size();
    int C[m+1][n+1];
    for(int i = 0; i <= m; i++){
        for(int j = 0; j <= n; j++){
            B[i][j] = ' ';
            C[i][j] = 0;
        }
    }
    for(int i = 0; i <= m; i++){
        C[i][0] = 0;
    }
    for(int i = 0; i<= n; i++){
        C[0][i] = 0;
    }
    for(int i = 1; i <= m; i++){
        for(int j = 1; j <= n; j++){
            if(x[i-1] == y[j-1]){
                C[i][j] = C[i-1][j-1]+1;
                B[i][j] = 'g';//left-up
            }
            else{
                if(C[i-1][j] >= C[i][j-1]){
                    C[i][j] = C[i-1][j];
                    B[i][j] = 'u';//up
                }
                else{
                    C[i][j] = C[i][j-1];
                    B[i][j] = 'l';//left
                }
            }
        }
    }
    backtrack(m, n, x, y);
    return C[m][n];
}

int main(int argc, const char * argv[]) {
    // insert code here...
    string x = "abcbdab";
    string y = "bdcaba";
    int result = LCS(x, y);
    cout << result << endl;
    return 0;
}

## Longest Common Substring

### Problem Statement

Given two sequences X = (x1, x2,..., xm) and Y = (y1, y2,...yn), find the longest common substring Z of X and Y, where a substring is a string that is contiguous (different from the subsequence).

### Example

Input:
    X = "abcdxyz"
    Y = "xyzabcd"
Output:
    Z = "abcd"
    
Input:
    X = “AGGTAB”
    Y = “GXTABYB”
Output:
    Z = "TAB”

### Optimal Substructure

The optimal substructure is similar to the previous problem (Longest Common Subsequence). Let c[i,j] represent the LCS of subsequence of Xi and Yj. If xi ≠ yj, we will not get any substring in this position. Therefore, the longest common substring will be 0. As the recursion below:

c[i,j] = 0 if i = 0 or j = 0  
c[i,j] = c[i-1,j-1]+1 if xi = yj
c[i,j] = 0 if xi ≠ yj

### Pseudo Code

function Longest-Common-Substring(X[], Y[]):
    m = x.length;
    n = y.length;
    max = inf;
    for i = 0 to m:
        for j = 0 to n:
            if i == 0 || j == 0:
                c[i,j] = 0;
            else if X[i] == Y[j]:
                c[i,j] = c[i-1, j-1]+1
                if max < c[i,j]:
                    pos = i
            else:
                c[i,j] = 0;
    print(X[max-pos:pos]);

Use max to record the current best answer and use pos to record the right-most position of the current longest common substring.

### Complexity

Complexity of dynamic programming will be the size of table times the operation complexity to fill a single cell of the table. Therefore, the overall complexity will be $O(mn)$.

### C++ Code

#include <iostream>
#include <stdio.h>
#include <string.h>
using namespace std;

string LCS(string x, string y){
    int max = 0;
    int max_i = 0;
    int DP[x.size()+1][y.size()+1];
    
    for(int i = 0; i <= x.size(); i++){
        for(int j = 0; j <= y.size(); j++){
            if(i == 0 || j == 0) DP[i][j] = 0;
            else if(x[i-1] == y[j-1]){
                DP[i][j] = DP[i-1][j-1]+1;
                if(max < DP[i][j]){
                    max_i = i;
                    max = DP[i][j];
                }
            }
            else DP[i][j] = 0;
        }
    }
    return x.substr(max_i-max, max);
}


int main(int argc, const char * argv[]) {
    string s = "AGGTAB";
    string s1 = "GXTABYB";
    cout << LCS(s, s1);
    return 0;
}

## Optimal Binary Search Tree

### Problem Statement

Given a sequence K = (k1, k2,...kn) of n distinct keys in sorted order. (ki < ki+1)

Build a binary search tree from these keys. For each key, we have probability pi that a search will be for ki.

Some searches may not be in sequence K, and we call these searches "dummy keys" d0,d1,... dn, which represents a value not in K.

In particular, d0 represents all values less than k1, dn represents all values greater than kn, the dummy key di represents all values between ki and ki+1. Also, for each key, we have probability qi that a search will be for di.

Each key is either successful (finding some key ki) or unsuccessful (finding some dummy key di).

Expected search time of a node ki is defined as the (depth(ki)+1)*pi. Our final goal is to minimize the total expected search cost of the tree, which is:

$E[cost] =\sum\limits_{i=1}^{n}(depth_T(k_i)+1)\cdot p_i + \sum\limits_{i=0}^{n}(depth_T(d_i)+1)\cdot q_i$

Note that $\sum\limits_{i=1}^{n}p_i + \sum\limits_{i=1}^{n}q_i = 1$.

Therefore, $E[cost] =\sum\limits_{i=1}^{n}(depth_T(k_i))\cdot p_i + \sum\limits_{i=0}^{n}(depth_T(d_i))\cdot q_i+1$

### Example

If we have the information as follow:

    i:      0      1      2      3      4      5
  p[]:          0.15   0.10   0.05   0.10   0.20
  q[]:   0.05   0.10   0.05   0.05   0.05   0.10
  
Two binary search tree as the following graph:

    Cost of graph 1:  2.80
    Cost of graph 2:  2.75 (better)


### Recursion
Consider that tree T is the optimal binary search tree, then the left subtree and the right subtree will also be the optimal substructure of tree T (skip the proof of optimal substructure).

Therefore, if we have a subsequence k' (k[i, j]) start from index i and end at index j of the original subsequence. We will need to find a key kr (i ≤ r ≤ j) such that the subtree will be the optimal substructure. Then we can write down the formula of the cost function as follow:

$E[i,j] =\sum\limits_{m=i}^{j}(depth_T(k_m)+1)\cdot p_m + \sum\limits_{m=i-1}^{j}(depth_T(d_m)+1)\cdot q_m$

Therefore, $E[i,j] =\sum\limits_{m=i}^{j}(depth_T(k_m))\cdot p_m + \sum\limits_{m=i-1}^{j}(depth_T(d_m))\cdot q_m+ \sum\limits_{m=i}^{j}p_m+\sum\limits_{m=i-1}^{j}q_m$

If we decide $k_r$ to be the root of the k[i,j], we will need to calculate the cost of left subtree and right subtree. The cost of left and right subtree are as follow:

$E[i,r-1] =\sum\limits_{m=i}^{r-1}(depth_{Tl}(k_m))\cdot p_m + \sum\limits_{m=i-1}^{r-1}(depth_{Tl}(d_m))\cdot q_m+ \sum\limits_{m=i}^{r-1}p_m+\sum\limits_{m=i-1}^{r-1}q_m$

$E[r+1,j] =\sum\limits_{m=r+1}^{j}(depth_{Tr}(k_m))\cdot p_m + \sum\limits_{m=r}^{j}(depth_{Tr}(d_m))\cdot q_m+ \sum\limits_{m=r+1}^{j}p_m+\sum\limits_{m=r}^{j}q_m$

Since $depth_T = depth_{Tr}+1 = depth_{T_l}+1$

Then, $E[i,j]=\sum\limits_{m=i}^{j}(depth_{Tl}(k_m))\cdot p_m + \sum\limits_{m=i-1}^{j}(depth_{Tl}(d_m))\cdot q_m+ 2\cdot \sum\limits_{m=i}^{j}p_m + 2\cdot \sum\limits_{m=i-1}^{j}q_m$

Finally, $E[i,j]-(E[i,r-1]+E[r+1,j]) =depth_{T_l}(k_r)\cdot p_r +2p_r+\sum\limits_{m=i}^{r-1}p_m + \sum\limits_{m=r+1}^{j}p_m+\sum\limits_{m=i-1}^{j}q_m$

Since $depth_{Tr}+1 = depth_{T}-1$, $E[i,j]-(E[i,r-1]+E[r+1,j]) = \sum\limits_{m=i}^{j}p_m+\sum\limits_{m=i-1}^{j}q_m$

Let $w[i,j] = \sum\limits_{m=i}^{j}p_m+\sum\limits_{m=i-1}^{j}q_m$

Therefore, we can write recursion relation as follow:

E[i,j] = E[i,r-1] + E[r+1, j] + w[i,j] 

Next, we will consider the initial condition of the recursion. Since i ≤ r ≤j , the boundary condition of the recursion will be r = i and r = j.

If r == i:
    E[i,j] = E[i, i-1] + E[i+1, j] + w[i,j]
Else if r== j:
    E[i,j] = E[i,j-1] + E[j+1, j] + w[i,j] 

E[i,i-1] might not be that intuitive; however, we will still have the dummy node as the left child of the tree. Therefore, E[i,j] will be defined as the probability q[i-1] of d[i-1].

On the other hand, E[j+1,j] will also be defined as the probability q[j] of d[j]. Note that the table should be (n+1)xn since we might need to calculate the value of E[n+1, n].

Therefore, the complete recursion relation will be as follow:

e[i,j] = q[i-1] if j = i-1;

e[i,j] = min{e[i,r-1]+e[r+1,j]+w[i,j]} for r = i, i+1,..., j.

### Pseudo Code

function Optimal-Binary-Search-Tree(p[],q[]):
    for i = 1 to n+1:
        e[i,i-1] = q[i-1];
        w[i,i-1] = q[i-1];
    for l = 1 to n:
        for i = 1 to n-l+1:
            j = i+1-1;
            e[i,j] = inf;
            w[i,j] = w[i,j-1]+p[j]+q[j];
            for r = i to j:
                t = e[i,r-1]+[r+1,j]+w[i,j];
                if t < e[i,j]:
                    e[i,j] = t;
                    root[i,j] = r;

The method of filling up the table is similar to that of the problem "Matrix-Chain Multiplication". Use root[i,j] to backtrack the structure of the tree.

### Complexity
Complexity of dynamic programming will be the size of table times the operation complexity to fill a single cell of the table. Therefore, the overall complexity will be $O(n^3)$.

Knuth offers a brilliant algorithm to find optimal binary search in $O(n^2)$!!!

### C++ Code

#include <iostream>
#include <string.h>
using namespace std;

double e[100][100];
double w[100][100];
int root[100][100];

void optimalBinarySearchTree(double p[], double q[], int N){
    for(int i = 1; i <= N+1; i++){
        for(int j = 1; j <= N; j++){
            e[i][j] = 0;
        }
    }
    for(int i = 1; i <= N+1; i++){
        e[i][i-1] = q[i-1];
        w[i][i-1] = q[i-1];
    }
    for(int l = 1; l <= N; l++){
        for(int i = 1; i <= N-l+1; i++){
            int j = l+i-1;
            e[i][j] = INT_MAX;
            w[i][j] = p[j] + q[j] + w[i][j-1];
            for(int r = i; r <= j; r++){
                double t = e[i][r-1] + e[r+1][j] + w[i][j];
                if(t < e[i][j]){
                    e[i][j] = t;
                    root[i][j] = r;
                }
            }
        }
    }
    cout << e[1][N] << endl;
}

int main(int argc, const char * argv[]) {
    // insert code here...
    double p[] = {0, 0.15, 0.1, 0.05, 0.1, 0.2};
    double q[] = {0.05, 0.1, 0.05, 0.05, 0.05, 0.1};
    int N = sizeof(p)/sizeof(p[1])-1;
    optimalBinarySearchTree(p, q, N);
    return 0;
}
//Output
//2.75