1. 程式人生 > >NumPy學習指南 學習筆記(二) 常用函式

NumPy學習指南 學習筆記(二) 常用函式

1.  檔案讀寫
通常情況下,資料是以檔案形式儲存的。學會讀寫檔案是深入學習Numpy的基礎。
1.1 建立單位矩陣,並存入txt檔案
i2 = np.eye(2)
i2
Out[84]: 
array([[ 1.,  0.],
       [ 0.,  1.]])

使用savetxt 函式將資料儲存到檔案中,當然需要指定檔名以及要保持的陣列。
np.savetxt("eye.txt", i2)

如儲存成.py 檔案,程式碼如下:
import numpy as np


i2 = np.eye(2)
print i2
np.savetxt("eye.txt", i2)

1.2 讀入CSV 檔案
Numpy 中的loadtxt 函式可以方便地讀取CSV檔案,自動切分欄位,並將資料載入NumPy陣列。
資料中以蘋果公司的歷史股價資料為例展開敘述。股價資料儲存在csv檔案中,第一列為股票程式碼以標識股票(蘋果公司股票程式碼為AAPL),第二列為dd-mm-yyyy格式的日期,第三列為空,隨後各列依次是開盤價、最高價、最低價和收盤價,最後一列為當日的成交量。下面為一行資料:
AAPL,28-01-2011, ,344.17,344.4,333.53,336.1,21144800
下面我們只關注股票的收盤價和成交量。上面的例項資料中,收盤價為336.1,成交量為22144800。將收盤價和成交量分別載入到兩個陣列中,如下所示:
c,v = np.loadtxt('data.csv', delimiter = ',', usecols = (6, 7), unpack = Ture)
可以看到,資料儲存在data.csv 檔案中,我們設定分隔符為,(英文標點逗號),usecols的引數為一個元組,以獲取第7欄位至第8欄位的資料,也就是股票的收盤價和成交量資料。unpack引數設定為True,意思是分拆儲存不同列的資料,即分別將兩個陣列賦值給變數c 和 v。

2. 計算成交量加權平均價格(VWAP)
VWAP(Volume-Weighted Average Price, 成交量加權平均價格)是一個非常重要的經濟學量,它代表著金融資產的“平均”價格。某個價格的成交量越高,該價格所佔的權重就越大。VWAP就是以成交量為權重計算出來的加權平均值,常用於演算法交易。
計算步驟如下:
   (1) 將資料讀入陣列。
   (2)計算VWAP。
import numpy as np


 c,v = np.loadtxt('data.csv', delimiter = ',', usecols = (6, 7), unpack = Ture)
 vwap = np.average(c, weights = v)
 print "VWAP = " VWAP
 
 
 The output is
 VWAP = 350.589549353

   2.1 算術平均值函式
Numpy 中mean函式很友好,可以計算陣列元素的算術平均值。用法如下:
 print "mean = ", np.mean(c)
 mean = 351.03766667
   2.2 時間加權平均價格
在經濟學中,TWAP(Time-Weighted Average Price,時間加權平均價格)是另一種“平均”價格的指標。
TWAP 是一個VWAP的變種,基本思想就是最近的價格重要性大一些,所以應該對近期的價格給以較高的權重。最簡單的方法就是用arange函式建立一個從0開始依次增長的自然數序列,自然數的個數即為收盤價的個數。當然,這並不一定是正確的計算TWAP的方式。事實上,關於股價分析的大部分示例都僅僅是為了說明問題。計算TWAP的程式碼如下:
 t = np.arange(len(c))
 print "twap = ", np.average(c, weights = t)
 輸出結果:
 twap = 352.428321839
例子中,TWAP的值甚至比算術平均值還要高。


   2.3 找最大值和最小值
min函式和max函式能夠滿足需求。我們按如下步驟來找最大值和最小值。
首先讀入資料,將每日最高價和最低價的資料載入陣列:
 h,l = np.loadtxt('data.csv', delimiter = ',', usecols = (4, 5), unpack = Ture)
下方程式碼即可獲取價格區間:
 print "highest = ", np.max(h)
 print "lowest = ",np.min(l)
 程式將返回如下結果:
 highest = 364.9
 lowest = 333.53
