A Deep Learning Black Box Code

I mentioned in a previous blog that for the basic machine learning tasks such as classification or regression, we could replace the traditional machine learning tools, such as SVM or random forest, by a fully-connected deep learning solution.  If we do such tasks often, we should be able to write a black-box solution based on a few simplifications: (1) all layers have the same number of neurons; (2) dropout is applied only to the last layer; (3) We use ReLu.  I have not really work on this, but below is an outlier of the Keras code, which automatically explores the optimal NN solution for the MNIST dataset.

from keras.datasets import mnist
from keras import models, layers
from keras.utils import to_categorical
import numpy as np

(train_images, train_labels), (test_images, test_labels) = mnist.load_data()
train_images = train_images.reshape(-1, 28*28).astype('float32')/255
validate_images = train_images[-10000:]
train_images = train_images[:-10000]
test_images = test_images.reshape(-1, 28*28).astype('float32')/255
train_labels = to_categorical(train_labels)
validate_labels = train_labels[-10000:]
train_labels = train_labels[:-10000]
test_labels = to_categorical(test_labels)

def model(input_shape, n_classes, n_layer=3, nodes_in_layer=64, dropout=0.5):
    net = models.Sequential()
    net.add(layers.Dense(nodes_in_layer, input_shape=input_shape, activation='relu'))
    for i in range(n_layer-1):
        net.add(layers.Dense(nodes_in_layer, activation='relu'))
    if dropout>0:
            net.add(layers.Dropout(dropout))
    net.add(layers.Dense(n_classes, activation='softmax'))
    return net

input_shape=(28*28,)

n_classes=10
n_layers=[1, 2, 3, 4]
nodes=[16, 32, 64, 128]
dropouts=[0, 0.2, 0.5]

n_search=10
best_acc=0.0
for i in range(n_search):
    n_layer = n_layers[np.random.randint(4)]
    nodes_in_layer = nodes[np.random.randint(4)]
    dropout = dropouts[np.random.randint(3)]
    print("INFO> Test new model: n_layer=%d, nodes=%d, dropout=%.2f" % (n_layer, nodes_in_layer, dropout))
    net = model(input_shape, n_classes, n_layer, nodes_in_layer, dropout)
    net.compile(optimizer='adam', loss='categorical_crossentropy', metrics=['accuracy'])
    history = net.fit(train_images, train_labels, epochs=5, batch_size=64, validation_data=(validate_images, validate_labels))
    loss,acc = net.evaluate(test_images, test_labels)
    if acc > best_acc:
        print('INFO> Better model found, accuracy=', acc)
        best_acc=acc
        net.save('my_best_model.h5')
    else:
        print('INFO> Keep existing model, new model has worse accuracy=', acc)

An example run of the above script gives the following output, where the optimal model is highlighted:

INFO> Test new model: n_layer=1, nodes=16, dropout=0.50
INFO> Better model found, accuracy= 0.9193
INFO> Test new model: n_layer=3, nodes=32, dropout=0.50
INFO> Better model found, accuracy= 0.9577
INFO> Test new model: n_layer=4, nodes=64, dropout=0.50
INFO> Better model found, accuracy= 0.9707
INFO> Test new model: n_layer=2, nodes=64, dropout=0.50
INFO> Keep existing model, new model has worse accuracy= 0.9674
INFO> Test new model: n_layer=3, nodes=64, dropout=0.50
INFO> Better model found, accuracy= 0.9729
INFO> Test new model: n_layer=1, nodes=64, dropout=0.00
INFO> Keep existing model, new model has worse accuracy= 0.9674
INFO> Test new model: n_layer=2, nodes=32, dropout=0.00
INFO> Keep existing model, new model has worse accuracy= 0.9597
INFO> Test new model: n_layer=2, nodes=16, dropout=0.20
INFO> Keep existing model, new model has worse accuracy= 0.9383
INFO> Test new model: n_layer=4, nodes=64, dropout=0.20
INFO> Keep existing model, new model has worse accuracy= 0.9672
INFO> Test new model: n_layer=2, nodes=32, dropout=0.00
INFO> Keep existing model, new model has worse accuracy= 0.9583

Recurrent Neural Network

Introduction

For normal FNN, the elements in the input vector $\mathbf{x}$ has no intrinsic order.  Therefore, if one uses FNN to solve a computer vision problems, it will not be able to take advantage of the texture patterns formed by neighboring pixels, as the "neighboring" relationship is lost.  CNN instead takes advantage of the spatial relationship of a neighborhood and use that to form features, therefore, it is superior for computer vision problems.  Similarly, when the input is a sequence data, which could be a temporal time series or positional data such as a text stream, we do not want to scramble the order of the data points, instead should try to recognize sequence patterns, which requires a different architecture than FNN.  You may argue why not use CNN to recognize the sequence pattern.  CNN has two limitations, first, CNN can only look at patterns of a fixed size (determined by the size of its filters), but sometimes we are interested in long-range patterns.  For instance, to predict weather, no longer we have to look for patterns within the past few days, but few weeks, few months, even few years, such long-range pattern cannot be captured well by local CNN filters.  Second, CNN patterns are translational invariant, which means they are not context sensitive.  Patterns embedded in sequence data are often interpreted differently, depending on its context.  For example, the word "Chinese" within a document could mean either Chinese people or Chinese language, its exact meaning depends on the context.  CNN does not capture context.

We so far have been dealing with inputs of fixed size.  This is true for all the traditional machine learning methods we discussed previously (such as SVM, random forest), as well as for CNN.  Therefore to classify ImageNet pictures, we need to make sure they are all of the same dimensions.  In reality, many inputs are of variable length. For example, to transcribe voice, the voice wave form lasts differently depending on what the contents are and the rate of speech.  When XBox Kinect recognizes a kick action, the whole kick video sequence comes at varying number of frames, as there are quick kicks and slow kicks.  To read an email and then classify it into spam or ham, the length of the emails vary.  Traditionally for the XBox Kinect problem, we apply Bayesian graph theory to construct generative models to explain the video sequence, which is very complicated and heavily rely on oversimplified assumptions.  Therefore, we need a new NN architecture to handle variable input length, as well as mechanisms to discover patterns of varying time spans.

When input elements are fed into our machine one piece at a time, the machine must have memory mechanism.  If the machine makes a decision simply based on the current input token, word "Chinese" will always be associated with either a language or people, preventing the true meaning from being assigned based on its context.  In translation, the output sequence should also be variable.  The new architecture invented to overcome the challenge of variable-length input/output sequence is called recurrent neural network (RNN).  In this blog, we try to use simple examples to convince you RNN theoretically is capable of capturing patterns in a sequence of arbitrary length.  RNN has applications in language translation, video classification, photo captioning, trend prediction (arguably stock price prediction, not sure if that can work in theory.  See Table 4 and conclusion in [1]), etc.

Architecture

The basic unit of RNN is a cell as depicted in Figure 1.  A cell has an internal state, which is a vector $\mathbf{h}$, so a cell can store multiple internal elements.  At time $t$, it reads in the input $\mathbf{x}_{t}$, mixed with its current hidden status $\mathbf{h}_{t-1}$ to transition into a new state $\mathbf{h}_{t}$; $\mathbf{h}_{t}$ can then be used to produce an output vector $\mathbf{y}_{t}$.  Imagine $\mathbf{h}_{0}$ is initialized as a zero vector, the cell is able to read in input $\mathbf{x}$, one at a time, and produces a sequence of hidden states $[\mathbf{h}_1,\mathbf{h}_2, \ldots, \mathbf{h}_n]$, for our purpose, we can assume $\mathbf{y} = \mathbf{h}$ for now.

Figure 1. An RNN cell rolled out in time sequence [2].

The update math can be described as the following:

$$\begin{equation}
\mathbf{h}_{t} = \phi(\mathbf{h}_{t-1}\mathbf{W}_{h}+ \mathbf{x}_{t}\mathbf{W}_x + \mathbf{b}),
\end{equation}$$

where $\mathbf{W}_x$, $\mathbf{W}_h$ and $\mathbf{b}$ are cell parameters to be determined by optimization.  $\mathbf{h}_{t}$ is the hidden state of the cell, i.e., it is the impression of all the data the RNN cell has viewed so far: $[\mathbf{x}_{1}, \mathbf{x}_{2}, \ldots, \mathbf{x}_{t}]$.  After we watched a movie, the snapshot of our brain is $\mathbf{h}_t$, which contains all the effects of the whole movie.  This vector is then used as the input to classify the movie into a good one or a bad one.

