Numpy习题

Generate matrices A, with random Gaussian entries, B, a Toeplitz matrix, where A ∈ R

n×m and B ∈ Rm×m,for n = 200, m = 500.

import numpy as np
from scipy.linalg import toeplitz

n_size = 200
m_size = 500
mu = 0.0
sigma = 0.1
np.random.seed(0)
A = np.random.normal(mu, sigma, (n_size, m_size))
print(A)
print(A.shape)
base = list(range(1, m_size + 1))
# print(base)
B = toeplitz(base, base)
print(B)
print(B.shape)

输出：

[[ 0.17640523  0.04001572  0.0978738  ..., -0.0514234  -0.10180419
  -0.00778548]
 [ 0.03827324 -0.00342423  0.10963468 ...,  0.00941923 -0.11476109
  -0.03581141]
 [ 0.05559627  0.08924739 -0.04223148 ..., -0.02467936 -0.10793432
  -0.01142256]
 ...,
 [-0.16058274 -0.035732   -0.09369679 ..., -0.1123346  -0.06918121
   0.05208148]
 [-0.00682909 -0.09105953  0.01875499 ..., -0.0513064   0.00768273
  -0.09786463]
 [ 0.02218512 -0.09420463 -0.0277842  ...,  0.0726182  -0.11678305
  -0.12852076]]
(200, 500)
[[  1   2   3 ..., 498 499 500]
 [  2   1   2 ..., 497 498 499]
 [  3   2   1 ..., 496 497 498]
 ...,
 [498 497 496 ...,   1   2   3]
 [499 498 497 ...,   2   1   2]
 [500 499 498 ...,   3   2   1]]
(500, 500)

Exercise 9.1: Matrix operations

Calculate A + A, AA^T, A^TA and AB. Write a function that computes A(B − λI) for any λ.

# 计算A + A
sum = A + A
print('sum:\n', sum)

# 计算A * （A^T）
res1 = np.dot(A, A.T)
print('res1:\n', res1)

# 计算（A^T) * A
res2 = np.dot(A.T, A)

# 计算A * B
res3 = np.dot(A, B)
print('res3:\n', res3)


def func(A, B, s, size):  # 计算A(B − λI)
    return np.dot(A, (B - s * np.identity(size)))


result = func(A, B, 0.5, m_size)
print('result:\n', result)

输出：

sum:
 [[ 0.35281047  0.08003144  0.1957476  ..., -0.10284679 -0.20360838
  -0.01557095]
 [ 0.07654649 -0.00684846  0.21926937 ...,  0.01883846 -0.22952219
  -0.07162282]
 [ 0.11119254  0.17849478 -0.08446296 ..., -0.04935871 -0.21586863
  -0.02284511]
 ...,
 [-0.32116548 -0.07146399 -0.18739357 ..., -0.22466921 -0.13836242
   0.10416295]
 [-0.01365817 -0.18211906  0.03750997 ..., -0.1026128   0.01536546
  -0.19572927]
 [ 0.04437024 -0.18840927 -0.05556839 ...,  0.1452364  -0.2335661
  -0.25704153]]
res1:
 [[ 4.984799   -0.19033257 -0.22672115 ...,  0.10982337  0.33523647
   0.2211457 ]
 [-0.19033257  4.77802726 -0.05226853 ...,  0.1097388   0.19826363
  -0.20766612]
 [-0.22672115 -0.05226853  4.66202048 ..., -0.24424636  0.10341813
   0.07524698]
 ...,
 [ 0.10982337  0.1097388  -0.24424636 ...,  4.78065551 -0.21885776
  -0.1114218 ]
 [ 0.33523647  0.19826363  0.10341813 ..., -0.21885776  5.18029085
  -0.29804843]
 [ 0.2211457  -0.20766612  0.07524698 ..., -0.1114218  -0.29804843
   4.28037337]]
