1.2.1. Aim of the analysis

This paper extends on the initial stages of our time series analysis on GSS (Green, Social and Sustainability) bonds in Dey (Reference Dey2024), specifically focussing on the daily values from the Standard and Poor’s (S&P) Green Bond Index and whether or not we can create accurate prediction models using neural networks. This paper builds on the initial paper, where we hope to develop a model that can assist with GSS bond index prediction, which will have wider applications such as index price modelling and investment portfolio analyses for actuaries and non-actuaries alike. For the purposes of this paper, we continue to look at predicting a rolling 1-day value of the index, based on the prior x days index value over the period 2013 to 2023 inclusive.

We retain the same Baseline model i.e. today’s value equals yesterday’s value over the course of the full date range of 31 January 2013 to 17 February 2023. Similarly to the initial paper, we aim to see if we can accurately create a time series model with non-traditional methods such as neural networks and more complex neural network architectures, using N-BEATS. Analyses using traditional stationarity techniques such as ARIMA and multivariate techniques have been deferred to later papers.

Towards the end of this paper, we extend our analysis by looking at using the prior x index values (window) to predict the next y days in the future (horizon), rerunning the initial neural network architectures analysed in Dey (Reference Dey2024). We restrict the input window to 2 and 5 days, whilst retaining a future horizon to 1 day in this paper. In Section 3, for N-BEATS, we train the model with an input window of 1 day in addition.

1.2.2. Data and method

We have retained the same underlying data as per Dey (Reference Dey2024) – the S&P Green Bond Index values between 31 January 2013 and 17 February 2023, splitting the data using 70%/20%/10% ratios for training/validation/test data sets. Section 3 of this paper extends our analysis by introducing a more complex model architecture, the N-BEATS model architecture.

The latter part of this paper, Section 4, revisits the models analysed in our initial paper but extends the input window of data to 2 days and 5 days respectively. These models can be categorised into the following model categories: Deep Neural Network (DNN), Convolutional Neural Network (CNN), Long Short-Term Memory (LSTM) and Gated Recurrent Unit (GRU) neural network architectures. This is discussed in Section 4 and Appendix 3 in more detail.

The loss function used during training the models was set to Mean Absolute Error (MAE) for all models, with an Adam optimiser used to update the weights to minimise the loss function as part of the training process and L2 regularisation to reduce any overfitting during the training process. Hyperparameter tuning of all models was completed via the open-source library Optuna, using the Bayesian optimisation algorithm Tree-structure Parzen Estimator (TPE). These techniques are discussed in more detail in Dey (Reference Dey2024).

1.2.3. Results and conclusions

The results of our analysis for N-BEATS (Section 3) were inconclusive: the models from each category produced comparable results to the Baseline model with differences in Mean Absolute Percentage Error (MAPE) of up to approximately +0.04% and hence with no material outperformance.

Using similar neural network architectures as per the initial paper (Section 4), if we extend the input window to 2 and 5 days respectively, there was no improvement in model performance when compared with the Baseline, with differences in MAPE of up to around +0.3%. Widening the window from 2 days to 5 days did result in an improvement in some models, though again the overall model performance was still worse than the Baseline. Hence, the overall results were inconclusive.

We are potentially not sharing sufficient correlated information e.g. as per a multivariate analysis, for our models to learn underlying material information and patterns in the data to result in a model that materially outperforms the Baseline. We aim to address this in future papers (please see below and Section 5 for more details).

1.2.4. Next steps

For future papers, we will expand our analysis to include the following:

1. 1. Expanding the analysis to general GSS bonds. The analysis in this paper is based on a single green bond index. We will look to expand our analysis to the wider GSS bond universe and over differing date ranges for the data to see if there are general underlying conclusions across different data sets and GSS bond indices.

2. 2. Expanding the analysis to include any potential relationships with the general market such as stock market and oil prices i.e. move to a multivariate analysis in subsequent papers. This is examined e.g. in Wang et al. (Reference Wang, Tang and Guo2021) which, when coupled with the CEEMDAN-LSTM model, seems to produce materially improved model predictions based on green bond time series data when compared to a baseline model.

3. 3. Extend our analysis to alternative neural network architecture types, such as Graph Neural Networks (GNNs).

1.2.5. Additional disclaimers

Please note the following:

1. a. Information within this paper is valid up to 31 May 2024. Hence, there may be updates beyond this date that are not reflected in this paper e.g. changes to any legislation mentioned or updates to any open-source libraries used.

2. b. This paper is not intended to be a comprehensive audit of models. Neither is this paper recommending or promoting one approach over another, nor promoting any of the sources or references stated in this paper. Any user of this paper should still reference the underlying legislation, reference any standard mentioned in this paper, and should there be any conflict, the underlying information in the relevant standard, reference or legislation supersedes any information presented in this paper.

3. c. Though the work in this paper does not fall under the Financial Reporting Council’s Technical Actuarial Standards, this paper has been reviewed both within the Working Party and by the Institute and Faculty of Actuaries’ Data Science Practice Board.

