## Proposed Interpolation Model
### Patch-Based NN Architecture  
When considering spatio-temporal fields, the application of the model defined by (16) and (17) should be considered with care to account for the underlying dimensionality, especially for the covariance model. For this reason, a global representation of the spatio-temporal field is most likely to fail due to computational limitations. Following our previous works on analog data assimilation [21; 22], we consider a patch-based representation as sketched in Figure 1 (A patch is a \(P\times P\) subregion of a 2D field with \(P\) the width and the height of the patch). This patch-based representation is fully embeddedin the considered NN architecture to make explicit both the extraction of the patches from a 2D field and the reconstruction of a 2D field from the collection of patches. The latter involves a reconstruction operator which is learnt from data. Regarding model \(\mathcal{F}\), the proposed architecture proceeds as follows:  
* At a given time \(t\), the first layer of the network, which is parameter-free in terms of training, comes to decompose an input field \(x_{t}\) into a collection of \(N_{p}\)\(P\times P\) patches \(x_{\mathcal{P}_{s},t}\), where \(P\) is the width and height of each patch and \(s\) the patch location in the global field. Each patch is decomposed onto an EOF basis \(\mathcal{B}\) according to: \[z_{\mathcal{P}_{s},t}=x_{\mathcal{P}_{s},t}\mathcal{B}^{T}\] (18) with \(z_{\mathcal{P}_{s},t}\) the EOF decomposition of the patch \(x_{\mathcal{P}_{s},t}\). The EOF decomposition matrix \(\mathcal{B}\) is trained offline as preprocessing step; For each \(s\in[1,...,N_{p}]\), we predict \(z_{\mathcal{P}_{s},t+1}\) using an EOF-patch-based model \(\mathcal{F}^{\mathcal{P}_{s}}\). This model is implemented based on a residual architecture to mimic a numerical integration scheme (typically, an Euler or 4 th-order Runge-Kutta scheme) of an approximate Ordinary Differential Equation (ODE) parametrized by the residual block of our residual network. By contrast to other neural networks models, This architecture grantee the physical interpretability of our dynamical model as stated in [27]. In order to enhance the modeling capabilities of our approximate model, The residual block is a classic Multilayer Perceptron (MLP) network with bilinear layers;
* The third layer is a reconstruction network \(\mathcal{F}_{r}\). It combines the predicted patches \(x_{\mathcal{P}_{s},t}=z_{\mathcal{P}_{s},t}\mathcal{B},s\in[1,...,N_{p}]\) to reconstruct the output field \(x_{t}\). This reconstruction network \(\mathcal{F}_{r}\) involves a convolution neural network [38]. The details of the considered parameterizations for the second and third layers are given in Section 4. To train the mean dynamical model \(\mathcal{F}\), we apply a two-step procedure. We first learn the local dynamical models \(\mathcal{F}^{\mathcal{P}_{s}},s\in[1,...,N_{p}]\) based on the minimization of the EOF-patch based forecasting error. The reconstruction network \(\mathcal{F}_{r}\) is then optimized using the same criterion over the global field. This training procedure allows the patch based models to be interpreted as local dynamical models and the reconstruction network as a post-processing operator. Other training configurations could be envisaged, we can for example train the all model according to a forecasting error over the global field. However, this results in inconsistent patch models \(\mathcal{F}^{\mathcal{P}_{s}}\) that can't be used in assimilation experiments for patch reconstruction issues. Regarding the covariance model \(\mathcal{F}_{\Sigma}\), we also consider a patch-based representation of the spatial domain. More precisely, a block-diagonal parameterization of the covariance model \(\mathcal{F}_{\Sigma}\) is addressed by training diagonal patch-level covariance models in the EOF space. It may be noted that a diagonal parameterization of the covariance in the EOF space forms a full covariance matrix in the original patch space.