BÀI THỰC HÀNH KIỂM TRA GIỮA KỲ MÔN QUẢN TRỊ RỦI RO TÀI CHÍNH

  1. Đọc dữ liệu từ link google sheet
library(googlesheets4)
gs4_deauth()
df <- read_sheet("https://docs.google.com/spreadsheets/d/10nZIVWn3i5djs0tRZbI--8nntNbfbGoNvYKX8yI7Eeo/edit?gid=1109581498#gid=1109581498", sheet = "hose3")
  1. 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>, …
vn <- df[, c("Date", "VNINDEX")]
vn 
# A tibble: 3,335 × 2
   Date                VNINDEX
   <dttm>                <dbl>
 1 2012-04-11 00:00:00    459.
 2 2012-04-12 00:00:00    465.
 3 2012-04-13 00:00:00    463.
 4 2012-04-16 00:00:00    468.
 5 2012-04-17 00:00:00    473.
 6 2012-04-18 00:00:00    472.
 7 2012-04-19 00:00:00    467.
 8 2012-04-20 00:00:00    466.
 9 2012-04-23 00:00:00    465.
10 2012-04-24 00:00:00    466.
# ℹ 3,325 more rows
  1. Kiểm tra dữ liệu
str(vn)
tibble [3,335 × 2] (S3: tbl_df/tbl/data.frame)
 $ Date   : POSIXct[1:3335], format: "2012-04-11" "2012-04-12" ...
 $ VNINDEX: num [1:3335] 459 465 463 468 473 ...
head(vn)
# A tibble: 6 × 2
  Date                VNINDEX
  <dttm>                <dbl>
1 2012-04-11 00:00:00    459.
2 2012-04-12 00:00:00    465.
3 2012-04-13 00:00:00    463.
4 2012-04-16 00:00:00    468.
5 2012-04-17 00:00:00    473.
6 2012-04-18 00:00:00    472.
tail(vn)
# A tibble: 6 × 2
  Date                VNINDEX
  <dttm>                <dbl>
1 2026-05-15 00:00:00   1922.
2 2026-05-18 00:00:00   1928.
3 2026-05-19 00:00:00   1913.
4 2026-05-20 00:00:00   1913.
5 2026-05-21 00:00:00   1897.
6 2026-05-22 00:00:00   1877.
summary(vn)
      Date                        VNINDEX      
 Min.   :2012-04-11 00:00:00   Min.   : 375.3  
 1st Qu.:2015-09-16 12:00:00   1st Qu.: 594.0  
 Median :2019-04-23 00:00:00   Median : 957.1  
 Mean   :2019-05-03 01:14:16   Mean   : 941.2  
 3rd Qu.:2022-12-27 12:00:00   3rd Qu.:1228.3  
 Max.   :2026-05-22 00:00:00   Max.   :1927.9  

4.1. Thống kê mô tả dữ liệu giá chứng khoán

library(fBasics)
basicStats(vn$VNINDEX) # cách 1
            X..vn.VNINDEX
nobs         3.335000e+03
NAs          0.000000e+00
Minimum      3.752600e+02
Maximum      1.927940e+03
1. Quartile  5.940250e+02
3. Quartile  1.228350e+03
Mean         9.411805e+02
Median       9.571400e+02
Sum          3.138837e+06
SE Mean      6.318758e+00
LCL Mean     9.287914e+02
UCL Mean     9.535695e+02
Variance     1.331555e+05
Stdev        3.649048e+02
Skewness     4.207220e-01
Kurtosis    -5.983190e-01
summary(vn$VNINDEX) # cách 2
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  375.3   594.0   957.1   941.2  1228.3  1927.9 

4.2. Thống kê mô tả dữ liệu return của giá chứng khoán

return_vn <- diff(log(vn$VNINDEX)) 
head(return_vn) 
[1]  0.014112788 -0.005906590  0.012333898  0.009733368 -0.001439154
[6] -0.010817362
tail(return_vn) 
[1] -0.0020067279  0.0032939030 -0.0078159773  0.0001568152 -0.0085772103
[6] -0.0104716883
basicStats(return_vn)
              return_vn
