NDArray介紹

機器學習處理的物件是資料，資料一般是由外部感測器（sensors）採集，經過數字化後儲存在計算機中，可能是文字、聲音，圖片、視訊等不同形式。
這些數字化的資料最終會載入到記憶體進行各種清洗，運算操作。
幾乎所有的機器學習演算法都涉及到對資料的各種數學運算，比如：加減、點乘、矩陣乘等。所以我們需要一個易用的、高效的、功能強大的工具來處理這些資料並組支援各種複雜的數學運算。

在C/C++中已經開發出來了很多高效的針對於向量、矩陣的運算庫，比如：OpenBLAS,Altlas,MKL等。

對於Python來說Numpy無疑是一個強大針對資料科學的工具包，它提供了一個強大的高維資料的陣列表示，以及支援Broadcasting的運算，並提供了線性代數、傅立葉變換、隨機數等功能強大的函式。

MXNet的NDArray與Numpy中的ndarray極為相似，NDAarray為MXNet中的各種數學計算提供了核心的資料結構，NDArray表示一個多維的、固定大小的陣列，並且支援異構計算。那為什麼不直接使用Numpy呢？MXNet的NDArray提供額外提供了兩個好處：

• 支援異構計算，資料可以在CPU,GPU，以及多GPU機器的硬體環境下高效的運算
• NDArray支援惰性求值，對於複雜的操作，可以在有多個計算單元的裝置上自動的並行運算。

NDArray的重要屬性

每個NDarray都具有以下重要的屬性，我們可以通過相應的api來訪問：

• ndarray.shape：陣列的維度。它返回了一個整數的元組，元組的長度等於陣列的維數，元組的每個元素對應了陣列在該維度上的長度。比如對於一個n行m列的矩陣，那麼它的形狀就是(n,m)。
• ndarray.dtype：陣列中所有元素的型別，它返回的是一個numpy.dtype的型別，它可以是int32/float32/float64等，預設是'float32'的。
• ndarray.size：陣列中元素的個數，它等於ndarray.shape的所有元素的乘積。
• ndarray.context：陣列的儲存裝置，比如：cpu()或gpu(1)

import mxnet as mx
import mxnet.ndarray as nd

a = nd.ones(shape=(2,3),dtype='int32',ctx=mx.gpu(1))
print(a.shape, a.dtype, a.size, a.context)
 

NDArray的建立

一般來常見有2種方法來建立NDarray陣列：

1. 使用ndarray.array直接將一個list或numpy.ndarray轉換為一個NDArray
2. 使用一些內建的函式zeros,ones以及一些隨機數模組ndarray.random建立NDArray，並預填充了一些資料。
3. 從一個一維的NDArray進行reshape

import numpy as np

l = [[1,2],[3,4]]
print(nd.array(l)) # 從List轉到NDArray
print(nd.array(np.array(l))) # 從np.array轉到NDArray

# 直接利用函式建立指定大小的NDArray
print (nd.zeros((3,4), dtype='float32'))
print (nd.ones((3,4), ctx=mx.gpu()))
# 從一個正態分佈的隨機數引擎生成了一個指定大小的NDArray，我們還可以指定分佈的引數，比如均值，標準差等
print (nd.random.normal(shape=(3,4))) 
print (nd.arange(18).reshape(3,2,3))

NDArray的檢視

一般情況下，我們可以通過直接使用print來檢視NDArray中的內容，我們也可以使用nd.asnumpy()函式，將一個NDArray轉換為一個numpy.ndarray來檢視。

a = nd.random.normal(0, 2, shape=(3,3))
print(a)
print(a.asnumpy())

基本的數學運算

NDArray之間可以進行加減乘除等一系列的數學運算，其中大部分的運算都是逐元素進行的。

shape=(3,4)
x = nd.ones(shape)
y = nd.random_normal(0, 1, shape=shape)
x + y # 逐元素相加
x * y # 逐元素相乘
nd.exp(y) # 每個元素取指數
nd.sin(y**2).T # 對y逐元素求平方，然後求sin，最後對整個NDArray轉置
nd.maximum(x,y) # x與y逐元素求最大值

