科学计算与数据可视化1

科学计算与数据可视化1

@(程序设计)

numpy模块

• Numpy最重要的一个特色就是其N维数组对象（即ndarray）该对象是一个快速而灵活的大数据集容器。
• 使用Numpy，开发人员能够执行如下操做：

一、数组的算数和逻辑运算。

二、傅立叶变换和用于图形操做的例程。

三、与线性代数有关的操做。

numpy.ndarray()

• ndarray是Numpy的数组类，其中的全部元素必须是相同的数据类型。ndarray类的重要对象属性有：

• 利用array函数，能够将序列类型的对象（元组、列表和其余数组）转换成数组类型ndarray。

  • ndarray.ndim：数组维度。两行及以上为二维
  • ndarray.shape：表示数组各维度大小的元组。
  • ndarray.size：数组元素的总个数，等于shape属性中元组元素的乘积。
  • ndarray.dtype：数组中元素的数据类型
  • numpy.ndarray()就是numpy的构造函数，咱们可使用这个函数建立一个ndarray对象。构造函数有以下几个可选参数：

np.savetxt()

• 将array保存到txt文件，并保持原格式

np.loadtxt()

• loadtxt(fname, dtype=<class 'float'>, comments='#', delimiter=None, converters=None, skiprows=0, usecols=None, unpack=False, ndmin=0)
• np.loadtxt()用于从文本加载数据。文本文件中的每一行必须含有相同的数据。

arange()

• numpy.arange([start, ]stop, [step, ]dtype=None)

arange函数是numpy内置的相似range的函数，其返回的是数组对象，而不是列表。

• range()与arrange()的区别

一、range()和np.arange()的返回类型不一样，range()返回的是range；
二、object，而np.arange()返回的是ndarray类型；
三、range()不支持步长为小数，而np.arange()支持步长(step)为小数；
四、range()和np.arange()均可用于迭代；
五、range()和np.arange()都有三个参数，以第一个参数为起点，第三个参数为步长，截止到第二个参数以前的不包括第二个参数的数据序列。
五、range()可用于迭代，而np.arange做用远不止于此，它是一个序列，可被当作向量使用。

array()

• array(...)
  array(object, dtype=None, copy=True, order='K', subok=False, ndmin=0)
• 数组的特色：

一、数组是相同数据类型的元素的集合。
二、数组中的各元素的存储是有前后顺序的，它们在内存中按照这个前后顺序连续存放在一块儿。
三、数组元素用整个数组的名字和它本身在数组中的顺序位置来表示。例如，a[0]表示名字为a的数组中的第一个元素，a[1]表明数组a的第二个元素，以此类推。

Array[0:]　　——>切片从前面序号“0”开始到结尾，包括“0”位
　　[2, 3, 9, 1, 4, 7, 6, 8]
Array[:-1]　　——>切片从后面序号“-1”到最前，不包括“-1”位
　　[2, 3, 9, 1, 4, 7, 6]
Array[3:-2]　　——>切从前面序号“3”开始(包括)到从后面序号“-2”结束(不包括)
　　[1, 4, 7]
Array[3::2]　　——>从前面序号“3”(包括)到最后，其中分隔为“2”
　　[1, 7, 8]
Array[::2]　　——>从整列表中切出，分隔为“2”
　　[2, 9, 4, 6]
Array[3::]　　——>从前面序号“3”开始到最后，没有分隔
　　[1, 4, 7, 6, 8]
Array[3::-2]　　——>从前面序号“3”开始，往回数第二个，由于分隔为“-2”
　　[1, 3]　　
Array[-1]　　　　——>此为切出最后一个
　　8
Array[::-1]　　　　——>此为倒序
　　[8, 6, 7, 4, 1, 9, 3, 2]

ndarray()

建立二维数组、访问数组对象属性

import numpy as np
list1=[5,6.5,9,2,3,7.8,5.6,4.9]
arr1=np.array(list1)
print(arr1)
print(arr1.dtype)   #数组中元素的数据类型
print(arr1.ndim)    #数组维度。两行及以上为二维
print(arr1.shape)   #表示数组各维度大小的元组。行*列
print(arr1.size)    #数组元素的总个数，等于shape属性中元组元素的乘积。
list2=[[1,2,3,4,5],[6,7,8,9,10]]
arr2=np.array(list2)
print(arr2)
print(arr2.dtype)   #数组中元素的数据类型
print(arr2.ndim)    #数组维度。两行及以上为二维
print(arr2.shape)   #表示数组各维度大小的元组。行*列
print(arr2.size)    #数组元素的总个数，等于shape属性中元组元素的乘积。