Now let us try to understand why this simple recurrent architecture works.  I have not read about any "proofs" as elegant as what we did for FNN in a previous blog, where it is rather convincing to see why FNN can model any arbitrary function.  For RNN, the best I can do is to present some special cases to convince you RNN is indeed quite capable.

Memory

RNN cells can memorize data.  Imagine the input is the time series of daily stock price, and output $y$ is the most-recent three day average, or any prediction of tomorrow's price based on the past three days.  In the prediction at time $t$, we only receive a one-day input $\mathbf{x}_{t}$, thus there needs to be a way for CNN to retain $\mathbf{x}_{t-1}$ and $\mathbf{x}_{t-2}$, i.e., it should have short-term memory capability in order to use three pieces of data to generate the output $y_t$.

Let us look at the special case, where stock price $x_t$ is a scalar:
$$\begin{equation} \mathbf{h}_0 = {[ 0, \; 0, \; 0]}^\intercal, \\ \mathbf{h}_t = \phi(\mathbf{h}_{t-1}\mathbf{W}_{h}+x_{t}\mathbf{W}_x +\mathbf{b}), \\ \phi(x) = x, \\ \mathbf{W}_h = \begin{bmatrix} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0 \end{bmatrix}, \\ \mathbf{W}_x = {[ 0, \; 0, \; 1 ]}^\intercal, \\ \mathbf{b} = {[ 0, \; 0, \; 0 ]}^\intercal. \\ \end{equation}$$

Our hidden state contains three elements, reserved to store the stock prices for the past three days.  We follow the time sequence, the following $\mathbf{h}$ sequence can be obtained based on Equation 2:

$$\begin{equation} \mathbf{h} = \left\{ \mathbf{h}_0=\begin{bmatrix} 0 \\ 0 \\ 0 \end{bmatrix}, \mathbf{h}_1=\begin{bmatrix} 0 \\ 0 \\ x_1 \end{bmatrix}, \mathbf{h}_2=\begin{bmatrix} 0 \\ x_1 \\ x_2 \end{bmatrix}, \mathbf{h}_3=\begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix}, \\ \mathbf{h}_4=\begin{bmatrix} x_2 \\ x_3 \\ x_4 \end{bmatrix}, \ldots, \mathbf{h}_t=\begin{bmatrix} x_{t-2} \\ x_{t-1} \\ x_t \end{bmatrix}, \ldots \right\}. \end{equation}$$

This way RNN can retain a short-term memory.  Although at each time point, we could use all stock prices for the past three days as the feature vector input in real applications, we here choose a scalar input form to illustrate the memory mechanism.

Running Average and Derivative

Let us consider another special case, where our hidden state has three elements, reserved to store the running average $h_{\Sigma}$, running derivative $h_{\nabla}$, and memorize the previous input $h_{mem}$, respectively:

$$\begin{equation} \mathbf{h}= {[ h_{\Sigma}, \; h_{\nabla}, \; h_{mem} ]}^\intercal, \\
\mathbf{h}_0 = {[ 0, \; 0, \; 0]}^\intercal, \\ \mathbf{h}_t = \phi(\mathbf{h}_{t-1}\mathbf{W}_{h}+x_{t}\mathbf{W}_x +\mathbf{b}), \\ \phi(x) = x, \\ \mathbf{W}_h = \begin{bmatrix} \beta_1 & 0 & 0 \\ 0 & \beta_2 & -(1-\beta_2) \\ 0 & 0 & 0 \end{bmatrix}, \\ \mathbf{W}_x = {[ 1-\beta_1, \; 1-\beta_2, \; 1 ]}^\intercal, \\ \mathbf{b} = {[ 0, \; 0, \; 0 ]}^\intercal. \\ \end{equation}$$

The recurrent equation can be simplified into:

$$\begin{equation}
\begin{bmatrix} h_{\Sigma,t} \\ h_{\nabla,t} \\ h_{mem,t} \end{bmatrix} =
\begin{bmatrix} \beta_1 h_{\Sigma,t-1} + (1-\beta_1)x_t \\
\beta_2 h_{\nabla,t-1} + (1-\beta_2) (x_t - x_{t-1}) \\ x_t \end{bmatrix}.
\end{equation}$$

We follow the time sequence, we obtain the following $\mathbf{h}$ sequence:

$$\begin{equation} \mathbf{h} = \left\{ \mathbf{h}_0=\begin{bmatrix} 0 \\ 0 \\ 0 \end{bmatrix}, \mathbf{h}_1=\begin{bmatrix} (1-\beta_1)x_1 \\ (1-\beta_2)x_1 \\ x_1 \end{bmatrix}, \\
\mathbf{h}_2=\begin{bmatrix} (1-\beta_1)\left(\beta_1  x_1+  x_2 \right) \\ (1-\beta_2)\left( \beta_2(x_1-0)+ (x_2-x_1)\right) \\ x_2 \end{bmatrix}, \ldots, \\
\mathbf{h}_t=\begin{bmatrix} (1-\beta_1)\left(\beta_1^{t-1}x_1+\beta_1^{t-2}x_2+\cdots + x_t)\right) \\  (1-\beta_2)\left(\beta_2^{t-1}(x_1-0)+\beta_2^{t-2}(x_2-x_1)+ \cdots + (x_t-x_{t-1})\right)  \\ x_t \end{bmatrix}\right\}. \end{equation}$$
Now we can see how the three elements in our hidden state work.  $h_{mem}$ simply memorizes the current input, so that $x_{t-1}$ is available to be used with $x_t$ at time point $t$.  $h_{\Sigma}$ is for calculating a running average, which is similar to the momentum optimization.  The current average is first decayed by a factor $\beta_1$, the new input is weighted by $(1-\beta_1)$ and used to update the running average.  $h_{\nabla}$ is for calculating a running first-order derivative.  The current derivative is first decayed by a factor $\beta_2$, the new derivative estimated based on $(x_{t}-x_{t-1})$ is weighted by $(1-\beta_2)$ and used to update the derivate memory.   We here assume $\beta_1$ and $\beta_2$ are decay rate within $[0, 1]$.  When $\beta$ approaches one, we see data far into the past, the new input has little weight to change the average or derivate too much.  When $\beta$ approaches one, we only see the effect of very recent inputs.  Imagine we create many RNN cells with different $beta$ values, we can have a range of cells that capture the temporal average or derivative of varying time windows.  If you stack such memory cells, say one $h_{\nabla}^{(2)}$ on top of the output of another $h_{\nabla}^{(1)}$, we can get higher-order derivatives with a multi-layer RNN architecture (Deep RNN).

Above thought experiments are far from any proof, they nevertheless offer us some confidence that with enough cells and enough layers, an RNN is able to capture many temporal patterns.  In fact, it was believed that RNN can be used to simulate arbitrary computer program [3]:

Thus RNNs can be thought of as essentially describing computer programs.  In fact, it has been shown that RNNs are turing complete (for more information refer to the article: On the Computational Power of Neural Nets, by H. T. Siegelmann and E. D. Sontag, proceeding of the fifth annual workshop on computational learning theory, ACM, 1992.) in the sense that given the proper weights, they can simulate arbitrary programs.

When the activation of its hidden state can encode different temporal patterns (just as feature neuron in CNN), they can certainly be used as powerful features to for downstream predictions.  If our thought experiment is representative, the memory exercise tells us we would need an RNN cell to contain at least $m$ hidden elements, if we want to predict the future based on the past $m$ time points.  Therefore, RNN cells should have reasonable complexity to be useful.  The second exercise tells us activities captured by a cell contains an exponential term on $\beta^t$.  During optimization, our $\beta$ may become greater than one, which will result in an explosion and make optimization impossible, as early steps make infinite contribution.  When $\beta$ becomes too small, the early steps have virtually no contribution, this is called vanishing gradient problem, which prevents RNN from looking far into the past.  Therefore this RNN architecture is actually very hard to train.   Also, it is inefficient if we would need large number of memory elements, when there is a need to use stock price in the past few years to predict tomorrow's price.  In practice this classical RNN architecture is now called SimpleRNN, the architecture we discussed above has been replaced by more efficient alternatives, such as LSTM (long short term memory) to be introduced later.

RNN Applications

RNN has a few variants that make it suitable to handle a number of very interesting use cases.

Case A, each input produces an output.  Say stock price prediction, each new input $x_t$ of today's price will enable us to combine it with the data it has seen so far (coded in $\mathbf{h}_{t-1}$) to predict tomorrow's price.  We can conveniently ignore the first few outputs as they were based on very few inputs, and start using the prediction only after it has seen enough historical data points.

