TIỂU LUẬN MÔN CÁC MÔ HÌNH NGẪU NHIÊN

---

library(xlsx) 

## Warning: package 'xlsx' was built under R version 4.3.3

library(quantmod)

## Warning: package 'quantmod' was built under R version 4.3.3

## Loading required package: xts

## Warning: package 'xts' was built under R version 4.3.3

## Loading required package: zoo

## Warning: package 'zoo' was built under R version 4.3.3

## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

## Loading required package: TTR

## Warning: package 'TTR' was built under R version 4.3.3

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

library(tidyverse)

## Warning: package 'ggplot2' was built under R version 4.3.3

## Warning: package 'tidyr' was built under R version 4.3.3

## Warning: package 'dplyr' was built under R version 4.3.3

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.0     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::first()  masks xts::first()
## ✖ dplyr::lag()    masks stats::lag()
## ✖ dplyr::last()   masks xts::last()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(xts)
library(PerformanceAnalytics)

## Warning: package 'PerformanceAnalytics' was built under R version 4.3.3

## 
## Attaching package: 'PerformanceAnalytics'
## 
## The following object is masked from 'package:graphics':
## 
##     legend

library(forecast)

## Warning: package 'forecast' was built under R version 4.3.3

library(rugarch)

## Warning: package 'rugarch' was built under R version 4.3.3

## Loading required package: parallel
## 
## Attaching package: 'rugarch'
## 
## The following object is masked from 'package:purrr':
## 
##     reduce
## 
## The following object is masked from 'package:stats':
## 
##     sigma

library(goftest)
library(VineCopula)

## Warning: package 'VineCopula' was built under R version 4.3.3

library(tseries)

## Warning: package 'tseries' was built under R version 4.3.3

library(FinTS)

## Warning: package 'FinTS' was built under R version 4.3.3

## 
## Attaching package: 'FinTS'
## 
## The following object is masked from 'package:forecast':
## 
##     Acf

library(ggplot2)
library(corrplot)

## Warning: package 'corrplot' was built under R version 4.3.3

## corrplot 0.92 loaded

library(ggcorrplot)

## Warning: package 'ggcorrplot' was built under R version 4.3.3

library(aTSA)

## 
## Attaching package: 'aTSA'
## 
## The following objects are masked from 'package:tseries':
## 
##     adf.test, kpss.test, pp.test
## 
## The following object is masked from 'package:forecast':
## 
##     forecast
## 
## The following object is masked from 'package:graphics':
## 
##     identify

library(ADGofTest)

## 
## Attaching package: 'ADGofTest'
## 
## The following object is masked from 'package:goftest':
## 
##     ad.test

library(kableExtra)

## Warning: package 'kableExtra' was built under R version 4.3.3

## 
## Attaching package: 'kableExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     group_rows

library(knitr)

#Nhập dữ liệu chỉ số giá VN Index của Việt Nam và chỉ số giá CSI300 của Trung Quốc vào lần lượt P1 và P2
P1 <- read.xlsx('C:/Users/HP/Downloads/VNIndexHistoricalData.xlsx', sheetIndex = 1, header = T)
P2 <- read.xlsx('C:/Users/HP/Downloads/ShanghaiCompositeHistoricalData1.xlsx', sheetIndex = 1, header = T)
##Khai báo dữ liệu chuỗi thời gian, sắp xếp theo thời gian thứ tự từ ngày 01/01/2018 đến ngày 01/01/2024
VNI <- xts(P1[,-1], order.by = P1$Date)
CSI300 <- xts(P2[,-1],order.by = P2$Date)
#Gộp hai chuỗi giá dựa trên những ngày trùng nhau
dt <- merge.xts(CSI300,VNI,join = 'inner')
#Tính tỷ suất lợi nhuận theo công thức logt - logt-1
rt <- CalculateReturns(dt,method = 'log')
rt <- rt[-1]
rt2 <- CalculateReturns(rt)
rt2 <- rt2[-1]
# ma hoa so lieu (an)
stt <- 1:nrow(dt)
dt <- data.frame(dt,stt)
stt.tsln <- 1:length(rt)
tsln.1 <- data.frame(rt,stt.tsln)

## Warning in data.frame(rt, stt.tsln): row names were found from a short variable
## and have been discarded

# Bieu do duong CSI300
dt %>% ggplot(aes(x=stt,y=CSI300)) + geom_line() + theme_bw()+
    labs(x="Date",y="Price",title = "Biểu đồ đường chỉ số CSI300")

# Bieu do duong CSI300
dt %>% ggplot(aes(x=stt,y=VNI)) + geom_line() + theme_bw()+
    labs(x="Date",y="Price",title = "Biểu đồ đường chỉ số VNI")

rt_df <- as.data.frame(rt)
a <- rt_df %>% summarise(Min = min(VNI),
                      Max = max(VNI),
                      Mean = mean(VNI),
                      StDev = sd(VNI),
                      Skewness = skewness(VNI),
                      Kurtosis = kurtosis(VNI))
b <- rt_df %>% summarise(Min = min(CSI300),
                      Max = max(CSI300),
                      Mean = mean(CSI300),
                      StDev = sd(CSI300),
                      Skewness = skewness(CSI300),
                      Kurtosis = kurtosis(CSI300))