[5. 6.5 9. 2. 3. 7.8 5.6 4.9]
float64
1
(8,)
8
[[ 1 2 3 4 5]
[ 6 7 8 9 10]]
int32
2
(2, 5)
10

建立指定数据类型的数组对象

import numpy as np
list2=[[1,2,3,4,5],[6,7,8,9,10]]
arr2=np.array(list2)
arr3=np.array([10,20,30,40],dtype=np.float64)
print(arr3)
arr4=arr2.astype(np.float64) #转换数据类型 float->int
print(arr4)
print(arr4.dtype)

[10. 20. 30. 40.]
[[ 1. 2. 3. 4. 5.]
[ 6. 7. 8. 9. 10.]]
float64

练习eg13_array1.py()

产生3个数组l一、l二、l3，分别为[0 1 2 3]、[0 2 4 6]、[0 3 6 9]，输出l一、l二、l3以及(l1,l2,l3)，并将(l1,l2,l3)的内容经过np.savetxt方法存入aa.txt文件中，而后用np.loadtxt读出数据并显示。程序保存为eg13_array1.py

import numpy as np
l1=np.arange(4)  #易错
l2,l3=l1*2,l1*3
print('l1：',l1)
print('l2：',l1)
print('l3：',l1)
print((l1,l2,l3))
np.savetxt('aa.txt',(l1,l2,l3)) #易错
aa=np.loadtxt('aa.txt')
print(aa)

l1 [0 1 2 3]
l2 [0 2 4 6]
l3 [0 3 6 9]
(array([0, 1, 2, 3]), array([0, 2, 4, 6]), array([0, 3, 6, 9]))
aa: [[0. 1. 2. 3.]
[0. 2. 4. 6.]
[0. 3. 6. 9.]]

思考：继续eg13_array1.py

import numpy as np
l1=np.arange(4)
l2,l3=l1*2,l1*3
print(l1)
print(l2)
print(l3)
aa=np.loadtxt('aa.txt')
print('aa的类型：',aa.dtype)
print('aa的维度：',aa.ndim)
print('aa的总个数：',aa.size)
print('aa的形状：',aa.shape) #reshape()数组对象中的方法，用于改变数组的形状
print('bb:',aa.reshape(4,3))  #reshape()改成一个四维数组

[0 1 2 3]
[0 2 4 6]
[0 3 6 9]
aa的类型： float64
aa的维度： 2
aa的总个数： 12
aa的形状： (3, 4)
bb: [[0. 1. 2.]
[3. 0. 2.]
[4. 6. 0.]
[3. 6. 9.]]

定义数组的类型

import numpy as np
d1=np.loadtxt('bb.txt')
print('d1:',d1,'sum=',sum(d1))
d2=np.loadtxt('bb.txt',dtype=int)
print('d2:',d2,'sum=',sum(d2))
d3=np.loadtxt('bb.txt',dtype=str)
print('d3:',d3,'sum=',sum([eval(i) for i in d3]))  #这边不是很明白

d1: [1. 2. 3. 4. 5.] sum= 15.0
d2: [1 2 3 4 5] sum= 15
d3: ['1' '2' '3' '4' '5'] sum= 15

astype()显示转换类型

import numpy as np
aa=np.array(['1.25','-9.6','42'],dtype=np.str)
bb=aa.astype(float)
print(bb)

[ 1.25 -9.6 42. ]

数组和标量之间的运算

import numpy as np
arr=np.array([[1,2,3],[4,5,6]])
print(arr)
print(arr*2)
print(1/arr)
print(arr-arr*2)

[[1 2 3]
[4 5 6]]
[[ 2 4 6]
[ 8 10 12]]
[[1. 0.5 0.33333333]
[0.25 0.2 0.16666667]]
[[-1 -2 -3]
[-4 -5 -6]]

数组的元素级运算与函数

• 大小相等的数组之间的任何算数运算都会应用到元素级

import numpy as np
arr1=np.arange(1,16).reshape(3,5)
print(arr1)
arr2=np.arange(1,30,2).reshape(3,5)
print(arr2)
print(arr1*arr2)
print(arr1/arr2)
print(arr1+arr2)
print(arr1-arr2)

[[ 1 2 3 4 5]
[ 6 7 8 9 10]
[11 12 13 14 15]]
[[ 1 3 5 7 9]
[11 13 15 17 19]
[21 23 25 27 29]]
[[ 1 6 15 28 45]
[ 66 91 120 153 190]
[231 276 325 378 435]]
[[1. 0.66666667 0.6 0.57142857 0.55555556]
[0.54545455 0.53846154 0.53333333 0.52941176 0.52631579]
[0.52380952 0.52173913 0.52 0.51851852 0.51724138]]
[[ 2 5 8 11 14]
[17 20 23 26 29]
[32 35 38 41 44]]
[[ 0 -1 -2 -3 -4]
[ -5 -6 -7 -8 -9]
[-10 -11 -12 -13 -14]]

