In manifold theory, in what sense is the derivative a first-order approximation? As I move on from the calculus definition of the derivative to the differential geometric definition in terms of tangent spaces, I am wondering how to recover the notion that the derivative of a function represents the first order approximation.
Clearly the definitions of the derivative traffic in first order concepts: tangent vectors can be thought of equivalence classes of curves agreeing up to first order; one characterization of the tangent space at $x \in M$ is the dual space to $I/I^2$, where $I$ is the ideal of $C^\infty$ functions that vanish at $x$, and so on. 
These definitions seem to be well-suited to a characterization of $Df$ as approximating $f$ up to first order at $x \in M$, but I am having trouble putting my finger on it and expressing it precisely. 
I know that when we're in $\mathbb{R}^n$ the "approximation" property of the derivative works very well because we can identify the tangent space at a given point with $\mathbb{R}^n$ itself. But what about in the more general setting? In what sense do we have "approximation up to first order" there? What is it exactly that's agreeing with $f$ up to $o(\|h\|)$?
 A: One way to phrase this is in terms of Taylor's theorem for maps between Euclidean spaces, after some "chart-nonsense". Maybe you're looking for something not involving coordinates, in which case this answer may not be satisfying.
Let's say that $f: M^m \to N^n$ is a map between manifolds. Given a fixed "basepoint" $x \in M$, choose coordinate charts $(\phi, U)$, where $x \in U \subset M$, and $(\psi, V)$, where $f(x) \in V \subset N$. By this I mean that $\phi: U \to \phi(U) \subset \mathbb{R}^m$ and $\psi: V \to \psi(V) \subset \mathbb{R}^n$ are homeomorphisms to open sets in Euclidean spaces. For notation's sake, let $\tilde{f} = \psi \circ f \circ \phi^{-1}: \phi(U) \to \psi(V)$ be the map between Euclidean spaces corresponding to $f$ in our coordinates, and let $\tilde{x} = \phi(x) \in \phi(U) \subset{\mathbb{R}^m}$. 
Taylor's theorem applied to $\tilde{f}$ says that
$$ \tilde{f}(\tilde{x} + h) = \tilde{f}(\tilde{x}) + D\tilde{f}_{\tilde{x}} \cdot h  + o(||h||).$$
Here $D\tilde{f}_{\tilde{x}}$ is the $n \times m$ matrix of partial derivatives of $\tilde{f}$ at $\tilde{x}$. We view $D\tilde{f}_{\tilde{x}}$ as a linear map $\mathbb{R}^m \to \mathbb{R}^n$.
$D\tilde{f}_{\tilde{x}}$ is precisely the matrix for $Df_x: T_x M \to T_{f(x)} N$ in our chosen coordinates, which give us bases for the tangent spaces $T_x M$ and  $T_{f(x)} N$.
Edit: If our manifolds are endowed with Riemannian metrics, then the exponential map gives us a canonical way to identify a neighborhood of the origin in the tangent space $T_x M$ with a neighborhood of $x \in M$ (and similarly for $N$). Then Euclidean coordinates on $T_x M$ give us geodesic normal coordinates on $M$, and the above discussion of course applies in those coordinates. Radial distances (i.e., distances to $x$) are preserved under the exponential map. Thus, with the notation from above, $||h||$ is precisely the distance from $x$ to $\phi(\tilde{x} + h) = \text{exp}_x(h)$, and you can think of "$o$" as being measured in terms of distances on $N$. In this setting one should think that we are approximating $f(y)$ (for $y$ near $x$) by $f(x)$ plus a first-order correction, which is given precisely by applying the derivative $Df_x$ to a tangent vector $h$ that measures the deviation of $y$ from $x$ (modulo translating via $\text{exp}_x$ and $\text{exp}_{f(x)}$ between the manifolds and the tangent spaces, where "plus" makes sense).