3.2.1. Basic building block

The basic building blocks generate a backcast and forecast as previously mentioned. The backcast can be viewed as iterative adjustments to the input vectors via residuals, whilst the forecast from each basic building block can be viewed as a partial forecast which is cumulated and form part of the overall global forecast.

Focusing on the l ^th basic building block, as per Oreshkin et al. (Reference Oreshkin, Carpov, Chapados and Bengio2019): each basic building block has a fork architecture which takes an input vector (or residual outputs from a prior basic building block) ${{\boldsymbol x}_l}$ and produces two outputs: $\widehat {{{\boldsymbol x}_l}}$ i.e. the block’s best estimate of ${{\boldsymbol x}_l}$ and the block’s forward forecast $\widehat {{{\boldsymbol y}_l}}$ . The input of the very first block in the model is the original data input. The remaining basic building block inputs are residual outputs from a prior basic building block, and can be thought of as running a sequential analysis of the input signal (Oreshkin et al., Reference Oreshkin, Carpov, Chapados and Bengio2019). See Section 3.2.2 for more details.

The internal structure of a basic building block broadly consists of two sections: a fully connected section and a basis layer section.

The fully connected layer produces backward (backcast) and forward (forecast) predictors, $\theta _l^b$ and $\theta _l^f$ respectively, of expansion coefficients. The predictors determine how much each basis function contributes to the overall approximation of the function. Please see Section 3.2.1.2 for more details on basis function.

The above two sections of the internal structure of a basic building block are discussed in more detail immediately below in Sections 3.2.1.1 and 3.2.1.2.

3.2.2. Doubly residual stacking

A unique feature of the N-BEATS model architecture is the concept of doubly residual stacking. Instead of a single residual branch of outputs, the architecture proposes two branches as mentioned: one for a backcast residual and one for a partial forecast. Mathematically, we represent this as:

$${{\boldsymbol x}_l} = {{\boldsymbol x}_{l - 1}} - \;\widehat {{{\boldsymbol x}_{l - 1}}}, \ \ \widehat {\boldsymbol y} =\mathop \sum \nolimits_l \widehat {{{\boldsymbol y}_l}},$$

where ${{\boldsymbol x}_l}$ and $\widehat {{{\boldsymbol x}_{l - 1}}}$ are defined earlier in Section 3.2.1, $\widehat {\boldsymbol y}\;$ is the global forecast i.e. overall model output which is the sum of individual $\widehat {{{\boldsymbol y}_l}}\;$ . Each $\widehat {{{\boldsymbol y}_l}}\;$ can be viewed as a partial forecast.

3.2.3. Interpretability

There are 2 types of configurations or design patterns for an N-BEATS architecture:

1. 1. Generic architecture, which does not assume any specific underlying time series model and relies on a generalised deep learning approach.

2. 2. Interpretable architecture, which is an augmented model that often includes a trend and seasonality model to the generalised N-BEATS architecture discussed earlier for a time series forecasting problem.

For the purposes of this paper, we have focussed on the first configuration i.e. generic architecture as the aim of this paper is to seek a model architecture that accurately predicts the time series future values without any further input by the modeller or prior knowledge of any underlying time series features from the modeller. For more details on the above, and other considerations, please see Oreshkin et al. (Reference Oreshkin, Carpov, Chapados and Bengio2019).

3.2.4. Comparison with other model architectures

In Table 2, we provide a high-level comparison of the N-BEATS architecture versus DNNs, and LSTMs/GRUs. Please note that these are generalised comments based on experience and the comparison table below may differ for different data sets, problems and scenarios.

Table 2. Comparison of N-BEATS model architecture.

3.3.1. Code background

The underlying basis of our code for the N-BEATS model is taken from Zero to Mastery TensorFlow for Deep Learning Book (Bourke, Reference Bourke2023) and has been adapted for our analysis, where we look at a general N-BEATS architecture. We vary the number of layers within each block and the number of blocks per stack but set the total number of stacks equal to 1 (i.e. M is 1 in Figure 2). This is in part to simplify the trained model architecture given that the time taken to train such a model was close to 4 hours. Further, in doing so, we do not believe that taking such an approach would materially impact any overall conclusion, as effectively adding another stack would repeat a similar iterative fine-tuning process and potentially the hyperparameter optimisation could compensate by increasing the number of blocks in the single stack.

The original authors of the N-BEATS model also extend their analysis to allow for an ensemble technique (using multiple different loss functions and multiple different lookback periods in Figure 2) to make predictions when testing on M4 dataset. To retain consistency with the approach taken with other models within our initial paper, we have not extended our analysis to allow for such ensemble techniques.

We have assumed the same underlying architecture per block within each stack e.g. the same number of hidden layers within each block. Similarly, we have assumed the same underlying architecture for each hidden layer within each block e.g. number of neurons, activation function and L2 loss regularisation parameter. The output layer from the block is a generalised dense layer, each with the same number of neurons in effect as the total combined length of the input window and output horizon. Similarly, we have assumed that all output layers have the same activation function, though this may differ to the hidden layers above.

