INTRODUCTION
Consider an equation $$ \dot{x} = Ax + \varphi(x), $$ where $x \in \mathbb{R}^n$, $A \neq 0$ is constant matrix, $\varphi(x)$ is a germ of smooth vector field in $0$ s.t. $$ \left\{ \begin{align*} &\varphi(0) = 0\\ &\frac{d\varphi}{dx}(0) = 0 \end{align*}\right. . $$
There are plenty of methods to prove that this equation is smoothly equivalent to its linearisation $\dot{y} = Ay$ in some neighbourhood of zero. For example, there is method of homological equation described in Arnold's book. Or we can just substitute some unknown function and derive an ODE for this function.
The first question is, what other methods are known and which of them can be generalised to the case of parametric families of vector fields?
MY APPROACH
I've also conducted an experiment, and came up with the following method. We can introduce a parameter $\varepsilon$: $$ \dot{x}(t,\varepsilon) = Ax(t,\varepsilon) + \varepsilon\varphi(x(t,\varepsilon)). $$ When $\varepsilon = 0$ we get the linearised system, and when $\varepsilon = 1$ we have the original system. Now denote $$ y = \frac{\partial x(t,\varepsilon)}{\partial \varepsilon}. $$ Then we can differentiate the whole equation w.r.t. $\varepsilon$: $$ \frac{\partial}{\partial \varepsilon}\dot{x}(t,\varepsilon) = \frac{\partial}{\partial \varepsilon}Ax(t,\varepsilon) + \frac{\partial}{\partial \varepsilon}\varepsilon\varphi(x(t,\varepsilon)), $$ $$ \dot{y}(t,\varepsilon) = A y(t,\varepsilon) + \varphi(x(t,\varepsilon)) + \varepsilon\varphi'(x(t,\varepsilon))y(t,\varepsilon). $$ We know that the solution $x(t,\varepsilon)$ exists, is unique and smooth in some neighbourhood of zero $U$ for all $\varepsilon \in [0,1]$, so everything is alright and we can define smooth functions $F(t,\varepsilon) = \varphi(x(t,\varepsilon))$, $f(t,\varepsilon) = \varepsilon\varphi'(x(t,\varepsilon)) + A$ in $U \times [0,1]$. Then our equation $$ \dot{y}(t,\varepsilon) = f(t,\varepsilon)y(t,\varepsilon) + F(t,\varepsilon) $$ has a unique smooth solution which satisfies a condition $y(0,\varepsilon) = 0$ maybe in some smaller neighbourhood of zero $W$ for all $\varepsilon \in [0,1]$.
Finally, if we have two systems $$ \dot{x_0}(t) = A x_0(t) $$ and $$ \dot{x_1}(t) = A x_1(t) + \varphi(x_1(t)), $$ then the substitution $$ x_1(t) = x_0(t) + \int_{0}^{1} y(t, \xi) \, d\xi $$ turns the second system into the first one in some neighbourhood of zero $W$.
Is my approach correct? If not, how can it be fixed? Can we improve it? Can we extend it to the case of parametric families of vector fields?