Shalini Sharma
email: shalinis@iiitd.ac.in.

Émilie Chouzenoux
email: emilie.chouzenoux@inria.fr

Víctor Elvira
email: victor.elvira@ed.ac.uk

Angshul Majumdar
email: angshul@iiitd.ac.in

A dynamical model for stock forecasting via deep recurrent dictionary learning

Keywords: Time series analysis, State space models, Deep nonnegative matrix factorization, Kalman filtering, Bayesian smoothing, EM algorithm, Stock forecasting, Stock trading

State-space models (SSM) and recurrent neural networks (RNN) are widely used approaches for dynamical system modeling. In the case of SSMs, they include explicit modeling of all components, including the noise characterization, and thus allow for interpretability and uncertainty quantification. However, the underlying dynamical model parameters must be specified, and closed-form inference is possible only in a few simple cases. RNNs, on the other hand, can learn, through supervised training, rather complex nonlinearities from the data but lack the aforementioned advantages of SSMs. In this work, we combine the benefits of both approaches by introducing a Gaussian SSM whose state and evolution operators can be learnt from the data. In order to deal with the

ill-posedness of this parameter estimation problem, we propose an innovative factorized form of both the state and observation operators reminiscent of deep nonnegative matrix factorization models. An expectation-maximization method combined with a block alternating strategy is introduced to estimate each of the involved positive latent factors, while jointly performing the probabilistic state inference. Our resulting formulation and inference tool is deep recurrent dictionary learning (DRDL). We then specialize DRDL for the problem of stock forecasting by proposing an online training strategy and a probabilistic assessment of the trading decision. Numerical experiments on a problem of stock market data inference show its superiority among several state-of-the-art dynamic modeling tools.

Introduction

Modeling dynamical systems has been a topic of interest to signal processing, machine learning and control engineering researchers for more than five decades.

Applications range in areas as diverse as financial market analysis to electric demand forecasting. We propose a new dynamical recurrent modeling technique that combines the advantages of state-of-the-art deep learning tools with those of traditional state-space models. The proposed tool is then particularized to processing stock market time series. In the following, we review the literature of dynamic modeling around this particular application and we describe the contributions of the paper.

State-of-the-art review

Modeling the stock market is a well-known challenging problem [START_REF] Yang | 10 challenging problems in data mining research[END_REF]. The difficulty lies in the non-stationary and nonlinearity of the underlying dynamical process. Moreover, financial markets are not only influenced by consumer behavior but also by a myriad of external factors like natural disasters, administrative policies, political decisions, international relations, etc., to name a few. Therefore developing reliable algorithmic models for stock trading still remains a challenging yet interesting topic from the point of view of both finance and machine learning/signal processing [START_REF] Fama | Efficient capital markets a review of theory and empirical work[END_REF]; [START_REF] Shah | Stock market analysis: A review and taxonomy of prediction techniques[END_REF]. Auto-regressive moving average (ARMA) models have been used to model stock market Zhao-yang (2010); [START_REF] Atsalakis | Surveying stock market forecasting techniques-part i: Conventional methods[END_REF]; [START_REF] Moon | Ar and arma model order selection for time-series modeling with imagenet classification[END_REF]; [START_REF] Nowicka-Zagrajek | Modeling electricity loads in california: Arma models with hyperbolic noise[END_REF]. ARMA assumes the stochastic process to be stationary; this turns out to be too simplistic and consequently unrealistic for the stock market. This limitation was partially overcome by autoregressive integrated moving average (ARIMA) [START_REF] Ariyo | Stock price prediction using the arima model[END_REF] models (also referred as Box-Jenkins model). ARIMA has been used in the past for stock forecasting and trading [START_REF] Devi | An effective time series analysis for stock trend prediction using arima model for nifty midcap-50[END_REF], [START_REF] Petrica | Limitation of arima models in financial and monetary economics[END_REF]. However, Box-Jenkins/ARIMA methods could not model non-smooth variations in time series [START_REF] Makridakis | Arma models and the box-jenkins methodology[END_REF]; [START_REF] O'donovan | Short term forecasting: An introduction to the boxjenkins approach[END_REF]. ARIMA with regressors were introduced to overcome the limitations, leading to ARIMAX Chadsuthi et al. (2012); [START_REF] Ababio | Comparative study of stock price forecasting using arima and arimax models[END_REF]. Unfortunately, ARIMAX introduced other problems such as over/under-fitting because of the handling of extra predictors and variables.

