資訊數學(下)講義範例 R 程式碼與執行結果

請務必注意 R 程式碼英文字大小寫有區別!

Page-01

計算小羊從頭到尾移動的效果並印出

u <- c(0,4)
v <- c(3,0)
w <- c(0,-2)
a <- c(-2,-4)
b <- c(3,0)
c <- c(-6,4)

(tot <- u+v+w+a+b+c)

## [1] -2  2

計算小羊如何由 \(G\) 點回到原點並印出

O <- c(0,0)
G <- c(-2,2)
(O-G)

## [1]  2 -2

計算小羊若由 \(D\) 點開始出發，最後將位居何處，並印出

D <- c(3,2)
(D+tot)

## [1] 1 4

Page-02

向量 \(u\) 與 \(v\) 的加總

u <- c(6,8)
v <- c(-4,14)
(u+v)

## [1]  2 22

(u+1/2*v)

## [1]  4 15

向量 \(u\) 的長度 (內積開根號)

sqrt(u%*%u)

##      [,1]
## [1,]   10

向量 \(v\) 的長度 (內積開根號)

sqrt(v%*%v)

##          [,1]
## [1,] 14.56022

\(ru+sv=c(0,8r+14s) \Rightarrow 6r=4s \Rightarrow r:s = 2:3 = 1:1.5\)

調整為大約 \(10:21.84033 \ \ (10: \sqrt{477}\) )

(1*u+1.5*v)

## [1]  0 29

(w <- 1.5*v)

## [1] -6 21

sqrt(w%*%w)

##          [,1]
## [1,] 21.84033

五位主管的評分加總 (甲錄取)

A <- c(88,95,92)
B <- c(96,95,91)
C <- c(91,84,88)
D <- c(85,90,89)
E <- c(95,88,91)

(A+B+C+D+E)

## [1] 455 452 451

加權後之加總總分 (乙錄取)

(3*A+2*B+C+D+E)

## [1] 727 737 726

Page-03

計算 \(A \rightarrow B \rightarrow C \rightarrow D\) 點的飛行總距離

A <- c(0,0,0)
B <- c(2,6,5)
C <- c(4,5,3)
D <- c(5,5,5)

u <- B-A  # 或直接計算 u <- c(2,6,5)
v <- C-B  # 或直接計算 v <- c(2,-1,-2)
w <- D-C  # 或直接計算 w <- c(1,0,2)

(tot <- sqrt(u%*%u)+ sqrt(v%*%v)+ sqrt(w%*%w))

##          [,1]
## [1,] 13.29833

計算直接飛到 \(D\) 點的飛行距離

s <- u+v+w
(sa<-sqrt(s%*%s))

##          [,1]
## [1,] 8.660254

計算節省的距離

tot-sa

##          [,1]
## [1,] 4.638072

老王、老張的幸福向量 \(u_w、u_c\)

uw <- c(0.88,0.75,0.94,0.84,0.93)
uc <- c(0.92,0.80,0.79,0.82,0.95)

計算老王、老張的幸福向量長度

(sqrt(uw%*%uw))

##          [,1]
## [1,] 1.947049

(sqrt(uc%*%uc))

##         [,1]
## [1,] 1.91974

計算老王、老張與幸福最高境界的差距

up <- c(1,1,1,1,1) # 幸福最高境界
dw <- up-uw 
dc <- up-uc

計算老王、老張與幸福最高境界的差距之長度

(sqrt(dw%*%dw))

##           [,1]
## [1,] 0.3331666

(sqrt(dc%*%dc))

##           [,1]
## [1,] 0.3541186

Page-04~05

建立 \(A、B、C\) 矩陣，計算矩陣相加結果 \(M\) 並印出

A <- matrix(c(20,60,40,12,30,20), nrow=3) # 亦可使用 A <- matrix(c(20,12,60,30,40,20), nrow=3, byrow=T)
B <- matrix(c(25,55,30,18,40,30), nrow=3)
C <- matrix(c(25,50,45,20,35,30), nrow=3)

(M <- A+B+C)

##      [,1] [,2]
## [1,]   70   50
## [2,]  165  105
## [3,]  115   80

建立 \(D\) 矩陣

D <- matrix(c(2200,1800,1000,850,1500,900),nrow=2)
D

##      [,1] [,2] [,3]
## [1,] 2200 1000 1500
## [2,] 1800  850  900

將矩陣 \(D\) 的第一列與矩陣 \(M\) 的第一行做內積 (視為兩向量)

D[1,]%*%M[,1]

##        [,1]
## [1,] 491500

將矩陣 \(D\) 的第二列與矩陣 \(M\) 的第二行做內積 (視為兩向量)

D[2,]%*%M[,2]

##        [,1]
## [1,] 251250

以下是將矩陣中的向量取出，再轉換為矩陣的做法

D1 <- as.matrix(D[1,])    # 將矩陣D的第一列視為矩陣D1

M1 <- as.matrix(M[,1])    # 將矩陣M的第一行視為矩陣M1

dim(D1)       # 注意D1為3x1矩陣

## [1] 3 1

class(D1)     # 確認D1之類型

## [1] "matrix" "array"

class(M1)     # 確認M1之類型

## [1] "matrix" "array"

t(D1)%*%M1    # 將D1轉置為1x3矩陣後再與矩陣M1(3x1)相乘

##        [,1]
## [1,] 491500

Page-06

A1 <- matrix(c(1,3,2,4),nrow=2)
A2 <- matrix(c(3,1,4,2),nrow=2)
(A1%*%A2)

##      [,1] [,2]
## [1,]    5    8
## [2,]   13   20

B1 <- matrix(c(2,5,-1,3,4,0,-1,0,2),nrow=3)
B2 <- matrix(c(0,2,-1,1,0,1),nrow=3)
(B1%*%B2)

##      [,1] [,2]
## [1,]    7    1
## [2,]    8    5
## [3,]   -2    1

C1 <- matrix(c(2,1,0,-1),nrow=2)
C2 <- matrix(c(3,-2,0,2),nrow=2)
C3 <- matrix(c(1,0,2,3),nrow=2)
(C1%*%C2%*%C3)

##      [,1] [,2]
## [1,]    6   12
## [2,]    5    4

D1 <- matrix(c(0,2,1,2),nrow=2)
D2 <- matrix(c(0,1,3,2,-1,0),nrow=2)
D3 <- matrix(c(0,2,1,1,-1,0),nrow=3)
(D1%*%D2%*%D3)

##      [,1] [,2]
## [1,]    4   -1
## [2,]   18   -8

MA <- matrix(c(2,1,0,-1),nrow=2)
MB <- matrix(c(3,2,0,-2),nrow=2)
MC <- matrix(c(1,0,2,3),nrow=2)
(MA+2*MB-3*MC)

##      [,1] [,2]
## [1,]    5   -6
## [2,]    5  -14

Page-07

建立 \(K1、K2、K3\) 三個轉換矩陣

K1 <- matrix(c(-1,0,0,1),nrow=2)
K2 <- matrix(c(0,1,1,0),nrow=2)
K3 <- matrix(c(0.5,0,0,0.5),nrow=2)

建立 \(A、B、C\) 三點座標構成的矩陣 \(MM\) 並印出

MM <- matrix(c(1,3,2,5,2.5,4.5),nrow=2)
MM

##      [,1] [,2] [,3]
## [1,]    1    2  2.5
## [2,]    3    5  4.5

顯示 \(K1\) 矩陣與 \(MM\) 座標矩陣透過 \(K1\) 矩陣轉換後的結果

K1

##      [,1] [,2]
## [1,]   -1    0
## [2,]    0    1

K1%*%MM

##      [,1] [,2] [,3]
## [1,]   -1   -2 -2.5
## [2,]    3    5  4.5

顯示 \(K2\) 矩陣與 \(MM\) 座標矩陣透過 \(K2\) 矩陣轉換後的結果

K2

##      [,1] [,2]
## [1,]    0    1
## [2,]    1    0

K2%*%MM

##      [,1] [,2] [,3]
## [1,]    3    5  4.5
## [2,]    1    2  2.5

顯示 \(K3\) 矩陣與 \(MM\) 座標矩陣透過 \(K3\) 矩陣轉換後的結果

