Phases of Atomic or Molecular Orbitals
Complex Numbers

As you point out, general a wavefunction if complex-valued.  A complex number and be described by two real numbers and hence are of depicted on the 2D complex plane (sometimes called the argand plane).

Numbers in the complex plane can be expressed in multiple forms.  The one usually taught first is as the sum of a real part and imaginary part:

$$\begin{align} z = x + iy\;{} &amp; ; x,y \in \mathbb{R}\\
\Re(z) = x\; {} &amp; ; \Im(z) = y \\
|z| = {} &amp; \sqrt{x^2+y^2}
\end{align}$$

Alternatively it can be described by a magnitude ($r$) and a phase ($\phi$):

$$\begin{align}
z = r e^{i\phi}\; {} &amp; ;r\in\mathbb{R}^{+0}, \phi\in[0:2\pi)\\
\Re(z) = r \cos(\phi) \; {} &amp;; \Im(z) = r \sin(\phi) \\
|z|{} &amp; = r
\end{align}$$



Two signals or functions with the same magnitude are said to be perfectly/exactly 'out of phase' if there phases differ by $\pi$ as $e^{i \pi}=-1, and perfectly 'in phase' if they have the same phase.

Orbitals

Now we can apply this to some orbitals.  Consider a cross-section of the 2p$_{\text{x}}$ orbital in the $x$-$y$ plane through $z=0$:

$$\psi_{2px}(x,y,z) = r(x,y,z) e^{i \phi(x,y,z)} \propto x \exp\left(-\sqrt{x^2+y^2+z^2}\right)$$



The intensity of the colour represents magnitude, $r(x,y,z)$, of the wavefunction and the colour itself represents the phase - red for $\phi=0$, blue for $\phi=\pi$, with a colourless node separating the two lobes.

As the probability density is equal to the square modulus of the wavefunction: $$\rho(x,y,z) = |\psi(x,y,z)|^2$$
we can multiply the wavefunction by a constant phase factor $e^{i\theta}$ and still represent the same particle:

Consider:
$$\begin{align}
\tilde{\psi} {} &amp; = \psi e^{i\theta} \\
\to \tilde{\rho} {} &amp; = |\tilde{\psi}|^2 \\
&amp; = |\psi|^2 e^{i\theta} e^{-i\theta} \\
&amp; = |\psi|^2 \\
&amp; = \rho
\end{align}$$
So both $\psi$ and $\tilde{\psi}$ are valid wavefunctions with either one is, and the wavefunction can only be defined up to multiplicative phase factor, so we can talk about absolute phase, but we can talk about relative phase and phase differences.

Consider the out of phase 2p$_{\text{x}}$ orbital:
$$\tilde{\psi}_{2px}(x,y,z) = \psi_{2px}(x,y,z)e^{i\pi} = - \psi_{2px}(x,y,z) $$


Molecular Orbitals

As we just showed they both have the same probability density, but when we begin to combine them, the phases dictate how much the wavefunctions constructively/destructively interfere.  We can combine linear combinations of atomic orbitals to make molecular orbitals.

Consider adding together two adjacent wavefunctions with the same phase on different atoms:



The adjacent regions of different phase cancel one another out.  When the probability density is examined:



we can see than the chance of the particle being between the two atoms is very low and it is localised around each nuclei.  This would be an anti-bonding molecular orbital.

On the other hand, adding together two adjacent wavefunctions with the opposite phase on different atoms:



The adjacent regions of similar phase add together.  When the probability density is examined:



we can see than the chance of the particle being between the two atoms is very much higher and it is localised in a bond.  This would be a bonding molecular orbital.

This is why we talk about in phase and out of phase orbital when talking about bonding - the regions of in phase overlap of wavefunctions increases the probability that a particle will be found in there, and out of phase overlap decreases the probability.

The physical result is the wide variety of bonding phenomena we see in nature


covalent bonds: local in phase overlap
anti-bonding orbitals: local out of phase overlap
conjugation: extended in phase overlap etc...