SSM is another powerful approach for modeling and analysing time-series [START_REF] Särkkä | Bayesian filtering and smoothing[END_REF]; [START_REF] Newman | State-space models for ecological time-series data: Practical model-fitting[END_REF]. Many studies used SSMs for stock forecasting, and analysis Östermark (1991); [START_REF] Saini | Forecasting volatility in indian stock market using state space models[END_REF]; [START_REF] Elvira | Graphical inference in linear-Gaussian state-space models[END_REF]. The celebrated Kalman filter is a solution to the inference of a linear SSM where the noise is assumed to be Gaussian [START_REF] Kalman | A new approach to linear filtering and prediction problems[END_REF]. The literature illustrates its minimal use in stock forecasting [START_REF] Choudhry | Forecasting the weekly time-varying beta of uk firms: Garch models vs. kalman filter method[END_REF] but found its application in other financial analyses, see for example [START_REF] Wells | The Kalman filter in finance[END_REF]. To overcome the restrictive linearity assumption, extended Kalman filter (EKF) [START_REF] Ljung | Asymptotic behavior of the extended kalman filter as a parameter estimator for linear systems[END_REF], unscented Kalman filter (UKF) [START_REF] Wan | The unscented kalman filter for nonlinear estimation[END_REF] were introduced. Particle filters [START_REF] Djuric | Particle filtering[END_REF]; [START_REF] Doucet | A tutorial on particle filtering and smoothing: Fifteen years later[END_REF]; [START_REF] Elvira | Elucidating the auxiliary particle filter via multiple importance sampling [lecture notes[END_REF]; [START_REF] Ntemi | A jump-diffusion particle filter for price prediction[END_REF] further relaxed the Gaussianity assumption. The advantage of the SSM approach is that it can model uncertainty in the estimate [START_REF] Doucet | A tutorial on particle filtering and smoothing: Fifteen years later[END_REF]; [START_REF] Bach | Learning graphical models for stationary time series[END_REF]. Uncertainty is crucial for financial markets since it gives a measure of the associated risk [START_REF] Rigotti | Uncertainty and risk in financial markets[END_REF].

The main drawback of the aforesaid signal processing-based forecasting approaches is that they need the model's specification. Unfortunately, specifying an underlying model for the stock market is difficult, if not impossible. Several works in the literature have thus investigated the learning of model parameters in SSMs. In the case of linear-Gaussian state-space models (LG-SSMs), see for instance the methods in [START_REF] Sharma | Blind kalman filtering for short-term load forecasting[END_REF]; [START_REF] Shumway | An approach to time series smoothing and forecasting using the em algorithm[END_REF]; [START_REF] Khan | An expectation-maximization algorithm based Kalman smoother approach for event-related desynchronization (ERD) estimation from EEG[END_REF], and (Särkkä, 2013, Chapter 12). All these works consider the observation and state operators to be unknown and estimated from data using expectation-maximization (EM) methods. However, the aforementioned works can only consider linear models and do not account for any prior knowledge on the involved operators. Inferring model parameters for non-linear SSMs has been explored in more generic algorithms, e.g., particle MCMC methods [START_REF] Andrieu | Particle markov chain monte carlo methods[END_REF], SMC 2 Chopin et al. (2013), and nested PFs [START_REF] Crisan | Nested particle filters for online parameter estimation in discrete-time state-space markov models[END_REF]. In all these cases, the inference is costly, as they use Monte-Carlo sampling methods and thus do not generally scale well. The problem of scalability in SSM model inference has been mildly explored. Let us mention our two recent works Chouzenoux & Elvira (2020); [START_REF] Sharma | Recurrent dictionary learning for state-space models with an application in stock forecasting[END_REF], both focusing on LG-SSMs. In Chouzenoux & Elvira (2020), a sparsity prior in introduced on the linear matrices to infer, providing an interpretable and compressible model. Though this method is promising, it does not allow explicit control of the final dimension of the model easily, and as such, still requires an increased computational time at inference. In [START_REF] Sharma | Recurrent dictionary learning for state-space models with an application in stock forecasting[END_REF], we proposed an online (still EM-based) estimation approach in the context of stock market time series processing. The online processing allows a reduced complexity and memory burden while being beneficial to capturing non-linear phenomena in such volatile time series. However, the parametric estimation step lacked of robustness and lack of sufficient imposed structure on the estimated factors.

Neural network (NN) models represent another family of approaches for time series modeling. By construction, these methods excel when the model specification is missing, as they implicitly learn the model from the data through the training phase. In particular, the approximation capability of recurrent neural networks (RNNs) for dynamical systems allows to learn the underlying phenomena given enough training data Gonzalez-Olvera & Tang (2010); [START_REF] Won | Identification of finite state automata with a class of recurrent neural networks[END_REF]; [START_REF] Yin | Graph-based stock correlation and prediction for high-frequency trading systems[END_REF]. RNN and its subsequent versions are used in several studies for stock price forecasting [START_REF] Saad | Comparative study of stock trend prediction using time delay, recurrent and probabilistic neural networks[END_REF]; [START_REF] Tino | Financial volatility trading using recurrent neural networks[END_REF].

Deeper neural network architectures are known to yield better results than their shallow counterparts [START_REF] Bengio | Greedy layerwise training of deep networks[END_REF]; [START_REF] Shao | Learning deep and wide: A spectral method for learning deep networks[END_REF]. They are engineered to approximate highly non-linear function in high-dimensional spaces and are supposed to be more suitable for challenging problems [START_REF] Cheridito | Efficient approximation of high-dimensional functions with neural networks[END_REF]; [START_REF] Bianchini | On the complexity of neural network classifiers: A comparison between shallow and deep architectures[END_REF]. 1-D CNN performs better when compared to LSTM and RNN owing to their ease of training. There are studies, such as [START_REF] Sezer | Algorithmic financial trading with deep convolutional neural networks: Time series to image conversion approach[END_REF], that use them as financial forecasting. However, 1-D CNNs cannot process streaming data. Hence some works have recently proposed the combination of RNN with 1-D CNN in order to model time-series signals [START_REF] Long | Deep learning-based feature engineering for stock price movement prediction[END_REF]. It must be noted that deep neural networks are computationally intensive [START_REF] Abbe | Provable limitations of deep learning[END_REF]. Furthermore, deep neural networks only provide data estimates for each time step and do not provide uncertainty quantification, while such information would be necessary for stock forecasting to assess risks. Therefore, recent works have been dedicated to combine probabilistic forecasting and deep learning techniques, so as to predict the probability distribution of future events in the time series given its past/historical recordings [START_REF] Salinas | Deepar: Probabilistic forecasting with autoregressive recurrent networks[END_REF]; [START_REF] Jiang | Applications of deep learning in stock market prediction: Recent progress[END_REF]. Deep factor model based on dropout-based heuristic and complex semantics have also been considered in [START_REF] Chauhan | Uncertainty-aware lookahead factor models for quantitative investing[END_REF]. These techniques provide probability distributions as outputs, thanks to specific learning strategies inherited from the Bayesian NN litterature. However, these works, up to our knowledge, do not mention any explicit strategy to estimate uncertainty/confidence score on future predictions/decisions that would help to assess their reliability. Standard (non deep) machine learning models have also been combined to statistical time series modeling tools. For instance, the work [START_REF] Pai | A hybrid arima and support vector machines model in stock price forecasting[END_REF] uses ARIMA and support vector machine.

Another work uses a combination of ARIMA and random forests for the same task [START_REF] Kumar | Forecasting stock index returns using arima-svm, arima-ann, and arima-random forest hybrid models[END_REF]. Finally, several works, such as [START_REF] Ding | Deep learning for event-driven stock prediction[END_REF]; [START_REF] Chong | Deep learning networks for stock market analysis and prediction: Methodology, data representations, and case studies[END_REF]; [START_REF] Fischer | Deep learning with long short-term memory networks for financial market predictions[END_REF], use information mined from news articles and blogs via natural language processing for stock forecasting [START_REF] Cheng | Financial time series forecasting with multi-modality graph neural network[END_REF]; [START_REF] Li | Dynamic graph structure learning for multivariate time series forecasting[END_REF]. Strictly speaking, this is not artificial intelligence since these are dependent on human cognizance.

Contributions compared to existing literature

As a summary, SSMs are valuable tools for probabilistic time series modeling and inference. But there is a crucial need for new strategies to cope with 2019), the main difference being in the priors imposed in the latent factors (positivity, in the case of deep NMF, low-rank/sparsity in the case of deep DL). In our work, deep NMF is neither used for unsupervised representation learning nor in an NN framework. In contrast, it is embedded into a Gaussian SSM to model, allowing to track and predict complex latent phenomena in time series. A novel algorithm is proposed, that learns the positive latent factors jointly with the probabilistic state inferential task induced by our SSM. We call this modelling and inference tool Deep Recurrent Dictionary Learning (DRDL).

We further specialize this tool, to make it practically efficient in the context of large and volatile time series arising in stock market data. In particular, we perused the online training strategy we previously introduced in Sharma et al.

(2021).

Contributions in a nutshell. In this work, we:

• Introduce an LG-SSM model involving deep positive latent factors;

• Propose a new EM-based inference method to jointly perform the time series prediction task and the deep linear positive factors estimation;

• Devise efficient implementation strategies for practical use of the method in the context of stock market time series analysis;

• Investigate through experiments and benchmark comparisons on real financial data the performance of the novel DRDL approach.

Paper organization

The rest of the paper is organized as follows. The proposed DRDL model and its inference is presented in the following Section 2. The practical implementation of DRDL in the context of stock market data analysis is discussed in Section 3. The experimental results and their analysis are described in Section 4. The conclusions of this work are discussed in section 5.