這裡需要注意的是*運算是兩個NDArray之間逐元素的乘法，要進行矩陣乘法，必須使用ndarray.dot函式進行矩陣乘

nd.dot(x, y.T)

索引與切片

MXNet NDArray提供了各種擷取的方法，其用法與Python中list的擷取操作以及Numpy.ndarray中的擷取操作基本一致。

x = nd.arange(0, 9).reshape((3,3))
x[1:3] # 擷取x的axis=0的第1和第2行
x[1:2,1:3] # 擷取x的axis=0的第1行，axis=1的第一行和第二行

儲存變化

在對NDArray進行演算法運算時，每個操作都會開闢新的記憶體來儲存運算的結果。例如：如果我們寫y = x + y，我們會把y從現在指向的例項轉到新建立的例項上去。我們可以把上面的運算看成兩步：z = x + y; y = z。

我們可以使用python的內建函式id()來驗證。id()返回一個物件的識別符號，當這個物件存在時，這個識別符號一定是惟一的，在CPython中這個識別符號實際上就是物件的地址。

x = nd.ones((3,4))
y = nd.ones((3,4))
before = id(y)
y = x + y
print(before, id(y))

在很多情況下，我們希望能夠在原地對陣列進行運算，那麼我們可以使用下面的一些語句：

y += x
print(id(y))

nd.elemwise_add(x, y, out=y)
print(id(y))

y[:] = x + y
print(id(y))

在NDArray中一般的賦值語句像y = x，y實際上只是x的一個別名而已，x和y是共享一份資料儲存空間的

x = nd.ones((2,2))
y = x
print(id(x))
print(id(y))

如果我們想得到一份x的真實拷貝，我們可以使用copy函式

y = x.copy()
print(id(y))

Broadcasting

廣播是一種強有力的機制，可以讓不同大小的NDArray在一起進行數學計算。我們常常會有一個小的矩陣和一個大的矩陣，然後我們會需要用小的矩陣對大的矩陣做一些計算。

舉個例子，如果我們想要把一個向量加到矩陣的每一行，我們可以這樣做

# 將v加到x的每一行中，並將結果儲存在y中
x = nd.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
v = nd.array([1, 0, 1])
y = nd.zeros_like(x)   # Create an empty matrix with the same shape as x

for i in range(4):
    y[i, :] = x[i, :] + v
print (y)

這樣是行得通的，但是當x矩陣非常大，利用迴圈來計算就會變得很慢很慢。我們可以換一種思路：

x = nd.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
v = nd.array([1, 0, 1])
vv = nd.tile(v, (4, 1))  # Stack 4 copies of v on top of each other
y = x + vv  # Add x and vv elementwise
print (y)
# 也可以通過broadcast_to來實現
vv = v.broadcast_to((4,3))
print(vv)

NDArray的廣播機制使得我們不用像上面那樣先建立vv，可以直接進行運算

x = nd.array([[1,2,3], [4,5,6], [7,8,9], [10, 11, 12]])
v = nd.array([1, 0, 1])
y = x + v
print(y)

對兩個陣列使用廣播機制要遵守下列規則：

1. 如果陣列的秩不同，使用1來將秩較小的陣列進行擴充套件，直到兩個陣列的尺寸的長度都一樣。
2. 如果兩個陣列在某個維度上的長度是一樣的，或者其中一個數組在該維度上長度為1，那麼我們就說這兩個陣列在該維度上是相容的。
3. 如果兩個陣列在所有維度上都是相容的，他們就能使用廣播。
4. 如果兩個輸入陣列的尺寸不同，那麼注意其中較大的那個尺寸。因為廣播之後，兩個陣列的尺寸將和那個較大的尺寸一樣。
5. 在任何一個維度上，如果一個數組的長度為1，另一個數組長度大於1，那麼在該維度上，就好像是對第一個陣列進行了複製。

在GPU上運算

