A permutation $π$ of a multiset is said to be a {\em quasi-Stirling} permutation if there does not exist four indices $i<j<k<\ell$ such that $π_i=π_k$ and $π_j=π_{\ell}$. For a multiset $\mathcal{M}$, denote by $\overline{\mathcal{Q}}_{\mathcal{M}}$ the set of quasi-Stirling permutations of $\mathcal{M}$. The {\em qusi-Stirling polynomial} on the multiset $\mathcal{M}$ is defined by $ \overline{Q}_{\mathcal{M}}(t)=\sum_{π\in \overline{\mathcal{Q}}_{\mathcal{M}}}t^{des(π)}$, where $des(π)$ denotes the number of descents of $π$. By employing generating function arguments, Elizalde derived an elegant identity involving quasi-Stirling polynomials on the multiset $\{1^2, 2^2, \ldots, n^2\}$, in analogy to the identity on Stirling polynomials. In this paper, we derive an identity involving quasi-Stirling polynomials $\overline{Q}_{\mathcal{M}}(t)$ for any multiset $\mathcal{M}$, which is a generalization of the identity on Eulerian polynomial and Elizalde's identity on quasi-Stirling polynomials on the multiset $\{1^2, 2^2, \ldots, n^2\}$. We provide a combinatorial proof the identity in terms of certain ordered labeled trees. Specializing $\mathcal{M}=\{1^2, 2^2, \ldots, n^2\}$ implies a combinatorial proof of Elizalde's identity in answer to the problem posed by Elizalde. As an application, our identity enables us to show that the quasi-Stirling polynomial $\overline{Q}_{\mathcal{M}}(t)$ has only real roots and the coefficients of $\overline{Q}_{\mathcal{M}}(t)$ are unimodal and log-concave for any multiset $\mathcal{M}$, in analogy to Brenti's result for Stirling polynomials on multisets.