之前介紹過將值放在參數型別上。

今天來介紹一下如何可以做到類似運算的效果。

1
2
3
4
5

struct A{T}
end

a1 = A{5}()
a2 = A{3}()

在一些應用場景上會希望將參數欄位上的值做運算，例如加總。

這時候我們可以這樣做：

1
2
3
4
5
6

import Base:+

function +(::A{T}, ::A{S}) where {T, S}
	x = T + S
	return A{x}()
end

如此一來，就可以簡單搞定囉！

1

println(a1 + a2)

1

A{8}()

留言與分享

這篇被選為 NeurIPS 2018 最佳論文，他將連續的概念帶入了神經網路架構中，並且善用以往解微分方程的方法來做逼近，可以做到跟原方法（倒傳遞）一樣好的程度，而參數使用複雜度卻是常數，更短的訓練時間。

核心觀念

概念上來說，就是將神經網路離散的層觀念打破，將他貫通成為連續的層的網路架構。

連續和離散的差別來自於倒傳遞的過程：

$$
\mathbb{y}_{t+1} = \mathbb{y}_t - \eta \nabla \mathcal{L}
$$

其中 $\nabla \mathcal{L}$ 就是梯度的部份，是向量的，然而我們把他簡化成純量來看的話，他不過就是

$$
\frac{d \mathcal{L}}{dt}
$$

廣義上來說，一個函數的微分，如果是離散的版本就會是

$$
\frac{dy}{dt} = \frac{y(t + \Delta) - y(t)}{\Delta}
$$

如此一來，所形成的方程式就會是差分方程，然而連續的版本就是

$$
\frac{dy}{dt} = \lim_{\Delta \rightarrow 0} \frac{y(t + \Delta) - y(t)}{\Delta}
$$

這個所形成的會是微分方程。

從離散到連續

我們可以從離散的版本

$$
\frac{dy}{dt} = \frac{y(t + \Delta) - y(t)}{\Delta}
$$

把他轉成以下的樣貌

$$
y(t + \Delta) = y(t) + \Delta \frac{dy}{dt}
$$

要將他貫通的話，我們就得由從神經網路的基礎開始，如果是一般的前回饋網路（feed-forward network）當中的隱藏層是像下列這個樣子：

$$
h_{t+1} = f(h_t, \theta)
$$

我們可以發現像是 ResNet 這類的網路有 skip connection 的設置，所以跟一般的前回饋網路不同

$$
h_{t+1} = h_t + f(h_t, \theta)
$$

而 RNN 等等有序列概念的模型也有類似的結構，就是會是前一層的結果加上通過 $f$ 運算後的結果，成為下一層的結果。

這樣的形式跟我們前面提到的形式不謀而合

$$
y(t + \Delta) = y(t) + \Delta \frac{dy}{dt}
$$

只要我們把 $\Delta = 1$ 代入，就成了

$$
y(t+1) = y(t) + \frac{dy}{dt}
$$

以下給大家比對一下

$$
h_{t+1} = h_t + f(h_t, \theta) \\
y(t+1) = y(t) + \frac{dy}{dt}
$$

也就是，我們可以讓

$$
\frac{dy}{dt} = f(h_t, \theta)
$$

神奇的事情就發生了！神經網路 $f$ 就可以被我們拿來計算微分 $\frac{dy}{dt}$！

比較精確的說法是，把神經網路的層 $f$ 拿來逼近微分項，或是說梯度。這樣我們後面就可以用數值方法來逼近解。

$$
y(t + \Delta) = y(t) + \Delta \frac{dy}{dt} \\
\downarrow \\
y(t + \Delta) = y(t) + \Delta f(t, h(t), \theta_t)
$$

要拉成連續的還有一個重要的手段，就是將不同的層 $t$ 從離散的變成連續的，所以作者將 $t$ 做了參數化，將他變成 $f$ 的參數之一，如此一來，就可以在任意的層中放入資料做運算。

最重要的概念導出了這樣的式子

$$
h(t) \rightarrow \frac{dy(t)}{dt} = f(h(t), t, \theta) \rightarrow y(t)
$$

神經網路作為一個系統的微分形式

在傳統科學或是工程領域，我們會以微分式來表達以及建構一個系統。

$$
\nu = \frac{dx}{dt} = t + 1
$$

其實在這邊是一樣的道理，整體來說，我們是換成用神經網路去描述一個微分式，其實本質上就是這樣。

原本的層的概念就是用數學函數來建立的，而層與層之間傳遞著計算的結果。

$$
\mathbb{h_1} = \sigma(W_1 \mathbb{x} + \mathbb{b_1}) \\
\mathbb{y} = \sigma(W_2 \mathbb{h_1} + \mathbb{b_2})
$$