Paper organization

The rest of the paper is organized as follows. The proposed DRDL model and its inference is presented in the following Section 2. The practical implementation of DRDL in the context of stock market data analysis is discussed in Section 3. The experimental results and their analysis are described in Section 4. The conclusions of this work are discussed in section 5.

Proposed Work

Considered model

Let us consider an observed sequence (x k ) 1≤k≤K of vectors of size N x ≥ 1.

We aim to estimate (z k ) 1≤k≤K , a sequence of unknown hidden/latent vectors of size N z ≥ 1. The DRDL approach relies on the following re-parametrized

LG-SSM:

For every k ∈ {1, . . . , K},    z k = D 0 D 1 D 2 z k-1 + v 1,k , x k = H 0 H 1 H 2 z k + v 2,k .
(1) where D 0 ∈ [0, +∞) Nz×Nz , D 1 ∈ [0, +∞) Nz×Nz , and D 2 ∈ [0, +∞) Nz×Nz are three positive-valued linear factors leading to a multi-linear state oper-

ator D 0 D 1 D 2 . Similarly, H 0 ∈ [0, +∞) Nx×Nz , H 1 ∈ [0, +∞) Nz×Nz , H 2 ∈
[0, +∞) Nz×Nz are three positive-valued linear factors yielding the multi-linear observation model H 0 H 1 H 2 .1 The process noise (v 1,k ) 1≤k≤K is assumed to have a Gaussian distribution with zero-mean and symmetric definite positive covariance matrix Q ∈ R Nz×Nz . The observation noise (v 2,k ) 1≤k≤K , is also assumed zero-mean Gaussian with symmetric definite positive covariance matrix R ∈ R Nx×Nx . We consider z 0 ∼ N (z 0 , P 0 ) as the initial state, with z0 ∈ R and P 0 ∈ R Nz×Nz defined as definite symmetric positive matrix. The model in Eq. ( 1) can be interpreted as a multi-linear Gaussian model involving a sequence of K hidden states represented by (z k ) 1≤k≤K . As discussed earlier, classical inference approaches for SSM require specifying every model parameters.

In the model above, this would mean setting the positive latent factor matrices matrices {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } involved in both state and observation models. In practical applications such as stock market analysis, these parameters are unknown and must be learnt from the observed data. The objective is thus to provide a point-wise estimate of the positive latent factor matrices

{D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } and a probabilistic estimate of sequence (z k ) 1≤k≤K ,
given the observed sequence (x k ) 1≤k≤K . This can be seen as solving jointly (i) two deep NMF problems, (ii) a filtering/smoothing problem.

Discussion on the model

We now discuss the main characteristics of the proposed DRDL method.

The model is mathematically described in Eq. ( 1) and displayed in Fig. 1. The top equation describes the hidden state evolution, assuming Markovianity between two consecutive hidden states. The second equation links the hidden and observed states. A first interesting aspect, inherited from the SSM paradigm, is that two Gaussian noise terms are explicitly introduced in DRDL to cope with model uncertainty, which is in contrast with most deep learning models for time series processing (e.g., LSTM). A second novel feature of (1) lies in using deep NMF models instead of generic matrices (in the linear case) or functions (in the non-linear case), as it is usually the case in SSMs Andrieu et al. ( 2010

Z k-2 Z k-1 X k-2 X k-1 Z k X k D 0 D 1 D 2 H 0 H 1 H 2 H 0 H 0 H 1 H 1 H 2 H 2 D 0 D 1 D 2 Figure 1: Schematic Diagram for Deep Recurrent Dictionary Learning

DRDL inference algorithm

Using SSM models for time series processing (e.g., for a prediction task)

amounts to solving the so-called smoothing/filtering problem, i.e., the probabilistic estimation of the hidden state (z k ) 1≤k≤K . In our context, as the deep

NMF factors {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } involved
in the construction of the state transition and the observation transition models are most often unknown, we must also infer them from the observed data, jointly with the hidden states (through the aforementioned filtering/smoothing procedure). To do so, we propose an expectation-maximization (EM) approach (see (Särkkä, 2013, chap.12) and [START_REF] Shumway | An approach to time series smoothing and forecasting using the em algorithm[END_REF]). The EM method alternates iteratively between the probabilistic inference of the state (z k ) 1≤k≤K , while the positive

factors {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } are fixed (E-step)
, and the update of these factors, assuming fixed state (M-step). More precisely, the E-step consists in fixing the linear operators obtained in the previous M-step and applying the classical Kalman/RTS recursions, for obtaining the filtered/smooth distributions p(z k |x 1:k ) and p(z k , x 1:K ), respectively. Then, the M-step updates the

operators {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } by maximizing an upper bound of: φ K (D 0 , D 1 , D 2 , H 0 , H 1 , H 2 ) = log p(x 1:K |D 0 , D 1 , D 2 , H 0 , H 1 , H 2 ). (2)
We explicit hereafter the construction of the i+1-th EM update, given estimates from the previous iteration i.

E-step: Kalman/RTS inference

At this step, we considered the factors

D [i] 0 , D [i] 1 , D [i] 2 , H [i] 0 , H [i] 1 , H [i]
2 to be fixed (either from the previous M-step or from the initialization at the first iteration), and the goal is the probabilistic estimation of the latent states. As we aforementioned, ( 1) is a multi-linear Gaussian model whose observation operator

H 0 H 1 H 2 , evolution/state operator D 0 D 1 D 2 ,
p(z k |x 1:k ) = N (z k ; zk , P k ). (3)
For every k, the mean zk and the covariance P k are given by the Kalman iterations:

For k = 1, . . . , K:

Predict state:

   z k|k-1 = D [i] 0 D [i] 1 D [i] 2 zk-1 , P k|k-1 = D [i] 0 D [i] 1 D [i] 2 P k-1 (D [i] 0 D [i] 1 D [i] 2 ) ⊤ + Q. (4)
Update state:

                     y k = x k -H [i] 0 H [i] 1 H [i] 2 z k|k-1 , S k = H [i] 0 H [i] 1 H [i] 2 P k|k-1 (H [i] 0 H [i] 1 H [i] 2 ) ⊤ + R, K k = P k|k-1 (H [i] 0 H [i] 1 H [i] 2 ) ⊤ S -1 k , zk = z k|k-1 + K k y k , P k = P k|k-1 -K k S k K ⊤ k .
(5)

Hereabove, y k represents the measurement pre-fit residual, S k represents the pre-fit covariance, K k represents Kalman gain, zk represents the updated (a posteriori) state estimate, P k represents the updated (a posteriori) covariance estimate. The backward recursion from the RTS smoother allow to build the smoothing distribution p(z k |x 1:K ).

For k = K, . . . , 1

Backward Recursion:

                     z - k+1 = D [i] 0 D [i] 1 D [i] 2 zk , P - k+1 = D [i] 0 D [i] 1 D [i] 2 P k (D [i] 0 D [i] 1 D [i] 2 ) ⊤ + Q, G k = P k (D [i] 0 D [i] 1 D [i] 2 ) ⊤ [P - k+1 ] -1 , z s k = zk + G k [z s k+1 -z - k+1 ], P s k = P k -G k [P s k+1 -P - k+1 ]G ⊤ k . (6)
Consequently, for every time step k ∈ {1, . . . , K}, the RTS smoother provides:

p(z k |x 1:K ) = N (z k ; z s k , P s k ). (7)

M-step: Evolution operators update

The M-step performs an optimization step to increase the likelihood of the positive latent factors {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 }, given the smoothing distribution obtained in the E-step. It proceeds by building the upper-bound:

φ k (D 0 , D 1 , D 2 , H 0 , H 1 , H 2 ) ≥ Q(D 0 , D 1 , D 2 , H 0 , H 1 , H 2 ; Θ [i] ). (8)
Hereabove, Θ [i] = {Σ [i] , Φ [i] , B [i] , C [i] , ∆ [i] } gathers five quantities defined from the outputs of the E-step described in Sec. 2.3.1):

