ViR: the Vision Reservoir

Based on CRJ, we argue that for an image classification task it is necessary to consider: 1) a simple fixed non-random reservoir topology–nearly full connection structure with self-loop connections; 2) the weight value $\textit{$r_{i},r_{j},r_{k}$}\in(0,1]$ for reservoir connections; 3) the signs of the weights in 2) are deterministically generated by a rule with a certain degree of randomness; 4) input neurons are randomly connected to the reservoir.

Fig. 3 depicts the nearly fully connected topology. We first generate an $N\times N$ reservoir connection matrix W matching the topology, in which all the weights have the same original value $r_{k}$ and then we set:

$$\textbf{{W}}_{1,\textit{N}}=\textbf{{W}}_{\textit{q}+1,\textit{q}}=\textit{}{r_{i}},\qquad\textit{q}=1,...,\textit{N}-1$$

(1)

$r_{j}$ is the jump weight. Consider the jump size $1<\ell<\textit{N}$. The first jump begins from neuron 1 to $1+\ell$, then from $1+\ell$ to $2\ell+1$, until $n\ell+1>N$, where n is the number of jump. These bi-directional jump weights are reset to:

		$\displaystyle\textbf{{W}}_{1,\ell+1}=\textbf{{W}}_{\ell+1,2\ell+1}=...$		(2)
		$\displaystyle=\textbf{{W}}_{\ell+1,1}=\textbf{{W}}_{2\ell+1,\ell+1}=...=\textit{$r_{j}$}$

Finally, we randomly set an element, with the symmetrical element, to 0, indicating that the two neurons are disconnected with each other. Hence, W mainly consists of $r_{k}$, and the amount of $r_{j},r_{c}$ is relatively small, with two 0 elements.

For the randomness, there exists a random matrix with elements e randomly generated in (0,1] with the same size as W. If e $<$ 0.5 (the threshold we set), then the sign of the weight in W, corresponding to the same position as e, will be $-$, otherwise $+$. The input weight matrix V, compared to W, is a sparse matrix of the given shape with uniformly distributed values.

For the stability of the reservoir, W is typically scaled as $\textbf{{W}}\leftarrow\alpha\textbf{{W}}/|\lambda_{max}|$, where $|\lambda_{max}|$ is the spectral radius of W and $0<\alpha<1$ is a scaling parameter [19].

For the training process, assuming that $\textit{{u}(t)}=(u_{1}(t),...,u_{K}(t))$ means the activation of $K$ input units at time step $t$ from flattened patches, and $\textit{{x}(t)}=(x_{1}(t),...,x_{M}(t))$ represents the $M$ internal units in the reservoir with the output $\textit{{y}(t)}=(y_{1}(t),...,y_{Q}(t))$ containing $Q$ output units. The $M\times K$ matrix V stores input weights. The $M\times(K+M+Q)$ matrix $\textbf{{W}}_{out}$ shows the connections from a reservoir to output units. Additionally, all the values mentioned above are real-valued. Apart from $\textbf{W}_{out}$, all the matrices are not trainable and their entries are initialized with a suitable distribution and then fixed. Fig. 2 shows the described structure.

The internal units are updated according to:

$$\textbf{{x}}(\textit{t}+1)=\textit{f }(\textbf{{Vu}}(\textit{t}+1)+\textbf{{Wx}}(\textit{t})+\textbf{{b}}(\textit{t}+1))$$

(3)

where f is the activation function (typically sigmoidal); b(t+1) is an optional uniform i.i.d. noise. The readout is computed by:

	$\displaystyle\textbf{{y}}(\textit{t}+1)=\textbf{{W}}_{out}$	$\displaystyle[\textbf{{u}}(\textit{t}+1);\textbf{{x}}(\textit{t}+1);\textbf{{y}}(\textit{t});$		(4)
		$\displaystyle\textbf{{u}}(\textit{t}+1)^{2};\textbf{{x}}(\textit{t}+1)^{2};\textbf{{y}}(\textit{t})^{2}]$

where ’;’ is a concatenation operation.

As shown in Eq. 4, the essence of network is to approximate output weights $\textit{{W}}_{out}$ through training samples, to obtain the predictive ability of certain tasks. Hence, while training, it is essential to collect and store the state of the reservoir x(t) and the corresponding outputs y(t) in a matrix X and Y after warm-up [19]. The classic calculation method is the Least Square Method (LSM) [33]:

$$\textit{{W}}_{out}=(\textit{{X}}^{T}\textit{{X}})^{-1}\textit{{X}}^{T}\textit{{Y}}$$

