Synaptic communication is based on a biological Molecular Communication (MC) system which may serve as a blueprint for the design of synthetic MC systems. However, the physical modeling of synaptic MC is complicated by the possible saturation of the molecular receiver caused by the competition of neurotransmitters (NTs) for postsynaptic receptors. Receiver saturation renders the system behavior nonlinear in the number of released NTs and is commonly neglected in existing analytical models. Furthermore, due to the ligands' competition for receptors (and vice versa), the individual binding events at the molecular receiver are in general not statistically independent and the commonly used binomial model for the statistics of the received signal does not apply. Hence, in this work, we propose a novel deterministic model for receptor saturation in terms of a state-space description based on an eigenfunction expansion of Fick's diffusion equation. The presented solution is numerically stable and computationally efficient. Employing the proposed deterministic model, we show that saturation at the molecular receiver effectively reduces the peak-value of the expected received signal and accelerates the clearance of NTs as compared to the case when receptor occupancy is neglected. We further derive a statistical model for the received signal in terms of the hypergeometric distribution which accounts for the competition of NTs for receptors and the competition of receptors for NTs. The proposed statistical model reveals how the signal statistics are shaped by the number of released NTs, the number of receptors, and the binding kinetics of the receptors, respectively, in the presence of competition. In particular, we show that the impact of these parameters on the signal variance is qualitatively different depending on the relative numbers of NTs and receptors. Finally, the accuracy of the proposed deterministic and statistical models is verified by particle-based computer simulations.