Q(D 0 , D 1 , D 2 , H 0 , H 1 , H 2 ; Θ [i] ) = - K 2 tr Q -1 Σ [i] -C [i] (D 0 D 1 D 2 ) ⊤ -D 0 D 1 D 2 (C [i] ) ⊤ +D 0 D 1 D 2 Φ [i] (D 0 D 1 D 2 ) ⊤ - K 2 tr R -1 ∆ [i] -B [i] (H 0 H 1 H 2 ) ⊤ -H 0 H 1 H 2 (B [i] ) ⊤ +H 0 H 1 H 2 Σ [i] (H 0 H 1 H 2 ) ⊤ , (9)
with:

Σ [i] = 1 K K k=1 P s k + z s k (z s k ) ⊤ , Φ [i] = 1 K K k=1 P s k-1 + z s k-1 (z s k-1 ) ⊤ , B [i] = 1 K K k=1 x k (z s k ) ⊤ , (10)
C [i] = 1 K K k=1 P s k G ⊤ k-1 + z s k (z s k-1 ) ⊤ , ∆ [i] = 1 K K k=1 x k x ⊤ k .
The updates {D

[i+1] 0 , D [i+1] 1 , D [i+1] 2 , H [i+1] 0 , H [i+1] 1 , H [i+1] 2
} given the knowledge of Θ [i] , amounts to maximizing Q(•; Θ [i] ) under positivity constraints on the factors. In contrast with the linear unconstrained model case studied in (Särkkä, 2013, Chapter 12), the maximization problem here does not have a closed-form solution. It is highly non-convex due to the multi-linearity of our model. Luckily, it happens to be convex with respect to each of the factors. We thus propose to resort to the following alternating maximization step:

D [i+1] 0 = argmax D 0 ≥0 Q(D0, D [i] 1 , D [i] 2 , H [i] 0 , H [i] 1 , H [i] 2 ; Θ [i] ) D [i+1] 1 = argmax D 1 ≥0 Q(D [i+1] 0 , D1, D [i] 2 , H [i] 0 , H [i] 1 , H [i] 2 ; Θ [i] ) D [i+1] 2 = argmax D 2 ≥0 Q(D [i+1] 0 , D [i+1] 1 , D2, H [i] 0 , H [i] 1 , H [i] 2 ; Θ [i] ) H [i+1] 0 = argmax H 0 ≥0 Q(D [i+1] 0 , D [i+1] 1 , D [i+1] 2 , H0, H [i] 1 , H [i] 2 ; Θ [i]
)

H [i+1] 1 = argmax H 1 ≥0 Q(D [i+1] 0 , D [i+1] 1 , D [i+1] 2 , H [i+1] 0 , H1, H [i] 2 ; Θ [i]
)

H [i+1] 2 = argmax H 2 ≥0 Q(D [i+1] 0 , D [i+1] 1 , D [i+1] 2 , H [i+1] 0 , H [i+1] 1 , H2; Θ [i] )
This approach ensures by construction the following inequality:

Q(D [i+1] 0 , D [i+1] 1 , D [i+1] 2 , H [i+1] 0 , H [i+1] 1 , H [i+1] 2 ; Θ [i] ) ≥ Q(D [i] 0 , D [i] 1 , D [i] 2 , H [i] 0 , H [i] 1 , H [i] 2 ; Θ [i] ), (11)
which is key to guarantee the convergence properties for the EM iteration. Indeed, the proposed updates yield an increase of the lower bound of the marginal likelihood, so as a consequence, an increase of the marginal log-likelihood itself.

The overall procedure is thus guaranteed to yield a monotonic increase of the marginal log-likelihood function φ K and classical results about majorizationminimization methods allow to ensure convergence guarantees to a stationary point of φ K [START_REF] Jacobson | An expanded theoretical treatment of iteration-dependent majorize-minimize algorithms[END_REF]. The six sub-problems are quadratic programming (convex) problems and can be solved through several available solvers.

We decided to use the simple and fast alternating least squares approach Ci- This yields the following analytic updates:

D [i+1] 0 = ReLu C [i] (D [i] 2 ) ⊤ (D [i] 1 ) ⊤ (D [i] 1 D [i] 2 Φ [i] (D [i] 2 ) ⊤ (D [i] 1 ) ⊤ ) † , D [i+1] 1 = ReLu ((D [i+1] 0 ) ⊤ Q -1 D [i+1] 0 ) † (D [i+1] 0 ) ⊤ Q -1 C [i] (D [i] 2 ) ⊤ ×(D [i] 2 Φ [i] (D [i] 2 ) ⊤ ) † , D [i+1] 2 = ReLu ((D [i+1] 1 ) ⊤ (D [i+1] 0 ) ⊤ Q -1 D [i+1] 0 D [i+1] 1 ) † (D [i+1] 1 ) ⊤ ×(D [i+1] 0 ) ⊤ Q -1 C [i] (Φ [i] ) -1 , H [i+1] 0 = ReLu B [i] (H [i] 2 ) ⊤ (H [i] 1 ) ⊤ (H [i] 1 H [i] 2 Σ [i] (H [i] 2 ) ⊤ (H [i] 1 ) ⊤ ) † , H [i+1] 1 = ReLu ((H [i+1] 0 ) ⊤ R -1 H [i+1] 0 ) † (H [i+1] 0 ) ⊤ R -1 B [i] (H [i] 2 ) ⊤ ×(H [i] 2 Σ [i] (H [i] 2 ) ⊤ ) † , H [i+1] 2 = ReLu ((H [i+1] 1 ) ⊤ (H [i+1] 0 ) ⊤ R -1 H [i+1] 0 H [i+1] 1 ) † ×(H [i+1] 1 ) ⊤ (H [i+1] 0 ) ⊤ R -1 B [i] (Σ [i] ) -1 . (12)
Hereabove, (•) † denotes the pseudo-inverse operator. Moreover, ReLu (•)

states for the rectified linear unit function, that projects each entry of its input to the positive orthant.

The DRDL algorithm summarized

