Model of an oscillatory neural network with multilevel neurons for pattern recognition

The article presents a study on a new type of an oscillatory neural network (ONN) that uses multilevel neurons for pattern recognition on the basis of high order synchronization effect. The feature of this network architecture is a single oscillator (neuron) at the output with multilevel variation of its synchronization value with the main oscillator thus allowing classifications of an input pattern into a set of classes. ONN model is realized on thermally-coupled VO2-oscillators. ONN training was performed by using the trial-and-error method for the network parameters selection. It is shown that ONN is capable to classify 512 visual patterns (as a cell array 3x3, distributed by symmetry into 102 classes) into a set of classes with the maximum number of elements up to P=11. Classification capability of a network is studied depending on the interior noise level and synchronization effectiveness parameter. The obtained results allow designing multilevel output cascades of neural networks with high net data throughput and it may be applied in ONNs with various coupling mechanisms and oscillators topology.


Introduction
Hypotheses of functional importance of synchronization for information processed by the brain were put forward long ago [1], [2] and its experimental discovery [3], [4] encouraged creation of neural networks models with oscillatory dynamics and of neuromorphic algorithms of image processing based on synchronization [5]- [7].
Research of oscillatory neural networks is mainly based of the models of phase oscillators with slow variation in oscillation amplitude in a narrowband spectrum.Such models, primarily Kuramoto model [8], appeared to be very useful for studying systems of various topology with a large number (N>10 3 ) of oscillators, where such modes as global and cluster synchronization, synchronization by middle field [9] and mode of quasisynchronization [10] are feasible.
Alternatively, there is another class of ONN based on relaxation oscillators that generate multiple pulses of short duration and fixed amplitude.Such pulses (spikes) can code information at pulse-repetition frequency.There is an important distinguishing feature of a pulse type ONN from classic spiking neural networks (SNN): it is a self-oscillating mode of ONN neuro-oscillators, that is not indicative for either SNN, or real (biological) neurons performance because they are characterized by a forced response through generation of a single spike or their group when the neuron threshold level is achieved.However, ONNs are fascinating due to simplicity of hardware implementation and because the developed micro-and nanoelectronic self-oscillators can ensure their compactability and energy efficiency.Also pulse type ONNs with multifrequent periodic oscillation spectrum feature a specific mode of multiple frequency synchronization or, in other words, high order synchronization effect [11].We have recently demonstrated this effect experimentally by the example of thermallycoupled VO2-oscillators [12].In relaxation oscillators with vanadium dioxide-based film elements electric self-oscillations are activated by the effect of electric switching governed by metalinsulator transition (MIT) [13]- [16].Such features as processing speed of switching VO2 devices amounting to (10 ns) [17] and manufacturing technology that allows switching elements to be created with high level of nanoscalability make VO2-switch based oscillators perfect objects for research on neuro-oscillators ensembles to solve cognitive technology problems [18]- [21].
It should be noted that relaxation oscillators with high order synchronization effect can be realized by using electric coupling [18]- [21] and also by using not only VO2-switches but any other switching elements such as thyristors [22], tunneling diodes [23], resistive memory cells [24], spin torque oscillator [25], etc.
In the present paper, we study ONN of thermally-coupled VO2-oscillators and present a general concept of visual patterns recognition based on high order synchronization effect.We use the property of this effectmultiplicity of synchronous states variantsto extract object classes by using a single output oscillator and to realize uniquely high network throughput rate.

