美文网首页
3.Python scipy笔记(进行科学计算)

3.Python scipy笔记(进行科学计算)

作者: Louis_Duan | 来源:发表于2018-12-04 10:11 被阅读0次

    1 简介

    1.1 定义

    scipy是一个以numpy为基础的用来进行科学计算的一个第三方库。

    1.2 构成

    子模块           描述
    cluster         聚类算法
    constants       物理数学常数
    fftpack         快速傅里叶变换
    integrate       积分和常微分方程求解
    interpolate     插值
    io              输入输出
    linalg          线性代数
    odr             正交距离回归
    optimize        优化和求根
    signal          信号处理
    sparse          稀疏矩阵
    spatial         空间数据结构和算法
    special         特殊方程
    stats           统计分布和函数
    weave           C/C++ 积分
    

    2 插值

    2.1 一维插值

    from scipy.interpolate import interp1d
    interp1d(x,y,kind='linear',...) 
    

    参数 x,y是一系列已知的点,kind参数是插值的类型

    候选值 作用
    ‘zero’、‘nearest’ 0阶插值、最近邻插值
    ‘slinear’、‘linear’ 线性插值,用一条直线连接所有取样点,相当于一阶B样条曲线
    ‘quadratic’、‘cubic’ 二次插值,更高阶的曲线可以直接使用整数值指定

    2.2 代码示例

    import numpy as np
    from scipy.interpolate import interp1d
    
    #创建待插值的数据
    x = np.linspace(0, 10*np.pi, 20)
    y = np.cos(x)
    
    # 分别用linear和quadratic插值
    fl = interp1d(x, y, kind='linear')
    fq = interp1d(x, y, kind='quadratic')
    
    #设置x的最大值和最小值以防止插值数据越界
    xint = np.linspace(x.min(), x.max(), 1000)
    yintl = fl(xint)
    yintq = fq(xint)
    
    作者:火锅侠
    链接:https://www.jianshu.com/p/b306095309db
    來源:简书
    简书著作权归作者所有,任何形式的转载都请联系作者获得授权并注明出处。
    

    2.3 径向基函数

    径向基函数,简单来说就是 x 处的函数值只依赖于 x 和 某点 c 的距离。

    2.3.1 常用的径向基(RBF)函数有:

    高斯函数

    x = np.linspace(-3,3,100)
    plt.plot(x,np.exp(-1*x**2))
    plt.title(''Gaussian'')
    
    Figure_1.png

    Mutiquadric函数(多元二次函数)

    plt.plot(x,np.sqrt(1+x**2))
    plt.title(''Multiquadric'')
    
    Figure_1.png

    Inverse quadratic函数(逆二次函数)

    plt.plot(x,1./np.sqrt(1+x**2))
    plt.title(''Inverse Multiquadric'')
    
    Figure_1.png

    2.3.2 径向基函数插值

    from scipy.interpolate.rbf import Rbf
    # 使用multiquadric的核
    cp_rbf = Rbf(x,y,function="multiquadric")
    # 使用gaussian的核
    cp_rbf=Rbf(x,y,function="gaussian")
    # 使用inverse_multiquadric的核
    cp_rbf=Rbf(x,y,function="inverse_multiquadric")
    

    2.3.3 高维RBF插值

    import numpy as np
    import pylab as pl
    from  mpl_toolkits.mplot3d import Axes3D
    from scipy.interpolate.rbf import Rbf
    
    # 3 维数据点
    x,y=np.mgrid[-np.pi/2:np.pi/2:5j, -np.pi/2:np.pi/2:5j]
    z = np.cos(np.sqrt(x**2+y**2))
    fig = pl.figure(figsize=(12,6))
    ax = fig.gca(projection="3d")
    ax.scatter(x,y,z)
    
    # 3 维RBF插值
    zz = Rbf(x,y,z)
    xx, yy = np.mgrid[-np.pi/2:np.pi/2:50j, -np.pi/2:np.pi/2:50j]
    fig = pl.figure(figsize=(12,6))
    ax = fig.gca(projection="3d")
    ax.plot_surface(xx,yy,zz(xx,yy),rstride=1,cstride=1,cmap=pl.cm.jet)
    
    pl.show()
    

    3 概率统计方法

    4 曲线拟合

    什么是曲线拟合:找一个方程来拟合数据。
    基础包

    import numpy as np
    import matplotlib as mpl
    import matplotlib.pyplot as plt
    from numpy import polyfit,poly1d
    # polyfit返回的是系数,poly1d返回的函数式
    

    4.1 多项式拟合

    从一阶到九阶拟合多项式拟合正弦函数

    import numpy as np
    import matplotlib.pyplot as plt
    from numpy import  polyfit,poly1d
    x = np.linspace(-np.pi,np.pi,100)
    y = np.sin(x)
    y1 = poly1d(polyfit(x,y,1))
    y3 = poly1d(polyfit(x,y,3))
    y5 = poly1d(polyfit(x,y,5))
    y7 = poly1d(polyfit(x,y,7))
    y9 = poly1d(polyfit(x,y,9))
    x = np.linspace(-3 * np.pi,3 * np.pi,100)
    
    p = plt.plot(x, np.sin(x), 'k')
    p = plt.plot(x, y1(x),label='y1')
    p = plt.plot(x, y3(x),label='y3')
    p = plt.plot(x, y5(x),label='y5')
    p = plt.plot(x, y7(x),label='y7')
    p = plt.plot(x, y9(x),label='y9')
    
    a = plt.axis([-3 * np.pi, 3 * np.pi, -1.25, 1.25])
    
    plt.legend()
    plt.show()
    
    Figure_1.png

    4.2 最小二乘拟合

    什么是最小二乘:已知多个近似分布于同一直线上的点,想拟合出一个直线方程:设该直线方程为y=kx+b,调整参数k和b,使得所有点到该直线的距离平方之和最小,设此时满足要求的k=k0,b=b0,则直线方程为y=k0x+b0。

    相关模块

    from scipy.linalg import lstsq#求最小二乘解
    from scipy.stats import linregress#线性回归
    from scipy.optimize import leastsq
    #(func,x0,args=())
    func是自己定义的一个计算误差的函数
    x0是计算的初始参数值
    args是指定func的其他参数
    

    用leastsq函数拟合数据,待拟合函数为kx+b

    ###最小二乘法试验###
    import numpy as np
    from scipy.optimize import leastsq
    
    ###采样点(Xi,Yi)###
    Xi=np.array([8.19,2.72,6.39,8.71,4.7,2.66,3.78])
    Yi=np.array([7.01,2.78,6.47,6.71,4.1,4.23,4.05])
    
    ###需要拟合的函数func及误差error###
    def func(p,x):
        k,b=p
        return k*x+b
    
    def error(p,x,y,s):
        print (s)
        return func(p,x)-y #x、y都是列表,故返回值也是个列表
    
    #TEST
    p0=[100,2]
    #print( error(p0,Xi,Yi) )
    
    ###主函数从此开始###
    s="Test the number of iteration" #试验最小二乘法函数leastsq得调用几次error函数才能找到使得均方误差之和最小的k、b
    Para=leastsq(error,p0,args=(Xi,Yi,s)) #把error函数中除了p以外的参数打包到args中
    k,b=Para[0]
    print("k=",k,'\n',"b=",b)
    
    ###绘图,看拟合效果###
    import matplotlib.pyplot as plt
    
    plt.figure(figsize=(8,6))
    plt.scatter(Xi,Yi,color="red",label="Sample Point",linewidth=3) #画样本点
    x=np.linspace(0,10,1000)
    y=k*x+b
    plt.plot(x,y,color="orange",label="Fitting Line",linewidth=2) #画拟合直线
    plt.legend()
    plt.show()
    
    Figure_1.png

    5 最小化函数

    ??

    6 积分

    ??

    7 解微分方程

    ??

    8 稀疏矩阵

    Scipy 提供了稀疏矩阵的支持(scipy.sparse)。稀疏矩阵主要使用 位置 + 值 的方法来存储矩阵的非零元素
    常用的存储格式有:

    • COO 格式在构建矩阵时比较高效
    • CSC 和 CSR 格式在乘法计算时比较高效
      例子:
    from scipy.sparse import *
    import numpy as np
    A = coo_matrix([[1,2,0],[0,0,3],[4,0,5]])
    print (A)
      (0, 0)    1
      (0, 1)    2
      (1, 2)    3
      (2, 0)    4
      (2, 2)    5
    前面是非零值的坐标,后面是值
    A.toarray()
    array([[1, 2, 0],
           [0, 0, 3],
           [4, 0, 5]])
    
    • sparse.find()
      返回一个三元组,表示稀疏矩阵中的非零元素的(行,列,值)
      例子:
    C 
      (0, 0)    3
      (0, 2)    1
      (1, 1)    2
      (3, 3)    1
    from scipy import sparse
    
    row, col, val = sparse.find(C)
    print (row, col, val)
    [0 0 1 3] [0 2 1 3] [3 1 2 1]
    
    • sparse.issparse()
      查看一个矩阵是否为稀疏矩阵,返回值为True或者False

    9 线性代数

    numpy 和 scipy 中,负责进行线性代数部分计算的模块叫做 linalg。

    • numpy.linalg(推荐使用)和scipy.linalg 区别:
      1.scipy.linalg包含scipy.linalg中所有函数,同时还包含了很多其没有的函数。
      2.scipy.linalg速度更快。

    线性代数的基本操作对象是矩阵

    numpy.matrix创建矩阵

    A = np.mat("[1, 2; 3, 4]")
    print (repr(A))
    matrix([[1, 2],
            [3, 4]])
    

    转置和逆

    A.I
    A.T
    

    矩阵乘法

    A*B
    

    2维numpy.ndarray创建矩阵

    A = np.array([[1,2], [3,4]])
    >> array([[1, 2],
           [3, 4]])
    

    linalg.inv(A)
    

    转置

    A.T
    

    矩阵乘法

    A.dot(b)
    

    普通乘法

    A*B
    

    求解线性方程组

    x+3y+5z=10
    2x+5y+z=8
    2x+3y+8z=3

    >>A = np.array([[1, 3, 5],
                  [2, 5, 1],
                  [2, 3, 8]])
    >>b = np.array([10, 8, 3])
    >>linalg.solve(A,b)
    [-9.28  5.16  0.76]
    

    计算行列式

    >>A = np.array([[1, 3, 5],
                  [2, 5, 1],
                  [2, 3, 8]])
    >>linalg.det(A)
    >>-25.0
    

    计算矩阵或向量的模

    x_norm=np.linalg.norm(x, ord=None, axis=None, keepdims=False)
    

    向量范数:

    参数 说明 计算方法
    ord=2 二范数 每个元素的平方的和再开方
    ord=1 一范数 每个元素的绝对值的和
    ord=np.inf 无穷范数 max(|xi|)
    • 矩阵的范数
      ord=1:列和范数,即所有矩阵列向量绝对值之和的最大值
      ord=2:谱范数,即ATA矩阵的最大特征值的开平方。
      ord=inf:行和范数,即所有矩阵行向量绝对值之和的最大值.
    • axis:处理类型
      axis=1表示按行向量处理,求多个行向量的范数
      axis=0表示按列向量处理,求多个列向量的范数
      axis=None表示矩阵范数。
    >>A = np.array([[1, 2],
                  [3, 4]])
    >>linalg.norm(A)
    import numpy as np
    x = np.array([
        [0, 3, 4],
        [1, 6, 4]])
    print linalg.norm(A)
    
    print linalg.norm(A,'fro') # frobenius norm 默认值
    
    print linalg.norm(A,1) # L1 norm 最大列和
    
    print linalg.norm(A,-1) # L -1 norm 最小列和
    
    print linalg.norm(A,np.inf) # L inf norm 最大行和
    

    伪逆

    linalg.pinv:使用最小二乘的解法求伪逆
    linalg.pinv2:使用求奇异值的算法求解
    

    特征值和特征向量

    linalg.eig(A):返回矩阵的特征值和特征向量
    linalg.eigvals(A):返回矩阵的特征值
    linalg.eig(A,B):?  +
    

    相关文章

      网友评论

          本文标题:3.Python scipy笔记(进行科学计算)

          本文链接:https://www.haomeiwen.com/subject/rarfqqtx.html