A surprising result about the product of Blaschke matrices I have verified analytically the conjecture described bellow up to $n=4$, but have had no success trying to prove it. Any help would be much appreciated. 

Setup
Let $\{\lambda_i\}_{i=1}^n$ be real numbers and $g_i:\mathbb R \to \mathbb R$ for all $i\in \mathbb N $. Consider the following recurrence
$$ g_{i+1}(x)=\frac{1+g_{i}(\lambda_{i})g_{i}(x)}{g_{i}(\lambda_{i})-g_{i}(x)}\frac{x-\lambda_{i}}{1-\lambda_{i}x}. $$
with 
$$ g_1(x) = a+bx .$$
Next, let
$$ \mathbf M_i(x) = \begin{bmatrix} 1 & 0 \\ 0 & \frac{x-\lambda_i}{1-\lambda_i x}\end{bmatrix}\begin{bmatrix} g_i(\lambda_i) & -1 \\ 1 & g_i(\lambda_i)\end{bmatrix}, $$
and
$$ S(x) = \begin{bmatrix} 0 & 1 \end{bmatrix}\begin{bmatrix} g_{n+1}(0) & -1 \\ 1 & g_{n+1}(0)\end{bmatrix}\left[\mathbf M_n(x)\mathbf M_{n-1}(x)\dots\mathbf M_{2}(x)\mathbf M_1(x)\right]\begin{bmatrix} 0 \\ 1\end{bmatrix}.$$

Conjecture

For all $n\ge 2$ and for all $i,j\in\{1,\dots,n\}$, $$\frac{S(\lambda_i)}{S(\lambda_j)} = \frac{\lambda_j\left(1+g_1(\lambda_j)g_1(1/\lambda_j)\right)}{\lambda_i\left(1+g_1(\lambda_i)g_1(1/\lambda_i)\right)}.$$


The reason I consider this surprising is that the analytical formulas for $S(x)$ become extremely complicated very quickly as $n$ increases, yet these ratios continue to satisfy this simple equation.

Background
The motivation for this question comes from the following. Suppose $\mathbf A (x)$ is a $2\times 2$ matrix such that $\det (\mathbf A (x))$ has $n+1$ roots (inside the unit circle): $\{\lambda_i\}_{i=1}^n$ and $0$. 
To solve some forecasting problems (in which $\mathbf A (x)$ represents the signal structure) it is useful to find a matrix $\mathbf B (x)$ such that $\det(\mathbf A (x)\mathbf B (x))$ has no roots (inside the unit circle).
A procedure to obtain the matrix $\mathbf B (x)$ is as follows. To remove the first root, $\lambda_1$, you multiply $\mathbf A (x)$ by an unitary matrix times the Blaschke matrix for this root, i.e. 
$$ \mathbf A^* (x) = \mathbf A (x) \mathbf W_1 \begin{bmatrix} 1 & 0 \\ 0 & \frac{x-\lambda_1}{1-\lambda_1 x}\end{bmatrix}$$
where $\mathbf W_1$ is a matrix whose columns are the left singular vectors of the of $\mathbf A(\lambda_1)$. To remove the second root, $\lambda_2$, you repeat this step using $\mathbf A^*(x)$ instead of $\mathbf A(x)$. Then, one can simply repeat this step for each of the other roots.
The matrices $\mathbf M_i(x)$ above are (up to a constant), the transpose of the matrices that are multiplied to $\mathbf A^*(x)$ in each step described above. One simplifying aspect of the problem I am working with is that the second row of $\mathbf A(x)$ is equal to $0$ when evaluated at any of the roots.
 A: Upon setting $\lambda_{n+1}=0$, we see that
