The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Introduction to FARS

Gian Pietro Bellocca

Introduction

The FARS package provides a comprehensive framework for modeling and forecasting economic scenarios based on the Multilevel Dynamic Factor Model (MLDFM).

Installation

# Install FARS from CRAN
#install.packages("FARS")

# Install FARS from GitHub
#install.packages("devtools")
#devtools::install_github("GPEBellocca/FARS")

library(FARS)
## 
## Attaching package: 'FARS'
## The following object is masked from 'package:stats':
## 
##     density

Input data

library(zoo)
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(readxl)
# Input dep variable
data_input_path <- system.file("extdata", "Data_IMF.xlsx", package = "FARS")
data_input <- read_excel(data_input_path, sheet = "data")
data_ts <- ts(data_input[, -1], frequency = 4)
data_diff <- diff(log(data_ts)) * 400
dep_variable <- data_diff[, 'United States', drop = FALSE] 
dep_variable <- dep_variable[2:60]

# Input data
data_df_path <- system.file("extdata", "DataBase.xlsx", package = "FARS")
data_df <- openxlsx::read.xlsx(data_df_path, sheet = "fulldata", cols = 2:625)
data_df <- data_df[,1:519]
data <- as.matrix(data_df)
dimnames(data) <- NULL

# Generate dates
quarters <- as.yearqtr(seq(from = as.yearqtr("2005 Q2"), by = 0.25, length.out = 59))
dates <- as.Date(quarters)

MLDFM Model

This section demonstrates how to use the mldfm function to fit the Multilevel Dynamic Factor Model using the input data.

# MULTI-LEVEL DYNAMIC FACTOR MODEL
# 3 blocks
n_blocks <- 3 # 63 248 208
block_ind <- c(63,311,519)
r <- c(1,0,1,0,1,1,1)

mldfm_result <- mldfm(data, 
                      r=r, 
                      blocks = n_blocks, 
                      block_ind = block_ind , 
                      tol = 1e-6, 
                      max_iter = 1000,
                      method = 0) 
## Summary of Multilevel Dynamic Factor Model (MLDFM)
## ===================================================
## Number of observations:  59 
## Total number of factors:  5 
## Number of nodes:  5 
## Initialization method:  CCA 
## Number of iterations to converge:  46 
## 
## Factors per node:
##  - 1-2-3 :  1 factor(s)
##  - 1-3 :  1 factor(s)
##  - 1 :  1 factor(s)
##  - 2 :  1 factor(s)
##  - 3 :  1 factor(s)
## 
## Residual sum of squares (RSS):  15464.0395 
## Average RSS per observation:  262.1024
# Plot factors
# plot(mldfm_result, dates = dates)

Subsampling

# Subsampling
n_samples <- 100
sample_size <- 0.9
mldfm_subsampling_result <- mldfm_subsampling(data, 
                                         r=r, 
                                         blocks = n_blocks, 
                                         block_ind = block_ind , 
                                         tol = 1e-6, 
                                         max_iter = 1000, 
                                         method = 0,
                                         n_samples = n_samples,
                                         sample_size = sample_size,
                                         seed = 123) 
## Generating 100 subsamples...
## 
## Subsampling completed.
## Number of subsamples generated: 100

Creating stressed scenarios

# Create stressed scenario
scenario <- create_scenario(model = mldfm_result,
                                  subsample = mldfm_subsampling_result,
                                  data = data, 
                                  block_ind = block_ind, 
                                  alpha=0.95)
## Constructing scenario using 100 subsamples and alpha = 0.95...
## Scenario construction completed.

FARS computation

# Compute Quantiles
fars_result <- compute_fars(dep_variable, 
                              mldfm_result$Factors, 
                              scenario = scenario, 
                              h = 1,   
                              edge = 0.05, 
                              min = TRUE) 
## Running Factor-Augmented Quantile Regressions (FARS)...
## Factor-Augmented Quantile Regressions (FARS)
## ===========================================
## Forecasted quantiles:
##  - Time periods:  59 
##  - Quantile levels:  0.05 0.25 0.50 0.75 0.95 
## 
## Stressed scenario quantiles: YES
# Plot quantiles
plot(fars_result,dates=dates)

Density Estimation

# Density 
density <- nl_density(fars_result$Quantiles,  
                             levels = fars_result$Levels,  
                             est_points = 512, 
                             random_samples = 100000,
                             seed = 123)
## Estimating skew-t densities from forecasted quantiles...
## FARS Density
## ====================
## Time observations  : 59 
## Estimation points  : 512 
## Random samples     : 100000 
## Range of x values  : [ -30 , 10 ]
# Plot density
# plot(density, time_index = dates)

Growth-at-Risk (GaR)

#GaR
GaR0.05 <- quantile_risk(density, QTAU = 0.05)
GaR0.50 <- quantile_risk(density, QTAU = 0.50)
GaR0.95 <- quantile_risk(density, QTAU = 0.95)

Scenario Density Estimation

# Scenario Density 
scenario_density <- nl_density(fars_result$Scenario_Quantiles,  
                             levels = fars_result$Levels,  
                             est_points = 512, 
                             random_samples = 100000,
                             seed = 123)
## Estimating skew-t densities from forecasted quantiles...
## FARS Density
## ====================
## Time observations  : 59 
## Estimation points  : 512 
## Random samples     : 100000 
## Range of x values  : [ -30 , 10 ]

Growth-in-Stress (GiS)

#GiS
GiS0.05 <- quantile_risk(scenario_density, QTAU = 0.05)
GiS0.25 <- quantile_risk(scenario_density, QTAU = 0.25)
GiS0.50 <- quantile_risk(scenario_density, QTAU = 0.50)
GiS0.75 <- quantile_risk(scenario_density, QTAU = 0.75)
GiS0.95 <- quantile_risk(scenario_density, QTAU = 0.95)

Final GaR and GiS Plot

# Final GaR and GiS plot
library(ggplot2)

time <- dates[-1]

MLGaRGiS <- data.frame(
  time = time,
  dep_variable = as.vector(dep_variable[2:59]),
  GaR.05 = as.vector(GaR0.05[1:58]),
  GaR.95 = as.vector(GaR0.95[1:58]),
  GiS.05 = as.vector(GiS0.05[1:58]), 
  GiS.25 = as.vector(GiS0.25[1:58]),
  GiS.75 = as.vector(GiS0.75[1:58]),
  GiS.95 = as.vector(GiS0.95[1:58])
)


p <- ggplot(MLGaRGiS,aes(x=time,y=dep_variable)) + 
  geom_line() + theme_bw()+
  geom_ribbon(aes(ymin=GiS.05,ymax=GiS.95), alpha=0.1, fill="red") +
  geom_ribbon(aes(ymin=GiS.25,ymax=GiS.75), alpha=0.1, fill="grey10") +
  geom_line(aes(x=time, GaR.05), linewidth=0.1, colour="black", linetype="dashed") +
  geom_line(aes(x=time, GaR.95), linewidth=0.1, colour="black", linetype="dashed") +
  scale_y_continuous("Growth") +
  scale_x_date(date_labels = "%Y", date_breaks = "2 years") +
  theme(axis.text.x = element_text(angle = 90))
fig <- plotly::ggplotly(p)
fig

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.