Let $A\in \mathbb{S}^{N\times N}$ be a symmetric, real and stable matrix, i.e., $\rho(A)<1$, where $\rho(A)$ stands for the spectral radius of $A$. Then, $$\sum\limits\_{i=0}^{\infty} A^{2i}=\left( I-A^2 \right)^{-1}.\tag{$\star$}\label{star}$$


Now, let $A$ be asymmetric. It is not clear whether the power series $$F(A):=\sum\_{i=0}^{\infty} A^i (A^i)^{\top}$$ admits a *compact closed-form* expression as in equation \eqref{star} (when it converges). If the matrix $A$ is normal — i.e., it commutes with its transpose — and stable, then we have that $F(A)=\left(I-AA^{\top}\right)^{-1}$. In general, that is not the case.



> 
> **Framework.** Let $A$ be a random matrix generated as follows. Let $\widetilde{A}\in \mathcal{G}(N,p)$ be the adjacency matrix associated with the realization of a random graph over $N$ nodes with probability $p$ of arrow drawing, i.e., the entries $\widetilde{A}\_{ij}$ of $\widetilde{A}$ are i.i.d. $\mathsf{Bernoulli}(p)$ for all $i\neq j$ — this is often referred to as a binomial random graph. Set $\widetilde{A}\_{ii}=0$ for all $i$. Let $A$ be obtained from $\widetilde{A}$ by: **i) [off-diagonal]** normalizing the entries of $\widetilde{A}$ as $A\_{ij}:=\alpha\_1 \widetilde{A}\_{ij}/d\_{\max}$, for all $i\neq j$, where $d\_{\max}$ is the maximum *in-flow* degree of the random graph and $0<\alpha\_1<1$ is some constant; **ii) [diagonal]** setting the diagonal terms as $A\_{ii}:=\alpha-\sum\_{k\neq i} \widetilde{A}\_{ik}$, where $0<\alpha\_1\leq\alpha<1$. In other words, the rows of $A$ sum to $\alpha$ and its support is given by the realization of the binomial random graph on $N$ nodes.
> 
> 
> 


**Question.** Is it true that $$\mathbb{P}\left(\lVert\sum\_{i=0}^{\infty}A^i \left(A^i\right)^{\top} - \left(I-AA^{\top}\right)^{-1}\rVert\_\text{max}>\varepsilon\right)\overset{N\rightarrow \infty}\longrightarrow 0?$$


Is there any known result in this direction? Any reference would be appreciated. The specifics of the random nature of $A$ are not important other than that of $A$ being stable and randomly generated (so that it is more natural to consider the limit $N\rightarrow\infty$).


Some numerical experiments are pointing in the direction of a *concentration* of $\sum\_{i=0}^{\infty} A^i \left(A^i\right)^{\top}$ about $\left(I-AA^{\top}\right)^{-1}$ as the dimension $N$ grows large (I may attach them later, if necessary) under the referred framework.


**Context.** $F(A)$ is the covariance matrix associated with the time-series reflecting the state-evolution of certain linear stochastic dynamical systems. Having a closed form expression for the covariance in terms of $A$ is helpful for model identification purposes.



(Asymmetric) matrix power series in closed form: $\sum_{i=0}^{\infty} A^i \left(A^i\right)^{\top}={?}$
Let me answer the first question:


**Q:** Does $F=\sum\_{i=0}^{\infty} A^i (A^i)^{\top}$ for non-symmetric $A$ have a closed-form solution analogous to the solution $\sum\_{i=0}^{\infty} A^{2i}=\left( I-A^2 \right)^{-1}$ for a symmetric $A$?


**A:** Yes, in terms of the [vectorization operation](https://en.wikipedia.org/wiki/Vectorization_(mathematics)):
$$\operatorname{vec}(F)=(I-{A} \otimes A)^{-1}\cdot\operatorname{vec}(I).\tag{$\ast$}\label{ast}$$
This does require the inversion of an $N^2\times N^2$ matrix, rather than the inversion of an $N\times N$ matrix as in the case of a symmetric $A$.


The vectorization formula follows because $F$ solves the [Lyapunov equation](https://en.wikipedia.org/wiki/Lyapunov_equation)
$$AFA^\top=F-I.$$


I do not have an answer to the second question, whether \eqref{ast} converges for large $N$ to $(I-AA^\top)^{-1}$.


