Quantitative Finance Notes
Quantitative Finance Notes

Applied Stochastic Calculus

เมื่อสิ่งที่เราสนใจเป็นฟังก์ชันของตัวแปรสุ่ม เราจะหาการเปลี่ยนแปลงของมันอย่างไร

“ตรรกะและสัญชาตญาณต่างมีบทบาทที่จำเป็น ทั้งสองสิ่งขาดไปไม่ได้”

“Thus logic and intuition have each their necessary rôle. Each is indispensable.”

— Henri Poincaré · The Value of Science, บท I §V

ในบทนี้เราจะเข้าใจ Stochastic Calculus อย่าพึงถอดใจเมื่อเจอสิ่งที่ยาก Stochastic แปลว่าสุ่ม Calculus คือการศึกษาการเปลี่ยนแปลง เรเพียงแต่จะศึกษาการเปลี่ยนแปลงของระบบที่สุ่ม

บท Transition Density Functions ตอบว่า “จากจุดเริ่มต้นหนึ่ง อนาคตกระจายไปอย่างไร” แต่ในการเงิน เรามักไม่ได้สนใจตัวแปรนั้นเพียงตัวเดียว เราอาจสนใจ log ของราคา ราคายกกำลังสอง หรือมูลค่าสัญญาที่ขึ้นกับราคาและเวลา

ถ้า S เปลี่ยนแบบสุ่ม แล้ว F(S,t) จะเปลี่ยนตามกฎใด? คำตอบคือ Itô’s lemma ซึ่งเป็นกฎลูกโซ่สำหรับกระบวนการแบบ Itô และเป็นเครื่องมือหลักของบทนี้

1. time step ย่อครึ่งหนึ่ง แต่ขนาดช็อกไม่ได้ย่อครึ่งหนึ่ง

เริ่มจากเหรียญยุติธรรมและผลการโยนที่เป็นอิสระ ให้ Rᵢ เป็น +1 หรือ −1 ด้วยโอกาสอย่างละ 1/2 และ Cₙ เป็นยอดสะสมหลัง n ครั้ง

Cn=∑i=1nRi,E[Ri]=0,E[Ri2]=1,E[RiRj]=0(i≠j).C_n=\sum_{i=1}^{n}R_i,\qquad \mathbb E[R_i]=0,\qquad \mathbb E[R_i^2]=1,\qquad \mathbb E[R_iR_j]=0\quad(i\ne j).

ก่อนเริ่มทดลอง E[Cₙ] = 0 และ Var(Cₙ) = n แต่เมื่อรู้ยอดปัจจุบัน Cₘ แล้ว ค่าคาดหมายของยอดในอนาคตคือ E[Cₙ | Cₘ] = Cₘ เมื่อ n ≥ m ประวัติที่รู้เปลี่ยนคำตอบของยอดรวม แม้จะไม่เปลี่ยนโอกาสหัวหรือก้อยของครั้งถัดไป

ถ้าบีบ n ครั้งให้อยู่ในเวลารวม T และลดขนาดแต่ละ step เป็น ±√(T/n) จะมี

Δt=Tn,ΔW=±Δt,Var⁡ ⁣(∑i=1nΔWi)=T.\Delta t=\frac{T}{n},\qquad \Delta W=\pm\sqrt{\Delta t},\qquad \operatorname{Var}\!\left(\sum_{i=1}^{n}\Delta W_i\right)=T.

การเลือก ขนาด step ตาม √Δt ทำให้ความสุ่มไม่หายไปเมื่อแบ่งเวลาละเอียดขึ้น ลิมิตของการเดินสุ่มที่ปรับสเกลเหมาะสมนี้ให้ Brownian motion ส่วนในแบบจำลอง Brownian โดยตรง step ในช่วงเวลาที่มีขนาดจำกัดมีการแจกแจง

ΔWi=Wti+1−Wti∼N(0,Δt),ΔWi=Δt Zi,Zi∼iidN(0,1).\Delta W_i=W_{t_{i+1}}-W_{t_i} \sim\mathcal N(0,\Delta t), \qquad \Delta W_i=\sqrt{\Delta t}\,Z_i,\quad Z_i\overset{\text{iid}}{\sim}\mathcal N(0,1).

Brownian motion มีเส้นทางต่อเนื่อง step ในช่วงเวลาที่ไม่ทับกันเป็นอิสระ และ SD ของ step เท่ากับ √Δt สัญลักษณ์ dW ใช้ในสมการเชิงสุ่ม ไม่ใช่อนุพันธ์ธรรมดาที่เรานำไปหาร dt แล้วได้ความเร็วของเส้นทางที่เรียบ

ถ้าเวลาเป็นปี W มีหน่วย √ปี จึงต้องดูหน่วยของสัมประสิทธิ์ที่คูณ dW ด้วยเสมอ

2. step เล็กลง แต่ผลรวมกำลังสองยังเหลืออยู่

ให้แบ่งช่วง [0,T] เป็น n ช่วงเท่ากัน แล้ววัด quadratic variation บนตารางนี้ด้วย