Similar to our initial paper, we have used the same underlying data and training/validation/test splits within the data to train our model. We have run our code in Google Colab, using similar libraries including Tensorfow and Keras as per our initial paper. For hyperparameter optimisation, we have again used Optuna. For an in-depth discussion, please see Dey (Reference Dey2024).

3.3.2. Hyperparameter optimisation

We have used Optuna for hyperparameter optimisation. Below we list some key search spaces used within our hyperparameter optimisation.

1. 1. The number of epochs was set to 400.

2. 2. The number of trials within Optuna was set to 30.

3. 3. The batch size for the input data was set to 128.

4. 4. The activation list was as follows: elu, gelu, linear, relu, sigmoid, swish and tanh.

5. 5. The specific search space for each N-BEATS block includes:

   1. a. Number of neurons in each hidden layer from 4 to 512 at steps of 4.

   2. b. Number of hidden layers from 1 to 4 at steps of 1. As mentioned earlier, we have assumed all hidden layers have the same underlying architecture e.g. the same number of neurons, same activation function and same L2 regularisation.

   3. c. L2 regularisation with a search space of 1e-5 to 1e-1.

   4. d. Output theta layer or basis layer is a generalised deep layer.

      1. i. The number of neurons is set equal to the theta size – in effect the combined length is the window size and horizon size for this analysis.

      2. ii. The activation function is from the activation list above. No L2 regularisation function was applied to this layer.

      3. iii. No further adjustments were made for seasonality and trend i.e. via basis functions, as previously mentioned.

6. 6. Specific search space for number of N-BEATS blocks per stack from 2 to 30 at steps of 2.

7. 7. Similar to Dey (Reference Dey2024), the model was compiled on a loss function of Mean Absolute Error (MAE) using an Adam optimiser with a learning rate search space of 1e-5 to 1e-1.

3.3.3. Training the model

Below is a loss history curve based on training one of the models, using the search space, number of epochs and number of trials in Optuna as mentioned in Section 3.3.2 (Figure 3).

As we would expect, the model loss decreases as the number of epochs increases, suggesting that the model is able to converge to a solution.

Similarly, we have used the open source Keras visualizer library to demonstrate the underlying architecture of the trained model in Figure 4. Please note that we have taken a snippet of the first few blocks and initial stacks, given the size of the overall final trained model.

Figure 3. Training and validation loss curves (output from Google Colab).

The visualisation tool in general gives an alternative simplified representation of any final neural network model architecture. Further, given that we have generalised our code and set the overall number of stacks equal to 1 as mentioned earlier, there are cosmetic differences in Figure 4 when compared with Figure 2.

In Figure 4, the Keras visualisation tool helps visualise basic building blocks (represented by “NBeatsBlock”) showing:

1. 1. Subtraction (represented by “Subtract” in Figure 4) to represent the process of residuals or backcasts.

2. 2. Additions (represented by “Add” in Figure 4) to represent the process of partial forecasts.

The final trained model (with an input window of 1 day) has over 4 million trainable parameters (i.e. weights and biases). As a comparison, Model 1 from Appendix 3 in Dey (Reference Dey2024), which is a DNN model, has 625 trainable parameters. Similarly, Model 6 from Appendix 3 in Dey (Reference Dey2024), which is an LSTM model, has 47,629 trainable parameters.

3.3.4. Model outputs and observations

The paper Oreshkin et al. (Reference Oreshkin, Carpov, Chapados and Bengio2019) mentions that the input length of the window is typically set to a multiple of the length of the output horizon H, between 2H and 7H. For the purpose of this paper, we have trained the model setting this input length to a multiple of one (1 times) to align with Dey (Reference Dey2024), two (2 times) and five (5 times) in order to adhere to the recommendation above and in Oreshkin et al. (Reference Oreshkin, Carpov, Chapados and Bengio2019), as well as our analysis in Section 4. The output horizon is as per the initial paper i.e. 1 day.

In Table 3, we compare the performance of the trained N-BEATS models, which have been trained on the same data but on input windows of 1, 2 and 5 days, and with the same output horizon of 1 day.

Table 3. Comparison of performance measure, between the baseline model and N-BEATS, based on MAE and MAPE, varying the input window between 1, 2 and 5 days.

Though the N-BEATS models performed relatively accurately, based on the data set, split of data between training/validation/test, and hyperparameter optimisation, there is no improvement to the Baseline for each of our N-BEATS models based on either MAE or MAPE measures on the test data set. The model trained on a 5-day input window performs marginally better than the model trained on a 2-day input window. However, the model trained on a 1-day input window outperforms the N-BEATS models trained on a 2-day and 5-day input window. All trained models perform marginally worse than the Baseline model. Hence, the overall results are inconclusive.

One potential area of future study is to incorporate ensemble techniques as per the original authors of N-BEATS as mentioned in Section 3.3.1, as well as allow for the number of stacks to vary during hyperparameter optimisation when training the model.