Eigenkets of degenerate perturbation theory Suppose the original Hamiltonian is $H$ and we perturb it by a small potential $V$.  The basis kets of the original hamiltonian $H$ contains some degeneracy.
Since there's some degeneracy, we take the set of degenerate eigenkets of the original hamiltonian $H$ and diagonalise it with respect to the small potential $V$.  The result is a set of eigenkets of the full Hamiltonian $H+V$.
Since these are exactly the eigenkets of the full Hamiltonian $H+V$, does it mean that there are no higher order perturbation terms?
 A: No, you diagonalize in the truncated degeneracy subspace D, so your eigenkets are definitely not eigenkets of the full hamiltonian! 
Most texts are needlessly formal and lose you in a  thicket of abstract formalism, and you never thought to illustrate them by a daft minimal example. Here is one.
Take 
$$
H_0+\lambda V=  \begin{pmatrix}
 1  &  0 & 0 \\
  0 &  1 & 0 \\
0 & 0 & 2
\end{pmatrix} +\lambda  \begin{pmatrix}
 0  &  1 & 2 \\
  1 &  0 & 3\\
2 & 3 & 7
\end{pmatrix} .$$
Consider 
$$|1\rangle=\begin{pmatrix}
1 \\
0 \\
0
\end{pmatrix}  ,\qquad |2\rangle=\begin{pmatrix}
0 \\
1 \\
0
\end{pmatrix}  , \qquad |3\rangle=\begin{pmatrix}
0 \\
0 \\
1
\end{pmatrix}  ,$$
$$ |a\rangle=  \frac{1}{\sqrt{2}} \begin{pmatrix}
1 \\
-1 \\
0
\end{pmatrix}  , \qquad |b\rangle= \frac{1}{\sqrt{2}}\begin{pmatrix}
1 \\
1 \\
0
\end{pmatrix}  ,$$
so the orthogonal states $|a\rangle, |b\rangle$ diagonalize the 1-2 subspace perturbation, but obviously not the entire hamiltonian.
You must check by direct calculation, and also pursuant to your text or WP, that, to lowest order in λ, that is ignoring its square, 
the  energy eigenvalues and  eigenstates are
$$
E_A= 1-\lambda; \qquad E_B= 1+\lambda; \qquad E_C=2+\lambda 7 ;\\
|A\rangle= |a\rangle + \frac{\lambda}{\sqrt{2}}  |3\rangle= \frac{1}{\sqrt{2}}\begin{pmatrix}
1 \\
- 1\\ \lambda 
\end{pmatrix} ~; \qquad  |B\rangle=  |b\rangle -\lambda  \frac{5} {\sqrt{2}}|3\rangle= \frac{1}{\sqrt{2}}\begin{pmatrix}
1 \\ 1\\
-5\lambda   
\end{pmatrix} ~; \\
|C\rangle=  |3\rangle   -\lambda  \frac{1} {\sqrt{2}}|a\rangle +\lambda  \frac{5} {\sqrt{2}}|b\rangle= \begin{pmatrix}
2\lambda \\ 3\lambda\\
1
\end{pmatrix} ~.
$$  
Recall normalizations, orthogonalities, etc... all work to $O(\lambda^2)$ only. 
In conclusion, $|a\rangle, |b\rangle$ are not exact hamiltonian eigenstates--and you have all perturbation orders survive, in principle. They merely avoided zeros in the denominator of the resolvent, having removed the degeneracy. That is, they ensured that $|A\rangle$ will have no $|b\rangle$ correction to first order, and, likewise, $|B\rangle$ will have no $|a\rangle$ correction, because this basis has diagonalized V, so it fails to connect $|a\rangle$ with  $|b\rangle$ and they do not impact each other, so far. But higher orders are another matter.


*

*Fine print (wonkish) footnote on the above "so far": avoid until 3rd reading. 


In actuality, an extra arbitrary piece $\propto \lambda |b\rangle$ in $|A\rangle$ and  $\propto \lambda |a\rangle$ in $|B\rangle$ has no impact (linkage) on the eigenvalue equation at the first order in λ discussed above, as you may readily check: the O(λ) potential has superselected them out to this order. However, as you may check in Courant-Hilbert  Methods of Mathematical Physics, v1 p 248 , and discussed in this question , you do get two innocuous (impact-less) pieces to this order, namely 
$$
|A\rangle=\hbox{ above } - \lambda \frac{5}{4} |b\rangle,\qquad  
|B\rangle=\hbox{ above } + \lambda \frac{5}{4} |a\rangle,\\
\langle b|A\rangle=\frac{\lambda}{2}\langle b|V|3\rangle \frac{1}{E^{(0)}_{a,b}-E^{(0)}_{3}}\langle 3|V|a\rangle=-5/4, \qquad \langle a|B\rangle=5/4.
$$
These pieces, bilinear in the potential!,  structurally  "fell off" the 2nd order calculation, where λ2 was discounted by a λ in a subspace denominator $ E^{(1)}_{b}-E^{(1)}_{a}=2\lambda$. The particular 5/4 coefficients are fixed by consistency with the next order energies. But you have already seen this subtlety above, in zeroth order:   $|a\rangle, |b\rangle$ of  $O(\lambda^0)$  actually  were specified by  $O(\lambda)$ potential considerations. In any case, for reassurance, compute 
$$
E_A= 1-\lambda -\lambda^2/2 + 7\lambda^3/8+...; \\ E_B= 1+\lambda-25\lambda^2/2+625\lambda^3/8+...; \\ E_C=2+7\lambda  +13\lambda^2-79\lambda^3+... 
$$