K3

##      [,1] [,2]
## [1,]  0.5  0.0
## [2,]  0.0  0.5

K3%*%MM

##      [,1] [,2] [,3]
## [1,]  0.5  1.0 1.25
## [2,]  1.5  2.5 2.25

合併以下兩步驟構成 \(K21\) 矩陣

以 \(Y\) 軸做對稱軸做轉換
以 \(Y=X\) 做對稱軸做轉換

並顯示 \(MM\) 座標矩陣透過 \(K21\) 矩陣轉換後的結果

(K21 <- K2%*%K1)

##      [,1] [,2]
## [1,]    0    1
## [2,]   -1    0

(K21%*%MM)

##      [,1] [,2] [,3]
## [1,]    3    5  4.5
## [2,]   -1   -2 -2.5

合併以下兩步驟構成 \(K12\) 矩陣

以 \(Y=X\) 做對稱軸做轉換
以 \(Y\) 軸做對稱軸做轉換

並顯示 \(MM\) 座標矩陣透過 \(K12\) 矩陣轉換後的結果

(K12 <- K1%*%K2)

##      [,1] [,2]
## [1,]    0   -1
## [2,]    1    0

(K12%*%MM)

##      [,1] [,2] [,3]
## [1,]   -3   -5 -4.5
## [2,]    1    2  2.5

以 \(Y\) 軸做對稱軸做兩次轉換

(K1%*%K1)

##      [,1] [,2]
## [1,]    1    0
## [2,]    0    1

以 \(Y=X\) 做對稱軸做兩次轉換

(K2%*%K2)

##      [,1] [,2]
## [1,]    1    0
## [2,]    0    1

Page-08

A <- matrix(c(8,4,4,7),nrow=2)
solve(A)

##        [,1] [,2]
## [1,]  0.175 -0.1
## [2,] -0.100  0.2

B <- matrix(c(0,3,2,1),nrow=2)
solve(B)

##            [,1]      [,2]
## [1,] -0.1666667 0.3333333
## [2,]  0.5000000 0.0000000

解 \(k^2+8k-8 = 0\)，注意多項式係數升冪順序

polyroot(c(-8,8,1))

## [1]  0.8989795+0i -8.8989795+0i

Page-09

計算 \(A\) (即前述 \(K12\)) 的逆矩陣 \(AI\) 並將其作用在點 \(b = (3,-1)\) 上

A <- matrix(c(0,-1,1,0),nrow=2)
AI <- solve(A) 
b <- matrix(c(3,-1),nrow=2)
AI%*%b

##      [,1]
## [1,]    1
## [2,]    3

求對稱於原點, 複製一長度為兩倍大的物件之轉換所需要的矩陣 \(D\) 及其逆矩陣 \(DI\)

(D <- matrix(c(-2,0,0,-2),nrow=2))

##      [,1] [,2]
## [1,]   -2    0
## [2,]    0   -2

(DI <- solve(D))

##      [,1] [,2]
## [1,] -0.5  0.0
## [2,]  0.0 -0.5

以下用三個點座標構成的矩陣 \(M\) 來驗證

(M <- matrix(c(2,8,4,7,4,5),nrow=2))

##      [,1] [,2] [,3]
## [1,]    2    4    4
## [2,]    8    7    5

三個點透過矩陣 \(D\) 轉換得新座標矩陣 \(MN\)

(MN <- D%*%M)

##      [,1] [,2] [,3]
## [1,]   -4   -8   -8
## [2,]  -16  -14  -10

新座標再透過 \(D\) 的逆矩陣 \(DI\) 轉換回原來的座標

(DI%*%MN)

##      [,1] [,2] [,3]
## [1,]    2    4    4
## [2,]    8    7    5

Page-11

利用行列式計算 \(x、y\) 並印出 (1) ，\(\mbox{\n}\) 是換行指令

A <- matrix(c(1,2,1,4),nrow=2)
AX <- matrix(c(10,28,1,4),nrow=2)
AY <- matrix(c(1,2,10,28),nrow=2)

cat("X=",det(AX)/det(A), "\n")

## X= 6

cat("Y=",det(AY)/det(A), "\n")

## Y= 4

利用行列式計算 \(x、y\) 並印出 (2)

B <- matrix(c(4, 6, 5, -3),nrow=2)
BX <- matrix(c(7, 21, 5, -3),nrow=2)
BY <- matrix(c(4,6, 7, 21),nrow=2)

cat("x=",det(BX)/det(B), "\n")

## x= 3

cat("y=",det(BY)/det(B), "\n")

## y= -1

利用行列式計算 \(x、y\) 並印出 (3)

C <- matrix(c(0.8,0.3,0.6,0.4),nrow=2)
CX <- matrix(c(4,3,0.6,0.4),nrow=2)
CY <- matrix(c(0.8,0.3,4,3),nrow=2)

cat("x=",det(CX)/det(C), "\n")

## x= -1.428571

cat("y=",det(CY)/det(C), "\n")

## y= 8.571429

利用行列式計算 \(x、y\) 並印出 (4)

D <- matrix(c(3/5,1/4,1/2,4/3),nrow=2)
DX <- matrix(c(2,4,1/2,4/3),nrow=2)
DY <- matrix(c(3/5,1/4,2,4),nrow=2)

cat("x=",det(DX)/det(D), "\n")

## x= 0.9876543

cat("y=",det(DY)/det(D), "\n")

## y= 2.814815

Page-12

將雞兔蝶同籠聯立方程式係數部分建立 \(D\) 矩陣，利用行列式計算 \(x、y、z\) 並印出

D <- matrix(c(1,2,2,1,4,0,1,6,4),nrow=3)
DX <- matrix(c(10,38,14,1,4,0,1,6,4),nrow=3)
DY <- matrix(c(1,2,2,10,38,14,1,6,4),nrow=3)
DZ <- matrix(c(1,2,2,1,4,0,10,38,14),nrow=3)

cat("x=",det(DX)/det(D), "\n")

## x= 3

cat("y=",det(DY)/det(D), "\n")

## y= 5

cat("z=",det(DZ)/det(D), "\n")

## z= 2

Page-13

將聯立方程式係數部分建立 \(M\) 矩陣，利用行列式計算 \(x、y、z\) 並印出

M <- matrix(c(3,4,2,8,1,0,5,3,2),nrow=3)
MX <- matrix(c(60,46,24,8,1,0,5,3,2),nrow=3)
MY <- matrix(c(3,4,2,60,46,24,5,3,2),nrow=3)
MZ <- matrix(c(3,4,2,8,1,0,60,46,24),nrow=3)

cat("x=",det(MX)/det(M), "\n")

## x= 8

cat("y=",det(MY)/det(M), "\n")

## y= 2

cat("z=",det(MZ)/det(M), "\n")

## z= 4

(另解) 計算 \(M\) 的反矩陣 \(MI\)，透過 \(MI\) 乘以係數矩陣 \(b\) 解出 \((x,y,z)\) 向量 (\(Mx=b \Longrightarrow x=M^{-1}b\))，並印出

b <- matrix(c(60,46,24),nrow=3)
MI <- solve(M)                      
cat("x y z=", MI%*%b, "\n")

## x y z= 8 2 4

Page-14~15 (僅供參考)

利用 cbind 指令將矩陣 \(A\) 與 \(b\) 合併(並列)

A <- matrix(c(1,2,2,1,4,0,1,6,4),nrow=3)
b <- matrix(c(10,38,14),nrow=3)
(Ab <-cbind(A,b))

##      [,1] [,2] [,3] [,4]
## [1,]    1    1    1   10
## [2,]    2    4    6   38
## [3,]    2    0    4   14

Ab[2,] <- Ab[1,]*(-2)+Ab[2,]    
# 將合併矩陣之第一列乘以-2加到第二列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    1    1   10
## [2,]    0    2    4   18
## [3,]    2    0    4   14

Ab[3,] <- Ab[1,]*(-2)+Ab[3,]    
# 將合併矩陣之第一列乘以-2加到第三列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    1    1   10
## [2,]    0    2    4   18
## [3,]    0   -2    2   -6

