P <- symbol("P")
Qd <- 100 -2*P
Qd\[100 - 2 P\]
พีชคณิต (Algebra) เป็นรากฐานทางคณิตศาสตร์ที่มีความสำคัญอย่างยิ่งต่อการศึกษาวิชาเศรษฐศาสตร์ โดยเฉพาะเมื่อนักเศรษฐศาสตร์ต้องการ วิเคราะห์ ปรับเปลี่ยน หรือพยากรณ์พฤติกรรมทางเศรษฐกิจ ผ่านการตั้งสมมุติฐานและแบบจำลองทางปริมาณ
วิธีทำ
P <- symbol("P")
Qd <- 100 -2*P
Qd\[100 - 2 P\]
แทนค่า P = 30 แล้วคำนวณ Q_d
subs(Qd, P, 30)\[40\]
วิธีทำ
Q <- symbol("Q")
TR <- 100 -2*Q
TR\[100 - 2 Q\]
แทนค่า Q = 200แล้วคำนวณ TR
subs(Qd, P, 30)\[40\]
วิธีทำ ตั้ง \(Q_d - Q_s=0\) แล้วแก้สมการเพื่อหา P
P <- symbol("P")
Eq <- 120 -3*P -(20+2*P)
Eq\[100 - 5 P\]
สามารถแก้สมการได้โดยใช้คำสั่ง solve_sys()
solve_sys(Eq,P)P = 20
หรือแก้สมการเชิงเส้นตัวแปรเดียวได้โดยใช้ R จึงจำเป็นต้องการแปลง ตัวแปร Eq ไปเป็นฟังก์ชันใน R แล้วใช้คำสั่ง uniroot() และควรต้องทราบช่วงที่คำตอบอยู่ด้วยโดยอาศัยการวาดกราฟช่วย
# แปลงเป็นไปฟังก์ชันใน R
Eq <- as_func(Eq)
# วาดกราฟ Qd-Qs เพื่อหาช่วงที่มีคำตอบของสมการ
curve(Eq, from = 10, to = 30, col = "red",
xlab = "P", ylab = "Qd-Qs")
grid()จากกราฟ ช่วงระหว่าง \(P\in [10,30]\) มีจุดที่ \(Q_d-Q_s =0\) อยู่ เราจึงสามารถแก้สมการด้วยคำสั่ง uniroot()
uniroot(Eq, lower = 10, upper = 30)$root[1] 20
lower และ upper คือขอบเขตของการค้นหาคำตอบ
วิธีทำ สมการเชิงเส้นของต้นทุนรวม (Total Cost Function) ซึ่งเขียนอยู่ในรูปแบบทั่วไปของต้นทุนคือ:
\[ TC = FC + VC \cdot Q \]
โดยที่
\(TC\) = ต้นทุนรวม
\(FC\) = ต้นทุนคงที่ (Fixed Cost)
\(VC\) = ต้นทุนแปรผันต่อหน่วย (Variable Cost per unit)
\(Q\) = ปริมาณที่ผลิต
จากโจทย์ที่ให้ \(TC = 500 + 25Q\)
\(500\) คือ ต้นทุนคงที่ (FC)
\(25\) คือ ต้นทุนแปรผันต่อหน่วย (VC)
คำถาม 1: เมื่อลงมือผลิต \(Q = 40\) หน่วย
คำนวณต้นทุนรวมโดยแทน Q ลงในสมการ:
Q <- symbol("Q")
TC <- 500 +25*Q
subs(TC, Q, 40)\[1500\]
คำถาม 2: ต้นทุนคงที่คือเท่าใด?
ดูจากสมการได้เลย: \[ FC = 500 \]
วิธีทำ
Y = รายได้ต่อวัน
W = ค่าจ้างรายชั่วโมง (บาท)
H = จำนวนชั่วโมงทำงาน (ชั่วโมง)
ในกรณีนี้ต้องแทน W และ H เท่ากับ 150 และ 8 ตามลำดับจึงใช้ฟังก์ชัน subs() ที่มี list ด้วย
W <- symbol("W")
H <- symbol("H")
Y <- W * H
# แทนค่า
subs(Y, list(W = 150, H = 8))\[1200\]
Q <- symbol("Q")
pi <- symbol("pi")
pi <- 70*Q -(30*Q +500)
pi\[40 Q - 500\]
# แทนค่า Q = 20
subs(pi,Q, 20)\[300\]
วิธีทำ คือ หาค่า \(x, y\) ที่เป็นจำนวนเต็มไม่ติดลบ \((x\geq 0, y\geq 0\) และทำให้สมการข้างต้นเป็นจริง เราสามารถ “แปลงสมการให้ y อยู่ในรูปของ x
\[ 10x + 20y = 200 \Rightarrow y = \frac{200 - 10x}{20} = 10 - 0.5x \]
เราสามารถแก้สมการโดยใช้ caracas ได้โดย
x <- symbol("x","positive" = TRUE)
y <- symbol("y","positive" = TRUE)
eq <- 10*x + 20*y -200
# แก้ y ในรูปของ x
Sol <- solve_sys(eq, y)
Sol[[1]]$y\[10 - \frac{x}{2}\]
ที่ใช้ Sol[[1]]$y เพื่อให้สามารถแสดงผลเอกสาร HTML หรือ PDF ได้
ดังนั้น เราสามารถใช้่สมการแทนค่า \(x =1,2,\cdots\) จะกระทั้งถึงค่าที่ \(y\leq 0\)
x <- 2*(1:10)
knitr::kable(data.frame(x = x, y =10-0.5*x))| x | y |
|---|---|
| 2 | 9 |
| 4 | 8 |
| 6 | 7 |
| 8 | 6 |
| 10 | 5 |
| 12 | 4 |
| 14 | 3 |
| 16 | 2 |
| 18 | 1 |
| 20 | 0 |
คำสั่ง kable() จากชุดคำสั่ง knitr ทำให้แสดงค่าออกอยู่รูปของตารางได้ในเอกสารแบบ HTML หรือ PDF และที่กำหนดให้ x เป็นเลขคู่ ก็เพราะต้องการใช้ y เป็นจำนวนเต็ม
วิธีทำ
\(Q\) = ผลผลิต (Output)
\(L\) = จำนวนแรงงาน (Labor)
\(K\) = จำนวนทุน (Capital)
L <- symbol("L")
K <- symbol("K")
Q <- 5*L + 10*K
# แทนค่า L = 4, K = 6
subs(Q, list(L = 4, K = 6))\[80\]
วิธีทำ
\(X\) = มูลค่าส่งออก
\(M\) = มูลค่านำเข้า
\(Y\) = รายได้แห่งชาติ (National Income)
\(NX = X - M\) = ดุลการค้า
ขั้นตอนที่ 1: เขียนสมการ NX
Y <- symbol("Y")
X <- 500 + 3*Y
M <- 200 + 2*Y
NX <- X - M # ได้สมการ NX = 300 + Y
NX\[Y + 300\]
ขั้นตอนที่ 2: เมื่อ \(Y = 100\)
# แทนค่า Y = 100
subs(NX, Y, 100)\[400\]
วิธีทำ สมการรายได้: \[ Y = C + I \]
โดยที่:
\(C = 100 + 0.8Y\) คือการบริโภค
\(I = 200\) คือการลงทุน
ต้องการหาค่า \(Y\) ที่ทำให้ระบบสมดุล
ขั้นที่ 1 แทนค่าลงในสมการ
Y <- symbol("Y")
C <- 100 + 0.8*Y
I <- 200
eq <- Y - (C + I) # ตั้งสมการ: Y = C + I
eq\[0.2 Y - 300\]
ขั้นที่ 2 แก้สมการ
solve_sys(eq, Y)Y = 1500.00000000000
ในทางเศรษฐศาสตร์ สมการไม่เชิงเส้น (Nonlinear Equations) พบได้บ่อยมาก โดยเฉพาะเมื่อพฤติกรรมเศรษฐกิจมีลักษณะซับซ้อน เช่น ดุลยภาพที่เปลี่ยนตามกำลังสอง, สมการอรรถประโยชน์, ฟังก์ชันการผลิตแบบ Cobb-Douglas, ฯลฯ
เราต้องแก้สมการ:
\[ 100 - 5P^2 = 20 \Rightarrow 5P^2 = 80 \Rightarrow P^2 = 16 \Rightarrow P = \pm 4 \]
แต่ทางเศรษฐศาสตร์ ราคา (P) ไม่ควรติดลบ
ดังนั้น \(P=4\)
P <- symbol("P", "positive" = TRUE)
Qd <- 100 -5*P^2
eq <- Qd -20 # ตั้ง Qd = 20
solve_sys(eq, P)P = 4
L <- symbol("L")
K <- symbol("K")
Q <- 10 * L^0.6 * K^0.4
Q_sub <- subs(Q, K, 25)
eq <- Q_sub - 200
solve_sys(eq, L)L = 17.2354775202551
คำตอบ: ค่าของ \(L\) ที่ทำให้ \(Q = 200\)
x <- symbol("x")
y <- symbol("y")
U <- x^0.5 * y^0.5
U_sub <- subs(U, x, 4)
eq <- U_sub - 10
solve_sys(eq, y)y = 25.0000000000000
คำตอบ: ค่าของ \(y\) ที่ทำให้ผู้บริโภคได้อรรถประโยชน์ 10
t <- symbol("t","positive" = TRUE)
Q <- 100 * exp(0.1 * t)
eq <- Q - 200
solve_sys(eq, t)t = 6.93147180559945
คำตอบ: เวลาที่ทำให้ผลผลิตเติบโตจาก 100 เป็น 200
เมตริกซ์ (Matrix) คือ การจัดเรียงตัวเลขหรือตัวแปรในรูปแบบตารางที่มี แถว (rows) และ คอลัมน์ (columns) โดยใช้สำหรับคำนวณทางคณิตศาสตร์ วิศวกรรม เศรษฐศาสตร์ และวิทยาศาสตร์ข้อมูล
เมตริกซ์ \(A\) ขนาด \(m \times n\) (m แถว, n คอลัมน์) เขียนได้เป็น:
\[ A = \begin{bmatrix} a_{11} & a_{12} & \cdots & a_{1n} \\ a_{21} & a_{22} & \cdots & a_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m1} & a_{m2} & \cdots & a_{mn} \end{bmatrix} \]
โดยที่ \(a_{ij}\) หมายถึงสมาชิกของเมตริกซ์ที่อยู่ใน แถวที่ i คอลัมน์ที่ j
ตัวอย่างเมตริกซ์
\[ B = \begin{bmatrix} 4 & 7 \\ 2 & 6 \end{bmatrix} \] เป็นเมตริกซ์ \(B\) ขนาด \(2 \times 2\)
หรือ
\[ A = \begin{bmatrix} 1 & 2 & 3 \\ 4 & 5 & 6 \end{bmatrix} \]
เมตริกซ์ \(A\) นี้มี 2 แถว และ 3 คอลัมน์
หรือเมตริกซ์ขนาด \(2 \times 3\)
| คำศัพท์ | ความหมาย |
|---|---|
| Elements | ตัวเลขหรือพจน์ที่อยู่ในเมตริกซ์ |
| Order หรือ Size | ขนาดของเมตริกซ์: แถว × คอลัมน์ |
| Square Matrix | เมตริกซ์ที่มีจำนวนแถว = คอลัมน์ |
| Row Vector | เมตริกซ์ที่มีเพียง 1 แถว |
| Column Vector | เมตริกซ์ที่มีเพียง 1 คอลัมน์ |
| Zero Matrix | เมตริกซ์ที่ทุกองค์ประกอบเป็น 0 |
| Identity Matrix | เมตริกซ์หน่วย เช่น \(\begin{bmatrix}1 & 0 \\ 0 & 1\end{bmatrix}\) |
แน่นอนครับ! ต่อไปนี้คือวิธีการทำ การดำเนินการทางคณิตศาสตร์ (arithmetic operations) สำหรับ เมตริกซ์ (matrix) ด้วยแพ็กเกจ caracas ใน R ซึ่งรองรับทั้งเมตริกซ์เชิงสัญลักษณ์ (symbolic matrix) และเชิงตัวเลข (numeric matrix)
# สร้างเมตริกซ์ A เชิงสัญลักษณ์ด้วย caracas
A <- matrix(c("a11", "a12", "a21", "a22"), nrow = 2, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}a_{11} & a_{12}\\a_{21} & a_{22}\end{matrix}\right]\]
# สร้างเมตริกซ์ B เชิงตัวเลขด้วย caracas
B <- matrix(c(1, 2, 3, 4), nrow = 2, byrow = TRUE)
B <- as_sym(B)
B\[\left[\begin{matrix}1 & 2\\3 & 4\end{matrix}\right]\]
# บวก
A_plus_B <- A + B
A_plus_B\[\left[\begin{matrix}a_{11} + 1 & a_{12} + 2\\a_{21} + 3 & a_{22} + 4\end{matrix}\right]\]
# ลบ
A_minus_B <- A - B
A_minus_B\[\left[\begin{matrix}a_{11} - 1 & a_{12} - 2\\a_{21} - 3 & a_{22} - 4\end{matrix}\right]\]
# คูณ A ด้วย 3
A_scaled <- 3 * A
A_scaled\[\left[\begin{matrix}3 a_{11} & 3 a_{12}\\3 a_{21} & 3 a_{22}\end{matrix}\right]\]
# A %*% B: เมตริกซ์คูณกันแบบมาตรฐาน
AB <- A %*% B
AB\[\left[\begin{matrix}a_{11} + 3 a_{12} & 2 a_{11} + 4 a_{12}\\a_{21} + 3 a_{22} & 2 a_{21} + 4 a_{22}\end{matrix}\right]\]
ใช้ %*% สำหรับ matrix multiplication แบบ linear algebra
# หาค่า A^2
A2 <- mat_pow(A, pow = 2)
A2\[\left[\begin{matrix}a_{11}^{2} + a_{12} a_{21} & a_{11} a_{12} + a_{12} a_{22}\\a_{11} a_{21} + a_{21} a_{22} & a_{12} a_{21} + a_{22}^{2}\end{matrix}\right]\]
A_t <- t(A)
A_t\[\left[\begin{matrix}a_{11} & a_{21}\\a_{12} & a_{22}\end{matrix}\right]\]
คำสั่ง t() สามารถใช้ได้กับ เมตริกซ์เชิงตัวเลขปกติจาก R
# หาดีเทอร์มิแนนต์
det_A <- det(A)
det_A\[a_{11} a_{22} - a_{12} a_{21}\]
# หาตัวผกผัน (ถ้า det ≠ 0)
inv_A <- inv(A)
inv_A\[\left[\begin{matrix}\frac{a_{22}}{a_{11} a_{22} - a_{12} a_{21}} & - \frac{a_{12}}{a_{11} a_{22} - a_{12} a_{21}}\\- \frac{a_{21}}{a_{11} a_{22} - a_{12} a_{21}} & \frac{a_{11}}{a_{11} a_{22} - a_{12} a_{21}}\end{matrix}\right]\]
คำสั่ง det() เป็นคำสั่งมาตราฐานสำหรับตัวแปรแบบเมตริกซ์ใน R สำหรับการหาเมตริกซ์ผกผันสำหรับตัวแปรแบบเมตริกซ์ใน R คือคำสั่ง solve()
ในหนังสือเล่มนี้สนใจเมทริกซ์ ที่มีจำนวนแถวเท่ากับจำนวนคอลัมน์ เพราะเป็นเมตริกซ์ที่สามารถใช้แก้ระบบสมการเชิงเส้นได้ (Linear Systems of Equations)
คุณสมบัติของเมตริกซ์ที่สำคัญ สำหรับการใช้ในการแก้ ระบบสมการเชิงเส้น (Linear Systems of Equations) ซึ่งเป็นพื้นฐานสำคัญในคณิตศาสตร์และเศรษฐศาสตร์
สมมติให้ระบบสมการเชิงเส้นมีรูปแบบ \[ A \cdot \mathbf{x} = \mathbf{b} \] โดยที่
\(A\) คือเมตริกซ์สัมประสิทธิ์ขนาด \(n \times n\)
\(\mathbf{x}\) คือเวกเตอร์ของตัวแปร
\(\mathbf{b}\) คือเวกเตอร์คำตอบ
| คุณสมบัติ | อธิบาย | เหตุผลที่สำคัญ |
|---|---|---|
| เป็นเมตริกซ์สี่เหลี่ยมจัตุรัส (Square Matrix) | \(A\) มีจำนวนแถว = จำนวนคอลัมน์ | เพื่อให้สามารถหาตัวผกผัน (Inverse) ได้ |
| ดีเทอร์มิแนนต์ไม่เป็นศูนย์ (Non-singular / \(\det(A) \neq 0\)) | เมตริกซ์ \(A\) ต้องมีดีเทอร์มิแนนต์ \(\neq\) 0 | บอกว่า \(A\) มีตัวผกผัน \(\rightarrow\) มีคำตอบที่แน่นอนเพียงหนึ่งชุด |
| มีอันดับ (Rank) เท่ากับจำนวนตัวแปร | Rank(\(A\)) = จำนวนตัวแปร | บอกว่าระบบมีคำตอบเฉพาะเจาะจง |
| สามารถหาเมตริกซ์ผกผันได้ (Invertible Matrix) | \(A^{-1}\) มีอยู่จริง | ใช้แก้สมการด้วยสูตร \(\mathbf{x} = A^{-1} \cdot \mathbf{b}\) |
| ไม่เกิดความสัมพันธ์เชิงเส้นระหว่างแถวหรือคอลัมน์ (Linearly Independent) | แถวหรือคอลัมน์ของ \(A\) ต้องไม่ซ้ำซ้อน | ป้องกันระบบสมการไม่สมบูรณ์หรือไม่มีคำตอบแน่นอน |
ตัวอย่าง เมตริกซ์ที่สามารถใช้แก้สมการได้
\[ A = \begin{bmatrix} 2 & 1 \\ 5 & 3 \end{bmatrix}, \quad \det(A) = 2(3) - 5(1) = 1 \neq 0 \]
เป็นเมตริกซ์ \(2 \times 2\)
มีดีเทอร์มิแนนต์ \(\neq\) 0 \(\rightarrow\) มีตัวผกผัน \(\rightarrow\) ใช้แก้สมการได้
สามารถใช้การแทนค่าด้วยคำสั่ง subs() เช่น ต้องการแทนค่า \(a_{11}\) ด้วย 0 ดังนั้น
A[1,1] <- 0
A\[\left[\begin{matrix}0 & a_{12}\\a_{21} & a_{22}\end{matrix}\right]\]
คำนวณค่า det(A) อีกครั้งจะได้
det(A)\[- a_{12} a_{21}\]
ถ้า แทนค่า \(a_{22}\) ด้วยตัวแปร z จะได้
A[2,2] <- symbol("z")
A\[\left[\begin{matrix}0 & a_{12}\\a_{21} & z\end{matrix}\right]\]
det(A)\[- a_{12} a_{21}\]
solve_lin()สามารถทำได้โดย
ตัวอย่าง \[ \begin{aligned} 3x + 2y =&~ 16 \\ x - y =&~ 4 \end{aligned} \]
เขียนในรูป \(A \cdot X = B\) คือ
\[ A = \begin{bmatrix} 3 & 2 \\ 1 & -1 \end{bmatrix}, \quad X = \begin{bmatrix} x \\ y \end{bmatrix}, \quad B = \begin{bmatrix} 16 \\ 4 \end{bmatrix} \]
# เมทริกซ์สัมประสิทธิ์
A <- matrix(c(3, 1, 2, -1), nrow = 2)
A <- as_sym(A)
A\[\left[\begin{matrix}3 & 2\\1 & -1\end{matrix}\right]\]
# เวกเตอร์คำตอบ
B <-matrix(c(16, 4), nrow = 2)
B <- as_sym(B)
B\[\left[\begin{matrix}16\\4\end{matrix}\right]\]
# แก้สมการ: A %*% X = B
X <- solve_lin(A, B)
X\[\left[\begin{matrix}\frac{24}{5}\\\frac{4}{5}\end{matrix}\right]\]
คำตอบคือเมตริกซ์ 2×1 ที่แสดง \(x\) และ \(y\) ในรูปสัญลักษณ์
ถ้าต้องการผลลัพธ์ออกมาเป็นตัวเลข สามารถใช้คำสั่ง N()
N(X, digits = 3)\[\left[\begin{matrix}4.8\\0.8\end{matrix}\right]\]
digits คือจำนวนทศนิยมสูงที่ต้องการแสดง แต่ในตัวอย่างเป็นตัวเลขทศนิยมลงตัวแสดงผลแค่เพียง ทศนิยม 1 ตำแหน่ง
มี 2 ปัจจัยการผลิต (แรงงาน \(L\), ทุน \(K\)) ที่ต้องใช้ในการผลิตสินค้า 2 ชนิด (X, Y) โดยต้องการให้ใช้ทรัพยากรทั้งหมด:
\(3X + 2Y = 18\) (ใช้แรงงานทั้งหมด)
\(4X + Y = 16\) (ใช้ทุนทั้งหมด)
จัดสมการอยู่ในรูปเมทริกซ์คือ
\[ \begin{bmatrix} 3&2\\4&1 \end{bmatrix} \begin{bmatrix} X\\Y \end{bmatrix} = \begin{bmatrix} 18\\16 \end{bmatrix}\rightarrow Ax=b \]
A <- matrix(c(3,2,4,1), nrow = 2, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}3 & 2\\4 & 1\end{matrix}\right]\]
b <- matrix(c(18,16), nrow = 2, byrow = TRUE)
b <- as_sym(b)
b\[\left[\begin{matrix}18\\16\end{matrix}\right]\]
solve_lin(A, b)\[\left[\begin{matrix}\frac{14}{5}\\\frac{24}{5}\end{matrix}\right]\]
กำหนเให้มี 2 สินค้า \(A\) และ \(B\) ที่ราคาของกันและกันมีผลต่ออุปสงค์
\[ \begin{cases} Q^d_A = 100 - 2P_A + P_B \\ Q^d_B = 80 - P_B + 0.5P_A \\ Q^s_A = -20 + 3P_A \\ Q^s_B = -10 + 2P_B \end{cases} \]
เงื่อนไขดุลยภาพ: \(Q^d_A = Q^s_A\), \(Q^d_B = Q^s_B\)
\[ \begin{aligned} 100 - 2P_A + P_B &= -20 + 3P_A \\ 80 - P_B + 0.5P_A &= -10 + 2P_B \end{aligned} \]
รูปเมตริกซ์ \[ \begin{bmatrix} -5 & 1 \\ 0.5 & -3 \end{bmatrix} \begin{bmatrix} P_A \\ P_B \end{bmatrix} = \begin{bmatrix} -120 \\ -90 \end{bmatrix} \]
ดังนั้น \[ \begin{bmatrix} P_A \\ P_B \end{bmatrix} = \begin{bmatrix} -5 & 1 \\ 0.5 & -3 \end{bmatrix}^{-1} \begin{bmatrix} -120 \\ -90 \end{bmatrix} \]
สามารถใช้ caracas สร้างเมตริกซ์เชิงสัญลักษณ์ แล้วหาเมตริกซ์ผกผัน เพื่อหาคำตอบได้่ดังนี้
A <- matrix(c(-5,1,0.5,-3),nrow = 2, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}-5 & 1\\0.5 & -3\end{matrix}\right]\]
B <- matrix(c(-120,-90),nrow = 2, byrow = TRUE)
B <- as_sym(B)
B\[\left[\begin{matrix}-120\\-90\end{matrix}\right]\]
ดังนั้น
solve_lin(A, B)\[\left[\begin{matrix}31.0344827586207\\35.1724137931035\end{matrix}\right]\]
หรือจะใช้
inv(A)%*%B\[\left[\begin{matrix}31.0344827586207\\35.1724137931035\end{matrix}\right]\]
โดยตรงก็ได้
ตัวอย่างเศรษฐกิจ 3 ภาค
การเกษตร, อุตสาหกรรม, บริการ
สมมุติว่าเมตริกซ์เทคโนโลยี \(A\) และเวกเตอร์ผลผลิตขั้นสุดท้าย \(\mathbf{d}\) เป็น: \[ A = \begin{bmatrix} 0.2 & 0.1 & 0.1 \\ 0.4 & 0.3 & 0.2 \\ 0.1 & 0.2 & 0.2 \end{bmatrix}, \quad \mathbf{d} = \begin{bmatrix} 100 \\ 200 \\ 150 \end{bmatrix} \]
ต้องการหาผลผลิตรวม \(\mathbf{x}\) โดยใช้ \[ (I - A)\mathbf{x} = \mathbf{d} \] ดังนั้น
\[ \mathbf{x} = (I - A)^{-1}\mathbf{d} \]
เราสามารถใช้ caracas สร้าง เมตริกซ์ \(I\) ขนาด \(3\times 3\) ได้ด้วย คำสั่ง diag() ใน R แบบปกติ ร่วมกับคำสั่ง as_sym()
I <- as_sym(diag(3))
I\[\left[\begin{matrix}1 & 0 & 0\\0 & 1 & 0\\0 & 0 & 1\end{matrix}\right]\]
ดังนั้น
A <- matrix(c(0.2,0.1,0.1,0.4,0.3,0.2,0.1,0.2,0.2),nrow = 3,byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}0.2 & 0.1 & 0.1\\0.4 & 0.3 & 0.2\\0.1 & 0.2 & 0.2\end{matrix}\right]\]
d <- matrix(c(100,200,150),nrow = 3,byrow = TRUE)
d <- as_sym(d)
d\[\left[\begin{matrix}100\\200\\150\end{matrix}\right]\]
solve_lin(I-A,d)\[\left[\begin{matrix}232.970027247957\\517.711171662127\\346.049046321526\end{matrix}\right]\]
สมมุติ รัฐมี 3 ประเภทการใช้จ่าย: ด้านสาธารณสุข \(G_1\), การศึกษา \(G_2\), โครงสร้างพื้นฐาน \(G_3\) มีข้อจำกัดงบประมาณ และข้อกำหนดนโยบาย เช่น
\(G_1 + G_2 + G_3 = 1000\)
\(G_1 = 1.5G_2\)
\(G_3 = 2G_2\)
เขียนในรูปสมการ \[ \begin{cases} G_1 + G_2 + G_3 = 1000 \\ G_1 - 1.5G_2 = 0 \\ G_3 - 2G_2 = 0 \end{cases} \quad \Rightarrow \quad \text{ระบบสมการ 3x3} \]
้น** คือ:
\[ A \cdot \mathbf{G} = \mathbf{b} \] รูปสมการเชิงเมตริกซ์:
\[ \underbrace{ \begin{bmatrix} 1 & 1 & 1 \\ 1 & -1.5 & 0 \\ 0 & -2 & 1 \end{bmatrix} }_{A} \cdot \underbrace{ \begin{bmatrix} G_1 \\ G_2 \\ G_3 \end{bmatrix} }_{\mathbf{G}} = \underbrace{ \begin{bmatrix} 1000 \\ 0 \\ 0 \end{bmatrix} }_{\mathbf{b}} \]
# สร้างเมตริกซ์ A
A <- matrix(c(1, 1, 1, 1, -1.5, 0, 0, -2, 1), nrow = 3, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}1 & 1 & 1\\1 & -1.5 & 0\\0 & -2 & 1\end{matrix}\right]\]
# เวกเตอร์ b
b <- matrix(c(1000, 0, 0), nrow = 3)
b <- as_sym(b)
b\[\left[\begin{matrix}1000\\0\\0\end{matrix}\right]\]
# แก้สมการ A * G = b
solve_lin(A, b)\[\left[\begin{matrix}333.333333333333\\222.222222222222\\444.444444444444\end{matrix}\right]\]
จะได้ค่าของ \(G_1, G_2, G_3\) ที่สอดคล้องกับเงื่อนไขงบประมาณและสัดส่วนการใช้จ่าย
แม้โดยทั่วไป IS-LM เป็นสมการไม่เชิงเส้น แต่สามารถ linearize ได้รอบจุดดุลยภาพ เช่น:
\[ \begin{cases} Y = C_0 + cY - bR + G \\ M/P = kY - hR \end{cases} \Rightarrow \text{จัดให้อยู่ในรูป } AX = B \]
แปลงให้อยู่ในรูปเชิงเส้น (linearized) ได้รอบจุดดุลยภาพ โดยใช้แบบจำลองง่าย คือรูปแบบ IS-LM แบบเชิงเส้น:
\[ \begin{cases} Y = C_0 + cY - bR + G \quad \text{(IS)} \\ \frac{M}{P} = kY - hR \quad \text{(LM)} \end{cases} \]
ซึ่งเป็นระบบสมการ 2 ตัวแปร: \(Y\) (รายได้), \(R\) (อัตราดอกเบี้ย) มีขั้นตอนการจัดรูปเป็นเมตริกซ์คือ
1. จัดสมการ IS: \[ Y = C_0 + cY - bR + G \Rightarrow (1 - c)Y + bR = C_0 + G \] 2. สมการ LM: \[ \frac{M}{P} = kY - hR \Rightarrow -kY + hR = -\frac{M}{P} \]
เขียนในรูปเมตริกซ์ \(A \cdot \begin{bmatrix} Y \\ R \end{bmatrix} = \mathbf{b}\)
\[ \underbrace{ \begin{bmatrix} 1 - c & b \\ -k & h \end{bmatrix} }_{A} \cdot \begin{bmatrix} Y \\ R \end{bmatrix} = \underbrace{ \begin{bmatrix} C_0 + G \\ \frac{M}{P} \end{bmatrix} }_{\mathbf{b}} \]
# สร้างเมตริกซ์ A
c <- 0.6; b <- 5; k <- 0.5; h <- 4
A <- matrix(c(1 - c, b, -k, h), nrow = 2, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}0.4 & 5\\-0.5 & 4\end{matrix}\right]\]
# สร้างเวกเตอร์ b
C0 <- 200; G <- 300; MP <- 800
B <- matrix(c(C0 + G, MP), nrow = 2)
B <- as_sym(B)
B\[\left[\begin{matrix}500\\800\end{matrix}\right]\]
# แก้ระบบสมการ
solve_lin(A, B)\[\left[\begin{matrix}-487.80487804878\\139.024390243902\end{matrix}\right]\]
เราจะได้คำตอบของ \(Y\) และ \(R\) ที่เป็นจุดดุลยภาพของระบบ IS-LM ภายใต้พารามิเตอร์ที่กำหนด
ประกอบไปด้วย 1. สมการไม่เชิงเส้น (Nonlinear Equation) คือ สมการที่ความสัมพันธ์ของตัวแปรไม่เป็นเส้นตรง เช่น
ตัวแปรถูกยกกำลัง (ยกกำลังมากกว่า 1 หรือเป็นทศนิยม)
ตัวแปรคูณกันเอง
มีฟังก์ชันไม่เชิงเส้น เช่น \(\sin\), \(\log\), \(e^x\)
สร้างสมการเชิงสัญลักษณ์
x <- symbol("x", "positive" = TRUE)
f <- x^2+3*x-5
f\[x^{2} + 3 x - 5\]
สามารถใช้คำสั่ง solve_sys() แก้สมการหาคำตอบได้
solution <-solve_sys(f, x)
solution[[1]]$x\[- \frac{3}{2} + \frac{\sqrt{29}}{2}\]
\(\log(x) + y = 3\)
\(xy = 10\)
\(e^x + x = 5\)
x <- symbol("x", "positive" = TRUE)
f <- exp(x)+x
f\[x + e^{x}\]
วาดกราฟ
fx<- as_func(f)
curve(fx,from = 0,to = 2, col= "blue")
curve(5+0*x, from = 0, to = 2, add = TRUE, col = "red")จากกราฟ คำตอบประมาณคือ 1.3
สามารถใช้คำสั่ง solve_sys() แก้สมการหาคำตอบได้
solution <- solve_sys(f-5, x)
N(solution[[1]]$x)c: 1.30655864103935
2. ระบบสมการไม่เชิงเส้น (Nonlinear System of Equations) คือ กลุ่มของสมการตั้งแต่ 2 สมการขึ้นไป ที่มีตัวแปรร่วมกัน และ อย่างน้อย 1 สมการเป็นสมการไม่เชิงเส้น
ลักษณะของระบบสมการไม่เชิงเส้น
มีตัวแปรตั้งแต่ 2 ตัวขึ้นไป
มีสมการตั้งแต่ 2 สมการขึ้นไป
ความสัมพันธ์ระหว่างตัวแปรมีลักษณะโค้งหรือซับซ้อน
การแก้สมการอาจมีคำตอบมากกว่า 1 คำตอบ หรือไม่มีเลย
ต้องใช้วิธีทางเชิงสัญลักษณ์หรือเชิงตัวเลขในการแก้
สร้างสมการเชิงสัญลักษณ์
x <- symbol("x")
y <- symbol("y")
f1 <- x^2+y^2 -25
f2 <- x*y-12\[x y - 12\]
หาคำตอบในเชิงตัวแปร y เพื่อใช้วาดกราฟ
solution.1 <- solve_sys(f1,y)
# คำตอบที่ 1
solution.1[[1]]$y\[- \sqrt{25 - x^{2}}\]
# คำตอบที่ 2
solution.1[[2]]$y\[\sqrt{25 - x^{2}}\]
solution.2 <- solve_sys(f2,y)
solution.2[[1]]$y\[\frac{12}{x}\]
วาดกราฟ
fx1.1 <- as_func(solution.1[[1]]$y)
fx1.2 <- as_func(solution.1[[2]]$y)
fx2 <- as_func(solution.2[[1]]$y)
curve(fx1.1, from = -5, to = 5, ylim = c(-6, 6), col = "red", xlim = c(-6, 6),
ylab = "y", xlab = "x")
curve(fx1.2, from = -5, to = 5, add =TRUE , col ="red")
curve(fx2, from = -5, to = -0.0001, add = TRUE, , col ="blue")
curve(fx2, from = 0.0001, to = 5, add = TRUE, , col ="blue")
grid()จะเห็นว่า มีคำตอบทั้่งหมด 4 คำตอบ
ในการใช้คำสั่งใน solve_sys() ใน caracas จะต้องทำการรวมสมการทั้งสองด้วยคำสั่ง cbind() ก่อน และอย่าลืม ใช้ list() จะสามารถหาคำตอบได้ดังนี้
solve_sys(cbind(f1, f2), list(x, y))Solution 1:
x = -4
y = -3
Solution 2:
x = -3
y = -4
Solution 3:
x = 3
y = 4
Solution 4:
x = 4
y = 3
สร้างสมการเชิงสัญลักษณ์
x <- symbol("x")
y <- symbol("y")
f1 <- log(x)+y-3
f1\[y + \log{\left(x \right)} - 3\]
f2 <- exp(y)+x-6
f2\[x + e^{y} - 6\]
หาคำตอบในเชิงตัวแปร y เพื่อใช้วาดกราฟ
solution.1 <- solve_sys(f1,y)
solution.1[[1]]$y\[3 - \log{\left(x \right)}\]
จะเห็นว่า x ต้องการกว่าศูนย์เท่านั้น
solution.2 <- solve_sys(f2,y)
solution.2[[1]]$y\[\log{\left(6 - x \right)}\]
คำตอบของสมการนี้ ค่า x ไม่เกิน 6 เท่านั้น ดังนั้น เมื่อนำช่วงทั้ง 2 มา intersection กันจะได้ \(x\in(0,6)\)
ลองวาดกราฟ
fx1 <- as_func(solution.1[[1]]$y)
fx2 <- as_func(solution.2[[1]]$y)
curve(fx1, from = 0.0001, to = 5.999, col = "red", ylab ="y", xlab ="x", ylim= c(-1,12))
curve(fx2, from = 0.0001, to = 5.999, col = "blue",add = TRUE)
grid()จากกราฟ จะพบว่า ระบบสมการนี้ไม่มีคำตอบ
สร้างสมการเชิงสัญลักษณ์
L <- symbol("L")
K <- symbol("K")
Q <- sqrt(L)*sqrt(K)-20
Q\[\sqrt{K} \sqrt{L} - 20\]
แทนค่า \(K=2L\) ลงไปจะได้
Q <- subs(Q, K, 2*L)
Q\[\sqrt{2} L - 20\]
ค่า \(L>0\) วาดกราฟ
fx <- as_func(Q)
# กราฟ Q
curve(fx, from = 0, to = 30, xlab = "L", ylab = "Q", col = "red")
# กราฟ y = 0
curve(0 + 0*x, from = 0, to = 30, xlab = "L", ylab = "Q", col = "blue", add = TRUE)
grid()จากกราฟ ระบบสมการนี้มีคำตอบ
solution <- solve_sys(Q,L)
solution[[1]]$L\[10 \sqrt{2}\]
สร้างสมการเชิงสัญลักษณ์
x <- symbol("x")
y <- symbol("y")
f1 <- sin(x)+y-1
f1\[y + \sin{\left(x \right)} - 1\]
f2 <- x^2+y^2 -2
f2\[x^{2} + y^{2} - 2\]
แก้สมการ f1 และ f2 เพื่อวาดกราฟ
solution.f1 <-solve_sys(f1,y)
solution.f1[[1]]$y\[1 - \sin{\left(x \right)}\]
f2 เป็นสมการกำลัง ต้องมี 2 คำตอบ
solution.f2 <-solve_sys(f2,y)
solution.f2[[1]]$y\[- \sqrt{2 - x^{2}}\]
solution.f2[[2]]$y\[\sqrt{2 - x^{2}}\]
fx1 <- as_func(solution.f1[[1]]$y)
fx2.1 <- as_func(solution.f2[[1]]$y)
fx2.2 <- as_func(solution.f2[[2]]$y)
curve(fx1, from = -sqrt(2), to = 1.5, col = "red", ylim = c(-2, 2), xlab = "x", ylab = "y")
curve(fx2.1, from = -sqrt(1.999999), to = sqrt(1.999999), col = "blue", add = TRUE)
curve(fx2.2, from = -sqrt(1.999999), to = sqrt(1.999999), col = "blue", add = TRUE)
grid()จากกราฟ จะเห็นว่า ระบบสมการนี้มี 2 คำตอบ
solve_sys(cbind(f1,f2), list(x, y))could not solve y - sin(sqrt(2 - y**2)) - 1
จะพบว่า ระบบสมการนี้ caracas แก้สมการนี้ไม่ได้ต้องใช้ต้องใช้คำสั่งอื่นๆในอาร์ ในหนังสือเล่มนี้ เลือกใช้ ชุดคำสั่ง nleqslv โดยนี้โค้ดดังนี้
library(nleqslv)
# นิยามระบบสมการ
f <- function(vars) {
x <- vars[1]
y <- vars[2]
c(sin(x) + y - 1,
x^2 + y^2 - 2)
}
# ค่าเริ่มต้นกราฟ เลือกที่ใกล้คำตอบ จุดแรก
start <- c(-.5, 1)
# แก้สมการ
sol <- nleqslv(start, f)
sol$x[1] -0.3727731 1.3641995
# ค่าเริ่มต้นกราฟ เลือกที่ใกล้คำตอบ จุดที่2
start <- c(1.2, 0)
# แก้สมการ
sol <- nleqslv(start, f)
sol$x[1] 1.41416057 0.01224232
การคำนวณเชิงสัญลักษณ์ด้วย caracas
ตัวอย่าง
โจทย์:
หาค่า eigenvalues และ eigenvectors ของเมทริกซ์ \[A = \begin{bmatrix}
2 & 1 \\
1 & 2
\end{bmatrix}\]
# สร้างเมทริกซ์ A
A <- matrix(c(2, 1,1, 2), nrow = 2, byrow = TRUE)
# ถ้าไม่กำหนด A เป็นตัวแปรจะไม่สามารถ คำนวณด้วยcaracas ได้
A <- as_sym(A)
A\[\left[\begin{matrix}2 & 1\\1 & 2\end{matrix}\right]\]
# นิยาม lambda (ค่าเฉพาะ)
lambda<- symbol("lambda_")
# หา Determinant ของ (A - lambda * I)
I <- diag(1, 2) # I คือ Identity Matrix 2x2
I <- as_sym(I)\[\left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]\]
char_eq <-det(A - lambda* I)
char_eq\[\lambda_{}^{2} - 4 \lambda_{} + 3\]
# แก้สมการ characteristic equation เพิ่อหา eigenvalues
eigenvals <- solve_sys(char_eq, lambda)
eigenvals[[1]]$lambda\[1\]
eigenvals[[2]]$lambda\[3\]
# สมมติเอาค่าเฉพาะค่าแรกมาแทนหา eigenvector
lambda <- eigenvals[[1]]$lambda # ดึงค่าแรก
# สร้างระบบสมการ (A - lambda1 * I) * v = 0
V <- as_sym(c("v1", "v2"))
# สมการคือ (A - lambda1*I) %*% v = 0
eqs <- (A - lambda*I) %*% V
# แก้สมการหาเวกเตอร์เฉพาะ
eigenvec1 <- solve_sys(eqs, V)
eigenvec1[[1]]$v1\[- v_{2}\]
ความหมายของคำตอบคือ ที่ค่าเฉพาะ \(\lambda=1\) จะได้ค่าเวคเตอร์คือ \[\begin{bmatrix}v_1\\v_2\end{bmatrix}=\begin{bmatrix}-v_2\\v_2\end{bmatrix}\] ถ้าให้ \(v_2=1\) คือ \[\begin{bmatrix}v_1\\v_2\end{bmatrix}=\begin{bmatrix}-1\\1\end{bmatrix}\]
# สมมติเอาค่าเฉพาะค่าที่สองมาแทนหา eigenvector
lambda <- eigenvals[[2]]$lambda # ดึงค่าแรก
# สร้างระบบสมการ (A - lambda1 * I) * v = 0
V <- as_sym(c("v1", "v2"))
# สมการคือ (A - lambda1*I) %*% v = 0
eqs <- (A - lambda*I) %*% V
# แก้สมการหาเวกเตอร์เฉพาะ
eigenvec2 <- solve_sys(eqs, V)
eigenvec2[[1]]$v1\[v_{2}\]
ความหมายของคำตอบคือ ที่ค่าเฉพาะ \(\lambda=1\) จะได้ค่าเวคเตอร์คือ \[\begin{bmatrix}v_1\\v_2\end{bmatrix}=\begin{bmatrix}v_2\\v_2\end{bmatrix}\] ถ้าให้ \(v_2=1\) คือ \[\begin{bmatrix}v_1\\v_2\end{bmatrix}=\begin{bmatrix}1\\1\end{bmatrix}\]
การคำนวณค่าเฉพาะและเวคเตอร์เจาะจงในตัวอย่างที่ผ่านมา จะการทดแทนการทำควณด้วยมือ เมื่อให้เกิดความเข้าใจในขั้นตอนการคำนวณ ถ้าเชี่ยวชาญดีแล้วสามารถใช้คำสั่งของ caracas โดยตรงเพื่อให้ได้ค่าเฉพาะและเวคเตอร์เจาะจงได้เลย คือ
eigenval() สำหรับหาค่าเฉพาะ
eigenvec() สำหรับหาค่าเวคเตอร์เจาะจง
เช่น จากตัวอย่างเมตริกซ์ A ที่ผ่านมา
eigenval(A)[[1]]
[[1]]$eigval
c: 3
[[1]]$eigmult
[1] 1
[[2]]
[[2]]$eigval
c: 1
[[2]]$eigmult
[1] 1
ความหมายของผลลัพธ์
[[1]]$eigval: Eigenvalue ตัวแรกคือ 3
[[1]]$eigmult: Eigenvalue ตัวนี้มีความซ้ำ (multiplicity) เท่ากับ 1
[[2]]$eigval: Eigenvalue ตัวที่สองคือ 1
[[2]]$eigmult: ตัวนี้ก็มี multiplicity = 1 เช่นกัน
eigenvec(A)[[1]]
[[1]]$eigval
c: 1
[[1]]$eigmult
[1] 1
[[1]]$eigvec
c: [-1 1]ᵀ
[[2]]
[[2]]$eigval
c: 3
[[2]]$eigmult
[1] 1
[[2]]$eigvec
c: [1 1]ᵀ
ความหมายของผลลัพธจากคำสั่ง eigenvec()
| ลำดับ | Eigenvalue (\(\lambda\)) | Eigenvector (\(v\)) | Multiplicity |
|---|---|---|---|
| 1 | 1 | \([-1,\ 1]^T\) | 1 |
| 2 | 3 | \([1,\ 1]^T\) | 1 |
ตัวอย่างการคำนวณเชิงสัญลักษณ์ กำหนดเมทริกซ์เทคนิคการผลิต: \[A = \begin{bmatrix} 0.3 & 0.2 \\ 0.4 & 0.1 \end{bmatrix}\]
ตรวจสอบว่าเศรษฐกิจนี้เสถียรหรือไม่ โดยดูจาก largest eigenvalue ของ A
# สร้างเมตริกซ์ A
A <- matrix(c(0.3, 0.2,
0.4, 0.1), nrow = 2, byrow = TRUE)
# เปลี่ยนเป็นตัวแปรเชิงสัญลักษณ์
A <- as_sym(A)
A\[\left[\begin{matrix}0.3 & 0.2\\0.4 & 0.1\end{matrix}\right]\]
# สร้างแปร lambda
lambda<- symbol("lambda_")
# สร้างเมตริกซ์เอกลักษณ์ I แบบสัญลักษณ์
I <- diag(1, 2)
I <- as_sym(I)
I\[\left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]\]
# สร้างสมการ (A-lambda*I) เชิงสัญลักษณ์
char_eq <- det(A - lambda*I)
char_eq\[\lambda_{}^{2} - 0.4 \lambda_{} - 0.05\]
## คำตอบมีแค่ 2 ค่า
factor_(char_eq)\[1.0 \left(1.0 \lambda_{} - 0.5\right) \left(1.0 \lambda_{} + 0.1\right)\]
# แก้สมการ
eigenvals <- solve_sys(char_eq, lambda)
eigenvals[[1]]$lambda\[-0.1\]
eigenvals[[2]]$lambda\[0.5\]
ถ้า largest eigenvalue \(\lambda_{max} < 1\) \(\rightarrow\) ระบบผลิตมีเสถียรภาพ
ถ้า \(\lambda_{max} \geq 1\) \(\rightarrow\) ระบบจะขยายตัวไม่สิ้นสุด (ไม่เสถียร)
จากสมการ ค่า \(\lambda_{max} =.5 <0 1\) ดังนั้น ระบบผลิตมีเสถียรภาพ
ตัวอย่างการคำนวณ พิจารณาเมทริกซ์การเปลี่ยนผ่านของระดับรายได้ \[P = \begin{bmatrix} 0.8 & 0.2 \\ 0.4 & 0.6 \end{bmatrix}\]
หาค่าสถานะสมดุล (stationary distribution) โดยแก้ \(\pi P = \pi\) ภายใต้เงื่อนไข \(\pi_1 + \pi_2 = 1\)
# สร้างตัวแปรเชิงสัญลักษณ์
P <- matrix(c(0.8, 0.2, 0.4, 0.6), nrow = 2, byrow = TRUE)
P <- as_sym(P)
P\[\left[\begin{matrix}0.8 & 0.2\\0.4 & 0.6\end{matrix}\right]\]
pi <- as_sym(c("pi1", "pi2"))
pi\[\left[\begin{matrix}\pi_{1}\\\pi_{2}\end{matrix}\right]\]
# สมการ: pi * P = pi $\rightarrow$ transpose ทั้งสองข้าง
eqs <- t(P) %*% pi - pi
eqs\[\left[\begin{matrix}- 0.2 \pi_{1} + 0.4 \pi_{2}\\0.2 \pi_{1} - 0.4 \pi_{2}\end{matrix}\right]\]
sol <-solve_sys(eqs, pi)
sol[[1]]$pi1\[2.0 \pi_{2}\]
แทนค่า \(\pi_1 =2\pi_2\)
subs(pi,pi[1], 2*pi[2])\[\left[\begin{matrix}2 \pi_{2}\\\pi_{2}\end{matrix}\right]\]
จากเงื่อนไข \(\pi_1+\pi_2 =1\) ดังนั้น
sum(subs(pi,pi[1], 2*pi[2]))-1\[3 \pi_{2} - 1\]
sol <- solve_sys(sum(subs(pi,pi[1], 2*pi[2]))-1, pi[2])
sol[[1]]$pi2\[\frac{1}{3}\]
จาก \(\pi_2=1/3\) ดังนั้น \(\pi_1 = 2/3\)
คำตอบคือ stationary distribution หรือความน่าจะเป็นของแต่ละสถานะในระยะยาว
rm(pi)ตรวจสอบ
pi[1] 3.141593
ตัวอย่างการคำนวณ ระบบเศรษฐกิจที่มีการเปลี่ยนแปลงแบบ: \[\frac{d}{dt} \begin{bmatrix} K \\ L \end{bmatrix} = \begin{bmatrix} 0.1 & 0.05 \\ 0.02 & 0.08 \end{bmatrix} \begin{bmatrix} K \\ L \end{bmatrix}\]
วิเคราะห์เสถียรภาพด้วย eigenvalue
# สร้างเมตริกซฺสัญลักษณ์ A
A <- matrix(c(0.1, 0.05, 0.02, 0.08), nrow = 2, byrow = TRUE)
A <- as_sym(A)
A\[\left[\begin{matrix}0.1 & 0.05\\0.02 & 0.08\end{matrix}\right]\]
# สร้างเมตริกซฺสัญลักษณ์ lambda
lambda<- symbol("lambda_")
# สร้างเมตริกซ์เอกลักษณ์
I <- diag(1, 2) |> as_sym()
I <- as_sym(I)
I\[\left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]\]
# สร้างสมการ det(A)-lambda*I =0
char_eq <- det(A - lambda* I)
char_eq\[\lambda_{}^{2} - 0.18 \lambda_{} + 0.007\]
## คำตอบมีแค่ 2 ค่า
factor_(char_eq)\[1.0 \left(1.0 \lambda_{}^{2} - 0.18 \lambda_{} + 0.007\right)\]
eigenvals <- solve_sys(char_eq, lambda)
eigenvals[[1]]$lambda\[0.056833752096446\]
eigenvals[[2]]$lambda\[0.123166247903554\]
ถ้า eigenvalues มีค่าบวก \(\rightarrow\) การเติบโตแบบ “unstable” (exponential growth)
ถ้าค่าลบ \(\rightarrow\) ระบบกลับสู่ดุลยภาพได้
ไม่มีค่าเฉพาะที่เป็นลบ แสดงว่า “unstable”
ตัวอย่างการคำนวณ พิจารณาเกมแบบ 2 ผู้เล่นที่มีเมทริกซ์ผลตอบแทน (payoff matrix): \[A = \begin{bmatrix} 2 & -1 \\ -1 & 2\end{bmatrix}\]
ต้องการวิเคราะห์ว่า กลยุทธ์เสถียร (stable strategy) อยู่ตรงไหน โดยพิจารณาทิศทางที่ระบบมีแนวโน้มจะเคลื่อนที่เข้าไปหาในระยะยาว
# สร้างเมตริกซ์สัญลักษณ์ A
A <- matrix(c(2, -1, -1, 2), nrow = 2, byrow = TRUE)
A<- as_sym(A)
A\[\left[\begin{matrix}2 & -1\\-1 & 2\end{matrix}\right]\]
# สร้างตัวแปร lambdaและเมติกซ์เอกลักษณ์
lambda<- symbol("lambda_")
I <- diag(1, 2) |> as_sym()
I <- as_sym(I)
I\[\left[\begin{matrix}1 & 0\\0 & 1\end{matrix}\right]\]
# สร้างสมการ det(A-lambda*I) =0
char_eq <- det(A - lambda* I)
char_eq\[\lambda_{}^{2} - 4 \lambda_{} + 3\]
## คำตอบมีแค่ 2 ค่า
factor_(char_eq)\[\left(\lambda_{} - 3\right) \left(\lambda_{} - 1\right)\]
eigenvals <- solve_sys(char_eq, lambda)
eigenvals[[1]]$lambda\[1\]
eigenvals[[2]]$lambda\[3\]
จากการคำนวณได้ eigenvalues คือ \(\lambda_1 = 3\), \(\lambda_2 = 1\)
ผลการวิเคราะห์:
\(\lambda_1 = 3 > 0\) \(\rightarrow\) บวก
\(\lambda_2 = 1 > 0\) \(\rightarrow\) บวก
ดังนั้น ระบบไม่เสถียร (Unstable)
แปลว่า: ถ้ามี perturbation เล็กๆ (การเบี่ยงเบนเล็กน้อย) ระบบจะเคลื่อน หนีออกจากดุลยภาพ ไม่ได้เคลื่อนกลับเข้าไปหาเอง
ตัวอย่างการคำนวณ พิจารณา network ของประเทศ 3 ประเทศที่มีความเชื่อมโยงทางการค้า (Trade Influence): สมมติความสัมพันธ์การค้า (Trade Influence)
ประเทศ A ส่งออก ไปยัง B และ C
ประเทศ B ส่งออก ไปยัง A และ C
ประเทศ C ส่งออก ไปยัง A และ B
ตารางแสดงการเชื่อมโยง
| A | B | C | |
|---|---|---|---|
| A | 0 | 1 | 1 |
| B | 1 | 0 | 1 |
| C | 1 | 1 | 0 |
เขียนเป็นเมทริกซ์ \(C\)
\[ C = \begin{bmatrix} 0 & 1 & 1 \\ 1 & 0 & 1 \\ 1 & 1 & 0 \end{bmatrix} \]
โดยที่
แถว 1 คือ ประเทศ A
แถว 2 คือ ประเทศ B
แถว 3 คือ ประเทศ C
และ
คอลัมน์ 1 คือ ประเทศ A
คอลัมน์ 2 คือ ประเทศ B
คอลัมน์ 3 คือ ประเทศ C
ต้องการหาว่า ประเทศใดมีความสำคัญที่สุด ในเครือข่าย โดยใช้ eigenvector centrality
# สร้างเมตริกซสัญลักษณ์ C
C <- matrix(c(0, 1, 1,
1, 0, 1,
1, 1, 0), nrow = 3, byrow = TRUE)
C <- as_sym(C)
C\[\left[\begin{matrix}0 & 1 & 1\\1 & 0 & 1\\1 & 1 & 0\end{matrix}\right]\]
# สร้างสัญลักษณ์ lambda และเมตริกซ์เอกลักษณ์ 3x3
lambda<- symbol("lambda_")
I <- diag(1, 3)
I <- as_sym(I)
I\[\left[\begin{matrix}1 & 0 & 0\\0 & 1 & 0\\0 & 0 & 1\end{matrix}\right]\]
# สร้างสมการ det(C)- lambda*I =0
char_eq <- det(C - lambda* I)
char_eq\[- \lambda_{}^{3} + 3 \lambda_{} + 2\]
## คำตอบมีแค่ 2 ค่า
factor_(char_eq)c: 2
-(λ - 2)⋅(λ + 1)
# หาค่า eigenvector
eigenvals <- solve_sys(char_eq, lambda)
eigenvals[[1]]$lambda\[-1\]
eigenvals[[2]]$lambda\[2\]
# เลือกค่า eigenvalue ใหญ่สุดมาหา eigenvector
V <- as_sym(paste0("v",1:3))
eqs <- (C - eigenvals[[2]]$lambda*I) %*% V
eqs\[\left[\begin{matrix}- 2 v_{1} + v_{2} + v_{3}\\v_{1} - 2 v_{2} + v_{3}\\v_{1} + v_{2} - 2 v_{3}\end{matrix}\right]\]
# แก้สมการหาค่า V
sol<- solve_sys(eqs, V)
sol[[1]]$v1\[v_{3}\]
sol[[1]]$v2\[v_{3}\]
จากผลการคำนวณพบว่า \(v_1=v_2=v_3\) แสดงว่าไม่มีใครมีอิพธิพลเหนือใคร
แน่นอนครับ! ต่อไปนี้คือแบบฝึกหัดพีชคณิตในเศรษฐศาสตร์ 20 ข้อเรียงเป็นรายการ (list) ที่ไม่ซ้ำกับในไฟล์ของคุณ:
ให้ใช้ caracasในการหาคำตอบ
สมการอุปทานคือ \(Q_s = 50 + 3P\) จงหาปริมาณอุปทานเมื่อราคาสินค้าเท่ากับ 10
สมการต้นทุนคือ \(C = 20Q + 150\) จงหาต้นทุนเมื่อผลิต 30 หน่วย
สมการรายได้รวมคือ \(TR = Q(120 - Q)\) จงหาค่ารายได้รวมเมื่อ \(Q = 40\)
ให้กำไร \(\Pi = TR - TC\) โดย \(TR = 100Q\), \(TC = 40Q + 300\) จงหากำไรเมื่อ \(Q = 20\)
สมการดุลยภาพ: \(Q_d = Q_s\), โดย \(Q_d = 120 - 2P\), \(Q_s = 30 + 4P\) จงหาราคาและปริมาณดุลยภาพ
รายได้เฉลี่ย (AR) คือ \(AR = \frac{TR}{Q}\), โดย \(TR = Q(80 - Q)\) จงหาค่า AR เมื่อ \(Q = 20\)
ถ้าต้นทุนเฉลี่ยคือ \(AC = \frac{TC}{Q}\), โดย \(TC = 5Q^2 + 100\) จงหาค่า AC เมื่อ \(Q = 10\)
ให้ \(TC = 60 + 10Q\), \(TR = 25Q\) จงเขียนสมการกำไร \(\Pi\) และหาค่ากำไรเมื่อ \(Q = 5\)
ผู้ผลิตมีฟังก์ชันอุปสงค์ \(Q = 200 - 5P\) จงหาค่าราคาเมื่อ \(Q = 100\)
ฟังก์ชันผลผลิตคือ \(Q = 5\sqrt{L} \cdot \sqrt{K}\) จงหาค่า \(Q\) เมื่อ \(L = 16\), \(K = 9\)
ให้สมการอรรถประโยชน์ \(U = XY\) จงหาค่า \(U\) เมื่อ \(X = 10\), \(Y = 8\)
พิจารณาสมการ \(TR = 80Q - 2Q^2\) จงหาค่ารายได้รวมเมื่อ \(Q = 30\)
ให้ต้นทุนเพิ่ม (MC) คือ \(MC = 6Q\) จงหาค่า MC เมื่อ \(Q = 12\)
ฟังก์ชันราคาคือ \(P = 100 - Q\) จงหาค่าราคาสินค้าเมื่อปริมาณขายคือ 30 หน่วย
ผู้บริโภคมีฟังก์ชันอรรถประโยชน์ \(U = \sqrt{X} \cdot \sqrt{Y}\) จงหาค่า \(U\) เมื่อ \(X = 25\), \(Y = 16\)
ให้ \(TC = Q^2 + 10Q + 100\) จงหาต้นทุนเมื่อ \(Q = 5\)
ฟังก์ชันการผลิตคือ \(Q = 10L - 0.5L^2\) จงหาค่า \(Q\) เมื่อ \(L = 10\)
พิจารณาสมการ \(R = P(200 - P)\) จงหาค่ารายได้เมื่อ \(P = 40\)
หาก \(TC = 2Q + 50\), \(TR = 5Q\) จงหาค่ากำไร \(\Pi\) เมื่อ \(Q = 10\)
หากฟังก์ชันอุปสงค์คือ \(Q = a - bP\), โดย \(a = 150\), \(b = 3\) จงเขียนสมการ Q แล้วหาค่า Q เมื่อ \(P = 20\)