Qn=∑i=0n−1(ΔWi)2.Q_n=\sum_{i=0}^{n-1}(\Delta W_i)^2.

ต้องแยกสองกรณีที่ดูคล้ายกันออกจากกัน

step ที่ใช้ กำลังสองของ step เดียว ผลรวม Qₙ
เหรียญ ±√Δt เท่ากับ Δt ทุกครั้ง เท่ากับ T พอดี แม้ n ยังจำกัด
step Gaussian ของ Brownian Δt Zᵢ² ยังเป็นตัวแปรสุ่ม ยังสุ่มที่ n จำกัด แต่เข้าใกล้ T ในลิมิต

สำหรับ Gaussian increments อิสระ เราใช้ E[Z²] = 1 และ E[Z⁴] = 3 จึงได้

E[Qn]=T,E[(Qn−T)2]=Var⁡(Qn)=2T2n⟶0.\mathbb E[Q_n]=T,\qquad \mathbb E[(Q_n-T)^2] =\operatorname{Var}(Q_n) =\frac{2T^2}{n}\longrightarrow0.

นี่คือการลู่เข้าแบบ mean square: ค่าเฉลี่ยของความคลาดเคลื่อนยกกำลังสองลดสู่ศูนย์ ไม่ได้หมายความว่าทุกเส้นทางต้องเข้าใกล้ T มากขึ้นในทุกครั้งที่เพิ่ม n

ผลนี้ย่อเป็นกฎจำของ stochastic calculus ว่า

(dW)2=dt,dt dW=0,(dt)2=0.(dW)^2=dt,\qquad dt\,dW=0,\qquad (dt)^2=0.

กฎนี้เป็นสัญลักษณ์ของการคำนวณในลิมิต ไม่ใช่สมการของ step Gaussian แต่ละ step ที่ Δt ยังไม่เป็นศูนย์ เราต้องใช้ (ΔW)² = Δt Z² ซึ่งไม่เท่ากับ Δt เสมอ

เทียบกับเส้นทางเรียบ เช่น x(t) = t จะได้ผลรวมกำลังสอง n(T/n)² = T²/n → 0 ความแตกต่างตรงนี้ทำให้แคลคูลัสของ Brownian ต้องเก็บพจน์อันดับสองไว้

3. ลองยกกำลังสอง แล้วดูสิ่งที่กฎลูกโซ่ธรรมดาทำหาย

เริ่มด้วยเอกลักษณ์พีชคณิตที่ถูกต้องสำหรับทุก step

Wti+12−Wti2=2WtiΔWi+(ΔWi)2.W_{t_{i+1}}^2-W_{t_i}^2 =2W_{t_i}\Delta W_i+(\Delta W_i)^2.

รวมทุก step และใช้ W₀ = 0 จะได้

WT2=2∑i=0n−1WtiΔWi+Qn.W_T^2=2\sum_{i=0}^{n-1}W_{t_i}\Delta W_i+Q_n.

ผลรวมตรงกลางใช้ค่า W ก่อน รับช็อกของช่วงถัดไป นี่เป็นจุดเริ่มต้นของ Itô integral ซึ่งใช้ข้อมูลที่มีอยู่แล้ว ไม่เอาค่าในอนาคตมาย้อนกำหนดน้ำหนักของ step ปัจจุบัน

เมื่อแบ่งเวลาละเอียดขึ้น ผลรวมฝั่งซ้ายของแต่ละช่วงให้

∫0TWt dWt=WT2−T2.\int_0^T W_t\,dW_t=\frac{W_T^2-T}{2}.

จึงได้กฎการเปลี่ยนแปลง

d(Wt2)=2Wt dWt+dt.\boxed{d(W_t^2)=2W_t\,dW_t+dt.}

ถ้าใช้กฎลูกโซ่ธรรมดาแล้วเขียนแค่ 2W dW เราจะทำ dt หายไป แต่ E[W_T²] = T จึงเห็นได้ว่าพจน์ที่หายมีความหมายจริง

สำหรับ integrand ที่ไม่มองอนาคตและมี E[∫₀ᵀ Hₜ² dt] < ∞ จะมี E[∫₀ᵀ Hₜ dWₜ] = 0 ไม่ควรสรุปว่าปริพันธ์เชิงสุ่มทุกแบบมีค่าคาดหมายศูนย์โดยไม่ดู integrand อ่านนิยามเพิ่มเติมใน Miranda Holmes-Cerfon, Lecture 7: Stochastic Integration

ทดลอง: ดูพจน์ที่สะสมจากกำลังสองของ step

ห้องทดลองใช้ Gaussian increments ของ Brownian path เดียวกันบนตารางละเอียด แล้วรวม step ติดกันเมื่อเลือก n ที่เล็กลง จึงเปรียบเทียบความละเอียดโดยไม่เปลี่ยนเส้นทางสุ่มพื้นฐาน

ลองเพิ่ม n แล้วเปลี่ยน seed ผลต่างของ Qₙ กับ T อาจขึ้นหรือลงในเส้นทางใดเส้นทางหนึ่ง แต่ค่าคลาดเคลื่อนแบบ RMS ตามทฤษฎีคือ T√(2/n) และผลรวมแบบ Itô บนตารางต้องตรงกับ (W_T² − Qₙ)/2 ก่อนพิจารณาลิมิต

