\section*{Introduction}
The world around us is constantly changing, which makes dynamics of the artificial and natural systems at different spatial scales is of great interest to the scientific community. Traditionally, to describe a dynamical system one explores a time-dependent behaviour of a relevant physical quantity. In certain situations the choice of such a quantity is obvious, for instance description of a planetary motion is fulfilled by analyzing the trajectories, time-varying coordinates of the planets. On the other hand, in the case of the systems out of equilibrium the search for a suitable quantity could be a non-trivial task. Here, a bright primer known from the school textbooks is a microparticle in a liquid environment characterized by Brownian motion \cite{Brown}. It was found that observations over the trajectory of an individual particle does not allow to characterize the random motion, the only time-dependent mean-square displacement \cite{Einstein, Langevin} or velocity autocorrelation function averaged over many particles trajectories matter. Importantly, this fundamental finding has determined a general strategy for the analysis of non-equilibrium motion in physical systems of one or a few bodies within the standard statistical approaches or advanced machine learning techniques\cite{Granik, ASCnano}.

\begin{figure}[!hb]
		\centering
		\includegraphics[width=0.45\textwidth]{fig1.png}
		\caption{\label{fig1} Pictorial representation of the protocol realizing the concept one process - one number. A physical process of our interest is recorded (A) and then is converted into a set of frames with one-dimensional digital representations (B) of the system at each time moment. (C) These frames are merged in a space-time grid. (D) By using the renormalization procedure with spatio-temporal filters discussed in the text we attribute the initial physical process to a single number (spatio-temporal dissimilarity) characterzing the dissimilarity of the patterns in the digital system's representation on the space-time grid.}
	\end{figure}

A wide diversity of non-equilibrium processes in classical systems stimulates permanent search for distinct measures \cite{compression1, compression2, Falk} that allow  theoretical characterization of the observed dynamics. Of particular difficulty is the description and prediction of a collective behavior of classical systems consisting of many interacting and simultaneously moving elements\cite{protein, colloid} or participants\cite{flocks, sheeps, pedestrian1,pedestrian2}. The standard choice here is the calculations of correlation functions that are non-local in both time and space for some physical quantities characterized the system in question. Depending on the situation these quantities can be taken directly from the experiment as they are (for instance, time-dependent positions of the particles obtained with a tracking system \cite{colloid,Scalliet}), or one could define some useful features \cite{Gnesotto} on the basis of available experimental data.  The problem is that even when defining the spatio-temporal correlators there is no a universal method for extracting a meaningful information from them, for instance, a signature of a hidden dynamical order. One of the possible solutions to simplify the consideration is to decouple the intertwined spatial and temporal correlations and explore them separately \cite{flocks, sheeps}. It means that potentially important information on dynamics of the system in question may be lost. In certain situations the choice of measures taken from information theory and related to entropy or complexity gives very promising results \cite{fallpaper1, fallpaper2, Mutual_info, Palmer} in solving the characterization task.

Non-equilibrium physics is not limited to the case of the classical systems. At the moment, significant theoretical and experimental efforts are concentrated around quantum matter and phenomena out of equilibrium for which none classical analogues exist. Exploration of thermalization dynamics \cite{thermalization1}, hydrodynamic transport \cite{transport},  quantum information scrambling \cite{scrabmling_Google}, non-equilibrium phases that are formed in periodically-driven quantum systems \cite{CNN_MBL, DTC_IBM, DTC_Google}, topological order emerge in out-of-equilibrium states of programmable Rydberg atom arrays \cite{Lukin1}, quantum chaos \cite{quantum_chaos} and other non-trivial findings play a crucial role in further developing quantum technologies, information theory and condensed matter physics. In this regard, quantum simulators  and quantum computers facilitate the modeling of the non-equilibrium quantum dynamics \cite{Google_Majorana} and represent a promising alternative to classical supercomputers whose using is fundamentally limited due to the exponential growth of the Hilbert space. Similar to the classical case, theoretical analysis of non-equilibrium quantum systems is mainly based on calculations of temporal autocorrelators that are local in space and the spatial correlation functions that are local in time. Such a separate consideration of the spatial and temporal correlations inevitably leads to loss of information on interplay between space and time degrees of freedom.

In this paper we propose an alternative approach for exploring the non-equilibrium dynamics of classical and quantum systems by considering space and time degrees of freedom simultaneously and, thus, without uncoupling spatio-temporal correlations that could contain important information on the state of the dynamical system. To this aim, we follow the idea one process - one number which is schematically visualized in Fig.\ref{fig1}. It includes recording a process at given external conditions (Fig.\ref{fig1} A), converting the recorded movie into digital representation as a spatio-temporal image (Figs.\ref{fig1} B and C) and estimating patterns of the image by means of the structural complexity (dissimilarity) measure recently introduced in Refs.\cite{PNASComplexity, Dissimilarity} (Fig.\ref{fig1} D). As we will show by classical and quantum examples such a consideration of the dynamical process eases the search for critical values of the external parameters at which there are global changes in the dynamical properties of the physical system in question. Generality and flexibility in exploring completely different non-equilibrium quantum or classical systems within the same framework are important features of the proposed scheme. Despite of the fact that there is no correlation function in our work, we perform a characterization of discrete time crystal - a complex out-of-equilibrium phase of the quantum matter. For that we use the spatio-temporal dissimilarity introduced in this investigation and some measures taken from the information theory. It shows that the theoretical description of non-equilibrium dynamical processes without calculating correlation functions is possible.

\section*{Results}
\subsection*{Spatio-temporal dissimilarity}

Our approach to the analysis of the quantum or classical dynamics out of equilibrium is formulated on the basis of the recently proposed renormalization group (RG) algorithm \cite{PNASComplexity} for estimating the structural complexity of an object. The essence of this method is to quantify dissimilarities between the patterns that occur at different spatial scales.  Such an approach was found extremely useful in context of defining equilibrium  phase transitions of various nature~\cite{PNASComplexity, Khatami} as well as verification of quantum wave functions \cite{Dissimilarity}. Since in quantum case this approach is hardly related to the complexity of a quantum state itself and operate with the array of single measurements (bitstrings) the term dissimilarity used to name the procedure\cite{Dissimilarity}.

In this paper, to characterize a dynamical process we propose the spatio-temporal dissimilarity measure. The idea is based on dealing with two-dimensional patterns that occur when one puts the digital representations of the system states at different steps of its time evolution together. Having a one-dimensional discrete system of size $L$ and its evolution over $T$ time steps we can construct a matrix $B^0_{L\times T}$. At every step $k>0$ of a coarse-graining procedure we create a matrix $B^k_{L\times T}$ of the same size, whose elements are calculated as follows
    \begin{equation}\label{Coarse-graining procedure}
        b^k_{ij} = \frac{1}{\Lambda^k_x\Lambda^k_t} \sum_{l} \sum_{m} b^{k-1}_{\Lambda^k_x i+m, \Lambda^k_t j+l}.
    \end{equation}
This means that the whole matrix is divided into blocks of $\Lambda^k_x\times\Lambda^k_t$ size (in this particular case, $k$ is a degree, not an index) and the elements within each block are substituted with its average value. In this work $\Lambda_{\rm x} = \Lambda_{\rm t} = 2$ space-time filters were used, since, as we will show below, they provide a complete description of the one-dimensional gas (classical case) and discrete time crystal (quantum case). However, for arbitrary non-equilibrium system one needs to test different spatio-temporal filters to be sure on extracting maximal useful information with dissimilarity measure. Since we are dealing with quantum or quantum-like systems, the initial matrix is filled with either $-1$ or $1$.

After that we compute the partial spatio-temporal dissimilarities for the patterns separated by one step of RG. In case of simple averaging scheme used in our calculations, this can be done using the following equation

   \begin{equation}\label{D_k}
		\mathcal{D}_{\rm xt}^k = \frac{1}{2LT}\left|\left(\sum_{i,j=1}^{L,T} \left[(b^{k+1}_{ij})^2 - (b^{k}_{ij})^2 \right]\right) \right|.
   \end{equation}

	Finally, given features emerging at every new scale, we obtain
	\begin{equation}\label{Complexity}
		\mathcal{D}_{\rm xt} = \sum_{k=1}^{N}\mathcal{D}_{\rm xt}^k ,
	\end{equation} where $N$ is the total number of RG steps.  Thus, within this protocol each dynamical process can be associated to a single number characterizing the dissimilarity of the system's patterns in the spatio-temporal grid.
 The described procedure is then technically identical to the calculation of the structural complexity of the two-dimensional grayscale image \cite{PNASComplexity}.

An independent choice of the initial filter sizes ($\Lambda_{\rm x}$ and $\Lambda_{\rm t}$) allows a selective quantifying either spatial or temporal patterns. For instance, if one chooses $\Lambda_{\rm x} = 1$ and $\Lambda_{\rm t} > 1$, then $\mathcal{D}_{\rm xt}$ boils down to $\mathcal{D}_{\rm t}$ that estimates temporal patterns individually for each spatial degree of freedom and then averages over them. On the other hand, for $\Lambda_{\rm x} > 1$ and $\Lambda_{\rm t} = 1$, one gets $\mathcal{D}_{\rm x}$ that probes spatial patterns in the non-equilibrium system.

Thus, the complete set of dissimilarity functions, $\mathcal{D}_{\rm xt}$, $\mathcal{D}_{\rm t}$ and $\mathcal{D}_{\rm x}$ including their partial contributions can be used for a low-level characterization of non-equilibrium dynamics in a quantum
 or classical system without searching for quantities that contain basic information on the system and calculating correlation functions for them.

\subsection*{1D lattice gas model}

\begin{figure}[t]
\centering
\includegraphics[width=0.8\columnwidth]{fig2.png}
\caption{\label{fig:clg_model} (Top) Schematic representation of the CLG model on a chain with periodic boundaries. At each time step one randomly chosen active particle (depicted in orange) moves to the unoccupied position. The dynamic stops when there are no active particles (t = 5). (Bottom) Examples of the time evolution of the L = 16 chain with N = 7 (adsorbing state) and N = 9 (active state) occupied sites.}
\end{figure}

\begin{figure}[t]
\centering
\includegraphics[width=1.0\columnwidth]{fig3.png}
\caption{\label{fig:1d_gas} Comparison of the $\mathcal{D}_{\rm xt}$ (red squares) and CID (blue circles) behaviour in case of CLG problem on 1D chains of different size.}
\end{figure}

To assess the viability of the spatio-temporal dissimilarity approach we first consider a simple conserved lattice gas (CLG) model in 1D~\cite{clg1, clg2} which is one of the sandpile models to study self-organised criticality. The initial state of the system in question is a random distribution of $N$ occupied sites on the $L\ge N$ chain. Particles, for which one of its neighbouring site is occupied are considered as active. Dynamics in the CLG model can be described in the following way: at each time step one should randomly choose an active particle and move it to the unoccupied neighbouring position. For $\rho = N/L \le 0.5$ the system ends up in a steady adsorbing state without active sites meanwhile for higher densities the dynamical fluctuations will take place at any time. Thus, we have a nonequilibrium phase transition with the order parameter $f_a$ denoting a fraction of active states on the final lattice configuration. Schematic representation of the CLG model as well as the examples of real time evolution of the system belonging to different phases are given in Fig.~\ref{fig:clg_model}.

The authors of Ref.~\onlinecite{clg_cid} have shown, that the phase transition between the adsorbing and active states in this model can be easily defined in the thermodynamic limit with using a computable information density (CID). This quantity is simply equal to the ratio between the sizes of the losslessly compressed and the initial data files, containing a final state of the chain. The idea is based on the intuitive assumption that the more structured the file, the less information is needed to reproduce its content, which is similar to the concept of Kolmogorov complexity~\cite{Kolmogorov}. It was demonstrated that the method facilitates accurate description of the thermodynamic phase transitions in equilibrium magnetic systems\cite{compression2}, reveals hidden dynamical ordering in classical systems out of equilibrium \cite{compression3} and helps to solve other interesting problems \cite{compression4, compression5}. Importantly, an effective implementation of the compression algorithm assumes the possibility to generate large amount of data and the data itself should fit a dictionary with discrete set of items.

\begin{figure}[t]
\centering
\includegraphics[width=1.0\columnwidth]{fig4.png}
\caption{\label{fig:d_avg} Comparison of the $\mathcal{D}_{\rm xt}$ (purple triangles), $\mathcal{D}_{\rm x}$ (orange squares) and $\mathcal{D}_{\rm t}$ (green circles) behaviour in case of CLG problem on 1D chain with $L=64$, $T=1024$. The results are averaged over 100 independent starting configurations at each density.}
\end{figure}

Here we compare $\mathcal{D}_{\rm xt}$ with CID in case of small chains with up to $L = 256$ sites. As can be seen from Fig.~\ref{fig:1d_gas}, the maximum value of dissimilarity is not sensitive to the size of the considered system. In turn, CID is not only exhibit strong fluctuations, but also lies in a significantly different range depending on $L$. Here we used LZMA compression algorithm as it is one of the extensions of LZ77 used in Ref.~\onlinecite{clg_cid} and shows the best results among other standard techniques. When decreasing the system size the amplitude of $\mathcal{D}_{\rm xt}$ drops a bit for $\rho > 0.5$. This is due to the rearrangement of the partial contributions to dissimilarity from different scales $\mathcal{D}_{\rm xt}^k$ (see Fig.~\ref{fig:partial_lg} in Appendix A) and the fact that we do not consider $k = 0$ step. It should be mentioned, that to get smooth graph with zero value at $\rho = 0.5$ we still have to consider sufficient number of time steps such that the dynamics at the critical point stops.

The deeper analysis of the dissimilarity (see Fig.~\ref{fig:d_avg}) shows that the behaviour of $\mathcal{D}_{\rm xt}$ is similar to the fully spatial case $\mathcal{D}_{\rm x}$. Moreover, the slight differences occur only in active phase when the temporal part is high enough. Therefore, spatial patterns in digital representation provide most of meaningful information about the system. Nevertheless, having only the $\mathcal{D}_{\rm t}$ one can still define a phase transition analysing the derivative behaviour, as it was proposed in Ref.~\onlinecite{PNASComplexity}.

\subsection*{Patterns of discrete time crystal}

The main focus in our paper is on the classification of dynamical phases realized in driven non-equilibrium systems of qubits. More specifically, we are interested in exploring  the properties of a periodically driven quantum system hosting the discrete time crystal (DTC) state. As it was shown in the previous works\cite{DTC1,DTC2,DTC3} DTC is formed as a result of  breaking discrete time translation symmetry in a Floquet system with time-periodic Hamiltonian and retaining memory of the initial state indefinitely. Experimentally, it can be realized in systems of superconducting qubits \cite{DTC_IBM, DTC_Google, Khemani}, but other platforms, for instance quantum-optical systems \cite{Fedorov} are also promising candidates for observing this non-trivial state. The main signature of DTC, the spatio-temporal order can be detected and probed with the combination of temporal autocorrelators and qubit-qubit spatial correlation functions. As it was demonstrated theoretically \cite{Khemani} and experimentally\cite{DTC_IBM, DTC_Google}, these correlation functions can be defined using the results (bitstrings) of the numerous projective measurements performed in a number of consecutive Floquet cycles during the evolution. In this work, by using the developed spatio-temporal dissimilarity measure we show that one can perform a complete low-level characterization of DTC with a limited number of bitstrings, directly without calculating correlation functions.

    \begin{figure}
		\centering
		\includegraphics[width=0.45\textwidth]{fig5.png}
		\caption{\label{fig:DTCcircuitInUniversalGates} Example of DTC quantum circuit for 4-qubit system where single-qubit rotational gates ($X_\epsilon=R_x((1-\epsilon)\pi)$, $Z_{J}=R_z(J_i)$, $Z_h=R_z(h_i)$) and CNOT-gate were used.}
	\end{figure}

Time evolution of the quantum system needed to stabilize non-equilibrium discrete time crystal state is described by the Floquet unitary $U_{\rm F}$ \cite{DTC_Google}:
\begin{multline}
\label{U_F}
        U_{\rm F} =
        e^{ \scalebox{1.0}{$-\frac{i}{2} \left(1-\epsilon\right) \pi \sum\limits_{i} X_i$} }   \\
        \cdot e^{\scalebox{1.0}{$-\frac{i}{2} \sum\limits_{i<j} J_{ij} Z_i Z_{j} -\frac{i}{2} \sum\limits_{i} h_i Z_i $ }}
\end{multline}
which is defined in term of two unitaries.
The first one corresponds to the imperfect global spin-flip (imperfect $\pi$-pulse). Here $X_i$ is the Pauli X-gate acting on the $i$-th qubit, $\epsilon\in[0;0.5]$ is the spin-flip imperfection parameter. The second unitary represents a combination of a random nearest-neighbor Ising interactions, $J_{ij}$ and a random longitudinal field, $h_i$. Here $Z_i$ is the Pauli Z-gate on the $i$-th qubit. Following to Ref.~\onlinecite{DTC_Google} we use a trotterized time evolution to simulate the DTC with a quantum simulator.

In Fig.~\ref{fig:DTCcircuitInUniversalGates} we show a fragment of the quantum circuit that corresponds to one trotter cycle for a 4-qubit system driven with Floquet unitary described above. It comprises rotation single-qubit gates and entangled CNOT gates. To construct a complete quantum circuit that imitates periodically driven quantum system for a given number of time steps we randomly choose  nearest neighbors interaction parameters $J_{ij}$ and local magnetic fields $h_i$ from the ranges $\left[-0.75\pi;-0.25\pi\right]$ and $\left[-\pi;\pi\right]$, respectively. These parameters are fixed for each quantum circuit instance independent on its depth. Wherein the imperfect spin-flip, $\epsilon$ plays a role of a global parameter whose value allows to choose a particular non-equilibrium phase. To be more concrete, $\epsilon$ is varied within [0, 0.5]. In our work we perform simulations of the DTC on the 16-qubit one-dimensional system with open boundaries.

\begin{figure}[!t]
	\centering
	\includegraphics[width=0.45\textwidth]{fig6.pdf}
	\caption{Examples of bitstring evolution within the time window [0, 63] for particular circuits constructed with five representative values for the $\epsilon$ parameter. \label{bitstrings}}
\end{figure}

While the previous works devoted to the DTC were mainly focused on the analysis of different forms of the spin polarization, here we show that the non-trivial state can be fully described at a more basic level of bitstrings. From
Fig.~\ref{bitstrings} it follows that there are two different types of bitstring temporal evolution depending on the parameter of imperfect rotation. At small values of $\epsilon$ one can observe period-doubled oscillations between initial bitstring and its flipped counterpart. At some times the projective measurement shows slightly modified bitstrings that differ from the initial one by a few bits. The latter means that quantum state becomes slightly delocalized in the Hilbert space with respect to the initial trivial state. At the same time, there are nearly perfect oscillations at the level spin polarization \cite{Khemani} within the DTC phase, when each spin is fully flipped. Thus, although from the point of view of structure, the quantum state changes during the evolution due to disorder, nevertheless it reproduces the period-doubled dynamics in terms of the correlation functions local in space, which is a signature of DTC.

As $\epsilon$ increases, bit-string patterns are becoming more and more diverse over time, which means that different basis functions are involved and the delocalization of the quantum state increases. For $\epsilon = 0.3$ we observe a complete loss of information about the initial configuration.

\begin{figure}[!b]
	\centering
	\includegraphics[width=0.40\textwidth]{fig7.pdf}
	\caption{ Probability distribution of the Hamming distances $P_n(d)$ computed between the initial state and bitstrings measured at time steps 62 (blue) and 63 (orange) of the Floquet evolution. The histograms are obtained using 512 independent circuits.}
	\label{Hamming}
\end{figure}

To quantify the bitstring dynamics described above we employ the Hamming distance that measures the number of bits by which one bitstring differs from another and is defined as $d$. Firstly, the probability distribution of Hamming distances, $P_{n} (d)$ between the initial and time-evolved states is calculated. However, in contrast to  Ref.~\onlinecite{Khemani} for our analysis we use a limited number of measurements, which doesn't allow to reconstruct the complete probability of the quantum system at each time step. Dependence of $P_{n}(d)$ on the parameter of imperfect rotation presented in Fig.~\ref{Hamming} unveils two boundary cases. The first one is observed at small $\epsilon$ and represents a signature of the DTC phase according to the results of Ref.~\onlinecite{Khemani}. The probability distributions of the Hamming distances calculated between the initial state and two consecutive states obtained after a moderate number of time-evolution steps are localized at $d = 0$ and $d = 16$.  Another distinct regime is revealed for $\epsilon \ge 0.23$. The calculated $P(d)$ is centered at $d = 8$ for even and odd time steps, which indicates the thermal phase. While DTC and thermal regimes of dynamics can be identified with the probability distribution of the Hamming distance, the transition between them on the level of the calculated $P(d)$ is smooth and doesn't show any feature that can be associated to the critical point. The presented distributions are not symmetrical with respect to $d = \frac{L}{2}$, since the number of bitstrings we used is much smaller than the total size of the Hilbert space for 16-qubit system. In Appendix B we complete the picture of the non-equilibrium dynamics in the considered quantum system by analyzing the behaviour of the Shannon and von Neumann entropies.

As it was shown previously~\cite{Khemani}, the transition point between DTC and thermal phases can be defined with the Edwards-Anderson spin glass order parameter \cite{Anderson}, which represents a total spin-spin correlation function averaged over several time steps and different instances of disorder. In our work to explore a dynamical order and define critical value of $\epsilon$ we use the developed $\mathcal{D}_{\rm xt}$ measure that allows to define the critical point directly from the spatio-temporal representation of bitstrings. The key moment of our scheme is that at each time step of the quantum evolution we perform only one measurement of all the qubits and the resulting bitstrings obtained from several consecutive time steps are used to form qubit-time grid. For each value of $\epsilon$ the dissimilarity $\mathcal{D}_{\rm xt}$ was averaged over 512 realization of disorder (circuit instances) each of which assumes a random initial state. We have found that the minimal number of the time cycles required for calculating spatio-temporal dissimilarity in the case of DTC is equal to 16. Moreover, an increase in this number does not change the $\mathcal{D}_{\rm xt}$ estimates. In Appendix C we discuss different averaging schemes with a fixed disorder or a fixed initial state.

The resulting dissimilarities are given in Fig.~\ref{dissimilarity}. One can clearly distinguish between discrete time crystal and thermal phases by analyzing the behaviour of the calculated $\mathcal{D}_{\rm xt} (\epsilon)$. In the DTC region the dissimilarity demonstrates a rapid growth with increasing $\epsilon$, which can be explained by the fact that the target wave function becomes more delocalized in the Hilbert space. In turn, in the thermal phase $\mathcal{D}_{\rm xt}$ is independent on this parameter and shows nearly zero slope. The critical value of the imperfect rotation parameter can be found as a crossing point of the corresponding linear fits \cite{glass1} of the results deep inside the DTC and thermal phases. Dependence of the $\epsilon_c$ on the system's size presented in Fig.\ref{dissimilarity} bottom shows the stabilization of the critical value for $L \geq 16$. The obtained $\epsilon_c\approx0.16$ is in good agreement with theoretical estimations based on the Edwards-Anderson order parameter\cite{Khemani} and lies within the range of experimental estimates~\cite{DTC_Google}.

The results presented in Fig.\ref{dissimilarity} were obtained for random initial states of zero entanglement, which is a standard setup widely used in the literature devoted to the DTC problem. In Appendix D we discuss the simulations of the DTC phase with different entangled initial states.

As it was shown in the classical case for one-dimensional CLG the exploration of the non-equilibrium critical behaviour with spatio-temporal dissimilarity can be successfully performed using much less data about the system in comparison with other measures. Thus, in the quantum case it is instructive to compare the total number of measurements required for characterization of DTC with $\mathcal{D}_{\rm xt}$ and standard order parameter of the spin glass state used in previous works \cite{Khemani}. To get one point on the graph presented in Fig.~\ref{dissimilarity} we use 8192 measurements (number of time steps $\times$ number of disorders). In turn, to estimate the glass order parameter for the given $\epsilon$ in the real experiment the authors of Ref.~\onlinecite{DTC_Google} have performed 1600000 measurements which is mainly caused by noise and gates imperfections. In contrast, theoretical results on the correlation functions presented in the same work for the system of 16 qubits were obtained with 40000 measurements, which is still considerably larger than  those needed to calculate the dissimilarity.

\begin{figure}
	\centering
	\includegraphics[width=0.45\textwidth]{fig8.pdf}
	\caption{ (Top) Dependence of the spatio-temporal dissimilarity on the critical imperfect rotation parameter. Blue squares and red circles indicate the internal regions of the DTC and thermal phases, respectively. Dashed blue and red lines show the linear fittings for the corresponding points. Green dash-dotted line denotes the value of $\epsilon_c$ corresponding to the DTC-thermal transition. (Bottom) Dependence of the $\epsilon_c$ on the system size $L$. These results were obtained with 4096 ($L$=8), 6144 ($L$=12), 8192 ($L$=16), 4000 ($L$=20) and 4800 ($L$ =24) measurements. } \label{dissimilarity}
\end{figure}

\section*{Discussion and Outlook}
We propose a flexible and simple approach for diagnosing non-equilibrium dynamics in classical and quantum systems. It allows identifying different dynamical regimes including hidden order out of equilibrium and accurate detecting critical parameters. Taking a time-resolved digital representation of a dynamical process as input, our method quantifies the variety or dissimilarity of spatio-temporal patterns by representing it as a single number. Thus, distinguishing different parts of the same process or different dynamical processes boils down to comparison of the corresponding numbers which considerably facilitates and accelerates characterization of non-equilibrium systems and constructing phase diagrams.

By using the proposed spatio-temporal dissimilarity measure we have analyzed non-equilibrium dynamics of one-dimensional classical and quantum models. In classical case we show that  $\mathcal{D}_{\rm xt}$ reveals a hidden order in one-dimensional conserved lattice gas model for much smaller system's sizes than in the case of the compression measure. In turn, non-trivial dynamically driven discrete time crystal phase in a non-equilibrium quantum system can be fully described on the level of bitstrings without calculating standard correlation functions, which also substantially reduces the total number of the measurements.

In this work to extract information on a non-equilibrium one-dimensional system we explore two-dimensional spatio-temporal grid formed with digital representations of the system in question measured at different times. To consider higher dimensional physical systems, two different strategies can be adopted. The first one assumes using one-dimensional representation for spatial degrees of freedom of the high-dimensional system at each time step. As it was shown previously, in Ref.~\onlinecite{Dissimilarity} such a dimensional reduction in space still allows an accurate description of equilibrium phase transitions in magnetic systems, however, it assumes testing different numeration schemes to order spatial degree of freedom and to mimic the topology of the initial system in one-dimension representation. Another possibility is to use high-dimensional spatio-temporal filters when calculating the dissimilarity, $\mathcal{D}_{\rm xt}$.

The dissimilarity method supports processing digital representations of completely different dynamical processes. For instance, in this work we are concentrated on the characterization of the non-equilibrium quantum systems for which measurement's outputs are given in the binary format, bitstrings that consist of 0 and 1 (-1 and 1 in our calculations). Vector fields with floating-point type components are another important example of a digital representation that can be handled with dissimilarity\cite{PNASComplexity}, which paves the way to characterize the dynamics of off- and in-lattice models. From the perspective of real experiments one can calculate spatio-temporal dissimilarity directly for a stream of experimental images in RGB or other formats without mapping onto a mathematical model. This can be rationalized in non-scientific applications, for instance as a pre-processing scheme in automatic video content classification systems.

\section*{Acknowledgements}
The research funding from the Ministry of Science and Higher Education of the Russian Federation (Ural Federal University Program of Development within the Priority-2030 Program) is gratefully acknowledged. Quantum simulations were performed on the Uran supercomputer at the IMM UB RAS.

\section*{Appendix A: Partial dissimilarity for CLG}\label{AppDK}
\begin{figure*}[t]
\centering
\includegraphics[width=.85\textwidth]{fig9.png}
\caption{\label{fig:partial_lg} The behavior of $\mathcal{D}_{\rm xt}^k$ for configurations corresponding to the maximum $\mathcal{D}_{\rm xt}$ in adsorbing (green squares) and active (orange circles) phases for $L=8$ (a) and $L=256$ (b). (c) Dependency of $\mathcal{D}_{\rm xt}^0$ on the density of particles in the system for $L=8$ (red diamonds) and $L=256$ (purple triangles).}
\end{figure*}

As one can see from Fig.~\ref{fig:1d_gas}, the total $\mathcal{D}_{\rm xt}$ for the CLG problem exhibits a symmetric behaviour close to the one presented in Ref.~\onlinecite{clg_cid} for the thermodynamic limit. However, unlike the CID, which is just a single value, we can still define the phase of a considered configuration by analyzing its partial spatio-temporal dissimilarities. According to Fig.~\ref{fig:partial_lg}, configurations corresponding to the active phase have monotonously decreasing $\mathcal{D}_{\rm xt}^k$ profiles while the ones from adsorbing phase have more complex $\mathcal{D}_{\rm xt}^k$ behaviour. Moreover, the $k=0$ step can be used as an indicator of how close we are to the phase transition point.

\section*{Appendix B: Analysis of the DTC entanglement entropy}
Hamming distance analysis as well as dissimilarity consideration presented in the main text are based on the results of a limited number of measurements, which is clearly not enough to reconstruct the true probabilities of the basis functions of the quantum system at each cycle of the Floquet unitary evolution. To complete the picture of non-equilibrium dynamics of the quantum system in question, we analyze time-dependent behaviour of the information theory measures, Shannon entropy and von Neumann entanglement entropy. For a given quantum state the Shannon entropy is defined as
\begin{eqnarray}
H(X) = -\sum_{n} p_n \log_2 p_n,
\end{eqnarray}
where $p_n$ is the probability of the $n$th basis function. The obtained time dependence of $H(X)$ are visualized in Fig.\ref{fig:EntangledStates} left. One can see that the Shannon entropy calculated at different $\epsilon$ demonstrates a fast growth and saturation for time steps $t > 20$. It corresponds to a delocalization of the quantum state in the Hilbert space. For instance, in the case of $\epsilon = 0.03$ such a delocalization can be qualified as moderate. However, from the Shannon entropy alone calculated in single $\sigma^z$ basis it is not clear whether the initial trivial state transforms into an entangled one during the evolution or it remains trivial. To clarify this point we have calculated the von Neumann entropy, $S_N = -{\rm Tr} \rho_{\rm A} \log_2 \rho_{\rm A}$, where $\rho_{\rm A} = {\rm Tr}_{\rm B} \rho_{\rm AB}$ is the reduced density matrix for a half-system biparition into parts A and B. It is found that for $\epsilon = 0.03$ the quantum system remains non-entangled during the whole evolution time. For larger values of $\epsilon$ non-equilibrium dynamics creates entanglement between spatial subsystems. The thermal phase is characterized by a very rapid growth of the von Neumann entropy with time cycles. For $\epsilon=0.30$ the Shannon entropy and the von Neumann entropy reach the values of 15.35 and 7.22, respectively. It means that at this $\epsilon$ the wave function of the system is close to the Haar random state that is characterized by the ''Page value'' $\log_2 2^{\frac{L}{2}} - \frac{1}{2}$ for von Neumann entropy \cite{Page} and $\log_2 2^{L} - 1 + \gamma$ for the Shannon entropy (here $\gamma$ is Euler constant) \cite{supremacy}. Similar to the Shannon entropy we cannot discriminate DTC-thermal transition analyzing the calculated time dependencies of $S_N$.

	\begin{figure}
		\centering
		\includegraphics[width=0.45\textwidth]{fig10.pdf}
		\caption{\label{fig:EntangledStates} Time evolution of the Shannon entropy $H(X)$ (left) and entanglement entropy $S_N$ (right) calculated at different values of the $\epsilon$ parameter.}
	\end{figure}

\section*{Appendix C: Averaging DTC data}

The results of the dissimilarity calculations used for constructing the phase diagram, Fig.\ref{dissimilarity} in the main text are averaged over initial states as well as disorder instances. In this section we perform a separate consideration of these two kinds of averaging. Fig.\ref{average} A gives the spatio-temporal dissimilarities calculated with fixed initial states of three types and only averaged over disorders. More specifically we consider the polarized state $\ket{\Psi_{\rm polar}} = \ket{0000..00}$, random wave function $\Psi_{\rm R}$ for which the state of each qubit is randomly chosen between 0 and 1, and the N\'eel state $\ket{\Psi_{\rm Neel}} = \ket{0101..01}$  One can see that each initial state results in a distinct dependence $\mathcal{D}_{\rm xt}(\epsilon)$. The critical $\epsilon$ parameters are slightly different within this series of experiments and vary from 0.15 to 0.20.

 	\begin{figure}
		\centering
		\includegraphics[width=0.45\textwidth]{fig11.pdf}
		\caption{\label{average} (A) Spatio-temporal dissimilarities calculated with different initial wave functions: N\'eel (green), random (orange) and polarized (blue) states. In each case $\mathcal{D}_{\rm xt}$ is averaged over 512 disorder instances. (B) Examples of the $\mathcal{D}_{\rm xt} (\epsilon)$ functions calculated for three different instances of disorder. Each curve is averaged over 512 trivial initial states.}
	\end{figure}
In turn, comparison of the dissimilarities obtained with different disorders and averaged over 1024 random initial states are presented in Fig.\ref{average} B. Variation of disorder gives rise to a slight difference in behaviour between the calculated dissimilarities, which also affect the estimation of the critical point in each case. Based on these results we can conclude that the spatio-temporal dissimilarity measure demonstrates a rather robust behaviour for different choices of the averaging data procedure.

\section*{Appendix D: Initial entangled states for DTC}
The use of the time-dependent spin polarization of qubits as indicator of the dynamical state of the system is not well defined if the initial state of the quantum system is non-trivial and characterized by non-zero entanglement. In this case the initial zero averaged magnetization of individual qubits can challenge probing dynamical properties of the quantum system. A possible solution of such a problem was recently demonstrated in Ref.~\onlinecite{DTC_Google}, where the authors used an ancilla qubit to confirm the DTC-MBL dynamics of the quantum systems initially prepared in the highly-entangled Haar-random state with zero magnetization.

Our dissimilarity approach is suitable for extracting useful information on the unitary dynamics of the zero-magnetization entangled quantum states without introducing additional ancilla-qubits. To demonstrate this we use Dicke and Schr\"odinger cat
wave functions of different entanglement as initial states for unitary dynamics with the DTC protocol. The Dicke quantum state is defined as
\begin{eqnarray}
|\Psi_D \rangle = \frac{1}{\sqrt{C^N_D}} \sum_{j} P_{j}(|0\rangle^{\otimes N-D} \otimes |1\rangle^{\otimes D}),
\label{eq:Dicke_wf}
\end{eqnarray}
where the sum goes over all possible permutations of qubits. The parameter $D$ controls the number of 1's in each basis function contributing to $\Psi_D$.

In turn, the Schr\"odinger cat state is controlled with the angle $\theta$,
\begin{eqnarray}
|\Psi_{\theta} \rangle = {\rm cos}(\frac{\theta}{2}) |0\rangle^{\otimes N} + {\rm sin}(\frac{\theta}{2}) |1\rangle^{\otimes N}.
\end{eqnarray}
For $\theta = 0$ or $\pi$ one deals with trivial states discussed in the main text.
Comparing the structures of the Dicke and cat states one can consider them as two limiting examples of the localized wave functions. Such a difference in the basis functions contributing to the quantum state probability manifests itself not only at the level of entanglement entropy, but also in the partial dissimilarities calculated in Ref.\onlinecite{Dissimilarity}.

As representative examples of the entangled quantum states that are characterized by zero local qubit polarization  we consider the Dicke state with $D= 8$ and Schr\"odinger cat with $\theta = \frac{\pi}{2}$ that are used as initial ones in the DTC protocol. Fig.\ref{entangled} illustrates the obtained $\mathcal{D}_{\rm xt}$ as functions of $\epsilon$. In contrast to non-equilibrium dynamics initialized with different trivial wave functions (presented in the main text), using the initial entangled states results in a more diverse behaviour of the $\mathcal{D}_{\rm xt} (\epsilon)$ functions  within the DTC part of the phase diagram.

\begin{figure}[!hb]
	\centering
	\includegraphics[width=0.45\textwidth]{fig12.pdf}
	\caption{\label{entangled} Comparison of the spatio-temporal dissimilarities calculated with different initial entangled states within DTC protocol. Each curve is obtained with averaging over 512 realizations of disorder. }
\end{figure}

The main peculiarity of the spatio-temporal dissimilarity for the Schr\"odinger cat state is a large standard deviation at small values of the imperfect rotation parameter. The deviation gradually decreases as approaching to $\epsilon$ = 0.25. Such strong fluctuations of the calculated dissimilarities are mainly due to smallness of the space-time grid (16$\times$16) used in our study. As the result, each of 1024 circuits produces the $B^{0}_{L \times T}$ matrix that contains a distinct sequence of 0000..00 and 1111..11 bitstrings. It can be rationalized by considering the bitstings examples presented in Fig.\ref{cat_bitstrings}. For $\epsilon = 0.01$ one could randomly choose a fragment of 16 consecutive bitstrings. Each of such fragments would be characterized by different number of ''0'' and ''1'' bits, which explains fluctuations of the polarization and dissimilarity.

\begin{figure}[!h]
	\centering
	\includegraphics[width=0.45\textwidth]{fig13.pdf}
	\caption{Examples of bit-string sequence obtained in a time window [0, 63] with Schr\"odinger cat state $\Psi_{\theta=\pi/2}$ as initial state of the qubits system. \label{cat_bitstrings}}
\end{figure}

In turn, the bitstrings corresponding to the Dicke state constitutes of the same number of the ''0'' and ''1'' bits, which together with the large number of the basis functions that contribute to $\Psi_{D}$ (Eq.\ref{eq:Dicke_wf}) results in nearly constant profile of the mean spatio-temporal dissimilarity and the corresponding standard deviation within the working area. A zoom view (inset of Fig.\ref{entangled}) revealed an oscillating dissimilarity that, nevertheless, reproduces the main trends of the $\mathcal{D}_{\rm xt} (\epsilon)$ function calculated for the trivial initial state.


\section{Introduction}
\vspace{-0.2cm}
Accurate recognition of elderly and disordered speech remains a challenging task \cite{yu2018development,xiong2020source,liu2021recent,ye2021development,deng21d_interspeech,geng2022speaker, yue2022acoustic}.
Speech disorders such as dysarthria are associated with neuro-motor conditions
and often found among elderly adults experiencing neurocognitive disorders such as Alzheimer’s disease (AD). ASR technologies tailored for their needs not only improve their quality of life, but also support large scale automatic early diagnosis of neurocognitive impairment
\cite{ferri2005global,vipperla2010ageing,rudzicz2014speech,zhou2016speech,mirheidari2019dementia,ye2021development}
\par
Elderly and disordered speech bring many challenges to current deep learning based ASR technologies predominantly targeting non-aged, healthy adult users. First, a large mismatch between such data and non-aged, healthy adult voice is often observed.
Second, the co-occurring disabilities and mobility limitations often found among elderly and disordered speakers lead to the difficulty in collecting large quantities of such data that are essential for data intensive ASR system development. To this end, a series of recent researches have been conducted to transfer knowledge from normal speech that is easier to collect and ASR systems trained on such out-of-domain data. These include, but not limited to, speed and temporal perturbation based data augmentation \cite{geng2022speaker, ye2021development}, adversarial data augmentation \cite{jin2022personalized}, auto-encoder features \cite{xie21b_interspeech, bhat22_interspeech}, domain adaptation of hybrid\cite{takashima2020two,yu2018development,xiong2020source,liu2021recent,ye2021development,deng21d_interspeech}, and end-to-end ASR systems \cite{wang2021improved, shor19_interspeech, green21_interspeech, tobin2022personalized, wang2022conformer}.

\par
An alternative solution to address the above data scarcity problem is to use self-supervised learning (SSL) models \cite{baevski2020wav2vec, chen2022wavlm, hsu2021hubert} that are pre-trained on large quantities of unlabelled data. The associated neural speech representations produced by these pre-trained ASR systems are also inherently robust to domain mismatch. Although they have been successfully applied to a range of normal speech processing tasks including speech recognition \cite{baevski2020wav2vec, conneau2020unsupervised, chen2022wavlm, hsu2021hubert}, speech emotion recognition \cite{pepino2021emotion} and speaker recognition \cite{vaessen2022fine}, very limited researches on SSL pre-trained models for disordered and elderly speech have been conducted \cite{violeta2022investigating, baskar2022speaker, hernandez2022cross}.
Among these, wav2vec2.0 and WavLM models were applied to Japanese electrolaryngeal and English dysarthric speech in \cite{violeta2022investigating}, where an overall WER of 51.8\% was reported on the benchmark UASpeech task. Speaker adaptation of wav2vec2.0 models using fMLLR and x-vectors for dysarthric speech recognition was investigated in \cite{baskar2022speaker}. Cross-lingual SSL based dysarthric speech recognition was studied in \cite{hernandez2022cross}.

\par
Efforts on applying SSL pre-trained speech models to dysarthric and elderly speech are also confronted with the same data scarcity problem. Direct fine-tuning of a large number of model parameters on limited impaired or elderly speech data rapidly leads to overfitting and poor generalisation. This issue is further exasperated when the limited training data provides insufficient coverage of words in the test data. For example, approximately 39\% of words in the benchmark UASpeech test set do not occur in the training data. This was found to produce a larger disparity in performance between seen and unseen words in dysarthric speech, and between impaired speakers with high and very low intelligibility, for end-to-end ASR systems including normal speech pre-trained models, than hybrid ASR systems that are  trained using phonetic targets \cite{liu2021recent,wang2021improved, violeta2022investigating}.
\par
To this end, this paper explores a series of approaches to integrate SSL pre-trained wav2vec2.0 models into hybrid TDNN \cite{peddinti2015time} and Conformer \cite{gulati2020conformer} ASR systems for dysarthric and  elderly speech recognition in order to exploit the diversity and complementarity among them, including
a) input feature fusion between standard acoustic front-ends and wav2vec2.0 speech representations; b) frame-level joint decoding \cite{swietojanski2013revisiting, liu2021recent, hu2022exploitinguti} of TDNN systems separately trained using standard acoustic features alone and with additional wav2vec2.0 representations; and c) cross system multi-pass decoding \cite{cui2022two} involving the TDNN or Conformer system N-best outputs to be rescored using domain adapted wav2vec2.0 models. Finally, wav2vec2.0 representations are utilized to facilitate neural network based acoustic-to-articulatory (A2A) inversion \cite{hu2022exploitinguti}. The resulting articulatory movement features are used to construct multi-modal dysarthric and elderly ASR systems.
\par
Experiments conducted on the UASpeech \cite{kim2008dysarthric} dysarthric and DementiaBank Pitt \cite{becker1994natural} elderly speech corpora suggest TDNN and Conformer ASR systems constructed by integrating domain adapted wav2vec2.0 models consistently outperform the standalone wav2vec2.0 models by statistically significant WER reductions of 8.22\% and 3.43\% absolute (26.71\% and 15.88\% relative) on the two tasks respectively.
The lowest published WERs of 22.56\% (52.53\% on very low intelligibility, 39.09\% on unseen words) and 18.17\% are obtained
on the UASpeech test set of 16 dysarthric speakers and the DementiaBank Pitt test set.