Oscillator circuit and method of ONN organization
The basic element in the studied ONN is an oscillator implemented on the circuit of a relaxation generator (Fig. 1a) based on a VO2-switch.We described in detail in [12] the process of fabrication and electric I-V characteristics of an electric VO2-switch.I-V characteristic is approximated well by a piecewise linear function that has two conduction states (low-resistance and high-resistance) and also a region of negative differential resistance (NDR).These switches may be used in the circuit of a relaxation generator with power supply Ip holding the operation point in the region NDR of I-V characteristic and with capacity C, parallel connected to the switch.In addition, there is a source of noise Un that models circuit interior or exterior noises, for instance, current noise of a switch [26].Oscillator output signal is current Isw(t) flowing through VO2-switch that is used to calculate synchronization level of two oscillators (see section 2.3).Besides, current signal directly conditions the effect of thermal coupling inside the network.We used thermal coupling to connect oscillators in a network [12], [27].Thermal coupling mechanism is based on mutual thermal effect of switches due to their Joule heating and dependence of gate threshold voltage Uth on temperature.The choice of this coupling type is determined by the simplicity of a computing model realization when oscillators are electrically separated from each other unlike in capacitive or resistive couplings [14], [15], [18], [21].An example of oscillators' interaction organization via thermal coupling of VO2 switches is shown in Fig. 1b.The model of thermal coupling is based on reduction of gate threshold voltage Uth by the value Δ at thermal effect of other switches.This effect, as we have already said, is determined by Joule heat generation on switches when pulse current is flowing through them at the moment of capacity discharge C.More detailed presentation of the mathematical model of thermally coupled relaxation oscillators is given in our previous paper [28].
It should be noted that we also use the concept of one-way thermal coupling of oscillators in the numerical model, it may be realized physically when a resistance heater is used as a connecting element in one of the circuits instead of a switch [27].
For numerical modeling we used the following parameters of I-V characteristic (Uth = 5 V, Uh = 1.5 V, Ubv = 0.82 V, Roff = 9.1 k, Ron = 615 ).In this circuit capacity is a constant parameter C=100 nF, its value significantly determines the frequency range of oscillator operation [29], and its natural frequency F0.In our case frequency range was 165 Hz≤F0≤1266 Hz at the range of feeding currents 550 µA≤Ip≤1061 µA.I sw (t)

ONN structure
The studied ONN consists of an input pattern that translates a pattern as 33 matrix on the layer of processing neurons including 10 oscillators, and of an output layer with only one oscillator (output neuron) №11 (see Fig. 2).The layer of processing neurons consists in its turn of a 33 oscillator matrix connected by similar couplings (i,j=j,i=mat, where i, j are the numbers of neighbor oscillators), connection between neighbor oscillators is made only by horizontal and vertical lines.Thus, the central oscillator №6 has four couplings in the matrix, the corner oscillators have two couplings, and the oscillators in the center of faces have three couplings (with the central oscillator and two corner ones).Besides, there is the basic oscillator №1 (first neuron) on the processing layer and the synchronization order of all other oscillators is measured against this oscillator.Oscillator №1 (Fig. 2) affects unidirectionaly all other oscillators in this matrix (1,i≠0=1, и i,1=0, где i=2…10) and they in their turn affect unidirectionaly oscillator №11 (output neuron) in the output layer (i,11≠0=11 и 11,i=0, where i=2…10).between the basic oscillator №1 and the oscillator in the output layer №11 serves as the control value.The values of feeding currents for these oscillators are fixed and may differ from currents in the matrix, therefore they are indicated as I1 and I11 (see.Fig. 3).

Method of synchronization order definition
To denote high-order synchronization value for N interacting oscillators in [30] we worked out the ratio of integers k1:k2:k3:..kN, where kN is harmonic order of N-th oscillator at the common frequency of the network synchronization Fs.In our case when we define synchronization order only between two oscillators (№1 и №2-11), it is possible to use the concept of subharmonics ratio as in [12], that is defined as a simple fraction SHRi,j=ki:kj ~ ki/kj , where i and j are the oscillators numbers between which the synchronization order is measured.It is obvious that to define the synchronization order for the first oscillator SHR1,1 in relation to itself is not reasonable because SHR1,1=1.
In our study the main technical problem was the problem of defining the synchronization order between the first (basic) oscillator №1 and the oscillator of the output layer №11 characterized by the value SHR1,11: It can be seen that the same value of SHR1,11 may be expressed in several ways: as a ratio, a simple fraction or a real number, for example, SHR1,11=10:3=10/3=3.33.Further, we will use this property to present the results more vividly.
The method of SHR definition for two oscillators is described in detail in our papers [12].
The presence or absence of synchronization was determined by a threshold method by using the value of synchronization effectiveness η so that at η ≥ ηth the two-oscillator system was considered as being synchronized.We varied the threshold value ηth during the numerical experiment.
Nevertheless, the commonly used value is ηth=90% in most cases (where it is not specified).
Current oscillograms Isw(t) of oscillators №1-11 were calculated simultaneously and contained ~250 000 points with time interval t=10µs.After that the oscillograms were automatically processed, synchronization order SHR1,i was determined, where i =2…11.

