shIf your matrices are "time-dependent" you need additional assumptions to bound the rate of convergence and it will usually be impossible to obtain something which depends only on the minimum of the second eigenvalues.
You may choose $a$ and $b$ at will in the example below, for example such that both matrices are doubly stochastic.
The product $\mathbf{A}\mathbf{B} \ = \ \mathbf{B}\mathbf{A}$ of the matrices
\begin{equation}
\mathbf{A} \ := \
\begin{pmatrix}
\frac{1}{3} + a^2 & \frac{1}{3} - a^2 & \frac{1}{3} \\
\frac{1}{3} - a^2 & \frac{1}{3} + a^2 & \frac{1}{3} \\
\frac{1}{3} & \frac{1}{3} & \frac{1}{3} \\
\
\end{pmatrix}
\end{equation}
and
\begin{equation}
\mathbf{B} \ := \
\begin{pmatrix}
\frac{1}{3} + \frac{1}{4} b^2 & \frac{1}{3} + \frac{1}{4} b^2 & \frac{1}{3} - \frac{1}{2} b^2 \\
\frac{1}{3} + \frac{1}{4} b^2 & \frac{1}{3} + \frac{1}{4} b^2 & \frac{1}{3} - \frac{1}{2} b^2 \\
\frac{1}{3} - \frac{1}{2} b^2 & \frac{1}{3} - \frac{1}{2} b^2 & \frac{1}{3} + b^2 \\
\
\end{pmatrix}
\end{equation}
is
\begin{equation}
\mathbf{E} \ := \
\begin{pmatrix}
\frac{1}{3} & \frac{1}{3} & \frac{1}{3} \\
\frac{1}{3} & \frac{1}{3} & \frac{1}{3} \\
\frac{1}{3} & \frac{1}{3} & \frac{1}{3} \\
\
\end{pmatrix} \, ,
\end{equation}
so whatever the values of $a$ and $b$ and independently of the initial conditions as long as all later matrices have row and column sums equal to 1, the process converges after just two iterations.
The matrices were constructed in the following way: Take an orthogonal basis of $\mathbb{R}^3$ containing the vector $(1,1,1)$, here $x = (1,1,1); y=(1,-1,0); z=(1/2,1/2,-1)$ and construct the matrices $A$ and $B$ as $k_x*x^tx + k_y*y^ty + k_z*z^tz$ by choosing the constants $k_x, k_y,k_z$ appropriately.
The eigenvalues of $A$ are $1; 2a^2; 0$, those of $B$ are $1; 1.5b^2; 0$, but while $A$ projects onto the orthogonal complement of $(1/2,1/2,-1)$, $B$ projects onto the orthogonal complement of $(1,-1,0)$, leaving only $k(1,1,1)^{(t)}$ afterwards, where $k$ is constant.
This is clearly a more general phenomenon, so without rather specific additional assumptions it will be very difficult/essentially impossible to progress, as the above result is independent of the second eigenvalues and there is no upper bound on the convergence rate in terms of them.
One might think about two things here:
1) There should be a lower bound on the convergence rate in terms of the supremum of the second eigenvalues if this is bounded away from 1 as in the time independent case because projecting will not increase the norm of the image.
2) If all the matrices $A(t)$ are invertible and their eigenvalue $\lambda_{small}$ of smallest modulus is bounded away from zero, there should be an estimate of the rate of convergence from above involving $\inf (|\lambda_{small}|)$.
Concerning the literature on Markov chains I can not really help you, but you may wish to consult
"Non-negative Matrices and Markov Chains" by E. Seneta, Springer, reprint 2006.
The example I just constructed from scratch.