m <- rbind(a,b)
rownames(m) <- c('VNI','CSI300')
kable(m, format = 'pandoc', caption = 'Bảng 1: Thống kê mô tả chuỗi TSLN', table.attr = "style='width:100%;'") %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed"))

## Warning in kable_styling(., bootstrap_options = c("striped", "hover",
## "condensed")): Please specify format in kable. kableExtra can customize either
## HTML or LaTeX outputs. See https://haozhu233.github.io/kableExtra/ for details.

Bảng 1: Thống kê mô tả chuỗi TSLN

	Min	Max	Mean	StDev	Skewness	Kurtosis
VNI	-0.0690762	0.0486002	0.0000040	0.0127318	-0.9710012	3.891447
CSI300	-0.0820870	0.0577737	-0.0001202	0.0124147	-0.3262891	3.145251

pivot_longer(rt_df, cols = everything(), names_to = "TSLN", values_to = "Value") %>% ggplot(aes(x = TSLN, y = Value)) +
  geom_boxplot(fill = "navy", color = "black") +
  theme_minimal() +
  labs(title = "Hình 7: Biểu đồ hộp của TSLN VNI và CSI300",
       x = "Chỉ số",
       y = "Tỷ suất sinh lợi")

# Bieu do duong rt
ggplot(tsln.1,aes(x=stt.tsln,y=CSI300)) + geom_line() + theme_bw()+
    labs(x="Date",y="%",title = "Biểu đồ đường tỷ suât lợi nhuận chỉ số CSI300")

# Bieu do duong tsln_vni
ggplot(tsln.1,aes(x=stt.tsln,y=VNI)) + geom_line() + theme_bw()+
    labs(x="Date",y="%",title = "Biểu đồ đường tỷ suât lợi nhuận chỉ số VNI")

#Kiểm định tính dừng
adf_vni <- adf.test(rt$VNI)

## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag  ADF p.value
## [1,]   0 36.1    0.99
## [2,]   1 49.9    0.99
## [3,]   2 60.7    0.99
## [4,]   3 69.9    0.99
## [5,]   4 77.4    0.99
## [6,]   5 83.9    0.99
## [7,]   6 89.6    0.99
## [8,]   7 94.8    0.99
## Type 2: with drift no trend 
##      lag  ADF p.value
## [1,]   0 36.1    0.99
## [2,]   1 49.9    0.99
## [3,]   2 60.7    0.99
## [4,]   3 69.9    0.99
## [5,]   4 77.4    0.99
## [6,]   5 83.8    0.99
## [7,]   6 89.6    0.99
## [8,]   7 94.8    0.99
## Type 3: with drift and trend 
##      lag  ADF p.value
## [1,]   0 36.0    0.99
## [2,]   1 49.9    0.99
## [3,]   2 60.7    0.99
## [4,]   3 69.9    0.99
## [5,]   4 77.3    0.99
## [6,]   5 83.8    0.99
## [7,]   6 89.6    0.99
## [8,]   7 94.8    0.99
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

adf_CSI300 <- adf.test(rt$CSI300)

## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0  38.1    0.99
## [2,]   1  54.0    0.99
## [3,]   2  65.2    0.99
## [4,]   3  75.5    0.99
## [5,]   4  85.0    0.99
## [6,]   5  94.4    0.99
## [7,]   6 102.5    0.99
## [8,]   7 110.1    0.99
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0  38.1    0.99
## [2,]   1  53.9    0.99
## [3,]   2  65.2    0.99
## [4,]   3  75.5    0.99
## [5,]   4  85.0    0.99
## [6,]   5  94.4    0.99
## [7,]   6 102.5    0.99
## [8,]   7 110.1    0.99
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0  38.1    0.99
## [2,]   1  54.0    0.99
## [3,]   2  65.2    0.99
## [4,]   3  75.5    0.99
## [5,]   4  85.0    0.99
## [6,]   5  94.5    0.99
## [7,]   6 102.6    0.99
## [8,]   7 110.3    0.99
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

#Kiểm định phân phối chuẩn
jq_vni <- jarque.bera.test(rt$VNI)
jq_CSI300 <- jarque.bera.test(rt$CSI300)
# Ước lượng và trích xuất phần dư từ mô hình ARMA tối ưu
arima_VNI <- autoarfima(rt$VNI,ar.max = 2, ma.max = 2, criterion = 'AIC', method = "full")

arima_CSI300 <- autoarfima(rt$CSI300,ar.max = 2, ma.max = 2, criterion = 'AIC', method = "full")

re_VNI <- arima_VNI$fit@fit$residuals
re_CSI300 <- arima_CSI300$fit@fit$residuals
#Kiểm định tương quan chuỗi bậc 2 cho phần dư
lj_vni <- Box.test(re_VNI,type = 'Ljung-Box', lag = 2)
lj_CSI300 <- Box.test(re_CSI300,type = 'Ljung-Box', lag = 2)

#Kiểm định tương quan chuỗi bậc 2 cho phần dư bình phương
lj_vni.2 <- Box.test(re_VNI^2,type = 'Ljung-Box', lag = 2)
lj_CSI300.2 <- Box.test(re_CSI300^2,type = 'Ljung-Box', lag = 2)

#Kiểm định hiệu ứng ARCH
ar_vni <- ArchTest(re_VNI, lags = 2)
ar_CSI300 <- ArchTest(re_CSI300, lags = 2)