4. Itô’s lemma: อ่าน drift และช็อกของฟังก์ชันใหม่

ให้ตัวแปร Y ตาม stochastic differential equation

dYt=a(Yt,t) dt+b(Yt,t) dWt.dY_t=a(Y_t,t)\,dt+b(Y_t,t)\,dW_t.

a เป็น drift ในหน่วย Y ต่อเวลา ส่วน b เป็นขนาดช็อกในหน่วย Y ต่อ √เวลา ค่า a(Yₜ,t) อาจสุ่มเพราะขึ้นกับสถานะ Yₜ คำว่า “ส่วนที่คาดการณ์ได้” จึงหมายถึงรู้สัมประสิทธิ์เมื่อกำหนดสถานะปัจจุบันแล้ว ไม่ได้หมายความว่าเส้นทางของ drift ทั้งหมดไม่สุ่ม

สำหรับฟังก์ชัน F(y,t) ที่หาอนุพันธ์ตามเวลาได้หนึ่งครั้งและตาม y ได้สองครั้ง โดยอนุพันธ์เหล่านี้ต่อเนื่อง หรือเรียกว่า C¹,² กฎ Itô ให้

dF(Yt,t)=(Ft+aFy+12b2Fyy)dt+bFy dWt.\boxed{ dF(Y_t,t)= \left(F_t+aF_y+\frac12 b^2F_{yy}\right)dt +bF_y\,dW_t. }

อนุพันธ์และสัมประสิทธิ์ทั้งหมดในสูตรประเมินที่ (Yₜ,t) สำหรับกระบวนการที่มีเงื่อนไขความเรียบและการมีอยู่ของคำตอบเหมาะสม

พจน์ สิ่งที่ทำให้ F เปลี่ยน
Fₜ dt เวลาเดินไป แม้ตรึงค่า Y
aFᵧ dt drift ของ Y ผ่านความชันของ F
½b²Fᵧᵧ dt ความผันผวนผ่านความโค้งของ F หรือ Itô correction
bFᵧ dW ช็อกที่เหลืออยู่ใน F

ที่มาของพจน์ใหม่คือการขยาย Taylor: เมื่อ dY มีขนาดช็อกตาม √dt พจน์ (dY)² จึงมีส่วน b²dt ที่ยังต้องเก็บไว้ หาก F เป็นเส้นตรงจะไม่มี Fᵧᵧ แต่เมื่อ F โค้ง พจน์แก้ไขนี้มักไม่เป็นศูนย์

เชื่อมกลับไปยังสมการการแพร่

ในบทก่อน dY = √2 c dW จึงมี a = 0 และ b = √2 c เมื่อเลือก F(Y) = Y² จะได้

d(Y2)=22 cY dW+2c2dt,E[YT2∣Yt=y]=y2+2c2(T−t).d(Y^2)=2\sqrt2\,cY\,dW+2c^2dt, \qquad \mathbb E[Y_T^2\mid Y_t=y]=y^2+2c^2(T-t).

เทอม c²Fᵧᵧ จึงเป็นตัวเดียวกับสัมประสิทธิ์ในสมการ Backward ของบทก่อน ส่วน 2c²(T−t) ตรงกับความแปรปรวนที่เราได้จาก Gaussian

ถ้า F เป็นมูลค่าสัญญา V(S,t) สูตรนี้ช่วยหาการเปลี่ยนแปลงของ V โดยใช้ความชันและความโค้ง แต่ Itô’s lemma เพียงอย่างเดียวยังไม่กำหนดราคา Option ต้องมีเงื่อนไข payoff, no arbitrage และแบบจำลองการคิดราคาที่สอดคล้องกันด้วย

5. ใช้กับ GBM: ทำไม log ของราคาจึงมี −½σ²

สมมติราคาตาม Geometric Brownian Motion ที่มี μ,σ คงที่ และ S₀ > 0

dSt=μSt dt+σSt dWt.dS_t=\mu S_t\,dt+\sigma S_t\,dW_t.

เมื่อเวลาเป็นปี μ มีหน่วยต่อปี และ σ มีหน่วยต่อ √ปี ในตัวอย่างจะใช้ μ = 0.10 และ σ = 0.20 ค่า μ เป็น drift ของผลตอบแทนแบบสัดส่วนในโมเดล ไม่ใช่ผลตอบแทนที่เส้นทางจริงต้องได้

เลือก F(S) = log S จะมี F_S = 1/S และ F_SS = −1/S² จึงได้

dlog⁡St=(μ−12σ2)dt+σ dWt.d\log S_t =\left(\mu-\frac12\sigma^2\right)dt+\sigma\,dW_t.

log เป็นฟังก์ชันเว้า อนุพันธ์อันดับสองจึงติดลบและทำให้ drift ของ log ต่ำกว่า μ เมื่อ σ ไม่เป็นศูนย์ นี่คือผลของ Itô correction ไม่ใช่ค่าธรรมเนียมหรือการหักผลตอบแทนตามใจ

