dg.differential geometry - Smooth Dependence of ODEs on Initial Conditions

Dear all,
The following is a theorem known to many, and is essential in elementary differential geometry. However, I have never seen its proof in Spivak or various other differential geometry books.

Let $t_0$ be real, and $x_0 in mathbb{R}^n$ and $a,b>0$. Let $f:[t_0-a,t_0 + a] times overline{B(x_0,b)}rightarrow mathbb{R}^n$ be $C^k$ for $kge 1$.

Then $f$ is Lipschitz continuous, with which it is easy, using the contraction mapping theorem of complete metric spaces, to prove that the ODE:
$dfrac{d}{dt}alpha(t,x)=f(t,alpha(t,x)),quad alpha(t_0,x)=x$

has a continuous solution in an open neighbourhood of $(t_0,x_0)$. In other words, the ODE
$x'(t)=f(t,x(t));x(t_0)=x_0$ has a family of solutions which depends continuously on the initial condition $x_0$.

The theorem that I'd like to prove is that, in fact, if $f$ is $C^k$, then $alpha$ is $C^k$, for any $kge 1$.

I'd like an "elementary" proof that needs no calculus on Banach spaces or any terribly hard theory such as that, but hopefully something elementary, such as the contraction mapping theorem. I currently have an attempt of a proof that looks at perturbations of linear ODEs, but it is incorrect (I think). The proof can be found on page 6 of http://people.maths.ox.ac.uk/hitchin/hitchinnotes/Differentiable_manifolds/Appendix.pdf. I believe that there is a typo in the claim:

"Apply the previous lemma and we get

$mathrm{sup}_{left| tright|leq epsilon}left|lambda(t,x)y-{alpha(t,x+y)+alpha(x)}right|=o(left|yright|).$"

but more importantly, what it should be replaced by is incorrect. What is needed is that $||A-B_y||=o(||y||)$ but I do not see why this is.

Thank you for your time and effort.