#Trình bày kết quả
test_result <- data.frame(Test = c("J-B","Q(2)","Q(2)^2", "ARCH(2)"),
                          P_value_VNI = c(jq_vni$p.value, lj_vni$p.value, lj_vni.2$p.value, ar_vni$p.value),
                          P_value_CSI300 = c(jq_CSI300$p.value, lj_CSI300$p.value, lj_CSI300.2$p.value, ar_CSI300$p.value))
kable(test_result, 
  caption = "Bảng 3: Kết quả các kiểm định", 
  label = 'Ghi chú: Q (2) and Q2 (2) lần lượt là kiểm định Ljung-Box Q2 cho tương quan chuỗi bậc 2 của phần dư và bình phương phần dư của lợi suất', 
  format = 'pandoc') %>% kable_styling(
            bootstrap_options = c("striped", "hover", "condensed"), 
            full_width = F)

## Warning in kable_styling(., bootstrap_options = c("striped", "hover",
## "condensed"), : Please specify format in kable. kableExtra can customize either
## HTML or LaTeX outputs. See https://haozhu233.github.io/kableExtra/ for details.

Bảng 3: Kết quả các kiểm định

Test	P_value_VNI	P_value_CSI300
J-B	0.0000000	0.0000000
Q(2)	0.9998155	0.9997269
Q(2)^2	0.0000000	0.0000217
ARCH(2)	0.0000000	0.0000534

pearson <- cor(rt$VNI,rt$CSI300, method="pearson")
spearman <- cor(rt$VNI,rt$CSI300, method="spearman")
kendall <- cor(rt$VNI,rt$CSI300, method="kendall")

#Trình bày kết quả
relat <- data.frame('Tương quan' = 'VNI-CSI300',
                    Pearson = pearson[1,1],
                    spearman = spearman[1,1],
                    Kendall = kendall[1,1])
kable(relat, 
      col.names = c("Phương pháp", "Pearson", "Spearman","Kendall"),
      caption = "Bảng 4: Kết quả hệ số tương quan",
      format = 'pandoc',
      align = c("l", "c", "c")) %>%
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)

## Warning in kable_styling(., bootstrap_options = c("striped", "hover",
## "condensed"), : Please specify format in kable. kableExtra can customize either
## HTML or LaTeX outputs. See https://haozhu233.github.io/kableExtra/ for details.

Bảng 4: Kết quả hệ số tương quan

Phương pháp	Pearson	Spearman	Kendall
VNI-CSI300	0.0401281	0.042613	0.0286111

par(mfrow = c(1,2))

corr <- cor(rt)
ggcorrplot(corr, hc.order = TRUE,
   outline.col = "white",
   ggtheme = ggplot2::theme_gray,
   colors = c("navy", "white", "gold"),
   lab = TRUE, lab_col = 'black',title = 'Hình 8: Trực quan hóa hệ số tương quan với phương pháp Pearson')

chart.Correlation(rt, histogram=TRUE, pch=19)

## Warning in par(usr): argument 1 does not name a graphical parameter

rt %>% ggplot(aes(VNI, CSI300)) +
              geom_point(col = 'gold',shape = TRUE) + 
              geom_smooth(method = 'lm',se = T, col = 'navy') + #thêm đường hồi quy với phương pháp hồi quy tuyến tính  
              labs(title = 'Hình 11: Biểu đồ Scatter TSLN của CSI300 và VNI', x = 'VNI', y = 'CSI300')

## `geom_smooth()` using formula = 'y ~ x'

#CSI300

arma_CSI300 <- autoarfima(rt$CSI300,ar.max = 2,ma.max = 2,
                          criterion = "AIC",
                          method = "full")
arma_CSI300$fit

## 
## *----------------------------------*
## *          ARFIMA Model Fit        *
## *----------------------------------*
## Mean Model   : ARFIMA(0,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##        Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.00012    0.000325 -0.36937  0.71185
## sigma   0.01241    0.000230 53.96314  0.00000
## 
## Robust Standard Errors:
##        Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.00012    0.000314 -0.38319  0.70158
## sigma   0.01241    0.000512 24.23342  0.00000
## 
## LogLikelihood : 4324.731 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -5.9378
## Bayes        -5.9306
## Shibata      -5.9378
## Hannan-Quinn -5.9351
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                  5.778e-05  0.9939
## Lag[2*(p+q)+(p+q)-1][2] 3.020e-04  0.9996
## Lag[4*(p+q)+(p+q)-1][5] 2.806e+00  0.4436
## 
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic   p-value
## Lag[1]                      16.46 4.979e-05
## Lag[2*(p+q)+(p+q)-1][2]     18.96 8.324e-06
## Lag[4*(p+q)+(p+q)-1][5]     28.52 8.241e-08
## 
## 
## ARCH LM Tests
## ------------------------------------
##              Statistic DoF   P-Value
## ARCH Lag[2]      19.67   2 5.344e-05
## ARCH Lag[5]      32.35   5 5.072e-06
## ARCH Lag[10]     51.43  10 1.452e-07
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.9296
## Individual Statistics:            
## mu    0.1506
## sigma 1.6347
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          0.61 0.749 1.07
## Individual Statistic:     0.35 0.47 0.75
## 
## 
## Elapsed time : 0.03497815

