Low-rank Approximation Research Articles

Data‐driven linear complexity low‐rank approximation of general kernel matrices: A geometric approach

SummaryA general, rectangular kernel matrix may be defined as where is a kernel function and where and are two sets of points. In this paper, we seek a low‐rank approximation to a kernel matrix where the sets of points and are large and are arbitrarily distributed, such as away from each other, “intermingled”, identical, and so forth. Such rectangular kernel matrices may arise, for example, in Gaussian process regression where corresponds to the training data and corresponds to the test data. In this case, the points are often high‐dimensional. Since the point sets are large, we must exploit the fact that the matrix arises from a kernel function, and avoid forming the matrix, and thus ruling out most algebraic techniques. In particular, we seek methods that can scale linearly or nearly linearly with respect to the size of data for a fixed approximation rank. The main idea in this paper is to geometrically select appropriate subsets of points to construct a low rank approximation. An analysis in this paper guides how this selection should be performed.

Hyperspectral image reconstruction based on low-rank coefficient tensor and global prior

ABSTRACT Compressive sensing (CS) has been employed to compress and store hyperspectral images (HSI) to transfer extensive data efficiently, and obtaining high-quality reconstructed images is crucial for subsequent applications. Model-based methods usually capture the priors of HSI, including sparsity, global and nonlocal low-rankness, nonlocal self-similarity and global correlation. However, it is challenging to reasonably integrate different regularizations in a united framework and exploit various priors mentioned to reconstruct high-quality images. In this paper, a tensor optimization model based on the low-rank coefficient tensor and global prior (LRCTGP) is proposed for HSI reconstruction. Complementary regularization terms are integrated into a united framework and effectively promote the reconstruction results. First, we apply the mode- tensor-matrix product to decompose the original HSI rather than employing the regularization to capture the spectral low-rankness and simultaneously constrain the original and feature HSI. Then, tensor ring decomposition is employed to constrain the coefficient tensor with fewer bands, which is more efficient in the low-rank approximation and has smaller calculation costs than applying it to the original HSI. Moreover, integrating BM3D as a regularizer is more efficient than other patch-based models. Finally, considering the smoothness and global correlation, spatial-spectral total variation (SSTV) is applied to compensate for the shortcomings after decomposing the original and feature HSI, which complements the overall structure and details of reconstructed images and improves the reconstruction quality. Alternating direction method of multipliers (ADMM) is used to optimize the proposed model. Experimental results of the LRCTGP model on different HSI datasets are better than existing state-of-the-art approaches, which proves the superiority.

Robust Low-Rank Graph Multi-View Clustering via Cauchy Norm Minimization

Graph-based multi-view clustering methods aim to explore the partition patterns by utilizing a similarity graph. However, many existing methods construct a consensus similarity graph based on the original multi-view space, which may result in the lack of information on the underlying low-dimensional space. Additionally, these methods often fail to effectively handle the noise present in the graph. To address these issues, a novel graph-based multi-view clustering method which combines spectral embedding, non-convex low-rank approximation and noise processing into a unit framework is proposed. In detail, the proposed method constructs a tensor by stacking the inner product of normalized spectral embedding matrices obtained from each similarity matrix. Then, the obtained tensor is decomposed into a low-rank tensor and a noise tensor. The low-rank tensor is constrained via nonconvex low-rank tensor approximation and a novel Cauchy norm with an upper bound is proposed to handle the noise. Finally, we derive the consensus similarity graph from the denoised low-rank tensor. The experiments on five datasets demonstrate that the proposed method outperforms other state-of-the-art methods on five datasets.

High-Throughput Condensed-Phase Hybrid Density Functional Theory for Large-Scale Finite-Gap Systems: The SeA Approach.

