Thực hành bài giữa kỳ
Nguyễn Hiếu Nghi
2026-06-09
- đọc dữ liệu từ link gg sheet
library(googlesheets4)
gs4_deauth()
df <- read_sheet("https://docs.google.com/spreadsheets/d/10nZIVWn3i5djs0tRZbI--8nntNbfbGoNvYKX8yI7Eeo/edit?gid=1870826978#gid=1870826978", sheet = "hose3")- Lấy dữ liệu mã cổ phiếu
df# A tibble: 3,335 × 241
Date AAA AAM ABT ACB ACC ACL ADP AGR ANV
<dttm> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 2012-04-11 00:00:00 4582 9448 10913 3583 2760 5507 422 7418 1115
2 2012-04-12 00:00:00 4560 9557 10972 3681 2851 5446 461 7418 1115
3 2012-04-13 00:00:00 4476 9557 10972 3681 2851 5446 505 7339 1094
4 2012-04-16 00:00:00 4624 9557 11300 3681 2967 5538 477 7655 1115
5 2012-04-17 00:00:00 4624 9703 11420 3667 2993 5660 521 7971 1126
6 2012-04-18 00:00:00 4560 9849 11330 3611 3110 5568 532 8208 1115
7 2012-04-19 00:00:00 4348 9849 11390 3555 3226 5844 581 7813 1061
8 2012-04-20 00:00:00 4370 9885 11420 3569 3071 6119 581 7497 1105
9 2012-04-23 00:00:00 4370 9885 11420 3611 3006 6241 581 7339 1159
10 2012-04-24 00:00:00 4433 9995 11539 3611 3110 6241 526 7418 1207
# ℹ 3,325 more rows
# ℹ 231 more variables: APG <dbl>, ASM <dbl>, ASP <dbl>, BCE <dbl>, BIC <dbl>,
# BMC <dbl>, BMI <dbl>, BMP <dbl>, BRC <dbl>, BSI <dbl>, BTP <dbl>,
# BTT <dbl>, BVH <dbl>, C32 <dbl>, C47 <dbl>, CCI <dbl>, CCL <dbl>,
# CDC <dbl>, CHP <dbl>, CIG <dbl>, CII <dbl>, CLC <dbl>, CLW <dbl>,
# CMG <dbl>, CMV <dbl>, CMX <dbl>, CNG <dbl>, COM <dbl>, CSM <dbl>,
# CTD <dbl>, CTG <dbl>, CTI <dbl>, CTS <dbl>, CVT <dbl>, D2D <dbl>, …tmt<-df[, c("Date", "TMT")]
tmt# A tibble: 3,335 × 2
Date TMT
<dttm> <dbl>
1 2012-04-11 00:00:00 3089
2 2012-04-12 00:00:00 3089
3 2012-04-13 00:00:00 3089
4 2012-04-16 00:00:00 3201
5 2012-04-17 00:00:00 3201
6 2012-04-18 00:00:00 3201
7 2012-04-19 00:00:00 3089
8 2012-04-20 00:00:00 3089
9 2012-04-23 00:00:00 3201
10 2012-04-24 00:00:00 3089
# ℹ 3,325 more rows- Kiểm tra dữ liệu
str(tmt)tibble [3,335 × 2] (S3: tbl_df/tbl/data.frame)
$ Date: POSIXct[1:3335], format: "2012-04-11" "2012-04-12" ...
$ TMT : num [1:3335] 3089 3089 3089 3201 3201 ...head(tmt)# A tibble: 6 × 2
Date TMT
<dttm> <dbl>
1 2012-04-11 00:00:00 3089
2 2012-04-12 00:00:00 3089
3 2012-04-13 00:00:00 3089
4 2012-04-16 00:00:00 3201
5 2012-04-17 00:00:00 3201
6 2012-04-18 00:00:00 3201tail(tmt)# A tibble: 6 × 2
Date TMT
<dttm> <dbl>
1 2026-05-15 00:00:00 11900
2 2026-05-18 00:00:00 11800
3 2026-05-19 00:00:00 11800
4 2026-05-20 00:00:00 11600
5 2026-05-21 00:00:00 11800
6 2026-05-22 00:00:00 11950summary(tmt) Date TMT
Min. :2012-04-11 00:00:00 Min. : 2415
1st Qu.:2015-09-16 12:00:00 1st Qu.: 5365
Median :2019-04-23 00:00:00 Median :10106
Mean :2019-05-03 01:14:16 Mean :12105
3rd Qu.:2022-12-27 12:00:00 3rd Qu.:15550
Max. :2026-05-22 00:00:00 Max. :51306 4.1 Thống kê mô tả dữ liệu giá chứng khoán
library(fBasics)
basicStats(tmt$TMT) #cách 1 X..tmt.TMT
nobs 3.335000e+03
NAs 0.000000e+00
Minimum 2.415000e+03
Maximum 5.130600e+04
1. Quartile 5.365000e+03
3. Quartile 1.555000e+04
Mean 1.210531e+04
Median 1.010600e+04
Sum 4.037120e+07
SE Mean 1.504014e+02
LCL Mean 1.181042e+04
UCL Mean 1.240019e+04
Variance 7.543968e+07
Stdev 8.685602e+03
Skewness 1.443408e+00
Kurtosis 1.971220e+00summary(tmt$TMT) #cách 2 Min. 1st Qu. Median Mean 3rd Qu. Max.
2415 5365 10106 12105 15550 51306 4.2 Thống kê mô tả dữ liệu return của giá chứng khoán
return_tmt <- diff(log(tmt$TMT))
head(return_tmt)[1] 0.00000000 0.00000000 0.03561585 0.00000000 0.00000000 -0.03561585tail(return_tmt)[1] -0.033060862 -0.008438869 0.000000000 -0.017094433 0.017094433
[6] 0.012631747basicStats(return_tmt) return_tmt
nobs 3334.000000
NAs 0.000000
Minimum -0.125174
Maximum 0.329148
1. Quartile -0.018422
3. Quartile 0.017912
Mean 0.000406
Median 0.000000
Sum 1.352884
SE Mean 0.000598
LCL Mean -0.000767
UCL Mean 0.001579
Variance 0.001194
Stdev 0.034548
Skewness 0.452047
Kurtosis 3.1545385.1 Vẽ đồ thị và hàm mật độ giá chứng khoán
ts.plot(tmt$TMT) 
plot(tmt, type = 'l') 
hist(tmt$TMT)
d1 = density(tmt$TMT)
plot(d1, type = "l")
5.2 Vẽ đồ thị và hàm mật độ của return giá chứng khoán
plot(return_tmt, type = 'l')
hist(return_tmt)
d2 = density(return_tmt)
plot(d2, type = "l")
- Kiểm tra tính dừng, kiểm tra phân phối chuẩn
library(tseries)
adf.test(tmt$TMT)
Augmented Dickey-Fuller Test
data: tmt$TMT
Dickey-Fuller = -1.8847, Lag order = 14, p-value = 0.6272
alternative hypothesis: stationaryadf.test(return_tmt)
Augmented Dickey-Fuller Test
data: return_tmt
Dickey-Fuller = -15.213, Lag order = 14, p-value = 0.01
alternative hypothesis: stationarykpss.test(return_tmt) #Ho: stationary
KPSS Test for Level Stationarity
data: return_tmt
KPSS Level = 0.1457, Truncation lag parameter = 9, p-value = 0.1normalTest(tmt$TMT, method = c("jb")) #H0: normal distribution
Title:
Jarque-Bera Normality Test
Test Results:
STATISTIC:
X-squared: 1700.667
P VALUE:
Asymptotic p Value: < 2.2e-16 normalTest(return_tmt,method = c("jb")) #H0: normal distribution
Title:
Jarque-Bera Normality Test
Test Results:
STATISTIC:
X-squared: 1499.2644
P VALUE:
Asymptotic p Value: < 2.2e-16 - Kiểm tra tự tương quan bằng đồ thị và bằng kiểm định
acf(return_tmt)
pacf(return_tmt)
Box.test(return_tmt, lag = 12) #h0: no autocorrelation
Box-Pierce test
data: return_tmt
X-squared = 24.514, df = 12, p-value = 0.0173Box.test(return_tmt, lag = 12, type = c("Ljung-Box")) #ho: no autocorrelation
Box-Ljung test
data: return_tmt
X-squared = 24.548, df = 12, p-value = 0.01711- Ước lượng mô hình ARIMA #8.1 Xây dựng mô hình AR
pacf(return_tmt)
m1 <- ar(return_tmt, method = c("mle"))
m1
Call:
ar(x = return_tmt, method = c("mle"))
Coefficients:
1
0.0702
Order selected 1 sigma^2 estimated as 0.001187m1$aic 0 1 2 3 4 5 6
14.4551208 0.0000000 1.9778631 2.1077406 0.9258144 2.7564381 4.2214368
7 8 9 10 11 12
6.1869164 8.1734228 9.2517810 9.8430512 10.6556452 12.6239812 #Ước lượng AR(p)
m21 <- arima(return_tmt, order = c(1,0,0))
m21
Call:
arima(x = return_tmt, order = c(1, 0, 0))
Coefficients:
ar1 intercept
0.0702 4e-04
s.e. 0.0173 6e-04
sigma^2 estimated as 0.001187: log likelihood = 6498.29, aic = -12990.59#Ước lượng AR(p)
m22<- arima(return_tmt, order = c(13,0,0))
m22
Call:
arima(x = return_tmt, order = c(13, 0, 0))
Coefficients:
ar1 ar2 ar3 ar4 ar5 ar6 ar7 ar8
0.0717 0.0008 0.0254 -0.0310 0.0089 -0.0137 -0.0025 -0.0032
s.e. 0.0173 0.0173 0.0173 0.0173 0.0173 0.0173 0.0173 0.0173
ar9 ar10 ar11 ar12 ar13 intercept
0.0167 -0.0209 0.0186 0.0058 -0.0382 4e-04
s.e. 0.0174 0.0173 0.0173 0.0173 0.0173 6e-04
sigma^2 estimated as 0.001182: log likelihood = 6505.42, aic = -12980.84AIC(m21,m22) df AIC
m21 3 -12990.59
m22 15 -12980.84# hàm tính p-value
coeftest <- function(model){
coef <- model$coef
se <- sqrt(diag(model$var.coef))
t <- coef / se
p <- 2 * (1 - pnorm(abs(t)))
data.frame(coef, se, t, p)
}
coeftest(m21) coef se t p
ar1 0.0701578105 0.0172726086 4.0617959 4.869665e-05
intercept 0.0004058068 0.0006420015 0.6320964 5.273239e-01coeftest(m22) coef se t p
ar1 0.0717338278 0.0173024979 4.14586543 3.385326e-05
ar2 0.0007899830 0.0173462914 0.04554189 9.636754e-01
ar3 0.0253758286 0.0173435638 1.46312655 1.434328e-01
ar4 -0.0309695476 0.0173435961 -1.78564742 7.415636e-02
ar5 0.0088759687 0.0173467907 0.51167785 6.088765e-01
ar6 -0.0137141356 0.0173498967 -0.79044480 4.292680e-01
ar7 -0.0024606920 0.0173498142 -0.14182815 8.872158e-01
ar8 -0.0032494737 0.0173490257 -0.18730007 8.514254e-01
ar9 0.0167493861 0.0173504017 0.96536014 3.343645e-01
ar10 -0.0208530400 0.0173447490 -1.20226819 2.292596e-01
ar11 0.0186430433 0.0173405387 1.07511327 2.823240e-01
ar12 0.0058397422 0.0173420068 0.33673970 7.363131e-01
ar13 -0.0382357766 0.0172956364 -2.21071812 2.705536e-02
intercept 0.0004054601 0.0006196464 0.65434115 5.128920e-01# Kiểm định phần dư sau khi đã chọn bậc cho mô hình
plot(m21$residuals, type = 'l')
qqplot(m21$residuals, return_tmt)
Box.test(m21$residuals, lag = 12)
Box-Pierce test
data: m21$residuals
X-squared = 9.307, df = 12, p-value = 0.67659.1 Kiểm định hiệu ứng ARCH (ARCH-LM test)
library(fGarch)
m4 <- garchFit(~ arma (1,0) + garch(1,0), data = return_tmt, trace = FALSE)
summary(m4)
Title:
GARCH Modelling
Call:
garchFit(formula = ~arma(1, 0) + garch(1, 0), data = return_tmt,
trace = FALSE)
Mean and Variance Equation:
data ~ arma(1, 0) + garch(1, 0)
<environment: 0x0000028c8c018550>
[data = return_tmt]
Conditional Distribution:
norm
Coefficient(s):
mu ar1 omega alpha1
-0.00123133 0.00302328 0.00083031 0.32178008
Std. Errors:
based on Hessian
Error Analysis:
Estimate Std. Error t value Pr(>|t|)
mu -1.231e-03 5.431e-04 -2.267 0.0234 *
ar1 3.023e-03 2.056e-02 0.147 0.8831
omega 8.303e-04 2.816e-05 29.482 <2e-16 ***
alpha1 3.218e-01 2.914e-02 11.042 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Log Likelihood:
6604.035 normalized: 1.980814
Description:
Tue Aug 25 01:19:58 2026 by user: DELL 5110
Standardised Residuals Tests:
Statistic p-Value
Jarque-Bera Test R Chi^2 2874.259158 0.000000000
Shapiro-Wilk Test R W 0.966452 0.000000000
Ljung-Box Test R Q(10) 10.474208 0.399919045
Ljung-Box Test R Q(15) 16.996324 0.319084221
Ljung-Box Test R Q(20) 23.602323 0.260195071
Ljung-Box Test R^2 Q(10) 25.792783 0.004029108
Ljung-Box Test R^2 Q(15) 30.899514 0.009060859
Ljung-Box Test R^2 Q(20) 41.960735 0.002798349
LM Arch Test R TR^2 25.590811 0.012258413
Information Criterion Statistics:
AIC BIC SIC HQIC
-3.959229 -3.951896 -3.959232 -3.956606 9.2 Ước lượng ARIMA+GARCH
m5 <-garchFit(~ arma(1,0) + garch(1,1), data = return_tmt, trace = FALSE)
summary(m5)
Title:
GARCH Modelling
Call:
garchFit(formula = ~arma(1, 0) + garch(1, 1), data = return_tmt,
trace = FALSE)
Mean and Variance Equation:
data ~ arma(1, 0) + garch(1, 1)
<environment: 0x0000028c8e1358c8>
[data = return_tmt]
Conditional Distribution:
norm
Coefficient(s):
mu ar1 omega alpha1 beta1
-0.00058777 0.02208209 0.00018643 0.17589328 0.67267809
Std. Errors:
based on Hessian
Error Analysis:
Estimate Std. Error t value Pr(>|t|)
mu -5.878e-04 5.362e-04 -1.096 0.273
ar1 2.208e-02 1.943e-02 1.136 0.256
omega 1.864e-04 3.546e-05 5.257 1.46e-07 ***
alpha1 1.759e-01 2.340e-02 7.516 5.62e-14 ***
beta1 6.727e-01 4.689e-02 14.345 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Log Likelihood:
6636.742 normalized: 1.990624
Description:
Tue Aug 25 01:19:59 2026 by user: DELL 5110
Standardised Residuals Tests:
Statistic p-Value
Jarque-Bera Test R Chi^2 2349.2907140 0.0000000
Shapiro-Wilk Test R W 0.9696294 0.0000000
Ljung-Box Test R Q(10) 7.4994685 0.6775991
Ljung-Box Test R Q(15) 14.9355910 0.4560666
Ljung-Box Test R Q(20) 22.7818593 0.2996322
Ljung-Box Test R^2 Q(10) 4.4677833 0.9237876
Ljung-Box Test R^2 Q(15) 8.0549284 0.9215617
Ljung-Box Test R^2 Q(20) 15.9890032 0.7173062
LM Arch Test R TR^2 7.7680079 0.8029869
Information Criterion Statistics:
AIC BIC SIC HQIC
-3.978249 -3.969083 -3.978254 -3.974970 predict(m5,5) meanForecast meanError standardDeviation
1 -0.0003088308 0.02662113 0.02662113
2 -0.0005945858 0.02807384 0.02806769
3 -0.0006008959 0.02924569 0.02923912
4 -0.0006010352 0.03020445 0.03019754
5 -0.0006010383 0.03099477 0.03098759- Tính và dự báo Var
#chạy hàm loss để tính Var
nreturn_tmt <- -return_tmt
#10.1 Ước tính Var theo RISK METRICS
#cách 1: tính thủ công
library(rugarch)
spec1 <- ugarchspec(mean.model = list(armaOrder=c(1,0)),
variance.model = list(model="iGARCH", garchOrder=c(1,1)))
m6 <- ugarchfit(spec=spec1,data = nreturn_tmt)
Box.test(residuals(m6, standardize=TRUE), lag=10, type="Ljung")
Box-Ljung test
data: residuals(m6, standardize = TRUE)
X-squared = 7.7834, df = 10, p-value = 0.65ugarchforecast(m6, n.ahead = 1)
*------------------------------------*
* GARCH Model Forecast *
*------------------------------------*
Model: iGARCH
Horizon: 1
Roll Steps: 0
Out of Sample: 0
0-roll forecast [T0=1979-02-17]:
Series Sigma
T+1 -8.529e-05 0.02521# Tính VaR (Phương pháp của JP Morgan xem mean = 0):
VaR1 <- 0 + qnorm(0.95) * 0.02521
# Cách 2: tính từ RMeasure.R
setwd("C:/Users/DELL 5110/Downloads")
source("RMeasure.R")
VaR2 <- RMeasure(0, 0.02521)
Risk Measures for selected probabilities:
prob VaR ES
[1,] 0.950 0.04146676 0.05200099
[2,] 0.990 0.05864723 0.06719005
[3,] 0.999 0.07790476 0.08488434# 10.2. Ước tính VaR theo ECONOMICTRICS
m7 <- garchFit(~ arma(1,0) + garch(1,1), data = nreturn_tmt, trace=F)
summary(m7)
Title:
GARCH Modelling
Call:
garchFit(formula = ~arma(1, 0) + garch(1, 1), data = nreturn_tmt,
trace = F)
Mean and Variance Equation:
data ~ arma(1, 0) + garch(1, 1)
<environment: 0x0000028c88aa1708>
[data = nreturn_tmt]
Conditional Distribution:
norm
Coefficient(s):
mu ar1 omega alpha1 beta1
0.00058776 0.02208206 0.00018643 0.17589327 0.67267811
Std. Errors:
based on Hessian
Error Analysis:
Estimate Std. Error t value Pr(>|t|)
mu 5.878e-04 5.362e-04 1.096 0.273
ar1 2.208e-02 1.943e-02 1.136 0.256
omega 1.864e-04 3.546e-05 5.257 1.46e-07 ***
alpha1 1.759e-01 2.340e-02 7.516 5.62e-14 ***
beta1 6.727e-01 4.689e-02 14.345 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Log Likelihood:
6636.742 normalized: 1.990624
Description:
Tue Aug 25 01:20:07 2026 by user: DELL 5110
Standardised Residuals Tests:
Statistic p-Value
Jarque-Bera Test R Chi^2 2349.2908510 0.0000000
Shapiro-Wilk Test R W 0.9696294 0.0000000
Ljung-Box Test R Q(10) 7.4994726 0.6775987
Ljung-Box Test R Q(15) 14.9355951 0.4560663
Ljung-Box Test R Q(20) 22.7818639 0.2996320
Ljung-Box Test R^2 Q(10) 4.4677827 0.9237876
Ljung-Box Test R^2 Q(15) 8.0549275 0.9215617
Ljung-Box Test R^2 Q(20) 15.9890030 0.7173063
LM Arch Test R TR^2 7.7680070 0.8029869
Information Criterion Statistics:
AIC BIC SIC HQIC
-3.978249 -3.969083 -3.978254 -3.974970 predict(m7,1) meanForecast meanError standardDeviation
1 0.00030883 0.02662113 0.02662113# Tính VaR:
VaR3 <- 0.00030883 + qnorm(0.95)*0.02662113
# 10.3. Ước lượng QUANTILE ESTIMATION
VaR4 <- quantile(nreturn_tmt, 0.95)
# 10.4. Ước lượng MONTE CARLO SIMULATION
mean_tmt <- mean(nreturn_tmt)
sd_tmt <- sd(nreturn_tmt)
# Tiến hành mô phỏng Monte Carlo
set.seed(42)
sim1 <- rnorm(10000, mean = mean_tmt, sd = sd_tmt)
ts.plot(sim1)
VaR5 <- quantile(sim1, 0.95)- Var theo số tiền
VaR1_100 = VaR1*100000
VaR1_100[1] 4146.676vaR3_100 = VaR3*100000
vaR3_100[1] 4409.669VaR4_100 = VaR4*100000
VaR4_100 95%
6262.562 VaR5_100 = VaR5*100000
VaR5_100 95%
5649.106