Number of matrices with weakly increasing rows and columns

I'm curious as to how many matrices there are of size $m \times n$ with elements of the set $\{1, \ldots , k\}$ such that each row and column is weakly increasing?

The answer should be expressable as a determinant.

I'm thinking that this could be solved by counting non-intersecting lattice paths somehow and using Lindström–Gessel–Viennot lemma, however I'm unsure of how to construct the matrix.

Thanks!

-
Hint: you can view a weakly increasing sequence of numbers $0\leq j_1\leq j_2\leq\ldots j_n$ as corresponding to a lattice path that has level $j_i$ above its starting point when it first attains the line $i$ places to the right of its starting point, as described in this answer. – Marc van Leeuwen Mar 15 '12 at 10:48
Why should the answer be expressible as a determinant? – Michael Hardy Mar 15 '12 at 12:51
Can you do any special cases, maybe see some patterns? – Gerry Myerson Mar 23 '12 at 11:33
@MarcvanLeeuwen could you please explain your idea more detailed? It is not clear how to keep constraint of increasing rows and columns while calculating the total number... – Pavel Podlipensky Mar 8 '13 at 0:24

You're asking, essentially, about the number of plane partitions inside $m\times n\times(k-1)$ box.
The answer is given by MacMahon formula, $$\prod_{i=1}^m\prod_{j=1}^n\frac{i+j+k-2}{i+j-1}$$ (sanity check: for $k=1$ this is $1$, for $k=2$ or $n=1$ this is a binomial coefficient).
This formula, indeed, can be derived using LGV to count the number of non-intersecting paths, say, from $s_i=(i-1,-n-i+1)$ to $t_j=(m+j-1,-j+1)$ or (equivalently) to $t'_j=(m+j-1,0)$ (where $1\leqslant i,j \leqslant k-1$); the corresponding determinant $$\det\Bigl(P(s_i \rightarrow t_j')\Bigr)=\det\Biggl(\binom{m+n+j-1}{n+i-1}\Biggr)=\det\begin{pmatrix} \binom{m+n}n & \binom{m+n+1}n& \ldots & \binom{m+n+k-2}n \\ \binom{m+n}{n+1} & \binom{m+n+1}{n+1} & \ldots & \binom{m+n+k-2}{n+1} \\ \ldots & \!\ldots & & \!\ldots \\ \binom{m+n}{n+k-2} &\binom{m+n+1}{n+k-2}& \ldots & \binom{m+n+k-2}{n+k-2} \end{pmatrix}$$ can be computed using Vandermonde determinant.