  Increasingly large and complex spatial datasets pose massive inferential
challenges due to high computational and storage costs. Our study is motivated
by the KAUST Competition on Large Spatial Datasets 2023, which tasked
participants with estimating spatial covariance-related parameters and
predicting values at testing sites, along with uncertainty estimates. We
compared various statistical and deep learning approaches through
cross-validation and ultimately selected the Vecchia approximation technique
for model fitting. To overcome the constraints in the R package GpGp, which
lacked support for fitting zero-mean Gaussian processes and direct uncertainty
estimation-two things that are necessary for the competition, we developed
additional \texttt{R} functions. Besides, we implemented certain
subsampling-based approximations and parametric smoothing for skewed sampling
distributions of the estimators. Our team DesiBoys secured victory in two out
of four sub-competitions, validating the effectiveness of our proposed
strategies. Moreover, we extended our evaluation to a large real spatial
satellite-derived dataset on total precipitable water, where we compared the
predictive performances of different models using multiple diagnostics.