$$S(x) = \frac{1}{x}\begin{bmatrix} 0 & 1 \end{bmatrix}\left[\mathbf M_{n+1}(x)\mathbf M_n(x)\mathbf M_{n-1}(x)\dots\mathbf M_{2}(x)\mathbf M_1(x)\right]\begin{bmatrix} 0 \\ 1\end{bmatrix}.$$ Thus (working with $n$ instead of $n+1$), we see that we have to control the lower right entry of $$\mathbf M(x)=\mathbf M_n(x)\mathbf M_{n-1}(x)\dots\mathbf M_{2}(x)\mathbf M_1(x).$$
By induction, one verifies that $g_i(\lambda_i)-g_i(x)$ is a scalar multiple (factor depending on $i$) of $\frac{x-\lambda_i}{1-x\lambda_{i-1}}$ for $i\ge 2$. Set $$\Gamma(x)=\prod_{i=1}^n\frac{x-\lambda_i}{1-x\lambda_i}\text{ and }g(x)=\prod_{i=1}^n(g_i(\lambda_i)-g_i(x)).$$ Thus
$$
g(x)=c\Gamma(x)(1-x\lambda_n)\tag{1}
$$
for a constant $c$.
In the following, we drop the argument $x$ if the meaning is clear. From the recursion for $g_i$ we obtain
$$\mathbf M_i\begin{bmatrix} 1 \\ g_i\end{bmatrix}=(g_i(\lambda_i)-g_i(x))\begin{bmatrix} 1 \\ g_{i+1}\end{bmatrix}.$$ Iterated application yields $$\mathbf M\begin{bmatrix} 1 \\ g_1\end{bmatrix}=g\begin{bmatrix} 1 \\ g_{n+1}\end{bmatrix}.\tag{2}$$
For a term $A$ in $x$ define $\tilde A$ by $\tilde A(x)=A(1/x)$. Then
$$\tilde{\mathbf M}_i^t\mathbf M_i=(g_i(\lambda_i)^2+1)\begin{bmatrix} 1 & 0 \\ 0 & 1\end{bmatrix},$$ hence $$\tilde{\mathbf M}^t\mathbf M=\gamma\begin{bmatrix} 1 & 0 \\ 0 & 1\end{bmatrix},$$ where $$\gamma=\prod_{i=1}^n(g_i(\lambda_i)^2+1).$$ Now we apply the tilde transformation to equation (2):
$$\tilde{\mathbf M}\begin{bmatrix} 1 \\ \tilde g_1\end{bmatrix}=\tilde g\begin{bmatrix} 1 \\ \tilde g_{n+1}\end{bmatrix}.$$ Transposing and multiplying from the right with $M$ yields
$$\gamma\begin{bmatrix} 1 & \tilde g_1 \end{bmatrix}=\begin{bmatrix} 1 & \tilde g_1 \end{bmatrix}\tilde{\mathbf M}^t\mathbf M=\tilde g\begin{bmatrix} 1 & \tilde g_{n+1} \end{bmatrix}\mathbf M.$$ Together with the previous relation (2) for $\mathbf M$ we obtain four linear equations for the four entries of $\mathbf M=\begin{bmatrix} m_1 & m_2\\ m_3 & m_4 \end{bmatrix}$ in terms of $\gamma,g_1,\tilde g_1,g_{n+1},\tilde g_{n+1},g,\tilde g$, namely
\begin{align*}
m_1+m_2g_1&=g\\
m_3+m_4g_1&=gg_{n+1}\\
m_1\tilde g+m_3\tilde g_{n+1}\tilde gt&=\gamma\\
m_2\tilde g+m_4\tilde g_{n+1}\tilde g&=\tilde g_1\gamma.
\end{align*}
The rank of the system is $3$. Use the first, second and fourth equation to express $m_1,m_2,m_3$ in terms of $m_4$. As expected, the third equation does not allow to solve for $m_4$, rather it collapses to
$$(g_{n+1}\tilde g_{n+1} + 1)g\tilde g=(g_1\tilde g_1+1)\gamma,\tag{3}$$
a relation we will need in a moment.
We can determine $m_4$ from
$$m_1m_4-m_2m_3=\det\mathbf M=\gamma\Gamma(x)$$ and obtain
$$
m_4(x)=\frac{(g\tilde g_1g_{n+1} + \tilde g\Gamma)\gamma}{(g_{n+1}\tilde g_{n+1} + 1)g\tilde g}.$$
Together with (3) this gives
$$m_4(x)=\frac{g\tilde g_1g_{n+1} + \tilde g\Gamma}{g_1\tilde g_1 + 1}.$$
For $1\le i\le n-1$ we have $g(\lambda_i)=0$, and $\lambda_i$ is not a pole of $\tilde g_1$, nor of $g_{n+1}$.
However, $\lambda_i$ is a root of $\Gamma$, and potentially a pole of $\tilde g$. So in order to evaluate $\tilde g\Gamma$ in $\lambda_i$ we apply the tilde operation to equation (1) and obtain $$\tilde g=c\tilde\Gamma\left(1-\frac{\lambda_n}{x}\right).$$
It follows from the definition of $\Gamma$ that $\tilde\Gamma\Gamma=1$, hence $$\tilde g\Gamma=c\frac{x-\lambda_n}{x}.$$
Now we can evaluate in $\lambda_i$ and $\lambda_j$, and obtain
$$
\frac{m_4(\lambda_i)}{m_4(\lambda_j)}
=\frac{\lambda_j\;(g_1(\lambda_j)g_1(1/\lambda_j) + 1)\;(\lambda_n - \lambda_i)}{\lambda_i\;(g_1(\lambda_i)g_1(1/\lambda_i) + 1)\;(\lambda_n - \lambda_j)}.$$
Thus, as explained in the beginning, with $\lambda_n=0$ we get the case $n-1$ of the question.
