EEG signals acquired at different electrodes can be modelled as Signals on Graph, where the graph structure reflects the underlying brain Functional Connectivity (FC), representing brain region interactions. FC gives crucial information to detect cognitive states, such as those involved in motor imagery(MI)-based brain-computer interface (BCI) tasks. FC estimators typically leverage similarity measures between pairs of signals recorded at EEG electrodes in order to build a graph. Here, we consider a multilayer network in order to integrate several graphs, corresponding to different FC measures. Specifically, each layer of the network corresponds to one FC estimator. We model the overall Laplacian as a Cartesian product of the single graphs with a path graph where the number of nodes is the number of layers of the multilayer network. Then, we investigate the mathematical formulation and properties of the covariance matrix of the Laplacian of the product graph, which enables us to write in a novel way the Jensen divergence between two conditions. In order to test our method, we perform study on real EEG data recorded during motor imagery experiments and we demonstrate the ability of our method to statistically distinguish observations under different MI tasks. Our findings opens new perspectives in the detection of connectivity states that can be modelled as multilayer network.