อินทิเกรตจาก 0 ถึง t โดยใช้ W₀ = 0:

log⁡StS0=(μ−12σ2)t+σWt,St=S0exp⁡ ⁣[(μ−12σ2)t+σWt].\log\frac{S_t}{S_0} =\left(\mu-\frac12\sigma^2\right)t+\sigma W_t, \qquad \boxed{S_t=S_0\exp\!\left[\left(\mu-\frac12\sigma^2\right)t+\sigma W_t\right].}

คำตอบนี้เป็นบวกทุกเวลาจำกัดเมื่อ S₀ > 0 และทำให้ log return เป็น Normal ส่วนราคาเป็น Lognormal ไม่ใช่เพราะ simple return ถูกสมมติให้เป็น Normal แล้วราคาจะเป็น Lognormal โดยทันที

สำหรับ σ > 0 และ t > 0:

log⁡StS0∼N ⁣((μ−12σ2)t,σ2t).\log\frac{S_t}{S_0} \sim\mathcal N\!\left(\left(\mu-\frac12\sigma^2\right)t,\sigma^2t\right).

เมื่อ S₀ = 100, μ = 0.10, σ = 0.20 และ t = 1 ปี จะมี

สิ่งที่วัด สูตร ค่า
ค่าคาดหมายของราคา S₀e^(μt) 110.5171
มัธยฐานของราคา S₀e^((μ−σ²/2)t) 108.3287
ค่าคาดหมายของ log return (μ−σ²/2)t 0.08
SD ของ log return σ√t 0.20

ค่าเฉลี่ยกับมัธยฐานต่างกันเพราะการแจกแจงราคาเบ้ ค่า 108.3287 จึงไม่ใช่ E[S₁] และไม่มีค่าใดเป็นการทำนายราคาจริง ณ สิ้นปี

อีกตัวอย่าง: ราคายกกำลังสอง

เลือก F(S) = S² มี F_S = 2S และ F_SS = 2 ดังนั้น

d(St2)=(2μ+σ2)St2 dt+2σSt2 dWt.d(S_t^2)=(2\mu+\sigma^2)S_t^2\,dt+2\sigma S_t^2\,dW_t.

คราวนี้ฟังก์ชันนูนทำให้ Itô correction เป็นบวก สังเกตเครื่องหมายเทียบกับกรณี log S แล้วจะเห็นบทบาทของความโค้งโดยไม่ต้องท่องสูตรแยกกัน

6. ใช้กับ mean reversion: กลับเข้าหาค่าเฉลี่ยได้โดยยังมีช็อก

GBM ให้ขนาดช็อกโตตามราคา แต่อีกกลุ่มหนึ่งของแบบจำลองมีแรงดึงกลับเข้าหาค่ากลางหรือ mean reversion ตัวอย่างคือ Ornstein–Uhlenbeck process หรือ Vasicek model เมื่อใช้กับอัตราดอกเบี้ยระยะสั้น

drt=κ(θ−rt) dt+σr dWt,κ>0.dr_t=\kappa(\theta-r_t)\,dt+\sigma_r\,dW_t, \qquad \kappa>0.

θ เป็นระดับระยะยาว κ เป็นความเร็วในการกลับเข้าหา θ และ σᵣ เป็นขนาดช็อกของระดับ r สมมติทั้งสามเป็นค่าคงที่สำหรับคำตอบในหัวข้อนี้ ในบทนี้เปลี่ยน γ และ r̄ ของต้นฉบับเป็น κ และ θ เพื่อให้อ่านง่าย

เมื่อ r ต่ำกว่า θ พจน์ drift เป็นบวก เมื่อ r สูงกว่า θ drift เป็นลบ แต่ช็อกยังอาจผลักให้เดินออกห่างได้ mean reversion เป็นแรงดึงใน drift ไม่ใช่คำสัญญาว่าทุก step ต้องเข้าใกล้ค่าเฉลี่ย

ตั้ง uₜ = rₜ − θ จะได้ du = −κu dt + σᵣ dW แล้วใช้ integrating factor e^(κt):

d(eκtut)=σreκt dWt.d(e^{\kappa t}u_t)=\sigma_r e^{\kappa t}\,dW_t.

จึงได้คำตอบ

rt=θ+(r0−θ)e−κt+σr∫0te−κ(t−s) dWs.r_t=\theta+(r_0-\theta)e^{-\kappa t} +\sigma_r\int_0^t e^{-\kappa(t-s)}\,dW_s.

สำหรับ r₀ ที่กำหนดแน่นอนและ σᵣ > 0 คำตอบมีการแจกแจง Normal โดย

E[rt∣r0]=θ+(r0−θ)e−κt,Var⁡(rt∣r0)=σr22κ(1−e−2κt).\begin{aligned} \mathbb E[r_t\mid r_0]&=\theta+(r_0-\theta)e^{-\kappa t},\\ \operatorname{Var}(r_t\mid r_0)&=\frac{\sigma_r^2}{2\kappa}\left(1-e^{-2\kappa t}\right). \end{aligned}

