第一節：基本知識(2)

但凡預測就會存在誤差，在這裡我們再定義一個『損失函數』\(diff()\)：

\[loss = diff(y, \hat{y})\] - 以簡單線性迴歸的損失函數為例，所求的值為殘差平方和，可將此式改寫為：

\[loss = diff(y, \hat{y}) = \sum \limits_{i=1}^{n} \left(y_{i} - \hat{y_{i}}\right)^{2}\]

由於『損失函數』內存在著『預測函數』的Output，因此我們可以簡單的將『預測函數』與『損失函數』進行合併：

\[loss = diff(y, f(x))\]

我們先以線性迴歸為例在上式進行全部展開，可以獲得下列結果

\[loss = diff(y, f(x)) = \sum \limits_{i=1}^{n} \left(y_{i} - \left(b_{0} + b_{1}x_{1,i}\right)\right)^{2}\]

第一節：基本知識(3)

顯而易見的是，我們會希望我們會希望得到一個『預測函數』使\(loss\)越小越好，因此我們可以將上述問題轉變為一個『優化問題』(求函數極值)：

\[min(loss)\]

所以，至此為止我們要了解到，其實我們的任務其實就是給定一組『預測函數』與『損失函數』，而其中『預測函數』將包含若干『待解參數』，而透過上述求解『優化問題』的過程，我們將可以獲得所有的『待解參數』，而在最後我們就能使用『預測函數』進行後續的預測。
這裡我們就不演示如何求解任何方程式了，有興趣的同學的以在學期中修習「機器學習與演算法」及「深度學習理論與實務」等兩個課程。

第二節：決策樹(1)

我們現在要介紹一個最簡單的機器學習法：決策樹，決策樹牽涉到了最少的數學，僅僅是使用電腦大量運算來進行分類。
決策樹的原理是，將所有樣本不斷地進行分類，透過樹狀擴展規則最終將樣本分類到某一個群落，從而決定其樣本結果

F01

決策樹的原理是，在大量的變項(X)中找出「一個」能對「該樣本」有「最大分類效果」的規則，並按規則持續把樣本分割，直到再也找不出規則為止！

F02

第二節：決策樹(2)

讓我們來試著用iris的資料實作決策樹吧，請到這裡下載這份資料集，你可以用這樣的方式把它讀取進來。

– 需要特別注意的是，對於資料科學實驗的流程，我們通常會把樣本分為3個部分，分別是：

訓練組(Training set, Development set)：負責用來建構一個預測模型，我們之前學的三大流程就是用在這上面的，樣本可以隨意調整
驗證組(Validation set, Tuning set)：不參與模型訓練，但會用來指導模型訓練，以及建構重要的資訊所用，樣本可以隨意調整
測試組(Testing set, Hold-out set, Validation set)：最終模型會在上面運行一次(只能一次)，已確定最終的準確度，原則上樣本選取必須符合未來使用條件

– 實務上，訓練組跟驗證組在不選擇模型的前提下可以合併，但測試組是必須的。

iris <- read.csv('data/iris.csv')

#Split data

set.seed(1)

Train.sample <- sample(1:150, 100, replace = FALSE)

Train.data <- iris[Train.sample,]
Test.data <- iris[-Train.sample,]

第二節：決策樹(3)

讓我們直接使用套件「party」進行分析：

library(party)

tree.model <- ctree(formula = Species ~ ., data = Train.data)
tree.model

## 
##   Conditional inference tree with 3 terminal nodes
## 
## Response:  Species 
## Inputs:  Sepal.Length, Sepal.Width, Petal.Length, Petal.Width 
## Number of observations:  100 
## 
## 1) Petal.Length <= 1.9; criterion = 1, statistic = 92.735
##   2)*  weights = 32 
## 1) Petal.Length > 1.9
##   3) Petal.Width <= 1.6; criterion = 1, statistic = 48.709
##     4)*  weights = 35 
##   3) Petal.Width > 1.6
##     5)*  weights = 33

接著我們來利用這個決策樹模型做預測，來看看結果：

pred.y <- predict(tree.model, Test.data[,1:4])
table(pred.y, Test.data[,5])

##             
## pred.y       setosa versicolor virginica
##   setosa         18          0         0
##   versicolor      0         14         3
##   virginica       0          0        15

使用內建函數畫圖看看，效果是不是很不錯呢?

plot(tree.model)

第二節：決策樹(4)

讓我們使用真實資料來進行實驗吧！請至這裡下載範例資料

dat <- read.csv("data/ECG_train.csv", header = TRUE, fileEncoding = 'CP950', stringsAsFactors = FALSE, na.strings = "")

這份資料共包含了5000人，應用目標是希望透過心電圖的參數去預測幾個重要的指標，包含了：

AMI：這是個類別變項描述心肌梗塞的狀態，包含STEMI、NSTEMI及not-AMI
K：這是一個連續變項描述鉀離子的濃度
LVD：這是一個二元類別變項：1代表left ventricular dysfunction，0則代表正常
time與death：這組變項描述病患隔多久後死亡與否，這用來做存活分析之用