plot(arma_CSI300[["fit"]]@fit[["residuals"]],type="l")

ggplot(tsln.1,aes(x=stt.tsln,y=CSI300)) + geom_line() + theme_bw()+
    labs(x="Date",y="%",title = "Biểu đồ đường tỷ suât lợi nhuận chỉ số CSI300")

# vni

arma_vni <- autoarfima(rt$VNI,ar.max = 2,ma.max = 2,
                       criterion = "AIC",
                       method = "full")
arma_vni$fit

## 
## *----------------------------------*
## *          ARFIMA Model Fit        *
## *----------------------------------*
## Mean Model   : ARFIMA(0,0,2)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##        Estimate  Std. Error  t value Pr(>|t|)
## ma1    0.053837    0.026176   2.0568 0.039710
## ma2    0.045130    0.026414   1.7085 0.087536
## sigma  0.012695    0.000235  53.9746 0.000000
## 
## Robust Standard Errors:
##        Estimate  Std. Error  t value Pr(>|t|)
## ma1    0.053837    0.031428   1.7130 0.086705
## ma2    0.045130    0.032955   1.3694 0.170860
## sigma  0.012695    0.000696  18.2503 0.000000
## 
## LogLikelihood : 4291.54 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -5.8909
## Bayes        -5.8800
## Shibata      -5.8909
## Hannan-Quinn -5.8868
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                  4.634e-05  0.9946
## Lag[2*(p+q)+(p+q)-1][5] 4.618e-01  1.0000
## Lag[4*(p+q)+(p+q)-1][9] 2.543e+00  0.9475
## 
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic   p-value
## Lag[1]                      52.78 3.736e-13
## Lag[2*(p+q)+(p+q)-1][2]     82.04 0.000e+00
## Lag[4*(p+q)+(p+q)-1][5]    151.05 0.000e+00
## 
## 
## ARCH LM Tests
## ------------------------------------
##              Statistic DoF P-Value
## ARCH Lag[2]      93.19   2       0
## ARCH Lag[5]     130.62   5       0
## ARCH Lag[10]    163.52  10       0
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  0.9205
## Individual Statistics:             
## ma1   0.05523
## ma2   0.42330
## sigma 0.46825
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          0.846 1.01 1.35
## Individual Statistic:     0.35 0.47 0.75
## 
## 
## Elapsed time : 0.0324738

plot(arma_CSI300[["fit"]]@fit[["residuals"]],type="l")

ggplot(tsln.1,aes(x=stt.tsln,y=VNI)) + geom_line() + theme_bw()+
    labs(x="Date",y="%",title = "Biểu đồ đường tỷ suât lợi nhuận chỉ số VNI")

# CSI300

# neu nhu arma co thanh phan nao = 0 thi se dung cau lenh nay de setup "fixed.pars = list(ar1 = 0, ma1 = 0)"
garch_CSI300norm <- ugarchspec(mean.model = list(armaOrder = c(2,0),include.mean = FALSE),
                 variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                 fixed.pars = list( ar1 = 0),
                 distribution.model = 'norm')

garch_CSI300std <- ugarchspec(mean.model = list(armaOrder = c(0,2),include.mean = TRUE),
                           variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                           fixed.pars = list(ar1 = 0, ma1 = 0),
                           distribution.model = 'std')

garch_CSI300sstd <- ugarchspec(mean.model = list(armaOrder = c(0,2),include.mean = TRUE),
                           variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                           fixed.pars = list(ar1 = 0, ma1 = 0),
                           distribution.model = 'sstd')

garch_CSI300ged <- ugarchspec(mean.model = list(armaOrder = c(0,2),include.mean = TRUE),
                           variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                           fixed.pars = list(ar1 = 0, ma1 = 0),
                           distribution.model = 'ged')

garch_CSI300sged <- ugarchspec(mean.model = list(armaOrder = c(3,3),include.mean = TRUE),
                           variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                           fixed.pars = list( ma1 = 0),
                           distribution.model = 'sged')


fCSI300norm <- ugarchfit(data = rt$CSI300, spec = garch_CSI300norm)
fCSI300norm

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(2,0,0)
## Distribution : norm 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error   t value Pr(>|t|)
## ar1     0.000000          NA        NA       NA
## ar2     0.000511    0.028049  0.018222 0.985462
## omega   0.000006    0.000001  4.549101 0.000005
## alpha1  0.089938    0.009050  9.937701 0.000000
## beta1   0.876468    0.012060 72.675935 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error   t value Pr(>|t|)
## ar1     0.000000          NA        NA       NA
## ar2     0.000511    0.027906  0.018315 0.985388
## omega   0.000006    0.000002  2.806961 0.005001
## alpha1  0.089938    0.016603  5.416942 0.000000
## beta1   0.876468    0.017861 49.071045 0.000000
## 
## LogLikelihood : 4391.963 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.0274
## Bayes        -6.0129
## Shibata      -6.0274
## Hannan-Quinn -6.0220
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.9293 0.33505
## Lag[2*(p+q)+(p+q)-1][5]    4.0869 0.05319
## Lag[4*(p+q)+(p+q)-1][9]    7.3210 0.10014
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      1.523  0.2171
## Lag[2*(p+q)+(p+q)-1][5]     3.164  0.3779
## Lag[4*(p+q)+(p+q)-1][9]     3.797  0.6234
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.2657 0.500 2.000  0.6062
## ARCH Lag[5]    0.6560 1.440 1.667  0.8368
## ARCH Lag[7]    1.0583 2.315 1.543  0.9039
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  2.7422
## Individual Statistics:              
## ar2    0.02704
## omega  1.20411
## alpha1 0.35941
## beta1  0.37964
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.07 1.24 1.6
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value    prob sig
## Sign Bias            1.086 0.27746    
## Negative Sign Bias   1.709 0.08771   *
## Positive Sign Bias   1.754 0.07961   *
## Joint Effect         7.092 0.06901   *
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     40.87     0.002513
## 2    30     51.09     0.006877
## 3    40     62.41     0.010048
## 4    50     72.98     0.014733
## 
## 
## Elapsed time : 0.3173919

