AI Cafe¹
© Bản quyền thuộc về Quantopian Inc.
© Bản quyền sửa đổi thuộc về AI Cafe¹
Giấy phép Creative Commons Attribution 4.0.
Khuyến cáo

Open in Colab¶

Ước lượng Maximum Likelihood (MLEs)¶

By Delaney Granizo-Mackenzie and Andrei Kirilenko developed as part of the Masters of Finance curriculum at MIT Sloan.

Edited by AI Cafe¹

Mục tiêu của phương pháp này là tìm ra giá trị của các tham số (của một mô hình thống kê) mà ở đó xác suất xuất hiện đồng thời các giá trị mẫu (gọi là likelihood) là lớn nhất. Trong bài này, chúng ta sẽ thực hiện:

  1. Tính toán MLE cho phân phối chuẩn.
  2. Tính toán MLE cho phân phối mũ.
  3. Gắn kết một phân phối chuẩn vào lợi suất tài sản bằng cách sử dụng MLE.
In [1]:
import math
import matplotlib.pyplot as plt
import numpy as np
import scipy
import scipy.stats

Phân phối chuẩn¶

Chúng ta sẽ bắt đầu bằng cách lấy mẫu dữ liệu từ một phân phối chuẩn.

In [2]:
TRUE_MEAN = 40
TRUE_STD = 10
X = np.random.normal(TRUE_MEAN, TRUE_STD, 1000)

Bây giờ chúng ta sẽ định nghĩa các hàm mà, dựa trên dữ liệu của chúng ta, sẽ tính toán MLE cho các tham số $\mu$ và $\sigma$ của phân phối chuẩn.

Hãy nhớ lại công thức

$$\hat\mu = \frac{1}{T}\sum_{t=1}^{T} x_t$$

$$\hat\sigma = \sqrt{\frac{1}{T}\sum_{t=1}^{T}{(x_t - \hat\mu)^2}}$$

In [3]:
def normal_mu_MLE(X):
    # Get the number of observations
    T = len(X)
    # Sum the observations
    s = sum(X)
    return 1.0/T * s

def normal_sigma_MLE(X):
    T = len(X)
    # Get the mu MLE
    mu = normal_mu_MLE(X)
    # Sum the square of the differences
    s = sum( np.power((X - mu), 2) )
    # Compute sigma^2
    sigma_squared = 1.0/T * s
    return math.sqrt(sigma_squared)

Bây giờ hãy thử các hàm của chúng ta trên dữ liệu mẫu và so sánh chúng với hàm np.mean và np.std có sẵn.

In [4]:
print("Mean Estimation")
print(normal_mu_MLE(X))
print(np.mean(X))
print("Standard Deviation Estimation")
print(normal_sigma_MLE(X))
print(np.std(X))

Mean Estimation
39.97339283238262
39.97339283238264
Standard Deviation Estimation
10.110058391084257
10.110058391084253

Bây giờ hãy ước lượng cả hai tham số cùng một lúc bằng cách sử dụng hàm fit() tích hợp sẵn của scipy.

In [5]:
mu, std = scipy.stats.norm.fit(X)
print("mu estimate:",  str(mu))
print("std estimate:", str(std))

mu estimate: 39.97339283238264
std estimate: 10.110058391084253

Bây giờ hãy vẽ phân phối PDF cùng với dữ liệu để xem nó phù hợp như thế nào. Chúng ta có thể làm điều đó bằng cách truy cập vào pdf được cung cấp trong scipy.stats.norm.pdf.

In [6]:
pdf = scipy.stats.norm.pdf
# We would like to plot our data along an x-axis ranging from 0-80 with 80 intervals
# (increments of 1)
x = np.linspace(0, 80, 80)
plt.hist(X, bins=x, density='true')
plt.plot(pdf(x, loc=mu, scale=std))
plt.xlabel('Value')
plt.ylabel('Observed Frequency')
plt.legend(['Fitted Distribution PDF', 'Observed Data', ]);
No description has been provided for this image

Phân phối mũ¶

Hãy làm điều tương tự, nhưng cho phân phối mũ. Chúng ta sẽ bắt đầu bằng cách lấy mẫu một số dữ liệu.