Pattern classifier implementation and problems definition
A "black-and white" pattern is used as an input test pattern presented as a 33 matrix (without gradation of gray color, 3 by 3 pixels) where each cell may take the value 0 (white color) or 1 (blue color).The total number of patterns (figures) in the input layer matrix is 2 9 =512.Presuming that the pattern processing layer together with the output layer has certain symmetry, a set of 512 figures may be divided into 102 classes (see Fig. 4).and also at rotation by 90 • (see the example of transformation for class №4 figures in Fig. 5).
where the class number Nc changes within all possible range Nc=1...102.
II.There is a set of classes S with the elements number P (P<102), S ={Class № S1, Class № S2… Class № SP}, for which there is a corresponding set of high order synchronizations SHR ={SHR (1)  1,11, SHR (2) 1,11… SHR (P) 1,11} with the same number of elements.The set SHR does not have the same elements, i.e. each class of figures from set S corresponds to a unique synchronization order SHR1,11.By analogy with (2) the problem may be expressed as: Actually, problems I-III have increasing degree of complexity and are subcases of problem II with different parameter P. Nevertheless, it should be noted that problem I is the simplest one which may be solved by a common neural network with one bistable output neuron (bistability means presence or absence of synchronization), although two variants of answers in neural networks are more often organized by using two output neurons [31].
Whereas the output neuron in problems II and III should have multilevel properties, thus being different from a common bistable neuron.All three problems may be set for the circuit shown in Fig. 2. We should note the most striking difference of this neural network from variants presented in the literature.The ONN circuit has only one output neuron, nevertheless, the effect of multilevel high order synchronization used here and characterized by value SHR1,11 allows input pattern classification into P classes within set S. This increases net data throughput of a single neuron and enables us to create multilevel output cascades of neural networks with high functionality.

Technique of ONN training
To solve problems I-III it is necessary to be able to train the network.As the ONN with the high order synchronization effect has not been studied before there are no specially developed methods of this network training.One of the obvious ways here is to use the trial-and-error method for the network parameters selection: currents (ION, IOFF, I1, I11), couplings strength (1, 11, mat), noise amplitude Un and synchronization effectiveness threshold ηth.
To realize direct searching of all parameters we should first determine the ranges and steps (number of grades) of their variations.In our numerical experiment current range was determined by the presence of a single oscillator generation effect.S, for the oscillator circuit described in section 2.1 the generation was within the feeding current range of 550 µA ≤ Ip ≤ 1061 µA.
Therefore, the currents (ION, IOFF, I1, I11) also varied in this range.We determined the variation steps as Ip = 1 µA.So there were only 512 current grades.
For coupling strength  variation range was determined by the threshold values of I-V characteristic of a switch.So the condition when integral thermal action on the neuron  should not reduce the threshold voltage of I-V characteristic of a switch Uth below the holding voltage must be met: On the basis of formula ( 4), values of threshold voltages (Uth=5 V, Uh=1.5 V) and couplings configurations (see Fig. 2) the limits of coupling strength variation are subject to the following conditions: 0 ≤ 1 < (3.5-4mat), 0 ≤ 11 < (3.5/9), 0 ≤ mat < ((3.5/4)-1),where units for values  are volts.The spacing was chosen as  = 0.1 V.
Range 20 µV ≤ Un ≤ 900 µV with the number of grades equal 12 was chosen for the noise amplitude.
For the synchronization effectiveness threshold ηth the variation range was 15 % ≤ ηth < 100%, with the number of grades equal 25 and minimal spacing ηth=1%.Here we should specify that ηth does not belong to the network parameters but rather to the parameters of the algorithm of synchronization order SHR1,11 calculation.Nevertheless, the value ηth strongly affects the result of synchronization definition and the results of pattern recognition as a whole, because it is a conditionally chosen characteristic.Identifying its optimal value for the recognition problems solution is an important step in a network tuning and training.
There are a lot of network parameters and each parameter variant should be calculated in 102 classes together with oscillograms and synchronizations values SHR1,11 calculations at each stage; it is obvious that full direct searching of all parameters variants will take a lot of time and computational resources.
Therefore, we developed software that enabled direct search of only one specific parameter, in particular, one of the currents (ION, IOFF, I1, I11), in this case all other parameters were fixed and selected on the basis of our experience of working with the model.
For the pattern recognition experiments we used a workstation (Intel Xeon quad core processor, 4x2 GHz, 8GB RAM) running 64 bit Windows Server 2008.CPU time cost by a single run on a one core for the direct search procedure of all 512 current grades took ~2.4 hours.
In conclusion, for this section we should note that direct search is not an effective training method and development of more efficient algorithms for this type of networks is a separate global problem for future research.