– 除了性別(GENDER)和年齡(AGE)外，心電圖的重要參數包含了8個連續變項特徵(Rate、PR、QRSd、QT、QTc、Axes_P、Axes_QRS、Axes_T)以及31個二元類別變項描述相對應的rhythm。

– rhythm依序為：abnormal T wave、atrial fibrillation、atrial flutter、atrial premature complex、complete AV block、complete left bundle branch block、complete right bundle branch block、first degree AV block、incomplete left bundle branch block、incomplete right bundle branch block、ischemia/infarction、junctional rhythm、left anterior fascicular block、left atrial enlargement、left axis deviation、left posterior fascicular block、left ventricular hypertrophy、low QRS voltage、pacemaker rhythm、prolonged QT interval、right atrial enlargement、right ventricular hypertrophy、second degree AV block、sinus bradycardia、sinus pause、sinus rhythm、sinus tachycardia、supraventricular tachycardia、ventricular premature complex、ventricular tachycardia、Wolff-Parkinson-White syndrome

– 在大部分的狀態下，我們會使用後面的幾個變項去預測前面的4組變項。

第二節：決策樹(5)

– 這裡給一些範例語法，讓大家知道怎樣用在不同的依變項屬性上。

這是連續變項的範例：

subdat <- dat[!(dat[,'K'] %in% NA) & !(dat[,'Rate'] %in% NA) & !(dat[,'AGE'] %in% NA), c('K', 'Rate', 'AGE')]

tree.model <- ctree(formula = K ~ ., data = subdat)
plot(tree.model)

這是二元分類的範例：

subdat <- dat[!(dat[,'LVD'] %in% NA) & !(dat[,'GENDER'] %in% NA) & !(dat[,'Rate'] %in% NA), c('LVD', 'GENDER', 'Rate')]
subdat[,'GENDER'] <- as.factor(subdat[,'GENDER'])
subdat[,'LVD'] <- as.factor(subdat[,'LVD'])

tree.model <- ctree(formula = LVD ~ ., data = subdat)
plot(tree.model)

這是多分類任務的範例：

subdat <- dat[!(dat[,'AMI'] %in% NA) & !(dat[,'GENDER'] %in% NA) & !(dat[,'AGE'] %in% NA), c('AMI', 'GENDER', 'AGE')]
subdat[,'GENDER'] <- as.factor(subdat[,'GENDER'])
subdat[,'AMI'] <- as.factor(subdat[,'AMI'])

tree.model <- ctree(formula = AMI ~ ., data = subdat)
plot(tree.model)

第二節：決策樹(6)

決策樹有非常多的參數可以讓我們設定，我們這裡舉幾個常用的。

參數【mincriterion】是顯著水準
參數【maxdepth】是樹的最大深度
參數【minsplit】是說樣本大於多少，才考慮繼續分類
參數【minbucket】是說每個分類至少要有幾個樣本

subdat <- dat[!(dat[,'LVD'] %in% NA) & !(dat[,'GENDER'] %in% NA) & !(dat[,'Rate'] %in% NA), c('LVD', 'GENDER', 'Rate')]
subdat[,'GENDER'] <- as.factor(subdat[,'GENDER'])
subdat[,'LVD'] <- as.factor(subdat[,'LVD'])

tree.model <- ctree(formula = LVD ~ ., data = subdat,
                    controls = ctree_control(mincriterion = 0.95, maxdepth = 2, minsplit = 20, minbucket = 7))
plot(tree.model)

– 如果你想獲得預測機率，可以用這種方式：

prob_list <- predict(tree.model, subdat, type = 'prob')
prob_y <- do.call('rbind', prob_list)[,1]

第三節：梯度提升機(1)

讓我們再跟大家介紹一個機器學習模型：梯度提升機(Gradient Boosting Machines)
梯度提升機是一個超級受歡迎的機器學習法，在許多領域上都有非常成功的表現，也是Kaggle競賽時常勝出的主要演算法之一。

– 為什麼他這麼強呢?決策樹的缺點是模型都是建立在單一預測模的基礎上，容易因為樣本決定一個特定的極端條件。

所謂的「Boosting」，主要概念就是將新模型「有順序、循序漸進」的加入集成學習。在每一次建立新模型時都會新增一個新的弱(weak)模型，針對到目前為止集成學習的誤差進行訓練。

F03

第三節：梯度提升機(2)

根據剛剛的圖你會發現，每個模型其實都是根據之前模型所無法預測的部分繼續加強，因此其實我們不需要要求第一個模型就很準了，這也是為什麼我們希望一個「弱」模型的原因。

– 有趣的是，這樣模型的組合將可以讓我們擁有一個非線性模型。

我們把他的數學式寫下來(參數\(\eta\)是收縮率，通常是0.3)：

