前言

  • 买了一本迪哥写的《人工智能数学基础》，在掌握了 Python 的基本用法后，着手学习其中的内容
  • 里面很多都是考研数学的内容，让我意识到数学没有想象中的那么“没用”
  • Python 有很多封装好的数学工具，正文中有大量的数学公式使用了 LateX 的语句，参考在线 LaTeX 公式编辑器-编辑器 (latexlive.com)，很遗憾没能在本科期间认识到它们
  • 有些 LateX 语句并不能正确地渲染, 先凑合着用吧

1 基础篇

2 高等数学基础

例 2.6 求limx→1x2−1x−1\lim_{x \to 1}\frac{x ^ 2 - 1}{x - 1} 的极限

1
2
3
4
5
6
import sympy
from sympy import oo
import numpy as np
x = sympy.Symbol('x')
f = (x ** 2 - 1) / (x - 1)
sympy.limit(f, x, 1)  # limit(数学表达式 expression, 变量 variable, 要趋向的值 value)

2\displaystyle 2

例 2.11 求y=arcsinsinxy = arcsin\sqrt{sinx}的导数

y′=11−sinx⋅12sinx⋅cosx=cosx2sinx−sin2xy'=\frac{1}{\sqrt{1-sinx} }\cdot\frac{1}{2\sqrt{sinx}}\cdot{cosx}=\frac{cosx}{2\sqrt{sinx-sin^2x}}

1
2
3
from sympy import *
from sympy.abc import x  # 导入变量 x
diff(asin(sqrt(sin(x))))  # diff 求导函数

cos(x)21−sin(x)sin(x)\displaystyle \frac{\cos{\left(x \right)}}{2 \sqrt{1 - \sin{\left(x \right)}} \sqrt{\sin{\left(x \right)}}}

例 2.12 求f(x,y)=x2+3xy+y2f(x,y)=x^2+3xy+y^2在点(1,2)(1,2)处的偏导数

fx(x,y)=2x+3yf_{x}(x,y)=2x+3y

fy(x,y)=3x+2yf _{y}(x,y)=3x+2y

$ f _{x}(1,2)=2x+3y| _{y=2}^{x=1}=8$

fy(1,2)=3x+2y∣y=2x=1=7f _{y}(1,2)=3x+2y| _{y=2}^{x=1}=7

1
2
3
4
from sympy import *
from sympy.abc import x, y, f
f = x ** 2 + 3 * x * y + y ** 2
diff(f, x)  # 对 x 求偏导

2x+3y\displaystyle 2 x + 3 y

1
diff(f, y)

3x+2y\displaystyle 3 x + 2 y

1
2
fx = diff(f, x)
fx.evalf(subs={x:1, y:2})  # 以字典的形式传入多个变量的值, 返回计算后的数学表达式。

8.0\displaystyle 8.0

1
2
fy = diff(f, y)
fy.evalf(subs={x:1, y:2})

7.0\displaystyle 7.0

2.6 方向导数 2.7 梯度

函数在某点的梯度是一个向量, 它的方向与取得最大方向导数的方向一致

$gradf(x,y)=\frac{\partial f}{\partial x}\textbf{i}+\frac{\partial f}{\partial y}\textbf{j} $

\frac{\partial f}{\partial l}=gradf\cdot \frac{\overrightarrow{l}}{\left | \overrightarrow{l} \right | }

用梯度gradf点乘l的单位向量就得到方向导数, 这是计算方向导数最简便的方法

例 2.15 用 Python 编程实现梯度下降法求解f(x,y)=x−y+2x2+2xy+y2f(x,y)=x-y+2x^2+2xy+y^2的最小值

使用常规方法:

z=x−y+2x2+2xy+y2z=x-y+2x^2+2xy+y^2

zx′=1+4x+2yz'_{x}=1+4x+2y

zy′=−1+2x+2yz'_{y}=-1+2x+2y

求得驻点(−1,32)(-1,\frac{3}{2})

A=zxx′′∣y=32x=−1=4A=z''_{xx}| _{y=\frac{3}{2}}^{x=-1}=4

B=zxy′′∣y=32x=−1=2B=z''_{xy}| _ {y=\frac{3}{2}}^{x=-1}=2

C=zyy′′∣y=32x=−1=2C=z''_{yy}| _ {y=\frac{3}{2}}^{x=-1}=2

