# 十，使用数组进行文件输入和输出
# Numpy 可以将数据以【文本】或【二进制】的形式存入硬盘，或从硬盘载入。
# 由于大部分用户更倾向于使用pandas等其他工具来载入文本或表格型数据，因此，这里只讲 Numpy 的内建二进制格式。
# 10.1，np.save(),np.savez()，np.savez_compressed()和 np.load() 是高效存取硬盘数据的工具函数。
# np.save(),数组在默认情况下是以未压缩的格式进行存储的，后缀名是 .npy
# np.savez(),用于保存单或多个数组，数组在默认情况下是以未压缩的格式进行存储的，后缀名是 .npz 必须以key-value对传入函数
# np.savez_compressed(),用于保存单个或多个数组，如果你的数据已经压缩好了，使用此方法将数据存入已经压缩好的文件中，后缀名是 .npz 必须以key-value对传入函数
# 示例1：保存单个数组 np.save(path,arr)，载入单个数组 np.load()：
arr=np.arange(10)
np.save('saved_arr',arr) #当文件的存放路径中没有写.npy时，后缀会被自动加上。
np.load('saved_arr.npy') #load方法时，不要漏写.npy
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
# 示例2：保存单个或多个数组 np.savez(path,key1=arr1,key2=arr2)，载入多个数组 np.load()：
# 用于在未压缩文件中保存单个或多个数组。
# np.savez('array_archive01',key=arr) #以字典型的对象存储 注意，key是变量的形式
# arch01=np.load('array_archive01.npz')
# arch01['key']
np.savez('array_archive02',key0=arr,key1=arr,key2=arr) #以字典型的对象存储 注意，key是变量的形式
arch02=np.load('array_archive02.npz')
arch02['key2']
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
# 示例3：np.savez_compressed(path,key1=arr1[,key2=arr2])，载入数组 np.load()：
# np.savez_compressed('array_archive_compressed01.npz',key=arr)
# arch_compressed01=np.load('array_archive_compressed01.npz')
# arch_compressed01['key']
np.savez_compressed('array_archive_compressed02.npz',key1=arr1,key2=arr2)
arch_compressed02=np.load('array_archive_compressed02.npz')
arch_compressed02['key1']
array([0, 1, 2, 3, 4, 5, 6, 7])
# 十一，线性代数
# 线性代数（比如矩阵的乘法，分解，行列式等方阵数学）是所有的数组类库的重要组成部分。
# 以下函数既可以作为数组实例的方法被调用，也可以作为Numpy的顶层函数使用。
# linear algebra functions 线性代数函数
# 1）nunpy命名空间中的线性代数函数：
# dot,vdot,inner,outer,matmul,tensordot,einsum,einsum_path and kron.
# 2）numpy.linalg中的线性代数函数：
# dig 将一个方阵的对角（或非对角）元素作为一维数组返回，或者将一维数组转换成一个方阵，并在非对角线上有零点。
# dot 这个应该是numpy的命名空间的函数，不属于numpy.linalg。
# 矩阵点乘 arr1.dot(arr2) 等效于 np.dot(arr1,arr2) 等效于 arr1 @ arr2
# trace 计算对角元素和
# det 计算矩阵的行列式
# eig 计算矩阵的特征值和特征向量
# inv 计算矩阵的逆矩阵
# pinv 计算矩阵的Moore_Penrose伪逆
# qr 计算QR分解
# svd 计算奇异值分解
# solve 求解 x的线性系统 Ax=b,其中A是方阵
# lstsq 计算Ax=b的最小二乘解
# 示例1 矩阵的点乘 xarr.dot(yarr) 等效于 np.dot(xarr,yarr) 等效于 xarr @ yarr
x=np.arange(1,7,dtype='float64').reshape((2,3)) # x=np.arange(1,7).reshape((2,3)) x.astype('float64')
y=np.array([[6.,23.],[-1,7],[8,9]])
z=np.ones(3)
x
array([[1., 2., 3.],
[4., 5., 6.]])
y
array([[ 6., 23.],
[-1., 7.],
[ 8., 9.]])
# x.dot(y)
np.dot(x,y)
# x @ y #或者使用 @
array([[ 28., 64.],
[ 67., 181.]])
# 一个二维数组与一个长度合适的一维数组之间的矩阵点乘，其结果是一个一维数组： （2，3）点乘（3，）
# z.shape # (3,1)
x.dot(z)
# x @ np.ones(3) #或者使用 @
array([ 6., 15.])
# 示例2 numpy.linalg 的标准函数集中 inv,qr的用法：
from numpy.linalg import inv,qr
X=np.random.randn(5,5)
mat=X.T.dot(X)
inv(mat)
array([[ 0.65012083, 0.39118584, 1.44000575, 0.62332612, 0.03844465],
[ 0.39118584, 0.96791675, 0.60758529, 0.72360998, -0.14917781],
[ 1.44000575, 0.60758529, 6.03307846, 1.93675845, 0.4891006 ],
[ 0.62332612, 0.72360998, 1.93675845, 1.28105523, 0.01168604],
[ 0.03844465, -0.14917781, 0.4891006 , 0.01168604, 0.15744843]])
mat.dot(inv(mat))
array([[ 1.00000000e+00, -4.82092301e-16, 2.23523512e-15,
7.50355993e-17, 1.37572727e-14],
[-3.01124445e-15, 1.00000000e+00, -1.53712953e-14,
-1.88230336e-15, 8.07520682e-15],
[-1.79688647e-15, -1.87859847e-15, 1.00000000e+00,
-3.97681915e-16, 1.45616874e-14],
[-5.93572733e-15, 1.04975882e-14, -3.67122121e-15,
1.00000000e+00, -3.47918402e-14],
[-3.17375114e-15, 4.27233867e-15, 2.03809481e-15,
-1.13216349e-15, 1.00000000e+00]])
q,r=qr(mat)
r
array([[ -5.77811439, 0.20538457, -12.87202449, 3.16352559,
3.16100464],
[ 0. , -7.02156113, -14.02550707, -0.01789824,
-1.75269181],
[ 0. , 0. , -1.34826783, -2.81469224,
0.54679839],
[ 0. , 0. , 0. , -3.68656824,
1.78405193],
[ 0. , 0. , 0. , 0. ,
0.09874386]])
# 十二，伪随机数生成
# numpy.random 模块填补了python模块的不足，可以高效的生成多种概率分布下的完整样本值数组。
# numpy.random 中的部分函数列表：
# seed 像随机数生成器传递随即状态种子
# permutation 返回一个序列的随机排列，或者返回一个乱序的整数范围序列
# shuffle 随机排列一个序列
# rand 从均匀分布中抽取样本
# randint 根据给定的由低到高的范围抽取随机整数 np.random.randint(start,end,size=n),n是想生成的一维随机数组的大小
# randn 从均值0方差1的正态分布中抽取样本 （MATLAB 型接口）
# binomial 从二项分布中抽取样本
# normal 从正态（高斯）分布中抽取样本
# beta 从beta分布中抽取样本
# chisquare 从卡方分布中抽取样本
# python 中random.randint()用法： random.randint(0,1) 表示从0和1随机选取一个整数，
# random.randint(0,5) 表示从0，1，2，3，4，5中随机选取一个整数
# 如，可以使用normal 来获得一个4x4的正态分布样本数组：
samples=np.random.normal(size=(4,4)) # 注意参数的写法
samples
array([[ 4.71435164e-01, -1.19097569e+00, 1.43270697e+00,
-3.12651896e-01],
[-7.20588733e-01, 8.87162940e-01, 8.59588414e-01,
-6.36523504e-01],
[ 1.56963721e-02, -2.24268495e+00, 1.15003572e+00,
9.91946022e-01],
[ 9.53324128e-01, -2.02125482e+00, -3.34077366e-01,
2.11836468e-03]])
# 我们称这些数为伪随机数，因为他们是由有确定性行为的算法根据随机数生成器中的随机种子生成的。
# 可以通过 np.random.seed 来更改 Numpy的随机数种子：
np.random.seed(1234) # 设置一个全局的随机数种子 1234，之后所有的随机数都与此种子有关系。
# 如果想要得到一个独立于上述设定的随机数种子，那么，需要使用 numpy.random.RandomState(1234) 创建一个独立的随机数生成器对象，
# 然后使用该对象的random()方法来创建随机数：
rng=numpy.random.RandomState(1234)
rng.random(10)
array([0.19151945, 0.62210877, 0.43772774, 0.78535858, 0.77997581,
0.27259261, 0.27646426, 0.80187218, 0.95813935, 0.87593263])
# 十三： 示例 随机漫步
# 11.1 模拟简单的随机漫步
# 通过这个例子来讲解如何使用数组的运算：
# 从0开始，步进为 1 和 -1，且两种步进发生的概率相等。
# 1）使用random 模块来利用纯python 实现一个1000步的随机漫步：
import random
position = 0
walk = [position]
steps = 1000
for i in range(steps):
# if random.randint(0,1):
# step = 1
# else:
# step = -1
step = 1 if random.randint(0,1) else -1 # python中三元表达式： x if condition else y
position += step
walk.append(position)
plt.plot(walk[:100]) #将前100步画出来
[<matplotlib.lines.Line2D at 0x1f53b25c308>]
# 2) 使用 数组的运算 来模拟一个1000步的随机漫步：
# 上面的walk 只是对随机步进的累积，并且可以通过一个数组表达式实现。因此，可以使用 np.random 一次性生成 1000 个随机 1，0数字的一维数组
nsteps=1000
draws=np.random.randint(0,2,size=nsteps) # 生成 元素为 0，1随机数的一维数组，长度1000。
steps=np.where(draws>0,1,-1) # 转化为1000个1，-1随机数
walk=steps.cumsum() # 得到一个长度与steps一致的一维数组， 非聚合函数 所有元素的累积和，第一个加数为0 见 7.3，数学和统计方法 聚合函数和其他函数
print(walk.min())
print(walk.max())
plt.plot(walk[:100]) #将前100步画出来
-2
49
<matplotlib.lines.Line2D at 0x1f53b32cec8>]
array([0, 1, 1, 1, 1, 1, 0, 0, 0, 1])
# 复杂一点的统计，随机漫步的某一步达到了某个特定值。这里假设我们想知道漫步中是何时连续朝某个方向走了10步。
# np.abs(walk)>=10给我们一个布尔值数组，用于表明漫步是否在同一个方向累计走了至少十步（距离是否超过10），
# 而我们想知道第一个超过10或-10的位置，此时可以使用argmax来计算，此函数可以用在数组（既可以是数值型数组，也可以是布尔类型数组）
# 中查找最大值，如果最大值重复，取第一个，用在布尔类型数组中，就是返回布尔值数组中True的第一个位置.
(np.abs(walk>=10)).argmax() # 随机漫步朝一个方向的最大距离所在的数组中第一个位置，也就是最早穿越10的位置，可以认为是时间
35
# 3) 一次性模拟多次随机漫步
# 如果我们的目标是模拟多次随机漫步（比如，5000个），我们可以稍微的修改上面的代码。
# 只需要给numpy.random函数的size参数传入一个二元元组，就可以产生一个二维数组，就可以一次性模拟5000个随机漫步过程。
nwalks=5000
nsteps=1000
draws=np.random.randint(0,2,size=(nwalks,nsteps)) # 生成 元素为 0，1随机数的二维数组，shape是（5000，1000）。
steps=np.where(draws>0,1,-1)
walks=steps.cumsum(1) # axis=1,表示按列的方向对行操作。
walks[:5]
array([[ -1, -2, -3, ..., 28, 29, 30],
[ 1, 0, 1, ..., 14, 15, 16],
[ 1, 0, -1, ..., -28, -27, -28],
[ -1, 0, -1, ..., -24, -23, -22],
[ 1, 0, -1, ..., -6, -5, -6]], dtype=int32)
print(walks.max()) #随机漫步朝一个方向的最大距离
print(walks.min()) #随机漫步朝一个方向的最小距离
131
-141
# 计算出这5000个随机漫步中，朝一个方向距离超过30的，最小穿越时间。
# 由于不是所有的随机漫步都能够有超过30的，所以可以使用any来检查, any()用法见 7.4，布尔值数组的方法 既可以作为数组实例的方法被调用，也可以作为顶级numpy函数使用。
hits30=(np.abs(walks)>=30).any(1) # 1表示对每一行进行查找
hits30[:5] # 前五行中有没有达到30的
array([ True, False, True, False, True])
hits30.sum() # 5000个随机漫步中，有多少个随机漫步至少有一次达到了30步及以上。
3427
# 将最早穿越30距离的时间放在数组里：
crossing_times=(np.abs(walks[hits30])>=30).argmax(1) # walks[hits30]表示walk中达到30次的随机漫步，布尔索引。 返回一维数组
crossing_times[:20]
array([925, 803, 395, 583, 881, 583, 705, 737, 445, 269, 427, 299, 961,
267, 619, 583, 859, 261, 319, 755], dtype=int64)
# 求个穿越30距离的平均时间：
crossing_times.mean()
502.4765100671141