然而變成連續之後，我們等於是用神經網路中的層去建立跟描繪微分形式。

$$
\frac{d h(t)}{dt} = \sigma(W(t) \mathbb{x}(t) + \mathbb{b(t)}) \\
\frac{d y(t)}{dt} = \sigma(W(t) \mathbb{h}(t) + \mathbb{b(t)})
$$

是不是跟如出一轍呢？

$$
\frac{dy(t)}{dt} = f(h(t), t, \theta)
$$

向前傳遞解微分式

我們可以來計算看看隱藏層是長什麼樣子的。在隱藏層的微分式中，也是利用隱藏層去計算出來的。

$$
\frac{dh(t)}{dt} = f(h(t), t, \theta)
$$

基本上，我們只要對上式做積分就可以了。

$$
h(t) = \int f(h(t), t, \theta) dt
$$

這是一個怎樣的概念呢？我們可以來看看下圖。

我們做積分這件事其實是用 $h(t_0)$ 來推斷 $h(t_1)$ 的，這跟神經網路的向前傳遞是一樣的行為。

$$
h(t_1) = F(h(t), t, \theta) \bigg|_{t=t_0}
$$

這樣的積分動作，我們可以用 $t_0$ 時間點的資訊來解 $h(t_1)$。

這樣的解法在程式上就會交由 ODE Solver 去處理。

$$
h(t_1) = ODESolve(h(t_0), t_0, t_1, \theta, f)
$$

反向傳遞解函數

$$
\mathcal{L}(t_0, t, \theta) = \mathcal{L}(ODESolve(\cdot))
$$

$$
\frac{\partial \mathcal{L}}{\partial h(t)} = -a(t)
$$

adjoint state

$$
a(t) = \int -a(t)^T \frac{\partial f}{\partial h} dt = - \frac{\partial \mathcal{L}}{\partial h(t)}
$$

$$
a(t) = \int_{t_1}^{t_0} -a(t)^T \frac{\partial f(h(t), t, \theta)}{\partial h(t)} dt
$$

擴充狀態（augmented state）

$\frac{d \theta}{dt} = 0$

$\frac{dt}{dt} = 1$

let $\begin{bmatrix}
h \\
\theta \\
t
\end{bmatrix}$ be a augmented state

augmented state function:

$$
f_{aug}(\begin{bmatrix}
h \\
\theta \\
t
\end{bmatrix}) =
\begin{bmatrix}
f(h(t), t, \theta) \\
0 \\
1
\end{bmatrix}
$$

augmented state dynamics:

$$
\frac{d}{dt}
\begin{bmatrix}
h \\
\theta \\
t
\end{bmatrix}

f_{aug}(
\begin{bmatrix}
h \\
\theta \\
t
\end{bmatrix})
$$

augmented adjoint state:

$$
\begin{bmatrix}
a \\
a_{\theta} \\
a_t
\end{bmatrix}
$$

$a = \frac{\partial \mathcal{L}}{\partial h}$

$a_{\theta} = \frac{\partial \mathcal{L}}{\partial \theta}$

$a_t = \frac{\partial \mathcal{L}}{\partial t}$

$$
\frac{d a_{aug}}{dt} = -
\begin{bmatrix}
a \frac{\partial f}{\partial h} \\
a \frac{\partial f}{\partial \theta} \\
a \frac{\partial f}{\partial t}
\end{bmatrix}
$$

留言與分享

參考從幾何向量空間到函數空間| 線代啟示錄。

1. 由 $\mathbb{R}^n$ 拓展到 $\mathbb{R}^{\infty}$ 所需俱備的條件是什麼？

由於一個向量 $\mathbb{v} \in \mathbb{R}^{\infty}$，在無限維度下我們需要考慮一個問題，就是 norm。

如果這個空間有定義 norm 的話，我們就要考慮他有沒有收斂，也就是 $||\mathbb{v}||^2$ 要存在。

所以條件就是

$$
||\mathbb{v}||^2 = \sum_{i=1}^{\infty} v_i^2
$$

要收斂。

2. 從 $\mathbb{R}^{\infty}$ 無限維度的向量空間再拓展到 $C^{\omega}$ 函數空間，所需要俱備的條件是什麼？

一個無限維度的向量是一個離散的版本，由剛剛的式子可以看的出來

$$
||\mathbb{v}||^2 = \sum_{i=1}^{\infty} v_i^2
$$

而一個（解析）函數則是連續的

$$
||f||^2 = \int f^2(x) dx
$$

除了以上的 norm 要收斂外，從離散到連續應該有些假設或是條件才是。

3. 函數的基底

Fourier series

$$
f(x) = a_0 + a_1 \cos x + b_1 \sin x + a_2 \cos 2x + b_2 \sin 2x + \cdots
$$

所以基底就是