(5)

where $(\ast)^{T}$ and $(\ast)^{-1}$ mean the transpose and the inverse of matrix operator, respectively. Usually, in case X is ill-posed and irreversible, the ridge regression method [45] is utilized (Eq. 6) instead of LSM. For other optimizing algorithms refer to [35, 39].

$$\textit{{W}}_{out}=(\textit{k}\textit{{I}}+\textit{{X}}^{T}\textit{{X}})^{-1}\textit{{X}}^{T}\textit{{Y}}$$

(6)

where k is a relatively small positive number called regular coefficient, and I is an identity matrix. In our work, we use a trainable linear layer (LL) to approximate $\textit{{W}}_{out}$ and Eq. 4 can be rewritten as:

	$\displaystyle\textbf{{y}}(\textit{t}+1)=\textbf{{LL}}$	$\displaystyle[\textbf{{u}}(\textit{t}+1);\textbf{{x}}(\textit{t}+1);\textbf{{y}}(\textit{t});$		(7)
		$\displaystyle\textbf{{u}}(\textit{t}+1)^{2};\textbf{{x}}(\textit{t}+1)^{2};\textbf{{y}}(\textit{t})^{2}]$

Fig. 2 depicts the proposed model for image classification, and the key component is the ViR core which is composed of a reservoir with the internal topology mentioned above and a residual block.

The typical reservoir receives a time-series sequence as input. Similar to the treatment in ViT, some pretreatments are used to convert raw images to patches as the time series input. We reshape the original image $\textbf{x}\in\mathbb{R}^{\textit{H}\times\textit{W}\times\textit{C}}$ into a series of patches $\textbf{x}_{\textit{p}}\in\mathbb{R}^{\textit{N}\times(\textit{P}^{2}\cdot\textit{C})}$, where $(\textit{H, W})$ and $(\textit{P, P})$ respectively represent the resolution of x and $\textbf{x}_{\textit{p}}$, and $\textit{N}=\textit{H}\textit{W}/\textit{P}^{2}$ is the number of patches, also serving as the input time steps. With a trainable linear projection E (Eq. 8), we flatten the patches and map to D dimensions to make the inputs suitable for the reservoir. We refer to the output of this projection as the patch coding. Therefore, the input at time step t+n can be written as:

$$\textbf{{u}}(\textit{t+n})=\textbf{x}^{n}_{\textit{p}}\textbf{{E}},\quad\textbf{{E}}\in\mathbb{R}^{(\textit{P}^{2}\cdot\textit{C})\times\textit{D}},\textit{n}=0,1,...,\textit{N}-1$$

(8)

Then the ViR core receives the data. It consists of the nearly fully connected reservoir (Fig. 3) and a residual block including a layer-norm (LN) operation and a feed-forward (FF) layer. The non-linear function in FF is a Gaussian Error Linear Unit (GELU) [16]. The processing of the reservoir has been shown in Eq. 3 and Eq. 4, and LN is applied to $\textbf{y}_{\textit{r}}$ (the output of the reservoir) with the feed-forward layer with a residual connection between them, shown as:

$$\textbf{{y}}_{\textit{c}}=\textbf{{FF}}(\textbf{{LN}}(\textbf{{y}}_{\textit{r}}))+\textbf{y}_{\textit{r}}$$

(9)

A classification operation is attached to $\textbf{{y}}_{\textit{c}}$ (the output of the ViR core ), which is implemented by a MLP, and we can obtain the image representation and the final classification results y by Eq. 10:

$$\textbf{{y}}=\textbf{{LN}}(\textbf{{y}}_{\textit{c}})$$

(10)

We further stack reservoirs to get the deep ViR to enhance the network performance.

The first one is the series reservoir consisting of L reservoirs depicted in Fig. 4. Similarly, the input matrices $\textbf{{{V}}}^{(l)}$ and reservoir connection matrices $\textbf{{W}}^{(l)}$ are constant and generated as mentioned above, l $\in$ $\{$1,2,…,L$\}$, with similar training process described in subsection 3.1.

$\textbf{{U}}^{(l)}$ is the output matrix. Units in reservoir l updating and the output of reservoir l are given by:

$$\textit{$\textbf{x}^{(l)}$}(\textit{t}+1)=\textit{f }(\textit{$\textbf{V}^{(l)}$}\textit{$\textbf{y}^{(l-1)}$}(\textit{t}+1)+\textit{$\textbf{W}^{(l)}\textbf{x}^{(l)}$}(\textit{t}))$$

