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")

待续.