Take a linear ordinary differential equation of the form : $$ \sum_{k=0}^n p_{n-k}(z)(z (z-1))^k \partial_k f = 0 $$ Where $p_i$ is a polynomial fraction of degree $i$, without zeros at $0$ or $1$, $p_{0} = 1$

Using the Frobenius method, we can find a basis of solution around $0$ (a regular singularity) of the form $f_i(z) = z^{\gamma_i} a_i(z)\ \ i \in \{1 \cdots l \}$, with $a_i$ a power serie at $0$. The $\gamma_i$ are solutions of a polynomial equation, and I'm assuming that no pair of them is separated by an integer (so the solutions are assured to be independent).

The same can be done around $1$ ($g_i(z)$) and around $\infty$ ($h_i(z)$), with the same assumption. So in the end I have three basis of solutions, of size $l$.

I would like to know **if there's a way to find the coefficients of the matrices linking those three basis** (without having to solve the equation).
Explicitely, I'm looking for $A$ and $B$ in :
$$
f_i = A_i^j g_j, \quad f_i = B_i^j h_j
$$

To give an example, if the equation is of order $2$, we obtain an hypergeometric equation, whose solutions are hypergeometric functions. The coefficients of the matrices $A$ and $B$ are then ratios of gamma functions. Generalized hypergeometric equations (a subclass of the equations above) behave in the same way (ratio of gamma function).

Any insight is appreciated, thanks in advance,