\par
The main contributions of the paper are summarized below:\\
\textbf{1)} To the best of our knowledge, this paper presents the first work to systematically investigate the generalization capability of SSL pre-trained models for dysarthric and elderly speech recognition tasks. In contrast, related previous works \cite{violeta2022investigating, baskar2022speaker, hernandez2022cross} analyzed their performance disparity across different speech pathology severity subsets but not their generalization to unseen words or sentences. \\
\textbf{2)} To address such generalization issues, this paper proposes novel approaches to integrate SSL pre-trained models into in-domain data constructed ASR systems. Drawing strengths from both, the lowest published WER of 22.56\% (52.53\% on very low intelligibility group, 39.09\% on unseen words) is obtained on the UASpeech test set of 16 dysarthric speakers. In contrast, prior researches focused on using stand-alone domain adapted pre-trained models \cite{violeta2022investigating, baskar2022speaker, hernandez2022cross}. \\
\textbf{3)} This work also presents the first use of SSL based speech representations to cross-domain A2A inversion. In contrast, prior researches on cross-domain A2A inversion were conducted predominantly using non-SSL based domain adaptation methods \cite{hu2022exploiting, hu2022exploitinguti}.
\vspace{-0.6cm}
\section{Pre-trained Wav2vec2.0 ASR Model}
\vspace{-0.2cm}
\subsection{model architecture}
\vspace{-0.2cm}
The wav2vec2.0 model is an SSL pre-trained model that jointly learns latent contextualized speech representations and an inventory of discretized speech units. It consists of three components (left part of Fig. 1), including a \textbf{speech feature encoder $f(\cdot)$}
containing multi-layer 1D CNN based convolution blocks, a \textbf{contextual transformer network $g(\cdot)$} and a \textbf{quantization module $h(\cdot)$} consisting of 2 codebooks with 320 entries each. The input of wav2vec2.0 model is raw speech audio $\mathcal{X}$, which is first encoded into a sequence of representations $\mathcal{Z}=f(\mathcal{X})$ with a stride of 20 ms and a receptive field of 25 ms by the CNN-based speech feature encoder. These CNN extracted features are then fed into transformer based contextual encoder to produce contextualized speech representations $\mathcal{C}=g(\mathcal{Z})$. A logit $l \in \mathcal{R}^{2 \times 320}$ is calculated from each vector $z_i \in \mathcal{Z}$ before fed into Gumbel-Softmax \cite{jang2016categorical} operators to choose the best code among all entries of each codebook. The selected codes are concatenated together before a linear transformation is applied to produce the quantized vector $q_i$ during pretraining.
\vspace{-0.4cm}
\subsection{Model pre-training and fine-tuning}
\vspace{-0.2cm}
Wav2vec2.0 model is pre-trained using a weighted interpolation between the contrastive and diversity loss functions. For the speech representation of time step $t$ that is randomly masked, the contrastive loss based training aims to reproduce the quantized vector $ q_t$ for the masked frame at the output of contextual transformer encoder, while discriminating against $K$ ``distractors''  representations ${\tilde q} \in {\tilde Q}$ that are randomly sampled from the other masked time steps within the same speech utterance. The contrastive loss function is defined as
\begin{small}
\begin{equation}
    \setlength{\abovedisplayskip}{2pt}
    \setlength{\belowdisplayskip}{2pt}
    L_{m} = - \log \frac{\exp({\rm sim}(c_t, q_t)/ \kappa)}{\sum_{\widetilde q \in \widetilde Q} \exp({\rm sim}(c_t, \widetilde q)/ \kappa)}
