$\renewcommand{\a}{\mathbf a}$The conclusion is indeed very intuitive. However, the proof of it is rather tedious, along the lines of the proof of the Picard–Lindelöf theorem.

Accordingly, assume the following standard conditions:

For each $j\in[n]:=\{1,\dots,n\}$, let $f_j$ be $L$-Lipshitz for some positive real $L$. Take any real $x$. Suppose that for some positive real $M$ and $b$ and for all $j\in[n]$ we have $|f_j|\le M$ on $[x-b,x+b]$.
Take any $T\in(0,\min(1/L,b/M))$.

By time rescaling, without loss of generality (wlog)
\begin{equation*}
\sum_{j\in[n]}a_j=1. \tag{1}\label{1}
\end{equation*}

Then the initial value problem
\begin{equation*}
z'(t)=f_\a(z(t)),\quad z(0)=x
\end{equation*}
has a unique solution on the interval $[0,T]$, where $f_\a:=\sum_{j\in[n]}a_j f_j$.
This problem can be rewritten in integral form, as follows:

\begin{equation*}
z(t)=x+\int_0^t f_\a(z(s))\,ds \tag{2}\label{2}
\end{equation*}
for $t\in[0,T]$.

For a natural $K$, suppose that
\begin{equation*}
0=t_0<\cdots<t_K=T.
\end{equation*}
For $k\in[K]$ and $j\in\{0\}\cup[n]$, let
\begin{equation*}
h_k:=t_k-t_{k-1},\quad t_{k,j}:=t_{k-1}+h_k\sum_{r=1}^j a_r,
\end{equation*}
so that $t_{k,0}=t_{k-1}$ and $t_{k,n}=t_n$, in view of \eqref{1}.

The role of $f_{\theta(t)}(x)$ will be played by
\begin{equation*}
g_t(x):=g_{\a,t}(x):=g_\a(t,x):=\sum_{k\in[K]}\sum_{j\in[n]} f_j(x)\,1(t_{k,j-1}\le t<t_{k,j});
\end{equation*}
so, $g_t=f_j$ on every interval of the form $[t_{k,j-1},t_{k,j})$.

Consider the initial value problem
\begin{equation*}
y'(t)=g_t(y(t)),\quad y(0)=x,
\end{equation*}
where $y'(t)$ is the right derivative of $y$ at $t$. This problem can be rewritten in integral form, as follows: for each $(k,j)\in[K]\times[n]$ and all $t\in[t_{k,j-1},t_{k,j}]$
\begin{equation*}
y(t)=y(t_{k,j-1})+\int_{t_{k,j-1}}^t f_j(y(s))\,ds. \tag{3}\label{3}
\end{equation*}

**The intuition:** So, within each (say small) time interval $[t_{k-1},t_k)$ of length $h_k$, we let the rate of change of $y$ be $f_j(y)$ for the fraction $a_j$ of the length $h_k$ of the time interval $[t_{k-1},t_k)$. Since $y$ will vary little within the small time interval $[t_{k-1},t_k)$, the effect of this switching between the $f_j$'s will be almost the same as the effect of using the average rate $f_\a=\sum_{j\in[n]}a_jf_j$.

Now the formalities:
Reasoning as in the proof of the Picard–Lindelöf theorem, we see that problem \eqref{3}
has a unique solution on the interval $[0,T]$.

Moreover, for all $(k,j)\in[K]\times[n]$ and all $t\in[t_{k,j-1},t_{k,j}]$
\begin{equation*}
|y(t)-y(t_{k,j-1})|\le Mh_k a_j\le Mh_k, \tag{3.5}\label{3.5}
\end{equation*}
so that
\begin{equation*}
|f_j(y(t))-f_j(y(t_{k,j-1}))|\le LMh_k,
\end{equation*}
so that
\begin{equation*}
|f_j(y(t_{k,j-1}))-f_j(y(t_{k-1}))|\le nLMh_k \tag{4}\label{4}
\end{equation*}
and
\begin{equation*}
|y(t)-y(t_{k,j-1})-f_j(y(t_{k,j-1}))(t-t_{k,j-1})| \\
=\Big|\int_{t_{k,j-1}}^t (f_j(y(s))-f_j(y(t_{k,j-1})))\,ds\Big|
\le LMh_k^2,
\end{equation*}
so that
\begin{equation*}
\Big|y(t_k)-y(t_{k-1})-\sum_{j\in[n]}f_j(y(t_{k,j-1}))h_k a_j\Big| \le nLMh_k^2,
\end{equation*}
so that, in view of \eqref{4},
\begin{equation*}
\Big|y(t_k)-y(t_{k-1})-\sum_{j\in[n]}f_j(y(t_{k-1}))h_k a_j\Big| \le nLMh_k^2+nLMh_k^2,
\end{equation*}
that is,
\begin{equation*}
\Big|y(t_k)-y(t_{k-1})-f_\a(y(t_{k-1}))h_k\Big| \le 2nLMh_k^2
\end{equation*}
for all $k\in[K]$.

Similarly, but a bit more simply, we get
\begin{equation*}
\Big|z(t_k)-z(t_{k-1})-f_\a(z(t_{k-1}))h_k\Big| \le 2LMh_k^2.
\end{equation*}
for all $k\in[K]$.

So, for all $k\in[K]$, letting $h_k:=h:=T/K$,

$D_k:=|y(t_k)-z(t_k)|$, and $E_k:=D_k+c$ for $c:=4nMh$, we get
\begin{equation}
|D_k-D_{k-1}|\le hLD_{k-1}+4nLMh^2,
\end{equation}
\begin{equation}
D_k\le (1+hL)D_{k-1}+4nLMh^2,
\end{equation}
\begin{equation}
E_k\le (1+hL)E_{k-1},
\end{equation}
\begin{equation}
E_k\le (1+hL)^K E_0=(1+hL)^K c\le e^{LT}c=e^{LT}4nMh,
\end{equation}
\begin{equation}
D_k\le E_k\le Ch,
\end{equation}
where $C:=e^{LT}4nM$.

Also, by \eqref{3.5}, for all $k\in[K]$ and all $t\in[t_{k-1},t_k]$ we have
$|y(t)-y(t_{k-1})|\le Mh$;
similarly, $|z(t)-z(t_{k-1})|\le Mh$.

Thus, $|y-z|\le(C+2M)h=(C+2M)T/K$ on $[0,T]$. Choosing now $K$ to be large enough, we make $|y-z|$ however small uniformly on $[0,T]$. $\quad\Box$