Results
As we have mentioned before, to solve problems I-III we should be able to train the network, and one of the obvious variants is direct search of parameters values (ION, IOFF, I1, I11), (1, 11, mat ), Un and ηth.
Solution of problem I.
For example, at ION=797µA (see Fig. 7a) we have a neural network that recognizes Class №102 with the corresponding synchronization order SHR1,11=10:9, in this case there is no synchronization for all other classes.It should be noted that incorrect solution of training, which do not correspond to the problem I conditions occur rather often, for example, when two or more classes correspond to the same value of SHR1,11 (see Fig. 7b).Such neural network parameters result in ambiguous recognition of figures and their classes.The same holds true for frequent cases of wrong training when there is no oscillatory synchronization for any input class (figure).

Solution of problem II.
If to analyze all results of direct search of current values ION, which arose in problem I solution with the same network parameters, then we will see that there will be solution variants for problem II as well.Fig. 7 demonstrates the variant of set S and SHR for P=2, and Fig. 7d demonstrates the variant of set S and SHR for P=5.
The condition of problem II requires that a certain class from set S should correspond to the unique synchronization order from set SHR, at the same time there is no output oscillator synchronization for other classes.
Table 1 presents a list of found current values ION, when the conditions for problem II are met for all emerging values P. It can be seen that value P reaches its maximum P=7.Naturally, there might be more than one solution with P≥1.For example, solution with P=2 is obtained both at currents ION= 864µA and ION= 813µA, see Fig. 7с, e, and also at some other current values.
Therefore, we have introduced quantity NP that is equal to the total number of solutions for the given value P during direct search of ION.Values of NP are also given in Table 1.It can be seen from the table that the most frequent solution appears for P=1 (Np =19) and the number of solutions reduces with the increase of P. We constructed a 3D graph shown in Fig. 9 to find the dependence of maximum possible P and Np on the noise amplitude value in network Un.
It can be seen that when the noise increases above Un = 400 µV the number of solutions NP and the value P sharply fall down.In this case the maximum value of P still does not exceed P = 7 and occurs in a wide noise range up to Un ≤ 200 µV.The basic quantity of solution variants is distributed in the range 1 ≤ P ≤ 4. It is interesting that the maximum value for Np corresponds to P = 2 at noise level Un = 100 µV.In this case the integral value ∑NP for all values P also has We constructed a 3D graph shown in Fig. 10 to find the dependence of maximum possible P and Np on the value of threshold synchronization effectiveness.
A general tendency for reduction of Np at reduction of ηth below 40% can be seen, in this case it is interesting to note that growth of ηth up to ηth = 99% in average does not change the values of P and Np.This seems to be caused by the fact that synchronization occurring during problems I and II solution has a high value of effectiveness η>99%.In this case reduction of just adds "extra" synchronous states which do not meet the problems conditions.
The presence of maximum possible P reaching P = 11 at ηth = 35% is also an interesting fact.In this case we can also observe some resonance for values of P. Its presence might be caused by reduction of maximum P both at growth of ηth, conditioned by solutions discarding and at reduction of ηth conditioned by emerging of "extra" synchronizations.