Ab[2,] <- Ab[2,]*(1/2)          
# 將合併矩陣之第二列乘以1/2
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    1    1   10
## [2,]    0    1    2    9
## [3,]    0   -2    2   -6

Ab[1,] <- Ab[2,]*(-1)+Ab[1,]    
# 將合併矩陣之第二列乘以-1加到第一列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    0   -1    1
## [2,]    0    1    2    9
## [3,]    0   -2    2   -6

Ab[3,] <- Ab[2,]*2+Ab[3,]       
# 將合併矩陣之第二列乘以2加到第三列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    0   -1    1
## [2,]    0    1    2    9
## [3,]    0    0    6   12

Ab[3,] <- Ab[3,]*(1/6)          
# 將合併矩陣之第三列乘以1/6
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    0   -1    1
## [2,]    0    1    2    9
## [3,]    0    0    1    2

Ab[1,] <- Ab[3,]*1+Ab[1,]       
# 將合併矩陣之第三列乘以1加到第一列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    0    0    3
## [2,]    0    1    2    9
## [3,]    0    0    1    2

Ab[2,] <- Ab[3,]*(-2)+Ab[2,]    
# 將合併矩陣之第三列乘以-2 加到第二列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    0    0    3
## [2,]    0    1    0    5
## [3,]    0    0    1    2

利用 cbind 指令將矩陣 \(A\) 與 \(I_3\) 合併(並列)

A <- matrix(c(1,2,2,1,4,0,1,6,4),nrow=3)
(AI3 <-cbind(A,diag(3)))

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    0    0
## [2,]    2    4    6    0    1    0
## [3,]    2    0    4    0    0    1

AI3[2,] <- AI3[1,]*(-2)+AI3[2,]    
# 將合併矩陣之第一列乘以-2加到第二列
AI3

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    0    0
## [2,]    0    2    4   -2    1    0
## [3,]    2    0    4    0    0    1

AI3[3,] <- AI3[1,]*(-2)+AI3[3,]    
# 將合併矩陣之第一列乘以-2加到第三列
AI3

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    0    0
## [2,]    0    2    4   -2    1    0
## [3,]    0   -2    2   -2    0    1

AI3[2,] <- AI3[2,]*(1/2)          
# 將合併矩陣之第二列乘以1/2
AI3

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1  0.0    0
## [2,]    0    1    2   -1  0.5    0
## [3,]    0   -2    2   -2  0.0    1

AI3[1,] <- AI3[2,]*(-1)+AI3[1,]    
# 將合併矩陣之第二列乘以-1加到第一列
AI3

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0   -1    2 -0.5    0
## [2,]    0    1    2   -1  0.5    0
## [3,]    0   -2    2   -2  0.0    1

AI3[3,] <- AI3[2,]*2+AI3[3,]       
# 將合併矩陣之第二列乘以2加到第三列
AI3

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0   -1    2 -0.5    0
## [2,]    0    1    2   -1  0.5    0
## [3,]    0    0    6   -4  1.0    1

AI3[3,] <- AI3[3,]*(1/6)          
# 將合併矩陣之第三列乘以1/6
AI3

##      [,1] [,2] [,3]       [,4]       [,5]      [,6]
## [1,]    1    0   -1  2.0000000 -0.5000000 0.0000000
## [2,]    0    1    2 -1.0000000  0.5000000 0.0000000
## [3,]    0    0    1 -0.6666667  0.1666667 0.1666667

AI3[1,] <- AI3[3,]*1+AI3[1,]       
# 將合併矩陣之第三列乘以1加到第一列
AI3

##      [,1] [,2] [,3]       [,4]       [,5]      [,6]
## [1,]    1    0    0  1.3333333 -0.3333333 0.1666667
## [2,]    0    1    2 -1.0000000  0.5000000 0.0000000
## [3,]    0    0    1 -0.6666667  0.1666667 0.1666667

AI3[2,] <- AI3[3,]*(-2)+AI3[2,]    
# 將合併矩陣之第三列乘以-2 加到第二列
AI3

##      [,1] [,2] [,3]       [,4]       [,5]       [,6]
## [1,]    1    0    0  1.3333333 -0.3333333  0.1666667
## [2,]    0    1    0  0.3333333  0.1666667 -0.3333333
## [3,]    0    0    1 -0.6666667  0.1666667  0.1666667

亦可直接利用 \(A\) 的逆矩陣 \(A^{-1}\) 求解

(AI <- AI3[,4:6])            # 合併矩陣之後三欄即為欲求之反矩陣

##            [,1]       [,2]       [,3]
## [1,]  1.3333333 -0.3333333  0.1666667
## [2,]  0.3333333  0.1666667 -0.3333333
## [3,] -0.6666667  0.1666667  0.1666667

cat("x y z=", AI%*%b, "\n")    # 印出(x,y,z)之解

## x y z= 3 5 2

Page-16 (僅供參考)

利用 cbind 指令將矩陣 \(A\) 與 \(b\) 合併(並列)

A <- matrix(c(1,2,4,5,4,-2,-2,3,1),nrow=3)
b <- matrix(c(21,17,9),nrow=3)
(Ab <-cbind(A,b))

##      [,1] [,2] [,3] [,4]
## [1,]    1    5   -2   21
## [2,]    2    4    3   17
## [3,]    4   -2    1    9

Ab[2,] <- Ab[1,]*(-2)+Ab[2,]    
# 將合併矩陣之第一列乘以-2加到第二列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    5   -2   21
## [2,]    0   -6    7  -25
## [3,]    4   -2    1    9

Ab[3,] <- Ab[1,]*(-4)+Ab[3,]    
# 將合併矩陣之第一列乘以-4加到第三列
Ab

##      [,1] [,2] [,3] [,4]
## [1,]    1    5   -2   21
## [2,]    0   -6    7  -25
## [3,]    0  -22    9  -75

Ab[2,] <- Ab[2,]*(-1/6)          
# 將合併矩陣之第二列乘以 -1/6
Ab

##      [,1] [,2]      [,3]       [,4]
## [1,]    1    5 -2.000000  21.000000
## [2,]    0    1 -1.166667   4.166667
## [3,]    0  -22  9.000000 -75.000000

Ab[1,] <- Ab[2,]*(-5)+Ab[1,]     
# 將合併矩陣之第二列乘以-5加到第一列
Ab

##      [,1] [,2]      [,3]        [,4]
## [1,]    1    0  3.833333   0.1666667
## [2,]    0    1 -1.166667   4.1666667
## [3,]    0  -22  9.000000 -75.0000000

Ab[3,] <- Ab[2,]*(22)+Ab[3,]    
# 將合併矩陣之第二列乘以22加到第三列
Ab

##      [,1] [,2]       [,3]       [,4]
## [1,]    1    0   3.833333  0.1666667
## [2,]    0    1  -1.166667  4.1666667
## [3,]    0    0 -16.666667 16.6666667

Ab[3,] <- Ab[3,]*(-1)/(16+2/3)    
# 將合併矩陣之第三列乘以-1/(16+2/3)
Ab

##      [,1] [,2]      [,3]       [,4]
## [1,]    1    0  3.833333  0.1666667
## [2,]    0    1 -1.166667  4.1666667
## [3,]    0    0  1.000000 -1.0000000

Ab[1,] <- Ab[3,]*(-(3+5/6))+Ab[1,]  
# 將合併矩陣之第三列乘以(-(3+5/6))加到第一列
Ab

##      [,1] [,2]          [,3]      [,4]
## [1,]    1    0 -4.440892e-16  4.000000
## [2,]    0    1 -1.166667e+00  4.166667
## [3,]    0    0  1.000000e+00 -1.000000

Ab[2,] <- Ab[3,]*(1+1/6)+Ab[2,]   
# 將合併矩陣之第三列乘以(1+1/6)加到第二列
Ab

##      [,1] [,2]          [,3] [,4]
## [1,]    1    0 -4.440892e-16    4
## [2,]    0    1  0.000000e+00    3
## [3,]    0    0  1.000000e+00   -1

利用 cbind 指令將矩陣 \(A\) 與 \(b1、b2\) 合併(並列)