先建立第一個模型：\(y = F_1(x) + \epsilon_1\)
根據第一個模型的殘差建立第二個模型：\(\epsilon_1 = F_2(x) + \epsilon_2\)
先將上面兩個式子合併，產生一個整合模型，並取得整合模型的殘差：\(y = F_1(x) + \eta F_2(x) + \tilde{\epsilon_2}\)
再根據這個殘差建立第三個模型：\(\tilde{\epsilon_2} = F_3(x) + \epsilon_3\)
再將三個式子合併，產生一個整合模型，並取得整合模型的殘差：\(y = F_1(x) + \eta F_2(x) + \eta F_3(x) + \tilde{\epsilon_3}\)
再根據這個殘差建立第四個模型：\(\tilde{\epsilon_3} = F_4(x) + \epsilon_4\)
再將三個式子合併，產生一個整合模型，並取得整合模型的殘差：\(y = F_1(x) + \eta F_2(x) + \eta F_3(x) + \eta F_4(x) + \tilde{\epsilon_4}\)
依此類推…

第三節：梯度提升機(3)

本週讓我們直接使用真實資料來進行實驗吧！我們就以二元分類預測LVD為例。

– 注意原始資料中存在「遺漏值」，一個比較簡單的方法是透過套件「mice」進行「多重插補法」

library(mice)

subdat <- dat[!(dat[,"LVD"] %in% NA), c(-1, -2, -4, -5)]

subdat[,'LVD'] <- as.factor(subdat[,'LVD'])
subdat[,'GENDER'] <- as.factor(subdat[,'GENDER'])
for (i in 1:31) {subdat[,paste0('rhythm.', i)] <- as.factor(subdat[,paste0('rhythm.', i)])}

used_dat.x <- subdat[,-1]
mice_dat <- mice(used_dat.x, m = 1, maxit = 5, meth = 'cart', seed = 123, printFlag = FALSE)
impute_dat.x <- mice:::complete(mice_dat, action = 1)

set.seed(0)
all_idx <- 1:nrow(subdat)

train_idx <- sample(all_idx, nrow(subdat) * 0.6)
valid_idx <- sample(all_idx[!all_idx %in% train_idx], nrow(subdat) * 0.2)
test_idx <- all_idx[!all_idx %in% c(train_idx, valid_idx)]

train_X <- impute_dat.x[train_idx,]
valid_X <- impute_dat.x[valid_idx,]
test_X <- impute_dat.x[test_idx,]

train_Y <- subdat[train_idx,"LVD"]
valid_Y <- subdat[valid_idx,"LVD"]
test_Y <- subdat[test_idx,"LVD"]

第三節：梯度提升機(4)

現在讓我們直接使用套件進行實驗吧，我們選用的套件是「xgboost」。

– 這個套件實現的是「eXtreme Gradient Boosting」，跟傳統的梯度提升機不完全一樣，但原理是類似的。

– 他必須吃矩陣格式的資料，並且還要轉換成自己的格式。

library(xgboost)

train_X_mat <- model.matrix(~ ., data = train_X)
xgb.data_train <- xgb.DMatrix(data = train_X_mat[,-1], label = as.integer(train_Y) - 1L)

valid_X_mat <- model.matrix(~ ., data = valid_X)
xgb.data_valid <- xgb.DMatrix(data = valid_X_mat[,-1], label = as.integer(valid_Y) - 1L)

xgb_fit <-  xgb.train(data = xgb.data_train, watchlist = list(eval = xgb.data_valid),
                      early_stopping_rounds = 10, eval_metric = 'auc', verbose = FALSE,
                      max.depth = 5, eta = 0.3,
                      nthread = 2, nrounds = 100, objective = "binary:logistic")

讓我們直接預測看看結果如何。

library(pROC)

pred_valid <- predict(xgb_fit, valid_X_mat[,-1])
roc_valid <- roc(valid_Y ~ pred_valid)
plot(roc_valid)
text(0.5, 0.5, paste0('AUC = ', formatC(roc_valid[['auc']], 4, format = 'f')), col = 'red')

第三節：梯度提升機(5)

你應該有注意到有許多「超參數」可以調整，你要試著找到一組最佳的參數解。

– 按照我們先前的了解，【eta】應該越小就需要越多的【nrounds】，我們來試試看

xgb_fit <-  xgb.train(data = xgb.data_train, watchlist = list(eval = xgb.data_valid),
                      early_stopping_rounds = 10, eval_metric = 'auc', verbose = FALSE,
                      max.depth = 5, eta = 0.1,
                      nthread = 2, nrounds = 300, objective = "binary:logistic")

讓我們再看看預測結果如何。

library(pROC)

pred_valid <- predict(xgb_fit, valid_X_mat[,-1])
roc_valid <- roc(valid_Y ~ pred_valid)
plot(roc_valid)
text(0.5, 0.5, paste0('AUC = ', formatC(roc_valid[['auc']], 4, format = 'f')), col = 'red')

第三節：梯度提升機(6)

結果更好對嗎?我們把最新的結果用在預測「測試組」上，這才是我們的最終結果，並且我們只能用一次：