AC−B2>0AC-B^2>0

A>0A>0

原函数在(−1,32)(-1,\frac{3}{2})处具有极小值

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np


def Fun(x, y):  # 原函数
    return x - y + 2 * x * x + 2 * x * y + y * y


def PxFun(x, y):  # 对 x 求偏导
    return 1 + 4 * x + 2 * y


def PyFun(x, y):  # 对 y 求偏导
    return -1 + 2 * x + 2 * y


if __name__ == "__main__":
    fig = plt.figure()  # 创建自定义图象
    # MatplotlibDeprecationWarning: Axes3D(fig) adding itself to the figure is deprecated since 3.4. 
    # Pass the keyword argument auto_add_to_figure=False and use fig.add_axes(ax) to suppress this warning. 
    # The default value of auto_add_to_figure will change to False in mpl3.5 and True values will no longer work in 3.6.
    ax = Axes3D(fig, auto_add_to_figure=False)  # 创建 3D 图象
    fig.add_axes(ax)
    # 在绘制三维图表时, 需要用到 mgrid 函数, 它会返回一个密集的多维网格
    # np.mgrid[开始坐标:结束坐标:步长(步长为复数表示点数，左闭右闭, 步长为实数表示间隔，左闭右开)]
    X, Y = np.mgrid[-2:2:40j, -2:2:40j]  
    Z = Fun(X, Y)
    # ax.plot_surface: https://blog.csdn.net/weixin_43584807/article/details/102331755
    # X, Y, Z 2D 数组形式的数据值
    # rstride 数组行距
    # cstride 数组列距
    # cmap 曲面块颜色映射
    ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap="rainbow")
    ax.set_xlabel('x')  # 设置坐标轴标签
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    # 梯度下降
    step = 0.008  # 取步长
    x = 0  # 初始值
    y = 0
    tag_x = [x]  # 绘制图像的列表
    tag_y = [y]
    tag_z = [Fun(x, y)]
    new_x = x
    new_y = y
    Over = False
    while not Over:
        new_x -= step * PxFun(x, y)  # 往最大方向减少
        new_y -= step * PyFun(x, y)
        if Fun(x, y) - Fun(new_x, new_y) < 7e-9:
            Over = True
        x = new_x  # 更新数据
        y = new_y
        tag_x.append(x)  # 添加用于绘制图象的点
        tag_y.append(y)
        tag_z.append(Fun(x, y))
    ax.plot(tag_x, tag_y, tag_z, 'r+]')  # 绘制点: 'r'表示红色(好像没有作用?)
    print('(x, y)~(' + str(x) + ',' + str(y) + ')')  # 输出结果
    plt.show()  # 显示图象
    
(x, y)~(-0.9993608094022046,1.498965767887478)

png

2.9.5 高手点拨: 求导的三种方式

已知f(x)=x5+2x4+3x2+5f(x)=x^5+2x^4+3x^2+5, 求f′(1)f'(1)

使用 Sympy 的 diff 函数
1
2
3
4
5
6
import sympy
from sympy.abc import x, f

f = x ** 5 + 2 * x ** 4 + 3 * x ** 2 + 5
fx = sympy.diff(f)
fx.evalf(subs={x:1})

19.0\displaystyle 19.0

使用 scipy.misc 模块下的 derivative 函数
1
2
3
4
5
6
7
8
import numpy as np
from scipy.misc import derivative


def f(x):
    return x ** 5 + 2 * x ** 4 + 3 * x ** 2 + 5

derivative(func=f, x0=1, dx=1e-6, n=1)  # 函数为 f, 取值为 1, 间距为 1e-6, 一阶导数
18.999999999991246
使用 NumPy 模块里的 poly1d 构造f(x)f(x)
1
2
3
4
import numpy as np

p = np.poly1d([1, 2, 0, 3, 0, 5])  # 构造多项式
np.polyder(p, 1)(1.0)  # 一阶导数在 x=1.0 时的取值
19.0

1
p.deriv(1)(1.0)  # 一阶导数在 x=1.0 时的取值
19.0

2.10 习题

limx→1sin(lnx)lim_{x\to1}sin(lnx)
1
2
3
4
5
import sympy
from sympy.abc import x, f

f = sympy.sin(sympy.log(x))
sympy.limit(f, x, 1)

0\displaystyle 0