A <- matrix(c(1,2,1,-1,3,4,2,1,-1),nrow=3)
b1 <- matrix(c(4,6,2),nrow=3)
b2 <- matrix(c(4,6,5),nrow=3)
(Ab12 <-cbind(A,b1,b2))

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1   -1    2    4    4
## [2,]    2    3    1    6    6
## [3,]    1    4   -1    2    5

Ab12[2,] <- Ab12[1,]*(-2)+Ab12[2,]  
# 將合併矩陣之第一列乘以-2加到第二列
Ab12

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1   -1    2    4    4
## [2,]    0    5   -3   -2   -2
## [3,]    1    4   -1    2    5

Ab12[3,] <- Ab12[1,]*(-1)+Ab12[3,]  
# 將合併矩陣之第一列乘以-1加到第三列
Ab12

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1   -1    2    4    4
## [2,]    0    5   -3   -2   -2
## [3,]    0    5   -3   -2    1

Ab12[2,] <- Ab12[2,]*(1/5)          
# 將合併矩陣之第二列乘以1/5
Ab12

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1   -1  2.0  4.0  4.0
## [2,]    0    1 -0.6 -0.4 -0.4
## [3,]    0    5 -3.0 -2.0  1.0

Ab12[1,] <- Ab12[2,]*(1)+Ab12[1,]   
# 將合併矩陣之第二列乘以1加到第一列
Ab12

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1    0  1.4  3.6  3.6
## [2,]    0    1 -0.6 -0.4 -0.4
## [3,]    0    5 -3.0 -2.0  1.0

Ab12[3,] <- Ab12[2,]*(-5)+Ab12[3,]    
# 將合併矩陣之第二列乘以-5加到第三列
Ab12

##      [,1] [,2]          [,3] [,4] [,5]
## [1,]    1    0  1.400000e+00  3.6  3.6
## [2,]    0    1 -6.000000e-01 -0.4 -0.4
## [3,]    0    0  4.440892e-16  0.0  3.0

# 對於b1，得 x1+1.4*x3=3.6 , x2-0.6*x2=-0.4 , 0*x1+0*x2+0*x3=0   => 無窮多解
# 對於b2，得 x1+1.4*x3=3.6 , x2-0.6*x2=-0.4 , 0*x1+0*x2+0*x3=3.0 => 無解

利用 cbind 指令將矩陣 \(A\) 與 \(I4\) 合併(並列)

A <- matrix(c(1,2,4,0,1,0,-1,1,0,1,3,2,-1,3,1,2),nrow=4)
(AI <-cbind(A,diag(4)))

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    1    0   -1    1    0    0    0
## [2,]    2    0    1    3    0    1    0    0
## [3,]    4   -1    3    1    0    0    1    0
## [4,]    0    1    2    2    0    0    0    1

AI[2,] <- AI[1,]*(-2)+AI[2,]    
# 將合併矩陣之第一列乘以-2加到第二列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    1    0   -1    1    0    0    0
## [2,]    0   -2    1    5   -2    1    0    0
## [3,]    4   -1    3    1    0    0    1    0
## [4,]    0    1    2    2    0    0    0    1

AI[3,] <- AI[1,]*(-4)+AI[3,]    
# 將合併矩陣之第一列乘以-4加到第三列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    1    0   -1    1    0    0    0
## [2,]    0   -2    1    5   -2    1    0    0
## [3,]    0   -5    3    5   -4    0    1    0
## [4,]    0    1    2    2    0    0    0    1

AI[2,] <- AI[2,]*(-1/2)         
# 將合併矩陣之第二列乘以-1/2
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    1  0.0 -1.0    1  0.0    0    0
## [2,]    0    1 -0.5 -2.5    1 -0.5    0    0
## [3,]    0   -5  3.0  5.0   -4  0.0    1    0
## [4,]    0    1  2.0  2.0    0  0.0    0    1

AI[1,] <- AI[2,]*(-1)+AI[1,]    
# 將合併矩陣之第二列乘以-1加到第一列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.5  1.5    0  0.5    0    0
## [2,]    0    1 -0.5 -2.5    1 -0.5    0    0
## [3,]    0   -5  3.0  5.0   -4  0.0    1    0
## [4,]    0    1  2.0  2.0    0  0.0    0    1

AI[3,] <- AI[2,]*(5)+AI[3,]     
# 將合併矩陣之第二列乘以5加到第三列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.5  1.5    0  0.5    0    0
## [2,]    0    1 -0.5 -2.5    1 -0.5    0    0
## [3,]    0    0  0.5 -7.5    1 -2.5    1    0
## [4,]    0    1  2.0  2.0    0  0.0    0    1

AI[4,] <- AI[2,]*(-1)+AI[4,]    
# 將合併矩陣之第二列乘以-1加到第四列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.5  1.5    0  0.5    0    0
## [2,]    0    1 -0.5 -2.5    1 -0.5    0    0
## [3,]    0    0  0.5 -7.5    1 -2.5    1    0
## [4,]    0    0  2.5  4.5   -1  0.5    0    1

AI[3,] <- AI[3,]*(2)            
# 將合併矩陣之第三列乘以2
AI

##      [,1] [,2] [,3]  [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.5   1.5    0  0.5    0    0
## [2,]    0    1 -0.5  -2.5    1 -0.5    0    0
## [3,]    0    0  1.0 -15.0    2 -5.0    2    0
## [4,]    0    0  2.5   4.5   -1  0.5    0    1

AI[1,] <- AI[3,]*(-0.5)+AI[1,]  
# 將合併矩陣之第三列乘以-0.5加到第一列
AI

##      [,1] [,2] [,3]  [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.0   9.0   -1  3.0   -1    0
## [2,]    0    1 -0.5  -2.5    1 -0.5    0    0
## [3,]    0    0  1.0 -15.0    2 -5.0    2    0
## [4,]    0    0  2.5   4.5   -1  0.5    0    1

AI[2,] <- AI[3,]*(0.5)+AI[2,]   
# 將合併矩陣之第三列乘以0.5加到第二列
AI

##      [,1] [,2] [,3]  [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0  0.0   9.0   -1  3.0   -1    0
## [2,]    0    1  0.0 -10.0    2 -3.0    1    0
## [3,]    0    0  1.0 -15.0    2 -5.0    2    0
## [4,]    0    0  2.5   4.5   -1  0.5    0    1

AI[4,] <- AI[3,]*(-2.5)+AI[4,]  
# 將合併矩陣之第三列乘以-2.5加到第四列
AI

##      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
## [1,]    1    0    0    9   -1    3   -1    0
## [2,]    0    1    0  -10    2   -3    1    0
## [3,]    0    0    1  -15    2   -5    2    0
## [4,]    0    0    0   42   -6   13   -5    1

AI[4,] <- AI[4,]*(1/42)         
# 將合併矩陣之第四列乘以1/42
AI

##      [,1] [,2] [,3] [,4]       [,5]       [,6]       [,7]       [,8]
## [1,]    1    0    0    9 -1.0000000  3.0000000 -1.0000000 0.00000000
## [2,]    0    1    0  -10  2.0000000 -3.0000000  1.0000000 0.00000000
## [3,]    0    0    1  -15  2.0000000 -5.0000000  2.0000000 0.00000000
## [4,]    0    0    0    1 -0.1428571  0.3095238 -0.1190476 0.02380952

AI[1,] <- AI[4,]*(-9)+AI[1,]    
# 將合併矩陣之第四列乘以-9加到第一列
AI

##      [,1] [,2] [,3] [,4]       [,5]       [,6]        [,7]        [,8]
## [1,]    1    0    0    0  0.2857143  0.2142857  0.07142857 -0.21428571
## [2,]    0    1    0  -10  2.0000000 -3.0000000  1.00000000  0.00000000
## [3,]    0    0    1  -15  2.0000000 -5.0000000  2.00000000  0.00000000
## [4,]    0    0    0    1 -0.1428571  0.3095238 -0.11904762  0.02380952

AI[2,] <- AI[4,]*(10)+AI[2,]    
# 將合併矩陣之第四列乘以10加到第二列
AI