High-throughput electronic structure calculations (often performed using density functional theory (DFT)) play a central role in screening existing and novel materials, sampling potential energy surfaces, and generating data for machine learning applications. By including a fraction of exact exchange (EXX), hybrid functionals reduce the self-interaction error in semilocal DFT and furnish a more accurate description of the underlying electronic structure, albeit at a computational cost that often prohibits such high-throughput applications. To address this challenge, we have constructed a robust, accurate, and computationally efficient framework for high-throughput condensed-phase hybrid DFT and implemented this approach in the PWSCF module of Quantum ESPRESSO (QE). The resulting SeA approach (SeA = SCDM + exx + ACE) combines and seamlessly integrates: (i) the selected columns of the density matrix method (SCDM, a robust noniterative orbital localization scheme that sidesteps system-dependent optimization protocols), (ii) a recently extended version of exx (a black-box linear-scaling EXX algorithm that exploits sparsity between localized orbitals in real space when evaluating the action of the standard/full-rank operator), and (iii) adaptively compressed exchange (ACE, a low-rank approximation). In doing so, SeA harnesses three levels of computational savings: pair selection and domain truncation from SCDM + exx (which only considers spatially overlapping orbitals on orbital-pair-specific and system-size-independent domains) and low-rank approximation from ACE (which reduces the number of calls to SCDM + exx during the self-consistent field (SCF) procedure). Across a diverse set of 200 nonequilibrium (H2O)64 configurations (with densities spanning 0.4-1.7 g/cm3), SeA provides a 1-2 order-of-magnitude speedup in the overall time-to-solution, i.e., ≈8-26× compared to the convolution-based PWSCF(ACE) implementation in QE and ≈78-247× compared to the conventional PWSCF(Full) approach, and yields energies, ionic forces, and other properties with high fidelity. As a proof-of-principle high-throughput application, we trained a deep neural network (DNN) potential for ambient liquid water at the hybrid DFT level using SeA via an actively learned data set with ≈8,700 (H2O)64 configurations. Using an out-of-sample set of (H2O)512 configurations (at nonambient conditions), we confirmed the accuracy of this SeA-trained potential and showcased the capabilities of SeA by computing the ground-truth ionic forces in this challenging system containing >1,500 atoms.

Optical coherence tomography image despeckling based on tensor singular value decomposition and fractional edge detection

Optical coherence tomography (OCT) imaging is a technique that is frequently used to diagnose medical conditions. However, coherent noise, sometimes referred to as speckle noise, can dramatically reduce the quality of OCT images, which has an adverse effect on how OCT images are used. In order to enhance the quality of OCT images, a speckle noise reduction technique is developed, and this method is modelled as a low-rank tensor approximation issue. The grouped 3D tensors are first transformed into the transform domain using tensor singular value decomposition (t-SVD). Then, to cut down on speckle noise, transform coefficients are thresholded. Finally, the inverse transform can be used to produce images with speckle suppression. To further enhance the despeckling results, a feature-guided thresholding approach based on fractional edge detection and an adaptive backward projection technique are also presented. Experimental results indicate that the presented algorithm outperforms several comparison methods in relation to speckle suppression, objective metrics, and edge preservation.

Hierarchical off-diagonal low-rank approximation of Hessians in inverse problems, with application to ice sheet model initialization

Obtaining lightweight and accurate approximations of discretized objective functional Hessians in inverse problems governed by partial differential equations (PDEs) is essential to make both deterministic and Bayesian statistical large-scale inverse problems computationally tractable. The cubic computational complexity of dense linear algebraic tasks, such as Cholesky factorization, that provide a means to sample Gaussian distributions and determine solutions of Newton linear systems is a computational bottleneck at large-scale. These tasks can be reduced to log-linear complexity by utilizing hierarchical off-diagonal low-rank (HODLR) matrix approximations. In this work, we show that a class of Hessians that arise from inverse problems governed by PDEs are well approximated by the HODLR matrix format. In particular, we study inverse problems governed by PDEs that model the instantaneous viscous flow of ice sheets. In these problems, we seek a spatially distributed basal sliding parameter field such that the flow predicted by the ice sheet model is consistent with ice sheet surface velocity observations. We demonstrate the use of HODLR Hessian approximation to efficiently sample the Laplace approximation of the posterior distribution with covariance further approximated by HODLR matrix compression. Computational studies are performed which illustrate ice sheet problem regimes for which the Gauss–Newton data-misfit Hessian is more efficiently approximated by the HODLR matrix format than the low-rank (LR) format. We then demonstrate that HODLR approximations can be favorable, when compared to global LR approximations, for large-scale problems by studying the data-misfit Hessian associated with inverse problems governed by the first-order Stokes flow model on the Humboldt glacier and Greenland ice sheet.

Compressing Transformers: Features Are Low-Rank, but Weights Are Not!

