Disclaimer: I like Tyberius' answer, but I'd like to go a little further.
$%
\newcommand{\ll}{\left\langle}\newcommand{\rr}{\right\rangle}
\newcommand{\lb}{\left|}\newcommand{\rb}{\right|}
\newcommand{\op}[1]{\mathbf{#1}}$
The beauty of Koopmanns' theorem is its simplicity. Quoting from the gold book:
Koopmans' theorem
Directly relates experimental ionization potentials with energy levels of molecular orbitals. The theorem states that the ionization potential required to remove an electron from the orbital $\Psi_i$ is given by the negative value of the energy of the orbital, $−\varepsilon_i$, as calculated within the Hartree–Fock approximation. The theorem is not applied to localized molecular orbitals, which are not eigenfunctions of the effective hamiltonian.
While Koopmans' text (in German: Physica 1934, 1, 104-113.) is a bit hard to comprehend, due to the ancient language I never learned, it basically boils down to two important facts.
There is a unitary transformation of the Lagrange multipliers which diagonalises them. These form a characteristic set of elementary wave functions, matrix, which in turn have a characteristic eigenvalue each. If all eigenvalues are different, then all elementary wave functions are specified through the Hamiltonian. (If there are equal eigenvalues, a freely chosen unitary transformation can be performed.)
In other words: There is one set of canonical orbitals that diagonalises the Fock matrix; each canonical orbital has an eigenvalue. For non-degenerate systems those orbitals are specified through the Hamiltonian.
This choice of functions has physical meaning, as the eigenvalue (disregarding a small correction) can be equated to the ionisation energies of the corresponding electron.
The second statement implies that there is no relaxation of the orbitals when removing an electron. This is oftentimes referred to as "frozen MO" approximation. Koopmanns himself appreciates that there should be a contraction of orbitals when removing an electron, however, he doesn't go into detail.
Nowadays we explain Koopmans' theorem in some simpler terms, where we remove the $k$th electron:
\begin{align}
&& E_N &= \sum_{i=1}^N H_{ii}
+ \frac12\sum_{i=1}^N\sum_{j=1}^N (J_{ij} - K_{ij})
+ V_{\mathrm{nuc}} \tag1\\
&& E_{N-1}^k &=\sum_{i=1}^{N-1} H_{ii}
+ \frac12\sum_{i=1}^{N-1}\sum_{j=1}^{N-1} (J_{ij} - K_{ij})
+ V_{\mathrm{nuc}} \tag2\\\hline
&& E_N - E_{N-1}^k &= H_{kk}
+ \frac12\sum_{i=1}^N (J_{ik}-K_{ik})
+ \frac12\sum_{i=1}^N (J_{kj}-K_{kj}) \tag{$1-2$}\\
\therefore&&
E_N - E_{N-1}^k &= H_{kk} + \sum_{i=1}^N (J_{ki}-K_{ki})\\
\therefore&&
E_N - E_{N-1}^k &= \varepsilon_k\\[2ex]
\text{with}&&
H_{ii} &= \langle \phi_i(\mathbf{x}_1)|
\mathbf{H}^\mathrm{c} | \phi_i(\mathbf{x}_1)\rangle\\
&&
J_{ij} &= \langle \phi_i(\mathbf{x}_1) \phi_j(\mathbf{x}_2) | r_{12}^{-1}
| \phi_i(\mathbf{x}_1) \phi_j(\mathbf{x}_2) \rangle\\
&&
K_{ij} &= \langle \phi_i(\mathbf{x}_1) \phi_j(\mathbf{x}_2) | r_{12}^{-1}
| \phi_j(\mathbf{x}_1) \phi_i(\mathbf{x}_2) \rangle\\
\end{align}
As you can see, Koopmans' theorem is not limited to the HOMO, but can be used for any occupied MO.
There is a very important point to be always be kept in mind when working with any Hartree-Fock based approaches:
Even though the equation
$$\op{F}_i\phi_i = \varepsilon_i\phi_i
\tag{3}\label{fock-pseudo}$$
suggests an eigenvalue problem, it is not. Remember the definition of the Fock operator and the operators contained
\begin{align}
&& \op{F}_i &= \op{H}^\mathrm{c} + \sum_j (\op{J}_j - \op{K}_j),\\
\text{with}&&
\op{J}_j\lb \phi_i\rr &=
\ll \phi_j(\op{x}_1) \rb r_{12}^{-1}
\lb \phi_j(\op{x}_1) \rr \lb \phi_i(\op{x}_2) \rr,\\
\text{and}&&
\op{K}_j\lb \phi_i\rr &=
\ll \phi_j(\op{x}_1) \rb r_{12}^{-1}
\lb \phi_i(\op{x}_1) \rr \lb \phi_j(\op{x}_2) \rr.
\end{align}
As you can see, the "one-electron" Fock operator depends on the solution of all "one-electron" Fock operators (ref. Szabó-Ostlund p. 115). The Hamilton operator is not the sum of all Fock operators, and the total HF energy is not the sum of all orbital energies. As a result of that, the canonical orbitals are actually unique solutions.
The Fock operator is associated with the whole wave function, and while unitary transformations will keep the wave function and its energy equivalent, everything else will fall apart.
Another important consideration is the fact that the Fock operator is only well defined for occupied MO. You sometimes find the statement
$$E_{N+1}^{l>N} - E_N = \varepsilon_l,$$
which can not be applied in the same way; it is very basis set dependent.
From all that above it is obvious, that Koopmans' theorem only works for systems where a single determinant approximation is reasonable. It also explains, why it only works for HF; although there are generalisations for DFT.
Is Koopmans' theorem still valid after a unitary transformation?
No. As stated above, the Fock operator is associated with the $N$-electron wave function. After a unitary transformation and then removal of an electron, the energy of the $N-1$-electron wave function is not conserved.
As an illustrative example:
\begin{align}
&& \sum_{i=1}^{N} H_{ii} &= \sum_{i=1}^{N} H_{ii}'\\[2ex]
\text{with}&&
H'_{ii} &= \langle \phi'_i(\mathbf{x}_1)|
\mathbf{H}^\mathrm{c} | \phi'_i(\mathbf{x}_1)\rangle\\
\text{and}&&
|\phi_i\rangle &\color{red}{\neq} |\phi'_i\rangle\\[2ex]
&& \sum_{i=1,i\neq k}^{N-1} H_{ii} + H_{kk}&=
\sum_{i=1,i\neq k}^{N-1} H_{ii}' + H_{kk}'\\
\text{for } H_{kk} = H_{kk}':
&& \sum_{i=1,i\neq k}^{N-1} H_{ii} &=
\sum_{i=1,i\neq k}^{N-1} H_{ii}'&&
\implies|\phi_i\rangle \color{red}{=} |\phi'_i\rangle\\
\end{align}
I have indicated the contradiction in red. Therefore it follows:
$$|\phi_i\rangle \neq |\phi'_i\rangle
\implies H_{kk} \neq H_{kk}'
\implies
\sum_{i=1,i\neq k}^{N-1} H_{ii} \neq \sum_{i=1,i\neq k}^{N-1} H_{ii}'
$$
You can follow that trough for the other terms and will see that the expectation value of energy of the $N-1$-electron wave function needs to be different from the expectation value of energy of the $N-1$-electron wave function after unitary transformation.
