It is true that if connection $\omega$ in a vector bundle is $\mathfrak{g}$-valued ($\mathfrak{g}$ being Lie Algebra of the structure Lie Group $G$) in a patch $U$, then it will be $\mathfrak{g}$-valued in all other patches due to the transformation law:
\begin{align} \omega_V = c^{-1}_{UV} \omega_U c_{UV} + c^{-1}_{UV} d c_{UV}. \end{align}
However, I'm trying to understand if it must take values in Lie algebra in order to be a connection?
On
I think I got it sorted out, at least motivationally.
Apart from the usual definition of a matrix of 1-forms with peculiar transformation properties, the connection could also be defined through identification of fibers at different points of base manifold. Unsurprisingly, it would make a lot of sense if these fibers are identified through multiplication by an element of Lie Group $G$. Now for the infinitesimal form of such identification we would naturally get Lie Algebra $\mathfrak{g}$.
Note that it is by no means required that the connection takes values in $\mathfrak{g}$ (unless we deal with Principal Bundle). The argument above might be considered an inspiration at best, not a mathematical treatment of the question.
On
Let the structure group be $\operatorname{O}(m),\operatorname{U}(m), \operatorname{SO}(m),\operatorname{SU}(m)$ (or $\operatorname{O}(k,l),...$), and the parallel transport preserve metric of Riemannian/Hermitian vector bundle $\sqcup_{x\in M} V_x$. $\newcommand{\op}[1]{\operatorname{#1}} \newcommand{\map}[5]{\begin{array}{cccc} #1 & #2 & \longrightarrow & #3 \\ & #4 & \longmapsto & #5 \end{array}} \newcommand{\K}{\mathbb{K}} \newcommand{\R}{\mathbb{R}} \newcommand{udb}[2]{\underbrace{#1}_{#2}} \newcommand{\e}[1]{\frac{\partial}{\partial #1}} \newcommand{\al}[1]{\begin{align*} #1 \end{align*}} \newcommand{\br}[1]{\left(#1\right)}$
Let $\gamma:[0,a]\to M$ be a smooth curve with initial value $\gamma(0)=p,\ \dot\gamma(0)=X$; Denote the parallel transport $P(\gamma)_t^0: V_{\gamma(t)}\to V_{\gamma(0)}$.
The definition of $G$ structure bundle chart ($G=\operatorname{SO}(m),...$) require that, the bundle chart vector frame $\phi^{-1}(x,e_j)$ is orthonormal. On Riemannian/Hermitian bundle, local $C^\infty$ orthonormal frame can be constructd by Gram-Schmidt procedure (So if we have this metric, we have this structure bundle chart), see https://math.stackexchange.com/a/3002983/634573. So, let $E_j|_x\in \sqcup_{x\in W} V_x$ be a local $C^\infty$ orthonormal frame around $p$. Compute the covariant derivative
$$ \nabla_X E_j=\lim_{t\to 0}\frac{P(\gamma)_t^0 E_j|_{\gamma(t)}-E_j|_{\gamma(0)}}{t}=A_X(p)_{kj} E_k|_p $$
By metric preserving, $P(\gamma)_t^0 E_j|_{\gamma(t)}$ is another orthonormal frame of $V_p$, so there is curve $t\xrightarrow{C^\infty} U(t)\in \operatorname{SO}(m),\cdots$ with $U(0)=I$, such that $$ P(\gamma)_t^0 E_j|_{\gamma(t)}= U(t)_{kj} E_k|_{p}=(E_1|_p,...,E_1|_p)U(t)\left(\array{ \vdots \\ \underbrace{1}_{\text{$j$-th}} \\ \vdots }\right) $$ Use notation $(w_1,...,w_n)=(v_1,..,v_n)\cdot \left(\array{ a_{11} &\cdots & a_{1n} \\ \vdots & \ddots &\vdots \\ a_{n1} & \cdots & a_{nn}}\right)=(a_{k1} v_k ,...,a_{kn} v_k )$. Now we have $$ \nabla_X E_j=(E_1|_p,...,E_1|_p)\cdot \underbrace{\lim_{t\to 0}\frac{U(t)-I}{t}}_{=\dot U(0)}\cdot\left(\array{ \vdots \\ \underbrace{1}_{\text{$j$-th}} \\ \vdots }\right)=\dot U(0)_{kj} E_k|_p $$ So we get Lie algebra valued $A_X(p)=\dot U(0)\in \mathsf{so}(m),\mathsf{su}(m),...$
Let $(F_1,...,F_m)$ be another local $C^\infty$ frame, and $B_X(x)$ the corresponding connection matrix; Let $(F_1,...,F_m)|_x=(E_1,...,E_m)|_x\cdot T(x)$. Then the the transformation rule:
$$ B=T^{-1}\cdot dT + T^{-1}\cdot A\cdot T $$
The proof is: write
$$(\nabla_X E^1,\ldots,\nabla_X E^m) (p)= (E^1,...,E^m)|_p\cdot A_X(p)$$
Let $(F^1,...,F^m)|_x=(E^1,...,E^m)|_x\cdot T(x)$. It can prove that
$$ \al{(\nabla_X F^1,..., \nabla_X F^m)&=(E^1,...,E^m)\cdot \partial_X T + (\nabla_X E_1,...,\nabla_XE^m)\cdot T \\ &= (E^1,...,E^m) \br{ \partial_X T + A_X\cdot T}} $$
On the other hand,
$$\al{(\nabla_X F^1, ..., \nabla_X F^m)&=(F^1, ..., F^m)\cdot \tilde A_X \\ &=(E^1,..., E^m)\cdot T\cdot \tilde A_X }$$
So $\tilde A_X=T^{-1}\cdot\partial_X T +T^{-1} A_X \cdot T$. Or proof via computation: $$ \al{\nabla_{X} \underbrace{F_j}_{\sum_{\alpha=1}^m T_{\alpha j} E_\alpha}&=\sum_{\alpha=1}^m\partial_X T_{\alpha j}E_{\alpha}+ \sum_{\alpha=1}^mT_{\alpha j} \underbrace{\nabla_X E_{\alpha}}_{\sum_{\beta=1}^m(A_X)_{\beta\alpha} E_\beta} \\ &= \sum_{\beta=1}^m\left(\partial_X T_{\beta j}+ \sum_{\alpha=1}^m T_{\alpha j} (A_X)_{\beta \alpha}\right)\underbrace{E_\beta}_{\sum_{k=1}^m(T^{-1})_{k\beta} F_k}} $$ Then
$$ \al{\underbrace{(B_X)_{kj}}_{\text{$k$ row, $j$ colume of $B_X$}}&=\underbrace{\sum_{\beta=1}^m (\partial_X T_{\beta j}) (T^{-1})_{k\beta}}_{\text{$k$ row, $j$ colume of $T^{-1}\cdot(\partial_X T)$}}+\underbrace{\sum_{\alpha,\beta=1}^m T_{\alpha j} (A_X)_{\beta\alpha} (T^{-1})_{k\beta}}_{\text{$k$ row, $j$ colume of $T^{-1}\cdot A_X\cdot T$}} \\ B_X&=T^{-1}\cdot(\partial_X T)+T^{-1}\cdot A_X\cdot T} $$
If $F_i$ is also an orthonormal frame, or equivalent $T(x)\in \operatorname{SO}(m),...$ then we also have $B_X(p)\in \mathsf{so}(m),...$ On the other hand, $$ T^{-1}\cdot(\partial_X T) +T^{-1}\cdot A_X\cdot T=\frac{d}{dt}|_{t=0} \left(T(p)^{-1}\cdot T\big(\gamma(t)\big) + T(p)^{-1}\cdot U(t)\cdot T(p) \right) $$
(1) If [structure group is matrix group $G\subset \operatorname{GL}(m,\mathbb K)$] require that [we select a family of basis $\{(E_1,...,E_m)|_x\cdot T\}_{x\in W \\ T\in G}$ of every fiber $V_p$], where $E_j|_x=\phi^{-1}(x,e_j)$ is bundle chart vector basis.
Another bundle chart vector basis $F_j|_x=\psi^{-1}(x,e_j)$ is also in this family, $$ \require{AMScd} \begin{CD} V_x @>{\text{identity map, $(F_1,...,F_m)|_x=(E_1,...,E_m)|_x\cdot g_{\psi\phi}(x)$}}>> V_x \\ @V{(W_\psi,\psi)}VV@VV{(W_\phi,\phi)}V\\ \mathbb K^n @>>{g_{\psi\phi}(x)}> \mathbb K^n \end{CD} $$
(2) Once again, assume that the parallel transport (or say it connection) maps between these basis, that is, $\{P(\gamma)_t^0 E_j|_{\gamma(t)}\}_{j=1}^m\in \{(E_1,...,E_m)|_{\gamma(0)}\cdot T\}_{T\in G}$. Or in other words, as OP's self answer say, we require that, the bundle chart representation of parallel transport $P(\gamma)_t^0:V_{\gamma(t)}\to V_{\gamma(0)}$ is an element of $G\subset \operatorname{GL}(m,\mathbb K)$.
By same argument, there is $t\xrightarrow{C^\infty} U(t)\in G$ with $U(0)=I$ such that $$ P(\gamma)_t^0 E_j|_{\gamma(t)}=(E_1|_p,...,E_1|_p)U(t)\left(\array{ \vdots \\ \underbrace{1}_{\text{$j$-th}} \\ \vdots }\right) $$ Then similarly we can get $\nabla_X E_j=\dot U(0)_{kj} E_k|_p$ and Lie algebra valued $A_X(p)=\dot U(0)\in \mathsf{g}$.
The converse is also ture. That is, if the connection matrix is Lie algebra valued, then connection preserve $G$ structure.
Let $\forall_X$ connection matrix $A_X$ be Lie algebra valued. Let $P(\gamma)_0^t E_\mu|_p= Z_\mu{}^k(t) E_k|_{\gamma(t)}$. Compute parallel transport equation $$ 0=\nabla_{\dot\gamma(t)} (Z_\mu{}^k(t)E|_{\gamma(t)})=\bigg(\dot Z_\mu{}^k(t)+\udb{\dot\gamma^i A_{i}}{=A_{\dot\gamma(t)}} \big(\gamma(t )\big)^k{}_{j} Z_\mu{}^j(t)\bigg) E_k|_{\gamma(t) } $$ We obtain the system of ODE equations $$ \al{\dot Z_\mu{}^1&= -A_{\dot\gamma(t)}\big(\gamma(t)\big)^1{}_j Z_\mu{}^j\\ &\cdots\\ \dot Z_\mu{}^m&= -A_{\dot\gamma(t)}\big(\gamma(t)\big)^m{}_j Z_\mu{}^j} \\ \text{initial value:}\ Z_\mu{}^k(0)= \delta_\mu{}^k $$ Written in vector form this is $$ \dot Z_\mu(t)=-A_{\dot \gamma(t)}\big(\gamma(t))\big) Z_\mu(t),\quad Z_\mu(0)=\br{\array{ \vdots \\ \underbrace{1}_{\text{$\mu$-th}} \\ \vdots }} $$ Write $Z_1,\ldots,Z_m$ in matrix form, that is $$ \dot Z(t)= -\big(A_{\dot\gamma(t)}\circ \gamma\big)(t)\cdot Z(t),\quad Z(0)=I $$ Because for each $t$, the matrix $A_{\dot \gamma(t)}\in \mathsf g$ is Lie algebra-valued. Then how to prove that $Z(t)\in G$? Let us abbreviate $-(A_{\dot{\gamma}(t)}\circ\gamma)(t)$ to $A(t)$.
We find that this ODE is a time-varying vector field on $\op{Matrix}(m,\K)\cong \K^{m^2}$ $$ \map{X:}{\R\times \K^{m^2}}{\K^{m^2}}{(t,B)}{A(t)B} $$ of the integral curve.
If we consider $M=\op{Matrix}(m,\K)$ as a real manifold (notation $M$ dose not mean the original base manifold), then the value domain should be $\sqcup_{B\in M}T_B M$. The matrix Lie group $G$ is a real sub-manifold of $M$. However, it’s troublesome to convert complex matrix to real matrix, via some linear isomorphism induced by $\op{GL}(m,\mathbb C),\op{GL}(2m,\R)$, and $(x_1+iy_1,...,x_m+iy_m)\longleftrightarrow \left\{\array{(x_1,y_1,...,x_m,y_m) \\ (x_1,...,x_m,y_1,...,y_m)}\right.$. When we write $\lim_{t\to 0}\frac{U(t)-I}{t}$ to get the tangent vector, we want this compatible in both real and complex case. An observation is that, $\frac{U(t)-I}{t}$ just use linearity, but not complex multiplication.
If $B\in G$, then $A(t)\cdot B\in T_{B} G$. This can be verified via curve: since $A(t)\in \mathsf g$ is Lie algebra valued, there exists a curve $s\mapsto \alpha(s)\in G$ such that $\alpha(0)=I,\ \dot \alpha(0)=A(t)$.
Then for the curve $s\mapsto \alpha(s) B\in G$, we have $$ \alpha(0)B=B\\\frac{d}{ds}\Big|_{s=0} \alpha(s)B=\udb{\dot\alpha (0)}{A(t)}B $$ Thus we get $A(t)B \in T_B G$.
Now the problem translates to: if $$ \map{X:}{\R\times M}{\sqcup_{p\in M}T_p M}{(t,p)}{X(t,p)} $$ is a time-varying tangent vector field on the $N$ dimension real manifold, and assume that, when $X$ is restricted to the sub-manifold $S$, it becomes a vector field on $S$. If the initial value of the integral curve $\phi(t)$ is also on $\phi(0)\in S$, does the whole integral curve stay inside the sub-manifold $S$ and not go out of it?
Recall that, the ODE on the (real) manifold is indeed well-defined and does not depend on the choice of the coordinate chart. The solution of the vector field $X$ — i.e. the integral curve –- is the map $\map{\phi:}{\R}{M}{t}{\phi(t)}$ satisfying $\frac{d}{dt} \phi(t)=X(t,\phi(t))$. When write down this equation, it already does not depend on the coordinates chart. But we still need to look at, how it look like in the coordinate chart $(U,x^1,... ,x^m)$. Since $\phi(t)=\phi^k(t)\e{x^k}$ and $X(t,p)=X^k(t,p)\e{x^k}$. The ODE should be $$ \frac{d}{dt}\phi^k(t)=X^k(t,\phi^1(t),... ,\phi^m(t)) $$ Both sides are contravariant in $\R^N$. Thus, for the case of sub-manifolds, we can choose a coordinate chart such that $x|_{S\cap U}(p)=(x^1,...,x^l ,0,... ,0)(p)$. Together with the assumption that $X|_{S}(t,p)\in T_p S$, we also have $X|_{S\cap U}(t,p)=(X^1,... ,X^l,0,... ,0)(t,p)$.
So the ODE in the coordinate chart can be solved first for the part of $\K^l$, and then construct a solution on $\K^m$ — let the $(l+1)\sim m$ component be zero. Then, by the uniqueness of the ODE, the integral curve will be restricted to $\K^l$ in the coordinate system. That is, the integral curve is restricted to the sub-manifold $S$.
Applying to our case $M=\op{Matrix}(m,\K)\cong \K^{m^2}$ and $S=G$, the ODE $\dot Z(t)=A(t)Z(t),\ Z(0)=I$, we get, the integral curve $Z(t)\in G$.
Yes, your connection needs to be valued in the Lie algebra $\mathfrak g$ of a Lie group $G$ which is the structure group of your bundle ($ G = GL(n)$ for a rank-$n$ vector bundle, $G = O(n)$ for a rank-$n$ real vector bundle with metric, etc). Here is the intuition behind that:
You choose some patch $U$ and use a local trivialization to write your vector bundle as $V \mid_U \cong U \times \mathbb R^n$. Via this local trivialization, you have identified all the fibers $V_x, x \in U$ with $\mathbb R^n$. So now you can write sections of $V$ over $U$ as $\mathbb R^n$-valued functions and multivariable calculus tells you how to take partial derivatives along given directions.
But you should not just naively start taking partial derivatives of the individual components of a section. There are many possible trivializations of $V$ over $U$, and if your friend chooses a different trivialization and you both take partial derivatives in the naive way, you will disagree with each other. This is why we need covariant derivatives.
So the connection enters: We recognize that in addition to the naive recipe for taking partial derivatives, we need to add on a "rotation" of the section as we take derivatives along a given direction. So say that at the point $x\in U$ we want to take a covariant derivative along the $X \in T_xM$ direction (and say $s$ is our section over $U$): We take "naive" multivariable calculus derivatives of the components of $s$ in the $X$ direction and then we ask our Lie algebra-valued connection 1-form $\omega$ how we should "rotate" $s$ in order to make the derivative covariant.
Now $\omega$ eats the tangent direction $X$ and gives us a Lie algebra element, $\omega(X) \in \mathfrak g$. Elements of $\mathfrak g$ are (intuitively speaking) infinitesimal "rotations," where "rotation" means "element of $G$." So we act on $s$ by the infinitesimal rotation $\omega(X)$ and we add the result to the "naive derivative" of $s$ along the $X$-direction.
(Above "let $\omega(X)$ act on $s$" meant implicitly that we have fixed some representation of $G$ on the fibers of $V$ so that an expression like $\omega(X) s$ makes sense.)
So the concise answer: $\omega$ has to be Lie-algebra valued because its job is to specify infinitesimal rotations that we add to the usual partial derivatives in order to make the derivative covariant.
I haven't attempted to convince you that this procedure will resolve the original difficulty of consistently defining derivatives of sections despite multiple ways of trivializing $V$ over $U$, but that's the essence of the transformation law you wrote down. If you'd like I could edit this answer with more detail along those lines.
I haven't read your reference by Frankel but I recommend that you additionally look at DuPont's Fibre Bundles and Chern-Weil Theory.