NumPy 中有一個ptp函式可以計算陣列的取值範圍。該函式返回的是陣列元素的最大值和最小值之間的差值。也就是說,返回值等於 max(array) - min(array)。呼叫ptp函式:
 print "Spread high price", np.ptp(h)
 print "Spread low price", np.ptp(l)
 
 我們將看到如下結果:
 
 Spread high price 24.86
 Spread low price 26.97

3. 統計分析
可以使用一些閥值來除去異常值,但其實有更好的方法,那就是中位數。將各個變數值按大小順序排列起來,形成一個數列,居於數列中間位置的那個數即為中位數。例如,我們有1、2、3、4、5 這5個數值,那麼中位數就是中間的數字3。

3.1 計算中位數
計算收盤價的中位數。使用median 函式找到中位數。
  c = np.loadtxt('data.csv', delimiter = ',', usecols = (6, ), unpack = Ture)
  print "median = ", np.median(c)
  
  輸出內容如下:
  
  median = 352.055
驗證函式正確性的話,可以使用msort 函式,對資料進行排序,然後去查詢中位數。
  sorted_close = np.msort(c)
  print "sorted = ", sorted_close

獲取位於中間的數字:
  N = len(c)
  print "middle = ", sorted[(N-1)/2]
  
  輸出結果:
  
  middle = 351.99

這個值和median 函式給出的值不一樣,這是因為所給樣本為偶數個,中位數應該為中間兩個數的平均值。
print "average middle = ", (sorted[N/2] + sorted[(N - 1)/2])/2
  
  輸出結果:
  
  average middle = 352.055

   3.2 方差的計算
  print "variance = ", np.var(c)
  
  結果:
  
  variance = 50.1265178889
方差:是指各個資料與所有資料算術平均數的離差平方和除以資料個數所得到的值。
一些書裡面說,應該用資料個數減去1去除離差平方和。-- 樣本個數減1(即 n-1)稱為自由度,減1可以保證樣本方差是一個無偏估計量。
  print "variance from definition = ", np.mean((c - c.mean()) ** 2)
  
  輸出結果:
  
  variance from definition = 50.1265178889
我們直接在陣列c上呼叫了mean方法c.mean(), ndarray 物件有mean方法,這將給你帶來便利。從現在開始,記住這種使用者是正確無誤的。
示例程式碼 simplestats.py 如下:
  import numpy as np
  
  c = np.loadtxt('data.csv', delimiter = ',', usecols = (6, ), unpack = True)
  print "median =", np.median(c)
  sorted "sorted =", sorted
  
  N = len(c)
  print "middle =", sorted[(N-1)/2]
  print "average middle =", (sorted[N/2] + sorted[(N-1)/2])/2
  
  print "variance =", np.var(c)
  print "variance from definition =", np.mean((c - c.mean()) ** 2)

4. 分析股票的收益率
在學術文獻中,收盤價的分析常常是基於股票收益率和對數收益率的。簡單收益率是指相鄰兩個價格之間的變化率,而對數收益率是指所有價格取對數後兩兩之間的差值。對數收益率可以用來衡量價格的變化率。由於收益率是一個比值,例如我們用美元除以美元(也可以是其他貨幣單位),因此它是無量綱的。投資者最感興趣的是收益率的方差或者標準差,因其代表著投資風險的大小。
   (1) 首先計算簡單收益率。NumPy中的diff函式可以返回一個由相鄰陣列元素的差值構成的陣列。這點類似與微積分中的微分。為了計算收益率,我們還需要用差值除以以前一天的價格。diff 返回的陣列比收盤價陣列少一個元素。經過仔細考慮,使用如下程式碼:
returns = np.diff( arr ) / arr[ : -1]
注意,我們沒有用收盤價陣列中的最後一個值做除數。接下來,用std函式計算標準差:
print "Standard deviation =", np.std(returns)
輸出結果如下:
Standard deviation = 0.0129221344368

   (2)  對收益率計算起來甚至更簡單一些。我們先用log函式得到每個收盤價的對數,再對結果使用diff 函式即可。
logreturns = np.diff(np.log(c) )
一般情況下,我們應檢查輸入陣列以確保其不含有零和負數。否則,將得到一個錯誤提示。
   (3)  我們很可能對哪些交易日的收益率為正值非常感興趣。在完成了前面的步驟之後,我們只需要用where函式就可以做到這一點。where函式可以根據指定的條件返回所有滿足條件的陣列元素的索引值。 輸入如下程式碼:
posretindices = np.where(return > 0)
print "Indices with positive returns", posretindices

