In this thesis mathematical models describing the growth of a solid tumour in the presence of an immune response are presented. Specifically, attention is focused on the interactions between cytotoxic T-lymphocytes (CTLs) and tumour cells in a small, avascular multicellular tumour. At this stage of the disease the CTLs and the tumour cells are considered to be in a state of dynamic equilibrium or cancer dormancy. The precise biochemical and cellular mechanisms by which CTLs can control a cancer and keep it in a dormant state are still not completely understood from a biological and immunological point of view. The mathematical models focus on the spatio-temporal dynamics of tumour cells, immune cells, chemokines and “chemo-repellors” in an immunogenic tumour. The CTLs and tumour cells are assumed to migrate and interact with each other in such a way that lymphocyte-tumour cell complexes are formed. These complexes result in either the death of the tumour cells (the normal situation) or the inactivation of the lymphocytes and consequently the survival of the tumour cells. In the latter case, we assume that each tumour cell which survives its “brief encounter” with the CTLs undergoes certain beneficial phenotypic changes. We explore the dynamics of the model under these assumptions and show that the process of the immuno-evasion can arise as a consequence of these encounters.Our computational simulations suggest that the proposed mechanism is able to mimic various dynamics of immunoevasion during the lifespan of a mouse. We also highlight the differential spatiotemporal contributions to evasion due, respectively, to: i) a decrease in the probability pi of being lethally hit; ii) a decrease in the probability, embedded in k+ i , that a tumour cell is recognized by a CTL. In particular, our model suggests that a decrease in the parameters pi is needed to produce evasion, which does not occur in the case where pi remains constant at its baseline level inferred from the experimental data. However, the role of the parameters k+ i is important since it can greatly accelerate the simulated process. Moreover, our computational simulations also show that the proposed mechanism can also deeply affect the spatial patterning of the tumour. In particular, our model suggests that to have a uniform invasion profile for the tumour cells necessitates also having a decrease in the recognition rate, embedded in the parameters k+ i . These parameters also differentially shape the spatial distribution of the various classes of tumour cells. Also in this thesis, we discuss mathematical models of the interactions between a tumour and both the innate and the cellular part of the adaptive immune system. We have developed and formulated spatiotemporal models of the interactions between macrophages, natural killer cells, cytotoxic T lymphocytes and tumour cells. In addition to presenting computational simulations of our ODE and PDE models, we investigate the linear stability analysis of steady states of the model and the effect of the initial conditions on the behaviour of the ODE solution. We show that limit cycle behaviour could be obtained by making some changes in the parameter values, which gave us oscillations in the solution of the ODE and PDE systems. We observe that there is a slowly damped oscillation in the behaviour of the tumour, natural killer and CTL cells. Also we note that the solution converges to the second steady state where the tumour size is small (dormant state).A model of cancer invasion and metastasis is also discussed in this thesis. This model attempts to describe the interactions between cancer cells, urokinase plasminogen activator (uPA), plasminogen activator inhibitor-1 (PAI-1), plasmin, extracellular matrix (ECM) and the immune response. The mathematical model focuses on the effect of the immune response on cancer invasion by assuming that there is some form of limit cycle behaviour between the cancer cells and the effector cells. The work we present in this chapter develops a mathematical model for tumour invasion with an immune response using a continuum model in 1 and 2 space dimensions. This model consists of a system of nonlinear partial differential equations and examines the effector cell response the tumour invasion. This model consists of effector cells, tumour cells, ECM, uPA, PAI-1, and plasmin. First, we set all spatial components of the model to zero and consider only the reaction kinetics in order to compare between the behaviour of our model and the original Chaplain and Lolas model (Chaplain and Lolas, 2005). The spatially homogeneous simulation shows the behaviour of solutions have regular oscillations because there is a closed orbit. Second, we present the computational results of the spatio-temporal model, and we note from these simulations that the tumour size of our model is smaller than the tumour size of the Chaplain and Lolas model because the immune cells are interacting with the tumour cells, and also the degradation of ECM is less than that in the Chaplain and Lolas model. In addition, the number of tumour cell clusters in our model is less than those in the Chaplain and Lolas model. Also we found the tumour clusters of the mathematical model which was discussed in this chapter to have the same range than the tumour clusters of the Chaplain and Lolas model. The final model presented in this thesis is a mathematical model of cancer cells and effector cells which exhibit standing-wave behaviour between them. We show tha the wave of invading cancer cells can be stopped by the wave of effector cells or ECM.This model also focuses on the effect of the mutation of cancer cells to another subpopulation which is more malignant and which has the ability to invade the ECM or the effector cells to occupy space. The numerical simulations discussed in this chapter are essentially associated with an initial model of two equations representing the effector cells and tumour cells, such that there is a standing wave between these species. We note that the solution of the mathematical model is a travelling wave and also has a standing wave solution (i.e. the wave of effector cells stops the wave of tumour cells when they meet). This phenomenon occurs when the two diffusion coefficients are the same. We calculate the wave speed to illustrate that the speed tends to zero when the two waves meet - a positive speed of tumour cells refers to an invading tumour, and a negative speed refers to the decreasing of effector cells. After this we modify the model by adding an equation for a second cancer cell population T2, which is a sub-population 2 of tumour cells. This is to reflect the fact that cancer is a progressive disease, and as such it becomes more malignant as the cancer cells undergo successive mutations. We show in this case how the new type of cancer cells start to invade the effector cells after the failure of the first type. The third model discussed in this chapter is arrived at by adding an ECM equation to the second model, and it explains how the standing wave arise from two types of equations - the first one contains diffusion, and the second one has no diffusion.All the mathematical models in this thesis use numerical analysis of nonlinear partial differential equations and computational simulations to obtain insight into the underlying biological systems. The systems of nonlinear partial differential equations were numerically solved by a PDE solver in MATLAB for 1D and COMSOL for 2D. We used the MATLAB PDE solver pdepe which uses the method described in Skeel and Berzins (1990) for the spatial discretisation and the MATLAB routine ode15s for the time integration.The numerical simulations demonstrate the existence of cell distributions that are quasi-stationary in time and heterogeneous in space.