เมื่อเวลายาวขึ้น อิทธิพลของ r₀ ลดลงและการแจกแจงเข้าใกล้ stationary distribution

r∞∼N ⁣(θ,σr22κ).r_\infty\sim\mathcal N\!\left(\theta,\frac{\sigma_r^2}{2\kappa}\right).

นี่หมายถึงการแจกแจงระยะยาวคงรูป ไม่ได้หมายความว่าเส้นทางแต่ละเส้นหยุดสุ่ม หากเริ่มจาก r₀ คงที่ กระบวนการยังไม่อยู่ใน stationary distribution ตั้งแต่เวลา 0

ตัวอย่างสมมติให้ r₀ = 5%, θ = 3%, κ = 1.5 ต่อปี และ σᵣ = 0.01 ต่อ √ปี หรือ 1 จุดเปอร์เซ็นต์ต่อ √ปี หลังหนึ่งปี ค่าคาดหมายประมาณ 3.4463% และ SD ประมาณ 0.5628 จุดเปอร์เซ็นต์ ส่วน SD ระยะยาวประมาณ 0.5774 จุดเปอร์เซ็นต์

OU/Vasicek มีหาง Gaussian บนเส้นจำนวนจริง จึงยอมให้อัตราดอกเบี้ยติดลบได้ การเลือกใช้ต้องพิจารณาว่าสมบัตินี้เหมาะกับสิ่งที่กำลังจำลองหรือไม่

7. SDE หนึ่งสมการ เชื่อมกับความหนาแน่นสองมุมมอง

กลับไปที่ตัวแปรทั่วไป dY = a(Y,t)dt + b(Y,t)dW และใช้สัญลักษณ์เดิมจากบทก่อน k(y,t;z,T) เมื่อมี transition density และความเรียบเพียงพอ สมการ Forward คือ

∂k∂T=−∂∂z[a(z,T)k]+12∂2∂z2[b(z,T)2k].\boxed{ \frac{\partial k}{\partial T} =-\frac{\partial}{\partial z}\big[a(z,T)k\big] +\frac12\frac{\partial^2}{\partial z^2}\big[b(z,T)^2k\big]. }

Forward ตรึงจุดเริ่มต้น y,t แล้วติดตามการแจกแจงปลายทาง สัมประสิทธิ์อยู่ภายในอนุพันธ์ เพราะอาจเปลี่ยนตามตำแหน่ง z ส่วน Backward ตรึงปลายทาง z,T และทำงานกับจุดเริ่มต้น:

∂k∂t+a(y,t)∂k∂y+12b(y,t)2∂2k∂y2=0.\boxed{ \frac{\partial k}{\partial t} +a(y,t)\frac{\partial k}{\partial y} +\frac12b(y,t)^2\frac{\partial^2k}{\partial y^2}=0. }

ถ้า a = 0 และ b = √2 c คงที่ ทั้งคู่จะกลับเป็นสมการในบท Transition Density ทันที แต่เมื่อสัมประสิทธิ์ขึ้นกับสถานะ จะสลับตัวแปรและเปลี่ยนเครื่องหมายอย่างเดียวไม่ได้

จาก log ของ GBM กลับเป็นความหนาแน่นของราคา

เพราะ log return เป็น Normal เราจึงแปลงกลับเป็นความหนาแน่นของราคาปลายทาง z > 0 ได้ สำหรับ τ = T − t > 0, ราคาเริ่มต้น s > 0 และ σ > 0:

k(s,t;z,T)=1zσ2πτexp⁡ ⁣[−(log⁡(z/s)−(μ−σ2/2)τ)22σ2τ].k(s,t;z,T) =\frac{1}{z\sigma\sqrt{2\pi\tau}} \exp\!\left[ -\frac{\left(\log(z/s)-(\mu-\sigma^2/2)\tau\right)^2}{2\sigma^2\tau} \right].

ตัวคูณ 1/z มาจากการเปลี่ยนตัวแปร log z กลับเป็น z จะลืมพจน์นี้ไม่ได้ เพราะพื้นที่ใต้ความหนาแน่นราคาต้องยังเป็น 1

สำหรับ OU, drift ดึงมวลกลับเข้าหา θ ขณะที่ diffusion กระจายมวลออก จึงมีการแจกแจงระยะยาวที่สมดุลกันได้ การมีสัมประสิทธิ์ที่ไม่ขึ้นกับเวลาเพียงอย่างเดียวไม่ได้รับประกันว่าจะมี stationary distribution สำหรับทุกแบบจำลอง

สมการเหล่านี้ต้องใช้ร่วมกับเงื่อนไขเริ่มต้น ขอบเขต และกฎความน่าจะเป็นที่เลือก การหาความหนาแน่นภายใต้แบบจำลองยังไม่เท่ากับการกำหนดราคาสัญญาภายใต้ risk-neutral measure อ่านเรื่องตัวดำเนินการและเงื่อนไขต่อได้ใน Lecture 10: Forward and Backward Equations for SDEs

8. จากสมการต่อเนื่อง สู่ขั้นตอนที่คอมพิวเตอร์คำนวณได้