##      [,1] [,2] [,3] [,4]       [,5]       [,6]        [,7]        [,8]
## [1,]    1    0    0    0  0.2857143  0.2142857  0.07142857 -0.21428571
## [2,]    0    1    0    0  0.5714286  0.0952381 -0.19047619  0.23809524
## [3,]    0    0    1  -15  2.0000000 -5.0000000  2.00000000  0.00000000
## [4,]    0    0    0    1 -0.1428571  0.3095238 -0.11904762  0.02380952

AI[3,] <- AI[4,]*(15)+AI[3,]    
# 將合併矩陣之第四列乘以15加到第三列
AI

##      [,1] [,2] [,3] [,4]       [,5]       [,6]        [,7]        [,8]
## [1,]    1    0    0    0  0.2857143  0.2142857  0.07142857 -0.21428571
## [2,]    0    1    0    0  0.5714286  0.0952381 -0.19047619  0.23809524
## [3,]    0    0    1    0 -0.1428571 -0.3571429  0.21428571  0.35714286
## [4,]    0    0    0    1 -0.1428571  0.3095238 -0.11904762  0.02380952

Page-18

計算 \(A\) 矩陣相對應的特徵值與特徵向量 (1)

A <- matrix(c(-1, 0, 0, 1),nrow=2)
eigen(A)

## eigen() decomposition
## $values
## [1]  1 -1
## 
## $vectors
##      [,1] [,2]
## [1,]    0   -1
## [2,]    1    0

計算 \(A\) 矩陣相對應的特徵值與特徵向量 (2)

A <- matrix(c(0.5, 0, 0, 0.5),nrow=2)
eigen(A)

## eigen() decomposition
## $values
## [1] 0.5 0.5
## 
## $vectors
##      [,1] [,2]
## [1,]    0   -1
## [2,]    1    0

計算 \(A\) 矩陣相對應的特徵值與特徵向量

A <- matrix(c(-4, 3, -6, 5),nrow=2)
eigen(A)

## eigen() decomposition
## $values
## [1]  2 -1
## 
## $vectors
##            [,1]       [,2]
## [1,]  0.7071068 -0.8944272
## [2,] -0.7071068  0.4472136

計算 \(A\) 矩陣相對應的特徵值與特徵向量

A <- matrix(c(5, 4, 2, 4, 5, 2, 2, 2, 2),nrow=3)
eigen(A)

## eigen() decomposition
## $values
## [1] 10  1  1
## 
## $vectors
##            [,1]       [,2]       [,3]
## [1,] -0.6666667  0.0000000  0.7453560
## [2,] -0.6666667 -0.4472136 -0.5962848
## [3,] -0.3333333  0.8944272 -0.2981424

Page-23~25 (僅供參考)

求矩陣 \(A\) 之 \(SVD\) 分解，設定 \(U\) 為 4x4， \(V\) 為 5x5

(A <- matrix(c(1,0,0,0,0,0,0,4,0,3,0,0,0,0,0,0,2,0,0,0),nrow=4))

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    1    0    0    0    2
## [2,]    0    0    3    0    0
## [3,]    0    0    0    0    0
## [4,]    0    4    0    0    0

svd(A,nu=4,nv=5)

## $d
## [1] 4.000000 3.000000 2.236068 0.000000
## 
## $u
##      [,1] [,2] [,3] [,4]
## [1,]    0    0    1    0
## [2,]    0    1    0    0
## [3,]    0    0    0    1
## [4,]   -1    0    0    0
## 
## $v
##      [,1] [,2]      [,3] [,4]       [,5]
## [1,]    0    0 0.4472136    0 -0.8944272
## [2,]   -1    0 0.0000000    0  0.0000000
## [3,]    0    1 0.0000000    0  0.0000000
## [4,]    0    0 0.0000000    1  0.0000000
## [5,]    0    0 0.8944272    0  0.4472136

求矩陣 \(A\) 之 \(SVD\) 分解，設定 \(U\) 為 4x4， \(V\) 為 3x3

(A <- matrix(c(1,0,0,3,0,0,4,0,0,1,0,2),nrow=4))

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    0    1
## [3,]    0    4    0
## [4,]    3    0    2

svd(A,nu=4,nv=3)

## $d
## [1] 4.000000 3.741657 1.000000
## 
## $u
##      [,1]       [,2]          [,3]       [,4]
## [1,]    0 -0.2223748  5.547002e-01 -0.8017837
## [2,]    0 -0.1482499 -8.320503e-01 -0.5345225
## [3,]    1  0.0000000  0.000000e+00  0.0000000
## [4,]    0 -0.9636241 -1.110223e-16  0.2672612
## 
## $v
##      [,1]       [,2]       [,3]
## [1,]    0 -0.8320503  0.5547002
## [2,]    1  0.0000000  0.0000000
## [3,]    0 -0.5547002 -0.8320503

自行構建一個 10x10 矩陣(影片例)

(A <- matrix(c(rep(0,20),rep(c(0,1,1,1,1,1,1,1,1,0),2),rep(c(0,1,1,0,0,0,0,1,1,0),2),rep(c(0,1,1,1,1,1,1,1,1,0),2),rep(0,20)),nrow=10))

##       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
##  [1,]    0    0    0    0    0    0    0    0    0     0
##  [2,]    0    0    1    1    1    1    1    1    0     0
##  [3,]    0    0    1    1    1    1    1    1    0     0
##  [4,]    0    0    1    1    0    0    1    1    0     0
##  [5,]    0    0    1    1    0    0    1    1    0     0
##  [6,]    0    0    1    1    0    0    1    1    0     0
##  [7,]    0    0    1    1    0    0    1    1    0     0
##  [8,]    0    0    1    1    1    1    1    1    0     0
##  [9,]    0    0    1    1    1    1    1    1    0     0
## [10,]    0    0    0    0    0    0    0    0    0     0

用 image 指令畫圖

image(A,col=topo.colors(256))

將 \(SVD\) 分解的結果存入 \(s\)

s <- svd(A,nu=10,nv=10)

取出奇異值做成矩陣

(D <- diag(s$d))

##           [,1]     [,2]         [,3]         [,4]       [,5]         [,6] [,7]
##  [1,] 6.040896 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
##  [2,] 0.000000 1.872853 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
##  [3,] 0.000000 0.000000 7.918465e-16 0.000000e+00 0.0000e+00 0.000000e+00    0
##  [4,] 0.000000 0.000000 0.000000e+00 1.275686e-16 0.0000e+00 0.000000e+00    0
##  [5,] 0.000000 0.000000 0.000000e+00 0.000000e+00 3.9807e-33 0.000000e+00    0
##  [6,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 1.670685e-49    0
##  [7,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
##  [8,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
##  [9,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
## [10,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.0000e+00 0.000000e+00    0
##       [,8] [,9] [,10]
##  [1,]    0    0     0
##  [2,]    0    0     0
##  [3,]    0    0     0
##  [4,]    0    0     0
##  [5,]    0    0     0
##  [6,]    0    0     0
##  [7,]    0    0     0
##  [8,]    0    0     0
##  [9,]    0    0     0
## [10,]    0    0     0

用 image 指令畫圖

image(D,col=topo.colors(256))

計算 \(A\) 的轉置矩陣 \(A^T\) 乘以 \(A\)

AA <- t(A)%*%A

將特徵值與特徵向量的計算結果存入 \(g\)

g  <- eigen(AA)

取出特徵值做成矩陣

(V <- diag(g$values))

##           [,1]     [,2]         [,3]        [,4]         [,5] [,6] [,7] [,8]
##  [1,] 36.49242 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##  [2,]  0.00000 3.507577 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##  [3,]  0.00000 0.000000 3.647876e-16 0.00000e+00 0.000000e+00    0    0    0
##  [4,]  0.00000 0.000000 0.000000e+00 1.97373e-16 0.000000e+00    0    0    0
##  [5,]  0.00000 0.000000 0.000000e+00 0.00000e+00 4.930381e-32    0    0    0
##  [6,]  0.00000 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##  [7,]  0.00000 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##  [8,]  0.00000 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##  [9,]  0.00000 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
## [10,]  0.00000 0.000000 0.000000e+00 0.00000e+00 0.000000e+00    0    0    0
##       [,9]         [,10]
##  [1,]    0  0.000000e+00
##  [2,]    0  0.000000e+00
##  [3,]    0  0.000000e+00
##  [4,]    0  0.000000e+00
##  [5,]    0  0.000000e+00
##  [6,]    0  0.000000e+00
##  [7,]    0  0.000000e+00
##  [8,]    0  0.000000e+00
##  [9,]    0  0.000000e+00
## [10,]    0 -4.749287e-16