二元函数及说明

一维数组的索引和切片（同列表）

• 数组切片是原始数组的视图，数据并不会被复制，即视图上的任何修改都会直接反映到源数组上。

读csv文件

• c,v=np.loadtxt('data.csv', delimiter=',', usecols=(6,7), unpack=True,skiprows=1)

  • fname：读取文件的文件名。例如C:/Dataset/iris.txt。
  • delimiter：数据之间的分隔符。如使用逗号","。
  • dtype：数据类型。如float，str等。
  • usecols：选取数据的列。
  • skiprows : list-like or integer, default None须要忽略的行数（从文件开始处算起），或须要跳过的行号列表（从0开始）。

平均值和加权平均值

加权平均值的概念

• np.mean()
• np.average()
• 对应的权值列表：weights=[]
• 将各数值乘以相应的权数，而后加总求和获得整体值，再除以总的单位数。

import numpy as np
a=(70,80,60)
print(np.mean(a)) #平均值
print(np.average(a,weights=[3,3,4]))  #加权平均值  average()

70.0
69.0

VWAP(Volume-Weighted Average Price)：成交量加权平均价格

• 表明金融资产的“平均”价格，某个价格的成交量越高，该价格所占的权重就越大
• VWAP就是以成交量为权重计算出来的加权平均值

TWAP(Time-Weighted Average Price)：时间加权平均价格

• 只是一个变种，基本的思想就是最近的价格重要性大一些,仅仅为了说明问题，并不必定很是正确
• 选用arange函数建立一个从0开始依次增加的天然数序列，天然数的个数即为收盘价的个数

例题：eg3_vwap_twap.py

利用data.csv文件，读出收盘价、成交量，计算成交量加权平均价格VWAP、算术平均值和时间加权平均价格TWAP，程序保存为eg13_vwap_twap3.py

import numpy as np
c,v=np.loadtxt('data.csv', delimiter=',', \
               usecols=(6,7), unpack=True,skiprows=1)   #usecols列数是从0开始的
print('收盘价：',c)
print('成交量：',v)
vwap=np.average(c,weights=v)
print('成交量加权平均价格VWAP=%f'%vwap)  #这边看不太懂？成交量加权平均价格
print('算术平均值mean1=',np.mean(c)) #算术平均值
print('算术平均值mean1=',c.mean())   #算术平均值
t=np.arange(len(c))
print('时间加权平均价格=',np.average(c,weights=t))  #时间加权平均价格

收盘价： [336.1 339.32 345.03 344.32 343.44 346.5 351.88 355.2 358.16 354.54
356.85 359.18 359.9 363.13 358.3 350.56 338.61 342.62 342.88 348.16
353.21 349.31 352.12 359.56 360. 355.36 355.76 352.47 346.67 351.99]
成交量： [21144800. 13473000. 15236800. 9242600. 14064100. 11494200. 17322100.

1. 1. 1. 1. 1. 1. 17184100.
2. 1. 1. 1. 1. 1. 14395400.
3. 1. 1. 1. 1. 1. 16192700.
4. 16824200.]

成交量加权平均价格VWAP=350.589549
算术平均值mean1= 351.0376666666667
算术平均值mean1= 351.0376666666667
时间加权平均价格= 352.4283218390804

ptp函数

• ptp函数能够计算数组的取值范围，返回数组元素的最大值和最小值之间的差值，即返回值=max(array)-min(array)

示例4：eg13_max_min4.py

利用data.csv，计算最高价中的最大值和最小值以及最高价和最低价中最大值和最小值之间的差值，程序保存为eg13_max_min4.py

import numpy as np
hst,lst=np.loadtxt('data.csv',delimiter=',',usecols=(4,5),\
                   unpack=True,skiprows=1)
hh=np.max(hst)
ll=np.min(hst)
print('最高价中最大值=',hh)
print('最高价中最小值=',ll)  #np.ptp()
print('最高价中最大值和最小值之间的差值=','%.2f'%np.ptp(hst))    
print('最低价中最大值和最小值之间的差值=','%.2f'%np.ptp(lst))

最高价中最大值= 364.9
最高价中最小值= 340.04
最高价中最大值和最小值之间的差值= 24.86
最低价中最大值和最小值之间的差值= 26.97

数组的基本统计分析函数及说明

argmax、argmin、max、min

import numpy as np
b=np.array([1,5,8,9,334])
print(np.argmax(b))  #np.argmax();np.argmin
print(np.argmin(b))    #reshape()的乘积恰好是arange()的列*行
print(np.arange(7,17).reshape(2,5))