\end{equation}
\end{small}
\par
\noindent
where ${\rm sim}(c_t, q_t)$ is the cosine similarity between the contextualized vector $c_t$ and the quantized representation $q_t$. $\kappa$ is the temperature and set $\widetilde Q$ consists of both the target quantization $q_t$ and $K$ distractors. The entropy based diversity loss is designed to regularize the use of the quantized codebook representations.
\par
During the fine-tuning stage, a randomly initialized linear projection layer is added on top of contextual transformer network $g(\cdot)$ into $V$ classes representing the vocabulary entries of the task considered. The model then is fine-tuned with standard Connectionist Temporal Classification (CTC) loss \cite{graves2006connectionist} in a supervised fashion.
\vspace{-0.4cm}
\subsection{Speech Feature Extraction}
\vspace{-0.2cm}
The SSL pre-trained wav2vec2.0 model can be used as a standalone speech recognition system after task fine-tuning. Alternatively their contextual speech representations can be incorporated into various back-end target domain data trained ASR systems via feature fusion based on, e.g. TDNN or Conformer architectures, as considered in this paper. To this end, a ``Bottleneck module'' (circled in blue, top middle in Fig. 1) is introduced immediately above the highest positioned contextual transformer block to produce more compact speech representations\footnote{Alternative approaches of extracting bottleneck features, e.g. either from the CNN feature encoder outputs, contextual transformer network outputs, or a weighted average of the hidden layer outputs of each transformer block as suggested in \cite{pepino2021emotion} led to performance degradation.}. The ``Bottleneck module'' contains a stack of four interleaving convolutional and feedforward layers: the first 1D transposed de-convolution CNN layer is used to change the stride length to 10ms to allow a frame rate synchronization with that of the back-end ASR systems; a fully connected (FC) block, which consists of a linear layer, rectified linear unit (ReLU) activation and dropout module, is used to change the dimension of extracted features from 1024 to 256\footnote{256 empirically outperforms other \#dim, e.g. 128, 512 or 768}. This is followed by a CNN layer and a final FC block which are used  to revert the stride back to 20ms and restore the dimensionality to 1024. The final wav2vec2.0 speech representations are extracted from the first FC layer outputs and concatenated with the standard acoustic front-ends via feature fusion in hybrid TDNN or Conformer E2E ASR systems ((a) or (b) in Fig. 1).
\vspace{-0.3cm}
\section{Wav2vec2.0 features based A2A inversion}
\vspace{-0.2cm}
In this paper, the out-of-domain non-aged normal speech of the TaL dataset \cite{ribeiro2021tal} containing parallel ultrasound, video data collected from 81 native speakers of English is used to construct mixture density networks (MDN) based neural A2A inversion models \cite{hu2022exploiting}. Due to the large acoustic domain mismatch, a direct cross-domain application of the A2A inversion model trained on normal, non-aged acoustic-articulatory parallel data to the elderly and dysarthric speech is problematic, as shown in the previous researches on cross-domain A2A inversion \cite{hu2022exploiting,hu2022exploitinguti}. To this end, the large acoustic domain mismatch is minimized using either baseline multi-level adaptive networks (MLAN) \cite{liu2021recent, hu2022exploiting, hu2022exploitinguti} or the above domain fine-tuned SSL wav2vec2.0 speech representations of Section 2.3.
\par
The SSL based domain adaptation includes the following steps: 1) The pre-trained wav2vec2.0 model is fine-tuned on both out-of-domain normal speech and in-domain dysarthric or elderly speech simultaneously; 2) the resulting fine-tuned model is then used to produce more domain invariant speech representations that will exhibit smaller mismatch between these types of data. The resulting domain adapted speech representations are used in A2A inversion model training and articulatory feature generation for dysarthric or elderly audio data. More details of the MDN based A2A inversion models can be found in \cite{hu2022exploiting}. The articulatory movement features generated from dysarthric or elderly speech are concatenated with the standard acoustic features via feature fusion to construct multi-modal TDNN or Conformer systems ((c) or (d) in Fig. 1).
\vspace{-0.3cm}
\section{integration with TDNN/Conformer system}
\vspace{-0.2cm}
In order to exploit the diversity and complementarity among the SSL pre-trained wav2vec2.0 model and back-end TDNN or Conformer ASR systems, and to improve the generalization performance on unseen words as well as dysarthric speech of very low intelligibility, a series of integration approaches are explored in this paper.
\par
\noindent
\textbf{1) input feature fusion} is used to fuse standard acoustic frontends and domain adapted wav2vec2.0 speech representations via either feature concatenation before fed into the TDNN/Conformer system input layer, {or optionally feature fusion at the TDNN hidden layer (connections (a) and (b) in Fig. 1, centre left).
\par
\noindent
\textbf{2) time-synchronous frame-level joint decoding} is used to combine two or more hybrid TDNN systems separately trained using either standard acoustic features alone, or with additional domain fine-tuned wav2vec2.0 representations, by a frame level linear interpolation of system specific acoustic scores (green dotted box (e) in Fig. 1, bottom right). Let $P_{a+w}(y_t|X_1^{1:T})$ and $ P_{a}(y_t|X_2^{1:T})$ denote the frame level acoustic scores at $t$-th frame of TDNN systems constructed with or without additional wav2vec2.0 representations. Then the final frame-level score is obtained by
\begin{small}
\begin{equation}
    \vspace{-0.1cm}
    \setlength{\abovedisplayskip}{2pt}
    \setlength{\belowdisplayskip}{2pt}
    {\hat P(y_t|X_1^{1:T}, X_2^{1:T})} = \alpha P_{a+w}(y_t|X_1^{1:T}) + \beta P_{a}(y_t|X_2^{1:T})