test_X_mat <- model.matrix(~ ., data = test_X)
pred_test <- predict(xgb_fit, test_X_mat[,-1])
roc_test <- roc(test_Y ~ pred_test)

plot(roc_test)
text(0.5, 0.5, paste0('AUC = ', formatC(roc_test[['auc']], 4, format = 'f')), col = 'red')

極限梯度提升機可以進行變項重要性分析。

importance_matrix <- xgb.importance(model = xgb_fit)
xgb.plot.importance(importance_matrix = importance_matrix)

第四節：人工神經網路(1)

接著我們再介紹一個強大的預測工具，人工神經網路。
大腦是神經系統的一部分，而他是由多個神經元互相結合而成的。(下圖為Ramón y Cajal在1905年所畫下的神經細胞樣貌)

F04

– 神經細胞的構造如下，不論是何種神經元皆可分成：接收區、觸發區、傳導區和輸出區。

F05

神經細胞傳導訊息的過程是透過電位的變化，在接收區會接收到一個電位改變的訊號，再交給突觸來處理。

F06

透過樹突(dendrite)能接收上一個神經元的訊息，而有些會在接收訊息後產生抑制性作用，有些會產生興奮性作用，然後這些訊號再透過神經元整合，之後再透過軸突(axon)將訊號傳導出去。
我們根據這樣的生物學知識，開始來用電腦模擬一個簡單的神經元。

第四節：人工神經網路(2)

第一代的人工神經網路是1958年由Frank Rosenblatt所發展的感知機(perceptron)，此時還稱不上是神經「網路」，他的結構如下：

F07

讓我們用數學語言來描述他：

\[ \begin{align} \mbox{weighted sum} & = w_{0} + w_{1}x_1 + w_{2}x_2 + \dots \\ \hat{y} & = step(\mbox{weighted sum}) \end{align} \]

仔細一看，他跟線性回歸實在是很像，差別只是在output經過了一個轉換。

第四節：人工神經網路(3)

顯然大腦的構成不是僅僅使用一個神經元，既然邏輯斯回歸的構造與神經元很像，所以我們試著把多個邏輯斯回歸結合在一起吧！

F08

我們同樣也需要用套件來進行運算，讓我們首先載入「mxnet」套件

library(mxnet)

在MxNet中，我們要依序把下面兩個部份準備好：

資料：對於不同的預測函數需要指定不同的結構，你只要記住你想要用的部分即可
模型結構：負責定義預測函數

後面我們依序學習如何編寫這2個部分，並了解怎樣把2個部分結合在一起運算。

第四節：人工神經網路(4)

原來的X矩陣是可以直接用的，但要經過一個轉置

– 在MxNet內，我們用最後一個維度描述樣本數：

TRAIN.X.array <- t(train_X_mat[,-1])
TRAIN.Y.array <- t(model.matrix(~ -1 + factor(train_Y)))

VALID.X.array <- t(valid_X_mat[,-1])

讓我們定義模型架構：

library(mxnet)

data = mx.symbol.Variable(name = 'data')
fc1 = mx.symbol.FullyConnected(data = data, num.hidden = 50, name = 'fc1')
act1 = mx.symbol.Activation(data = fc1, act.type = 'relu', name = 'act1')
fc2 = mx.symbol.FullyConnected(data = act1, num.hidden = 2, name = 'fc2')
out_layer = mx.symbol.SoftmaxOutput(data = fc2, name = 'out_layer')

第四節：人工神經網路(5)

開始訓練，另外我們定義一下模型評估指標：

my.eval.metric.mlogloss <- mx.metric.custom(
  name = "m-logloss", 
  function(real, pred) {
    real1 = as.numeric(as.array(real))
    pred1 = as.numeric(as.array(pred))
    pred1[pred1 <= 1e-6] = 1e-6
    pred1[pred1 >= 1 - 1e-6] = 1 - 1e-6
    return(-mean(real1 * log(pred1), na.rm = TRUE))
  }
)

mx.set.seed(0)

mlp_fit = mx.model.FeedForward.create(symbol = out_layer,
                                      X = TRAIN.X.array, y = TRAIN.Y.array,
                                      optimizer = "adam", learning.rate = 0.001, beta1 = 0.9, beta2 = 0.999,
                                      array.batch.size = 20, num.round = 20,
                                      ctx = mx.cpu(), 
                                      eval.metric = my.eval.metric.mlogloss)

再讓我們進行預測吧！

library(pROC)

pred_valid <- predict(mlp_fit, VALID.X.array, array.layout = "colmajor")[2,]
roc_valid <- roc(valid_Y ~ pred_valid)
plot(roc_valid)
text(0.5, 0.5, paste0('AUC = ', formatC(roc_valid[['auc']], 4, format = 'f')), col = 'red')

覺得準嗎?

第五節：卷積神經網路(1)

神經網路最初的構想是利用電腦模擬大腦運作的過程，而前面那個神經網路中的「全連接層」所模擬的是腦細胞傳遞訊息的思考過程。