In [7]:
TRUE_LAMBDA = 5
X = np.random.exponential(TRUE_LAMBDA, 1000)

numpy định nghĩa phân phối mũ là $$\frac{1}{\lambda}e^{-\frac{x}{\lambda}}$$

Vì vậy, chúng ta cần sử dụng giá trị nghịch đảo của MLE. Tại đó, nó được thể hiện là

$$\hat\lambda = \frac{T}{\sum_{t=1}^{T} x_t}$$

Ở đây chỉ là nghịch đảo, vì vậy

$$\hat\lambda = \frac{\sum_{t=1}^{T} x_t}{T}$$

In [8]:
def exp_lamda_MLE(X):
    T = len(X)
    s = sum(X)
    return s/T
In [9]:
print("lambda estimate:", str(exp_lamda_MLE(X)))

lambda estimate: 5.19533144751805
In [10]:
# The scipy version of the exponential distribution has a location parameter
# that can skew the distribution. We ignore this by fixing the location
# parameter to 0 with floc=0
_, l = scipy.stats.expon.fit(X, floc=0)
In [11]:
pdf = scipy.stats.expon.pdf
x = range(0, 80)
plt.hist(X, bins=x, density='true')
plt.plot(pdf(x, scale=l))
plt.xlabel('Value')
plt.ylabel('Observed Frequency')
plt.legend(['Fitted Distribution PDF', 'Observed Data', ]);
No description has been provided for this image

MLE cho Lợi suất (returns) của Tài sản¶

Bây giờ chúng ta sẽ lấy một số lợi suất thực tế và cố gắng ép chúng vào một phân phối chuẩn bằng cách sử dụng MLE.

In [12]:
# Cài thư viện lấy dữ liệu miễn phí
!curl -fsSLO https://raw.githubusercontent.com/algo-stocks/data/master/data.py
In [13]:
from data import get_prices

start = '2022-01-01'
end = '2025-01-01'

closes = get_prices('FPT', start_date=start, end_date=end)
returns = closes.pct_change()[1:]

Hãy sử dụng hàm fit của scipy để lấy các giá trị MLEs cho $\mu$ và $\sigma$.

In [15]:
mu, std = scipy.stats.norm.fit(returns)
pdf = scipy.stats.norm.pdf
x = np.linspace(-.1,.1, num=100)
h = plt.hist(returns, bins=x, density='true')
l = plt.plot(x, pdf(x, loc=mu, scale=std))
No description has been provided for this image

Tất nhiên, việc ước lượng này sẽ không có ý nghĩa nếu chúng ta chưa kiểm tra xem chúng có tuân theo phân phối chuẩn hay không. Chúng ta có thể kiểm tra điều này bằng cách sử dụng kiểm định Jarque-Bera. Kiểm định Jarque-Bera sẽ bác bỏ giả thuyết về phân phối chuẩn nếu giá trị p-value nhỏ hơn một ngưỡng.

In [16]:
from statsmodels.stats.stattools import jarque_bera
jarque_bera(returns)
Out[16]:
(array([306.70352218]),
 array([2.51290672e-67]),
 array([0.10809672]),
 array([6.13164315]))

Giá trị p-value quá nhỏ (2.51290672e-67) nên có thể kết luận rằng tập mẫu không tuân theo quy luật phân phối chuẩn. Hay xem kết quả bài kiểm định Jarque-Bera với dữ liệu được lấy mẫu từ phân phối chuẩn sau:

In [17]:
jarque_bera(np.random.normal(0, 1, 100))
Out[17]:
(np.float64(0.2612005890317266),
 np.float64(0.8775684732297053),
 np.float64(0.049438133164470335),
 np.float64(3.2300252710561663))

Lưu ý: Tài liệu này mang tính chất chia sẻ kiến thức khoa học dữ liệu. Các mã chứng khoán xuất hiện trong tài liệu chỉ mang tính chất minh hoạ, không khuyến nghị đầu tư. Chúng tôi không chịu trách nhiệm trước mọi rủi ro nào do sử dụng tài liệu này.