The question is a good one and has been asked also by Ehrenfest (1932): "Einige die Quantenmechanik betreffende Erkundigungsfragen". The answer was given by Pauli (1933): "Einige die Quantenmechanik betreffenden Erkundigungsfragen". Unfortunately I'm not aware of a english translation of these two publications. However one can find a slightly different form of the answer also in Pauli's book "General Principles of Quantum Mechanics" p.16. In that book Pauli writes
a single real function is not sufficient in order to construct from wavefunctions of the form (3.1) a non-negative probability function that is constant in time when integrated over the whole space.
I will try to summarize his arguments here:
A wave packet to describe a single particle (basically deBroglie's idea) can be written generally like $$ u(x,t) = \int U(k) e^{i(kx-\omega t)} dk = \int U(k) e^{ikx} dk \, e^{-i\omega t} $$ where $U(k)$ is the Fourier transform of $u(x,0)$. The complex conjugate of this wave packet is $$ u^*(x,t) = \int U^*(k) e^{-i(kx-\omega t)} dk = \int U^*(k) e^{-ikx} dk \, e^{i\omega t} $$ One can also define such wave packets in Electrodynamics. But in Quantum mechanics we have an additional condition, namely that the probability $P(x,t)$ to find a particle must always be positive and the total probability to find a single particle somewhere must be one, so $$ P(x,t) \geq 1 \\ \int P(x,t)\, dx = 1 $$ Pauli argues that the simplest ansatz to construct such a function $P(x,t)$ from $u(x,t)$ is a definite quadratic form from the functions $u$ and $u^*$, that means $$ P(x,t) = a u^2 + b {u^*}^2 + c u u^* $$ Now from the form of $u(x,t)$ and $u^*(x,t)$ we see that $$ u^2 \sim e^{-2i\omega t}\ \text{and}\ {u^*}^2 \sim e^{2i\omega t} $$ and an integral over space over these two functions can never be time independent. So the constants $a$ and $b$ must be zero in the ansatz for $P(x,t)$. Only the product of a wave packet and it's complex conjugate will yield a time independent total probability: $$ 1 = \int P(x,t)\, dx = \int uu^*\, dx = \iiint U(k)U^*(k') e^{i(kx-k'x)} \, e^{-i\omega t} e^{i\omega t} dk dx \\ = \iint U(k)U^*(k') \delta(k-k') dk dk' = \int \left|{U(k)}\right|^2 dk = \int P(k)\, dk $$ Since the product $uu^* = Re[u]^2 + Im[u]^2$ it follows that - as Pauli said - in order to compute a meaningful probability from wave packets of the form $u(x,t)$ one needs the real and the imaginary part of $u(x,t)$ and the wave function in quantum mechanics must be complex.