– 但回到我們的手寫數字分類問題，當我們看到這些手寫數字時，我們一眼就能認出他們了，但從「圖片」到「概念」的過程真的這麼簡單嗎?

我們是否有個疑問，當光影映射至視網膜時，這些訊號是如何向後傳遞的?

– 現在我們面對的是視覺問題，看來除了模擬大腦思考運作的過程之外，我們還需要模擬眼睛的作用！

第五節：卷積神經網路(2)

1962年時David H. Hubel與Torsten Wiesel共同發表了一篇研究：Receptive fields, binocular interaction and functional architecture in the cat’s visual cortex，這篇研究旨在探討生物視覺系統的運作方式，並獲得了1981年的Nobel prize

F09

– 他們的研究發現，貓咪在受到不同形狀的圖像刺激時，感受野的腦部細胞會產生不同反應

F10

這樣的研究給了我們什麼啟示呢?也就是在最開始接觸光影訊號的腦神經細胞，每個細胞會辨認一種特定的簡單特徵，而且這些簡單特徵「無論出現在視野範圍的哪裡」，他們都會被激活。

第五節：卷積神經網路(3)

電腦學家受到了上面生物研究的啟發，因此在多層感知器(網路中的全連接層)之前增加了一個構造，叫做「卷積器」(convolution filter)

– 卷積器模擬了感受野最初的細胞，他們負責用來辨認特定特徵，他們的數學模式如下：

F11

簡單來說，卷積器透過一個「濾鏡」對圖片進行全局的搜索，而卷積後的新圖片稱為「特徵圖」

– 「特徵圖」的意義是什麼呢?卷積器就像是最初級的視覺細胞，他們專門辨認某一種簡單特徵，那這個「特徵圖」上面數字越大的，就代表那個地方越符合該細胞所負責的特徵。

F12

第五節：卷積神經網路(4)

那假定我們有3個卷積器，分別負責「＼」、「╳」、「／」，那將會產生3個特徵圖，分別如下：

F13

獲得特徵圖之後，還記得我們為了增加神經網路的數學複雜性，會添加一些非線性函數做轉換，因此卷積神經網路在經過卷積層後的特徵圖會再經過非線性轉換。
接著，由於連續卷積的特徵圖造成了訊息的重複，這時候我們經常會使用「池化層」(pooling layer)進行圖片降維，事實上他等同於把圖片的解析度調低，這同時也能節省計算量。

F14

第五節：卷積神經網路(5)

接著我們思考一下，連續的堆疊「卷積器」會產生什麼效果?

– 我們想像有一張人的圖片，假定第一個卷積器是辨認眼睛的特徵，第二個卷積器是在辨認鼻子的特徵，第三個卷積器是在辨認耳朵的特徵，第四個卷積器是在辨認手掌的特徵，第五個卷積器是在辨認手臂的特徵

– 第1.2.3張特徵圖中數值越高的地方，就分別代表眼睛、鼻子、耳朵最有可能在的位置，那將這3張特徵圖合在一起看再一次卷積，是否就能辨認出人臉的位置?

– 第4.5張特徵圖中數值越高的地方，就分別代表手掌、手臂最有可能在的位置，那將這2張特徵圖合在一起看再一次卷積，是否就能辨認出手的位置?

– 第4.5張特徵圖對人臉辨識同樣能起到作用，因為人臉不包含手掌、手臂，因此如果有個卷積器想要辨認人臉，他必須對第1.2.3張特徵圖做正向加權，而對第4.5張特徵圖做負向加權

所以，我們可以在神經網路中堆疊卷積器，而越淺層的卷積器所辨認的特徵越簡單，越深層的卷積器所辨認的特徵越複雜。

F15

第六節：利用卷積神經網路做手寫數字辨識(1)

讓我們一起來實現1989年Yann LeCun所發展的世界上第一個成功的卷積神經網路(Convolutional Neural Networks，CNN)，這個CNN被命名為LeNet。

F16

第六節：利用卷積神經網路做手寫數字辨識(2)

讓我們把剛剛MNIST的手寫數字資料讀取近來

– 請在這裡下載MNIST的手寫數字資料，並讓我們了解一下這筆資料的結構

library(data.table)

mnist = fread("data/MNIST.csv", data.table = FALSE)
mnist = data.matrix(mnist)

n.sample = dim(mnist)[1]

X = mnist[,-1]
X = t(X)
dim(X) = c(28, 28, 1, dim(mnist)[1])
for (i in 1:dim(X)[4]) {
  X[,,,i] = t(X[,,,i] )
}

Y = array(t(model.matrix(~ -1 + factor(mnist[,1]))), dim = c(10, n.sample))

再看一次隨便一張圖以及他的label吧：

library(OpenImageR)

imageShow(X[,,,25])

Y[,25]

##  [1] 0 0 1 0 0 0 0 0 0 0

第六節：利用卷積神經網路做手寫數字辨識(3)

把資料分成Training set以及Testing set：

– 另外，我們需要對X進行標準化(平均值33.4，標準差78.7)，有些人可以嘗試不先標準化的結果，你會發現根本無法訓練！