res3:
 [[ -770.16766139  -768.54712896  -766.84656508 ...,   137.33964977
    136.29110713   135.03895612]
 [  -21.2917197    -17.95722443   -14.62957762 ..., -1604.79749053
  -1607.75429431 -1610.94062028]
 [  107.89619671   106.54753873   105.37737552 ...,   620.30765296
    622.00621723   623.48891285]
 ...,
 [  270.00991989   267.0670624    264.05274091 ...,  1038.2843583
   1040.94024978  1043.45777884]
 [ 1299.92652164  1295.17208978  1290.23553886 ...,  1065.34345608
   1070.26459357  1075.20109653]
 [ -558.78922472  -557.96853555  -557.33625563 ...,   170.65842728
    170.37271597   169.85343857]]
result:
 [[ -770.25586401  -768.56713682  -766.89550198 ...,   137.36536147
    136.34200923   135.04284885]
 [  -21.31085632   -17.95551232   -14.68439496 ..., -1604.80220015
  -1607.69691377 -1610.92271458]
 [  107.86839858   106.50291503   105.39849126 ...,   620.31999264
    622.06018438   623.49462413]
 ...,
 [  270.09021126   267.08492839   264.0995893  ...,  1038.3405256
   1040.97484038  1043.4317381 ]
 [ 1299.92993618  1295.21761955  1290.22616136 ...,  1065.36910928
   1070.26075221  1075.25002885]
 [ -558.80031728  -557.92143323  -557.32236354 ...,   170.62211818
    170.4311075    169.91769895]]

Exercise 9.2: Solving a linear system
Generate a vector b with m entries and solve Bx = b.

vector_b = np.random.random((500, 1))
# print('vector_b:\n', vector_b)
sloution = np.linalg.solve(B, vector_b)
print('sloution:\n', sloution)

输出：

