Gaussian process (GP) regression provides a flexible, nonparametric framework for probabilistic modeling, yet remains computationally demanding in large-scale applications. For one-dimensional data, state space (SS) models achieve linear-time inference by reformulating GPs as stochastic differential equations (SDEs). However, SS approaches are confined to gridded inputs and cannot handle multi-dimensional scattered data. We propose a new framework based on kernel packet (KP), which overcomes these limitations while retaining exactness and scalability. A KP is a compactly supported function defined as a linear combination of the GP covariance functions. In this article, we prove that KPs can be identified via the forward and backward SS representations. We also show that the KP approach enables exact inference with linear-time training and logarithmic or constant-time prediction, and extends naturally to multi-dimensional gridded or scattered data without low-rank approximations. Numerical experiments on large-scale additive and product-form GPs with millions of samples demonstrate that KPs achieve exact, memory-efficient inference where SDE-based and low-rank GP methods fail.