In this work we review the application of the theory of Gaussian processes to the modeling of noise in pulsar-timing data analysis, and we derive various useful and optimized representations for the likelihood expressions that are needed in Bayesian inference on pulsar-timing-array datasets. The resulting viewpoint and formalism lead us to two improved parameter-sampling schemes inspired by Gibbs sampling. The new schemes have vastly lower chain autocorrelation lengths than the Markov Chain Monte Carlo methods currently used in pulsar-timing data analysis, potentially speeding up Bayesian inference by orders of magnitude. The new schemes can be used for a full-noise-model analysis of the large datasets assembled by the International Pulsar Timing Array collaboration, which present a serious computational challenge to existing methods.