limx→8x3−2x−8lim_{x\to8}\frac{\sqrt[3]{x}-2}{x-8}
1
2
3
4
5
import sympy
from sympy.abc import x, f

f = (x ** (1/3) - 2) / (x - 8)
sympy.limit(f, x, 8)

112\displaystyle \frac{1}{12}

求y=x4−2x3+5sinx+ln3y=x^4-2x^3+5sinx+ln3的导数
1
2
3
4
5
import sympy
from sympy.abc import x, y

y = x ** 4 - 2 * x ** 3 + 5 * sympy.sin(x) + sympy.log(3)
diff(y)

4x3−6x2+5cos(x)\displaystyle 4 x^{3} - 6 x^{2} + 5 \cos{\left(x \right)}

求z=(3x2+y2)4x+2yz=(3x^2+y^2)^{4x+2y}在点(1,2)(1,2)的偏导数

lnz=(4x+2y)⋅ln(3x2+y2)lnz=(4x+2y)\cdot ln(3x^2+y^2)

1zdz=[4ln(3x2+y2)+24x2+12xy3x2+y2]dx+[2ln(3x2+y2)+8xy+2y23x2+y2]dy\frac1zdz=\left[4ln(3x^2+y^2)+\frac{24x^2+12xy}{3x^2+y^2}\right]dx+\left[2ln(3x^2+y^2)+\frac{8xy+2y^2}{3x^2+y^2}\right]dy

1
2
3
4
5
6
import sympy
from sympy.abc import x, y, z

z = (3 * x ** 2 + y ** 2) ** (4 * x + 2 * y)
zx = sympy.diff(z, x)
zx

(3x2+y2)4x+2y(6x(4x+2y)3x2+y2+4log(3x2+y2))\displaystyle \left(3 x^{2} + y^{2}\right)^{4 x + 2 y} \left(\frac{6 x \left(4 x + 2 y\right)}{3 x^{2} + y^{2}} + 4 \log{\left(3 x^{2} + y^{2} \right)}\right)

1
zx.evalf(subs={x:1, y:2})

84401203.0927369\displaystyle 84401203.0927369

1
2
zy = sympy.diff(z, y)
zy

(3x2+y2)4x+2y(2y(4x+2y)3x2+y2+2log(3x2+y2))\displaystyle \left(3 x^{2} + y^{2}\right)^{4 x + 2 y} \left(\frac{2 y \left(4 x + 2 y\right)}{3 x^{2} + y^{2}} + 2 \log{\left(3 x^{2} + y^{2} \right)}\right)

1
zy.evalf(subs={x:1, y:2})

48788945.5463684\displaystyle 48788945.5463684

求方向导数和梯度

求函数z=x2+y2z=x^2+y^2在点(1,2)(1,2)处沿点(1,2)(1,2)到点(2,2+3)(2,2+\sqrt{3})方向的方向导数, 以及在点(1,2)(1,2)的梯度

1
2
3
4
5
6
7
8
import sympy as sp
import numpy as np
from sympy.abc import x, y, z

z = x ** 2 + y ** 2
zx, zy = z.diff(x), z.diff(y)
gradz = np.array([zx.evalf(subs={x:1, y:2}), zy.evalf(subs={x:1, y:2})], dtype=float)  # 求梯度
gradz
array([2., 4.])

1
2
3
A = np.array([1, 2])
B = np.array([2, 2 + 3 ** 0.5], dtype=float)
gradz * (B - A) / np.linalg.norm(B - A)  # 计算方向导数
array([1.        , 3.46410162])

3 微积分

例 3.6 定积分

应用 SciPy 科学计算库求∫03cos2(ex)dx\int_{0}^{3} cos^2(e^x)dx

1
2
3
4
5
import numpy as np
from scipy.integrate import quad

func = lambda x:np.cos(np.exp(x)) ** 2  # 定义被积分函数
quad(func, 0, 3)  # 调用 quad 积分函数
(1.296467785724373, 1.397797133112089e-09)

输出结果(积分值, 误差)

例 3.7 二重积分

求∬De−xx−y2dxdy\iint_{D}e^{-x^x-y^2}dxdy, 其中D={(x,y)∣0≤x≤10,0≤y≤10}D=\left \{ (x,y) | 0 \le x \le 10,0 \le y \le 10\right \}

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
import numpy as np
from scipy.integrate import dblquad  # 二重积分