Transformer and its variants achieve excellent results in various computer vision and natural language processing tasks, but high computational costs and reliance on large training datasets restrict their deployment in resource-constrained settings. Low-rank approximation of model weights has been effective in compressing CNN models, but its application to transformers has been less explored and is less effective. Existing methods require the complete dataset to fine-tune compressed models, which are both time-consuming and data-hungry. This paper reveals that the features (i.e., activations) are low-rank, but model weights are surprisingly not low-rank. Hence, AAFM is proposed, which adaptively determines the compressed model structure and locally compresses each linear layer's output features rather than the model weights. A second stage, GFM, optimizes the entire compressed network holistically. Both AAFM and GFM only use few training samples without labels, that is, they are few-shot, unsupervised, fast and effective. For example, with only 2K images without labels, 33% of the parameters are removed in DeiT-B with 18.8% relative throughput increase, but only a 0.23% accuracy loss for ImageNet recognition. The proposed methods are successfully applied to the language modeling task in NLP, too. Besides, the few-shot compressed models generalize well in downstream tasks.

Spectral Feature Augmentation for Graph Contrastive Learning and Beyond

Although augmentations (e.g., perturbation of graph edges, image crops) boost the efficiency of Contrastive Learning (CL), feature level augmentation is another plausible, complementary yet not well researched strategy. Thus, we present a novel spectral feature argumentation for contrastive learning on graphs (and images). To this end, for each data view, we estimate a low-rank approximation per feature map and subtract that approximation from the map to obtain its complement. This is achieved by the proposed herein incomplete power iteration, a non-standard power iteration regime which enjoys two valuable byproducts (under mere one or two iterations): (i) it partially balances spectrum of the feature map, and (ii) it injects the noise into rebalanced singular values of the feature map (spectral augmentation). For two views, we align these rebalanced feature maps as such an improved alignment step can focus more on less dominant singular values of matrices of both views, whereas the spectral augmentation does not affect the spectral angle alignment (singular vectors are not perturbed). We derive the analytical form for: (i) the incomplete power iteration to capture its spectrum-balancing effect, and (ii) the variance of singular values augmented implicitly by the noise. We also show that the spectral augmentation improves the generalization bound. Experiments on graph/image datasets show that our spectral feature augmentation outperforms baselines, and is complementary with other augmentation strategies and compatible with various contrastive losses.

Beyond Graph Convolutional Network: An Interpretable Regularizer-Centered Optimization Framework

Graph convolutional networks (GCNs) have been attracting widespread attentions due to their encouraging performance and powerful generalizations. However, few work provide a general view to interpret various GCNs and guide GCNs' designs. In this paper, by revisiting the original GCN, we induce an interpretable regularizer-centerd optimization framework, in which by building appropriate regularizers we can interpret most GCNs, such as APPNP, JKNet, DAGNN, and GNN-LF/HF. Further, under the proposed framework, we devise a dual-regularizer graph convolutional network (dubbed tsGCN) to capture topological and semantic structures from graph data. Since the derived learning rule for tsGCN contains an inverse of a large matrix and thus is time-consuming, we leverage the Woodbury matrix identity and low-rank approximation tricks to successfully decrease the high computational complexity of computing infinite-order graph convolutions. Extensive experiments on eight public datasets demonstrate that tsGCN achieves superior performance against quite a few state-of-the-art competitors w.r.t. classification tasks.

Synthesis of sparse planar arrays with multiple patterns: An automatic coordinate pairing approach

An innovative algorithm to synthesize a non-uniform sparse planar array with as few elements as possible is proposed in this paper. The algorithm is aimed at multiple different beampatterns, and can directly estimate the exact positions of the planar array elements, avoiding the problem of coordinate pairing between X-axis and Y-axis. To achieve this, the proposed method utilizes the rotational invariance caused by the translation invariant structure of the steering vector. It derives the array element position from the data matrix based on unitary ESPRIT and employs a low-rank approximation operation to construct a matrix containing the position information of the reconstructed sparse array elements. Automatic pairing is then realized by assigning the real and imaginary parts of the matrix eigenvalues to the abscissa and ordinate of elements one-to-one. The paper presents two different simulation scenarios to demonstrate the performance of the proposed algorithm, with the results showing a reduction in the number of array elements by 20%-30% while ensuring low reconstruction error, low time, and space complexity. Overall, the proposed algorithm provides an effective solution for synthesizing non-uniform sparse planar arrays with fewer elements, which can significantly reduce costs and complexity in various applications.

