NdNiO3 array fabrication
NdNiO3 thin film (50 nm) was deposited on a 2 inch LaAlO3 (LAO) (100) single-crystal wafer (1 mm thickness) at room temperature using radio-frequency (RF) sputtering. The deposition process utilized co-sputtering of Ni and Nd targets with 75 W DC and 125 W RF power, respectively. Deposition was performed at 5 mTorr pressure maintained by argon/oxygen gas mixture with flow rates of 40 sccm and 10 sccm, respectively. After deposition, the wafer underwent post-annealing at 550 °C for 24 h in 25 sccm O2 gas flow to enhance crystallinity. A bilayer lift-off process (LOR5B and AZ1512) was used to define Pd electrode patterns on the NNO/LAO substrate. Pd electrodes (50 nm) were subsequently deposited via sputtering. Ti/Au (10 nm/100 nm) electrodes were then patterned and deposited using the same bilayer lift-off process. During the Au deposition step, a smaller-area Au pad was also deposited on top of each Pd electrode to improve electrical contact during device measurements. After Au and Pd metal deposition, the samples were placed in a tube furnace and annealed at 115 °C for 20 min in a hydrogen/argon (5%/95%) atmosphere at a flow rate of 35 sccm. This process led to hydrogenation of the NdNiO3 (H-NNO) thin films beneath the Pd electrodes, forming the proton-modulated nickelate junctions used in the device array.
KPFM
We have used an atomic force microscopy (AFM) instrument (NT-MDT Solver Next AFM) for the AFM and KPFM experiments. The instrument consists of an HA_NC tip made of monocrystalline silicon. The curvature of the tip was less than 10 nm, and it has a spring constant of 2.8 N m−1. AC modulation and frequency values were 0.5 V and 131 KHz. During a two-step scan process, topography of the samples was recorded using AFM in tapping mode. In the next step, contact potential difference (VCPD) between the tip (which was retracted by Δz = 10 nm, where z denotes the out-of-plane (vertical) distance between the AFM tip and the sample surface) and the sample was recorded and mapped spatially.
Electrochemical impedance spectroscopy
A Gamry Reference 3000 potentiostat was used to perform electrochemical impedance spectroscopy studies on the devices. Since our H-NNO devices are solid-state devices with two electrodes, we have connected one of them to the working electrode of the potentiostat, and the other electrode to the reference and counter electrode of the potentiostat. The measurements were done in the frequency range of 10 Hz to 1 MHz. The DC voltage applied was of 1.5 V amplitude, and a perturbation AC signal of 100 mV was applied.
Electrical measurement of NdNiO3 array
The electrical characterization of the NdNiO3 array was performed in a FormFactor Summit probe station using a Keithley 4200-SCS semiconductor analyser. A 4225-PMU module was used to generate electric pulses ranging from 500 ns to 10 μs. For real-time, high-precision current measurements during voltage pulse application, a 4225-RPM remote amplifier was utilized. To monitor the relaxation process following a voltage pulse, a constant 0.1 V voltage bias was applied and the current was measured in real time at 100 ns intervals. All electric pulses were applied to the Pd electrode, while the Au electrode was kept ground. The device resistance was extracted by fitting the current–voltage curve within the linear low voltage regime (−0.1 to +0.1 V). All measurements were conducted in air at room temperature.
Equivalent circuit modelling of Pd–Au and Pd–Pd devices
The circuit simulation based on Pd–Au and Pd–Pd devices was performed using Cadence. To emulate the hydrogen cloud movement beneath the Pd electrode, a Verilog-A compact model was developed, treating hydrogen migration as a resistive switching memory device with its resistance state modulated by applied voltage (positive RESET, negative SET). A voltage-driven resistance evolution mechanism was incorporated to capture the nonlinear I–V behaviour during large voltage sweeps, incorporating a scaling factor to adjust the resistance dynamics. For Pd–Au devices with only one hydrogen cloud per device, the model used a single resistive switching memory element to reflect unidirectional hydrogen migration. By contrast, Pd–Pd devices were modelled using two elements with opposite polarities connected in series to represent simultaneous hydrogen cloud expansion and shrinkage. The model, integrated with capacitive effects, was incorporated into a SPICE circuit simulation framework, enabling equivalent circuit analysis of both Pd–Au and Pd–Pd devices based on experimental pulse measurement data.
COMSOL simulation for spatial interactions in 2 × 3 Pd–Pd arrays
The spatial distribution of electric potential in the 2 × 3 Pd–Pd array was simulated using the COMSOL AC/DC solver module. The entire array was simplified as a 2D surface, with Pd electrodes modelled as conductive metal regions. The resistivity values for NNO and H-NNO were set to 2.5 × 10−6 Ω m and 8.85 Ω m, respectively. The simulated geometry was based on experimental dimensions, where each Pd pad was a square of 120 μm in length, with a 10 μm interval between adjacent Pd pads. The hydrogen cloud beneath each Pd electrode was assumed to be uniformly distributed and extended 3.5 μm beyond each side of the electrode, resulting in a total hydrogen cloud length of 17 μm. By applying different voltage configurations to each Pd electrode, the potential distribution across the nickelate film surface was computed using the two-dimensional model.
Large-scale Pd–Pd spatiotemporal processing layer modelling
For the large-scale neural network simulation of spoken digit recognition tasks, a 128-node Pd–Pd spatiotemporal processing array was implemented in Python. Each Pd node was initialized with a random hydrogen cloud thickness (x) between 2 μm and 2.5 μm, and its conductance (G) was calculated based on x using the equation provided in Supplementary Text 5. The spatiotemporal state evolution was updated every time step (500 ns) through three key processes: voltage potential evolution, hydrogen cloud evolution and current readout. The voltage potential across the NNO film was influenced by capacitive charging and discharging effects, which depend on the number of Pd nodes receiving a voltage spike at a given time step. This process introduced both nonlinear temporal characteristics and spatial interactions between nodes. Simultaneously, the hydrogen cloud thickness at each node was dynamically updated based on the local electric field, which is determined by the voltage difference between the Pd node and the NNO film, as well as the current cloud thickness. The combined effects of capacitive behaviour and hydrogen migration result in the final current accumulation and relaxation behaviours, which define the short-term memory and spatial interaction properties of the spatiotemporal processing layer. In the end, the system’s spatiotemporal response was captured by measuring the output current It under a small read voltage 0.1 V applied to half of the nodes. The resulting spatiotemporal processing layer output of each spoken digit signal had a size of (1, 128 × Nsample), where Nsample represents the number of times each signal is sampled.
Large-scale Pd–Au output layer training
A single-layer linear regression model is used as the output layer for spoken digit classification, implemented in Python. The trained weight values are directly linear mapped to the conductance states of the Pd–Au output layer array for hardware realization. The total spatiotemporal processing layer output X has a size of (Ntrain, 128 × Nsample), where Ntrain represents the total number of training samples and Nsample represents the number of times each signal is sampled. The target matrix Y has a size of (Ntrain, 10), where each row is a one-hot encoded representation of the corresponding digit label.
To evaluate the spatiotemporal processing layer’s performance, we use fivefold cross-validation. The dataset is split into training (Xtrain, Ytrain) and testing (Xtest, Ytest) sets using an 80–20 split. This process is repeated five times, ensuring that each subset serves as the test set once while the remaining four are used for training32.
For linear regression, the relationship between input and output is modelled using weights W of size (128 × Nsample, 10) and a bias vector b of size (1, 10). The predicted class probabilities Ŷ with size (Ntrain, 10) are computed as
$$\hat{Y}={XW}+{\textbf{b}}.$$
(1)
The final class label ŷ with size (Ntrain, 1) is determined using the argmax operation:
$$\,\hat{y}=\arg \mathop{\max }\limits_{i\,\in 1,2,3,4\ldots .10}\hat{{Y}_{i}}.$$
(2)
The optimal weights W and bias b are obtained by minimizing the mean squared error (MSE) loss, defined as
$$L=\frac{1}{{N}_{\mathrm{train}}}\times \,\mathop{\sum }\limits_{i=1}^{{N}_{\mathrm{train}}}{||{Y}_{i}-{\hat{Y}}_{i}||}^{2},$$
(3)
where Ntrain represents the number of training samples.
Seizure detection
CHB-MIT Scalp EEG Database (version 1.0.0) available from PhysioNet37,38,39 including EEG recordings sampled at 256 Hz from 23 electrodes placed according to the International 10–20 system was used. Our training dataset included EEG data from 5 patients, each contributing 20 seizure and 20 non-seizure 10 s EEG segments (a total of 200 signal segments). For clear comparison, both normalized seizure and normal signal clips from channels 13 to 15 are presented in Supplementary Fig. 17. The whole network set-up remained consistent with that used in the spoken digit recognition task, except for a reduced pad array size (from 128 to 46) to accommodate the smaller number of EEG input channels (23 versus 64 in the spoken digit dataset). The classification accuracy among different threshold set-ups is shown in Supplementary Fig. 18, where a 10 s EEG clip was sampled at 10 evenly spaced time points (that is, sampling interval of 1 s per point). As observed, there exists an optimal threshold setting that maximizes classification accuracy, which is expected: if the threshold is set too low, both seizure and non-seizure signals may generate excessive activations, introducing noise and reducing accuracy. Conversely, if the threshold is too high, useful information that differentiates the two classes is filtered out, also degrading performance.
