Fluid antenna systems (FAS) enable unprecedented spatial diversity within a compact form factor by flexibly switching among high-density antenna ports. To activate this capability, channel state information (CSI) over the ports is required, which implies high estimation overhead because the number of ports is usually very large. Conventional estimation schemes tend to first estimate the CSI for a small number of ports and then infer the CSI for the remaining antenna ports by interpolation exploiting correlation characteristics. However, existing correlation-based techniques lack generalization ability, and the fundamental limits of interpolating the CSI from sparse observations remain poorly understood. This paper adopts a generative modeling framework for characterizing the channel correlation among the FAS ports that departs fundamentally from covariance-descriptive models. Specifically, we represent the spatially sampled channel as a $p$th-order autoregressive (AR) Gauss-Markov process, which provides a principled and tunable tradeoff between model complexity and approximation accuracy via the AR order. In so doing, we can characterize the limits of channel interpolation by deriving the globally optimal minimum mean-square error (MMSE) estimator and establishing a tight lower bound on the minimum number of observations required to meet a prescribed reconstruction error. To reduce the complexity of MMSE estimation, we then exploit the state-space structure due to the ${\rm AR}(p)$ model and develop a Kalman filtering/smoothing-based interpolation algorithm. The resulting method attains the optimal MMSE performance with strictly linear complexity $\mathcal{O}(N)$ with $N$ denoting the number of ports, resulting in a scalable, efficient, and theoretically grounded framework for practical FAS channel reconstruction.