(11)

$$\textit{$\textbf{y}^{(l)}$}(\textit{t}+1)=\textit{$\textbf{U}^{(l)}$}[\textit{$\textbf{y}^{(l)}$}(\textit{t});\textit{$\textbf{x}^{(l)}$}(\textit{t}+1);\textit{$\textbf{y}^{(l)}$}(\textit{t})^{2};\textit{$\textbf{x}^{(l)}$}(\textit{t}+1)^{2}]$$

(12)

where l $\in$ $\{$1,2,…,L$\}$ and $\textbf{{y}}^{(0)}(\textit{t})=\textbf{{u}}(\textit{t})$, $\textbf{{y}}^{(L)}(\textit{t})=\textbf{{y}}(\textit{t})$.

The other one is a parallel reservoir, depicted in Fig. 4.

L reservoirs simultaneously receive the same input sequence u(t) and all L reservoirs can be trained simultaneously referring to the above process. The units in L reservoirs and the output units are updated as follows:

$$\textit{$\textbf{x}^{(l)}$}(\textit{t}+1)=\textit{f }(\textit{$\textbf{V}^{(l)}$}\textbf{{u}}(\textit{t}+1)+\textit{$\textbf{W}^{(l)}\textbf{x}^{(l)}$}(\textit{t}))$$

(13)

$$\textit{$\textbf{y}^{(l)}$}(\textit{t}+1)=\textit{$\textbf{U}^{(l)}$}[\textit{$\textbf{y}^{(l)}$}(\textit{t});\textit{$\textbf{x}^{(l)}$}(\textit{t}+1);\textit{$\textbf{y}^{(l)}$}(\textit{t})^{2};\textit{$\textbf{x}^{(l)}$}(\textit{t}+1)^{2}]$$

(14)

The final output of a parallel reservoir is the arithmetic mean of L reservoir outputs, given by:

$$\textit{{y}}(\textit{t})=\sum_{l=1}^{L}\textit{$\textbf{y}^{(l)}$}(\textit{t})/\textit{L}$$

(15)

Usually, we evaluate the quality of a reservoir through three indicators: the Small-World (SW) characteristics, the Largest Lyapunov Exponent (LLE), and the Memory Capacity (MC) [37]. The reservoir with SW characteristics has a stronger information processing ability and faster information dissemination speed than the common reservoir [21, 1]. LLE of a reservoir represents the dynamics and stability, and usually the value of LLE should be approximately equal to 0 for great performance. MC is the property of a reservoir to keep the previous input information, related to the data processing capability.

Small-World Characteristics: RC is a brain-inspired neural network, which mimics the connections of cerebral neurons, cortical neural connectivity has been shown to exhibit a small-world (SW) network topology [58], which has a shorter average path length and larger clustering coefficient than regular networks. In a reservoir, a small clustering degree means that dynamic information flow through the reservoir nodes is not ‘too cluttered’. Also, a small average path length can allow for the representation of a variety of dynamic time scales.

Due to the bi-directional connections, we can view the interconnection topology as an undirected graph G=(J, E), J and E representing neurons and the connections in the reservoir, respectively. This is similar to the kernel connections in the ViT. Define $\textbf{{l}}_{\textbf{{G}}}$ to be the average path length between vertex pairs in J, which is calculated by [38]:

$$\textbf{{l}}_{\textbf{{G}}}=\frac{2}{M(M+1)}\sum_{\begin{subarray}{c}i\geq j\end{subarray}}d_{ij}$$

(16)

where $d_{ij}$ is the geodesic distance from vertex i to vertex j, and M is the number of nodes (equal to the number of neurons in the reservoir). If $i=j$, the distance is 0 [38]. Moreover, if there are no self-loop connections, Eq. 16 should be rewritten as:

$$\textbf{{l}}_{\textbf{{G}}}=\frac{2}{M(M-1)}\sum_{\begin{subarray}{c}i>j\end{subarray}}d_{ij}$$

(17)

Assuming that the neighborhood $|\textit{M}_{i}|$ of the node $\textit{v}_{i}$ is the nodes next to $\textit{v}_{i}$, $\textit{M}_{i}$ is the number of nodes in $|\textit{M}_{i}|$. The local clustering coefficient $\textbf{{C}}_{\textbf{{i}}}$ of $\textit{v}_{i}$ is given by a ratio, i.e., the connecting edges between nodes in the neighborhood divided by the number of possible connecting edges between them, shown in Eq. 19 [58]:

$$\textbf{{M}}_{\textbf{{i}}}={\textit{v}_{j}:\textit{e}_{ij}\in\textbf{{E}}\vee\textit{e}}_{ji}\in\textbf{{E}}$$

(18)

$$\textbf{{C}}_{\textbf{{i}}}=\frac{2|{\textit{e}_{jk}:\textit{v}_{j},\textit{v}_{k}\in\textbf{{M}}_{i},\textit{e}}_{jk}\in\textbf{{E}}|}{M(M-1)}$$

(19)

The overall clustering level $\bar{\textit{C}}$ is $\bar{\textbf{{C}}}=\frac{1}{M}\sum_{\begin{subarray}{c}i=1\end{subarray}}^{\begin{subarray}{c}M\end{subarray}}C_{i}$.

The degree of SW, named small-worldness, is given as:

$$\delta=\frac{C}{C_{0}}/\frac{l_{G}}{l_{0}}$$

(20)

where $\textit{C}_{0}$ and $\textit{l}_{0}$ represent the average path length and clustering level of a regular network with similar size. If $\delta>1$, the network is a SW network [17].

Lyapunov Exponent: The Largest Lyapunov Exponent (LLE) [59] represents the stability of reservoir dynamics, expounding the sensitivity of initial conditions of a system to small perturbations, defined as [51]:

$$\textbf{$\lambda$}=\lim_{k\to\infty}\frac{1}{k}ln(\frac{\gamma_{k}}{\gamma_{0}})$$

(21)

where $\gamma_{0}$ is the initial distance between the perturbed and the unperturbed trajectory, $\gamma_{k}$ is the distance at time k. For ordered dynamic systems, $\lambda<0$ and for chaotic systems $\lambda>0$. At $\lambda\approx 0$, a phase transition occurs (called the critical point, or edge of chaos) with the best computational capability [3, 4].

Memory Capacity: Memory Capacity (MC) [20] is another metric to measure the learning performance of our model. With receiving a random input u(t) at a time t, the reservoir is trained to generate the desired output $\textbf{y}_{d}(\textit{t})=u(t-\tau)$. The output y is learned from u(t) which was $\tau$ steps earlier. Then the MC is given as:

$$\textbf{MC}_{\tau}=\frac{cov^{2}(\textbf{u}_{\tau},\textbf{y})}{v^{2}(\textbf{u}_{\tau})v^{2}(\textbf{y})},\qquad\textbf{MC}=\sum_{\tau=1}^{T}\textbf{MC}_{\tau}$$

(22)

where $cov^{2}(\textbf{u}_{\tau},\textbf{y})$ indicates the covariance between the true value $\textbf{u}_{\tau}$ and the predicted value y. $v^{2}(\ast)$ means the variance of a series $\ast$. T is the maximal time delay we set.

We have experimentally demonstrated the SW characteristics, largest Lyapunov exponents, and memory capacities of our model, in Section 4.1.

Table 2 shows the normalized average path length $\frac{l_{G}}{l_{0}}$, clustering coefficient $\frac{C}{C_{0}}$, and small-worldness $\delta$, $\frac{C}{C_{0}}/\frac{l_{G}}{l_{0}}$. The baseline is the classic regular network and random network detailed in [23].

Our model has a lower clustering level and a smaller average path than the regular one. It still exhibits the SW property, referring to the fact that our topology has great potential for information calculation and processing. Similar to the brain [42], a small average path indicates low wiring cost for physical implementations. By applying a standard algorithm given in [46], inputs are sampled from Gaussian distribution [54] for the memory capacity and the largest Lyapunov exponents.

Fig. 5(a) shows we get the peak Memory Capacity (MC) when the spectral radius $\rho$(W) is around 1. Also, Fig. 5(b) shows that the largest Lyapunov Exponent (LLE, defined in subsection 3.4) is close to 0 when $\rho$(W) is around 1$\sim$1.5. Meanwhile, the influences of the input scaling (IS) to MC and LLE are also shown in Fig. 5. And the results indicate that we should constrain the $\rho$(W) and the IS to about 1 for image tasks.

Without any pre-training, we compare our models ViR-1, ViR-3, ViR-6, and ViR-12 with ViT-1, ViT-3, ViT-6, and ViT-12 by performing image classification tasks on MNIST, CIFAR10, and CIFAR100. Table 3 shows the classification accuracy and the number of parameters. It can be seen that ViR-1 has more competitive results than ViT-12 and therefore ViR models have a promising advantage in saving parameters. According to the aforementioned analysis, the reservoir with a nearly fully connected topology has great information processing capability and rich dynamics. Hence, our model is suitable for handling image tasks.