print(mean(X))

## [1] 33.40891

print(sd(X))

## [1] 78.67774

X = (X - mean(X))/sd(X)

set.seed(0)
Train.sample = sample(1:n.sample, n.sample*0.6, replace = FALSE)

Train.X = X[,,,Train.sample]
dim(Train.X) = c(28, 28, 1, n.sample * 0.6)
Train.Y = Y[,Train.sample]

Test.X = X[,,,-Train.sample]
dim(Test.X) = c(28, 28, 1, n.sample * 0.4)
Test.Y = Y[,-Train.sample]

我們先把Training/Test的資料寫出，這樣下次我們就可以比較快的讀進來了：

fwrite(x = data.table(mnist[Train.sample,]),
       file = 'data/train_data.csv',
       col.names = FALSE, row.names = FALSE)

fwrite(x = data.table(mnist[-Train.sample,]),
       file = 'data/test_data.csv',
       col.names = FALSE, row.names = FALSE)

第六節：利用卷積神經網路做手寫數字辨識(4)

準備好資料後，讓我們定義模型結構：

# input
data <- mx.symbol.Variable('data')

# first conv
conv1 <- mx.symbol.Convolution(data=data, kernel=c(5,5), num_filter=10, name = 'conv1')
relu1 <- mx.symbol.Activation(data=conv1, act_type="relu")
pool1 <- mx.symbol.Pooling(data=relu1, pool_type="max",
                          kernel=c(2,2), stride=c(2,2))

# second conv
conv2 <- mx.symbol.Convolution(data=pool1, kernel=c(5,5), num_filter=20, name = 'conv2')
relu2 <- mx.symbol.Activation(data=conv2, act_type="relu")
pool2 <- mx.symbol.Pooling(data=relu2, pool_type="max",
                          kernel=c(2,2), stride=c(2,2))

# first fullc
flatten <- mx.symbol.Flatten(data=pool2)
fc1 <- mx.symbol.FullyConnected(data=flatten, num_hidden=150, name = 'fc1')
relu3 <- mx.symbol.Activation(data=fc1, act_type="relu")

# second fullc
fc2 <- mx.symbol.FullyConnected(data=relu3, num_hidden=10, name = 'fc2')

# Softmax
lenet <- mx.symbol.SoftmaxOutput(data = fc2, name = 'lenet')

這個網路結構如下：

原始圖片(28x28x1)要先經過10個5x5的「卷積器」(5x5x1x10)處理，將使圖片變成10張「一階特徵圖」(24x24x10)
接著這10張「一階特徵圖」(24x24x10)會經過ReLU，產生10張「轉換後的一階特徵圖」(24x24x10)
接著這10張「轉換後的一階特徵圖」(24x24x10)再經過2x2「池化器」(2x2)處理，將使圖片變成10張「降維後的一階特徵圖」(12x12x10)

– 第二層卷積組合

再將10張「降維後的一階特徵圖」(12x12x10)經過20個5x5的「卷積器」(5x5x10x20)處理，將使圖片變成20張「二階特徵圖」(8x8x20)
接著這20張「二階特徵圖」(8x8x20)會經過ReLU，產生20張「轉換後的二階特徵圖」(8x8x20)
接著這20張「轉換後的二階特徵圖」(8x8x20)再經過2x2「池化器」(2x2)處理，將使圖片變成20張「降維後的二階特徵圖」(4x4x20)

– 全連接層

將「降維後的二階特徵圖」(4x4x20)重新排列，壓製成「一階高級特徵」(320)
讓「一階高級特徵」(320)進入「隱藏層」，輸出「二階高級特徵」(150)
「二階高級特徵」(150)經過ReLU，輸出「轉換後的二階高級特徵」(150)
「轉換後的二階高級特徵」(150)進入「輸出層」，產生「原始輸出」(10)
「原始輸出」(10)經過Softmax函數轉換，判斷圖片是哪個類別

第六節：利用卷積神經網路做手寫數字辨識(5)

我們使用剛剛那個模型評估指標(同樣多是多分類任務)並開始訓練(這個要跑一陣子)：

my.eval.metric.mlogloss <- mx.metric.custom(
  name = "m-logloss", 
  function(real, pred) {
    real1 = as.numeric(as.array(real))
    pred1 = as.numeric(as.array(pred))
    pred1[pred1 <= 1e-6] = 1e-6
    pred1[pred1 >= 1 - 1e-6] = 1 - 1e-6
    return(-mean(real1 * log(pred1), na.rm = TRUE))
  }
)

mx.set.seed(0)

lenet_model = mx.model.FeedForward.create(symbol = lenet,
                                          X = Train.X, y = Train.Y,
                                          optimizer = "sgd", learning.rate = 0.05, momentum = 0.9,
                                          array.batch.size = 100, num.round = 20,
                                          ctx = mx.gpu(1),
                                          eval.metric = my.eval.metric.mlogloss)

看看預測準確度：