$$
<\beta> = <1, \cos x, \sin x, \cos 2x, \sin 2x, \cdots>
$$

4. 非週期性函數基底

Legendre polynomial

5. Least square problem

$$
(A^TA)\hat{y} = A^Tb
$$

留言與分享

我在使用的時候有注意到在參數化型別中使用值的方式與傳統封裝的方式有效能上的差異。

所以我就做了一些測試。

在參數化型別中使用值：

1
2
3

struct A{T}

end

傳統型別封裝：

1
2
3

struct B
    x::Int64
end

全文出現的程式碼為實際測試程式碼

因為 Julia 有提供好用的 @code_llvm 及 @code_native 來觀察一行程式碼實際被轉換成 LLVM 或是組合語言的時候會產生多少行的程式碼，藉此我們可以用低階程式碼來評估是否有效率。程式碼的行數愈少是越有效率的。

建立

我們來測試一個物件被建立需要多少行的程式碼。

A - LLVM

1

@code_llvm A{5}()

1
2
3
4
5
6
7

;  @ REPL[1]:3 within `Type'
define nonnull %jl_value_t addrspace(10)* @japi1_Type_12238(%jl_value_t addrspace(10)*, %jl_value_t addrspace(10)**, i32) #0 {
top:
  %3 = alloca %jl_value_t addrspace(10)**, align 8
  store volatile %jl_value_t addrspace(10)** %1, %jl_value_t addrspace(10)*** %3, align 8
  ret %jl_value_t addrspace(10)* addrspacecast (%jl_value_t* inttoptr (i64 140407726014496 to %jl_value_t*) to %jl_value_t addrspace(10)*)
}

B - LLVM

1

@code_llvm B(5)

1
2
3
4
5
6

;  @ REPL[1]:2 within `Type'
define { i64 } @julia_Type_12221(%jl_value_t addrspace(10)*, i64) {
top:
  %.fca.0.insert = insertvalue { i64 } undef, i64 %1, 0
  ret { i64 } %.fca.0.insert
}

A - Assembly

1

@code_native A{5}()

1
2
3
4
5

; ┌ @ REPL[1]:3 within `Type'
	movq	%rsi, -8(%rsp)
	movabsq	$140407726014496, %rax  # imm = 0x7FB338A20020
	retq
; └

B - Assembly

1

@code_native B(5)

1
2
3
4
5

; ┌ @ REPL[1]:2 within `Type'
	movq	%rsi, %rax
	retq
	nopw	%cs:(%rax,%rax)
; └

取值

接著測試從物件當中取值出來的效能。

定義取值的方法：

1
2

get_value(::A{T}) where {T} = T
get_value(b::B) = b.x

事先建立好物件：

1
2

a = A{5}()
b = B(5)

A - LLVM

1

@code_llvm get_value(a)

1
2
3
4
5

