Excel có 2 công cụ tuyệt vời là solver và goal seek. Với những dự án nhỏ, hoặc dữ liệu nhỏ, mô hình ít thì 2 công cụ này cũng đủ dùng. Nhưng khi dữ liệu nhiều, số lượng mô hình tăng, excel không còn phù hợp. Có nhiều công cụ thay thế excel giải quyết những vẫn đề này, trong đó có R là 1 công cụ rất mạnh.
Trong phần này tôi sử dụng hàm optim của gói stats. Các bước như sau:
Bước 1: Viết hàm sử dụng cho mục đích tối ưu, tìm min hoặc max của 1 object
Bước 2: Sử dụng hàm optim để tìm tham số phù hợp
Viết hàm tìm tham số phù hợp với dạng phân phối.
Đầu vào bao gồm:
#' @title Optimze every formula
#' @param param list or vector contains list of parameters
#' @param actual_pd vector
#' @param formula text type
#' @example fnc_optimize(c(1,1,1), c(0.1, 0.3, 0.4, 0.9), 'dgamma(x = x, shape = a, scale = b) * scalar')
fnc_optimize <- function(param, actual_pd, formula){
stopifnot(length(actual_pd) >= 1)
x <- 1:length(actual_pd)
a <- param[1]
b <- param[2]
scalar <- param[3]
fitted_pd <- eval(parse(text = formula))
sse <- sum((actual_pd - fitted_pd)^2)
return(sse)
}Ví dụ sử dụng phân phối gamma để fit
gamma_formula <- 'dgamma(x = x, shape = a, scale = b) * scalar'
pd <- c(0.1, 0.3, 0.4, 0.9)
res <- optim(c(1, 1, 1), fnc_optimize, actual_pd = pd, formula = gamma_formula)
res## $par
## [1] 3.336554 5.127654 46.730584
##
## $value
## [1] 0.02294065
##
## $counts
## function gradient
## 340 NA
##
## $convergence
## [1] 0
##
## $message
## NULLĐầu ra có 2 giá trị quan trọng nhất cần quan tâm:
#' @title Export parameters and pitpd
#' @param period_len length of period will predicted
#' @example fnc_predict_pitpd(period_len = 10, actual_pd = c(0.1, 0.3, 0.4, 0.9), formula = 'dgamma(x=x, shape=a, scale=b)*scalar')
fnc_predict_pitpd <- function(period_len, actual_pd, formula) {
res <- optim(c(1, 1, 1), fnc_optimize, actual_pd = actual_pd, formula = formula)
a <- res$par[1]
b <- res$par[2]
scalar <- res$par[3]
x <- seq(period_len)
pred <- eval(parse(text = formula))
sse <- res$value
return(data.frame(a = a,
b = b,
scalar = scalar,
sse = sse,
period = t(pred)))
}Ví dụ: dự báo cho 10 kỳ
fnc_predict_pitpd(period_len = 10, actual_pd = pd, formula = gamma_formula) %>%
t() %>%
knitr::kable()| a | 3.3365536 |
| b | 5.1276544 |
| scalar | 46.7305841 |
| sse | 0.0229406 |
| period.1 | 0.0590203 |
| period.2 | 0.2452887 |
| period.3 | 0.5205083 |
| period.4 | 0.8387976 |
| period.5 | 1.1625085 |
| period.6 | 1.4645723 |
| period.7 | 1.7275836 |
| period.8 | 1.9419759 |
| period.9 | 2.1041075 |
| period.10 | 2.2145538 |
Viết vòng lặp cho 5 loại phân phối và danh sách các giá trị pd
#' @title Scale pit pd
#' @param period_len length of period will predicted
#' @param input_pd data.frame with columns are actual_pd
#' @param input_formula list formual to optimize
#' @example fnc_predict_pitpd_five_formula(20, input_pd = cohort_fit_distribution)
fnc_predict_pitpd_five_formula <-
function(period_len, input_pd) {
# list formula to fit
list_formula <-
list(
gamma = 'dgamma(x = x, shape = a, scale = b) * scalar',
cauchy = 'scalar/(pi*a*(1+((x-b)/a)^2))',
exponential = 'a*exp(-a*x)*scalar',
loglogit = 'scalar*(a/b)*(x/b)^(a-1)/(1+(x/b)^a)^2',
lognorm = 'scalar*dlnorm(x, a, b)'
)
# list object to loop (same name with ...)
ll <- cross(list(actual_pd = input_pd, formula = list_formula)) %>%
transpose()
# id name
cname <- cross2(names(input_pd), names(list_formula)) %>%
transpose() %>%
pmap_dfr(data.frame) %>%
set_names(c('Segment', 'Distribution'))
pred <- ll %>%
pmap_dfr(fnc_predict_pitpd,
period_len = period_len,
.id = 'id')
out <- cname %>% bind_cols(pred)
return(out)
}Thực hành với dữ liệu: Đầu vào là dữ liệu PD thực tế theo năm của từng danh mục
head(df)## Autoloan_G1 CC_G1 Mortage_G5
## y1 0.015753652 0.015340336 0.618131868
## y2 0.023808917 0.017209285 0.049450549
## y3 0.010171047 0.009240085 0.008241758
## y4 0.003696245 0.003850036 0.000000000tictoc::tic()
output <- fnc_predict_pitpd_five_formula(10, input_pd = df)
tictoc::toc()## 0.23 sec elapsedKết quả
output %>%
t() %>%
knitr::kable()| Segment | Autoloan_G1 | CC_G1 | Mortage_G5 | Autoloan_G1 | CC_G1 | Mortage_G5 | Autoloan_G1 | CC_G1 | Mortage_G5 | Autoloan_G1 | CC_G1 | Mortage_G5 | Autoloan_G1 | CC_G1 | Mortage_G5 |
| Distribution | gamma | gamma | gamma | cauchy | cauchy | cauchy | exponential | exponential | exponential | loglogit | loglogit | loglogit | lognorm | lognorm | lognorm |
| id | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 | 14 | 15 |
| a | 5.1691860 | 3.5074529 | 1.0181969 | 0.9989865 | 1.3060123 | 0.1539415 | 0.2953774 | 0.3207815 | 2.5118029 | 3.3649808 | 2.6932034 | 2.9437420 | 0.6762375 | 0.7073124 | -0.9938052 |
| b | 0.4021953 | 0.6149202 | 0.3953014 | 1.7732236 | 1.6207018 | 0.6712425 | 1.5389001 | 5.0228347 | -1.2275540 | 1.9819620 | 2.0473477 | 0.4797391 | 0.4547107 | 0.5541649 | 0.7126139 |
| scalar | 0.05318037 | 0.04755569 | 2.98382250 | 0.07869961 | 0.07693377 | 1.66241104 | 0.09094084 | 0.07534348 | 3.03352871 | 0.05655741 | 0.05142470 | 2.26929566 | 0.05426241 | 0.04807505 | 2.91981266 |
| sse | 8.391286e-07 | 1.634227e-08 | 1.826390e-05 | 2.916338e-07 | 3.918079e-07 | 1.142901e-04 | 1.191815e-04 | 3.305539e-05 | 1.801130e-05 | 8.325210e-09 | 2.138549e-07 | 1.448987e-05 | 5.165187e-08 | 1.253483e-07 | 3.271995e-06 |
| period.1 | 0.01579231 | 0.01536379 | 0.61799129 | 0.01568161 | 0.01529583 | 0.61815369 | 0.01999196 | 0.01753644 | 0.61811820 | 0.01573678 | 0.01533181 | 0.61822726 | 0.01575460 | 0.01532635 | 0.61813640 |
| period.2 | 0.02364209 | 0.01718127 | 0.04986764 | 0.02384738 | 0.01729226 | 0.04552631 | 0.01487903 | 0.01272412 | 0.05014289 | 0.02378380 | 0.01729497 | 0.04849238 | 0.02378722 | 0.01729892 | 0.04960103 |
| period.3 | 0.010666775 | 0.009339730 | 0.004002966 | 0.009998375 | 0.008864046 | 0.014955540 | 0.011073729 | 0.009232389 | 0.004067684 | 0.010097909 | 0.008954495 | 0.010004152 | 0.010308401 | 0.008990890 | 0.007313613 |
| period.4 | 0.0029452084 | 0.0037787730 | 0.0003206373 | 0.0042013675 | 0.0043414993 | 0.0073358789 | 0.0082416301 | 0.0066988532 | 0.0003299781 | 0.0037416372 | 0.0042034161 | 0.0032336745 | 0.0035164848 | 0.0040845999 | 0.0015452759 |
| period.5 | 6.213581e-04 | 1.300429e-03 | 2.565288e-05 | 2.193282e-03 | 2.436713e-03 | 4.341791e-03 | 6.133839e-03 | 4.860566e-03 | 2.676843e-05 | 1.550390e-03 | 2.103884e-03 | 1.343755e-03 | 1.159040e-03 | 1.839786e-03 | 4.136145e-04 |
| period.6 | 1.105741e-04 | 4.039878e-04 | 2.050859e-06 | 1.326653e-03 | 1.531448e-03 | 2.866356e-03 | 4.565114e-03 | 3.526738e-03 | 2.171505e-06 | 7.276591e-04 | 1.145398e-03 | 6.552589e-04 | 3.914076e-04 | 8.501211e-04 | 1.310061e-04 |
| period.7 | 1.749674e-05 | 1.169426e-04 | 1.638750e-07 | 8.837564e-04 | 1.043734e-03 | 2.032593e-03 | 3.397589e-03 | 2.558937e-03 | 1.761565e-07 | 3.784464e-04 | 6.718936e-04 | 3.569264e-04 | 1.378907e-04 | 4.067475e-04 | 4.709306e-05 |
| period.8 | 2.540542e-06 | 3.214572e-05 | 1.308961e-08 | 6.292437e-04 | 7.542879e-04 | 1.515975e-03 | 2.528658e-03 | 1.856718e-03 | 1.429014e-08 | 2.134621e-04 | 4.191880e-04 | 2.108526e-04 | 5.089824e-05 | 2.017579e-04 | 1.869091e-05 |
| period.9 | 3.454465e-07 | 8.494233e-06 | 1.045240e-09 | 4.701892e-04 | 5.694945e-04 | 1.173912e-03 | 1.881955e-03 | 1.347201e-03 | 1.159243e-09 | 1.284158e-04 | 2.750265e-04 | 1.325290e-04 | 1.967118e-05 | 1.035899e-04 | 8.034852e-06 |
| period.10 | 4.460078e-08 | 2.175692e-06 | 8.344635e-11 | 3.643895e-04 | 4.447080e-04 | 9.357901e-04 | 1.400646e-03 | 9.775043e-04 | 9.403993e-11 | 8.137173e-05 | 1.880606e-04 | 8.747754e-05 | 7.939844e-06 | 5.491799e-05 | 3.689290e-06 |