predict_Y = predict(lenet_model, Test.X)
predict_cat = max.col(t(predict_Y))
real_cat = max.col(t(Test.Y))

confusion_table = table(predict_cat, real_cat)
cat("Testing accuracy rate =", sum(diag(confusion_table))/sum(confusion_table))

## Testing accuracy rate = 0.983869

print(confusion_table)

##            real_cat
## predict_cat    1    2    3    4    5    6    7    8    9   10
##          1  1656    0    4    0    3    3    7    0    4   10
##          2     0 1833    3    1    3    0    1    5    4    5
##          3     1    3 1634    4    1    0    3   10    9    1
##          4     0    0    6 1716    0    2    1    3    3    8
##          5     0    8    5    0 1588    1    8    1    2   14
##          6     1    1    1   14    0 1536    6    0   16   10
##          7     5    1    0    0    0    1 1632    0    5    0
##          8     0    4    3    2    3    3    0 1731    0   16
##          9     0    1    0    4    1    4    3    1 1631    6
##          10    0    0    0    1    7    1    0    2    1 1572

你是否覺得非常驚訝，準確度居然有98%以上！而這就是卷積神經網路的威力。

– 我們可以把模型儲存/讀取模型：

#Save model
mx.model.save(lenet_model, "model/lenet", iteration = 0)

#Load model
lenet_model = mx.model.load("model/lenet", iteration = 0)

第七節：利用經典模型進行圖像預測(1)

在深度學習領域，演算法只能算是引擎，而數據可以說是石油。沒有石油的加入，引擎是沒辦法體現它的性能的！
即使到了現在，一個單純蒐集數據的研究案仍然不太可能取得計劃支持，但這項重要的工作總得有人做。史丹佛大學的李飛飛從2005年開始堅持下來努力做這件事情，並於2007年創辦ImageNet，他收集大量帶有標註信息的圖片數據供電腦視覺模型訓練，而每年將會利用其資料進行ILSVRC圖像識別競賽。

– 她在2009年CVPR所釋出的資料可以說是這一波深度學習突破中最重要的基石！

F17

– 我們在這裡提供了一個已經預訓練好的模型列表

第七節：利用經典模型進行圖像預測(2)

讓我們隨便上網找一張圖，使用
有了這個東西後，我們就能隨便找張圖片進行測試，先讓我們把圖片顯示出來：

library(OpenImageR)

img<- readImage('test.jpg') 
resized_img <- resizeImage(img, 224, 224, method = 'bilinear')

imageShow(resized_img)

第七節：利用經典模型進行圖像預測(3)

讓我們讀取模型進行預測吧：

library(mxnet)

res_model <- mx.model.load("model/resnet-50", 0)

dim(resized_img) <- c(dim(resized_img), 1)
pred_prob  <- predict(res_model, resized_img)
which.max(pred_prob)

## [1] 3

模型認為這張圖是最有可能是第3個類別，讓我們把標籤讀取進來做比對，並顯示出前5個預測的結果：

synsets <- readLines('model/chinese synset.txt', encoding = 'UTF-8')

pred_prob <- as.numeric(pred_prob)
names(pred_prob) <- synsets
pred_prob <- sort(pred_prob, decreasing = TRUE)
pred_prob <- formatC(pred_prob, 4, format = 'f')

head(pred_prob, 5)

##             n01484850 大白鯊               n01491361 虎鯊 
##                     "0.9971"                     "0.0027" 
##             n01494475 鎚頭鯊 n02071294 殺人鯨,逆戟鯨,虎鯨 
##                     "0.0001"                     "0.0001" 
##               n02066245 灰鯨 
##                     "0.0000"

你覺得準嗎，上網隨便找張圖測試他的能力！

第八節：訓練一個圖像分類模型並分享他的預測功能(1)

請你點擊這裡下載一個我們已經寫好的網頁應用程式

– 稍微試一下他的能力，你可以隨便下載一張圖片，這跟剛剛的一樣

在MxNet裡面要基於已經訓練完的模型做新的任務也非常簡單，我們可以用這個程式碼抓取倒數第二層的Layer，而抓到這個Layer後我們就能隨意的添加我們想要的結構。

– 舉例來說，我們可以定義我們要整個resnet-50除了最後一個全連接層外的所有結構，只把最後一層的FC從分1000類轉變成分2類：

library(magrittr)
library(mxnet)

res_model <- mx.model.load("model/resnet-50", 0)

all_layers <- res_model$symbol$get.internals()
flatten0_output <- which(all_layers$outputs == 'flatten0_output') %>% all_layers$get.output()

fc1 <- mx.symbol.FullyConnected(data = flatten0_output, num_hidden = 2, name = 'fc1')
softmax <- mx.symbol.SoftmaxOutput(data = fc1, name = 'softmax')

這樣就定義好了新的結構，至於你自己想做的任務是幾個類別就自行決定！

第八節：訓練一個圖像分類模型並分享他的預測功能(2)

這裡我們要把這個模型拿來做貓狗分類任務，我們需要用到第一節課教到的技能讀取貓狗資料：