nobs        3334.000000
NAs            0.000000
Minimum       -0.069102
Maximum        0.077067
1. Quartile   -0.004576
3. Quartile    0.006742
Mean           0.000423
Median         0.001194
Sum            1.409016
SE Mean        0.000208
LCL Mean       0.000015
UCL Mean       0.000830
Variance       0.000144
Stdev          0.012011
Skewness      -0.753109
Kurtosis       4.949390

5.1. Vẽ đồ thị và hàm mật độ giá chứng khoán

ts.plot(vn$VNINDEX) # cách 1: đồ thị giá

plot(vn, type = 'l') # cách 2: đồ thị giá

hist(vn$VNINDEX)

d1 = density(vn$VNINDEX)
plot(d1, type = "l")

5.2. Vẽ đồ thị và hàm mật độ của return giá chứng khoán

plot(return_vn, type = 'l')

hist(return_vn)

d2 = density(return_vn)
plot(d2, type = "l")

  1. Kiểm tra tính dừng, kiểm tra phân phối chuẩn
library(tseries)
adf.test(vn$VNINDEX)

    Augmented Dickey-Fuller Test

data:  vn$VNINDEX
Dickey-Fuller = -2.265, Lag order = 14, p-value = 0.4661
alternative hypothesis: stationary
adf.test(return_vn)

    Augmented Dickey-Fuller Test

data:  return_vn
Dickey-Fuller = -14.706, Lag order = 14, p-value = 0.01
alternative hypothesis: stationary
kpss.test(return_vn) # H0: stationary

    KPSS Test for Level Stationarity

data:  return_vn
KPSS Level = 0.042477, Truncation lag parameter = 9, p-value = 0.1
normalTest(vn$VNINDEX, method = c("jb")) # H0: normal distribution

Title:
 Jarque-Bera Normality Test

Test Results:
  STATISTIC:
    X-squared: 147.9806
  P VALUE:
    Asymptotic p Value: < 2.2e-16 
normalTest(return_vn, method = c("jb")) # H0: normal distribution

Title:
 Jarque-Bera Normality Test

Test Results:
  STATISTIC:
    X-squared: 3724.973
  P VALUE:
    Asymptotic p Value: < 2.2e-16 
  1. Kiểm tra tự tương quan bằng đồ thị và bằng kiểm định
acf(return_vn) 

pacf(return_vn) 

Box.test(return_vn, lag = 12) # H0: no autocorrelation

    Box-Pierce test

data:  return_vn
X-squared = 37.152, df = 12, p-value = 0.000211
Box.test(return_vn, lag = 12, type = c("Ljung-Box")) # H0: no autocorrelation

    Box-Ljung test

data:  return_vn
X-squared = 37.24, df = 12, p-value = 0.0002042
  1. Ước lượng mô hình ARIMA
# 8.1. Xây dựng mô hình AR
pacf(return_vn)

m1 <- ar(return_vn, method = c("mle"))
m1

Call:
ar(x = return_vn, method = c("mle"))

Coefficients:
      1        2        3        4        5        6        7        8  
 0.0632   0.0336  -0.0038  -0.0162   0.0286  -0.0156   0.0057   0.0098  
      9       10       11       12  
-0.0322   0.0340   0.0232  -0.0558  

Order selected 12  sigma^2 estimated as  0.0001425
m1$aic
        0         1         2         3         4         5         6         7 
15.325157  4.306591  3.082695  5.078960  6.117447  5.341304  6.663063  8.624042 
        8         9        10        11        12 
10.291214  9.409725  7.676903  8.383212  0.000000 
# Ước lượng AR(p)
m21 <- arima(return_vn, order = c(1, 0, 0))
m21

Call:
arima(x = return_vn, order = c(1, 0, 0))

Coefficients:
         ar1  intercept
      0.0624      4e-04
s.e.  0.0173      2e-04

sigma^2 estimated as 0.0001437:  log likelihood = 10019.08,  aic = -20032.16
m22 <- arima(return_vn, order = c(12, 0, 0))
m22

Call:
arima(x = return_vn, order = c(12, 0, 0))

