Functional data analysis is the study of random phenomena which live on a continuum and are probabilistically modelled as random functions. Through studying the statistical properties of the underlying stochastic process which generates the functional data, one can build a model for the prediction, forecasting, and uncertainty quantification of future realizations of the random process. Heteroscedastic models play a key role in improving uncertainty quantification, and the study of such models in functional data analysis is a relatively new area. This thesis presents four different heteroscedastic functional data models and shows how each model can be used to improve uncertainty estimates in different applications. Chapter 2 examines the functional stochastic volatility model. It formulates an analogue of the discrete-time stochastic volatility model for functional time series and examines its Bayesian inference. It then demonstrates in an application to SPX option surfaces how the functional stochastic volatility model can improve quantile estimates for portfolio risk management compared to a constant volatility model. Chapter 3 examines a generalized autoregressive conditional heteroscedasticity (GARCH)-based functional time series model which is used to improve forecasts of maximum daily temperature for 111 cities across the United States. The functional model uncovers interpretable spatial effects that group geographical regions of the country, while the time series aspect exploits time dependence of weather forecast errors to de-bias the forecasts and reduce variance. Chapter 4 models phase uncertainty in functional data with Dirichlet warping functions. The data set consists of a collection of curves drawn from an aerospace dynamical system which exhibit both horizontal phase variation and vertical amplitude variation. The proposed model decouples these two sources of variation for a parsimonious probabilistic model that can be used to generate new sample curves and bypass experimental costs associated with sampling curve data. Chapter 5 studies a prior distribution on stationary covariance functions built from Lévy processes. This model is particularly suited towards modelling quasiperiodic functional data where the Fourier spectrum exhibits localized peaks that can be captured by linear combinations of a parameterized basis function. The resulting model improves the uncertainty quantification of long-range quasiperiodic forecasts compared to a Gaussian process with a fixed covariance function.