วิธี Euler–Maruyama ประมาณการเปลี่ยนแปลงโดยประเมินสัมประสิทธิ์ที่ต้นช่วง แล้วใช้ช็อกใหม่ในแต่ละ step:

Yi+1=Yi+a(Yi,ti)Δt+b(Yi,ti)Δt Zi,Zi∼iidN(0,1).Y_{i+1}=Y_i+a(Y_i,t_i)\Delta t +b(Y_i,t_i)\sqrt{\Delta t}\,Z_i, \qquad Z_i\overset{\text{iid}}{\sim}\mathcal N(0,1).

ช็อกในสูตรนี้คูณด้วย √Δt หากใช้ Δt แทน ความแปรปรวนรวมจะยุบลงเมื่อแบ่งเวลาละเอียดขึ้นและได้คนละแบบจำลอง

สำหรับ GBM สูตร Euler เป็น

Si+1Euler=SiEuler(1+μΔt+σΔt Zi).S_{i+1}^{\mathrm{Euler}} =S_i^{\mathrm{Euler}}\left(1+\mu\Delta t+\sigma\sqrt{\Delta t}\,Z_i\right).

ใช้ค่าจากตัวอย่างในเอกสาร: S₀ = 100, μ = 0.1, σ = 0.2, Δt = 0.01 ปี และช็อกสองค่า 0.12, −0.25 จะได้

step การคำนวณ ราคาหลังจบ step
1 100 × (1 + 0.1×0.01 + 0.2×0.1×0.12) 100.34
2 100.34 × (1 + 0.1×0.01 + 0.2×0.1×(−0.25)) 99.93864

ต้นฉบับปัดค่าที่สองเป็น 99.94 สำหรับการแสดงผล แต่การคำนวณต่อควรเก็บความละเอียดไว้

เมื่อมีคำตอบ exact ให้ใช้เป็นตัวเทียบ

GBM ที่ μ,σ คงที่มีสูตรเปลี่ยนสถานะที่แน่นอนบนจุดเวลาที่เลือก:

Si+1exact=Siexactexp⁡ ⁣[(μ−12σ2)Δt+σΔt Zi].S_{i+1}^{\mathrm{exact}} =S_i^{\mathrm{exact}} \exp\!\left[\left(\mu-\frac12\sigma^2\right)\Delta t +\sigma\sqrt{\Delta t}\,Z_i\right].

Euler กับ exact ควรใช้ ช็อกเดียวกัน เพื่อดูความคลาดเคลื่อนจากวิธีคำนวณ หากสุ่มคนละเส้นทาง ผลต่างจะปนทั้งความสุ่มและความคลาดเคลื่อนของวิธี

คำว่า exact ในที่นี้หมายถึงคำตอบของ GBM ที่สมมติไว้ ณ จุดเวลาที่คำนวณ เส้นเชื่อมจุดบนกราฟเป็นภาพช่วยอ่าน ไม่ใช่การแสดงรายละเอียดทุก step ของเส้นทางต่อเนื่อง และคำตอบ exact ของแบบจำลองไม่ได้รับประกันความถูกต้องของแบบจำลองต่อโลกจริง

Euler อาจให้ราคาศูนย์หรือติดลบถ้าตัวคูณ 1 + μΔt + σ√Δt Z ไม่เป็นบวก แม้คำตอบ GBM จริงจะเป็นบวก ตัวอย่างสมมติ μ = 0.1, σ = 1, Δt = 1 และ Z = −2 ให้ตัวคูณ Euler = −0.9 แต่ตัวคูณ exact = e^(−2.4) > 0 การตัดค่าติดลบให้เป็นศูนย์จะเปลี่ยนวิธีและการแจกแจง จึงไม่ควรทำโดยไม่อธิบาย

เมื่อย่อ Δt ภายใต้เงื่อนไขที่เหมาะสม Euler–Maruyama ลู่เข้าในความหมายทางความน่าจะเป็นที่กำหนด แต่ error ของเส้นทางเดียวไม่จำเป็นต้องลดทุกครั้ง การตรวจความแม่นควรใช้หลายเส้นทางร่วมกับช็อกที่จับคู่กัน

แหล่งความคลาดเคลื่อน เกิดจากอะไร สิ่งที่ช่วยลด
การแบ่งเวลา แทน SDE ด้วย step ที่ยังมีขนาดจำกัด ย่อ Δt หรือเลือกวิธีเหมาะสม
Monte Carlo sampling error ประมาณค่าคาดหมายด้วยเส้นทางจำนวนจำกัด เพิ่มจำนวนเส้นทางและวัด standard error
Model error สมมติฐานไม่ตรงกับสิ่งที่ต้องการอธิบาย ตรวจแบบจำลองกับข้อมูลและโจทย์จริง

9. ตัวสุ่มที่มีค่าเฉลี่ยและ SD ถูก ยังอาจมีหางผิด

สองเอกสารพูดถึงการสร้างตัวสุ่ม Normal สำหรับจำลอง SDE แนวคิด inverse transform คือ ถ้า U กระจายสม่ำเสมอบน (0,1) แล้ว Z = Φ⁻¹(U) จะมีการแจกแจง Standard Normal โดย Φ เป็น CDF ของ Normal