Accuracy: All the ViR models except Series-ViR perform better than ViT models with the same depths. Depths are the layers of the ViR or the encoders of the ViT. Interestingly, the shallow model ViR-1 usually rivals ViT-1, ViT-3, and ViT-12. This indicates that shallow ViR models have a competitive performance compared with the deep ViT models, but with lower time costs and fewer parameters. The great MC of the reservoir donates the high accuracy to some degree.

Convergence Speed: As depicted in Fig. 1, the initial test accuracy of the ViR is already close to the best results of the ViT. Meanwhile, the ViT-12 suffers over-fitting. The convergence speed of most ViR models is much faster than ViT. This reflects the SW characteristics and rich dynamics from LLE.

The Number of Parameters: Due to the fact that training a reservoir is limited to the output layer, $\textbf{{W}}_{out}$, much fewer parameters are trained compared to other models. This is one of the attractive characteristics of RC. In this aspect, the ViR has the promising advantage compared with ViT.

Table 3 shows the trainable parameters of the ViR and the ViT on different datasets with different layers. The number of parameters of the ViR is about 15% of the ViT with the same depth. But the comparison between shallow ViR models and deep ViT models shows that the number of parameters of the ViR could be about 5% of the ViT while guaranteeing an acceptable accuracy.

Table 3 also indicates that, with the same depth, the series reservoir usually performs worse than the parallel one. One possible explanation is that: the series one always randomizes the input to the next reservoir, having no substantial help for tasks. However, the parallel one randomizes the inputs only once and then gets the arithmetic mean results. As a result, the enhancement of accuracy is predictable. With the increase of reservoir layers, the gains from deep structures can’t offset the influence of time and computational costs compared to shallow ones. This reveals that it is not necessary to stack too many reservoirs.

Considering the potential crash risk from memory overloading, it is necessary to evaluate and reduce its usage in a single task while ensuring acceptable results. The memory footprint of the ViR, compared with ViT models, on MNIST and CIFAR100 datasets, is shown in Fig. 6.

For small patch size, we observed that the memory footprint of the ViT models is about 2.5 to 5 times to ours, with the same depth. Although the large patch size of the ViT could reduce the memory footprint, our models (with small patch size) still have a slight superiority in such cases. It is worth noticing that, for both the series ViR and the parallel ViR the memory footprints within the optimal number of layers are nearly equal.

We evaluate the robustness of the ViR from two aspects: i.e. corruptions of input images and perturbations of system hyperparameters.

Corruptions of input images. The input images are selected from the CIFAR10-C dataset, which consists of imposed corruptions from noise, blur, weather and digital categories on the original CIFAR10 dataset [15].

To achieve the comprehensive evaluation of the robustness to a given type of corruption, we score the classification performance across five corruption severity levels, denoted by s ($1\leq s\leq 5$), and aggregate these scores according to:

$$\textbf{CE}_{c}=\sum_{s=1}^{5}\textbf{E}_{s,c}-\textbf{E}_{clean},\textbf{CE}=\frac{\textbf{CE}_{c}}{n}$$

(23)

where $\textbf{E}_{clean}$ is the top-1 error rate on clean CIFAR10 dataset; $\textbf{E}_{s,c}$ is the top-1 error rate on CIFAR10-C, c is the type of corruptions. CE is arithmetic mean of all types of $\textbf{CE}_{c}$, and n represents the number of corruption types. It should be pointed out that, the smaller value of $\textbf{CE}_{c}$ or CE represents the stronger robustness.

Table 4 compares the robustness performance by means of $\textbf{CE}_{c}$ and CE. The results show that the number of layers within parallel-ViR has little influence on robustness. With the increasing encoding numbers of the ViT, it has increasingly robustness, which is contrary to in the cases of the series ViR model.

Perturbations of system hyperparameters. From Fig. 1, we can see that the number of layers has a smaller impact on our models, compared with the ViT, perhaps because the spectral radius of W in the ViR is constrained to be around 1 when training [37]. Moreover, the influence of IS is also reduced because of the constrained spectral radius, shown in Fig. 7. It has the same conclusion for other hyperparameters of the ViR, such as the influence of the number of neurons N, $r_{i},r_{j},r_{k}$ and the jump size. SW characteristics support our experimental results, which make RC have strong robustness. However, the hyperparameters in ViT, such as the number of heads, patch size, and MLP size etc., will significantly impact the classification accuracy. This indicates that our model has a certain degree of robustness to the perturbations of system hyperparameters.