\end{equation}
\end{small}
\par
\noindent
where $\alpha$ and $\beta$ are the weights assigned to each system.
\par
\noindent
\textbf{3) cross system multi-pass decoding} involves a first decoding pass TDNN or Conformer system produced N-best outputs to be rescored using domain adapted wav2vec2.0 models.  Consider a sequence of acoustic features $X$ and its corresponding TDNN or Conformer (cfm) decoded N-best recognition hypotheses in \{$ Y_1, Y_2, ..., Y_n$\}. Let ${\textbf{s}}^{w2v} = [s_1^{w2v}, s_2^{w2v}, ..., s_n^{w2v}]^T$ and ${\textbf{s}}^{tdnn/cfm} = [s_1^{tdnn/cfm}, s_2^{tdnn/cfm}, ..., s_n^{tdnn/cfm}]^T$ denote the N-best score vectors produced by the first pass TDNN/Conformer and the second pass wav2vec2.0 models. The final 1-best output produced in the second pass wav2vec2.0 rescoring is obtained by interpolating the system specific N-best entry level scores as
\begin{small}
\begin{equation}
    \setlength{\abovedisplayskip}{2pt}
    \setlength{\belowdisplayskip}{2pt}
    {\hat Y}_{best} = \mathop{\arg}\min_{i} \{\alpha s_i^{w2v} + \beta s_i^{tdnn/cfm}\}
    \vspace{-0.1cm}
