Why need non-degeneracy in this theorem? In Sakurai's book Modern Quantum Mechanics, p276 chapter 4, near equation(4.4.59), there is a theorem saying:

Therorem. Suppose the Hamiltonian is invariant under time reversal and
the energy eigenket $|n>$ is nondegenerate; then the corresponding
energy eigenfunction is real or, more generally, a real function times
a phase factor independent of x).

I don't know why the requirement of "nondegenerate" is a must. If $H$ is invariant under time reversal operator $\Theta$, then $H\psi_n=E_n\psi_n$ leads to $H\Theta\psi_n = E_n\Theta\psi_n$. Since $\Theta\psi_n=\psi_n^*$, we can always choose $\phi_n=\psi_n+\psi_n^*$ and $\varphi_n=i(\psi_n-\psi_n^*)$ to be the new real eigenfunctions, even if they are degenerate (if no degenearcy, $\phi_n$ and $\varphi_n$ are linear dependent).
The above discussion is within the spinless Hamiltonian, as implied by the context of the theorem in the book. Of course, feel free to generalize to spinful cases.
 A: The condition on the theorem is necessary.
If you have two independent real eigenfunctions $\psi_1$ and $\psi_2$, then the eigenfunction $\psi_1 + \mathrm{e}^{\mathrm{i}\alpha}\psi_2$ for $\alpha\in\mathbb{R}$ suitably chosen is neither real nor real times an overall phase factor.
Note that the theorem is not meant to be a claim about existence of real eigenfunctions (you've shown how to construct real eigenfunctions in the degenerate case), but a claim about all eigenfunctions being either real or real times an overall phase factor.
A: A Hamiltonian is invariant under time reversals if $[H , T ] = 0$. Remember that $T$ is an anti-unitary operator. Now, consider the state $|n\rangle$ which satisfies
$$
H|n\rangle = E_n |n\rangle
$$
Now, consider the state $T|n\rangle$. Since $H$ commutes with $T$ and $E_n$ is real, it can easily be shown that the state $T|n\rangle$ also has energy $E_n$ (check this yourself). Now, if the state $|n\rangle$ is non-degenerate, then we must have
$$
T |n\rangle = e^{i\theta_n}|n\rangle , \qquad \theta_n \in [0,2\pi) \tag{1}
$$
$\theta_n$ is some phase. Now, lets look at the wave-function. Act on both sides of the equation above with $\langle x|$. Since time reversals does not act on the operator $x$, we have $T| x \rangle  = | x \rangle$. Thus,
$$
\langle x | T |n\rangle = e^{i\theta_n}\langle x |n\rangle \implies e^{i\theta_n} \psi_n(x) = \langle n | T^\dagger|x\rangle = \langle n | T^{-1} |x\rangle = \langle n |x\rangle  = \langle x |n \rangle ^* = \psi_n(x)^*
$$
Here, we have used the fact that for anti-linear operators, the adjoint is defined as $\langle n | T^\dagger x\rangle = \langle x | T n\rangle$ (instead of $\langle n | L^\dagger x\rangle = \langle L n | x \rangle$ which is the definition for linear operators). We have also used that the time-reversal operator is unitary so $T^\dagger = T^{-1}$.
The final result is
$$
\psi_n(x)^* = e^{i\theta_n} \psi_n(x). \tag{2}
$$
Here, $\theta_n$ does not itself depend on $x$. We now write
$$
\psi_n(x) = e^{-i\theta_n/2} f_n(x)
$$
Substituting this into equation (2), we find
$$
[ e^{-i\theta_n/2} f_n(x) ] ^* = e^{i\theta_n} [ e^{-i\theta_n/2} f_n(x) ] \implies e^{i\theta_n/2} f_n(x)^* = e^{i\theta_n}  e^{-i\theta_n/2} f_n(x)  \implies f_n(x)^* = f_n(x) .
$$
Thus, we find that we can write
$$
\psi_n(x) = e^{-i\theta_n/2} f_n(x) \quad \text{where} \quad f_n(x)^* = f_n(x) .
$$
QED.

Let me also discuss what happens if we relax the non-degeneracy requirement. In this case, there is no relation between the states $|n\rangle$ and $T|n\rangle$. Now recall that for all states in the Hilbert space, $T^2$ acts as $+1$ or $-1$ (depending on the spin of the particle). We consider two cases:

*

*If $T^2 |n\rangle = |n\rangle$, we define two states $|n\rangle_\pm = |n\rangle \pm T |n\rangle$. These states satisfy $T |n\rangle_\pm = \pm |n\rangle_\pm$ which is precisely equation (1) with phase $1$ or $-1$. The rest of the proof follows exactly showing that the wave-function $\psi_{\pm,n}(x)$ is a phase times a real function.


*If $T^2 |n\rangle = - |n\rangle$, it is not possible to define any linear combination such that (1) is true. The theorem completely fails in this case.

Let me make a comment about the phase $\theta_n$. This phase freedom is actually not physical. This is because, we can rescale the state by an unphysical phase $|n\rangle \to e^{i\alpha_n} |n\rangle$. In this case, (1) is modified to
$$
T (  e^{i\alpha_n} |n\rangle )  = e^{i\theta_n} (  e^{i\alpha_n} |n\rangle ) \implies  e^{-i\alpha_n} T   |n\rangle  = e^{i( \theta_n + \alpha_n)} |n\rangle \implies T   |n\rangle  = e^{i( \theta_n + 2\alpha_n)} |n\rangle
$$
We used the fact that $T$ is an anti-linear operator so $T(c |n\rangle) = c^* T |n\rangle$. We can now see that if we choose $\alpha_n = - \frac{1}{2} \theta_n$, we can remove the phase in equation (1). In terms of the new phase rescaled state, the equation is simply $T|n\rangle = |n\rangle$ which then implies that $\psi_n(x)$ is real.
