Claim: Let $\underset{(n \times n)}{\textbf{B}}$ be a positive definite matrix with eigenvalues $\lambda_1 \geq \lambda_2 \geq ... \geq \lambda_n \geq 0$ and associated normalized eigenvectors $\textbf{e}_1, \textbf{e}_2, ... \textbf{e}_n$. Then, $$ \underset{\textbf{x} \neq \textbf{0}}{max} \hspace{6pt} \frac{\textbf{x}^T\textbf{B}\textbf{x}}{\textbf{x}^T\textbf{x}}=\lambda_1 \hspace{10pt} (\text{when }\textbf{x}=\textbf{e}_1).$$
Here is my proof (so far??) (which I gained most of my insight from Applied Multivariate Statistical Analysis (6th Edition) by Johnson and Wichern).
Proof: Let $\underset{(n \times n)}{\textbf{P}}$ be the orthogonal matrix whose columns are the eigenvectors $\textbf{e}_1,\textbf{e}_2,...\textbf{e}_n$ and $\Lambda$ be the diagonal matrix with eigenvalues $\lambda_1,\lambda_2,...,\lambda_n$. Let $\textbf{B}^{1/2}=\textbf{P}\Lambda^{1/2}\textbf{P}^T$ and $\underset{(n \times 1)}{\textbf{y}}=\underset{(n \times n)}{\textbf{P}^T}\underset{(n \times 1)}{\textbf{x}}$ where $\textbf{x}\neq \textbf{0}$, and thus $\textbf{y} \neq \textbf{0}$. Now notice,
\begin{eqnarray*}
\frac{\textbf{x}^T\textbf{B}\textbf{x}}{\textbf{x}^T\textbf{x}} &=& \frac{\textbf{x}^T\textbf{B}^{1/2}\textbf{B}^{1/2}\textbf{x}}{\textbf{x}^T\textbf{P}\textbf{P}^T\textbf{x}} \hspace{50pt} (\text{by definition of a square root matrix})\\
&=& \frac{\textbf{x}^T\textbf{P}\Lambda^{1/2}\textbf{P}^T\textbf{P}\Lambda^{1/2}\textbf{P}^T\textbf{x}}{\textbf{x}^T\textbf{P}\textbf{P}^T\textbf{x}} \hspace{20pt} (\text{by definition of }\textbf{B}^{1/2})\\
&=& \frac{\textbf{x}^T\textbf{P}\Lambda^{1/2}\textbf{P}^T\textbf{P}\Lambda^{1/2}\textbf{P}^T\textbf{x}}{\textbf{y}^T\textbf{y}} \hspace{20pt} (\text{by definition of }\textbf{y})\\
\end{eqnarray*}
Now let $\textbf{x}=\textbf{e}_1$. Thus,
by the definition of spectral decomposition, that $\textbf{e}_k^T\textbf{e}_j=1$ if $k=j$ and $\textbf{e}_k^T\textbf{e}_j=0$ if $k \neq j$. Thus, now,
$$\frac{\textbf{y}^T \Lambda^{1/2}\textbf{I}\Lambda^{1/2}\textbf{y}}{\textbf{y}^T\textbf{y}}= \frac{\textbf{y}^T\Lambda\textbf{y}}{\textbf{y}^T\textbf{y}}=\frac{\lambda_1}{1}=\lambda_1.$$
However, my question is, how do I use the max part in the claim?
Using the given eigenbasis, any vector can be written as
$$\mathbf x=x_1\mathbf e_1+x_2\mathbf e_2+\cdots+x_n\mathbf e_n,$$
and by definition of eigenvalue/vector, $\mathbf B\mathbf e_i=\lambda_i\mathbf e_i$, so
$$\mathbf B\mathbf x=x_1\lambda_1\mathbf e_1+x_2\lambda_2\mathbf e_2+\cdots+x_n\lambda_n\mathbf e_n.$$
The eigenvectors are orthonormal; $\mathbf e_1^T\mathbf e_1=1,\;\mathbf e_1^T\mathbf e_2=0$, etc; so the quadratic form is
$$\mathbf x^T\mathbf B\mathbf x=(x_1\mathbf e_1+x_2\mathbf e_2+\cdots+x_n\mathbf e_n)^T(x_1\lambda_1\mathbf e_1+x_2\lambda_2\mathbf e_2+\cdots+x_n\lambda_n\mathbf e_n)$$
$$=x_1^2\lambda_1+x_2^2\lambda_2+\cdots+x_n^2\lambda_n,$$
and because $\lambda_1$ is the largest eigenvalue, and $x_i^2\geq0$, this is
$$\leq x_1^2\lambda_1+x_2^2\lambda_1+\cdots+x_n^2\lambda_1$$
$$=(x_1^2+x_2^2+\cdots+x_n^2)\lambda_1$$
$$=(\mathbf x^T\mathbf x)\lambda_1.$$
Therefore, for any $\mathbf x\neq\mathbf0$,
$$\frac{\mathbf x^T\mathbf B\mathbf x}{\mathbf x^T\mathbf x}\leq\lambda_1.$$
This alone doesn't mean that $\lambda_1$ is the maximum; it only means that $\lambda_1$ is an upper bound, which might never be reached. But we can see that it is the maximum; it is reached when $x_1\neq x_2=\cdots=x_n=0$, or $\mathbf x\propto\mathbf e_1$.