Coefficients:
         ar1     ar2      ar3      ar4     ar5      ar6     ar7     ar8
      0.0632  0.0336  -0.0038  -0.0161  0.0286  -0.0157  0.0057  0.0098
s.e.  0.0173  0.0173   0.0173   0.0173  0.0173   0.0173  0.0173  0.0173
          ar9    ar10    ar11     ar12  intercept
      -0.0321  0.0339  0.0232  -0.0557      4e-04
s.e.   0.0173  0.0173  0.0173   0.0173      2e-04

sigma^2 estimated as 0.0001425:  log likelihood = 10032.22,  aic = -20036.43
AIC(m21, m22)
    df       AIC
m21  3 -20032.16
m22 14 -20036.43
# 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.062437853 0.0172868313 3.611874 0.0003039925
intercept 0.000422271 0.0002220688 1.901532 0.0572323850
coeftest(m22)
                   coef           se          t            p
ar1        0.0631558558 0.0172937242  3.6519523 0.0002602542
ar2        0.0336015173 0.0173247513  1.9395094 0.0524393461
ar3       -0.0037952099 0.0173247080 -0.2190634 0.8266006396
ar4       -0.0161266241 0.0173153689 -0.9313474 0.3516738805
ar5        0.0286019207 0.0173129367  1.6520548 0.0985233738
ar6       -0.0156638994 0.0173230800 -0.9042214 0.3658780252
ar7        0.0056791668 0.0173218005  0.3278624 0.7430157176
ar8        0.0098392809 0.0173123573  0.5683386 0.5698050757
ar9       -0.0321138522 0.0173092381 -1.8553013 0.0635532103
ar10       0.0339422163 0.0173214056  1.9595532 0.0500480323
ar11       0.0231749045 0.0173209618  1.3379687 0.1809066464
ar12      -0.0557080012 0.0172915663 -3.2216862 0.0012743861
intercept  0.0004215796 0.0002240997  1.8812147 0.0599427210
# 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_vn)

Box.test(m21$residuals, lag = 10)

    Box-Pierce test

data:  m21$residuals
X-squared = 14.949, df = 10, p-value = 0.1339
# 8.2. Xây dựng mô hình MA (chỉ là ví dụ, vì mô hình AR(1) phần dư đã hết tự tương quan)
acf(return_vn)

m3 <- arima(return_vn, order = c(1, 0, 1))
m3

Call:
arima(x = return_vn, order = c(1, 0, 1))

Coefficients:
         ar1      ma1  intercept
      0.3921  -0.3292      4e-04
s.e.  0.1940   0.1988      2e-04

sigma^2 estimated as 0.0001435:  log likelihood = 10020.34,  aic = -20032.67
coeftest(m3)
                   coef           se         t          p
ar1        0.3920670065 0.1940188137  2.020768 0.04330379
ma1       -0.3291660365 0.1987908526 -1.655841 0.09775403
intercept  0.0004165141 0.0002296079  1.814023 0.06967414
# Kiểm định phần dư
plot(m3$residuals, type = 'l')

qqplot(m3$residuals, return_vn)

Box.test(m3$residuals, lag = 10)

    Box-Pierce test

data:  m3$residuals
X-squared = 12.578, df = 10, p-value = 0.2482

9.1. Kiểm định hiệu ứng ARCH (ARCH-LM test)

library(fGarch)
m4 <- garchFit(~ arma(1,0) + garch(1,0), data = return_vn, trace = FALSE)
summary(m4)

Title:
 GARCH Modelling 

Call:
 garchFit(formula = ~arma(1, 0) + garch(1, 0), data = return_vn, 
    trace = FALSE) 

Mean and Variance Equation:
 data ~ arma(1, 0) + garch(1, 0)
<environment: 0x306baa4e0>
 [data = return_vn]

Conditional Distribution:
 norm 

Coefficient(s):
        mu         ar1       omega      alpha1  
0.00054347  0.03762936  0.00010369  0.31379950  

Std. Errors:
 based on Hessian 

Error Analysis:
        Estimate  Std. Error  t value Pr(>|t|)    