用 image 指令畫圖

image(V,col=topo.colors(256))

擷取矩陣的一部分

(RA <- A[3:8,3:8])

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    1    1
## [2,]    1    1    0    0    1    1
## [3,]    1    1    0    0    1    1
## [4,]    1    1    0    0    1    1
## [5,]    1    1    0    0    1    1
## [6,]    1    1    1    1    1    1

用 image 指令畫圖

image(RA,col=topo.colors(256))

將 \(SVD\) 分解的結果存入 \(s\)

s <- svd(RA,nu=6,nv=6)

取出奇異值做成矩陣

(D <- diag(s$d))

##          [,1]     [,2]         [,3]         [,4]         [,5] [,6]
## [1,] 5.048676 0.000000 0.000000e+00 0.000000e+00 0.000000e+00    0
## [2,] 0.000000 1.584574 0.000000e+00 0.000000e+00 0.000000e+00    0
## [3,] 0.000000 0.000000 6.613243e-16 0.000000e+00 0.000000e+00    0
## [4,] 0.000000 0.000000 0.000000e+00 1.695558e-17 0.000000e+00    0
## [5,] 0.000000 0.000000 0.000000e+00 0.000000e+00 5.908398e-34    0
## [6,] 0.000000 0.000000 0.000000e+00 0.000000e+00 0.000000e+00    0

用 image 指令畫圖

image(D,col=topo.colors(256))

透過函數繪製影像圖 (修改自 R 範例)

x <- y <- seq(-4*pi, 4*pi, len = 500)
r <- sqrt(outer(x^2, y^2, "+"))
image(z = z <- cos(r^2)*exp(-r/6), col  = gray((0:256)/256))

image(z, axes = FALSE, main = "Math can be beautiful ...",
      xlab = expression(cos(r^2) * e^{-r/6}))

image(z, axes = FALSE, main = "Math can be beautiful ...",
      xlab = expression(cos(r^2) * e^{-r/6}),col=cm.colors(16))

Page-37 (僅供參考)

模擬醉漢走路，路寬左右各十米處有懸崖

x<-0 ; count<-0
repeat{x<- x + sample(c(-1,1),1) ; print(x); count<-count+1; if(x>=10 || x<=-10) {print(count);break}}

## [1] -1
## [1] 0
## [1] 1
## [1] 2
## [1] 3
## [1] 2
## [1] 1
## [1] 0
## [1] 1
## [1] 0
## [1] 1
## [1] 2
## [1] 1
## [1] 0
## [1] -1
## [1] -2
## [1] -3
## [1] -2
## [1] -1
## [1] -2
## [1] -3
## [1] -4
## [1] -3
## [1] -4
## [1] -5
## [1] -4
## [1] -5
## [1] -6
## [1] -5
## [1] -4
## [1] -3
## [1] -2
## [1] -3
## [1] -2
## [1] -1
## [1] 0
## [1] -1
## [1] 0
## [1] -1
## [1] -2
## [1] -1
## [1] 0
## [1] -1
## [1] -2
## [1] -3
## [1] -2
## [1] -1
## [1] 0
## [1] 1
## [1] 2
## [1] 1
## [1] 2
## [1] 3
## [1] 2
## [1] 1
## [1] 2
## [1] 1
## [1] 0
## [1] -1
## [1] 0
## [1] 1
## [1] 0
## [1] 1
## [1] 0
## [1] 1
## [1] 2
## [1] 3
## [1] 4
## [1] 5
## [1] 6
## [1] 7
## [1] 8
## [1] 7
## [1] 6
## [1] 7
## [1] 8
## [1] 9
## [1] 8
## [1] 9
## [1] 10
## [1] 80

另一種寫法，紀錄過程

x<-c() ; count<-0
repeat{x<- rbind(x, sample(c(-1,1),1)) ; count<-count+1; if(sum(x)>=10 || sum(x)<=-10) {print(x); print(count); break}}

##       [,1]
##  [1,]   -1
##  [2,]    1
##  [3,]    1
##  [4,]    1
##  [5,]    1
##  [6,]   -1
##  [7,]   -1
##  [8,]    1
##  [9,]    1
## [10,]   -1
## [11,]    1
## [12,]   -1
## [13,]    1
## [14,]   -1
## [15,]   -1
## [16,]    1
## [17,]    1
## [18,]    1
## [19,]   -1
## [20,]   -1
## [21,]    1
## [22,]    1
## [23,]    1
## [24,]   -1
## [25,]    1
## [26,]    1
## [27,]   -1
## [28,]   -1
## [29,]    1
## [30,]    1
## [31,]    1
## [32,]   -1
## [33,]   -1
## [34,]   -1
## [35,]   -1
## [36,]    1
## [37,]    1
## [38,]    1
## [39,]    1
## [40,]   -1
## [41,]   -1
## [42,]    1
## [43,]    1
## [44,]    1
## [45,]   -1
## [46,]    1
## [47,]    1
## [48,]    1
## [1] 48

為什麼你通常會輸光? (模擬100次，賭徒與莊家財力 10:100)

for ( i in 1:100) {
  x<-0 ; count<-0
  repeat{x<- x + sample(c(-1,1),1) ; count<-count+1; if(x>=100 || x<=-10) {cat("x=",x,"count=",count, "\n"); break}} }

## x= -10 count= 106 
## x= 100 count= 3324 
## x= -10 count= 208 
## x= -10 count= 18 
## x= -10 count= 88 
## x= -10 count= 240 
## x= -10 count= 1004 
## x= -10 count= 280 
## x= -10 count= 1062 
## x= 100 count= 5698 
## x= -10 count= 118 
## x= -10 count= 1014 
## x= -10 count= 968 
## x= -10 count= 276 
## x= -10 count= 726 
## x= -10 count= 34 
## x= -10 count= 268 
## x= -10 count= 50 
## x= -10 count= 738 
## x= -10 count= 30 
## x= -10 count= 72 
## x= -10 count= 50 
## x= -10 count= 250 
## x= -10 count= 168 
## x= -10 count= 480 
## x= -10 count= 416 
## x= -10 count= 162 
## x= -10 count= 300 
## x= -10 count= 48 
## x= -10 count= 194 
## x= -10 count= 2572 
## x= -10 count= 1090 
## x= 100 count= 4884 
## x= -10 count= 174 
## x= -10 count= 34 
## x= -10 count= 328 
## x= -10 count= 6706 
## x= -10 count= 18 
## x= -10 count= 1336 
## x= -10 count= 46 
## x= 100 count= 6514 
## x= -10 count= 230 
## x= -10 count= 72 
## x= -10 count= 84 
## x= 100 count= 2528 
## x= 100 count= 5362 
## x= -10 count= 180 
## x= -10 count= 18 
## x= -10 count= 100 
## x= -10 count= 124 
## x= -10 count= 66 
## x= -10 count= 116 
## x= -10 count= 108 
## x= -10 count= 972 
## x= -10 count= 138 
## x= -10 count= 36 
## x= -10 count= 112 
## x= -10 count= 38 
## x= -10 count= 1322 
## x= -10 count= 698 
## x= -10 count= 136 
## x= -10 count= 740 
## x= -10 count= 172 
## x= 100 count= 3762 
## x= -10 count= 1938 
## x= -10 count= 26 
## x= 100 count= 3096 
## x= -10 count= 92 
## x= -10 count= 120 
## x= -10 count= 1084 
## x= -10 count= 32 
## x= -10 count= 246 
## x= -10 count= 30 
## x= -10 count= 50 
## x= -10 count= 36 
## x= -10 count= 2146 
## x= -10 count= 52 
## x= -10 count= 850 
## x= -10 count= 166 
## x= -10 count= 4730 
## x= -10 count= 1088 
## x= 100 count= 4468 
## x= -10 count= 20 
## x= 100 count= 4844 
## x= -10 count= 32 
## x= -10 count= 50 
## x= -10 count= 236 
## x= -10 count= 112 
## x= -10 count= 60 
## x= -10 count= 1386 
## x= -10 count= 98 
## x= -10 count= 494 
## x= -10 count= 80 
## x= -10 count= 196 
## x= -10 count= 2686 
## x= -10 count= 224 
## x= -10 count= 998 
## x= -10 count= 410 
## x= -10 count= 58 
## x= -10 count= 698