sloution:
 [[-0.19933145]
 [ 0.46112091]
 [-0.28651194]
 [-0.36007313]
 [ 0.45699715]
 [-0.1136803 ]
 [ 0.06858123]
 [-0.02439619]
 [ 0.01277973]
 [ 0.07296403]
 [ 0.10677244]
 [-0.53755885]
 [ 0.70831459]
 [-0.50730103]
 [-0.12370938]
 ...

Exercise 9.3: Norms

Compute the Frobenius norm of A: ||A||F and the infinity norm of B: ||B||∞. Also find the largest and smallest singular values of B.

Frobenius_norm_A = np.linalg.norm(A, ord='fro')
print('Frobenius_norm_A:\n', Frobenius_norm_A)
infinity_norm_B = np.linalg.norm(B, np.inf)
print("infinity_norm_B:\n", infinity_norm_B)
U, sigma_singular, VT = np.linalg.svd(B)
max_singular_values = max(sigma_singular)
min_singular_values = min(sigma_singular)
print('max_singular_values:', max_singular_values)
print('min_singular_values:', min_singular_values)

输出：

Frobenius_norm_A:
 31.5386895558
infinity_norm_B:
 125250.0
max_singular_values: 87334.5204564
min_singular_values: 0.500004934835

Exercise 9.4: Power iteration
Generate a matrix Z, n × n, with Gaussian entries, and use the power iteration to find the largest eigenvalue and corresponding eigenvector of Z. How many iterations are needed till convergence?

Optional: use the time.clock() method to compare computation time when varying n.

根据《数值分析与算法》一书中所提供的幂法，得到下面的做法：

import numpy as np
from scipy.linalg import toeplitz
import time


def power_iteration(A, n_size):
    u = np.random.randint(0, 10, size=n_size)
    num_iter = 0
    start = time.clock()
    next_lam = 1
    pre_lam = 0
    while abs(next_lam - pre_lam) > 0.00001:
        pre_lam = next_lam
        num_iter += 1
        v = np.dot(A, u)
        next_lam = max(abs(v))
        u = v/next_lam
    end = time.clock()
    return u, next_lam, num_iter, end - start


n_size = 200
m_size = 500
mu = 0.0
sigma = 0.1
np.random.seed(0)
Z = np.random.normal(mu, sigma, (n_size, n_size))
eigenvector, lam, num_iters, time = power_iteration(Z, n_size)
print('the largest eigenvalue: ', lam)
print('the corresponding eigenvector: ', eigenvector)
print('the number of iterations: ', num_iters)
print('the time of iterations: ', time)

这里取n = 200输出：

the largest eigenvalue:  1.38374435975
the corresponding eigenvector:  [-0.23448063 -0.23857682 -0.45558981 -0.22375566  0.49945187 -0.04983987
  0.22586259  0.27561021 -0.25905579  0.27833055  0.10392899  0.10268637
  0.27160415  0.21686688 -0.41957299  0.17116327 -0.17859623  0.16854047
 -0.74401348 -0.41213827  0.42046522  0.34716834  0.05462962  0.68213181
 -0.26967417  0.22156731 -0.05764227  0.15715657 -0.34579995  0.22889528
 -0.43437992 -0.01709412  0.27750828 -0.32423244 -0.29912263 -0.33907246
 -0.13649076 -0.33804399 -0.12328962  0.07250919 -0.0715865   0.2687078
  0.12319253  0.07354623 -0.09890419  0.1051427  -0.37148188 -0.67210815
  0.21043618  0.04296881  0.35635652  0.28666593  0.05953488 -0.52401179
  0.10210595  0.58416896  0.49189862  0.01899168  0.05358757  0.28626702
 -0.31010217 -0.75048177  0.12335981 -0.80049872  0.65508207  0.00988996
  0.92204464 -0.41418936 -0.16885457 -0.16748709 -0.11395704  0.16557866
 -0.22648951 -0.09870179  0.57416871 -0.07121615 -0.62246095 -0.20429331
 -0.49580195  0.30371582 -0.75810448 -0.00355056  0.41848717  0.07971469
  0.28304665  0.48246567 -0.28948918  0.57802151  0.17559049 -0.14606612
 -0.00850067  0.43262172  0.04759938  0.52080095  0.06836337 -0.14321346
  0.20378579 -0.05444155 -0.49864496 -0.11220402  0.03517885  0.53004052
  0.15553566  0.25167114  0.00179509 -0.62592081 -0.31307944 -0.12592512
  0.05764811  0.10770587 -0.26141666  0.32830107 -0.23000649  0.32579796
 -0.42460092 -0.54333748 -0.07546899 -0.1202275  -0.22148577  0.46604526
 -0.37403106 -0.19443784 -0.02189891  0.3196339   0.45715056 -0.95319163
 -0.11534788  0.3458295  -0.07180602 -0.18107832 -0.4378873  -0.70818384
  0.06307144 -0.75049885  0.15457617 -0.65519509 -0.24892478  0.04379539
 -0.05674392 -0.45369933 -0.60605137  0.14269135 -0.15789244  0.71103203
  0.10365828 -0.50814752 -0.3803296  -0.3921604   0.34154827  0.03067514
  0.36761996 -0.01638594 -0.47892651 -0.62457798  0.62175791  0.60974948
 -0.34869693 -0.15907353  0.47225017  0.357265   -0.01573852 -0.36644187
  1.          0.11222906  0.14525716 -0.24414313  0.84366931 -0.77558439
  0.61776754  0.08202249 -0.27274454  0.05338878 -0.58895049  0.43344831
  0.32846749 -0.06703431  0.06770655  0.06771246  0.28469391 -0.04991912
 -0.03054763 -0.07803063 -0.23599849  0.26720273  0.05086394 -0.05258738
 -0.48809871  0.40606355 -0.34869432  0.07958705  0.04265011 -0.32885641
 -0.35524689  0.68314129 -0.14964527 -0.222974    0.09612197 -0.42326915
  0.11732494  0.15071393]
the number of iterations:  1813
the time of iterations:  0.13679134993804806

n = 100时：

the largest eigenvalue:  1.00512510268
the corresponding eigenvector:  [ 0.25042316  0.45906092 -0.43155535  0.17084169 -0.31571496 -0.32947339
  0.09475778 -0.07051941 -0.47606022 -0.44952443  1.         -0.2734709
  0.83039778 -0.07502614 -0.30195033  0.43850008 -0.28705981 -0.34490665
 -0.8302928  -0.14746054 -0.20792767  0.46113797  0.11940885 -0.43873384
 -0.4026181   0.18259404 -0.29167808 -0.57137324 -0.26113824  0.0447006
  0.6175856  -0.20932656  0.1497286   0.40798432 -0.48649966  0.05361968
 -0.15613001 -0.63377699 -0.19008996 -0.50509828  0.20765657 -0.67042046
 -0.04948563  0.1925169   0.30199847 -0.05921492  0.08473466  0.40901562
  0.56484163  0.60907853  0.52669863  0.50567972 -0.04177342  0.231218
 -0.87164059 -0.18402144  0.25972169 -0.06419379  0.64224447 -0.09154909
 -0.26810684  0.01138896 -0.40211867 -0.11467527  0.52627956  0.08511445
  0.06288489 -0.06932288  0.19850499 -0.0042502  -0.44447338  0.00977573
 -0.26873374  0.02494805  0.4358977  -0.4499939   0.42533057 -0.00221055
 -0.3709814  -0.44295519  0.0393151   0.24926072  0.26211296  0.29795194
  0.06267313  0.02638521 -0.34356405  0.40541127 -0.04932551 -0.53103175
 -0.07837982  0.02486878 -0.0533286  -0.31930375 -0.20309002 -0.00346277
  0.37197201  0.34570481  0.2211697   0.48039722]
the number of iterations:  54074
the time of iterations:  1.5376810345121839

Exercise 9.5: Singular values
Generate an n × n matrix, denoted by C, where each entry is 1 with probability p and 0 otherwise. Use the linear algebra library of Scipy to compute the singular values of C. What can you say about the relationship between n, p and the largest singular value?

import numpy as np

def one_or_zero(p):
    if np.random.random() > 1 - p:
        return 1
    else:
        return 0

def largest_singular_value(A):
    U, S, V = np.linalg.svd(A)
    return max(S)

for n in range(2, 200, 10):
    C = np.zeros((n, n))
    for k in range(1, 10):
        p = k / 10
        for i in range(n):
            for j in range(n):
                C[i][j] = one_or_zero(0.5)
        value = largest_singular_value(C)
        print('the largest singular value: ', value, 'p: ', p, 'n: ', n)

输出：

the largest singular value:  1.41421356237 p:  0.1 n:  2
the largest singular value:  1.0 p:  0.2 n:  2
the largest singular value:  1.41421356237 p:  0.3 n:  2
the largest singular value:  1.61803398875 p:  0.4 n:  2
the largest singular value:  1.0 p:  0.5 n:  2
the largest singular value:  1.0 p:  0.6 n:  2
the largest singular value:  1.0 p:  0.7 n:  2
the largest singular value:  1.0 p:  0.8 n:  2
the largest singular value:  1.61803398875 p:  0.9 n:  2
the largest singular value:  6.09803062927 p:  0.1 n:  12
the largest singular value:  5.36168216768 p:  0.2 n:  12
the largest singular value:  6.63331122699 p:  0.3 n:  12
the largest singular value:  6.13206436997 p:  0.4 n:  12
the largest singular value:  6.71450698537 p:  0.5 n:  12
the largest singular value:  6.43079837273 p:  0.6 n:  12
the largest singular value:  6.57675764838 p:  0.7 n:  12
the largest singular value:  6.72227602896 p:  0.8 n:  12
the largest singular value:  6.43344014997 p:  0.9 n:  12
the largest singular value:  11.5611945294 p:  0.1 n:  22
the largest singular value:  10.861829775 p:  0.2 n:  22
the largest singular value:  11.2604419494 p:  0.3 n:  22
the largest singular value:  10.7399483624 p:  0.4 n:  22
the largest singular value:  11.427691543 p:  0.5 n:  22
the largest singular value:  11.3774277544 p:  0.6 n:  22
the largest singular value:  11.1564436536 p:  0.7 n:  22
the largest singular value:  12.663423727 p:  0.8 n:  22
the largest singular value:  10.861248278 p:  0.9 n:  22
...
the largest singular value:  87.0336062627 p:  0.1 n:  172
the largest singular value:  86.8825137409 p:  0.2 n:  172
the largest singular value:  86.1686471393 p:  0.3 n:  172
the largest singular value:  86.7796838102 p:  0.4 n:  172
the largest singular value:  87.3335809673 p:  0.5 n:  172
the largest singular value:  86.6432805456 p:  0.6 n:  172
the largest singular value:  86.3703671617 p:  0.7 n:  172
the largest singular value:  87.0565486864 p:  0.8 n:  172
the largest singular value:  86.0200230692 p:  0.9 n:  172
the largest singular value:  90.7639134334 p:  0.1 n:  182
the largest singular value:  91.5999226806 p:  0.2 n:  182
the largest singular value:  91.913131263 p:  0.3 n:  182
the largest singular value:  92.0368133674 p:  0.4 n:  182
the largest singular value:  91.29022559 p:  0.5 n:  182
the largest singular value:  91.089694189 p:  0.6 n:  182
the largest singular value:  91.2739383084 p:  0.7 n:  182
the largest singular value:  91.5252808991 p:  0.8 n:  182
the largest singular value:  91.9294664348 p:  0.9 n:  182
the largest singular value:  95.9877956052 p:  0.1 n:  192
the largest singular value:  96.3236948637 p:  0.2 n:  192
the largest singular value:  96.1700156123 p:  0.3 n:  192
the largest singular value:  97.0419055764 p:  0.4 n:  192
the largest singular value:  96.3372471599 p:  0.5 n:  192
the largest singular value:  96.6318842389 p:  0.6 n:  192
the largest singular value:  96.2941907085 p:  0.7 n:  192
the largest singular value:  96.8051842424 p:  0.8 n:  192
the largest singular value:  96.2303337938 p:  0.9 n:  192

可以看到the largest singular value在n比较小时其值受到p的影响还是比较大，但在当n比较大时，其值并没有太大的变化，而且n比较时的数值并没有太大的代表性，所以可以只看n比较大时the largest singular value的变化，而对于p相同时，n增加10，the largest singular value的变化还是很明显的，并且随着n的增大而增大。

Exercise 9.6: Nearest neighbor
Write a function that takes a value z and an array A and finds the element in A that is closest to z. The function should return the closest value, not index.
Hint: Use the built-in functionality of Numpy rather than writing code to find this value manually. In particular, use brackets and argmin.

import numpy as np

def get_nearest_neighbor(z, A):
    B = abs(A - z)
    return A[np.argmin(B)]

z = np.random.rand() * 10
A = np.random.rand(100) * 10
print('z: ', z)
print('A: ', A)
print('the nearest neighbor: ', get_nearest_neighbor(z, A))

输出：

z:  5.140756589959544
A:  [ 2.38548496  9.49087376  1.28884664  0.63702924  6.67374531  9.51206904
  5.84934177  8.86422373  4.18638771  4.99405837  3.96835962  7.17309549
  9.3484444   1.04547294  1.88207947  4.10086577  5.2077526   3.7403351
  0.49729221  2.69092716  4.49822518  2.48918472  9.79170684  1.99704282
  6.90493967  1.35927367  0.58062983  3.58336034  3.34607267  1.63150145
  7.13385139  5.0682374   2.87122973  7.39368304  6.20506582  7.67409263
  2.50880758  2.15953104  3.44241519  5.26641213  7.80353381  9.44048273
  8.34729033  1.8021764   6.61065728  6.02415552  2.81593119  2.72401768
  5.28152751  3.18988028  3.57815635  8.04498665  6.87776589  6.5423706
  4.53766417  4.68319582  7.91227892  2.09053791  6.13504116  5.27837593
  3.90897699  6.74414173  4.44737472  4.3569735   9.65642109  6.08779349
  5.45481335  3.04881375  2.23628713  0.22330483  4.68740936  1.07713896
  0.3546105   1.53122082  7.82766266  8.33500204  8.45982052  5.94027026
  5.80184327  8.2297121   3.49990858  8.6759586   5.70309204  7.70432551
  8.30971419  9.86307761  6.39693714  1.44703552  9.81276596  4.73014594
  0.61624376  3.24323494  7.41945733  8.60669928  9.49122329  6.63813434
  5.43821111  3.24123646  4.52327842  0.33085725]
the nearest neighbor:  5.20775260306