;  @ REPL[8]:1 within `get_value'
define i64 @julia_get_value_12274() {
top:
  ret i64 5
}

B - LLVM

1

@code_llvm get_value(b)

1
2
3
4
5
6
7
8
9

;  @ REPL[5]:1 within `get_value'
define i64 @julia_get_value_12630({ i64 } addrspace(11)* nocapture nonnull readonly dereferenceable(8)) {
top:
; ┌ @ sysimg.jl:18 within `getproperty'
   %1 = getelementptr inbounds { i64 }, { i64 } addrspace(11)* %0, i64 0, i32 0
; └
  %2 = load i64, i64 addrspace(11)* %1, align 8
  ret i64 %2
}

A - Assembly

1

@code_native get_value(a)

1
2
3
4
5

; ┌ @ REPL[8]:1 within `get_value'
	movl	$5, %eax
	retq
	nopw	%cs:(%rax,%rax)
; └

B - Assembly

1

@code_native get_value(b)

1
2
3
4
5

; ┌ @ REPL[5]:1 within `get_value'
	movq	(%rdi), %rax
	retq
	nopw	%cs:(%rax,%rax)
; └

給大家參考。

留言與分享

應該不少人看到這個標題會摸不著頭緒到底要做什麼，但是看完 Julia 中常見的程式碼你就會了解了。

1

Array{Any, 2}

有沒有曾經納悶過那個數字 2 到底是怎麼進到參數的位置上的呢？

參數的位置不是只能放型別（type）嗎？

這同時也是我困惑已久的問題，就搜尋了一下，果不其然被我找到了方法：

1
2

struct A{T}
end

1

A{5}()

原來這麼簡單就可以完成了！語法上並沒有限定一定要是型別，要放型別以外的東西似乎是可以的。

我目前測試了可以的有：Int64、Float64、Complex、Char、Bool、Symbol，所以估計數字應該都是可以的。

不行的有：String、Array，估計物件或是陣列都是不行的。

定義範圍

不過使用上並沒有任何限制會有點危險，所以還是定義一下範圍會比較好，像是：

1
2
3
4
5
6

struct A{I}
    function A{I}() where {I}
        isa(I,Integer) || error("bad parameter")
        new()
    end
end

這樣就可以限制參數要是整數的範圍。

從參數取值

那我們能不能從型別的參數當中取值呢？

可以。

1

get_value(::A{I}) where A{I} = I

如此一來，我們就可以從型別中拿到值了。

好處？

這麼做有什麼好處？

當你把值的資訊放到型別當中，型別就多了一些資訊可以提供編譯器處理，這對於要自己設計型別階層可是非常好用的。

例如像是你可以將陣列的長度資訊儲存到型別上，這樣編譯器就可以處理陣列的長度資訊了。

這樣的程式風格會跟 dependent type language 有些相似了。

大家可以玩玩看。

留言與分享

20th century Cantor:

All mathematical objects can be defined as sets.

Fundamentals:

• numbers
• permutations
• partitions
• matrices
• sets
• functions
• relations

Geometry:

• hexagons
• points
• lines
• triangles
• circles
• spheres
• polyhedra
• topological space
• manifolds

Algebra:

• groups
• rings
• fields
• lattices

留言與分享

受到其他文章的啟發，我也來寫一篇為什麼我踏入生物資訊領域好了。

受到啟發應該算是從高中的時候說起，高中的時候喜歡數學、物理跟生物。對於數學，喜歡他的抽象及純粹，而物理可以解釋這個世界的法則，對於生物則是一直以來隱隱約約有些感覺的。小時候對於生命現象一直很好奇，對於生物的多樣性感到驚奇，但到了高中卻成了考卷上的考題，我不認為那是我要的。

還記得高中生物上到下視丘的時候會講到很多不同種的激素調控，我突然覺得這一切的背後似乎有著什麼，我對「調控」產生了興趣。接著到了高中快結束，終於上到近代的生物技術以及 DNA 分子的轉錄轉譯，雖然對當時的我來說有點複雜，但是我喜歡挑戰理解這種複雜的事物，我將他轉化成比較好理解的「設計圖」解釋。DNA 就像是一台車子的整體設計圖，RNA 就是將設計圖的一部份零件複製一份出來，並且製造出蛋白質，也就是真實的零件。理解了這些讓我非常開心。

大學念了醫學檢驗生物技術，但卻不是我的第一志願，不過我確定我對生物技術是有興趣的，我也非常認真對待我的選擇。在傳統的生物醫學研究都是花了十幾年的時間在研究一個蛋白或是一個基因的功能或是交互作用。

我大三的某天在逛維基百科（你沒看錯，我會去逛維基百科）被我發現了系統生物學這個領域，看到頁面的當下非常震驚，可以以一個系統的觀點切入生物的議題，那麼就可以不用那麼辛苦的一個基因一個基因研究了。而且系統的概念直接串起了在生化中學到的調控，他不只是 pathway，而是一個複雜的網路，可以藉由網路的調控或是反應機制，讓生物體做出特定的行為。生物體就是個巨大的機械，但是複雜度卻遠高於機械，也不像機構那樣那麼容易理解，很多事情是人類目前還不知道的。

因為這樣的未知，這樣的複雜，這樣的調控系統，讓我決心研究所要往這個方向走。

大四的時候有進階生物技術，接觸到定序技術、生物資訊、序列處理的議題。同時雙主修資訊工程，我更享受在資工系的課程當中，雖然他講的是程式、作業系統等等，但是對於（建造）系統的概念始終是保留的。我最有興趣的大概是離散數學、演算法跟訊號與系統了，離散數學中的圖論可以說是非常神妙，而圖論的用途也超級廣，可以拿來 model 很多不同的事物。演算法則是去證明一件事情可以被如何的完成是最快的，這些魔法都來自於數學。訊號與系統講述了如何去探知或是解析一個系統，我們怎麼從一個系統的行為當中去反推這個系統的架構。

到研究所真正接觸了生物資訊與我的認知相去不遠，不過還是少了點什麼，看了看課程發現了機器學習的課程，也詢問了學長關於這個領域，聽說還蠻推薦的，但是受限於開課時間，就乾脆自己去找了 coursera 上林軒田老師的機器學習課程看，大概一個月左右就把他看完了。看的當下非常開心，學到了跟演算法非常相似的技術，而當時大數據剛開始紅，所以就以這個技術為主軸開始了我的研究。

殊不知，後來的深度學習的崛起，AI 的爆紅，讓機器學習變得異常的熱門。不過我還是希望繼續做系統生物學，應該說是計算生物學。來把這迷樣的生物系統 model 出來吧！

留言與分享