為什麼你通常會輸光? (模擬100次，調高莊家的勝率，賭徒與莊家財力 10:100)

for ( i in 1:100) {
  x<-0 ; count<-0
  repeat{x<- x + sample(c(-1,1,-1,1,-1,1,-1,1,-1,1,-1),1) ; count<-count+1; if(x>=100 || x<=-10) {cat("x=",x,"count=",count, "\n"); break}} }

## x= -10 count= 42 
## x= -10 count= 36 
## x= -10 count= 136 
## x= -10 count= 54 
## x= -10 count= 368 
## x= -10 count= 62 
## x= -10 count= 50 
## x= -10 count= 22 
## x= -10 count= 74 
## x= -10 count= 34 
## x= -10 count= 224 
## x= -10 count= 128 
## x= -10 count= 108 
## x= -10 count= 56 
## x= -10 count= 84 
## x= -10 count= 162 
## x= -10 count= 386 
## x= -10 count= 128 
## x= -10 count= 34 
## x= -10 count= 176 
## x= -10 count= 88 
## x= -10 count= 60 
## x= -10 count= 70 
## x= -10 count= 70 
## x= -10 count= 112 
## x= -10 count= 24 
## x= -10 count= 116 
## x= -10 count= 44 
## x= -10 count= 40 
## x= -10 count= 40 
## x= -10 count= 240 
## x= -10 count= 220 
## x= -10 count= 70 
## x= -10 count= 48 
## x= -10 count= 56 
## x= -10 count= 36 
## x= -10 count= 160 
## x= -10 count= 90 
## x= -10 count= 88 
## x= -10 count= 100 
## x= -10 count= 104 
## x= -10 count= 62 
## x= -10 count= 74 
## x= -10 count= 140 
## x= -10 count= 122 
## x= -10 count= 70 
## x= -10 count= 48 
## x= -10 count= 116 
## x= -10 count= 94 
## x= -10 count= 64 
## x= -10 count= 144 
## x= -10 count= 106 
## x= -10 count= 32 
## x= -10 count= 42 
## x= -10 count= 88 
## x= -10 count= 104 
## x= -10 count= 32 
## x= -10 count= 56 
## x= -10 count= 202 
## x= -10 count= 416 
## x= -10 count= 36 
## x= -10 count= 56 
## x= -10 count= 72 
## x= -10 count= 160 
## x= -10 count= 28 
## x= -10 count= 44 
## x= -10 count= 54 
## x= -10 count= 104 
## x= -10 count= 126 
## x= -10 count= 50 
## x= -10 count= 34 
## x= -10 count= 98 
## x= -10 count= 170 
## x= -10 count= 16 
## x= -10 count= 18 
## x= -10 count= 36 
## x= -10 count= 128 
## x= -10 count= 32 
## x= -10 count= 52 
## x= -10 count= 90 
## x= -10 count= 46 
## x= -10 count= 748 
## x= -10 count= 24 
## x= -10 count= 84 
## x= -10 count= 72 
## x= -10 count= 42 
## x= -10 count= 48 
## x= -10 count= 588 
## x= -10 count= 104 
## x= -10 count= 26 
## x= -10 count= 52 
## x= -10 count= 62 
## x= -10 count= 148 
## x= -10 count= 86 
## x= -10 count= 30 
## x= -10 count= 44 
## x= -10 count= 64 
## x= -10 count= 26 
## x= -10 count= 98 
## x= -10 count= 52

為什麼你通常會輸光? (模擬100次，調整莊家的勝率至 501/1000，賭徒與莊家財力 10:100)

for ( i in 1:100) {
  x<-0 ; count<-0
  repeat{x<- x + sample(c(rep(1,499),rep(-1,501)),1) ; count<-count+1; if(x>=100 || x<=-10) {cat("x=",x,"count=",count, "\n"); break}} }

## x= 100 count= 15570 
## x= -10 count= 110 
## x= -10 count= 420 
## x= -10 count= 48 
## x= -10 count= 1062 
## x= -10 count= 44 
## x= -10 count= 1314 
## x= -10 count= 278 
## x= -10 count= 4454 
## x= -10 count= 1122 
## x= -10 count= 830 
## x= 100 count= 4046 
## x= -10 count= 98 
## x= -10 count= 182 
## x= -10 count= 156 
## x= 100 count= 2790 
## x= -10 count= 244 
## x= -10 count= 42 
## x= -10 count= 370 
## x= -10 count= 78 
## x= -10 count= 842 
## x= -10 count= 386 
## x= -10 count= 36 
## x= -10 count= 7740 
## x= -10 count= 140 
## x= -10 count= 58 
## x= -10 count= 62 
## x= -10 count= 488 
## x= 100 count= 1412 
## x= 100 count= 8102 
## x= -10 count= 72 
## x= -10 count= 150 
## x= -10 count= 114 
## x= -10 count= 684 
## x= -10 count= 44 
## x= -10 count= 282 
## x= -10 count= 996 
## x= -10 count= 34 
## x= -10 count= 26 
## x= 100 count= 4294 
## x= -10 count= 1000 
## x= -10 count= 1122 
## x= -10 count= 160 
## x= -10 count= 184 
## x= -10 count= 330 
## x= -10 count= 76 
## x= 100 count= 1724 
## x= -10 count= 88 
## x= -10 count= 50 
## x= -10 count= 30 
## x= -10 count= 226 
## x= -10 count= 1934 
## x= -10 count= 40 
## x= -10 count= 40 
## x= -10 count= 88 
## x= -10 count= 106 
## x= -10 count= 172 
## x= -10 count= 156 
## x= -10 count= 180 
## x= -10 count= 300 
## x= -10 count= 206 
## x= 100 count= 1526 
## x= -10 count= 86 
## x= -10 count= 36 
## x= -10 count= 94 
## x= -10 count= 38 
## x= -10 count= 686 
## x= -10 count= 92 
## x= -10 count= 32 
## x= -10 count= 80 
## x= -10 count= 4210 
## x= -10 count= 86 
## x= -10 count= 40 
## x= -10 count= 170 
## x= -10 count= 44 
## x= -10 count= 22 
## x= -10 count= 46 
## x= -10 count= 22 
## x= -10 count= 120 
## x= -10 count= 170 
## x= -10 count= 294 
## x= 100 count= 5390 
## x= -10 count= 152 
## x= -10 count= 46 
## x= -10 count= 1398 
## x= -10 count= 142 
## x= -10 count= 3028 
## x= -10 count= 278 
## x= -10 count= 138 
## x= -10 count= 38 
## x= 100 count= 1664 
## x= 100 count= 3586 
## x= 100 count= 596 
## x= -10 count= 608 
## x= -10 count= 106 
## x= -10 count= 120 
## x= -10 count= 210 
## x= -10 count= 1096 
## x= -10 count= 72 
## x= -10 count= 60

Page-38

一步轉移機率矩陣 (one-step transition probability matrix)

A<-matrix(c(.5,.4,.2,.3,.3,.4,.2,.3,.4),nrow=3, dimnames=list(c("晴","陰","雨"),c("晴","陰","雨")))
A

##     晴  陰  雨
## 晴 0.5 0.3 0.2
## 陰 0.4 0.3 0.3
## 雨 0.2 0.4 0.4

兩步轉移機率矩陣 (two-steps transition probability matrix)

A%*%A

##      晴   陰   雨
## 晴 0.41 0.32 0.27
## 陰 0.38 0.33 0.29
## 雨 0.34 0.34 0.32

# Ex: 0.41 = 0.5*0.5+0.3*0.4+0.2*0.2

三步轉移機率矩陣 (three-steps transition probability matrix)

A%*%A%*%A

##       晴    陰    雨
## 晴 0.387 0.327 0.286
## 陰 0.380 0.329 0.291
## 雨 0.370 0.332 0.298