Discussion and Conclusion
This paper presents a new model of ONN with high order harmonics synchronization that classifies (recognizes) two-dimensional visual patterns by unique synchronous states, i.e. in accordance with their definite symmetry class.The most outstanding feature of this neural network in comparison with variants published is that neurons possess not only bistable properties (presence or absence of synchronization with the basic neuron) but have multilevel synchronization.The mode considered here has only one output neuron, nevertheless, variation of the value of high order synchronization SHR1,11 allows its usage to classify the input pattern into P classes with the set S.
Actually parameter SHR1,11 belongs to the properties of an output neuron while the first (basic) neuron may be considered as a master generator in regard to which we calculate synchronization of other network neurons.
Multilevel high order synchronization increases net data throughput of a single neuron and enables creating multilevel output cascades of neural networks with high functionality.
It is shown that the problem of a network training may be solved by direct search of the values of oscillators' feeding currents for the variant of incomplete classification when the number of patterns is less than the maximum number of classes (P<102).We have found solution for problem II with the maximum value P=11.
Solution of the problem of complete classification (problem III) when ideal training is realized and each class out of 102 possible ones corresponds to its own unique value of synchronization order SHR1,11 has not been found yet.The ways of its solution are related to more expanded direct search of possible ONN parameters.For example, we may vary the parameters of I-V characteristic, the configuration of the network, distribution of coupling strength or the method of pattern translation into a network.Besides, an additional argument for searching new and more effective methods of pattern translation into ONN is that SHR1,11 varies only in a very narrow range near the average value (see Fig. 8).Addition of alternative methods of pattern translation into a network that vary, for example, current I11 may significantly expand the range of variations and thus the maximal attainable value P.
Besides, the technique of SHR determination may affect the result and future research aimed at its improvement also may lead to positive shift in this field.
In conclusion, we should note again that development of ONN models with high order synchronization effect hold the significant potential for increasing informational capacity and efficiency of artificial intelligence networks and development of their training techniques is an important direction for further research.

Fig. 1 .
Fig.1.Oscillator circuit (a) and example of oscillators' interaction via thermal coupling of VO2 switches (b).Isw(t)current signal in a VO2 switch which is an output signal; Unsource of noise.

Fig. 2 .Fig. 3 .
Fig. 2. ONN organization circuit for pattern recognition as a matrix of 3×3 elements.Digits indicate the sequence numbers of oscillators.

Fig. 7 .
Fig. 7. Training results at various values of current ION: examples of correct (а) and incorrect (b) solutions of problem I; examples of correct solutions of problem II (с)-(e).

Fig. 9 .
Fig. 9. Dependence of solution number NP for various values of P on noise level Un.

Fig. 10 .
Fig. 10.Dependence of solution number for various values of P on the threshold synchronization effectiveness ηth.Solution of problem III Solution of problem III when ideal training is realized and each class out of 102 possible ones corresponds to its own unique value of synchronization order SHR1,11 has not been found yet (see section 4).

Table 1 .
Variants of current values ION, at which problems I -II are solved for various P and the corresponding number of their solutions NP.Figure8demonstrates statistics of values SHR1,11 variation for various values P, here the synchronization order is given as a real number (see section 2.3).In other words, this figure shows all possible values SHR1,11 for all values ION from Table1.It can be seen that SHR1,11 has an average value ~ 1.079, in this case the range of variation is small with maximum deviation from the average value ~ 0.103.We suppose that weak dependence of SHR1,11 on currents ION and IOFF is conditioned by the permanence of currents I1=750µA, and I11=750µA that determine natural frequencies of outermost oscillators, and in this particular case the currents are equal.Equality of natural frequencies of oscillators №1 and №11 sets the variation area of SHR1,11 near SHR1,11=1:1=1.Thus, the matrix of processing layer oscillators just deviates the value SHR1,11 from some average one.Fig. 8. Values SHR1,11 at various P for all variants ION from Table 1.Dashed line indicates the level of average value.