Denoising for airborne transient electromagnetic data using noise-whitening-based weighted nuclear norm minimization

Abstract Recovering anomalous information covered under noise in late gates can enhance airborne transient electromagnetic (ATEM) detection. Conventional denoising mainly comprises filtering and gate correlation-based decomposition algorithms; the former fails to extract anomalies contaminated by noise and the latter relies on the correlation between gates, which may yield false late gate anomalies caused by early large-amplitude anomalies in early gates. In ATEM profiles, the correlation between anomalies in adjacent gates makes the anomalies to be measured to have low-rank characteristics relative to the noise-contaminated profiles; the noise is uniformly distributed in the profiles, which have nonlocal self-similarity. Therefore, the low-rank matrix approximation algorithm is applicable to ATEM data denoising. In this study, an algorithm—noise-whitening-based weighted nuclear norm minimization (NW-WNNM)—is designed to remove ATEM profile noise. First, we analyze the influence of patch size in block matching on anomalous and noisy patches and estimate the profile patch size adaptively. Then, we combine the estimation of noise variance in weighted nuclear norm minimization (WNNM) with the noise whitening of similar patch matrices to reduce the noise interference on the nuclear norm and add a whitening factor in the weight vector to make the soft-thresholding function applicable to the low-rank reconstruction of the whitened matrix. By analyzing the reconstructed low-rank matrix and its feature distribution, compared with WNNM, NW-WNNM can detect the feature information more accurately and eliminate the influence of noise on the nuclear norm. Simulation and field profile results indicate that NW-WNNM is superior to comparison denoising methods.

HIPPYlib-MUQ: A Bayesian Inference Software Framework for Integration of Data with Complex Predictive Models under Uncertainty

Bayesian inference provides a systematic framework for integration of data with mathematical models to quantify the uncertainty in the solution of the inverse problem. However, the solution of Bayesian inverse problems governed by complex forward models described bypartial differential equations (PDEs)remains prohibitive with black-boxMarkov chain Monte Carlo (MCMC)methods. We present hIPPYlib-MUQ, an extensible and scalable software framework that contains implementations of state-of-the art algorithms aimed to overcome the challenges of high-dimensional, PDE-constrained Bayesian inverse problems. These algorithms accelerate MCMC sampling by exploiting the geometry and intrinsic low-dimensionality of parameter space via derivative information and low rank approximation. The software integrates two complementary open-source software packages, hIPPYlib and MUQ. hIPPYlib solves PDE-constrained inverse problems using automatically-generated adjoint-based derivatives, but it lacks full Bayesian capabilities. MUQ provides a spectrum of powerful Bayesian inversion models and algorithms, but expects forward models to come equipped with gradients and Hessians to permit large-scale solution. By combining these two complementary libraries, we created a robust, scalable, and efficient software framework that realizes the benefits of each and allows us to tackle complex large-scale Bayesian inverse problems across a broad spectrum of scientific and engineering disciplines. To illustrate the capabilities of hIPPYlib-MUQ, we present a comparison of a number of MCMC methods available in the integrated software on several high-dimensional Bayesian inverse problems. These include problems characterized by both linear and nonlinear PDEs, various noise models, and different parameter dimensions. The results demonstrate that large (∼ 50×) speedups over conventional black box and gradient-based MCMC algorithms can be obtained by exploiting Hessian information (from the log-posterior), underscoring the power of the integrated hIPPYlib-MUQ framework.

An Unsupervised Image Denoising Method Using a Nonconvex Low-Rank Model with TV Regularization

In real-world scenarios, images may be affected by additional noise during compression and transmission, which interferes with postprocessing such as image segmentation and feature extraction. Image noise can also be induced by environmental variables and imperfections in the imaging equipment. Robust principal component analysis (RPCA), one of the traditional approaches for denoising images, suffers from a failure to efficiently use the background’s low-rank prior information, which lowers its effectiveness under complex noise backgrounds. In this paper, we propose a robust PCA method based on a nonconvex low-rank approximation and total variational regularization (TV) to model the image denoising problem in order to improve the denoising performance. Firstly, we use a nonconvex γ-norm to address the issue that the traditional nuclear norm penalizes large singular values excessively. The rank approximation is more accurate than the nuclear norm thanks to the elimination of matrix elements with substantial approximation errors to reduce the sparsity error. The method’s robustness is improved by utilizing the low sensitivity of the γ-norm to outliers. Secondly, we use the l1-norm to increase the sparsity of the foreground noise. The TV norm is used to improve the smoothness of the graph structure in accordance with the sparsity of the image in the gradient domain. The denoising effectiveness of the model is increased by employing the alternating direction multiplier strategy to locate the global optimal solution. It is important to note that our method does not require any labeled images, and its unsupervised denoising principle enables the generalization of the method to different scenarios for application. Our method can perform denoising experiments on images with different types of noise. Extensive experiments show that our method can fully preserve the edge structure information of the image, preserve important features of the image, and maintain excellent visual effects in terms of brightness smoothing.