We summarize in Alg. 1 the DRDL algorithm, for the probabilistic inference of the sequence of hidden state (z k ) 1≤k≤K , jointly with the point-wise estimation of the latent factors {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 }, assuming the data follows the DRDL model (1). In practice, DRDL algorithm is ran for a maximum number of iterations i max , set so as to reach stabilisation of the latent factors.

Application to Stock Trading

We now particularize the DRDL inference algorithm to the stock trading applications. In particular, we address the forecasting/trading tasks given a set of K daily (i.e., k is a day index) observations of stock market data.

Algorithm 1. DRDL (3 layers) inference algorithm.

Inputs. Prior parameters (z 0 , P 0 ) ; model noise covariance matrices

Q, R ; set of observations {x k } 1≤k≤K .
Initialization. Set positive latent factors

{D (0) 0 , D (0) 1 , D (0)
2 , H

0 , H

1 , H

2 }.

Recursive step. For i = 0, 1, . . . , i max :

(E step) Run the Kalman filter (4)-( 5) and RTS smoother (6) using latent factors {D

(i) 0 , D (i) 1 , D (i) 2 , H (i) 0 , H (i) 1 , H (i) 2 }. Calculate (Σ (i) , Φ (i) , B (i) , C (i) , ∆ (i) ) using (10). (M step) Compute {D (i+1) 0 , D (i+1) 1
, D

, H

, H

, H

} using (12).

Output. State filtering/smoothing pdfs (3) and (7) along with pointwise estimates of the latent factor from (12).

Online implementation

First, in order to better cope with high volatility of stock market quantities and allow immediate feedback to the users for on-the-fly trading, we propose here an online implementation of our DRDL approach. We make use of sliding windows of size of τ ∈ {1, . . . , K} time steps. The model parameters are inferred for every k ∈ {0, . . . , K -τ } using the last τ data points observed in the window, i.e. X k = {x j } k+τ j=k+1 , then followed by the EM approach described in detail above. The sliding window approach leverages two advantages. First, it helps in faster processing of the sequence as one can choose a small number τ of data points in the window. Second, it might provide better modeling, since the constant linear factors assumption is likely to better model the time series if τ is small. However, a too small τ might also degrades the inference capabilities.

Hence, it is essential to find a tradeoff in finding an optimal τ , as we will show in our experiments. When implementing the online strategy, a warm start approach is employed for the Kalman filter initialization. The observation/state factors are set to their most recent values, and the mean/covariance of the state for processing X k+1 are initialized using the results of the processing of X k . Let us note that, when we set τ = K, we retrieve the offline version of the algorithm, where the EM inference tool is ran only once.

Modeling and post-processing for stock market analysis tasks

Stock market data processing typically amounts to solving two distinct applicative problems, namely daily stock price forecasting and stock trading decision (among 3 options: buy/hold/sell) estimation. We hereafter explain how to post-process DRDL results to tackle both above-stated problems.

Stock forecasting

Let us first specify the observation model in stock forecasting. For a given window size τ > 0, for each k ∈ {0, . . . , K-τ }, we observe (x j ) k+1≤j≤k+τ ∈ R 15 , gathering 14 technical indicators 2 as well as the adjusted close price. Running our DRDL on the considered window yields the following mean estimate of the 15 quantities for the next time step indexed as k + τ + 1:

xk+τ+1 = H 0 H 1 H 2 z k+τ |k+τ -1 , (13)
2 We retained the relative strength index (RSI), the William percentage range, the absolute price oscillator (APO), the commodity channel index, the Chande momentum oscillator (CMO), the directional movement Indicator (DMI), the ultimator oscillator, the WMA, the exponential moving average (EMA), the Simple Moving Average (SMA), the triple EMA, the moving average convergence (MAC), the percentage price oscillator, the rate of change (ROC). Detailed definitions can be found in https://www.investopedia.com/terms/t/technicalindicator.asp with the covariance matrix

S k+τ +1 = H 0 H 1 H 2 (D 0 D 1 D 2 P k+τ (D 0 D 1 D 2 ) ⊤ + Q) + R. (14)
Hereabove, {D 0 , D 1 , D 2 , H 0 , H 1 , H 2 } are the factors estimated during the Mstep of our EM-based inference method, and (z k+τ |k+τ -1 , P k+τ ) are byproducts of the Kalman prediction step (4)-( 5), computed during the E-step of the EM.

The proposed methodology aims at predicting the entire 15-dimensional vector.

However, stock forecasting is typically focused on the prediction of a single quantity such as the adjusted close price.

Stock trading

In stock trading, a different set of inputs are passed to the model. For each window index k ∈ {0, . . . , K -τ }, the observed data points (x j ) k+1≤j≤k+τ ∈ R 17 , x j [i], for i ∈ {1, . . . , 14}, are the same 14 technical indicators as in stock forecating. Additionally, [x j [15], x j [16], x j [17]] ∈ {0, 1} 3 gathers the decisions "hold", "buy", or "sell", which are calculated for each stocks for every day so as to maximize the annualized returns. The labels are further turned into soft hot encoded vectors as explained in (Sharma et al., 2021, Sec. 3.3.2). The mean and covariance of these 17 quantities can be estimated for next day following ( 13) and ( 14). We then define our class label for next time step as

ℓ k+τ +1 = argmax i∈{1,2,3} x k+τ +1 [i + 14]. (15)

Probabilistic assessment of stock trading decision

We now describe the procedure to assess the uncertainty quantification associated to the DRDL predictions. Let k ∈ {0, . . . , K -τ } be the window index on which Algorithm 1 has been run. The probabilistic estimation of the quantities of interest for the next time step (i.e., one-day ahead prediction) x k+τ +1 conditioned to the data observed in the window x k:k+τ , reads as a multivariate Gaussian distribution

p(x k+τ +1 |x k:k+τ ) = N (x k+τ +1 ; xk+τ+1 , S k+τ +1 ) , (16)
with mean and covariance (x k+τ +1 , S k+τ +1 ), given by ( 13) and ( 14), respectively. Equation ( 16) assigns a probability score to any decision (e.g., trading)

based on the prediction output of the DRDL method. Let us focus on the particular example of assessing the uncertainty of the stock trading decision at time index k + τ + 1, given observations at indexes j ∈ {k + 1, . . . , k + τ }. The trading decision relies on the discrete maximization step (15). Let us express the probability mass function (pmf) of this decision, from the gaussian predictive probability density function (pdf) of the observed data points in Eq. ( 16).

The pmf can here be summarized as p k+τ +1 ∈ [0, 1] 3 where each p k+τ +1 [i],

i ∈ {1, 2, 3} is a probability, and

3 i=1 p k+τ +1 [i] = 1. Each p k+τ +1 [i] represents
the probability inferred by DRDL that the true value x k+τ +1 [i + 14] is greater than x k+τ +1 [j + 14], for j = {1, 2, 3} \ i. According to Eqs. ( 16) and ( 15), p k+τ +1 can be obtained through with

