# CRIMS2S

This page describes our submission for the 2021 S2S AI Challenge. The first section, **Submission information**, describes what data was used and our safeguards to prevent overfitting and improve reproductibility. The second section, **Methods**, describes our methods in more detail.

## Submission information

### Data used

* Multiple fields from the ECMWF hindcasts/forecasts. See the `crims2s/conf/fields/both_easy.yaml` file for a complete list.
* `t2m` and `tp` from NCEP and ECCC, as provided on Renku by the challenge organizers.
* `orog` for ECMWF, which was not downloaded through Renku, but we provide the file in our repo.

### Safeguards to prevent overfitting

We took the following measures to prevent overfitting.

* [x] We did not use 2020 observations in training (explicit overfitting and cheating)
* We also remove any training data for which the validity time is in 2020 (i.e. the last part of the late 2019 forecasts).
* [x] We did not repeatedly verify our model on 2020 observations and incrementally improved my RPSS (implicit overfitting)
* The train/validation data and the test data lived in two separate directories for the duration of our development.
* [x] We provide RPSS scores for the training period with script skill_by_year, as provided in the template submission.
* [x] We tried our best to prevent data leakage.
* [x] We honor the train-validate-test split principle. This means that the hindcast data is split into train and validate, whereas test is withheld.
- We use the train-validate split for all but our final experiments. In the final push to improve our results, we use the validation set for training. That is, all hyperparamater tuning and model building for calibrating the individual dynamical model forecast was done using the training data only, which was then validated using the validation data. After this is done, and in our final experiment, we include the validation data for training in order to increase sample size. The test set is in a separate folder and was never loaded for training purposes.
* [x] We did not use test explicitly in training or implicitly in incrementally adjusting parameters.
* [x] We considered cross-validation.
- Our validation strategy was not a cross validation, because the amount of time required to train our models make it impractical. Since we use data from all three centers (ECMWF, ECCC and NCEP), we used years 2010, 2017 and 2019 as validation years. That way, we have at least one validation year where the data is available for each center.

### Safeguards for Reproducibility
Notebook/code must be independently reproducible from scratch by the organizers (after the competition), if not possible: no prize
* [x] All training data is publicly available (no pre-trained private neural networks, as they are not reproducible for us)
- All the data was downloaded from renku, except the orography which we provide in the repo.
* [x] Code is well documented, readable and reproducible.
* Full instructions to execute the code are provided. Reproductibility hiccups are possible, as in any experimental code. We remain available if any support is required.
* Some parts of the source code are on the configurable end of the readability/configurability spectrum. This was done to maximize our ability to improve the performance metrics. We provide pointers to the most important parts of the source code in our **Methods** section. We remain available if more explanations are necessary.
* [x] Code to reproduce training and predictions is preferred to run within a day on the described architecture. If the training takes longer than a day, please justify why this is needed. Please do not submit training piplelines, which take weeks to train.
- Our models comprises deep neural architecture which typically train over many days depending on the GPU power used. In our case, the model training lasts about two days on a GTX 1080 ti. These training times are not out of the ordinary for a convolutional neural network in our experience.
- The data preparation step lasts about an hour for us, but is done using about 53 jobs in parallel on a computing cluster. That means that it could take more than 24 hours if only consumer-grade equipment is used. The usage of so many data is justified by a common principle in deep machine learning where a larger data volume typically improve model performance. We can make our version of the pre-processed data available upon request.

## Methods
The goal of our contributing model is to make global, 3-category, probabilistic forecasts of temperature and precipitation at the weeks 3-4 and 5-6 lead times. For this, we have developed a so-called "opportunistic mixture model" that blends together five different predictions. The predictions are based on dynamical ensemble forecasts from three global circulation models, as well as the observed climatology. Generally speaking, the blending of the different predictions is according to historical skill, hence the name "opportunistic". Since the observed climatology is used as one of the prediction models, our opportunistic model can be interpreted as a form of Bayesian Model Averaging (Raftery et al., 1997). The model architecture is displayed in the following schematic:


Briefly, the predictions include:

1. Climatology -- i.e. the observed climatological probability for a given category (1/3 by definition);
2. Dynamical ensemble forecasts obtained from models from three forecasting centres (ECCC, NCEP, and ECMWF) that have been calibrated using ensemble model output statistics (EMOS; Gneiting et al., 2007);
3. Dynamical ensemble forecasts from ECMWF that have been calibrated using a convolutional neural network (CNN).

These predictions are then combined using an additional CNN that weights each prediction based on historical skill. We now describe the main steps of our opportunistic model, which includes a detailed description of each prediction method, including implementation details and references to the source code.