4
0
[[ 7 8 9 10 11]
[12 13 14 15 16]]

axis

axis=0 跨行
axis=1 跨列

import numpy as np
a=np.arange(7,17).reshape(2,5)
print(a)
print(np.argmin(a,axis=0))
print(np.argmin(a,axis=1))

[[ 7 8 9 10 11]
[12 13 14 15 16]]
[0 0 0 0 0]
[0 0]

NumPy经常使用函数：计算中位数和方差

思考1：si13_var1.py

• 中位数：np.median()
• 方差：np.var()

import numpy as np
price=np.loadtxt('data.csv',delimiter=',',usecols=(6,),unpack=True,skiprows=1)
print(price)
print('中位数：',np.median(price))
print('方差：', np.var(price))

[336.1 339.32 345.03 344.32 343.44 346.5 351.88 355.2 358.16 354.54
356.85 359.18 359.9 363.13 358.3 350.56 338.61 342.62 342.88 348.16
353.21 349.31 352.12 359.56 360. 355.36 355.76 352.47 346.67 351.99]
中位数： 352.055
方差： 50.126517888888884

import numpy as np
price=np.loadtxt('data.csv',delimiter=',',usecols=(6,),unpack=True,skiprows=1)
print(price)
print('中位数',np.median(price))
#验证刚才求的中位数是否正确
#一、将数组进行排序
sorted = np.msort(price)
print(sorted)
#二、计算数组的元素个数
n = len(sorted)
print(n)
#print('middle', '=', sorted[n - 1]//2)  #n为奇数
print('中位数：', (sorted[n // 2] + sorted[(n - 1) // 2]) / 2)    #n为偶数
#方差等于各个数据与平均数之差的平方和的平均数,用来度量随机变量和其数学指望（即均值）之间的偏离程度
#((x1 - a)^2 + (x2 - a)^2 + (x3 - a)^2 ... + (xn - a)^2) / n
print('方差：', np.var(price))

[336.1 339.32 345.03 344.32 343.44 346.5 351.88 355.2 358.16 354.54
356.85 359.18 359.9 363.13 358.3 350.56 338.61 342.62 342.88 348.16
353.21 349.31 352.12 359.56 360. 355.36 355.76 352.47 346.67 351.99]
中位数 352.055
[336.1 338.61 339.32 342.62 342.88 343.44 344.32 345.03 346.5 346.67
348.16 349.31 350.56 351.88 351.99 352.12 352.47 353.21 354.54 355.2
355.36 355.76 356.85 358.16 358.3 359.18 359.56 359.9 360. 363.13]
30
中位数： 352.055
方差： 50.126517888888884

简单收益率和对数收益率

• 简单收益率：相邻两个价格之间的变化率
• 对数收益率：全部价格取对数后两两之间的差值，也能够用来衡量价格的变化率

计算历史波动率（如年波动率和月波动率）时，须要用到对数收益率。

* 若是a的x次方等于N（a>0，且a不等于1），那么数x叫作以a为底N的对数（logarithm），记做x=logaN。其中，a叫作对数的底数，N叫作真数。

• 年波动率 = 对数收益率的标准差std / 其均值mean，再除以交易日倒数的平方根。一般交易日取252天。

diff函数

• 返回一个由相邻数组元素的差值构成的数组

import numpy as np
a=np.array([1,2,5,4,3,7,8,38])
print(np.diff(a))

[ 1 3 -1 -1 4 1 30]

numpy.where（暂时不看）

where(condition, [x, y])：

-- x, y不为空: condition为True, 返回x; False, 返回y
-- x, y为空： 返回condition为True的数组下标

np.where(关系表达式)：数组中知足关系表达式的元素的下标数组

numpy.take（暂时不看）

eg13_diff5.py

np.take(数组，下标数组)：数组中由下标数组所表示的元素集合
利用data.csv，对于收盘价，计算简单收益率及标准差、对数收益率及标准差、年波动率和月波动率

import numpy as np
c=np.loadtxt('data.csv',delimiter=',',usecols=(6,),\
             unpack=True,skiprows=1)
returns=np.diff(c)/c[:-1]
print('简单收益率',returns)
print('简单收益率的标准差:',np.std(returns))
logreturns=np.diff(np.log(c))
print('对数收益率：',logreturns)
print('对数收益率的标准差：',np.std(logreturns))
posretindices=np.where(returns>0)
print('全部正值元素的元素值:',np.take(returns.posretindices)) #有问题
annual_volatility=np.std(logreturns)/np.mean(logreturns)
annual_volatility=annual_volatility/np.sqrt(1/252)
print('年波动率',annual_volatility)
print('月波动率',annual_volatility*np.sqrt(1/12))