(∀i ∈ {1, 2, 3}) p k+τ +1 [i] = Yi N (
Y i = y ∈ R 3 | y[i] ≥ y[j], j = {1, 2, 3} \ i . ( 18
)
Due to the intricate form of the constrained set in ( 18), the integral in ( 17) is intractable. It can be easily approximated with high precision by direct simulation. In practice, we sampled 10 4 three-dimensional sample from a normal standard distribution. The samples can be re-used for all time steps using coloring and shifting according to the covariance and mean, respectively. Thanks to this procedure, we can infer p k+τ +1 for every k, and then assess the next day stock trading outcome by using the standard cross-entropy loss:

log-loss = 1 K -τ + 1 K-τ k=0 3 i=1 -(L k+τ +1 [i] log(p k+τ +1 [i])), ( 19
)
where the true labels are denoted L k+τ +1 ∈ {0, 1} 3 for each time k + τ + 1

(hereagain, we use soft hot encoding representation).

Summarized pipeline

We provide in Alg. 2 the summary of our proposed pipeline for applying DRDL, in Algorithm 1, in the context of stock forecasting (steps a-b with N x = 15) and trading (steps a-b-c-d with N x = 17).

Algorithm 2. DRDL (3 layers) method for stock forecasting and trading.

Inputs. Prior parameters (z 0 , P 0 ) ; model noise covariance matrices Q, R ; set of observations {x k } 1≤k≤K ; windows size τ .

Initialization. Set positive latent factors

{D (0) 0 , D (0) 1 , D (0)
2 , H

0 , H

1 , H

2 }.

Window processing. For k = 0, 1, . . . , K -τ : a. Run DRDL algorithm 1 on sequence (x j ) k+1≤k≤k+τ , initialized with estimates from k -1th window (warm start).

b. Calculate one-step ahead predicted mean x k+τ +1 and its covariance S k+τ +1 using (13)-( 14).

c. Compute one-step ahead predicted label ℓ k+τ +1 using (15). (2019). The diversification also allows to assess the model robustness to various trends [START_REF] Kumar | Expanding the role of marketing: from customer equity to market capitalization[END_REF]. From the knowledge of the close prices, we build two observation sequences associated to the resolution of two specific problems, namely stock forecasting and stock trading, as described in Sec.

D 0 = D 1 = D 2 = Id and {H 0 , H 1 , H 2 } are estimated.
Note that, ignoring the positivity constraint, DRDL (1 layer) would identify with our previously published method RDL Sharma et al. (2021). We implement the sliding window approach described in Sec. 3.1, for various choices of τ described hereafter. In all experiments, the initial value are set as ; z0 is an all zero vector; P 0 =10 -7 Id, Q=10 -2 Id and R=10 -2 Id, where Id states for the identity matrix. Moreover, we set the dimension of the state as N z = 14, which also corresponds to the number of measured technical indicators, we observed better performance of the model. The entries of the linear factors to estimate are initialized at time 0 using independent realizations of a uniform distribution on

Compared methods

The proposed DRDL approach is analyzed by comparing with state-of-the- Sezer & Ozbayoglu (2018) and ARIMA [START_REF] Ariyo | Stock price prediction using the arima model[END_REF]. We select the stateof-the-art methods for each task for a fair comparison. For stock forecasting, the comparison is done with ARIMA and LSTM. The ARIMA parameter value are set to (p, d, q) = (5, 1, 5). The LSTM is customized from its original version to carry out regression tasks by replacing the softmax output layer with an affine layer. The Adam optimizer is used with learning rate 10 -4 and 200 epochs.

We used a mini-batch strategy where batch-size is fixed to 16 to reduce the objective function's mean square error (MSE). The evaluation of the methods is done using metrics like mean absolute error (MAE), root means square error (RMSE), SMAPE (Symmetric mean absolute percentage error), and Pearson correlation factor.

For the stock trading task, the comparison is made with CNN-TA, Multifilter neural network, and LSTM implemented with their original architecture The choice of window size is an essential aspect as it can enhance and limit the methodology's potential. To understand the model behavior, we present The challenging task with the DRDL approach is to preserve a balance between the computational time and optimal predictions. We experimented with various window sizes to analyze and preserve the best parameter for our future experiments. We present state that the performance of the approach increases as it feeds more data to the model for better understanding. We can preserve a balance parameter τ = 650,

In Table 5, we also present the computational times (train and test) for conducting trading simulations for our dataset. Hereagain, LSTM is the more demanding method at training. All methods have rather comparable test times, despite DRDL is implemented on CPU only. In particular, besides its probabilistic output, DRDL is not more costly than its competitors. Adding more layers to DRDL slightly increases its train time, but does not affect much the test time.

Annualized Returns

Stock market aims to analyze and evaluate the return on investment for a given stock. Every trader is indeed interested in evaluating his investment returns and taking risks accordingly. We simulate market scenarios [START_REF] Sezer | Algorithmic financial trading with deep convolutional neural networks: Time series to image conversion approach[END_REF] by evaluating the annualized returns by the predicted stock trading decisions provided by DRDL using 1 to 3 layers as well as the decisions from from the benchmark models. Table 9 presents a detailed study of nine stocks for DRDL methodology and state-of-the-art methods. We display only empirical values for nine stocks and the average results over the 180 stocks. To make it easy for readers we have highlighted best annualized returns in bold.

It is clearly evident that the DRDL approach yields higher returns when tested for a duration of 10 years when compared to annualized returns obtained from deep learning state-of-the-art methods predictions.

Portfolio diversification

Many researchers and traders believe that it is essential to know the associated sentiments associated with each stock to understand the stock market.

Traders support and recommend having a mix of stock sentiments in one's portfolio. The market is very well divided into three types of stock sentiments:

small-cap, mid-cap, and large-cap. To read about them in detail please refer to Sharma et al. (2021)[section 4.4.4]. To evaluate this sentiment using the predicted signals from the proposed approach, we calculated probabilistic quan- score to have a well-diversified portfolio. The score provides a confidence score that helps the investor decide where to invest in the market to have a balance of market sentiments and maximize returns.

To understand further, we present Table 10 which provides a log-loss score.

The log-loss score provides the confidence score in terms of its volatility nature, where the smaller value is considered, the better and less volatile. We evaluated the confidence score for the proposed approach for different configuration. The market capitalization of these stocks can be found 6 . The log-loss value provides the probabilistic inference for the predictions. The inference tries to penalize the events for which the method assigns a low probability. We observed that the log-loss value reached a meager value which indicated good prediction accuracy in large-cap stocks, which are expected to be least volatile. In contrast, we achieved a higher log-loss value for predictions associated with small-cap stocks as they are highly unstable and new to the market.

Acknowledgment

The CNRS-CEFIPRA project supported this work under grant NextGenBP PRC2017. E.C. acknowledges support from the European Research Council

Starting Grant MAJORIS ERC-2019-STG-850925.

Conclusion

In our approach, time-series sequences are modeled with a flexible Gaussian SSM.The transition matrices (state and observation models) are unknown, and are estimated thanks to an expectation-minimization strategy, assuming a particular deep NMF structure. The DRDL approach inherits advantages from sophisticated modeling techniques while quantifying the uncertainty in the predictions.We have then adapted the DRDL approach to deal with a challenging large scale financial time series problem, to target stock forecasting and trading tasks. In particular, the method is able to successfully operate in an online processing manner, allowing to capture piece-wise linear characteristics in the data. The results show that the proposed method outperforms the state-of-theart techniques. Given these promising results, we plan as future work to delve deeper into the area of financial forecasting, including the application of our technique in forecasting derivatives.

  the curse of dimensionality in learning SSM model parameters. In this work, we propose to impose a structured prior on the observation/state operators involved in an LG-SSM. We introduce deep nonnegative matrix factorized (deep NMF) models for both operators. Deep NMF De Handschutter et al. (2021) is a generalized form of NMF Cichocki et al. (2009) that models latent representations from complex data through a product of a (usually small) number of linear operators (called latent factors) satisfying positivity constraints. Deep NMF has been employed with success on various unsupervised machine learning tasks Yu et al. (2015); Trigeorgis et al. (2016); Xue et al. (2017); Chen et al. (2021); Flenner & Hunter (2017). When embedded into an NN structure, it leads to the so-called deep ReLu networks Liu & Liang (2021); Daubechies et al. (2022); Chen et al. (2019). Deep NMF shares connections with the recently introduced deep dictionary learning (deep DL) Tariyal et al. (2016); Mahdizadehaghdam et al. (

  );[START_REF] Chopin | SMC2: an efficient algorithm for sequential analysis of state space models[END_REF];[START_REF] Crisan | Nested particle filters for online parameter estimation in discrete-time state-space markov models[END_REF], taking advantage on the acknowledged representation power of deepNMF De Handschutter et al. (2021). One important benefit of the proposed approach w.r.t. most existing methods in the literature is that we avoid Monte Carlo simulation or complex optimization procedure, which is known to suffer more severely the curse of dimensionality. In our method, each latent factor can be understood as representations in abstract spaces of the phenomena occurring between both pairs of variables. Third, in contrast with the typical usage of deep NMF in machine learning, relying on backpropagation for their model training Chen et al. (2021); Flenner & Hunter (2017), DRDL model allows the construction of an handcrafted training strategy (see the next section), which benefits from a low computational cost, sound optimality guarantees (in terms of Bayesian estimator), and enables uncertainty quantification.

  chocki et al. (2009), reminiscent from the literature of deep nonnegative matrix factorization Chen et al. (2021), and the deep ReLu neural networks models Daubechies et al. (2022), both showing a satisfactory behavior in preliminary experiments. We start by computing each subproblem solution ignoring the positivity constraints, and then capped the negative entries of the obtained factors.

  3.2. The data is scaled by normalising the 14 technical indicator values. For both problems, we will compare DRDL and several state-of-the-art methods arising from signal processing and machine learning literature. In all experiments, each of the 180 observed time series is split into two parts, namely a train phase made of the first recorded 2546 days, and a test phase made of the next 2882 days. The train phase is used to learn the models parameters (for instance, the linear factors involved in DRDL), while the test phase is used to evaluate the performance of the learnt models, their parameters being fixed. More details about DRDL and the retained benchmark methods setting are provided in the next subsection.

Figure 2 Figure 2 :

 22 Figure 2 displays the evolution of 4 of the 14 technical indicators used as input of the inference tools, during the test phase. One can notice the high volatility in the observed data.

[ 0 ,

 0 10 -1 ]. As mentioned in Sec. 3.1, to initialize the next processed windows, we used warm start strategy. The estimation of the linear factors (i.e. M-step of the EM method) is only conducted during the training phase. A maximum of 50 iterations of the EM loop are used in Alg. 1, which was observed sufficient to reach stability of the estimated factors. During the test phase, the linear latent factors are fixed, and only the Kalman/RTS inference is ran (i.e., we inhibit M-step in Alg. 1). The scores shown are arithmetic means of 10 random trials, and are computed only during the test phase.

  art-methods deep learning namely Multi-filter neural network Long et al. (2019), Long short term memory Fischer & Krauss (2018), 2-D deep CNN (CNN-TA)

Figure 3 :

 3 Figure 3: Ground truth and inferred adjusted close price on the test phase for four different stocks, using DRDL with 1 to 3 layers, LSTM or ARIMA.

  d), LSTM approach failed to reach satisfying results which might be due to vanishing gradient issues. In cases (a-b), ARIMA performs quite good when compared to its performance in other cases (c-d). In contrast, DRDL (3 layers) reaches stable and satisfactory outcomes. DRDL (2 layers) outperforms DRDL (1 layer) and both benchmark methods but remains lower quality than its 3layers variant. 4.4. Numerical results for the stock trading problem 4.4.1. Influence of the window size

  We curated data for 180 stocks which comprises stocks from USA, UK, India and China. The data is prepared by scrapping daily adjusted close prices, open price, volume, high price, low price for a span of twenty years (i.e., 01/01/1998 to 01/10/2019) using yahoo finance API for Python. Having stocks from different market cap is always advisable by the traders, as it gives them breadth while investing in advanced as well as emerging markets Fawaz et al.

	The finance dataset used for experiments is curated from Yahoo finance
	repository. 3
	Output. Forecasting/trading predictions and log-loss value (19).