Low tensor train and low multilinear rank approximations of 3D tensors for compression and de-speckling of optical coherence tomography images

Objective. Many methods for compression and/or de-speckling of 3D optical coherence tomography (OCT) images operate on a slice-by-slice basis and, consequently, ignore spatial relations between the B-scans. Thus, we develop compression ratio (CR)-constrained low tensor train (TT)—and low multilinear (ML) rank approximations of 3D tensors for compression and de-speckling of 3D OCT images. Due to inherent denoising mechanism of low-rank approximation, compressed image is often even of better quality than the raw image it is based on. Approach. We formulate CR-constrained low rank approximations of 3D tensor as parallel non-convex non-smooth optimization problems implemented by alternating direction method of multipliers of unfolded tensors. In contrast to patch- and sparsity-based OCT image compression methods, proposed approach does not require clean images for dictionary learning, enables CR as high as 60:1, and it is fast. In contrast to deep networks based OCT image compression, proposed approach is training free and does not require any supervised data pre-processing. Main results. Proposed methodology is evaluated on twenty four images of a retina acquired on Topcon 3D OCT-1000 scanner, and twenty images of a retina acquired on Big Vision BV1000 3D OCT scanner. For the first dataset, statistical significance analysis shows that for CR ≤ 35, all low ML rank approximations and Schatten-0 (S 0) norm constrained low TT rank approximation can be useful for machine learning-based diagnostics by using segmented retina layers. Also for CR ≤ 35, S 0-constrained ML rank approximation and S 0-constrained low TT rank approximation can be useful for visual inspection-based diagnostics. For the second dataset, statistical significance analysis shows that for CR ≤ 60 all low ML rank approximations as well as S 0 and S 1/2 low TT ranks approximations can be useful for machine learning-based diagnostics by using segmented retina layers. Also, for CR ≤ 60, low ML rank approximations constrained with S p , p ∊ {0, 1/2, 2/3} and one surrogate of S 0 can be useful for visual inspection-based diagnostics. That is also true for low TT rank approximations constrained with S p , p ∊ {0, 1/2, 2/3} for CR ≤ 20. Significance. Studies conducted on datasets acquired by two different types of scanners confirmed capabilities of proposed framework that, for a wide range of CRs, yields de-speckled 3D OCT images suitable for clinical data archiving and remote consultation, for visual inspection-based diagnosis and for machine learning-based diagnosis by using segmented retina layers.

Randomized Low-Rank Approximation of Monotone Matrix Functions

Randomized Low-Rank Approximation of Monotone Matrix Functions

On Best Low Rank Approximation of Positive Definite Tensors

Tensors, or multi-indexed arrays, play an important role in machine learning and signal processing. These higher-order generalizations of matrices allow for preservation of higher-order structure in data, and low rank decompositions of tensors allow for recovery of underlying information. In many cases, the underlying tensor of interest is positive (semi)definite, i.e., the homogeneous polynomial associated to the tensor has nonnegative evaluation on all inputs. An archetypal problem is that one has a noisy measurement of some low rank signal tensor of interest. This measurement itself does not have low rank, so one must compute a best low rank approximation of the measured tensor to recover component information. As it turns out, the set of tensors of rank less than or equal to is, in general, not closed when , and, as a consequence, best low rank approximations can fail to exist. In the case when a best low rank approximation does not exist, near optimal low rank approximations suffer numerical issues and cannot be used to reliably approximate underlying component information. As a consequence, existence guarantees for best low rank tensor approximations are of great practical and theoretical interest. This article develops deterministic guarantees for the existence of best low rank approximations of tensors which are positive semidefinite. In particular, we show that the set of low rank positive semidefinite tensors is relatively closed as a subset of the set of positive semidefinite tensors. We use this fact to give a deterministic bound which may be used to guarantee the existence of a best low rank approximation of a noisy low rank positive semidefinite tensor. Furthermore, we show that this condition guarantees uniqueness of the canonical polyadic decomposition for best low rank approximations. In addition, for order three tensors, we prove that our bound is sharp and that it can be computed using semidefinite programming. The main results of this article are illustrated through numerical experiments, which show that our bound is highly predictive of numerical issues when attempting to recover underlying component information from a noisy low rank positive semidefinite tensor.