fCSI300std <- ugarchfit(data = rt$CSI300, spec = garch_CSI300std)
fCSI300std

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,2)
## Distribution : std 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000036    0.000274 -0.13059 0.896102
## ma1     0.000000          NA       NA       NA
## ma2    -0.008914    0.026927 -0.33105 0.740608
## omega   0.000005    0.000002  1.87711 0.060503
## alpha1  0.058247    0.005326 10.93568 0.000000
## beta1   0.912778    0.011176 81.67083 0.000000
## shape   6.098477    0.958917  6.35975 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000036    0.000272 -0.13124 0.895588
## ma1     0.000000          NA       NA       NA
## ma2    -0.008914    0.027525 -0.32385 0.746049
## omega   0.000005    0.000006  0.72569 0.468027
## alpha1  0.058247    0.029641  1.96511 0.049401
## beta1   0.912778    0.012022 75.92457 0.000000
## shape   6.098477    1.007772  6.05145 0.000000
## 
## LogLikelihood : 4427.999 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.0742
## Bayes        -6.0524
## Shibata      -6.0742
## Hannan-Quinn -6.0661
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                     0.8627 0.35300
## Lag[2*(p+q)+(p+q)-1][5]    4.4370 0.02000
## Lag[4*(p+q)+(p+q)-1][9]    7.7042 0.07372
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      4.815 0.02822
## Lag[2*(p+q)+(p+q)-1][5]     5.749 0.10302
## Lag[4*(p+q)+(p+q)-1][9]     6.392 0.25564
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]     1.030 0.500 2.000  0.3101
## ARCH Lag[5]     1.073 1.440 1.667  0.7114
## ARCH Lag[7]     1.518 2.315 1.543  0.8177
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.3011
## Individual Statistics:              
## mu     0.46747
## ma2    0.04994
## omega  0.19540
## alpha1 0.27554
## beta1  0.28594
## shape  0.10857
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.49 1.68 2.12
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value    prob sig
## Sign Bias            1.046 0.29589    
## Negative Sign Bias   2.562 0.01052  **
## Positive Sign Bias   1.104 0.26991    
## Joint Effect         9.224 0.02645  **
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     22.65      0.25293
## 2    30     24.67      0.69510
## 3    40     38.45      0.49475
## 4    50     62.06      0.09961
## 
## 
## Elapsed time : 0.5501051

fCSI300sstd <- ugarchfit(data = rt$CSI300, spec = garch_CSI300sstd)
fCSI300sstd

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,2)
## Distribution : sstd 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000081    0.000291 -0.28014 0.779370
## ma1     0.000000          NA       NA       NA
## ma2    -0.009773    0.027041 -0.36143 0.717778
## omega   0.000005    0.000003  1.74394 0.081169
## alpha1  0.058285    0.006050  9.63397 0.000000
## beta1   0.912896    0.011116 82.12684 0.000000
## skew    0.983515    0.035662 27.57871 0.000000
## shape   6.143533    0.975268  6.29933 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000081    0.000287 -0.28391 0.776482
## ma1     0.000000          NA       NA       NA
## ma2    -0.009773    0.027733 -0.35241 0.724533
## omega   0.000005    0.000007  0.63218 0.527267
## alpha1  0.058285    0.032797  1.77715 0.075544
## beta1   0.912896    0.013238 68.95840 0.000000
## skew    0.983515    0.034492 28.51405 0.000000
## shape   6.143533    1.028746  5.97187 0.000000
## 
## LogLikelihood : 4428.103 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.0729
## Bayes        -6.0475
## Shibata      -6.0730
## Hannan-Quinn -6.0635
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      0.865 0.35234
## Lag[2*(p+q)+(p+q)-1][5]     4.486 0.01728
## Lag[4*(p+q)+(p+q)-1][9]     7.756 0.07067
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      4.793 0.02858
## Lag[2*(p+q)+(p+q)-1][5]     5.724 0.10440
## Lag[4*(p+q)+(p+q)-1][9]     6.356 0.25935
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]     1.014 0.500 2.000  0.3141
## ARCH Lag[5]     1.057 1.440 1.667  0.7161
## ARCH Lag[7]     1.494 2.315 1.543  0.8225
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  1.5451
## Individual Statistics:              
## mu     0.46441
## ma2    0.04753
## omega  0.18041
## alpha1 0.28017
## beta1  0.29102
## skew   0.47432
## shape  0.11022
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.69 1.9 2.35
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value    prob sig
## Sign Bias           0.9462 0.34421    
## Negative Sign Bias  2.5185 0.01189  **
## Positive Sign Bias  1.0512 0.29335    
## Joint Effect        9.0812 0.02823  **
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     23.73       0.2069
## 2    30     23.81       0.7384
## 3    40     37.96       0.5174
## 4    50     57.94       0.1788
## 
## 
## Elapsed time : 0.7603269

