Stapledon decomposition conjecture for open order polytopes

Let $O_{\Pi}$ be an order polytope of dimension $d$. Let $p_{\Pi}(z)$ be the numerator polynomial of its open Ehrhart series, defined by

$$\sum_{n=1}^{\infty}\Omega^{\circ}_{\Pi}(n)z^n=\frac{p_{\Pi}(z)}{(1-z)^{d+1}}.$$

Stapledon decomposition conjecture. The polynomial $p_{\Pi}(z)$ can be decomposed as

$$p_{\Pi}(z)=a(z)+zb(z),$$

where $a(z)$ and $b(z)$ are symmetric and satisfy

$$a(z)=z^{d}a\left(\frac{1}{z}\right),\qquad b(z)=z^{d-1}b\left(\frac{1}{z}\right),$$

and $b(z)$ and $-a(z)$ have nonnegative coefficients.

This is the order-polytope version of the graph-polynomial decomposition conjecture, and the latter follows by summing over the order polytopes associated with acyclic orientations. The paper reports computer experimentation supporting both conjectures but does not establish this one.