mu     5.435e-04   1.896e-04    2.866  0.00416 ** 
ar1    3.763e-02   2.271e-02    1.657  0.09753 .  
omega  1.037e-04   3.413e-06   30.376  < 2e-16 ***
alpha1 3.138e-01   3.315e-02    9.465  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Log Likelihood:
 10148.65    normalized:  3.043986 

Description:
 Tue Jun  9 12:22:16 2026 by user:  



Standardised Residuals Tests:
                                  Statistic      p-Value
 Jarque-Bera Test   R    Chi^2  3333.392168 0.0000000000
 Shapiro-Wilk Test  R    W         0.937498 0.0000000000
 Ljung-Box Test     R    Q(10)    24.669549 0.0060083299
 Ljung-Box Test     R    Q(15)    39.082444 0.0006229258
 Ljung-Box Test     R    Q(20)    41.750664 0.0029800275
 Ljung-Box Test     R^2  Q(10)   207.983937 0.0000000000
 Ljung-Box Test     R^2  Q(15)   267.981576 0.0000000000
 Ljung-Box Test     R^2  Q(20)   277.183072 0.0000000000
 LM Arch Test       R    TR^2    165.246026 0.0000000000

Information Criterion Statistics:
      AIC       BIC       SIC      HQIC 
-6.085573 -6.078240 -6.085576 -6.082949 

9.2. Ước lượng ARIMA+GARCH

m5 <-garchFit(~ arma(1,0) + garch(1,1), data = return_vn, trace = FALSE)
summary(m5)

Title:
 GARCH Modelling 

Call:
 garchFit(formula = ~arma(1, 0) + garch(1, 1), data = return_vn, 
    trace = FALSE) 

Mean and Variance Equation:
 data ~ arma(1, 0) + garch(1, 1)
<environment: 0x30588d4a0>
 [data = return_vn]

Conditional Distribution:
 norm 

Coefficient(s):
        mu         ar1       omega      alpha1       beta1  
5.9490e-04  7.4157e-02  6.1102e-06  1.2476e-01  8.3884e-01  

Std. Errors:
 based on Hessian 

Error Analysis:
        Estimate  Std. Error  t value Pr(>|t|)    
mu     5.949e-04   1.712e-04    3.475 0.000510 ***
ar1    7.416e-02   1.968e-02    3.768 0.000164 ***
omega  6.110e-06   1.050e-06    5.821 5.84e-09 ***
alpha1 1.248e-01   1.325e-02    9.418  < 2e-16 ***
beta1  8.388e-01   1.602e-02   52.367  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Log Likelihood:
 10340.51    normalized:  3.101533 

Description:
 Tue Jun  9 12:22:16 2026 by user:  



Standardised Residuals Tests:
                                   Statistic    p-Value
 Jarque-Bera Test   R    Chi^2  4492.1725549 0.00000000
 Shapiro-Wilk Test  R    W         0.9456598 0.00000000
 Ljung-Box Test     R    Q(10)    15.2766879 0.12229521
 Ljung-Box Test     R    Q(15)    26.4199217 0.03383384
 Ljung-Box Test     R    Q(20)    28.0285591 0.10872504
 Ljung-Box Test     R^2  Q(10)     9.2064302 0.51263121
 Ljung-Box Test     R^2  Q(15)    10.1713349 0.80882626
 Ljung-Box Test     R^2  Q(20)    13.1578268 0.87051012
 LM Arch Test       R    TR^2     10.0226301 0.61397524

Information Criterion Statistics:
      AIC       BIC       SIC      HQIC 
-6.200066 -6.190900 -6.200070 -6.196786 
predict(m5, 5)
   meanForecast   meanError standardDeviation
1 -0.0001816455 0.009490514       0.009490514
2  0.0005814312 0.009664193       0.009638532
3  0.0006380185 0.009805269       0.009779043
4  0.0006422149 0.009939188       0.009912555
5  0.0006425260 0.010066547       0.010039527
  1. Tính và dự báo VaR
# Chạy hàm loss để tính VaR
nreturn_vn <- -return_vn
# 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_vn)
Box.test(residuals(m6, standardize=TRUE), lag=10, type="Ljung")

    Box-Ljung test