fCSI300ged <- ugarchfit(data = rt$CSI300, spec = garch_CSI300ged)
fCSI300ged

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,2)
## Distribution : ged 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000086    0.000221 -0.39024 0.696361
## ma1     0.000000          NA       NA       NA
## ma2    -0.023492    0.026904 -0.87317 0.382571
## omega   0.000005    0.000002  3.20791 0.001337
## alpha1  0.068428    0.007209  9.49174 0.000000
## beta1   0.898666    0.012497 71.91135 0.000000
## shape   1.343409    0.067212 19.98752 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## mu     -0.000086    0.000186 -0.46252 0.643710
## ma1     0.000000          NA       NA       NA
## ma2    -0.023492    0.028548 -0.82288 0.410574
## omega   0.000005    0.000002  2.10913 0.034933
## alpha1  0.068428    0.013557  5.04750 0.000000
## beta1   0.898666    0.010707 83.93632 0.000000
## shape   1.343409    0.078166 17.18653 0.000000
## 
## LogLikelihood : 4424.856 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.0699
## Bayes        -6.0481
## Shibata      -6.0699
## Hannan-Quinn -6.0617
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic  p-value
## Lag[1]                     0.9497 0.329790
## Lag[2*(p+q)+(p+q)-1][5]    5.5571 0.000443
## Lag[4*(p+q)+(p+q)-1][9]    8.9475 0.025064
## d.o.f=2
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      3.322 0.06835
## Lag[2*(p+q)+(p+q)-1][5]     4.384 0.20989
## Lag[4*(p+q)+(p+q)-1][9]     4.951 0.43660
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.6595 0.500 2.000  0.4167
## ARCH Lag[5]    0.7628 1.440 1.667  0.8044
## ARCH Lag[7]    1.1752 2.315 1.543  0.8835
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  2.2522
## Individual Statistics:              
## mu     0.49608
## ma2    0.02386
## omega  0.71209
## alpha1 0.32890
## beta1  0.34651
## shape  0.21824
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.49 1.68 2.12
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value    prob sig
## Sign Bias            1.011 0.31220    
## Negative Sign Bias   2.247 0.02476  **
## Positive Sign Bias   1.312 0.18958    
## Joint Effect         8.185 0.04234  **
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     16.58       0.6181
## 2    30     28.26       0.5041
## 3    40     30.65       0.8279
## 4    50     52.59       0.3370
## 
## 
## Elapsed time : 0.783998

fCSI300sged <- ugarchfit(data = rt$CSI300, spec = garch_CSI300sged)

## Warning in arima(data, order = c(modelinc[2], 0, modelinc[3]), include.mean =
## modelinc[1], : possible convergence problem: optim gave code = 1

fCSI300sged

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(3,0,3)
## Distribution : sged 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error   t value Pr(>|t|)
## mu     -0.000165    0.000324  -0.50788 0.611538
## ar1    -0.011627    0.028078  -0.41408 0.678816
## ar2    -0.924121    0.023691 -39.00738 0.000000
## ar3    -0.067023    0.189074  -0.35448 0.722981
## ma1     0.000000          NA        NA       NA
## ma2     0.914926    0.040845  22.40001 0.000000
## ma3     0.071187    0.210391   0.33836 0.735094
## omega   0.000005    0.000002   3.12513 0.001777
## alpha1  0.068340    0.007258   9.41570 0.000000
## beta1   0.898900    0.012508  71.86375 0.000000
## skew    0.984190    0.039129  25.15258 0.000000
## shape   1.356402    0.069037  19.64748 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error   t value Pr(>|t|)
## mu     -0.000165    0.000511  -0.32239 0.747156
## ar1    -0.011627    0.051311  -0.22659 0.820741
## ar2    -0.924121    0.026718 -34.58819 0.000000
## ar3    -0.067023    0.577031  -0.11615 0.907533
## ma1     0.000000          NA        NA       NA
## ma2     0.914926    0.104017   8.79592 0.000000
## ma3     0.071187    0.652867   0.10904 0.913172
## omega   0.000005    0.000003   1.94617 0.051634
## alpha1  0.068340    0.014150   4.82965 0.000001
## beta1   0.898900    0.011239  79.98399 0.000000
## skew    0.984190    0.077503  12.69867 0.000000
## shape   1.356402    0.082990  16.34416 0.000000
## 
## LogLikelihood : 4426.585 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.0654
## Bayes        -6.0254
## Shibata      -6.0655
## Hannan-Quinn -6.0505
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                          statistic p-value
## Lag[1]                       1.799 0.17986
## Lag[2*(p+q)+(p+q)-1][17]    10.279 0.01877
## Lag[4*(p+q)+(p+q)-1][29]    13.942 0.62195
## d.o.f=6
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                      3.590 0.05812
## Lag[2*(p+q)+(p+q)-1][5]     4.561 0.19194
## Lag[4*(p+q)+(p+q)-1][9]     5.101 0.41476
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.5457 0.500 2.000  0.4601
## ARCH Lag[5]    0.6388 1.440 1.667  0.8420
## ARCH Lag[7]    1.0646 2.315 1.543  0.9028
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  3.0517
## Individual Statistics:              
## mu     0.43568
## ar1    0.04773
## ar2    0.22416
## ar3    0.06763
## ma2    0.18256
## ma3    0.06810
## omega  0.68783
## alpha1 0.32857
## beta1  0.36617
## skew   0.44371
## shape  0.23899
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          2.49 2.75 3.27
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value    prob sig
## Sign Bias            1.213 0.22536    
## Negative Sign Bias   2.423 0.01553  **
## Positive Sign Bias   1.347 0.17831    
## Joint Effect         8.760 0.03267  **
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     19.27       0.4393
## 2    30     27.64       0.5372
## 3    40     39.99       0.4260
## 4    50     45.30       0.6237
## 
## 
## Elapsed time : 2.990064

