# Asymptotic behaviour of a two-dimensional recurrence relation

This problem comes out of a research in models of firm growth.

The model is simple: A firm has two parameters which are its size (number of employees) and job vacancies. A firm of size $n$ will produce job vacancies at rate $n\mu$ and $k$ job vacancies will be filled at rate $k\beta$. On the other hand, a firm is likely to quit at rate $\delta$, and new firms come into the market at rate $\alpha$ (new firm starts with 0 employee and 1 vacancy, but one can also model it as 1 employee and 0 vacancy if necessary).

So in terms of math, the equation describing the relation is: $D_t s_{n,k,t}=-(n\mu+k\beta+\delta)s_{n,k,t}+(k+1)\beta s_{n-1,k+1,t}+n\mu s_{n,k-1,t}$, when (n,k) is not in the boundary, where $s_{n,k,t}$ is the probability of a firm having size $n$ and $k$ vacancies.

I am particularly interested in the stationary distribution which is described by:
$s_{0,1}=\frac{\alpha}{\delta+\beta}$,
$s_{n,k}=0$ for $n<0$ or $k<0$,
and $(n\mu+k\beta+\delta)s_{n,k}=(k+1)\beta s_{n-1,k+1}+n\mu s_{n,k-1}$.

Moreover, the generating function $M(x,y)=\sum_{n=0}^{\infty}{\sum_{k=0}^{\infty}{s_{n,k}x^ny^k}}$ satisfies the PDE $\delta M=\beta(x-y)M_y+\mu x(y-1)M_x+\alpha y$.

What I would like to know is the asymptotic behavior of $s_n=\sum_{k=0}^{\infty}{s_{n,k}}$, or more explicitly the tail $\sum_{n=m}^{\infty}{s_n}$ when m is very large. The empirical Data on firm size suggests that it might appear as a power law (~ $n^{-\zeta}$, where $\zeta$ is a bit larger than 1), also it may not have a finite variance ($M_{xx}(1,1)=\infty$). The one-dimensional problem, in which we only consider the size (that is, vacancies are filled immediately), is completely characterized by Yule process and also suggests the power law.

Could anyone give me some insight on how to approach this problem? I've been considering this problem for a month, but since I am not familiar with highly technical tools in analysis or probability I get stuck completely.

Thank you!

-
The obvious approach would be to compute separately each generating series s_n(x)=sum(s_{n,k}x^k,k=0..infty). It seems your relations yield something like beta.x.s'_n(x)+(n.mu(1-x)+delta)s_n(x)=beta.s'_{n-1}(x). Starting from s_0(x)=s_{0,1}x one gets recursively s_n(x) as a function of s_{n-1}(x). I wonder if you tried this or if the result is too complicated to be estimated efficiently or what. – Did Aug 10 '11 at 23:47
Thank you for your reply. But indeed the solution is very complicated. It doesn't even has a closed form. – epsilon Aug 11 '11 at 4:35
To help people help you, you could show how you know this (that the solution is complicated and so on). – Did Aug 11 '11 at 8:27
I didn't remember what I got since s_1 was already complicated. But for example, using the equation you gave, mathematica gives s1[x] -> E^((n u x)/b - (d Log[x])/b - (n u Log[x])/b) C[1] - E^((n u x)/b - (d Log[x])/b - (n u Log[x])/b) s0 x^((d + n u)/ b) ((n u x)/b)^(-((d + n u)/b)) Gamma[(d + n u)/b, (n u x)/b], where b=beta, u=mu, d=delta, s0=s_{0,1}. I don't think it is possible to solve for s_2. – epsilon Aug 12 '11 at 0:24
What I was mainly thinking about was to find out the behaviour of M_xx(x,1) when x->1 (because M_xx doesn't have a finite variance), and then use some versions of tauberian theorem to get the tail. It actually worked for a modified model where the PDE for M was soluble. (I used results from this paper: epubs.siam.org/sidma/resource/1/sjdmec/v3/i2/p216_s1). But that model was kind of "artificial". – epsilon Aug 12 '11 at 0:28