d. Compute p k+τ +1 using (17)-(

18

).

Table 1

 1 which provides detailed information on the performance of the model on varying window sizes. The table offers an analysis of various metrics like Pearson

	correlation (r), RMSE (Root Mean Square Error), MAE (Mean Absolute Error),
	and SMAPE (Symmetric mean absolute percentage error) for different window
	sizes τ . The model's performance improves as the window size increases till
	a stabilization point. We can see that a balanced choice is τ = 650 to reach
	stabilized performance on this particular task and dataset. We further use this
	value in upcoming experiments.
	4.3.2. Comparison with benchmark models
	To understand better, we present table 2 which provides comprehensive anal-
	ysis on performance estimation on the stock forecasting problem using DRDL,
	LSTM, ARIMA, DeepAR Salinas et al. (2020), Nbeats Oreshkin et al. (2019),
	and TFT Lim et al. (2019). Table 2 presents comparison in terms of Pearson

Table 1 :

 1 Results of DRDL (3 layers) on stock forecasting problem for different window size. Scores averaged on test phase, on all the 180 stocks.

	Model	r	RMSE↓	MAE(%)	SMAPE(%)
	ARIMA	0.13	78.6 (1.89)	1.23(0.56)	65.5
	LSTM	0.24	297.5 (2.64)	6.12 (0.65)	47
	DeepAR	0.40	73.89 (1.85)	0.43(0.34)	58.13
	Nbeats	0.38	95.52 (1.99)	0.57 (0.12)	63.75
	TFT	0.52	35.79(1.56)	0.35(0.019)	38.25
	DRDL (1 layer)	0.65	23.24 (1.64)	0.19(0.02)	35.3
	DRDL (2 layers)	0.69	14.2 (1.43)	0.15(0.006)	23.2
	DRDL (3 layers) 0.71 13.35(0.37) 0.11(0.003)	18.4

