第十三周作业

今天老师讲了 Scipy 的用法，并要求我们做 Scipy 的课后练习。

第一题

随机生成一个 $m \times n (m \gt n)$ 的矩阵 $A$ 和一个 $m$ 维向量 $b$ ，求 $x = \arg \min_x ||Ax - b|| _2$

这时候就可以用 numpy 的 linalg ，专门处理线性代数问题的库。只要用 numpy.linalg.lstsq 函数，即可解出 $x$ 。

import numpy as np

if __name__ == "__main__":
    A = np.random.random((6, 4))
    b = np.random.random(6)
    x = np.linalg.lstsq(A, b)[0]
    print("A = ")
    print(A)
    print("")
    print("B = ")
    print(b)
    print("")
    print("x = ")
    print(x)

运行结果

A = 
[[0.94108162 0.34912602 0.01343064 0.25470686]
 [0.42072768 0.03906845 0.40904322 0.58866869]
 [0.95366445 0.92325452 0.03476652 0.53395704]
 [0.55189519 0.34647224 0.73896788 0.1849428 ]
 [0.32289582 0.13617389 0.59604818 0.42770833]
 [0.11844025 0.54501551 0.60428737 0.32446275]]

B = 
[0.04127154 0.4027248  0.9696257  0.81466285 0.11660395 0.91161533]

x = 
[-0.30132808  1.20957557  0.49407607  0.12972759]

第二题

求 $f(x) = sin^2(x - 2)e^{-x^2}$ 的最大值。

我们可以用scipy.optimize.minimize 函数来求函数的最小值。定义函数时，我们将 $f(x)$ 取相反数。

scipy.optimize.minimize中，必选参数有两个：fun, x0，分别表示函数和初值。

初值我设置为 0 ，由 Scipy 迭代出解。

import numpy as np
import scipy.optimize

if __name__ == '__main__':
    fun = lambda x : - np.power(np.sin(x - 2), 2) * np.exp(- np.power(x, 2))
    x0 = 0
    xmin = scipy.optimize.minimize(fun, x0)
    print(xmin)
    print("x = ")
    print(xmin.x)
    print("y = ")
    print(-xmin.fun)

运行结果

      fun: -0.9116854118471545
 hess_inv: array([[0.2680098]])
      jac: array([-1.49011612e-08])
  message: 'Optimization terminated successfully.'
     nfev: 18
      nit: 4
     njev: 6
   status: 0
  success: True
        x: array([0.21624132])
x = 
[0.21624132]
y = 
0.9116854118471545

可以看到，函数在 $x = 0.21624132$ 取得最大值 $0.9116854118471545$ ，结合上周作业画出的图像，容易确认这个解是最大值。

第三题

随机生成 n n <script type="math/tex" id="MathJax-Element-474">n</script> 个向量，求两两的距离。

我们可以使用 scipy.spatial.distance.euclidean 函数来求两个向量的欧几里得距离。

import numpy
import scipy.spatial.distance

if __name__ == '__main__':
    A = numpy.random.random((6, 4))
    B = numpy.zeros((6, 6))
    for i in range(0, A.shape[0]):
        for j in range(0, A.shape[0]):
            B[i][j] = scipy.spatial.distance.euclidean(A[i], A[j])
    print("A = ")
    print(A)
    print("")
    print("B = ")
    print(B)

运行结果

A = 
[[0.7099972  0.93121583 0.7913237  0.6047409 ]
 [0.27980001 0.32672134 0.07378547 0.28111213]
 [0.4690468  0.05131542 0.72802694 0.15770375]
 [0.96567853 0.67670334 0.39120523 0.61704902]
 [0.78609342 0.41034175 0.07176074 0.4306803 ]
 [0.73727238 0.50185982 0.21008417 0.31372443]]

B = 
[[0.         1.08170231 1.01790496 0.53888387 0.90838779 0.77950015]
 [1.08170231 0.         0.74493241 0.89806773 0.53450924 0.50950501]
 [1.01790496 0.74493241 0.         0.98091798 0.85709975 0.75335409]
 [0.53888387 0.89806773 0.98091798 0.         0.48987487 0.45557686]
 [0.90838779 0.53450924 0.85709975 0.48987487 0.         0.20873691]
 [0.77950015 0.50950501 0.75335409 0.45557686 0.20873691 0.        ]]