data:  residuals(m6, standardize = TRUE)
X-squared = 14.848, df = 10, p-value = 0.1377
ugarchforecast(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 0.0001809 0.009661
# Tính VaR (Phương pháp của JP Morgan xem mean = 0): 
VaR1 <- 0 + qnorm(0.95) * 0.009661
# Cách 2: tính từ RMeasure.R
setwd("~/Library/CloudStorage/OneDrive-Personal/0. Bai giang/QTRRTC/R")
source("RMeasure.R")
VaR2 <- RMeasure(0, 0.009661)

 Risk Measures for selected probabilities: 
      prob        VaR         ES
[1,] 0.950 0.01589093 0.01992787
[2,] 0.990 0.02247485 0.02574863
[3,] 0.999 0.02985473 0.03252946
# 10.2. Ước tính VaR theo ECONOMICTRICS
m7 <- garchFit(~ arma(1,0) + garch(1,1), data = nreturn_vn, trace=F)
summary(m7)

Title:
 GARCH Modelling 

Call:
 garchFit(formula = ~arma(1, 0) + garch(1, 1), data = nreturn_vn, 
    trace = F) 

Mean and Variance Equation:
 data ~ arma(1, 0) + garch(1, 1)
<environment: 0x30e5b9580>
 [data = nreturn_vn]

Conditional Distribution:
 norm 

Coefficient(s):
         mu          ar1        omega       alpha1        beta1  
-5.9490e-04   7.4157e-02   6.1102e-06   1.2476e-01   8.3884e-01  

Std. Errors:
 based on Hessian 

Error Analysis:
         Estimate  Std. Error  t value Pr(>|t|)    
mu     -5.949e-04   1.712e-04   -3.475 0.000510 ***
ar1     7.416e-02   1.968e-02    3.768 0.000164 ***
omega   6.110e-06   1.050e-06    5.821 5.84e-09 ***
alpha1  1.248e-01   1.325e-02    9.418  < 2e-16 ***
beta1   8.388e-01   1.602e-02   52.367  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Log Likelihood:
 10340.51    normalized:  3.101533 

Description:
 Tue Jun  9 12:22:16 2026 by user:  



Standardised Residuals Tests:
                                   Statistic    p-Value
 Jarque-Bera Test   R    Chi^2  4492.1724364 0.00000000
 Shapiro-Wilk Test  R    W         0.9456598 0.00000000
 Ljung-Box Test     R    Q(10)    15.2766878 0.12229521
 Ljung-Box Test     R    Q(15)    26.4199217 0.03383384
 Ljung-Box Test     R    Q(20)    28.0285592 0.10872504
 Ljung-Box Test     R^2  Q(10)     9.2064300 0.51263122
 Ljung-Box Test     R^2  Q(15)    10.1713347 0.80882627
 Ljung-Box Test     R^2  Q(20)    13.1578265 0.87051013
 LM Arch Test       R    TR^2     10.0226299 0.61397526

Information Criterion Statistics:
      AIC       BIC       SIC      HQIC 
-6.200066 -6.190900 -6.200070 -6.196786 
predict(m7, 1)
  meanForecast   meanError standardDeviation
1 0.0001816456 0.009490514       0.009490514
# Tính VaR: 
VaR3 <- 0.0001816456 + qnorm(0.95)*0.009490514
# 10.3. Ước lượng QUANTILE ESTIMATION
VaR4 <- quantile(nreturn_vn, 0.95)
# 10.4. Ước lượng MONTE CARLO SIMULATION
mean_vn <- mean(nreturn_vn)
sd_vn <- sd(nreturn_vn)
# Tiến hành mô phỏng Monte Carlo
set.seed(42)
sim1 <- rnorm(10000, mean = mean_vn, sd = sd_vn)
ts.plot(sim1)

VaR5 <- quantile(sim1, 0.95)
  1. VaR theo số tiền
VaR1_100 = VaR1*100000
VaR1_100
[1] 1589.093
VaR3_100 = VaR3*100000
VaR3_100
[1] 1579.215
VaR4_100 = VaR4*100000
VaR4_100
    95% 
1893.13 
VaR5_100 = VaR5*100000
VaR5_100
     95% 
1935.795