Let $S = S_{\mathbf{Q}} = M_{13}(\Gamma_1(N),\mathbf{Q})$, and $S_{\mathbf{C}} = S \otimes \mathbf{C}$ denote the corresponding space of modular forms over $\mathbf{C}$.
Let $V \subset S \times S$ be the subspace cut out by pairs of forms $(A,B)$ satisfying the following equation:
$$A \cdot E_{12} = B \cdot \Delta $$
As equations in the Fourier coefficients of $A$ and $B$ these are linear equations with coefficients in $\mathbf{Q}$. Since, by the $q$-expansion principle, a modular form can be recovered from some finite number of Fourier coefficients, $V$ is determined by the null space of some finite matrix with coefficients in $\mathbf{Q}$. Since a linear system over $\mathbf{Q}$ has the same rank over $\mathbf{C}$, it follows that $V_{\mathbf{C}} = V \otimes \mathbf{C}$, where $V_{\mathbf{C}}$ is the set of solutions in $S_{\mathbf{C}} \times S_{\mathbf{C}}$ of the same equations.
On the other hand, there is an isomorphism:
$$V_{\mathbf{C}} \rightarrow M_{1}(\Gamma_1(N),\mathbf{C})$$
given by
$$(A,B) \mapsto \frac{A}{\Delta} = \frac{B}{E_{12}}$$
The point is that $E_{12}$ and $\Delta$ do not have any common zeros, so the image of this map clearly consists of holomorphic forms. Hence the map is well defined. However, if $F$ has weight one, then $(A,B) = (F \cdot \Delta, F \cdot E_{12})$ maps to $F$, so the map is surjective. It is clearly injective, so it is an isomorphism.
It follows that the image of $V$ under this map gives a rational basis for $M_1(\Gamma_1(N))$. Since $V$ is then preserved by Hecke operators (as is obvious on $q$-expansions), the result follows.
Best Answer
More simply, the eigenvalues are real because the Hecke operators are Hermitian for the Peterson scalar product (a fact which can be checked by a straightforward computation). See for example the introduction by Serre on modular forms in Cours d'arithmétique.