def integrand(x,y):
    return np.exp(-x ** 2 - y ** 2)

x_a = 0
x_b = 10
y_a = 0
y_b = 10
"""
scipy.integrate.dblquad(func, a, b, gfun, hfun, args=(), epsabs=1.49e-08, epsrel=1.49e-08)
参数：
func：可调用的
至少有两个变量的 Python 函数或方法：y 必须是第一个参数，x 必须是第二个参数。

a, b：浮点数
x:a < b 的积分极限

gfun：可调用或浮点数
y 中的下边界曲线，它是一个函数，采用单个浮点参数 (x) 并返回浮点结果或表示恒定边界曲线的浮点数。

hfun：可调用或浮点数
y 中的上边界曲线(与 gfun 的要求相同)。

args：顺序，可选
传递给 func 的额外参数。

epsabs：浮点数，可选
绝对容差直接传递到内部一维正交积分。默认值为 1.49e-8。dblquad 尝试获得 abs(i-result) <= max(epsabs, epsrel*abs(i)) 的精度，其中 i = 从 gfun(x) 到 hfun(x) 的 func(y, x) 的内部积分，而 result 是数值近似值。请参阅下面的 epsrel。

epsrel：浮点数，可选
内部一维积分的相对容差。默认值为 1.49e-8。如果 epsabs <= 0, epsrel 必须大于 5e-29 和 50 * (machine epsilon).看易胜宝更多。

返回：
y：浮点数
结果积分。

abserr：浮点数
误差的估计。
"""
dblquad(integrand, x_a, x_b, lambda x:y_a, lambda x:y_b)
(0.7853981633974476, 1.375309851021853e-08)

例 3.8 定积分近似求解

用定义法求∫03cos2(ex)dx\int_{0}^{3} cos^2(e^x)dx的近似解

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
from numpy import *

a, b = 0, 3


def f(x):
    return cos(exp(x)) ** 2


def trape(n):
    h = (b - a) / n  # 将区间分成 n 等分后的长度
    x = a
    sum = 0
    for i in range(1, n):
        x2 = a + i * h
        sum += (f(x) + f(x2)) * h / 2  # 小梯形的值(上底加下底的和乘高除以二)
        x = x2
    return sum
1
trape(10)
0.944822326405313

1
trape(100)
1.2843391540917448

1
trape(1000)
1.2960750567338157

不定积分

$ \int lnx $

1
2
3
4
5
from sympy import *
from sympy.abc import x

expr = log(x)
integrate(expr,x)

xlog(x)−x\displaystyle x \log{\left(x \right)} - x

3.8 习题

∫12x2+1x4dx\int _{1}^{2} x^2 + \frac {1}{x^4}dx

$ =(\frac {1}{3}x^3 - \frac {1}{3}x^{-3}) | _{1}^{2}$

$ = \frac {21}{8}$

1
2
3
4
5
import numpy as np
from scipy.integrate import quad

func = lambda x: x ** 2 + x ** (-4)
quad(func, 1, 2)
(2.625, 2.914335439641036e-14)

$ \int _{-1}^{0}\frac {3x4+3x2+1}{x^2+1}dx$

$ =\int _{-1}{0}3x2+\frac {1}{1+x^2}dx$

$ = (x^3 + arctanx)| _{-1}^{0} $

$ = 1 + \frac{\pi}{4} $

1
2
3
4
5
import numpy as np
from scipy.integrate import quad

func = lambda x: (3 * x ** 4 + 3 * x ** 2 + 1) / (1 + x ** 2)
quad(func, -1, 0)
(1.7853981633974483, 1.9821901491273144e-14)

利用定积分的定义计算极限:

$ lim _{n \to \infty}\frac {1{p}+2{p}+…+n{p}}{n{p+1}}$

$ = lim {n \to \infty}\frac {\sum{i=1}{n}ip}{n^{p+1}} $

$ = \frac{1}{n} lim {n \to \infty}\sum{i=1}{n}(\frac{i}{n})p $

$ = \int _{0}{1}xpdx $

$ = \frac {1}{p+1} $

4 泰勒公式与拉格朗日乘子法

指对连, 三角断,
三角对数隔一换, 三角指数有感叹

例 4.7 根据exe^x的nn次泰勒多项式展开式, 计算ee的近似值

