Consider $A$ as a matrix, that when normalized represents an finite-state time-homogeneous Markov chain $M$ with entries $0\leq p_{i,j}\leq 1$, where each row sums up to $1$. We can also assume that $M$ is irreducible and aperiodic, hence it has an unique steady-state distribution $\pi$ ($\pi=\pi P$).

Now consider a parameter $0<\delta\leq 1$ is added to the diagonal elements of the matrix $A$ Markov chain. So $a'_{i,j}=a_{i,j}+\delta$, from where we get $A'$. All rows of $A'$ are then normalized to sum up to one. From it we get new transition matrix $P'$. We again compute the steady-state distribution of this newly constructed chain $M'$ and call it $\pi'$.

If we observe what is happening with the elements of $\pi'$ if $\delta$ is increased, we can see that each value in the distribution monitonically increases or decreases.

Now my question. How could I prove that this is true for all $M$ with above properties (or does anyone have an counter example?)?

It also seems that $\pi'$ converges to some $\pi_C$ if I set $\delta$ to some arbitrarily large number (not compute the limit). Is there a direct formula to compute such $\pi_C$?

For example, lets take this matrix $A$ for out initial matrix:

$$A = \begin{bmatrix} 0.056&0.084&0.242&0.255\\ 0.071&0.056&0.249&0.210\\ 0.086&0.095&0.056&0.080\\ 0.115&0.102&0.101&0.056 \end{bmatrix}$$

Matrix $P$ for our Markov chain $M$, that we get by normalizing rows of $A$ (note that values are rounded to three decimal places):

$$P = \begin{bmatrix} 0.087&0.132&0.380&0.401\\ 0.121&0.095&0.426&0.358\\ 0.271&0.300&0.175&0.254\\ 0.309&0.272&0.270&0.149 \end{bmatrix}$$

Steady state distribution for chain $M$ is $$\pi= \begin{bmatrix} 0.211&0.213&0.298&0.278 \end{bmatrix}$$

$P'$ for $\delta=0.5$ (diagonal elements of $A$ are increased by $0.5$ and $A'$ is normalized by rows)

$$P' = \begin{bmatrix} 0.489&0.074&0.213&0.224\\ 0.065&0.512&0.229&0.193\\ 0.105&0.116&0.681&0.098\\ 0.132&0.117&0.116&0.636 \end{bmatrix}$$

Corresponding steady-state distribution for $M'$ (with $\delta=0.5$) is: $$\pi'=\begin{bmatrix}0.172 & 0.180 & 0.351 & 0.296\end{bmatrix}$$ And also I computed $\pi_C$ (with $\delta=100$) which is: $$\pi_C= \begin{bmatrix} 0.139&0.153&0.395&0.312 \end{bmatrix}$$

Below is a plot of steady-state distribution with increasing $\delta$. We can observe that each value in $\pi$ either monotonically increases or monotonically decreases.

And a little background... I am constructing a new multi-criteria decision making method, which is in some sense simmilar to PageRank algorithm. We are trying to evaluate given alternatives, based on user preference. I discovered that increasing $\delta$ in some sense increases the separation between the alternatives.

-edit comment- Sorry for the not so clear explanation and the error that I made. Hope we can now understand eachother.

-
How is $0<\delta\leq1$ compatible with $\delta\to\infty$? If $\delta\to\infty$, then $M'$ should converge to the identity matrix which does not have a unique stationary distribution. –  Eckhard Jan 3 '13 at 10:40
Yes it is not compatible. In my method I am using $0< \delta \leq 1$, but I noticed that if I increase $\delta$ to some larger values (e. g. 100), the values converge to above mentioned $\pi_C$ and not to $[0.25,0.25,0.25,0.25]$. Sorry for the confusion, I am a decision support grad students and mathematics was just a part of my undergrad program. –  Nejc Jan 3 '13 at 10:58
I'm a bit confused: if I'm not mistaken, the new matrix has to have the form $$P_\delta := \frac{1}{1+\delta}(P+\delta I)$$ however, $\pi_\delta P_\delta = \pi_\delta$ if and only if $\pi_\delta P = \pi_\delta$ which implies that the stationary distribution is independent of $\delta$. I guess, it's not the case though. –  Ilya Jan 3 '13 at 11:44
@Ilya: That observation is actually compatible with a few short calculations I did. –  Eckhard Jan 3 '13 at 11:47
@Eckhard: I think, that using the average return time characterization of the steady-state distribution, one obtains that $\pi$ is independent on the diagonal elements and depends only on ratios of the off-diagonal ones. Thus, my previous guess was incorrect and indeed, $\pi_\delta \equiv \pi$ for $\delta>0$ whereas $\pi_\infty$ in undefined as $P_\delta \to I$. –  Ilya Jan 3 '13 at 11:53

So far I can give only a partial answer. Let $A$ be the original matrix and let $$\alpha(i):=\sum_j a_{ij}$$ be the sum of its rows. For simplicity assume that all $a_{ij}$ are positive - just to be sure that the correspondent stochastic matrix is irreducible and aperiodic. Define $$K:=\mathrm{diag}(k(1),\dots,k(n))$$ with $k(i) = 1/\alpha(i)$ to be the normalizing matrix, so that $P:=KA$ is a stochastic matrix. Denote the correspondent steady-state distribution by $\pi$. Now let $\delta>0$ be the perturbation of $A$ and define $$A_\delta:=A+\delta I$$ where $I$ is the identity matrix, so $\alpha_\delta(i) = \alpha(i)+\delta$ and for the normalizing constants $$k_\delta(i) = \frac{k(i)}{1+\delta k(i)}.$$ The new stochastic matrix is then $P_\delta = K_\delta A_\delta$, so for the new steady-state distribution we have $$\pi_\delta(K_\delta A_\delta) = \pi_\delta.$$ If we denote $\nu_\delta:=\pi_\delta - \pi$, by opening all the brackets we obtain the following equation on $\nu_\delta$: $$\tag{1} \nu_\delta(P_\delta - I) = \pi(P-P_\delta).$$ Although equation $(1)$ does not have the unique solution as $P_\delta$ admits the eigenvalue $1$ and thus $(P_\delta - I)$ is not invertible, we also know that $\sum_i \nu_\delta(i) = 0$ since $\nu_\delta = \pi_\delta-\pi$. Using this last condition, we can obtain the desired solution of $(1)$, however I am not sure whether it is possible to have the solution in a neat form. As a special case, when $A$ is itself stochastic then $\nu_\delta\equiv 0$ and thus $\pi_\delta\equiv \pi$.
I suppose you meant $K_\delta$ instead of $N_\delta$... Anyway, i do not understand what you meant by $\sum \nu_\delta = 0$. Isn't $\nu_\delta$ fixed for some $\delta$. How can we sum them? –  Nejc Jan 4 '13 at 10:49
@Nejc: I updated the notation a bit. I meant, that the sum of elements of the vector $\nu_\delta$ is always zero as $\nu_\delta$ is a difference between two probability distributions. –  Ilya Jan 4 '13 at 11:00
OK, now I understand everything you have written. Can I now get the monotonicity per position in distribution (for instance for $i$-th element), if I rewrite the formula for index $i$? Instead of $\nu_\delta$ I would look at $\nu_\delta(i)$. –  Nejc Jan 4 '13 at 12:57