อีกวิธีในเอกสารคือผลรวม uniform 12 ตัวลบ 6 ซึ่งมีค่าเฉลี่ย 0 และความแปรปรวน 1 จริง แต่มีค่าจำกัดอยู่ในช่วง [−6,6] จึงเป็นเพียง Normal approximation และไม่สามารถแทนหาง Normal ได้ทุกระดับ

ห้องทดลองในบทใช้ตัวสร้างเลขสุ่มที่ระบุ seed และแปลงเป็น Normal เพื่อให้ทดลองซ้ำได้ ส่วน Notebook ระบุวิธีและ seed ของตนเอง การใช้ seed เดียวกันต่างภาษาไม่ได้รับประกันว่าจะได้เลขชุดเดียวกัน ถ้าอัลกอริทึมสุ่มต่างกัน

การตรวจตัวสุ่มต้องดูทั้งชนิดการแจกแจงและหาง ควบคู่กับค่าเฉลี่ยและ SD

10. จากช็อกเดียว สู่ช็อกสองตัวที่สัมพันธ์กัน

ถ้าจำลองสองสินทรัพย์โดยให้ช็อกเป็นอิสระเสมอ เรากำลังกำหนดโครงสร้างความเสี่ยงบางอย่างไว้แล้ว การสร้าง correlated increments เริ่มจากช็อก Standard Normal ที่มี correlation ρ โดยใช้ Z₁,Z₂ ซึ่งเป็น อิสระต่อกัน และต่างก็เป็น Standard Normal:

ϕ1=Z1,ϕ2=ρZ1+1−ρ2 Z2,−1≤ρ≤1.\phi_1=Z_1,\qquad \phi_2=\rho Z_1+\sqrt{1-\rho^2}\,Z_2, \qquad -1\le\rho\le1.

จะได้ E[φ₁] = E[φ₂] = 0, Var(φ₁) = Var(φ₂) = 1 และ

Cov⁡(ϕ1,ϕ2)=E[ϕ1ϕ2]=ρ.\operatorname{Cov}(\phi_1,\phi_2) =\mathbb E[\phi_1\phi_2]=\rho.

ตัวอย่าง ρ = −0.6, Z₁ = 1, Z₂ = −0.5 ให้ φ₁ = 1 และ φ₂ = −0.6×1 + 0.8×(−0.5) = −1 นี่เป็นเพียงช็อกหนึ่งคู่ ค่า correlation เป็นคุณสมบัติของการแจกแจงร่วม ต้องดูหลายคู่จึงประมาณจากตัวอย่างได้

เมื่อ ρ = 1 จะมี φ₂ = φ₁ และเมื่อ ρ = −1 จะมี φ₂ = −φ₁ ทุกคู่ ส่วน ρ = 0 ให้ช็อกสองตัวอิสระในโครงสร้าง Gaussian นี้ การเลื่อน ρ ใช้ Z₁,Z₂ ชุดเดิมเพื่อให้เห็นผลของพารามิเตอร์

ในการจำลองตามเวลา ให้สุ่มคู่ Z₁,Z₂ ใหม่อย่างอิสระจากคู่ในช่วงก่อนหน้า เมื่อแปลงเป็น step Brownian ด้วย ΔW₁ = √Δt φ₁ และ ΔW₂ = √Δt φ₂ จะมี Cov(ΔW₁,ΔW₂) = ρΔt หรือเขียนย่อว่า

dW1 dW2=ρ dt.dW_1\,dW_2=\rho\,dt.

ρ นี้เป็น correlation ของ ช็อกในช่วงเวลาเดียวกัน ไม่ใช่ correlation ของระดับราคาที่วัดตามเวลาบนเส้นทางหนึ่งโดยอัตโนมัติ และ sample correlation จากจำนวนคู่จำกัดไม่จำเป็นต้องเท่ากับค่าเป้าหมายพอดี

เมื่อบวกตัวแปร Normal ต้องดู covariance ด้วย

ผลรวมเชิงเส้นของตัวแปรที่เป็น jointly Gaussian ยังเป็น Normal หากตัวแปรอิสระและแต่ละตัวเป็น Normal จะเข้าเงื่อนไขนี้ ส่วนการรู้ว่าแต่ละตัวมี marginal Normal อย่างเดียวไม่พอ

Var⁡ ⁣(∑iwiXi)=∑iwi2Var⁡(Xi)+2∑i<jwiwjCov⁡(Xi,Xj).\operatorname{Var}\!\left(\sum_i w_iX_i\right) =\sum_i w_i^2\operatorname{Var}(X_i) +2\sum_{i<j}w_iw_j\operatorname{Cov}(X_i,X_j).

จะตัดพจน์ covariance ทิ้งได้เมื่อเป็นศูนย์ ในสูตรสร้าง φ₁,φ₂ เราเริ่มจาก Z₁,Z₂ อิสระ แล้วจึงสร้างความสัมพันธ์ขึ้นอย่างควบคุมได้ ข้อสมมตินี้เติมให้ชัดจากสูตรสรุปในเอกสาร