即可輸出該陣列中所有正值元素的索引。
Indices with positive returns (array([ 0, 1, 4, 5, 6, 7, 9, 10, 11, 12, 16, 17, 18, 19, 21, 22, 23, 25, 28]), )

   (4)  在投資學中,波動率(volatility)是對價格變動的一種度量。歷史波動率可以根據歷史價格資料計算得出。計算曆史波動率(如年波動率或月波動率)時,需要用對數收益率。年波動率等於對數收益率的標準差除以其均值,再除以交易日倒數的平方根,通常交易日取252 天。 我們用std 和mean 函式來計算,程式碼如下所示:
  annual_volatility = np.std(logreturns)/np.mean(logreturns)
  annual_volatility = annual_volatility / np.sqrt(1. /252.)
  print annual_volatility

   (5) 請注意sqrt函式中的除法運算。在python中,整數的除法和浮點數的除法運算機制不同,必須使用浮點數才能得出正確的結果。與計算年波動率的方法類似,計算月波動率如下:
print "Monthly volatility", annual_volatility * np.sqrt(1. /12.)

5. 日期分析
首先,要讀入收盤價資料。隨後根據星期幾來切分收盤價資料,並分別計算平均價格。最後,我們將找出一週中哪一天的平均收盤價最高,哪一天的最低。
程式設計師不喜歡日期,因為處理日期總是很繁瑣。NumPy是面向浮點數運算的,因此需要對日期做一些專門處理。請自行嘗試如下程式碼,單獨編寫指令碼或使用本書附帶的程式碼檔案:
dates, close = np.loadtxt('data.csv', delimiter = ',', usecols = (1, 6), unpack = True)
執行以上程式碼後會得到一個錯誤提示:
ValueError:invalid literal for float(): 28-01-2011
按如下步驟處理日期。
   (1)顯然,NumPy 嘗試把日期轉換成浮點數。我們需要做的是顯式地告訴NumPy怎樣來轉換日期,而這需要用到loadtxt 函式中的一個特定的引數。這個引數就是converters, 它是一本書籍列和轉換函式之間進行對映的字典。
為此,我們必須寫出轉換函式:
# 星期一 0
# 星期二 1
# 星期三 2
# 星期四 3
# 星期五 4
# 星期六 5
# 星期日 6
import datetime

def datestr2num(s):
    return datetime.datetime.strptime(s, "%d-%m-%Y").date().weekday()

datestr2num('12-05-2017')
Out[100]: 4
我們將日期作為字串傳給datestr2num函式,如“12-05-2017”。這個字串首先會按照指定的形式“%d-%m-%Y"轉換成一個datetime物件。補充一點,這是由Python標準庫提供的功能,與NumPy無關。隨後,datetime 物件被轉換為date物件。最後,呼叫weekday方法返回一個數字。如同你在註釋中看到的,這個數字可以是0到6的整數,0代表星期一,6代表星期天。當然,具體的數字並不重要,只是用作標識。
   (2) 接下來,我們將日期轉換函式掛接上去,這樣就可以讀入資料了。
dates, close = np.loadtxt('data.csv', delimiter = ',', usecols = (1, 6), converters{1:datestr2num}, unpack = True)
print "Dates = ", dates
輸出結果如下:
Dates = [4. 0. 1. 2. 3. 4. 0. 1. 2. 3. 4. 0. 1. 2. 3. 4. 0. 1. 2. 3. 4. 0. 1. 2. 3. 4. 0. 1. 2. 3. 4.] 

如你所見,沒有出現星期六和星期天。股票交易在週末是休市的。
   (3)我們來建立一個包含5個元素的陣列,分別代表一週的5個工作日。陣列元素將初始化為0.
averages = np.zeros(5)
這個陣列將用於儲存各工作日的平均收盤價。
   (4) 我們已經知道,where函式會根據指定的條件返回所有滿足條件的陣列元素的索引值。take函式可以按照這些索引值從陣列中取出相應的元素。我們將用take函式來獲取各個工作日的收盤價。在下面的迴圈體中,我們將遍歷0到4的日期標識,或者說是遍歷星期一到星期五,然後用where函式得到各工作日的索引值並存儲在indices陣列中。再用take函式獲取這些索引值相應的元素值。最後,我們對每個工作日計算出平均值存放在averages陣列中。程式碼如下:
  for i in range(5):
      indices = np.where(dates == i)
      prices = np.take(close, indices)
      avg = np.mean(prices)
      price "Day", i, "prices", prices, "Average", avg
      