.

Table 2 :

 2 Comparative analysis of DRDL against state-of-the-art methods for stock forecasting problem: Pearson correlation score (r), MAE, RMSE , SMAPE scores and their respective std. deviation on the estimation of next time step adjusted close price, averaged over the data in the test set and the stocks.

	correlation factor r, RMSE, MAE and SMAPE. We can see that DRDL (3 lay-
	ers) architecture outperforms DRDL (2 layers), DRDL (1 layer) as well as the
	other benchmarks models. We also notice that the average performance of TFT
	is comparable to DRDL (1 layer) architecture. Fig A.4 in appendix section

Table 3 :

 3 Averaged time over 10 random runs for processing the dataset (train(hrs) and test(min)), for DRDL and its competitors.

	Method	Sharpe Ratio T-test
	DRDL (3 layers)	2.14	0.63
	DRDL (2 layers)	1.99	0.78
	DRDL (1 layer)	1.84	0.83
	ARIMA	1.03	0.89
	LSTM	0.88	1.34
	DeepAR	2.06	0.58
	TFT	1.22	0.67
	Nbeats	1.33	0.54

Table 4 :

 4 Comparison of Sharpe ratio and T-test score, for DRDL and its com-

	petitors averaged over 180 stocks.
	Appendix A displays the Pearson correlation analysis between ground truth
	daily adjusted close price time series and predicted ones along test phase, using
	DRDL (3 layers) for four representive stock cases. Table 4 presents the statisti-
	cal test (t-test) and stock market simulation (Sharpe Ratio) on the forecasting
	results. We observe that the proposed method with 3 layers architecture gives
	better performance as its average score for 185 stocks is smallest as compared to
	other state-of-the-art methods hence we can conclude that more similarity exists

Table 5 :

 5 Comparison of classification scores of different methods on the stock trading problem. All scores are averaged over 180 stocks and over the days of the test phase.

			F1 Score			Precision			Recall
	Window size (τ )
		Hold	Sell	Buy	Hold	Sell	Buy	Hold	Sell	Buy
	250	0.61	0.23	0.19	0.91	0.15	0.12	0.46	0.51	0.50
	300	0.61	0.23	0.19	0.91	0.15	0.12	0.46	0.50	0.51
	350	0.61	0.26	0.27	0.89	0.18	0.19	0.47	0.53	0.52
	500	0.61	0.26	0.27	0.89	0.18	0.19	0.47	0.52	0.53
	650	0.61	0.30 0.29 0.85 0.26 0.29 0.51 0.54 0.54
	700	0.61	0.30	0.29	0.87	0.21	0.20	0.48	0.53	0.54

Table 6 :

 6 Classification scores of DRDL (3 layers) for varying window size.

Table 3

 3

	presents the computational time for forecasting the next day closing
	price for our dataset. We distinguish the time required to train the methods
	(on the first ten years) and to test them (on the next ten years) using the walk-
	forward method described in (Sharma et al., 2021, Section 4.2.1). We observed

Table 7 :

 7 Classification scores of DRDL (2 layers) with varying window size.

			F1 Score			Precision			Recall
	Window size (τ )
		Hold	Sell	Buy	Hold	Sell	Buy	Hold	Sell	Buy
	250	0.84	0.10	0.10	0.85	0.10	0.10	0.85	0.12	0.12
	300	0.80	0.10	0.12	0.85	0.10	0.10	0.74	0.17	0.14
	350	0.68	0.15	0.15	0.82	0.10	0.10	0.62	0.30	0.31
	500	0.59	0.18	0.22	0.86	0.15	0.14	0.46	0.51	0.51
	650	0.59	0.19 0.23 0.88 0.15 0.12 0.45 0.52 0.51
	700	0.59	0.24	0.22	0.90	0.16	0.14	0.46	0.51	0.52

Table 8 :

 8 Classification scores of DRDL (1 layer) with varying window size.

	the highest computational time with the LSTM approach. The other methods
	have rather similar computational time, DRDL (1 layer) being the fastest. The
	computational time of DRDL (3 layers), reaching the best performance metrics,
	stays reasonable, and is comparable with the one of DeepAR and TFT. Here, we
	must recall that, in contrast with most of its competitors (except ARIMA), our
	implementation of DRDL method (for both train/test phases) does not exploit
	GPU facilities such as PyTorch. Complexity reductions could certainly occur if
	this was the case.
	The stock forecasting results are presented in Fig. 3. The comparison is
	carried out between the proposed DRDL for different layer number, LSTM,
	ARIMA, DeepAR, TFT, Nbeats method. We observed that in some cases (c-

Table 9 :

 9 Annualized returns resulting from the stock trading decisions of different methods during the test phase.

Table 6

 6

	which depicts the experimental performance

Table 10 :

 10 tification, as explained in sec. 3.3. The practitioner uses this quantification Comparative analysis of uncertainty quantification provided by DRDL using 1 to 3 layers. The quantification is listed for stocks with market capitalization categories. The log-loss is computed over the test phase.

	Small-cap	ALOKTEXT.BO ALKYLAMINE.BO ZEEMEDIA6.BO	1.04 1.17 0.89	1.20 1.19 0.99	1.09 1.34 0.23
		PVP.BO	1.34	1.98	2.78
		IOC.BO	1.02	1.05	0.87
	Mid-cap	TATACHEM.BO SPICEJET.BO	0.76 0.34	0.45 0.65	0.94 1.20
		BHEL.BO	0.20	0.51	1.15
		AAPL	1.13	0.98	1.11
	Large-cap	AMZN HINDZINC.BO ONGC.BO	0.11 0.03 0.20	0.41 0.65 0.13	0.43 0.45 0.09
		SIEMENS.NS	0.12	0.02	0.11

Throughout the paper, we consider three-terms factorizations, for the sake of readibility.The 3-layers modeling and inference methodology has the great advantage of being generic enough to be straightforwardly extended to any number, greater or equals to one, of factors.

https://yahoo.finance.com

http://ta-lib.org

http://www.ta4j.org

https://finance.yahoo.com/screener

which indicates stabilized performance. We further use this value in upcoming experiments.

Classification metrics

To explain the empirical analyses of the trading process (classification), we present confusion matrices. The trading process involves classifying the signal into three classes, namely "Buy," "hold," and "sell" classes. The summarized performance for 180 stocks by DRDL and other state-of-the-art methods is presented in Fig. A.5 in Appendix A. Among the three classes, we see the prediction of hold class is captured efficiently when compared to the other classes. The LSTM approach predicts the best score over the other state-of-the-art methods when compared to false negatives scores. However, in LSTM and CNN-TA approaches are highlighted many false positives for the "hold" class. It can be noted that these deep learning techniques have labeled most signals as hold class, jeopardizing the model behavior for the other classes ("buy," "sell"). The nature of the finance market is highly volatile and non-linear; hence we get to see a highly imbalanced dataset. However, we noticed that the DRDL approach handles it by imposing an activation function on the operators. These operators are expected to evolute continuously as we grow deeper with time sequence.

Table. 5 and