# vni 
garch_vni <- ugarchspec(mean.model = list(armaOrder = c(0,0),include.mean = TRUE),
                           variance.model = list(model = "sGARCH",garchOrder = c(1,1)),
                           distribution.model = 'ged')
fvni <- ugarchfit(data = tsln.1$VNI, spec = garch_vni)
fvni

## 
## *---------------------------------*
## *          GARCH Model Fit        *
## *---------------------------------*
## 
## Conditional Variance Dynamics    
## -----------------------------------
## GARCH Model  : sGARCH(1,1)
## Mean Model   : ARFIMA(0,0,0)
## Distribution : ged 
## 
## Optimal Parameters
## ------------------------------------
##         Estimate  Std. Error  t value Pr(>|t|)
## mu      0.000397    0.000219   1.8174 0.069152
## omega   0.000005    0.000004   1.3272 0.184434
## alpha1  0.125209    0.002771  45.1776 0.000000
## beta1   0.844694    0.015303  55.1964 0.000000
## shape   1.080865    0.032809  32.9438 0.000000
## 
## Robust Standard Errors:
##         Estimate  Std. Error  t value Pr(>|t|)
## mu      0.000397    0.000522  0.76074 0.446813
## omega   0.000005    0.000019  0.29174 0.770486
## alpha1  0.125209    0.064761  1.93341 0.053186
## beta1   0.844694    0.092934  9.08917 0.000000
## shape   1.080865    0.095245 11.34831 0.000000
## 
## LogLikelihood : 9084.779 
## 
## Information Criteria
## ------------------------------------
##                     
## Akaike       -6.2361
## Bayes        -6.2258
## Shibata      -6.2361
## Hannan-Quinn -6.2324
## 
## Weighted Ljung-Box Test on Standardized Residuals
## ------------------------------------
##                         statistic   p-value
## Lag[1]                      11.17 0.0008335
## Lag[2*(p+q)+(p+q)-1][2]     12.32 0.0004672
## Lag[4*(p+q)+(p+q)-1][5]     14.56 0.0005840
## d.o.f=0
## H0 : No serial correlation
## 
## Weighted Ljung-Box Test on Standardized Squared Residuals
## ------------------------------------
##                         statistic p-value
## Lag[1]                   0.004578  0.9461
## Lag[2*(p+q)+(p+q)-1][5]  0.208392  0.9919
## Lag[4*(p+q)+(p+q)-1][9]  0.398738  0.9993
## d.o.f=2
## 
## Weighted ARCH LM Tests
## ------------------------------------
##             Statistic Shape Scale P-Value
## ARCH Lag[3]    0.2128 0.500 2.000  0.6446
## ARCH Lag[5]    0.3622 1.440 1.667  0.9233
## ARCH Lag[7]    0.3832 2.315 1.543  0.9877
## 
## Nyblom stability test
## ------------------------------------
## Joint Statistic:  9.32
## Individual Statistics:              
## mu     5.45552
## omega  0.58299
## alpha1 0.04967
## beta1  0.05092
## shape  0.07576
## 
## Asymptotic Critical Values (10% 5% 1%)
## Joint Statistic:          1.28 1.47 1.88
## Individual Statistic:     0.35 0.47 0.75
## 
## Sign Bias Test
## ------------------------------------
##                    t-value      prob sig
## Sign Bias           3.5976 3.266e-04 ***
## Negative Sign Bias  0.3066 7.592e-01    
## Positive Sign Bias  1.6889 9.135e-02   *
## Joint Effect       38.8957 1.826e-08 ***
## 
## 
## Adjusted Pearson Goodness-of-Fit Test:
## ------------------------------------
##   group statistic p-value(g-1)
## 1    20     125.4    1.075e-17
## 2    30     157.5    1.317e-19
## 3    40     183.0    1.593e-20
## 4    50     224.6    2.581e-24
## 
## 
## Elapsed time : 0.9870448

# Sau khi lựa chọn mô hình tốt thông qua AIC thì sẽ đặt tên chốt mô hình
# VD mô hình của CSI300 có phân phối sged là tốt nhất thì sẽ chốt mô hình

bfCSI300 <- fCSI300sged

bfvni <-  fvni