\end{equation}
\end{small}
\par
\noindent
where $\alpha$ and $\beta$ are the weights assigned to the second pass wav2vec2.0 and the first pass TDNN or Conformer systems.
An example of two-pass rescoring based TDNN/Conformer and wav2vec2.0 system combination is shown in Fig. 1 (top right, red box (f)).

\begin{table*}[htbp]
\centering
\caption{The performance of wav2vec2.0 based ASR,  TDNN/Conformer based systems constructed without or with wav2vec2.0 (w2v) representations and optionally using the cross-domain inverted UTI articulatory features on 1) the UASpeech test set; AND 2) the DementiaBank Pitt development (Dev) and evaluation (Eval) sets; ``VL/L/M/H'' refer to intelligibility subgroups. ``INV'' and ``PAR'' refer to clinical investigator and elderly participant. ``ada. w2v'' and ``ada. BN'' are the abbreviations of cross-domain adapted wav2vec2.0 representations and bottleneck features from MLAN\cite{hu2022exploitinguti} respectively. ``+'' and ``$\rightarrow$'' stand for frame-level joint decoding and multi-pass decoding. $\dag$ and $\ast$ denote statistical significant (MAPSSWE,$\alpha$ = 0.05) differences obtained against the baseline wav2vec2.0 and Conformer systems (Sys. 4,12).}
\scalebox{0.77}{
\begin{tabular}{c|c|c|c|cc|cccc|c|cc|cc|c}
\hline\hline
\multirow{3}{*}{Sys} & \multirow{3}{*}{model}  & \multirow{3}{*}{\begin{tabular}[c]{@{}c@{}}acoustic \\ feature\end{tabular}} & \multirow{3}{*}{A2A input} & \multicolumn{7}{c|}{UASpeech WER (\%)}                                                                                & \multicolumn{5}{c}{DBank WER(\%)}                                                                                \\
\cline{5-16}
                     &                         &                                                                              &                            & \multirow{2}{*}{unseen} & \multirow{2}{*}{seen} & \multirow{2}{*}{VL}   & \multirow{2}{*}{L}    & \multirow{2}{*}{M}    & \multirow{2}{*}{H}    & \multirow{2}{*}{All} & \multicolumn{2}{c|}{Dev.}                     & \multicolumn{2}{c|}{Eval.}                    & \multirow{2}{*}{All}  \\
\cline{12-15}
                     &                         &                                                                              &                            &                         &                       &                       &                       &                       &                       &                       & PAR.                  & INV.                  & PAR.                  & INV.                  &                        \\
\hline\hline
$1$                  & \multirow{2}{*}{TDNN-f} & \multirow{2}{*}{FBK}                                                         & \ding{55}                            & 51.62                   & 16.98                 & 62.53                 & 31.92                 & 23.12                 & 13.67                 & 30.56                 & 47.93                 & 19.91                 & 36.66                 & 19.76                 & 33.80                  \\
$2$                  &                         &                                                                              & ada. BN                    & 54.95                   & 17.64                 & 63.35                 & 33.14                 & 25.76                 & 15.78                 & 32.27                 & 45.82                 & 19.21                 & 34.89                 & 18.42                 & 32.35                  \\
\hline
$3$                  & +LHUC                   & FBK                                                                          &   \ding{55}                         & 47.73                   & 15.49                 & 61.12                 & 28.95                 & 19.08                 & 11.94                 & 28.13                 & 45.49                 & 19.26                 & 35.44                 & 18.42                 & 32.33                  \\
\hline\hline
$4$                  & wav2vec2                & raw audio                                                                    &  \ding{55}                          & 55.75                   & 14.67                 & 62.76                 & 36.91                 & 25.24                 & 9.11                  & 30.78                 & 29.71                 & 14.29                 & 21.27                 & 15.32                 & 21.60                  \\
\hline
$5$                  & \multirow{2}{*}{TDNN-f} & \multirow{2}{*}{FBK+w2v}                                                     & \ding{55}                            & 55.01          & 14.65                 & 60.59$^\dag$ & 34.81$^\dag$ & 23.80$^\dag$ & 12.01                 & 30.48        & 27.21$^\dag$ & 13.20$^\dag$ & 19.13$^\dag$ & 12.87                 & 19.73$^\dag$  \\
$6$                  &                         &                                                                              & ada. w2v                   & \textbf{53.90$^\dag$}   & \textbf{14.17}        & \textbf{60.07$^\dag$} & \textbf{33.58$^\dag$} & \textbf{23.49$^\dag$} & 11.33                 & \textbf{29.75$^\dag$} & \textbf{26.43$^\dag$} & \textbf{12.77$^\dag$} & \textbf{18.29$^\dag$} & \textbf{12.10$^\dag$} & \textbf{19.08$^\dag$}  \\
\hline
$7$                  & Sys.3 + 5               & -                                                                            & -                          & 43.77$^\dag$   & 11.72$^\dag$ & 53.19$^\dag$ & 26.33$^\dag$ & 16.47$^\dag$ & 9.02                  & 24.29$^\dag$ & 26.98$^\dag$ & 12.82$^\dag$ & 18.06$^\dag$ & 12.21                 & 19.29$^\dag$  \\
$8$                  & Sys.3 + 6               & -                                                                            & -                          & 43.68$^\dag$   & 11.71$^\dag$ & 53.39$^\dag$ & 26.17$^\dag$ & 16.57$^\dag$ & 8.84                  & 24.25$^\dag$ & 26.11$^\dag$ & 12.60$^\dag$ & 17.80$^\dag$ & 11.32$^\dag$ & 18.78$^\dag$  \\
$9$                  & Sys.3+5+6               & -                                                                            & -                          & \textbf{43.54$^\dag$}   & \textbf{11.68$^\dag$} & \textbf{53.26$^\dag$} & \textbf{26.13$^\dag$} & \textbf{16.45$^\dag$} & 8.81                  & \textbf{24.18$^\dag$} & \textbf{26.05$^\dag$} & \textbf{12.52$^\dag$} & \textbf{17.59$^\dag$} & \textbf{11.43$^\dag$} & \textbf{18.69$^\dag$}  \\
\hline
$10$                 & Sys.3$\rightarrow$4     & -                                                                            & -                          & \textbf{40.06$^\dag$}   & \textbf{11.72$^\dag$} & \textbf{52.53$^\dag$} & \textbf{25.00$^\dag$} & \textbf{14.80$^\dag$} & \textbf{7.08$^\dag$}  & \textbf{22.83$^\dag$} & 34.17                 & 15.41                 & 25.49                 & 15.98                 & 24.55                  \\
$11$                 & Sys.9$\rightarrow$4     & -                                                                            & -                          & \textbf{39.09$^\dag$}   & \textbf{11.90$^\dag$} & \textbf{53.12$^\dag$} & \textbf{25.03$^\dag$} & \textbf{14.04$^\dag$} & \textbf{6.32$^\dag$}  & \textbf{22.56$^\dag$} & \textbf{25.27$^\dag$} & \textbf{12.28$^\dag$} & \textbf{16.84$^\dag$} & \textbf{11.54$^\dag$} & \textbf{18.17$^\dag$}  \\
\hline\hline
$12$                 & \multirow{2}{*}{Conformer}               & FBK                                                                          & \ding{55}                           & 99.30                   & 18.24                 & 66.77                 & 49.39                 & 46.47                 & 42.02                 & 50.03                 & 48.71                 & 20.97                 & 36.93                 & 19.42                 & 34.57                  \\
$13$                 &               & FBK+w2v                                                                      & ada. w2v                   & \textbf{72.80$^\ast$}   & \textbf{14.87$^\ast$} & \textbf{64.24$^\ast$} & \textbf{43.09$^\ast$} & \textbf{36.41$^\ast$} & \textbf{17.28$^\ast$} & \textbf{37.59$^\ast$} & \textbf{28.38$^\ast$} & \textbf{14.53$^\ast$} & \textbf{19.40$^\ast$} & \textbf{13.10$^\ast$} & \textbf{20.79$^\ast$}  \\
\hline
$14$                 & Sys.13$\rightarrow$4    & -                                                                            & -                          & \textbf{71.73$^\ast$}   & \textbf{14.50$^\ast$} & \textbf{63.46$^\ast$} & \textbf{42.75$^\ast$} & \textbf{35.69$^\ast$} & \textbf{16.53$^\ast$} & \textbf{36.95$^\ast$} & \textbf{27.89$^\ast$} & \textbf{14.28$^\ast$} & \textbf{19.15$^\ast$} & \textbf{13.10$^\ast$} & \textbf{20.45$^\ast$}  \\
\hline\hline
\end{tabular}
}
\vspace{-0.6cm}
\end{table*}
\vspace{-0.3cm}
\section{Experiments}
\vspace{-0.2cm}
\noindent
\textbf{Task Description:} The UASpeech corpus is the largest publicly available dysarthric speech corpus \cite{kim2008dysarthric} consisting of 16 dysarthric and 13 control speakers, and contains 155 common words and 300 uncommon words. The entire corpus is further divided into 3 subset blocks per speaker. The same set of 155 common words is used in all three blocks, while the uncommon words in each block are different. The data from Block 1 and 3 of all the 29 speakers are used as the training set, while the data of Block 2 collected from all the 16 dysarthric speakers serves as the test set. After removing excessive silence and speed perturbation based data augmentation \cite{geng2020investigation}, a total of 130.1 hours of audio data is used as the training set, while 9 hours of speech is used for performance evaluation. The DementiaBank Pitt \cite{becker1994natural} corpus contains about 33-hour audio recorded over interviews between the 292 elderly participants and clinical investigators.
After silence stripping and data augmentation, the duration of training data is increased to 58.9 hours, while the development and evaluation sets contain 2.5 hours and 0.6 hours of audio respectively. The test data word coverage rates of UASpeech and DementiaBank are 61\% and 98.7\%.
\\
\noindent
\textbf{Experimental Setup:} The pre-trained wav2vec2.0 models on UASpeech and DementiaBank corpora are the Large Wav2Vec2.0 model\footnote{https://huggingface.co/facebook/wav2vec2-large-960h-lv60} and Conformer based Wav2Vec2.0 with relative position embeddings\footnote{https://huggingface.co/facebook/wav2vec2-conformer-rel-pos-large-960h-ft} respectively which are pretrained with 60k hours of Libri-light data and fine-tuned with 960 hours of Librispeech data. During domain adaptation to UASpeech data, the wav2vec2.0 model is first fine-tuned with the non-augmented data to update all parameters except those of speech feature encoder, then further fine-tuned with the 130.1h augmented data to only update the parameters of the first 12 transformer blocks to ensure generalization to unseen words. The numbers of training epochs for the two fine-tuning stages are both 20. For the DementiaBank data, all the wav2vec2.0 parameters except those of the speech feature encoder are fine-tuned for 30 epochs on the 58.9h augmented data. The initial learning rate is 1e-5 and then linearly decays. As section 2.3 mentioned, the dimensionality of bottleneck wav2vec2.0 speech representations is 256 and that of ultrasound based articulatory features is 144. Experiment setups of back-end TDNN and Conformer systems follow \cite{jin2022personalized}. \\
\noindent
\textbf{Results on UASpeech:} \textbf{1)}
Although there is no significant overall WER difference between the baseline TDNN with standard 40-dim filter-bank (FBK) features and the domain fine-tuned wav2vec2.0 model in Table 1 (Sys.1 vs. 4), the wav2vec2.0 model degrades the performance by 4.13\% absolute on the unseen subset.
\textbf{2)} Fusing FBK front-ends with the fine-tuned wav2vec2.0 features (Fig. 1(a)-(b)) as TDNN system inputs (Sys. 5), or optionally further incorporating cross-domain A2A inverted features (Fig. 1(c)-(d)) for multimodal TDNN system training (Sys. 6), consistently reduce the overall WER by 0.3\% and 1.03\% absolute over the wav2vec2.0 model (Sys.5, 6  vs. 4). \textbf{3)} Large overall WER reductions of 6.49-6.60\% absolute (11.98-12.21\% on unseen words) over the  wav2vec2.0 model are obtained by frame-level joint decoding between the LHUC-SAT TDNN system trained using FBK features alone, and those trained with additional wav2vec2.0 representations and inverted articulatory features (Sys.7-9 vs. 4)\footnote{System weights empirically set as 3:2, 3:2 and 9:1:5 for Sys. 7, 8 and 9}. \textbf{4)} Further performance improvements are obtained by cross system multi-pass decoding (Sys.10, 11), which uses the wav2vec2.0 model to rescore the joint decoded TDNN systems outputs. The lowest published WER of 22.56\% on 16 dysarthric speakers (Sys.10-11, 52.53\% on very low intelligibility, 39.09\% on unseen words) is obtained\footnote{System weights empirically set as 2:9 and 4:15 for Sys.10 and 11}. \textbf{5)} Similar experiments on Conformer systems also demonstrate the effectiveness of the wav2vec2.0 features via feature fusion (Sys.13 vs. 12), where smaller system combination improvements are obtained by multi-pass decoding (Sys.14 vs. 13), due to architecture similarity between the wav2vec2.0 and Conformer models leading to reduced cross system complementarity. Finally, the performance of our best systems (Sys.10, 11, Table 1) are contrasted against published state-of-the-art results on UASpeech in Table 2.