Fast truncated SVD of sparse and dense matrices on graphics processors

We investigate the solution of low-rank matrix approximation problems using the truncated singular value decomposition (SVD). For this purpose, we develop and optimize graphics processing unit (GPU) implementations for the randomized SVD and a blocked variant of the Lanczos approach. Our work takes advantage of the fact that the two methods are composed of very similar linear algebra building blocks, which can be assembled using numerical kernels from existing high-performance linear algebra libraries. Furthermore, the experiments with several sparse matrices arising in representative real-world applications and synthetic dense test matrices reveal a performance advantage of the block Lanczos algorithm when targeting the same approximation accuracy.

Efficient algorithms for computing rank‐revealing factorizations on a GPU

AbstractStandard rank‐revealing factorizations such as the singular value decomposition (SVD) and column pivoted QR factorization are challenging to implement efficiently on a GPU. A major difficulty in this regard is the inability of standard algorithms to cast most operations in terms of the Level‐3 BLAS. This article presents two alternative algorithms for computing a rank‐revealing factorization of the form , where and are orthogonal and is trapezoidal (or triangular if is square). Both algorithms use randomized projection techniques to cast most of the flops in terms of matrix‐matrix multiplication, which is exceptionally efficient on the GPU. Numerical experiments illustrate that these algorithms achieve significant acceleration over finely tuned GPU implementations of the SVD while providing low rank approximation errors close to that of the SVD.

Robust Matrix Completion Method Based on TNNR and Total Row Difference for Recovering Optical Image

Matrix completion aims to recover a matrix from an incomplete matrix with many unknown elements and has wide applications in optical image recovery and machine learning, in which the popular method is to formulate it as a general low-rank matrix approximation problem. However, the traditional optimization model for matrix completion is less robust. This paper proposes a robust matrix completion method in which the truncated nuclear norm regularization (TNNR) is used as the approximation of the rank function and the sum of absolute values of the row difference, which is called as the total row difference, is used to constrain the oscillations of the missing matrix. By minimizing the value of total row difference in the objective, the proposed model controls the oscillation and reduces the impact of missing part in the process of matrix completion continuously. Furthermore, we propose a two-step iterative algorithm framework and design an ADMM algorithm for subproblem model that includes minimizing the total row difference. Experiments show that the proposed algorithm has more stable performance and better recovery effect and obviously reduces the sensitivity of the traditional TNNR models to the truncated rank parameter.

Multistep-Ahead Chaotic Time Series Prediction Based on Hierarchical Echo State Network With Augmented Random Features

Multistep-ahead chaotic time series prediction is a kind of highly nonlinear problem, which puts forward higher requirements both for the dynamical memory and nonlinearity of the model. Echo state network (ESN) is frequently employed in the realm of chaotic time series modeling and prediction, but the basic ESN has been proved to have an antagonistic trade-off between nonlinear transformation and memory capacity. To overcome this trade-off, a new architecture named hierarchical echo state network with augmented random features (HESN-ARF) is proposed. On the basis of traditional linear random projection, the proposed HESN-ARF further leverages nonlinear kernel transformation to construct augmented random features, which can enable the linear and nonlinear properties to be fully represented. Moreover, the HESN-ARF utilizes low-rank kernel approximation to further reduce the computational cost, preserving the advantage of efficient modeling as much as possible while ensuring the capacities of nonlinear transformation and dynamical memory simultaneously. The proposed HESN-ARF can mine and learn the latent evolution patterns hidden in the dynamic system layer by layer through the hierarchical strategy, and achieves excellent performance in multi-step-ahead chaotic time series prediction, as demonstrated by experimental findings on two synthetic chaotic systems and a real-world meteorological dataset.