Case B, all outputs are ignored, only the last output is used.  An RNN reads a movie review from beginning to the end, its last hidden state $\mathbf{h}_{last}$ encodes content for the whole review, therefore, this can be used as the feature representation of the review and be passed onto a FNN to classify the review into positive/negative or to be regressed into movie ratings.

Case C, only one input and the RNN generates a sequence of output.  Give the RNN a photo $\mathbf{x}_0$, RNN represents the photo into a hidden state, then this hidden state can recurrently output a sequence of words to be used as the capture.  This sounds magic, but it actually has been shown to work reasonably well, please watch the nice TED talk by Feifei Li [4], if you have not.

Case D, RNN takes an input sequence, then it outputs a new sequence.  This is what behind language translocation, where an English sentence is sent to an RNN, its hidden state at the end encodes the meaning of the sentence, it then spits out a new sentence in Chinese.  Case D is actually the combination of Case B and Case C.

Figure 2. Four different RNN use cases [5].

LSTM

We mentioned above that the SimpleRNN architecture has many practical issues, especially it cannot accommodate too many time steps.  When there are many time steps, the network suffers from either vanishing or exploding gradient problem, so RNN cannot have longer term memory.  LSTM is a technique that addresses these challenges by enabling longer short-term memories, so it is called Long Short-Term Memory cells.

Figure 3. Anatomy of an LSTM cell [6].

Figure 4. The update rules of LSTM [7].

LSTM cell is still quite a magic to me, so my confidence is not very high for what I write next.

The example I am using is inspired by a tool called LSTMVis [8].  Let us consider an LSTM cell to parse the following string, where numeric tokens are associated with different depth defined by parentheses (highlighted in three colors, one per depth), and the ultimate goal of the LSTM is to sum up all numbers associated with depth two and beyond (blue and green numbers).

( 4.5 ( 2.5 ( 3.5 5.5 ) 3.2 ( 1.0 ) ) 2.8 ( 6.8 9.0 ) )

The encounter of each left parenthesis increases the depth by one and a right parenthesis decreases the depth by one.  Our cell vector $\mathbf{c}$ naturally contains two elements, one for counting the depth and another for storing the sum, i.e., ${[ c_{depth}, c_{sum} ]}^\intercal$.  Our input vector $\mathbf{x}$ has four elements: ${[ x_l, x_r, x_n, x_v ]}^\intercal$, where $x_1$, $x_r$ and $x_n$ are one-hot encoder indicating the input is a left parenthesis, a right parenthesis, or a number, respectively.  When it is a number, $x_v$ stores the numerical input value.  Therefore string "(2.5 4.2)" is represented as four inputs:

$$\begin{equation} \begin{bmatrix}1\\0\\0\\0\end{bmatrix},  \begin{bmatrix}0\\0\\1\\2.5\end{bmatrix}, \begin{bmatrix}0\\0\\1\\4.2\end{bmatrix}, \begin{bmatrix}0\\1\\0\\0\end{bmatrix}. \end{equation}$$

Here we only look at the special case where output gate is always open, i.e., $\mathbf{o} = {[1, 1]}^\intercal$, and output state $\mathbf{h} = \mathbf{c}$.

Let us first consider the case where there is no feedback loop, i.e., $\mathbf{h}_{t-1}$ is not linked to the input at time $t$.  The input gate is determined by:
$$\begin{equation} \mathbf{i} = \rm{sign} \left(\begin{bmatrix} 1 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \end{bmatrix}\cdot \mathbf{x}\right). \end{equation}$$
Apply Equation 8 to Equation 7, we obtain a sequence of input gates:

$$\begin{equation} \begin{bmatrix}1\\0\end{bmatrix},  \begin{bmatrix}0\\1\end{bmatrix}, \begin{bmatrix}0\\1\end{bmatrix}, \begin{bmatrix}1\\0\end{bmatrix}. \end{equation}$$

The input gate is ${[1, 0]}^\intercal$, when it sees either "(" or ")", therefore, the input gate is open for $c_{depth}$ in order for it to be updated and the input gate for $c_{sum}$ is closed as no update on that element is needed.  We here replace $\rm{tahn}(\cdot)$ by a simple $\rm{sign}(\cdot)$ function :
$$\begin{equation} \rm{sign}(x) = \begin{cases} 1, \; \rm{if} \; x \gt 0 \\ 0,  \; \rm{if} \;x \le 0 \end{cases} \end{equation}$$
On the other hand, when it seems a number, the input gate will be ${[0, 1]}^\intercal$, which prepares $c_{sum}$ for an update and protects $c_{depth}$ from being modified.

The forget gate can be kept open all the time, i.e., ${[1, 1]}^\intercal$ for our simply example.

The new data is calculated by:
$$\begin{equation} \mathbf{g} =\rm{sign} \left( \begin{bmatrix} 1 & -1 & 0 & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix}\cdot \mathbf{x} \right). \end{equation}$$

Apply Equation 11 to Equation 7, we get the sequence:

$$\begin{equation} \begin{bmatrix}1\\0\end{bmatrix},  \begin{bmatrix}0\\2.5\end{bmatrix}, \begin{bmatrix}0\\4.2\end{bmatrix}, \begin{bmatrix}-1\\0\end{bmatrix}. \end{equation}$$

This means $\mathbf{g}$ is ${[1, 0]}^\intercal$ for "(", ${[-1, 0]}^\intercal$ for ")", and ${[0, x_v]}^\intercal$ for number inputs.