一步轉移機率矩陣 (\(RR,RN,NR,NN\))

P <- matrix(c(.7,0,.5,0,.3,0,.5,0,0,.4,0,.2,0,.6,0,.8),nrow=4)
row.names(P) <- c("RR","RN","NR","NN")
colnames(P) <- c("RR","RN","NR","NN")
P

##     RR  RN  NR  NN
## RR 0.7 0.3 0.0 0.0
## RN 0.0 0.0 0.4 0.6
## NR 0.5 0.5 0.0 0.0
## NN 0.0 0.0 0.2 0.8

\(RRNRNNR \Rightarrow RR \rightarrow RN \rightarrow NR \rightarrow RN \rightarrow NN \rightarrow NR\)

(P1 <- 0.3*0.4*0.5*0.6*0.2)

## [1] 0.0072

\(RRRNNRN \Rightarrow RR \rightarrow RR \rightarrow RN \rightarrow NN \rightarrow NR \rightarrow RN\)

(P2 <- 0.7*0.3*0.6*0.2*0.5)

## [1] 0.0126

Page-45

Entropy 的計算可用定義函數的方式簡化，例如：

ent <- function(x) {sum(-x*log2(x))}
x <- c(rep(1/16,8),rep(1/48,24))
ent(x)

## [1] 4.792481

問題不複雜也可直接輸入算式計算：

(H <- -(8*(1/16)*log2(1/16)+24*(1/48)*log2(1/48)))

## [1] 4.792481

若 \(32\) 隊奪冠機會完全相同，則由哪對奪冠不確定性高，此時 Entropy 達極大值

(H=32*-(1/32)*log2(1/32))

## [1] 5

Page-47

計算 \(H(X)\)

(HX <- -((1/2)*log2(1/2)+(1/4)*log2(1/4)+2*(1/8)*log2(1/8)))

## [1] 1.75

計算 \(H(Y)\)

(HY <- -(4*(1/4)*log2(1/4)))

## [1] 2

計算 \(H(Y|X)\)

(HYIX <- -((1/8)*log2(1/4)+2*(1/16)*log2(1/8)+(1/4)*log2(1/2)+(1/8)*log2(1/2)+2*(1/16)*log2(1/4)+2*(1/32)*log2(1/4)+(1/16)*log2(1/2)+2*(1/32)*log2(1/4)+(1/16)*log2(1/2)))

## [1] 1.625

計算 \(H(X|Y)\)

(HXIY <- -((1/8)*log2(1/2)+(1/16)*log2(1/4)+2*(1/32)*log2(1/8)+(1/16)*log2(1/4)+(1/8)*log2(1/2)+2*(1/32)*log2(1/8)+4*(1/16)*log2(1/4)+(1/4)*log2(1)))

## [1] 1.375

計算 \(I(X;Y)\)

(IXY <- ((1/8)*log2(1)+(1/16)*log2(1)+2*(1/32)*log2(1)+(1/16)*log2(1/2)+(1/8)*log2(2)+2*(1/32)*log2(1)+(1/16)*log2(1/2)+(1/16)*log2(1)+2*(1/16)*log2(2)+(1/4)*log2(2)))

## [1] 0.375

另一作法

(IXY <- HX-HXIY)

## [1] 0.375

(IYX <- HY-HYIX)

## [1] 0.375

可看出 \(I(X;Y) = I(Y;X)\)

(IXY==IYX)

## [1] TRUE

計算 \(H(X,Y)\)

(HXY <- -(2*(1/8)*log2(1/8)+(1/4)*log2(1/4)+6*(1/16)*log2(1/16)+4*(1/32)*log2(1/32)))

## [1] 3.375

可看出 \(I(X;Y) = H(X)+H(Y)-H(X,Y)\)

(IXY==HX+HY-HXY)

## [1] TRUE

計算 \(I(X;X)\)

(IXX <- ((1/2)*log2((1/2)/(1/4))+(1/4)*log2((1/4)/(1/16))+2*(1/8)*log2((1/8)/(1/64))))

## [1] 1.75

可看出 \(I(X;X) = H(X)\)

(IXX==HX)

## [1] TRUE

Page-48 (僅供參考)

用微分說明 Entropy \(H(x)=-(p \log_2 p + (1-p) \log_2(1-p))\) maximized when \(p=\frac{1}{2}\)

(fx <- deriv(~ -(x*log2(x)+(1-x)*log2(1-x)),"x"))

## expression({
##     .expr1 <- log2(x)
##     .expr3 <- 1 - x
##     .expr4 <- log2(.expr3)
##     .expr8 <- log(2)
##     .value <- -(x * .expr1 + .expr3 * .expr4)
##     .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
##     .grad[, "x"] <- -(.expr1 + x * (1/(x * .expr8)) - (.expr3 * 
##         (1/(.expr3 * .expr8)) + .expr4))
##     attr(.value, "gradient") <- .grad
##     .value
## })

x <- 0.5
eval(fx)

## [1] 1
## attr(,"gradient")
##      x
## [1,] 0

Page-50 (僅供參考)

將 \(A、B\) 矩陣相乘並印出結果

A <- matrix(c(1,3,1,2,3,4,0,5,2,1,2,3,2,5,-1,3),nrow=4)
B <- matrix(c(2,1,2,3,1,3,1,4,0,1,-2,2,4,0,0,-1),nrow=4)
(AB <- A%*%B)

##      [,1] [,2] [,3] [,4]
## [1,]   15   20    3    2
## [2,]   27   36   12    7
## [3,]    3   -1   -6    5
## [4,]   24   32    5    5

將 \(A\) 拆成 4 個小矩陣，分別乘以 \(B\)

A1B <- A[1,]%*%B                   
A2B <- A[2,]%*%B
A3B <- A[3,]%*%B
A4B <- A[4,]%*%B

將個四乘出來的部分合併

(JA <- rbind(A1B,A2B,A3B,A4B))

##      [,1] [,2] [,3] [,4]
## [1,]   15   20    3    2
## [2,]   27   36   12    7
## [3,]    3   -1   -6    5
## [4,]   24   32    5    5

與 \(A、B\) 矩陣相乘結果比較，完全一致

(JA==AB)

##      [,1] [,2] [,3] [,4]
## [1,] TRUE TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE TRUE
## [4,] TRUE TRUE TRUE TRUE

MapReduce 分治

(A23B <- A[2,1]*B[1,3]+A[2,2]*B[2,3]+A[2,3]*B[3,3]+A[2,4]*B[4,3])

## [1] 12

與 \(A、B\) 矩陣相乘結果之第二列第三行元素，完全一致

(A23B==AB[2,3])

## [1] TRUE

用 \(K\) 做 Encryption

K <- matrix(c(1,2,2,5),nrow=2)
MATH <- matrix(c(13,1,20,8),nrow=2)
ENC <- K%*%MATH

\(KI\) As a key to Decryption

KI <- solve(K)                     
DEC <- KI%*%ENC

判斷解碼後結果與原始文字所代表矩陣是否一致

(DEC==MATH)

##      [,1] [,2]
## [1,] TRUE TRUE
## [2,] TRUE TRUE

測試 MapReduce 分治算法時間

A <- matrix(rnorm(1000000),nrow=1000)
B <- matrix(rnorm(1000000),nrow=1000)

run1<-function(){
  AB <- A%*%B ; print(AB[3,2])
}
run2<-function(){
  s<-0
  {for (i in 1:nrow(A))  {s <- s + A[3,i]*B[i,2] } 
  } 
  print(s)
}
system.time(run1())

## [1] -7.385027

## 使用者   系統   流逝 
##   0.83   0.00   0.84

system.time(run2())

## [1] -7.385027

## 使用者   系統   流逝 
##      0      0      0

######

資訊數學(下)講義範例 R 程式碼與執行結果

gcchen

2025/01/31

Page-01

Page-02

Page-03

Page-04~05

Page-06

Page-07

Page-08

Page-09

Page-11

Page-12

Page-13

Page-14~15 (僅供參考)

Page-16 (僅供參考)

Page-18

Page-23~25 (僅供參考)

Page-37 (僅供參考)

Page-38

Page-45

Page-47

Page-48 (僅供參考)

Page-50 (僅供參考)