$ e\approx 1+1+\frac{1}{2!}+\frac{1}{3!}+…+\frac{1}{n!}$

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
import numpy as np
import pandas as pd


def f(n):
    sum1 = 1
    if n == 0:
        sum1 = 1
    else:
        m = n + 1
        for i in range(1, m):
            sum2 = 1.0
            k = i + 1
            for j in range(1, k):
                sum2 = sum2 * j
            sum1 = sum1 + 1.0 / sum2
    return sum1


num = 10
pd.DataFrame(np.array([[i,f(i)] for i in range(1, num + 1)]), columns=['n', 'e'])
n e
0 1.0 2.000000
1 2.0 2.500000
2 3.0 2.666667
3 4.0 2.708333
4 5.0 2.716667
5 6.0 2.718056
6 7.0 2.718254
7 8.0 2.718279
8 9.0 2.718282
9 10.0 2.718282

例 4.13 sinxsinx的nn阶泰勒多项式

sinx=x−x33!+x55!−x77!+...+(−1)m−1x2m−1(2m−1)!+R2m(x)sinx = x - \frac {x^3}{3!} + \frac {x^5}{5!} - \frac {x^7}{7!} + ...+ (-1)^{m-1}\frac {x^{2m-1}}{(2m-1)!}+R_{2m}(x)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
import numpy as np
import pandas as pd


def fsin(x):
    m = 20  # 项数
    sum = 0.0
    for i in range(1, m+1):
        n = 2 * i - 1  # 2m-1
        temp1, temp2, temp3 = 1, 1, 1
        for j in range(1, i): # (-1)^(m-1)
            temp1 = -temp1
        for j in range(1, n + 1):
            temp2 = temp2 * x  # x^(2m-1)
            temp3 = temp3 * j  # (2m-1)!
        sum += temp1 * temp2 / temp3
    return sum


pd.DataFrame({'np.sin(x)': np.array([np.sin(x) for x in range(-20, 1)]),
              'fsin(x)': np.array([fsin(x) for x in range(-20, 1)]), 
              'error': np.array([fsin(x) - np.sin(x) for x in range(-20, 1)])},
              index=np.array([x for x in range(-20, 1)]))
np.sin(x) fsin(x) error
-20 -0.912945 5364.411846 5.365325e+03
-19 -0.149877 666.994385 6.671443e+02
-18 0.750987 74.739042 7.398806e+01
-17 0.961397 8.185042 7.223645e+00
-16 0.287903 0.899283 6.113793e-01
-15 -0.650288 -0.606249 4.403901e-02
-14 -0.990607 -0.987967 2.640574e-03
-13 -0.420167 -0.420039 1.282664e-04
-12 0.536573 0.536578 4.880595e-06
-11 0.999990 0.999990 1.394300e-07
-10 0.544021 0.544021 2.831679e-09
-9 -0.412118 -0.412118 3.790196e-11
-8 -0.989358 -0.989358 2.858824e-13
-7 -0.656987 -0.656987 3.441691e-15
-6 0.279415 0.279415 3.497203e-15
-5 0.958924 0.958924 -2.775558e-15
-4 0.756802 0.756802 -6.661338e-16
-3 -0.141120 -0.141120 1.110223e-16
-2 -0.909297 -0.909297 0.000000e+00
-1 -0.841471 -0.841471 0.000000e+00
0 0.000000 0.000000 0.000000e+00

4.9 习题

limx→0(sinx−xcosxsin3x)lim _{x \to 0}(\frac {sinx-xcosx}{sin^3x})

=limx→0x−16x3−x(1−12x2)x3=lim _{x \to 0}\frac {x-\frac{1}{6}x^3-x(1-\frac{1}{2}x^2)}{x^3}

=limx→0−16x3+12x3x3=lim _{x \to 0}\frac {-\frac{1}{6}x^3+\frac{1}{2}x^3}{x^3}

=13=\frac{1}{3}

1
2
3
4
5
6
import sympy
from sympy import sin, cos, limit

x = sympy.Symbol('x')
f = (sin(x) - x * cos(x)) / (sin(x) ** 3)
limit(f, x, 0)

13\displaystyle \frac{1}{3}

Icon
Icon
Icon
Icon
上一篇：Python-人工智能数学基础(5-6)
下一篇：Python-莫烦python学习笔记（multithreading）