It is a standard fact that the (linear) Hales--Jewett theorem tensorizes to yield its multidimensional version. The same observation applies to the polynomial version.
I will describe the idea for the two-dimensional case: assume that $n=2k$ is an even integer, $A$ is a finite alphabet and $c:A^{[n]^2}\to [r]$ is an $r$-coloring. Split the discrete interval $[n]$ into two successive subintervals $I_0:=[k]$ and $I_1:=\{k+1,\dots,2k\}$; this induces a partition of the index set $[n]^2$ into four squares $I_0\times I_0$, $I_0\times I_1$, $I_1\times I_0$ and $I_1\times I_1$ that are all isomorphic to $[k]^2$. Therefore, the $r$-coloring $c$ also induces an $r$-coloring of $(A^4)^{[k]^2}$, which I shall denote by $C$. Then observe that a monochromatic, with respect to $C$, polynomial combinatorial line of $(A^4)^{[k]^2}$ unfolds to a two-dimensional polynomial combinatorial subspace of $A^{[n]^2}$ that is monochromatic with respect to $c$. Note that this subspace is quite special: its wildcard sets are successive and have the same size.
This argument yields the estimate
$$ MPHJ(d,m,k,r) \leq PHJ(d,k^{m^d},r), $$
where
- $PHJ(d,k,r)$ denotes the polynomial Hales--Jewett number, that is, the least positive integer $N$ such that for every integer $n\geq N$ and every alphabet $A$ with $k$ letters, an arbitrary $r$-coloring $c:A^{[n]^d}\to [r]$ has a monochromatic polynomial combinatorial line, and
- $MPHJ(d,m,k,r)$ denotes the corresponding multidimensional polynomial Hales--Jewett number, that is, the least positive integer $N$ such that for every integer $n\geq N$ and every alphabet $A$ with $k$ letters, an arbitrary $r$-coloring $c:A^{[n]^d}\to [r]$ has a monochromatic $m$-dimensional polynomial combinatorial subspace.
The paper of Shelah you cited has an ingenious proof that the shows that, for every fixed dimension $d$, the function $(k,r) \mapsto PHJ(d,k,r)$ is upper bounded by a primitive recursive function that belongs to the class $\mathcal{E}_{7+3d}$ of Grzegorczyk’s hierarchy. (Here, the constants are most likely not optimal, but the linear dependence on the dimension $d$ is necessary.) The same paper also has a proof that claims that the polynomial Hales--Jewett numbers are upper bounded by a primitive recursive function (namely, we have a uniform control of their growth as the dimension $d$ increases) but that proof is problematic (see Fact 4.7 in the first preprint).
Summing up, for every fixed dimension $d$, the function $(m,k,r)\mapsto MPHJ(d,m,k,r)$ is upper bounded by a primitive recursive function that belongs to the class $\mathcal{E}_{7+3d}$ of Grzegorczyk’s hierarchy. From this, one can easily extract bounds for the multidimensional polynomial Van der Waerden theorem.