Therefore, the final update rules for $\mathbf{c}$ is the following (recall that forget gate simply just pass $\mathbf{c}_{t-1}$ through according to Figure 4:

$$\begin{equation} \begin{aligned}
c_{depth,t} &= \begin{cases} c_{depth,t-1} + 1, \rm{if \; input \; is \; (} \\
c_{depth,t-1} - 1, \rm{if \; input \; is \; )} \\
c_{depth,t-1}, \rm{if \; input \; is \; numeric}\end{cases}, \\[2ex]
c_{sum, t} &= \begin{cases} c_{sum, t-1}, \rm{if \; input \; is \; ( \; or )} \\
c_{sum, t-1}+x_v, \rm{if \; input \; is \; numeric} \end{cases}. \end{aligned} \end{equation}$$

So far, we are summing up all numbers read from the input, regardless of their depths, but our ultimate goal is to only sum the numbers when its depth is greater than one.  Since now the input gate for $c_{sum}$ should only be opened, when $c_{depth} \gt 1$, we link our output $\mathbf{h}$ back to the next time step to influence our input gate calculation.  Our new input gate calculation formula is updated from Equation 8 into the following:

$$\begin{equation} \mathbf{i} = \rm{sign} \left(\begin{bmatrix} 1 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \end{bmatrix}\cdot \mathbf{x} + \begin{bmatrix} 0 & 0 \\ 1 & 0 \end{bmatrix}\begin{bmatrix}c_{depth} \\ c_{sum}\end{bmatrix} - \begin{bmatrix} 0 \\ 2\end{bmatrix} \right). \end{equation}$$

The input gate for $c_{depth}$ behaves the same as before, as the two added terms have no effect on the top output element.  When the input is a number, the input gate for $c_{sum}$ is controlled by the value of $c_{depth}-1$, i.e., it only opens for update when $c_{depth} \gt 1$, exactly what we want.  The linkage of output to the next input brings the context (memory), so that the same input $\mathbf{x}_t$ produces different gate values and our cell is updated in a context-sensitive manner.

With the capability of longer-term memory, LSTM network can be used to deal with very long sequences.  Why LSTM can avoid the vanishing/exploding gradient problem that troubles the SimpleRNN?  A good explanation is offered by a blog here [9].  The main nested product involved in back-propagation is $\frac{\nabla c_{t}}{\nabla c_{t-1}}$.  Based on Figure 4, the addition term containing $\mathbf{i}$ and $\mathbf{g}$ has no contribution, the derivative is basically $\mathbf{f}$ - the forget gate.  $\mathbf{f}$ value is never greater than one, therefore, there is no exploding problem.  When memory needs to be retained, $\mathbf{f}$ values should be close to one, so the vanishing gradient problem is also diminished.  The gradient is terminated when old memory needs to be wiped out ($\mathbf{f}$ is set to zero).

Through the few virtual experiments discussed here, hopefully we now have more faith in the recurrent network.  As LSTM is so successful, when people say they use RNN, they actually refer to LSTM by default.  The traditional RNN is called SimpleRNN instead.  Another RNN variant call Gated Recurrent Unit (GRU) has shown to be simpler than LSTM and perform equally well, but we do not dig into that for now.  For readers speaking Chinese, a nice video on RNN is available here [10].

I am very excited by the RNN architecture, because it seems to explain how we dreams.  In the video here (20:00 - 31:00), it was shown that we can feed an RNN with Shakespeare, with C-code, or with mathematical articles and the RNN seems to be able to spit out shakespeare, C-code, and math articles.  The RNN architecture is outline in Figure 5.  Where one input character triggers the output of the next character, then the next character, etc.  Starting with a character, the RNN generates a sequence of characters.  If we feed the input with Shakespeare, the output generated by the RNN is then compared to the desirable output and the cross-entropy loss function can be used to train the weights in the RNN.  With enough iterations, the network starts to generate text that share the many syntax features on the training sequence.

Figure 5. The input is one character at a time. The character is one-hot encoded.  The hidden state is updated accordingly and then used to generate logits for the output characters (softmaxed into probabilities).  One output character spits out as the result.  This is a sequence-to-sequence application depicted in Figure 2.A. (credit: [11])

Figure 6 & 7 shows some sample outputs where RNN generate Shakespeare and Linux C code.

Figure 6.  With sufficient training samples, the RNN can generate something syntactically looks like Shakespeare. 

Figure 7.  With sufficient training samples, the RNN can generate something syntactically looks like C programming code. (Credit [11])

The generated output seems sensible at their local structures, but non-sensible at the longer range.  This conveniently explains how we dream, fragments of dreams seem quite logical and movie-like, however, these fragments are loosely connected, they do not serve the purpose of one theme.  When the fragments accidentally make some sense, you are extremely impressed by the power of dreaming.  It is like how people chat, moving from one topic to the next, logical within each transition, but accomplishing nothing as the whole.

When the elements within the hidden states were examined, some elements appear to represent interesting sequence pattern, similar to the depth and summation elements in our previous example.  Figure 8 shows a few examples of such neurons, this reminds us that RNN is able to train meaningful neurons just like CNN neurons can represent complex patterns, such as faces or texts.

Figure 8. Three example elements in the hidden state, where each captures a distinct sequence pattern.  A. A neuron that is activated when it encounters a left quote, then deactivates when right quote is encountered, very similar to our previous $C_depth$ neuron. B. A neuron appears to count the number of characters it reads and gets reset when line break is encountered.  C. A neuron that counts the indention depth of the C code.  (Credit: [11])

Reference

1. http://lup.lub.lu.se/luur/download?func=downloadFile&recordOId=8911069&fileOId=8911070
2. Aurelien Geron, Hands-on machine learning with scikit-learn & tensorflow. p373.
3. Antonio Gulli, Deep Learning with Keras: Implementing deep learning models and neural networks with the power of Python, page 214.
4. https://www.ted.com/talks/fei_fei_li_how_we_re_teaching_computers_to_understand_pictures
5. Aurelien Geron, Hands-on machine learning with scikit-learn & tensorflow. p385.
6. Aurelien Geron, Hands-on machine learning with scikit-learn & tensorflow. p403.
7. Aurelien Geron, Hands-on machine learning with scikit-learn & tensorflow. p405.
8. http://lstm.seas.harvard.edu/
9. http://harinisuresh.com/2016/10/09/lstms/
10. https://www.youtube.com/watch?v=xCGidAeyS4M
11. http://cs231n.stanford.edu/slides/2016/winter1516_lecture10.pdf

Motif Discovery for DNA- and RNA-binding Proteins by Deep Learning

Bioinformatics Application of Deep Learning Technology

Resource

We will discuss the following paper: Alipanahi et al., Predicting the sequence specificities of DNA- and RNA-binding proteins by deep learning.  The paper is available here, the supplementary method is available here (subscription required).

All figures shown here should be credited to the paper, unless the source is explicitly provided.

Introduction

Genes are encoded by DNA sequences.  Whether a gene is expressed and how much it is expressed can be regulated by proteins called transcription factors (TF).  A TF protein can bind to certain DNA motifs, then recruit other factors to initiate the transcription of a piece of DNA sequence into an mRNA sequence.  In eukaryotic cells, mRNA consists of multiple exons separated by introns.  Introns are spliced out and exons are concatenated into a matured mRNA, which will then be translated into a protein.  RNA-binding proteins can regulate the splicing process and lead to the inclusion or exclusion of certain exons, therefore, producing different splicing isoforms, eventually different protein products of different functions.  Mutations in the genomic sequence can perturb the binding of proteins, then lead to undesirable under/over transcription of genes, or undesirable splicing forms, which can induce disease.  Therefore, it is of interest to understand what nucleotide pattern (motifs) a DNA/RNA-binding protein can recognize and how mutations impact such molecular recognition.

Let us describe two popular experimental approaches for identifying sequences enriched for protein binding sites.  First, protein binding microarrays (PBM).  The microarray contains many different short nucleotide sequences called probes at different array locations.  For example, an array can contains probes representing many variants of 10-nucleotide patterns (10-mers).  Even for a given 10-mer, it can include many probes where the same 10-mer is embedded within different sequence context of the longer probe sequence.  A TF protein is hybridized over the array and it will be stuck to the places where binding occurs and produce fluorescent readouts.  Therefore, PBM gives readout of probe sequences and their intensity readout indicates the strength of the binding to a specific TF.  A deep learning (DL) system aims to take these probe sequences as input $X$ and predict the binding score $y$ (i.e., intensity).  The natural choice of the loss function would be least-square error, as this is a regression task.  Second, ChIP-seq technology.  When a TF is cross-linked to genomics materials, DNA were chopped and all DNA fragments bound to the TF will be extracted through chromatin immunoprecipitation (ChIP) process, bound DNA/RNA fragments are NGS-sequenced to produce readout (NGS stands for next-generation sequencing).  ChIP-seq data only provides positive sequence fragments, one needs to added negative sequences (either taken from genome background, or the authors here generated negative sequence by permuting positive sequences), a DL system aims to classify an input sequence $X$ into binary $y$ of 1/0, i.e., bind or not bind.  Since this is a classification problem, cross-entropy loss function is a natural choice.  The paper also studied datasets generated using a HT-SELEX technology, where the data format is similar to that of ChIP-seq.

The motif discovery problem is a well-studied bioinformatics problem for the past twenty years (I also contributed an algorithm called GEMS).  The positive sequences are collected and algorithms look for patterns that are over-represented among them, compared to a pool of negative background sequences.  Alipanahi et al. used a DL approach, more specifically, a CNN + FNN architecture that archived better performance than existing tools.  It is a very clean application of CNN to bioinformatics.

Results

Architecture

In the diagram below, the system processes five input sequences in parallel.  This is rather a technical detail in order to make full use of the GPU resource.  For our purpose, we only need to focus on one sequence.  The sequence is first cast into two sequences, the original one and its reverse complement (e.g., ATGG is cast into CCAT as well).  This is because a TF could theoretically recognize either one of the DNA strand, or both in rare cases (Figure 4.e in the paper is such a rare case).  A motif is a sequence pattern represented by a 1-D CNN convolution filter.  More precisely, the input sequence is actually a 4-channel 1-D vector, where 4 channels encode the base composition A, T, C, G and the vector is a 1/0 binary vector one-hot encoding the input sequence.  For example, a short sequence "ATGG" can be represented as:


$$S_{ATGG} = \begin{bmatrix} A & C & G & T\\
. 25  & . 25 & . 25 & . 25\\
. 25 & . 25 & . 25 & . 25 \\
1 & 0 & 0 & 0 \\
0 & 0 & 0 & 1 \\
0 & 0 & 1 & 0 \\
0 & 0 & 1 & 0 \\
. 25 & . 25 & . 25 & . 25 \\
. 25 & . 25 & . 25 & . 25 \end{bmatrix}. $$

Those rows containing 0.25 are paddings (represent "any" nucleotide) to extend the original sequence, so that convolution will not miss motifs partially overlapping with the ends.  If the CNN is designed to capture motifs of length $m$ (it will capture motifs length $\lt m$ as well), there will be $m-1$ extra 0.25- rows on both ends.  A motif is represented as a convolution filter $M$ that is known as position-weighted matrix (PWM) in bioinformatics, it basically specifies the weighting of each nucleotide at each location.  The output of the convolution layer is:

$$\begin{equation} \begin{aligned}
 X_{i,k} &= \sum_{j=1}^{m} \sum_{l=1}^{4} S_{i+j,l}M_{k,j,l}, \\
Y_{i,k} &= \max(0, X_{i,k}-b_k),
\end{aligned} \end{equation}$$

where $k$ is the index of a particular motif (a motif here is called a feature in CNN), $i$ is a location on the input sequence, $m$ is the length of the motif or the convolution window, the convolution therefore is taken over all positions $j$ and over all four nucleotide channels $l$.  The sum is then passed through a ReLU function as shown in the second step (weights can be negative, so neurons with convolution sum less than $b_k$ will not be activated).  The result is a feature map representing where a given motif is detected along the input.  A pooling layer is followed.  The authors used a max-pooling layer for DNA-binding proteins and both max- and average-pooling for RNA-binding proteins.  In the latter case, $k$ features results in $2k$ pooling maps.  The pooled output will then serve as the input for a fully-connected neural network (FNN) to generate target scores, one per input sequence.

This architecture is very simple, as it only contains one CNN convolution layer and a pooling layer.  The usage of CNN here is because the motif can be located at arbitrary positions in the input sequence, just like a dog can appear in any location of an image, CNN, beings translation-invariant, is a natural option for this problem.

Performance

The CNN is then trained with a set of experimental data (either PBM, ChIP-seq, or HT-SELEX) and its performance was shown to be consistently superior than established bioinformatics motif discovery tools.  Without getting into the details of biology applications, the figure below shows their CNN solution, named DeepBind,was able to out-perform other methods on PBM, ChiP-seq, SELEX datasets according to various performance metrics (AUC or Pearson correlation scores).

Why CNN is Better?

Unfortunately, the paper did not discuss why DeepBind performed better.  Although it shows many examples of how the motifs identified makes biological sense, but that is not something unique, as traditional motif discovery tools can produce those results as well (the author did not say traditional tools failed in doing that).  It was only mentioned in the beginning of the paper:

"The sequence specificities of a protein are most commonly characterized using position weight matrices (PWMs), which are easy to interpret and can be scanned over a genomic sequence to detect potential binding sites. However, growing evidence indicates that sequence specificities can be more accurately captured by more complex techniques".

This hints the motivation of their study is to describe motifs with a language more sophisticated than PWM.  However, motifs in CNN are described by filters, which is basically the same as PWM, what is the strength of this approach then?

I think what happened probably is DeepBind contains $k$ motif sub-elements.  In previous methods, only one PWM is used to explain biological measurements (for that specific transcription factor or RNA-binding protein).  However, in DeepBind, multiple elements are combined through a FNN for activity prediction.  Although each motif element is equivalent to a PWM, their combination could non-linear motif descriptor cannot be captured by PWM.  Hypothetically, imagine there are two motif elements in a real system: AGGGT and TGGGA.  That is "GGG" in the middle, when "A" appears on the left, "T" is favored on the right, similarity a left "T" pairs with a right "A".  However, AGGGA and TGGGT are not valid patterns.  However, if we only allow this to be represented using one PWM, the result is an "averaged" motif looks like "A/T-GGG-A/T", which would encode wron patterns such as AGGGA and TGGGT.  We therefore hypothesize the CNN system provides a way to capture motifs in a non-linear way, more powerful than PWM.


Another possible reason is NN allows a motif consisting of wide-spacing sub elements, which are very hard to be discovered by traditional algorithms due to the vast search space.  But wide-spacing motifs are important particularly for RNA-binding proteins, as they bind to RNA structures in non-trivial manner.  As shown below, unlike DNA-binding, RNA-binding proteins may not recognize only one binding site, but may actually recognize multiple binding sites on the surface of RNA.  Unlike double-strand DNA fragments (the last structure in the figure below), single-strand RNA can form complex secondary structures, oftentimes, the RNA structure contributes significantly to the binding no less than the specific oligonucleotide bases.  The authors commented: "RNA recognition motif, the most abundant RNA binding domain, is highly flexible, and usually the assembly of multiple domains are needed for proper binding. Third, “indirect readout” in which the structure of nucleotides guide the binding process, is more prominent in RBPs; Gupta and Gribskov analyzed 211 protein-RNA structural complexes and reported that nearly three-quarters of protein-RNA interactions involve the RNA backbone, not just the bases."  A RNA-binding protein could recognize multiple binding sites on the surface in contact, therefore, this is better modeled by multiple motif features in DeepBind, where they are combined through FNN to predict the binding affinity.  Traditional tools identify motifs as one continuous binding element, without synergistic combinations of underlying elements.  I suspect this might be the main reason DeepBind performs well, at least for the more difficult RNA-binding proteins.

Protein-binding interface on RNA (DNA, the last example) (source).  Purple: binding sites, yellow: non-binding sites.

We here only see one convolution layer, instead of multiple repeats of convolution plus pooling layers.  Hierarchical motif patterns may not be necessary for motif discovery problems, because although shorted motif elements, say "ATC" might appear in multiple longer motifs, the frequency of A, T, and C in the longer motif might be very different compared to their frequencies in the short elements.  So unlikely eyes, noses and mouths made up of faces, short motif PWMs may not be very effective in concatenated to form larger PWMs.

How to Visualize CNN?

A PWM generated by traditional bioinformatics tools can be easily visualized as a motif logo, where the height of a character encodes the frequency or information content of a base in its PWM.  CNN here contains multiple motifs and combines them in a non-linear way, therefore, it is conceptually incorrect to visualize each convolution filter alone.  How to visualize DeepBind is a rather challenging topics, the authors did not address this in my view.  Instead they adopted a much naive approach -- basically visualizing individual feature.  For each sequence in the test set, if the sequence contains a positive output signal for a given motif feature map, the corresponding input sequence fragment (gives the maximum activation) is retained.  Then all such CNN activating fragments are aligned and base counts are compiled to construct a PWM and subsequently generating a motif logo.  Some example logos are shown in the figure below, where it shows DeepBind model indeed reproduced many known motifs.  As we mentioned the downside of this approach is it really does not visualize the non-linearity introduced in the FNN component as we hypothesized, why DeepBind performs better is really not answered in this simplified representation.

Examples of known motifs discovered by DeepBind.

Mutation Map

An alternative visualization that is conceptually more powerful is called mutation map.  Given an input sequence $s$, DeepBind predicts a binding score $p(s)$.  We can then mutation the nucleotide at each position into a different bases, i.e., resulting in a new sequence $\hat{s}$ and predict a new binding score $p(\hat{s})$.  The importance of a particular single-base mutation, therefore, can be quantified by:

$$\begin{equation} \Delta S_{ij} = p(\hat{s})-p(s). \end{equation}$$

The nice thing is the mutation score is based on the whole CNN plus FNN system, therefore, include all the non-linear effects.  However, I guess the above equation did not work well as it is.  As the predicted scores do not necessarily represents the change in the unit of binding energy; The scores and binding energy might be related through some non-linear function.  The authors empirically used the following formula instead, where the effect of a mutation is magnified by the binding strength of the mutation position:

$$\begin{equation} \Delta S_{ij} = (p(\hat{s})-p(s)) \cdot \max(0, p(\hat{s}), p(s)). \end{equation}$$

Two examples of mutation maps are shown below.  On the left, it shows the promoter region of LDL-R gene; LDL-R is the LDL receptor gene and LDL is known as the "bad cholesterol".  The wild type genome contains based C in the middle of promoter, which can be recognized by SP1 protein, which then drives the expression of the LDL-R, and that helps recognize excess levels of LDL and regulates its removal from the blood.  Patients with this base mutated to G seriously decrease the binding score (see the red dot at the bottom), as the mutant version is no longer recognized by SP1.  This leads to the deleterious condition known as hypercholesterolemia (high level of cholesterol) due to the low expression of LDL-R.  On the right, a mutation in the globin gene, from A to G, creates a new binding site TGATAA attracting a TF protein called GATA1.   The binding of GATA1 blocked the genome region originally reserved to serve as the promoter of globin.  As the result, globin is no longer expressed at sufficient level, which leads to a disease called $\alpha$ thalassemia (insufficient amount of hemoglobin to carry oxygen in blood).

How DeepBind is Optimized?

The DeepBind system contains many hyperparameters, including motif length $m$, number of motifs $k$, learning rate, momentum, batch size, weight decay, dropout rate, etc.  The training process first randomly sample 30 sets of such hyperparameters (see Figure below), using 3-fold cross validation to evaluate their performance.  The best performing set of hyperparameter is then used to repeat the DeepBind training process using the whole training dataset three times.  This is to avoid statistical fluctuations during the optimization process, the best performing set of DeepBind model is then chosen as the final model for the given protein and used to predict the test dataset and produce performance metrics.  One DeepBind model was trained per protein per experimental data type.

Conclusion

CNN is a DL architecture suitable to discover translation-invariant patterns in inputs.  DeepBind system uses CNN to capture protein-binding motifs on DNA/RNA, where NN architecture enables it to capture some non-linearity effects and/or aggregation effects of multiple underlying motif patterns, as well as aggregate motifs that are widely spaced.  The study showed DeepBind surpassed traditional bioinformatics motif discovery tools and was able to explain a number of disease causing mutations.  Our visual estimation shows the improvement, compared to traditional approach, is at about a few percentage, similar to what we generally see from DL applications in bioinformatics.  However, the architecture of DeepMind is a lot simpler compared to the many sophisticated motif discovery algorithms developed over the passed two decades.  The authors did not offer explanations on why DeepBind performed better.

Convolutional Neural Network (CNN)

Disclaimer: most contents in this blog are my own speculations on how human vision works; I did not research if they are supported by literature.  I am enthusiastic about CNN; these speculations have greatly helped boost my own understanding of CNN, therefore, I hope they could be useful to you as well, even if they might be nonsense. Bear in mind that the purpose of this blog series is try to check the "why" aspect of deep learning.  Nowadays, you for sure can find many "how-to" tutorials on CNN elsewhere.

Introduction

Given a cat photo, we can recognize the cat absolutely effortlessly.  How does it work?

First, we do not think the brain relies on an object model approach.  That is it has a cat model consisting of an assembly of some geometric shapes representing eyes, ears, whiskers, tails, paws, etc.  In fact the object model approach is intractable, as you would need a model for a sitting cat, for a sleeping cat, for a cat with paws hiding, ..., just seem too complicated to be feasible [1].  Anyway, in this approach, the brain will also need to contain models of millions of other objectives, then apply all these models to the input image to determine which model is the closest match.  This modelling approach does not really explain how we see cats in the photos below, as they would not match any reasonable cat model well.

 Figure 1. Cat photos that do not match any cat models well.

Second, computer vision experts have spent decades to construct image features that are robust against translation, rotation and scaling [2], as our brains seem to have no problem in recognizing an object under such transformations.  Do our brains use similar tricks to represent objects with some mathematically-elegant transformations?  It does not seem so, at least, our vision recognition system is not very capable of correcting rotation transformations.  The two faces below are $180^\circ$ rotations of each other.  If our brains transform them into some rotation-invariant representation, we would expect to see the second face by staring at the first, but we do not.

Figure 2.  We cannot see the other face, if we look at one, as our brains do not transform the image into rotation-invariant representation (source).

Third, we do not seem to see details, but only see concepts instead.  Look at the scenery photo below for two seconds, now cover it up.  I assume you can tell it is about a beach, and you might probably saw some other objects such as beach, mountain, ocean, bridge, sky.  Now are you able to tell if there are many people on the beach? Is there a lagoon on the left? Are there wild flowers on the hill slope? are there people in the ocean? are there clouds in the sky?  All those pixels have been sent to our brains as input, yet we cannot answer these questions, as we "looked at" these pixels but we did not "see" them.  This means the output of our vision system probably is a few abstract concepts and lots of details were not really passed to the memory or reasoning regions of the brain.  This is great, as we immediately see beach, without the need of prioritizing roads, people, waves, lagoons, bridges, and then try to follow some complicated decision-making rules to come to the conclusion that this is a beach.  We simply can see the meaning of each image without reasoning.  It only happens in movies, when an agent can name everything in a room after a few seconds; photo memory does not seem to be what our normal people see things, otherwise, the brain will be drained by too many objects and have a hard time capture a key event.  This means in computer vision, details are often not as important as you might have expected.

Figure 3.  Beautiful Torrey Pines beach at Del Mar, San Diego (source).

Computer vision is the field where deep learning really shines, as it provides a feasible approach to see objects the way no previous machine learning technologies was capable of.  Personally, I am convinced this might actually be how our vision system works more or less, as it could help explain many vision phenomena, including the examples mentioned above and some examples to be discussed at the end of the blog.  Here, we will look at how a specialized neural network (NN) called convolutional neural network (CNN) works.

The Rationale for CNN

Vision recognition for ImageNet is basically a function $f(x)$, where $x$ is the input image and output is the object id(s) representing the concepts in the image.  We mentioned in the previous chapter that a fully-connected NN (FNN) can be used to model any function, such as the FNN in Figure 4 (left).  In this example, not all pixels on the input image are relevant, therefore the NN can be simplified into a localized NN (where weights for non-dog pixel inputs are zeroed out, Figure 4, right) with the dog at its center of the input.

Figure 4.  FNN for image object recognition is a localized NN, as only the pixels belonging to the object contribute to the recognition.

An image can contain multiple instances of dogs, to recognize all the instances, we only need to slide the localized NN over the input image, whenever there is a dog, the NN produces a response signal.  This operation of sliding a function over an image to obtain a new response image (Figure 5) is known as convolution in computer vision.  For the ImageNet task, we only care about whether the image contains dogs or not without needing to know how many and where they are.  Therefore, our NN has to be translation invariant, i.e., $f(x - \delta) = f(x)$.  This can be easily implemented by applying the $max$ operator onto the output convoluted image.

Figure 5.  Apply the same localized NN over the input image allows it to recognize multiple instances, where it produces a signal at each matched location.

In traditional computer vision, we construct a template (a dog model) and slide that template over the input image, then the output convoluted image shows bright spots where template matches the input.  However, it is unrealistic to have a dog template (model) that can be used to recognize all dogs, as not only dogs look very different, maybe only a partial of its body is exposed, or it is a highly distorted cartoon dog, etc.  A more robust way is to model a dog as a collection of its body parts that are loosely connected.  If we have a function $\mathbf{g}(x)$ that has multiple dimensions aiming to recognize dog eyes, nose, ears, spot pattern, tails, respectively.   Each dimension is a convolution function mentioned above that produces an output image, we call a feature map.  Each dimension still uses a convolution filter, which is basically a weight matrix describing what pattern this feature is looking for, and the intensity of the output feature map represents the likelihood (not strictly in the probability sense, but just monotonically reflecting the likelihood) of the presence of that particular body parts at a corresponding underlying location.  I.e., we assign weights to each dimension and produce a final output to indicate if the dog is present, i.e., $f(x) = \mathbf{w(\delta)}\cdot \mathbf{g}(x-\delta)$.  $\delta$ here accounts for the relative location shift among the signals from different feature maps.

What is important here is this weighting (i.e., convolution) operation does not need to be applied to the raw image pixels, it can be applied to previous feature maps, where pixels encoding the presence of smaller body parts.  To understand this, say we aim to construct a feature representing a dog face, which consists of preliminary parts such as eyes, nose, ears and are arranged in some lose geometric positions.  This means the feature maps of eyes, nose and ears have already been generated and they will become the input images for the detection of the dog face feature; each input feature map (nose or eye) is a separate input image channel (imagine the raw image has three color channels); they are convoluted to produce a new dog face feature map.  The convolution relies on a weight matrix, which can implement loose logic, such as the eyes are within a fuzzy region and the nose can be in another fuzzy location.  Such filters sitting on top of lower-level feature maps can therefore encode a very flexible face model, unlike the traditional template-based match (Figure 6).

Figure 6.  One convolution weight matrix (left) can encode very different faces.  For the weight matrix, red and orange indicates the possible locations of eyes, blue for nose and green for mouth.  This matrix can match the three very different faces on the right equally well.  This is not something the traditional template-based method is capable of.

  Similarly, the feature maps of dog face, legs, bodies will need to be further combined to produce an even higher-level feature map say represents a sitting dog.  In the other direction of the feature hierarchy, the feature map of dog eyes can be decomposed into more simpler imaging feature maps of circles, different color patches, lines, curves, corners, etc.  Therefore, at the end we are looking at a multi-layer hierarchical convolutional neural network, where each layer relies on a small input image section (either the raw input image or the output feature maps of a previous layer), it performs convolution and produce the output feature map encoding whether the feature it is looking for is present at specific locations.  At the top layer, we form a collection of feature maps, where each tell us whether certain high-level body parts are observed with high probability, the weighted collection of this information allows us to make a decision on whether the object of interest is present.  This way, a dog can be detected even if a significant portion of its body parts are missing, or the relative geometry of its parts are highly distorted.

Another practical requirements for object recognition is to be able to scale down an image.  Despite there are many differences in details for the faces in Figure 7, those details are irrelevant for our purpose of recognizing them all as face.  So if we scale down the image into low-resolution representations, where the irrelevant details are smoothed out, it will be much easier for us to match the smaller image with some convolution filter encoding a face.  The process of taking a higher-resolution image and scaling it down into a lower-resolution one is called pixel pooling, i.e., pool multiple pixels into one pixel.  Pooling offers a few advantages.  First, pooling get rid of details, so that a face filter (weights) only requires two eyes on top followed by a nose below and a mouth below, but it does not mandate the exact looks of the eyes, nose or mouth, there are lots of flexibility in such filter elements.  Second, pooling brings the underlying part objects closer in space, so that the filter to convolute these objects into higher-level object can be a lot smaller.  Third, it is computationally much more efficient to compute and store without compromising our goal of object recognition, as fewer pixels are involved.  The most-often used pooling mode in CNN is max pooling, where the maximum value of the several original pixels is retained as the pooling output.

Figure 7. Scaling down a complex image helps remove irrelevant details for a recognition task.

In the traditional template matching strategy in image analysis, the recognition of individual parts and the relative geometric location of these parts are simultaneously encoded within the template at the pixel level.  This seriously restricted the variation of the objects a template can represent.  In the CNN approach, each hierarchical feature allows room for variation within individual features (e.g., there can be multiple eye features of different styles) as well as the spatial variations among them, and the pooling further provides even more rooms for the "disfiguring", therefore, CNN no longer recognizes objects at the pixel level, but start to incorporate tolerance and somewhat represent objects at the conceptual level.  I think this might be the main reason behind the superb performance of CNN in real-life object recognition.

CNN for Face Recognition

With convolution and pooling layers as building elements, we are ready to perform what used to be the most difficult computer vision challenge, classify photos according to the object(s) it contains (the ImageNet challenge).  Let us look at how CNN can be applied to recognize a face.

Figure 8. Architecture of a CNN that recognize a $16 \times 16$ image containing a face.

As shown in Figure 8, our input (A) is a black-and-white face of size $16 \times 16$.  The first hidden layer (B) actually is a group of layers of the same size $16 \times 16$ each, and each is the result of $2 \times 2$ convolution across the input image (this means each layer is the heatmap resulted from convoluting a particular neuron across the original image).  The weight matrix $w$ for one of the feature (blue) reads (the matrix is shown at the lower-left corner of the feature maps):

$$w = \begin{bmatrix} 0 & 0 \\ 1 & 1 \end{bmatrix},$$

where 1 stands for black and 0 for white.  It basically detects an edge (black object on white background) facing top.  The resultant heatmap is depicted below the layers associated with group B in Figure 8.  The convolution sum is passed through a non-linear function, typically ReLU function.  ReLU has no effect in our specific example, as we only consider positive weights, but weights can be negative in real applications and there is a threshold parameter $b$ to suppress signals that are too weak (remember ReLU = $\max(0, signal)$.  Similarly we can build three other features to detect the bottom-facing (orange), right-facing (green) and left-facing (purple) edges, the resultant heatmaps are shown in Figure 8 as well.  There will be other feature maps in this group capturing other shapes, however, let us focus on these four feature maps in this particular example.

The image is big, so let us down sample it.  The next step C is a max-pooling step.  The example here is to shrink a $2 \times 2$ area into one pixel, where the value of the output pixel is the maximum of the four original pixels.  This will produce four new feature layers in the second hidden layer group C, each feature map is now only one forth of the size of those in the previous layer B, accounting for a 75% saving on computational resource.  There is no cross-talking among the feature maps in this pooling step, all what we do is to cut the resolution in half feature by feature, which also brings the response signal 50% closer to each other.

Then it uses another convolution layer D.  The difference here is the convolution is a 3D convolution, i.e., the weights not only covers the width and height of a sub-region in the previous layer C, but also the depth, i.e., it mixes signals from all the feature maps in the previous layer under that sub 2D region.  A convolution filter corresponding to the $2 \times 2 \times 4$ weight pattern as shown in Figure 8 produces an new output feature map in the D group representing an eye feature.  Where it looks for an activated pattern consisting of 2-pixel top, a 2-pixel bottom, a 2-pixel left and a 2-pixel right -facing patterns in the previous layer C.  This filter means the neuron is detecting a black square, it produces a strong activation signals on its output heatmap, where an eye (represented by a black square) is present.  Notice the color for the right eye is a bit lighter, as it is not a perfect match (but still a good match) to the eye filter (therefore the activation is less strong).

Figure 9. Weights (filter) associated with the eye feature.

Similarly, layer D also contains a feature map representing where noses are and another representing where mouths are.  We then do another max-pooling, which also pulls the detected eyes, nose, and the mouth closer to each other, actually now they all fall within a $3 \times 3$ neighborhood.  Layer F contains a new features representing a face, whose $3 \times 3 \times 3$ filter consists of the convolution of signals from eye, nose and mouth feature maps.  It generates an activated pixel in the output layer G, where we know a face is detected.  Here we see an important role of pooling is to bring lower-level features closer together, so that a higher-level feature can be constructed with a much smaller convolution windows.  If our CNN does not contain any pooling layers, a face filter would require a convolution matrix of size $16 \times 16$, which would use too many parameters and contains unnecessary details.

The face filter is a 3D weight matrix, which can be written as:

$$w_{face} = \begin{bmatrix} w_e \\ w_n \\ w_m \end{bmatrix}, \\
w_e = \begin{bmatrix} 1 & 0 & 1  \\ 0 & 0 & 0 \\ 0&0&0 \end{bmatrix}, \\
w_n = \begin{bmatrix} 0 & 0 & 0  \\ 0 & 1 & 0 \\ 0&0&0 \end{bmatrix}, \\
w_m = \begin{bmatrix} 0 & 0 & 0  \\ 0 & 0 & 0 \\ 0&1&0 \end{bmatrix}.
$$

So now we use convolution and max-pooling to build a CNN that is capable of producing an activated output, when we show it a face.  If the input image is twice larger and contains a tiling of four faces, our output layer will contains four activated pixels on the face feature map, each signaling a face.

Now imagine there are many other feature maps within each hidden layer, capturing all kinds of interesting features, from simple features such as edges, diagonal lines, to medium features such as rectangles, circles, to more sophisticated features such as wheels, front of a car, then to rather complex features representing face, car, bike, dogs, etc.  We can see how single features at the earlier hidden layers can be combined (depth-convoluted) to produce median complex features and then further combined to create more and more complex features.  The last output layer of CNN can contains thousands of high-level features and their activation status.  Those features then serve as the input for a FNN, which typically ends with a softmax layer to produce probabilities for each object class.  The top scoring classes represents the objects we see in the image (we use softmax, as this is a multi-class classification problem, see previous blog).

In practice, we certainly do not hand construct the feature maps ourselves.  The parameters of the CNN filters would be initialized with random values, we then simply feed the CNN with lots of input images and they true labels.  The loss function will be a cross-entropy term to quantify the prediction accuracy.  The result of training and parameter optimization will enable the CNN (including the FNN layers) to determine the optimal filters, i.e., features, for each layer.  As the CNN are trained by many faces: long, round, man, woman, etc., there will be many neurons representing different faces, therefore, we do not relying on one face feature to recognize all faces.

Although the features found by CNN are rarely as clean as our face-detection example, they generally resembles such concepts.  An extremely nice videos is available here, it should convince you some CNN neurons (representing a feature) indeed represent rather meaningful objects, e.g., an edge, a face, or text.

VGGNet for ImageNet

The huge success of CNN is reflected it its ability to classify images provided by ImageNet into 1000 classes, this is the well-known ImageNet competition. The current state-of-the-art CNN classification already surpassed human performance of 5.1% error rate [3], therefore, this task is considered a solved problem.  VGGNet (Figure 10) is the champion of 2014 ImageNet competition [4], we discuss it because it has the nice architecture very similar to our face-recognition example, i.e., very easy to understand.

Figure 10. VGG16 architecture. (source)

It takes the input image of size $224 \times 224 \times 3$, do two $3 \times 3$ convolutions followed by a max pooling. Then repeat this a few times.  It gradually doubles the number of features, from 64 to 128, then to 256 and finally to 512.  Let us walk through the statistics in Figure 10, the input image has a size of 224 pixels, the size is preserved after convolution, then halved by POOL2 into 112, continued the cycle of convolution and pooling until it gets $7 \times 7 \times 512$ at the last pool layer, that is we identify 512 features, each corresponds to a likelihood heatmap of size $7 \times 7$.  All these neurons become the input features for a 3-layer FNN to finally generate 1000 probability values for the 1000 object classes.   Why would one uses two $3 \times 3$ convolution layers sequentially instead of just one convolution layer?  Two of such $3 \times 3$ filters basically play the role of a single $5 \times 5$ filter, however, the former has only 20 parameters while the latter has 26 parameters, plus two layers can model more non-linearity than a single layer.  The saving is more significant if we replace a $7 \times 7$ filter with three $3 \times 3$ filters.  So CNN tends to use multiple small filters than a large filter.  Strictly speaking, we are actually doing 3D convolution, as the first filter is $3 \times 3 \times 3$, given a three-channel color input image.  The second filter is $3 \times 3 \times 64$ sitting after the first hidden layer.

In addition, we see there are a total of 138 million parameters in VGGNet!  Among them, 102 million are used for the first fully connecting layer.  So CNN layers are very efficient in terms of the number of parameters it consumes, as it uses small filters.  On the other hand, the fully connected layers are expensive for optimization.  Microsoft's ResNet [5] does a better job in extracting CNN features, and it only requires a single fully connected layer at the end for the classification, therefore, it is way more efficient.  ResNet was the champion of 2015 competition, with an error rate of 3.57%!

Machines such as VGGNet took a long time to train.  In our own projects, we often cannot offer such long training time and the CPU/GPU resource, instead we take advantage of a technique called transfer learning.  E.g., if we are tasked to classify a subset of ImageNet photos into either dog or cat, instead of the original 1000 classes.  VGGNet was trained to classify 1000 classes, therefore, is not optimized for our simpler problem.  Why?  Say bone is one of the 1000 classes, VGGNet cannot use bone features to help dog recognition, as it was also tasked to distinguish bone from dog.  Now with only two classes, each can take advantage of features representing other classes to help boost its confidence.  But to train the VGGNet from beginning would be too costly.  We know almost all the key features required to distinguish a dog from a cat are already captured in the CNN portion of the VGGNet, therefore, we can freeze that part. Instead, we only replace the last three fully connected layers by one fully connect layer producing a single output, representing the probability of the dog class.  The new CNN can be optimized much quicker.

Miscellaneous Topics

Here we discuss some non-technical topics of interest and see how CNN can help us better understand some recognition activities.

Objects of Different Scales

The CNN described above can deal with object translocation, but it does not deal with scaling.  If our CNN was only trained with faces of the same scale, then we feed it with faces of other scales during its application, it is likely to fail.  Because if the eyes are way bigger than what it saw, CNN is still looking at a pupil at its last layer, there will be no activation for the eye feature map in the earlier layer.  For the same argument, the face neuron will not be activated as well.  When the input face is too small, it already converges into a tiny $2 \times 2$ eye, nose, mouth-like features at a layer earlier than where the these features are expected, therefore it will be missed as well.  To overcome this, some people use multi-scale CNN [6] consisting of multiple independent CNNs, each trained with input training images at different scales.  Their outputs are combined to be used in the final classification.  You can picture this approach as taking the input image, shrink it and magnify it into multiple copies, and then sent them to different people specialized at looking at large face or small face, respectively, one of them will probably be able to recognize the object.

Instead of using multiple CNNs, we could probably train one CNN using image augmentation technique.  That is we turn one training image into multiple images of varying resolutions, then use all of them to train one single CNN.  With the CNN being trained by objects of different scales, it will automatically construct features of varying scales and be able to deal with scaling challenge.  This certainly implies we need a lot more feature maps in this CNN, a large eye feature and a small eye feature, ...  Human brain contains 86 billion neurons, comparable to 300 billion stars in our the Milky Way [7].  Therefore, our brain can afford lots of neurons to form a rather gigantic network to cope with scaling issue.  Considering the pair of human eyes of a three-year old kid already consumed hundreds of millions of images, human CNN does not lack of training data [1], so this could work for real.

Why a CNN trained above could recognize both a tiny face and a large face? A tiny face might be recognized by a small-face feature map in one of the early layers of our CNN trained by data augmentation, how will this activation survive the multiple convolution and pooling layers afterward?  This is at least mathematically possible.  First, max pooling is different from average pooling, an activated signal on a heatmap does not disappear after pooling.  Second, an activation can survive future convolution, if the convolution filter happens to be a trivial filter, for instance:

$$w = \begin{bmatrix} 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \end{bmatrix}$$

is a $5 \times 5$ filter that will preserve the signals on a feature map in its previous layer.  Thus CNN trained by small faces will construct a neuron capturing a tiny face at one of the earlier layer, it then  build an identity neuron in each of the following layers to simply propagate the activation face signal all the way to the output layer.  Human brain likely consists of CNNs of mixed layering structure, so the tiny-face neuron could be wired directly to the output layer.  ResNet introduce shortcut wiring structure, probably can also enable it to be efficient in relaying an early-detected features to the back.

Illusion - Empirical Activation

In the CNN model, the last feature layer contains thousands of features, where each one gives an activity representing how strong the feature is present in the input.  Our brain only see some concepts, which might be just an unordered or loosely-ordered collection of a list of such activated features.  E.g., the concept of a beach scene is a list of features such as beach, sea, sky.  In the left side of Figure 11, we easily see a smiling president.  This might be because we have a "smiling" feature neuron activated due to its underlying features such as smiling mouth, smiling cheeks are activated.  A "face" neuron is also activated.  Although our CNN does not rotate the image $180^\circ$, we have seen enough upside down faces, so that we have an upside-down face neuron for such cases.  With both neurons activated, we see the concept of a "smiling"+"face" equals "smiling face".  Nothing is alarming to our brain.  On the contrary, the right side photo looks weird, although we all know its is simply the left photo.  Here, our "upside up face" neuron is activated.  However, none of the typically facial expression neuron is strong activated, at least not those facial expression neurons normally associated with the president.  So we see an "eerie face".  We recognize the president, as there are sufficient amount of president-specific features, but we have reservations, as some features that expected to be activated go silent.

Figure 11. Smiles and face are recognized independently.

Figure 12 shows a study, where researchers found out the time it took to recognize a person is shortest with the normal photo, slower with one eye moved up, and the slowest with two eyes moved up [8].  This seems rather straightforward to explain with CNN.  With the normal photo, many feature neurons fire up and collectively let us see Elvis or JFK easily.  When one eye is moved, the half face with a moved eye does not activate some feature neurons, so our concept of Elvis/JFK does not contain all the necessary feature members, we are unsure whom we see.  Then our eyes explore the photo and focus on the normal half of the face, it immediately activate all the necessary feature neuron needed and the Elvis/JFK concept is seen.   With both eyes moved, our eyes wonder around,  Elvis/JFK face neuron is still not firing.  Then we focus on the top portion and ignore the bottom, then focus on the bottom half and ignore the top, the Elvis/JFK concept is probably the closest concept to match either output, therefore, we go along with the decision.  I can feel I am doing a decision making (reasoning) step here, unlike no decision making is needed for the normal photo.

Figure 12. Research by Cooper and Wojan (source)

In Figure 13, they sure look like two roads extending into two different directions, but in fact, these are identical photos.  What might have happened is we have seen many scenes in our life similar to the Flatiron building in the New York City (Figure 14), where roads branching out in a "V" shape.  If there is a "V" shape feature neuron in our CNN, that same neuron fires strongly when the two photos in Figure 13 are viewed at once.  Our brain has no way to suppress/ignore this "V" firing, this feature contribute strongly to the concepts that we are looking at two different roads.

Figure 12. Two identical photos look very different.

Figure 13. The Flatiron building in New York City.

CNN is also the hero behind many most well known AI systems, such as Alpha GO/Zero.  For bioinformatics applications, CNN can also be applied to recognize one-dimensional patterns, such as to binding motifs in DNA sequences.  For medical applications, it can be applied to identify 3D tumor patterns in MRI or CT scans.

Let us look at a bioinformatics application next, where CNN is used to identify motifs for DNA/RNA-binding proteins (next blog).

Reference

1. https://www.ted.com/talks/fei_fei_li_how_we_re_teaching_computers_to_understand_pictures#t-370742
2. https://en.wikipedia.org/wiki/Scale-invariant_feature_transform
3. http://karpathy.github.io/2014/09/02/what-i-learned-from-competing-against-a-convnet-on-imagenet/
4. http://cs231n.stanford.edu/slides/2017/cs231n_2017_lecture9.pdf
5. https://arxiv.org/pdf/1512.03385v1.pdf
6. https://academic.oup.com/bioinformatics/article-lookup/doi/10.1093/bioinformatics/btx069
7. https://www.nature.com/scitable/blog/brain-metrics/are_there_really_as_many
8. http://www.public.iastate.edu/~ecooper/Facepaper.html