NDArray支援陣列在GPU裝置上運算，這是MXNet NDArray和Numpy的ndarray最大的不同。預設情況下NDArray的所有操作都是在CPU上執行的，我們可以通過ndarray.context來查詢陣列所在裝置。在有GPU支援的環境上，我們可以指定NDArray在gpu裝置上。

gpu_device = mx.gpu(0)
def f():
    a = mx.nd.ones((100,100))
    b = mx.nd.ones((100,100), ctx=mx.cpu())
    c = a + b.as_in_context(a.context)
    print(c)

f() # 在CPU上運算

# 在GPU上運算
with mx.Context(gpu_device):
    f()

上面語句中使用了with來構造了一個gpu環境的上下文，在上下文中的所有語句，如果沒有顯式的指定context，則會使用wtih語句指定的context。
當前版本的NDArray要求進行相互運算的陣列的context必須一致。我們可以使用as_in_context來進行NDArray context的切換。

NDArray的序列化

有兩種方法可以對NDArray物件進行序列化後儲存在磁碟，第一種方法是使用pickle，就像我們序列化其他python物件一樣。

import pickle

a = nd.ones((2,3))
data = pickle.dumps(a) # 將NDArray直接序列化為記憶體中的bytes
b = pickle.loads(data) # 從記憶體中的bytes反序列化為NDArray

pickle.dump(a, open('tmp.pickle', 'wb')) # 將NDArray直接序列化為檔案
b = pickle.load(open('tmp.pickle', 'rb')) # 從檔案反序列化為NDArray

在NDArray模組中，提供了更優秀的介面用於陣列與磁碟檔案（分散式儲存系統）之間進行資料轉換

a = mx.nd.ones((2,3))
b = mx.nd.ones((5,6))
nd.save("temp.ndarray", [a, b]) # 寫入與讀取的路徑支援Amzzon S3以及Hadoop HDFS等。
c = nd.load("temp.ndarray")

惰性求值與自動並行化

MXNet使用了惰性求值來追求最佳的效能。當我們在Python中執行a = b + 1時，Python執行緒只是將運算Push到了後端的執行引擎，然後就返回了。這樣做有下面兩個好處：

1. 當操作被push到後端後，Python的主執行緒可以繼續執行下面的語句，這對於Python這樣的解釋性的語言在執行計算型任務時特別有幫助。
2. 後端引擎可以對執行的語句進行優化，比如進行自動並行化處理。

後端引擎必須要解決的問題就是資料依賴和合理的排程。但這些操作對於前端的使用者來說是完全透明的。我們可以使用wait_to_read來等侍後端對於NDArray操作的完成。在NDArray模組一類將資料拷貝到其他模組的操作，內部已經使用了wait_to_read，比如asnumpy()。

import time
def do(x, n):
    """push computation into the backend engine"""
    return [mx.nd.dot(x,x) for i in range(n)]
def wait(x):
    """wait until all results are available"""
    for y in x:
        y.wait_to_read()

tic = time.time()
a = mx.nd.ones((1000,1000))
b = do(a, 50)
print('time for all computations are pushed into the backend engine:\n %f sec' % (time.time() - tic))
wait(b)
print('time for all computations are finished:\n %f sec' % (time.time() - tic))

除了分析資料的讀寫依賴外，後端的引擎還能夠將沒有彼此依賴的操作語句進行並行化排程。比如下面的程式碼第二行和第三行可以被並行的執行。

a = mx.nd.ones((2,3))
b = a + 1
c = a + 2
d = b * c

下面的程式碼演示了在不同裝置上並行排程

n = 10
a = mx.nd.ones((1000,1000))
b = mx.nd.ones((6000,6000), gpu_device)
tic = time.time()
c = do(a, n)
wait(c)
print('Time to finish the CPU workload: %f sec' % (time.time() - tic))
d = do(b, n)
wait(d)
print('Time to finish both CPU/GPU workloads: %f sec' % (time.time() - tic))

tic = time.time()
c = do(a, n) 
d = do(b, n) #上面兩條語句可以同時執行，一條在CPU上運算，一條在GPU上運算
wait(c)
wait(d)
print('Both as finished in: %f sec' % (time.time() - tic))

參考資源