輸出結果如下:
Day 0 prices [[ 339.32 351.88 359.18 353.21 355.36]] Average 351.79
Day 1 prices [[ 345.03 355.2 359.9 338.61 349.31 355.76]] Average 350.635
Day 2 prices [[ 344.32 358.16 363.13 342.62 352.12 352.47]] Average 352.136666667
Day 3 prices [[ 343.44 354.54 358.3 342.88 359.56 346.67]] Average 350.898333333
Day 4 prices [[ 336.1 346.5 356.85 350.56 348.16 360. 351.99]] Average 350.022857143

   (5) 如果你願意,還可以找出哪個工作日的平均收盤價是最高的,哪個是最低的。用max和min 函式即可,程式碼如下:
top = np.max(averages)
print "Highest average", top
print "Top day of the week", np.argmax(averages)
bottom = np.min(averages)
print "Lowest average", bottom
print "Bottom day of the week", np.argmin(averages)
輸出結果如下:
Highest average 352.136666667
Top day of the week 2
Lowest average 350.022857143
Bottom day of the week 4
argmin 函式返回的是averages 陣列中最小元素的索引值,這裡是4,也就是星期五。而argmax函式返回是averages陣列中最大元素的索引值,這裡是2,也就是星期三。示例程式碼如下weekdays.py檔案。
import numpy as np
from datetime import datetime
# 星期一 0
# 星期二 1
# 星期三 2
# 星期四 3
# 星期五 4
# 星期六 5
# 星期日 6
def datestr2num(s):
return datetime.strptime(s, "%d-%m-%Y").date().weekday()
dates, close=np.loadtxt('data.csv', delimiter=',', usecols=(1,6), converters={1:
datestr2num}, unpack=True)
print "Dates =", dates
averages = np.zeros(5)
for i in range(5):
    indices = np.where(dates == i)
    prices = np.take(close, indices)
    avg = np.mean(prices)
    print "Day", i, "prices", prices, "Average", avg
    averages[i] = avg

top = np.max(averages)
print "Highest average", top
print "Top day of the week", np.argmax(averages)
bottom = np.min(averages)
print "Lowest average", bottom
print "Bottom day of the week", np.argmin(averages)

6. 周彙總
我們要彙總整個交易週中從週一到週五的所有資料。資料覆蓋的時間段內有一個節假日:2月21日是總統紀念日。這天是星期一,美國股市休市,因此在我們的示例資料中沒有這一天的資料記錄。資料中的第一天為星期五,處理起來不太方便。按照如下步驟來彙總資料。
   (1)為了簡單起見,我們只考慮前3周的資料,這樣就避免了節假日造成的資料缺失。你可以稍後嘗試對其進行拓展。
close = close[: 16]
dates = dates[: 16]
   (2) 首先我們來找到示例資料中的第一個星期一。回憶一下,在Python中星期一對應的編碼是0,這可以作為where函式的條件。接著,我們要取出陣列中的首個元素,其索引值為0。但where函式返回的結果是一個多維陣列,因此要用ravel函式將其展平。
# 找到第一個星期一
first_monday = np.ravel(np.where(dates == 0))[0]
print "The first Monday idex is ", first_monday

輸出結果如下:

The first Monday index is 1

   (3)下面要做的是找到示例資料的最後一個星期五,方法和找第一個星期一類似。星期五相對應的編碼是4.此外,我們用-1作為索引值來定位陣列的最後一個元素。
# 找到最後一個星期五
last_friday = np.ravel(np.where(dates == 4))[-2]
print "The last Friday index is ", last_friday

輸出結果如下:

The last Friday index is 15

接下來建立一個數組,使用者儲存三週內每一天的索引值。
weeks_indices = np.arange(first_monday, last_friday + 1)
print "Weeks indices initial", weeks_indices

   (4) 按照每個子陣列5個元素,用split函式切分陣列:

weeks_indices = np.split(weeks_indices, 5)
print "Weeks indices after split", weeks_indices

輸出結果如下:

Weeks indices  after split [array([1, 2, 3, 4, 5]), array([6, 7, 8, 9, 10]), array([11, 12, 13, 14, 15])]

  (5) 在NumPy中,陣列的維度被稱作軸。現在我們來熟悉一下apply_along_axis 函式。這個函式會盜用另外一個由我們給出的函式,作用於每個陣列元素上。目前我們的陣列中有3個元素,分別對應於示例資料中的3個星期,元素中的索引值對應於示例資料中的1天。在呼叫apply_along_axis 時提供我們自定義的函式名summarize, 並指定要作用的軸或維度的編號(如取1)、目標陣列以及可變數量的summarize函式的引數。

weeksummary = np.apply_alone_axis(summarize, 1, weeks_indices, open, high, low, close)
print "Week summary", weeksummary

   (6) 編寫summarize函式。該函式將為每一週的資料返回一個元組,包含這一週的開盤價、最高價、最低價和收盤價,類似於每天的盤後資料。
def summarize(a, o, h, l, c):
    monday_open = o[a[0]]
    week_high = np.max(np.take(h, a))
    week_low = np.min(np.take(l, a))
    friday_close = c[a[-1]]
    
    return ("APPL", monday_open, week_high, week_low, friday_close)
注意,我們用take函式來根據索引值獲取陣列元素的值,並用max和min函式輕鬆計算出一週的最高股價和最低股價。一週的開盤價即為週一的開盤價,而一週的收盤價即為週五的收盤價。
Week summary [['APPL' '335.8' '346.7' '334.3' '346.5']
['APPL' '347.89' '360.0' '347.64' '356.85']
['APPL' '356.79' '364.9' '349.52' '350.56']]
   (7) 使用NumPy中的savetxt 函式,將資料儲存至檔案。
np.savetxt("weeksummary.csv", weeksummary, delimiter = ",", fmt = "%s")
如程式碼中所示,我們指定了檔名、需要儲存的陣列名、分隔符(在這個例子中為英文標點逗號)以及儲存浮點數的格式。
格式字串以一個百分號開始。接下來是一個可選的標誌字元:- 表示結果左對齊, 0 表示左端補0, + 表示輸出符號(正號+或負號-)。第三部分為可選的輸出寬度引數,表示輸出的最小位數。第四部分是精度格式符,以"." 開頭,後面跟一個表示精度的整數。最後是一個型別指定字元,在我們的例子中指定為字串型別。
------------------------------------------------
字元編碼                         含 義
------------------------------------------------
   c                         單個字元
   d或i                      十進位制有符號整數
   e或E                      科學記數法表示的浮點數
   f                         浮點數
   g或G                      自動在e、E和f中選擇合適的表示法
   o                         八進位制有符號整數
   s                         字串
   u                         十進位制無符號整數
   x或X                      十六進位制無符號整數
------------------------------------------------
我們剛剛完成的事情在其他程式語言中幾乎是無法完成的。我們定義了一個函式summarize, 把它作為引數傳給了apply_along_axis函式, 而summarize 的引數也通過apply_along_axis 的引數列表便捷地傳遞了過去。示例程式碼見weeksummary.py檔案。
import numpy as np
from datetime import datetime
# 星期一 0
# 星期二 1
# 星期三 2
# 星期四 3
# 星期五 4
# 星期六 5
# 星期日 6
def datestr2num(s):
    return datetime.strptime(s, "%d-%m-%Y").date().weekday()
dates, open, high, low, close=np.loadtxt('data.csv', delimiter=',', usecols=(1, 3, 4, 5, 6), converters={1: datestr2num}, unpack=True)
close = close[:16]
dates = dates[:16]
# get first Monday
first_monday = np.ravel(np.where(dates == 0))[0]
print "The first Monday index is", first_monday
# get last Friday
last_friday = np.ravel(np.where(dates == 4))[-1]
print "The last Friday index is", last_friday
weeks_indices = np.arange(first_monday, last_friday + 1)
print "Weeks indices initial", weeks_indices
weeks_indices = np.split(weeks_indices, 3)
print "Weeks indices after split", weeks_indices
def summarize(a, o, h, l, c):
    monday_open = o[a[0]]
    week_high = np.max( np.take(h, a) )
    week_low = np.min( np.take(l, a) )
    friday_close = c[a[-1]]
    return("APPL", monday_open, week_high, week_low, friday_close)
    
weeksummary = np.apply_along_axis(summarize, 1, weeks_indices, open, high, low, close)
print "Week summary", weeksummary
np.savetxt("weeksummary.csv", weeksummary, delimiter=",", fmt="%s")


待續.