### EMOS Model

Our model uses an ensemble of EMOS corrected forecasts as its basis. For its parameters $\mu$ and $\sigma$ (where $\mu$ is the ensemble mean and $\sigma$ is the sample standard deviation), our EMOS model uses a rolling window inside a model bank.


The rolling window has the effect of regularizing the EMOS parameters and avoiding overfitting. Indeed, the parameters for a specific week also have to be good for other weeks neighboring weeks of the year. Consequently, they cannot overfit on the data of a single week. We use a rolling window of 20 weeks for our final submission.

Unfortunately, the source code of this EMOS model is heavily factored and perhaps less easily understood. This was done to maximize the configurability of it during the competition. This way, it was possible to try more variants in less time, and optimize the metrics more efficiently. The downside is that the source code could probabily benefit from a rewrite using a fixed configuration, to make it more easily readable.

Nevertheless, we point to important parts of its source code.
The rolling window logic is implemented in the `PytorchRolllingWindowMultiplexer` and `WeeklyRollingWindowMultiplexer` in the `crims2s/training/model/util.py` file. The linear model logic is implemented in `crims2s/training/model/emos.py` file in the `LinearModel` class. The `MultiplexedNormalEMOSModel` and `RollingWindowNormalEMOSModel` classes interface the rolling window logic with the rest of our training framework. `TempPrecipEMOS` joins two EMOS models together so that we can forecast temperature and precipitation at the same time.

### Convolutional post-processing

One of the 5 predictions we mix is the ECMWF forecast post-processed using a convolutional neural network. The convolutional architecture that performs the post-processing is similar to the one of the weight model described below, in that it uses ConvBlocks and skip connections. It is different in that it outputs a 3-category probabilistic forecast instead of relative weights. We refer to the schematics of the weight model for an overview of the architecture.

For a more precise description, the source code of the convolutional post-processing model is found in the `ConvPostProcessingJG` class of the `crims2s/training/model/conv.py` file. Its configuration file is `crims2s/training/conf/experiment/model/conv_JG.yaml`.

### Mixing model

The mixing model brings all our sub-models together to produce forecasts. The schematic below shows its architecture. The "forward pass" and "backpropagation" arrows show how we go from model input (in blue) to the training loss (in red), and back. The backpropagation allows the "Conv. weight model" (upper middle box) to receive information from all 5 predictions (5 lower left blue boxes).


The source code of the mixing model is in the `MultiCenterBayesianUpdateModel` class in the`crims2s/training/model/bayes.py` file. Notably, the weighted averaging is done using `einsum` at line 255:

# Meaning of the keys: Batch, Model (prior or ecmwf or eccc), Category, Lead time, lAtitude, lOngitude
t2m = torch.einsum("bmclao,bmlao->bclao", t2m_forecasts, t2m_weights)
tp = torch.einsum("bmclao,bmlao->bclao", tp_forecasts, tp_weights)

The source code of the training logic is in the `S2STercilesModule` class.
It is located in the `crims2s/training/lightning.py` file. It uses the Pytoch Lightning interface to define training procedures in PyTorch. The actual training is launched using PyTorch Lightning at line 206 of the `crims2s/training/train.py` file.

The configuration file for the training module is `crims2s/training/conf/experiment/module/bayes.yaml`. The configuration file for the overall experiment is `crims2s/training/conf/experiment/bayes_multi_noval_jg.yaml`.

#### Weight model

The weight model has the task of outputting a relative weight for each model, for each variable (temperature and precipitation), for every initialization time, for every lead time, and at each grid cell. It does so using a convolutional architecture.

In the core of the convolutional architecture are two branches. The local branch assigns weights to every grid cell. The global branch captures some of the possibly non-local and larger scale features. It summarizes these features in a 64 dimensions vector. We postulate this vector could summarize important teleconnections such as the Madden-Julian Oscillation. This global context is then re-added to the local representations of the tiles. Then, two prediction heads go from this representation to a set of weights for each variable.

This schematic shows the architecture of the weight model, and then shows the architecture of its various components.


The weight model is defined in the `BiheadedWeightModel` class in the `crims2s/training/model/bayes.py` file. It is configured in `crims2s/training/conf/model/bayes_multi_jg_noval.yaml`, under the `weight_model` key.

## Github

This was an overview of the key aspects of our submission for the S2S AI challenge.
To take this overview further and execute the code, instructions are available in the [repository README](https://github.com/crim-ca/crims2s#reproducing).