Train_img <- array(0, dim = c(224, 224, 3, 200))
Train.y <- array(0, dim = c(2, 200))

for (i in 1:100) {
  
  # Cat
  img <- readImage(paste0('Dogs vs. Cats/cat.', i, '.jpg'))
  resized_img <- resizeImage(img, 224, 224, method = 'bilinear')
  Train_img[,,,2*i-1] <- resized_img
  Train.y[1,2*i-1] <- 1
  
  # Dog
  img <- readImage(paste0('Dogs vs. Cats/dog.', i, '.jpg'))
  resized_img <- resizeImage(img, 224, 224, method = 'bilinear')
  Train_img[,,,2*i] <- resized_img
  Train.y[2,2*i] <- 1
  
}

第八節：訓練一個圖像分類模型並分享他的預測功能(3)

– 接著，我們在開始訓練之前需要取得模型權重參數，我們可以將最後一層以外的部分填入resnet-50的參數，並以這為基礎開始訓練任務：

new_arg <- mxnet:::mx.model.init.params(symbol = softmax,
                                        input.shape = list(data = c(224, 224, 3, 32)),
                                        output.shape = NULL,
                                        initializer = mxnet:::mx.init.uniform(0.01),
                                        ctx = mx.cpu())

for (i in 1:length(new_arg$arg.params)) {
  pos <- which(names(res_model$arg.params) == names(new_arg$arg.params)[i])
  if (length(pos) == 1) {
    if (all.equal(dim(res_model$arg.params[[pos]]), dim(new_arg$arg.params[[i]])) == TRUE) {
      new_arg$arg.params[[i]] <- res_model$arg.params[[pos]]
    }
  }
}

for (i in 1:length(new_arg$aux.params)) {
  pos <- which(names(res_model$aux.params) == names(new_arg$aux.params)[i])
  if (length(pos) == 1) {
    if (all.equal(dim(res_model$aux.params[[pos]]), dim(new_arg$aux.params[[i]])) == TRUE) {
      new_arg$aux.params[[i]] <- res_model$aux.params[[pos]]
    }
  }
}

第八節：訓練一個圖像分類模型並分享他的預測功能(4)

可以開始訓練了：

my.eval.metric.mlogloss <- mx.metric.custom(
  name = "m-logloss", 
  function(real, pred) {
    real1 = as.numeric(as.array(real))
    pred1 = as.numeric(as.array(pred))
    pred1[pred1 <= 1e-6] = 1e-6
    pred1[pred1 >= 1 - 1e-6] = 1 - 1e-6
    return(-mean(real1 * log(pred1), na.rm = TRUE))
  }
)

mx.set.seed(0)

my_model <- mx.model.FeedForward.create(symbol = softmax,
                                       X = Train_img, y = Train.y,
                                       optimizer = "sgd", learning.rate = 0.001, momentum = 0.9,
                                       array.batch.size = 20, num.round = 20,
                                       arg.params = new_arg$arg.params, aux.params = new_arg$aux.params,
                                       ctx = mx.gpu(1),
                                       eval.metric = my.eval.metric.mlogloss)

第八節：訓練一個圖像分類模型並分享他的預測功能(5)

讓我們讀張新的圖看看這個模型準不準：

img <- readImage('Dogs vs. Cats/test_cat.3.jpg')
resized_img <- resizeImage(img, 224, 224, method = 'bilinear')

imageShow(resized_img)

看看模型的預測吧：

dim(resized_img) = c(224, 224, 3, 1)

pred_prob <- predict(my_model, resized_img)
pred_prob

##             [,1]
## [1,] 0.996811688
## [2,] 0.003188348

第八節：訓練一個圖像分類模型並分享他的預測功能(6)

再看看狗狗的預測結果吧：

img <- readImage('Dogs vs. Cats/test_dog.3.jpg')
resized_img <- resizeImage(img, 224, 224, method = 'bilinear')

imageShow(resized_img)

看看模型的預測吧：

dim(resized_img) <- c(224, 224, 3, 1)

pred_prob <- predict(my_model, resized_img)
pred_prob

##            [,1]
## [1,] 0.05498257
## [2,] 0.94501746

換多點貓狗圖片試試看吧！

結語

訓練自己的深度學習模型好玩嗎?你要不要試著把剛剛的貓狗識別模型保留下來，並且把他放入剛剛那個App內，看看他是否work。

– 如果你想把自己訓練的模型放到剛剛的App內，還要記得改變chinese synset.txt檔案喔！

今天的課程主要教授大家如何簡單的做出圖片分類任務，但其中有非常多細節都沒有提到，事實上要成功訓練一個深度學習模型是很難的，各位如果有興趣的話可能還需要再更深入的學習！

– 另外，深度學習的潛力可不僅僅只有圖像辨識任務，其他包含物件識別、物件分割、圖像生成，甚至是語言模型，都是可以做出來的，這也值得有興趣的人再進一步學習！

深度學期快速入門

第一節：基本知識(1)