ลองเชื่อมเครื่องมือเข้าด้วยกัน

  1. ถ้า Δt = 0.01 และ Z = 2 step ΔW และกำลังสองของ step เป็นเท่าไร กฎ (dW)² = dt บอกว่าค่านี้ต้องเท่ากับ 0.01 หรือไม่
  2. ถ้า dY = 0.3dt + 0.4dW และ F(Y) = Y² จงหา dF
  3. GBM มี μ = 0.10 และ σ = 0.20 drift ของ log S เป็นเท่าไร และเท่ากับ drift ของผลตอบแทน dS/S หรือไม่
  4. ใน OU ถ้า r สูงกว่า θ จะสรุปได้หรือไม่ว่า step ถัดไป r ต้องลด
  5. ถ้าจะเปรียบเทียบ Euler กับ exact ควรใช้ช็อกชุดเดียวกันหรือคนละชุด และเพราะเหตุใด
  6. ถ้าเพิ่มจำนวนเส้นทาง Monte Carlo อย่างเดียว ความผิดพลาดจาก Δt ที่หยาบและสมมติฐานโมเดลจะหายหรือไม่
เปิดเฉลยพร้อมเหตุผล
  1. ΔW = √0.01 × 2 = 0.2 และ (ΔW)² = 0.04 กฎเชิง differential ไม่ใช่การบังคับให้กำลังสองของ step ที่มีขนาดจำกัดเท่ากับ Δt
  2. F_Y = 2Y, F_YY = 2 จึงได้ dF = (0.6Y + 0.16)dt + 0.8Y dW พจน์ 0.16 มาจาก Itô correction
  3. 0.08 ต่อปี เทียบกับ drift 0.10 ต่อปีของผลตอบแทนเชิง differential dS/S จึงไม่เท่ากัน
  4. สรุปไม่ได้ drift เป็นลบ แต่ช็อกอาจเป็นบวกและมีขนาดใหญ่กว่าแรงดึงกลับ
  5. ชุดเดียวกัน เพื่อแยกความต่างของวิธีคำนวณออกจากความต่างของเส้นทางสุ่ม
  6. ไม่หาย การเพิ่มเส้นทางลด sampling error แต่ยังต้องตรวจ discretization error และ model error แยกกัน

ดาวน์โหลด Notebook ของบท Applied Stochastic Calculus เพื่อรันเอกลักษณ์ Itô ตรวจ moments ของ GBM/OU และทดลองการแบ่งเวลากับช็อกที่สัมพันธ์กัน

เปิด JavaScript เพื่อใช้ห้องทดลอง เนื้อหา สมการ ตาราง และเฉลยยังอ่านได้โดยไม่เปิด JavaScript

เครื่องมือแต่ละชิ้นมีหน้าที่ต่อกัน: quadratic variation อธิบายพจน์ที่กฎธรรมดาทิ้งไม่ได้, Itô’s lemma แปลงฟังก์ชันของกระบวนการ, Kolmogorov อธิบายการแจกแจง และการจำลองช่วยคำนวณเมื่อไม่สะดวกใช้คำตอบวิเคราะห์

อ้างอิงและขอบเขตการเรียบเรียง

  • เนื้อหาเรียบเรียงภาษาไทยใหม่จากสองไฟล์ เปลี่ยน X ของ Brownian เป็น W และเปลี่ยน γ,r̄ ของ Vasicek เป็น κ,θ เพื่อรักษาความต่อเนื่องกับบทก่อน ไม่เผยแพร่ PDF หรือภาพสไลด์ต้นฉบับ
  • ส่วนที่ขยายให้ชัด ได้แก่ การแยก coin increments จาก Gaussian increments, ความหมายของ (dW)² = dt, เงื่อนไข Itô integral และความเรียบของ F, ค่าเฉลี่ยกับมัธยฐานของ GBM, moments ระหว่างทางของ OU, ความคลาดเคลื่อนสามประเภท และข้อสมมติ jointly Gaussian / covariance ตัวอย่างและห้องทดลองคำนวณขึ้นใหม่ ยกเว้นตัวอย่าง Euler สอง step ที่ระบุว่ามาจากเอกสาร
  • Miranda Holmes-Cerfon, Applied Stochastic Analysis (Spring 2022): Lecture 7 — Stochastic Integration, Lecture 8 — Stochastic Differential Equations, Lecture 10 — Forward and Backward Equations ใช้ตรวจนิยาม เงื่อนไข และการเชื่อม SDE กับสมการความหนาแน่น
ภาพการ์ตูน Nuthdanai Wangpratham

ผู้เขียน · QuantCorner

Nuthdanai Wangpratham

แบ่งปันความรู้ด้าน Quantitative Finance และ AI เพื่อให้การเรียนรู้เรื่องการลงทุนเริ่มต้นได้จากตัวอย่างที่ลองทำตามได้

อ่านบทนี้แล้วอยากเรียนรู้ต่อ? ติดตามผู้เขียนและ QuantCorner ได้ที่นี่

ค้นหาในสมุดบันทึก