# lay phan du của mô hình
rs_CSI300 <- residuals(bfCSI300)/sigma(bfCSI300)
rs_vni <- residuals(bfvni)/sigma(bfvni)

# bieu do histogram của phan du CSI300
ggplot(rs_CSI300, aes(x=rs_CSI300)) +
    geom_histogram(aes(y = ..density..),binwidth = 0.05, fill = "navy", alpha = 0.5) +
    geom_density(color = "red") +
    labs(title = "Histogram with Density Plot", x = "Value", y = "Density") +
    theme_minimal()

## Warning: The dot-dot notation (`..density..`) was deprecated in ggplot2 3.4.0.
## ℹ Please use `after_stat(density)` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## Don't know how to automatically pick scale for object of type <xts/zoo>.
## Defaulting to continuous.

# bieu do histogram của phan du vni 
ggplot(rs_vni, aes(x=rs_vni)) +
    geom_histogram(aes(y = ..density..),binwidth = 0.05, fill = "navy", alpha = 0.5) +
    geom_density(color = "red") +
    labs(title = "Histogram with Density Plot", x = "Value", y = "Density") +
    theme_minimal()

## Don't know how to automatically pick scale for object of type <xts/zoo>.
## Defaulting to continuous.

# kiểm định xem phân dư có phải là phân phối đó không(distribution)
# Nếu như mà pvalue lon hơn 0.05 thì đúng là phân phối đó

fitdist(distribution = "sged", rs_CSI300, control = list())

## $pars
##         mu      sigma       skew      shape 
## -0.0051555  0.9971507  0.9818874  1.3585151 
## 
## $convergence
## [1] 0
## 
## $values
## [1] 2067.993 2034.143 2034.143
## 
## $lagrange
## [1] 0
## 
## $hessian
##            [,1]       [,2]       [,3]      [,4]
## [1,] 1661.75486  157.83797 -450.98834 -20.36565
## [2,]  157.83797 1990.16671  -56.07376 154.84044
## [3,] -450.98834  -56.07376  979.81256  77.29642
## [4,]  -20.36565  154.84044   77.29642 239.36588
## 
## $ineqx0
## NULL
## 
## $nfuneval
## [1] 90
## 
## $outer.iter
## [1] 2
## 
## $elapsed
## Time difference of 0.296145 secs
## 
## $vscale
## [1] 1 1 1 1 1

ppCSI300 <- pdist(distribution = "sged", q = rs_CSI300, mu = -0.01939748 , sigma = 1.00046879,
           skew =0.90303842,
           shape = 1.32455112) 
ad.test(ppCSI300,punif) 

## 
##  Anderson-Darling GoF Test
## 
## data:  ppCSI300  and  punif
## AD = 1.1851, p-value = 0.2733
## alternative hypothesis: NA

fitdist(distribution = "sged", rs_vni, control = list())

## $pars
##          mu       sigma        skew       shape 
## -0.03945386  0.98548017  0.82952145  1.15719938 
## 
## $convergence
## [1] 0
## 
## $values
## [1] 4159.575 3926.430 3926.430
## 
## $lagrange
## [1] 0
## 
## $hessian
##            [,1]      [,2]       [,3]     [,4]
## [1,]  6275.5848 1953.1649 -4675.3526 134.3469
## [2,]  1953.1649 3946.4692  -927.5763 551.7278
## [3,] -4675.3526 -927.5763  8339.2922 293.9490
## [4,]   134.3469  551.7278   293.9490 716.0586
## 
## $ineqx0
## NULL
## 
## $nfuneval
## [1] 103
## 
## $outer.iter
## [1] 2
## 
## $elapsed
## Time difference of 0.728766 secs
## 
## $vscale
## [1] 1 1 1 1 1

ppvni <- pdist(distribution = "sged", q = rs_vni, mu = -0.09274464 , sigma = 0.98803583,
               skew =0.84413499,
               shape = 1.09437326)
ad.test(ppvni,punif) 

## 
##  Anderson-Darling GoF Test
## 
## data:  ppvni  and  punif
## AD = 15.546, p-value = 2.06e-07
## alternative hypothesis: NA

# Giả sử u1 dài hơn u2, cắt u1 cho khớp với u2
if (length(ppCSI300) > length(ppvni)) {
    ppCSI300 <- ppCSI300[1:length(ppvni)]
} else if (length(ppvni) > length(ppCSI300)) {
    ppvni <- ppvni[1:length(ppCSI300)]
}
# ước lượng mô hình copula

BiCopSelect(ppCSI300, ppvni, familyset= 1:9, selectioncrit="AIC",indeptest = FALSE, level = 0.05)

## Bivariate copula: Gaussian (par = 0.05, tau = 0.03)

Stu <- BiCopEst(ppCSI300, ppvni, family = 14, method = "mle", se = T, max.df = 10) 
summary(Stu) 

## Family
## ------ 
## No:    14
## Name:  Survival Gumbel
## 
## Parameter(s)
## ------------
## par:  1.02  (SE = 0.01)
## 
## Dependence measures
## -------------------
## Kendall's tau:    0.02 (empirical = 0.03, p value = 0.15)
## Upper TD:         0 
## Lower TD:         0.02 
## 
## Fit statistics
## --------------
## logLik:  0.84 
## AIC:    0.32 
## BIC:    5.6