\begin{table}[htbp]
\vspace{-0.55cm}
\centering
\caption{WER (\%) of published and our best systems on UASpeech. ``15 spk'' stands for only 15 speakers in the test set, but the same 4 speakers are used in ``very low'' (VL) intelligibility group.}
\scalebox{0.72}{
\begin{tabular}{c|c|c}
\hline\hline
System                                                                   & VL      & All     \\
\hline\hline
Sheffield-2020 Fine-tuning CNN-TDNN speaker adaptation \cite{xiong2020source} (15 spk)                & 68.24 & 30.76  \\
CUHK-2021 NAS  DNN  +  Data  Aug.  +  LHUC-SAT  +  AV  fusion \cite{liu2021recent}        & 60.30 & 25.21  \\
CUHK-2022 DNN + Data Aug. + LHUC-SAT + AUV fusion \cite{hu2022exploiting}                    & 60.14 & 24.82  \\
CUHK-2022 DNN + Data Aug. + SBE Adapt + LHUC-SAT \cite{geng2022speaker}                             & 59.30 & 25.05  \\
CUHK-2022 TDNN + spectral basic GAN + LHUC-SAT \cite{jin2022personalized}                        & 59.18 & 27.85  \\
BUT-2022 Wav2vec2 + fMLLR + xvectors \cite{baskar2022speaker} (15 spk)                                 & 57.72 & 22.83  \\
\textbf{TDNN + wav2vec2 fea. + sys. combination (sys. 10 Table 1, ours) }& \textbf{52.53} & 22.83  \\
\textbf{TDNN + wav2vec2 fea. + sys. combination (sys. 11 Table 1, ours)} & 53.12 & \textbf{22.56}  \\
\hline\hline
\end{tabular}
}
\vspace{-0.3cm}
\end{table}
\noindent
\textbf{Results on DementiaBank:} \textbf{1)} Compared with the baseline TDNN system trained with 40-dim FBK features, the domain fine-tuned wav2vec2.0 model produces consistent and large WER reduction of 12.20\% absolute  (Sys.4 vs. 1). \textbf{2)} Similar to the UASpeech experiments, the incorporation of wav2vec2.0 features alone, or plus their corresponding cross-domain articulatory features, produces significant WER reduction up to 2.52\% (Sys.6 vs. 4). \textbf{3)} Frame-level joint decoding between the LHUC-SAT TDNN system using FBK features (Sys. 3), and the TDNN systems using additional wav2vec2.0 features (Sys. 5) and articulatory features (sys. 6) consistently outperforms the adapted wav2vec2.0 model (sys. 4) by a significant WER reduction of 2.91\% absolute (13.47\% relative, Sys. 9 vs. 4). \textbf{4)} After the cross system multi-pass decoding, the lowest published WER of 18.17\% on DementiaBank Pitt test set is obtained (Sys.11). Using its final recognition outputs on evaluation set to extract textual features, a state-of-the-art speech recognition based AD detection \cite{wang22l_interspeech} accuracy of 91.7\% is also obtained. A comparison of WERs and AD detection accuracy on DementiaBank between previously published ASR systems and our system is shown in Table 3.
\vspace{-0.6cm}
\begin{table}[htbp]
\centering
\caption{WER (\%) and AD detection accuracy (\%) between published systems on DementiaBank and our system.}
\scalebox{0.66}{
\begin{tabular}{c|c|c|c}
\hline\hline
\multirow{2}{*}{System}                                         & \multicolumn{2}{c|}{WER (\%)} & \multirow{2}{*}{Acc}  \\
\cline{2-3}
& All    & Eval PAR.      &                       \\
\hline\hline
TDNN + Data Aug. + domain adapt. + LHUC-SAT \cite{ye2021development}                           & 29.90 & 33.17        & 87.5                \\
TDNN + Data Aug. + SBE Adapt + LHUC-SAT \cite{geng2022speaker}                                & 29.16 & 31.89        & -                     \\
TDNN + spectral basic GAN + LHUC-SAT \cite{jin2022personalized}                           & 31.60 & 34.08        & -                     \\
Conformer + NAS + domain \& spk adapt  + sys. combination \cite{wang2022conformer}        & 24.2  & 25.5         & 87.5                \\
\textbf{TDNN + wav2vec2 fea. + sys. combination (sys. 11 table 1, ours)} & \textbf{18.17} & \textbf{16.84}        & \textbf{91.7}                \\
\hline\hline
\end{tabular}
}
\end{table}
\vspace{-0.6cm}
\section{CONCLUSION}
\vspace{-0.2cm}

This paper explores a series of approaches to integrate cross-domain adapted SSL pre-trained models into TDNN and Conformer ASR systems for dysarthric and elderly speech including input feature fusion, frame-level joint decoding and multi-pass decoding. In addition, domain adapted wav2vec2.0 representations are utilized in acoustic-to-articulatory (A2A) inversion to construct multi-modal dysarthric and elderly ASR systems. The TDNN and Conformer ASR systems integrated domain adapted wav2vec2.0 models consistently outperform the standalone domain fine-tuned wav2vec2.0 models by statistically significant WER reductions of 8.22\% and 3.43\% absolute (26.71\% and  15.88\% relative) on the benchmark UASpeech and DementiaBank Pitt test set respectively, while producing the lowest published WERs of 22.56\% and 18.17\% on the two tasks. Future research will study the personalization of pre-trained ASR models for dysarthric and elderly speakers.
\vspace{-0.3cm}

\clearpage
\bibliographystyle{IEEEtran}